黃潔鈺，南 哲，商學(xué)棽，王曉燕，杜 伊 (首都師范大學(xué)資源環(huán)境與旅游學(xué)院，北京 100048)

隨著人口的增長，人們對糧食需求增大，施肥量不斷升高，農(nóng)業(yè)污染已成為水環(huán)境污染的主要來源[1]。盡管氮輸入對糧食可持續(xù)發(fā)展十分必要，但由于過量的氮投入會增加環(huán)境的氮損失，從而對地表水、地下水以及土壤產(chǎn)生負面影響[2]，引起非點源污染和水體富營養(yǎng)化，對農(nóng)業(yè)養(yǎng)分管理是一個很大的挑戰(zhàn)[3]。東北地區(qū)作為我國重要的糧食生產(chǎn)基地，是我國非點源污染的熱點研究區(qū)[4]。2005—2015年間東北地區(qū)氮肥施用量呈逐年上升趨勢[5]，近年來雖有減緩，但仍處于一個較高水平，農(nóng)業(yè)活動帶來的非點源污染威脅不可忽視。

為了制定針對較低環(huán)境風(fēng)險的管理辦法，農(nóng)業(yè)系統(tǒng)養(yǎng)分管理或環(huán)境管理人員需要相應(yīng)的方法或工具評估風(fēng)險[6]，如指數(shù)法[7]和輸出系數(shù)法[8]等。LEMUNYON等[9]最先提出了對農(nóng)田磷素流失風(fēng)險進行分級評估的指數(shù)方法；而后學(xué)者們也發(fā)展了氮磷流失風(fēng)險評估指數(shù)[7,10]。這類指數(shù)方法主要考慮污染源因子和遷移因子，用簡單的數(shù)學(xué)方法評估各種污染源對水體氮、磷污染的風(fēng)險水平，具有數(shù)據(jù)需求低、操作簡便和運行周期短等特點，可實現(xiàn)對大型流域養(yǎng)分流失風(fēng)險的高精度評估[7]，已在國內(nèi)被應(yīng)用于分析流域非點源氮磷流失風(fēng)險的空間分布[11-12]，常見的評價指標為土壤全氮、化肥施用量、土壤侵蝕、降雨徑流、距河距離和坡度坡長等。由于不同地區(qū)的自然地理條件不同，如何根據(jù)區(qū)域具體特點對指標進行修改是需要注意的問題[13]。

農(nóng)業(yè)系統(tǒng)的氮流失是受多種因素影響的復(fù)雜時空變化過程[14-15]，而目前指標體系多應(yīng)用于評估空間分異，對時間分異的考慮較少。降雨的季節(jié)或年際差異對污染物輸出的主要途徑[16-17]以及河流氮通量[14,18]存在一定影響，從而影響氮流失風(fēng)險評價。目前，我國關(guān)于降水時空分異的研究多為降雨量較大的南方地區(qū)[11,14]。東北地區(qū)屬于溫帶季風(fēng)氣候區(qū)，降水時空分布不均勻[19]，年際差異大，降水年型變化對其非點源污染評估的影響值得關(guān)注。另一方面，農(nóng)用地類型或作物類型也是影響河流氮輸出的重要因素，不同作物施肥量和污染物輸出途經(jīng)存在較大差異[20-21]。近年來，為確保糧食安全，東北許多地區(qū)將旱地改為水田。由于水田在水肥管理和地面覆蓋等方面與旱地有較大差異，使得水田的非點源污染物輸出負荷大于旱地[21]，導(dǎo)致該地區(qū)農(nóng)業(yè)氮污染更為嚴峻。因此，在評估農(nóng)業(yè)種植氮流失風(fēng)險時，提高對這些控制因素的認識，對種植業(yè)氮流失治理具有一定指導(dǎo)意義。

筆者以東北3省(黑龍江、吉林、遼寧)為研究區(qū)，針對其種植結(jié)構(gòu)和氮污染來源特征，考慮降水年型差異和不同農(nóng)田水媒氮流失途徑，構(gòu)建符合區(qū)域特點的氮流失風(fēng)險評估體系，分析該地區(qū)氮流失風(fēng)險分布特征及主導(dǎo)影響因子，重點探究作物類型和降水年型對氮流失的影響，以期為該地區(qū)農(nóng)田養(yǎng)分管理和非點源污染防治提供參考。

1 材料與方法

1.1 研究區(qū)概況

東北3省位于我國東北部(圖1)，介于北緯38°43′～53°33′、東經(jīng)118°53′ ～135°05′之間，土地總面積約為7.87×105km2，平均海拔約為200 m。東北3省屬于溫帶季風(fēng)氣候區(qū)，年均降水量為300～800 mm，時空分布不均勻，降水多集中在6—9月，占全年降水的70%以上，且多以集中降雨形式出現(xiàn)。根據(jù)中國科學(xué)院地理科學(xué)與資源研究所的土地利用數(shù)據(jù)，2015年東北三省耕地面積為3.04×105km2，墾殖率為38.58%。水田和旱地面積比例約為1∶4.75。該地區(qū)農(nóng)業(yè)生產(chǎn)活動較多，主要種植制度為大田一熟，主要農(nóng)作物為玉米(50.80%，按播種面積計，后同)、水稻(18.82%)和大豆(14.89%)等。

1.2 數(shù)據(jù)來源及處理

該研究主要涉及的自然地理數(shù)據(jù)和農(nóng)業(yè)管理數(shù)據(jù)的來源和精度見表1[22-24]。采用Excel對數(shù)據(jù)進行預(yù)處理和分析，采用ArcGIS 10.2進行指數(shù)計算和數(shù)據(jù)空間分析，采用SPSS和RStudio平臺進行相關(guān)性分析、聚類分析和單因素方差分析等。

1.3 氮流失風(fēng)險評價體系

參考BECHMANN等[10]和CECCHI等[25]的研究，構(gòu)建適合該研究區(qū)的氮磷指數(shù)評價體系，評價因子分為源因子和遷移因子。源因子包括土壤全氮(SN，NS)和施氮強度(FI，IF)；遷移因子包括降雨徑流(R)、坡度坡長(LS，SL)、距河距離(D)、農(nóng)田徑流流失(RL，LR)和地下淋溶/側(cè)滲流失(LL，LL)。將各因子歸一化后，賦予不同權(quán)重，得到最終指數(shù)得分。早期體系中權(quán)重確定通常具有一定主觀性[7]，因此，有研究者結(jié)合指標的數(shù)值信息和數(shù)學(xué)方法進行改進，如改進的理想解法(TOPSIS)，相對常見的專家打分法更為客觀[26]。因子歸一化及TOPSIS求權(quán)重算法參考張汪壽等[26]的研究。最終構(gòu)建的氮指數(shù)計算公式為

(1)

式(1)中，I為氮指數(shù)；Si為歸一化的第i個源因子；Wi為第i個源因子的權(quán)重；Tj為歸一化的第j個遷移因子；Wj為第j個遷移因子權(quán)重。根據(jù)上述計算方法，將各歸一化因子進行加權(quán)疊加計算得到東北3省氮流失風(fēng)險分布，并運用ArcGIS中的自然斷點法進行分級[12]。自然斷點法通過尋求使每個分組與分組平均值的平均偏差最小化，同時使每個分組與其他分組的平均偏差最大化，以確定最佳分組[27]。其中，豐水年和平水年按照枯水年的劃分閾值進行劃分，確保分級采用統(tǒng)一閾值。

審圖號：GS(2022)1949號

1.4 源因子

1.4.1土壤全氮

土壤中氮元素含量高低直接影響其向水體傳輸?shù)牧縖28]。研究表明土壤全氮、堿解氮和速效氮等與包含無機氮的指標間具有空間相關(guān)性[29]，同時考慮到部分有機氮可能通過礦化作用轉(zhuǎn)化為無機氮，因此，以土壤全氮分布體現(xiàn)研究區(qū)土壤屬性差異分布對氮流失的影響。通過中國土壤數(shù)據(jù)庫獲得東北3省各樣點土壤類型和土壤全氮含量數(shù)據(jù)；對于同一地級市的樣點，按土壤亞類分類求取土壤全氮平均值，并賦值到土壤分類地圖中，得到土壤全氮空間分布。

1.4.2施氮強度

外源肥料投入是土壤中氮素的重要來源。通過查閱東北3省統(tǒng)計年鑒，獲得各地級市2014—2017年氮肥施用折純量，對于缺失數(shù)據(jù)通過其他年鑒數(shù)據(jù)或已有數(shù)據(jù)進行估算。按照土地利用分類，擬區(qū)分水田和旱地的氮肥施用強度。由于年鑒中施肥數(shù)據(jù)無法區(qū)分水田和旱地，將根據(jù)主要作物播種面積進行折算，公式如下：

(2)

式(2)中，F(xiàn)net為年鑒所示地級市所在省年氮素折純施用量(kg，以N計)表示，IF,p和IF,u分別為水田和旱地施氮強度，kg·hm-2；Ap和Au分別為水田和旱地作物播種面積，hm2；下標i表示第i種主要旱地作物。地級市年氮素折純施用量(Fnet)和作物播種面積(A)均來自各省統(tǒng)計年鑒，其中，復(fù)合肥按通用復(fù)合肥氮磷鉀比例1∶1∶1進行折純。

由于旱地作物類型較多，僅考慮東北3省主要作物(玉米、大豆、花生和蔬菜等)，其單位面積施氮量根據(jù)文獻[30-32]獲得。為了將總施肥量分配給不同作物，引入作物施肥系數(shù)(K)，K為某旱地作物單位面積施氮量與水稻單位面積施氮量的比值。令水稻K值為1，則其他作物i的K可表示為

Ki=IF，u,i/IF，p。

(3)

式(3)中，Ki為作物i的作物施肥系數(shù)。將式(3)代入式(2)，則水田和旱地施氮強度可分別表示為

(4)

IF,u=(Fnet-IF,p×Ap)/Au。

(5)

1.5 遷移因子

1.5.1降雨徑流

東北3省種植季為每年4—10月，涵蓋該區(qū)域全年90%的降雨，而種植季降雨的年際差異較大，2002—2018年東北3省種植季最大年降雨量約為最小年降雨量的1.5倍。采用Wischmeier經(jīng)驗公式[33]計算降雨量和徑流影響因子(R,MJ·mm·hm-2·h-1·a-1)，月降雨量取4—10月東北地區(qū)農(nóng)作物種植期降雨數(shù)據(jù)：

(6)

式(6)中，P為年平均降雨量，mm；Pi為月降雨量，mm。

為了研究降雨年際差異對東北3省種植季氮流失特征的影響，采用干旱指數(shù)(DI，ID)進行降水年型劃分[34]。干旱指數(shù)公式如下：

ID=(P-M)/σ。

(7)

式(7)中，ID為干旱指數(shù)；P為種植季降雨總量，mm；M為所衡量指標的多年平均值，mm；σ為多年降雨標準差，mm。將干旱指數(shù)ID≥0.35的年份劃分為豐水年，將-0.35≤ID<0.35的年份劃分為平水年，將ID<-0.35的年份劃分為枯水年。

1.5.2坡度坡長

地形因子主要由坡度和坡長決定。坡度體現(xiàn)某區(qū)域地形起伏狀況，影響產(chǎn)流能力；坡長則通過影響坡面降雨徑流的流速與流量來影響水流挾沙力，兩者會聯(lián)合影響土壤侵蝕強度。計算公式[35]為

SL=(L/22.1)m(65.41×sin2θ+4.56×sinθ+0.065)。

(8)

式(8)中，SL為坡度坡長因子；L為坡長，m；θ為坡度角，(°)；m為坡長指數(shù)。

1.5.3距河距離

距離因子反映某種土地利用單元距離匯水路線的距離，潛在的氮磷污染源距離河流越近，越容易發(fā)生遷移且流失風(fēng)險性越高。采用經(jīng)驗公式[25-26]進行計算，其表達式為

Di=e-0.090 533×di。

(9)

式(9)中，Di為像元i的距離因子；di為像元i到匯水路線的柵格距離；-0.090 533為經(jīng)驗指數(shù)。

1.5.4農(nóng)田徑流流失和地下淋溶/側(cè)滲流失

一般地，旱田氮流失的主要途徑為地表徑流和地下淋溶[5,16]。而在水田集中區(qū)，農(nóng)業(yè)面源污染主要通過暴雨徑流、人為排水和側(cè)滲3種途徑進入渠系[21]。側(cè)滲是指水田內(nèi)污染物的水平遷移流失，由于水田犁底層的存在，滲漏率較低，側(cè)滲是渠系化水田區(qū)面源污染的主要輸出方式。人為排水雖會產(chǎn)生一定氮排放，但持續(xù)時間和排放次數(shù)有限[20]，因此，在計算水田徑流流失時根據(jù)人工排水氮流失進行修正。綜上，設(shè)置氮徑流流失及地下流失(淋溶或側(cè)滲)因子。根據(jù)《全國農(nóng)田面源污染排放系數(shù)手冊》[36]及?；莸萚20]的研究，參考孫鋮等[5]的加權(quán)平均法，計算農(nóng)田氮徑流、旱田地下淋溶以及水田側(cè)滲流失系數(shù)(表2[20,36])，采用面積加權(quán)法，各地級市氮徑流(LR)或地下淋溶/側(cè)滲流失(LL)計算公式為

(10)

(11)

式(10)～(11)中，μi和λi分別為第i類種植模式農(nóng)田氮徑流和淋溶/側(cè)滲流失系數(shù)，kg·hm-2；Si為第i類種植模式播種面積占水田或旱田總播種面積的比例，%。

表2 農(nóng)田氮流失系數(shù)[20,36]

2 結(jié)果與分析

2.1 氮指數(shù)及其因子分布

圖2顯示，土壤全氮因子從整體來看以黑龍江為最高，三江平原及松嫩平原土壤肥沃，吉林東南部地區(qū)以林地為主，雖然耕地分布較少，但地力水平較西北部更高。施氮強度空間分布差異較大，遼寧、吉林和黑龍江平均施氮強度分別為192.89、185.59和76.55 kg·hm-2。施氮強度及其分布與文獻[6]相符。研究區(qū)遷移因子空間分布見圖3。

審圖號：GS(2022)1949號

審圖號：GS(2022)1949號

如圖3所示，坡度坡長因子在空間分布上與地形地貌關(guān)系密切，長白山和小興安嶺山地地勢起伏相對較大。在農(nóng)田流失因子方面，農(nóng)田徑流流失較高的地區(qū)集中在黑龍江省三江平原、吉林省中部以及遼寧省遼河平原。三江平原和遼河平原水稻種植分布較多，水田側(cè)滲流失明顯。淋溶強度較高的地區(qū)位于松嫩平原和遼河平原，與已有研究[5]相符。降雨徑流因子在空間分布上呈現(xiàn)東南部最高，中部松嫩平原中等，東北部地區(qū)較低。從平均值來看，豐、平、枯水年降雨徑流因子分別為199.5、185.3和171.4 MJ·mm·hm-2·h-1·a-1，呈現(xiàn)豐水年>平水年>枯水年。

采用TOPSIS法計算的權(quán)重結(jié)果見表3。在源因子中，豐水年與枯水年施氮強度權(quán)值相等，土壤全氮權(quán)值相等。不同降水年型遷移因子權(quán)值計算結(jié)果存在一定差異，與平水年和枯水年相比，豐水年降雨徑流因子分布差異較高，所占比例略有增加。

表3 TOPSIS法權(quán)值計算結(jié)果

運用ArcGIS中的自然斷點法將枯水年氮指數(shù)結(jié)果分為4個等級[12]。為便于比較，按照枯水年分級閾值劃分，得到不同降水年型氮流失風(fēng)險等級分布(圖4)。氮流失較高區(qū)域較為集中，極高、高風(fēng)險區(qū)面積約占總面積的22%，以耕地為主，集中在吉林省中部、遼寧省中部和南部以及黑龍江省中部和東部。

極高、高、中和低分別指氮指數(shù)為>0.131～0.385、>0.077～0.131、>0.035～0.077和0～0.035。 審圖號：GS(2022)1949號

2.2 氮流失風(fēng)險的空間特征區(qū)劃

為了便于區(qū)劃管理，根據(jù)東北3省各市氮指數(shù)及因子特征進行聚類分析(圖5)，其中，一級控制區(qū)氮指數(shù)最高，同時，具備高施氮強度和高農(nóng)田流失因子。二級控制區(qū)為氮指數(shù)在4組中處于中高水平，其施氮強度、農(nóng)田流失因子處于中高水平，而降雨徑流因子、坡度坡長因子較高。三級控制區(qū)氮指數(shù)處于中低水平，各因子顯示為中等水平。四級控制區(qū)具有較高的降雨徑流和土壤全氮因子，但因其施氮強度和農(nóng)田流失均較低，其氮指數(shù)也較低。

2.3 氮流失風(fēng)險與作物類型

研究區(qū)超過80%的水田屬于極高、高風(fēng)險區(qū)，約35%的旱地屬于極高、高風(fēng)險區(qū)。可見，研究區(qū)旱地氮流失風(fēng)險占主導(dǎo)地位，而水田是高風(fēng)險的關(guān)鍵區(qū)域。在不同作物的種植面積比例(表4)方面，按2.2節(jié)中對地級市的分組，在一級控制區(qū)，玉米種植面積比例最高，水稻面積占比則位列第2，除盤錦市以水稻種植為主(約占75%)外，其他5市均以玉米種植為主(占比均不低于55%)；在4個級別控制區(qū)中，二級控制區(qū)水稻面積占比最高。在二級控制區(qū)，玉米種植面積占比最高，氮指數(shù)與玉米面積占比呈現(xiàn)協(xié)同增加趨勢，這體現(xiàn)了玉米作為旱地作物在東北3省氮流失風(fēng)險中的主導(dǎo)地位。此外，值得注意的是，大豆種植面積比例越高，氮指數(shù)則越低，這可能是大豆的農(nóng)田氮流失系數(shù)較小所致。

審圖號：GS(2022)1949號

表4 東北3省城市聚類分組氮指數(shù)和種植面積比例

2.4 不同降水年型氮流失風(fēng)險分布

為探究降水年型的影響，統(tǒng)計研究區(qū)各縣域降雨徑流因子和氮指數(shù)進行單因素方差分析(表5)。結(jié)果顯示，在整個研究區(qū)和各控制分區(qū)內(nèi)，不同水文年降雨徑流因子間存在顯著差異(P<0.05)。而對于氮指數(shù)，僅在整個研究區(qū)和二級控制區(qū)內(nèi)不同水文年氮指數(shù)間存在顯著差異(P<0.05)，而在一、三和四級控制區(qū)則無顯著差異。這表明降水年型差異會造成降雨徑流因子和氮流失風(fēng)險的差異，而氮流失風(fēng)險的影響在二級控制區(qū)較為顯著。多重比較中，豐水年與平水年、豐水年與枯水年間差異較明顯，枯水年與平水年間差異不明顯，這說明豐水年對氮流失風(fēng)險分布的影響較大。

3 討論

3.1 作物類型對氮流失風(fēng)險的影響

作物類型對流域氮肥投入和全氮損失強度有較大影響[15]。多項研究表明，水稻作為主要水田作物，其氮素地下淋溶流失高于旱地作物[37]，可達到旱地作物的5倍以上[21]。為增加糧食產(chǎn)量，黑龍江省開展了“兩江一湖”土地整理重大項目，進行旱改水灌溉工程建設(shè)，三江平原分布的大量水田是該地區(qū)氮流失因素的主要來源。遼河平原氮流失風(fēng)險較高，這與遼河沿岸大量的水稻種植有關(guān)。

由于東北地區(qū)水田和旱地面積占比差異較大，旱地氮流失風(fēng)險仍占主導(dǎo)地位。玉米是東北地區(qū)種植面積最大的喜氮旱地作物，其單位經(jīng)濟產(chǎn)量需氮量較水稻、小麥等糧食作物更高，而其化肥利用率較低[38]，這意味著其氮流失(氣態(tài)和溶解態(tài))風(fēng)險更高。吉林省為玉米種植大省，玉米播種面積和產(chǎn)量也居于3省首位，這使吉林總施肥量遠高于遼寧和黑龍江。此外，蔬菜種植對氮流失風(fēng)險的影響也值得注意。研究表明，除化肥施用量高于其他旱地作物和水稻以外，蔬菜種植系統(tǒng)氮流失率也高于其他系統(tǒng)[39]。筆者研究中如遼寧省沈陽、鞍山、營口和朝陽，其蔬菜種植面積占比相對較高(>12%)，當?shù)厥卟朔N植的化肥施用策略亟需受到重視。

表5 降水年型對降雨徑流和氮指數(shù)影響的顯著性

3.2 降水年型對氮流失風(fēng)險的影響

降水作為驅(qū)動力，是影響氮流失的一個重要因素[40]，地表徑流和地下淋溶是降水影響農(nóng)田氮流失的主要途徑[5,16]。研究表明，在一定徑流量范圍內(nèi)硝酸鹽濃度隨徑流量的增加而增加[18]。HUANG等[14]對我國浙江省的研究發(fā)現(xiàn)，枯水年流域硝態(tài)氮輸出為豐水年的10%～50%。筆者研究中，降水年型對降雨徑流因子和氮指數(shù)均有顯著影響，氮流失風(fēng)險分布對降雨年型變化的響應(yīng)，體現(xiàn)在極高風(fēng)險區(qū)面積增加，與已有研究結(jié)果[14,18]相符。而有研究表明降雨不僅對氮流失有影響[14]，而且對作物產(chǎn)量、水分利用效率和氮素利用效率變化也存在一定影響[34]。因此，根據(jù)不同的降水年型制定不同的施肥計劃，對于防控農(nóng)業(yè)非點源污染尤為重要；而氮流失隨降雨季節(jié)變化而變化的狀況也應(yīng)受到關(guān)注，以便在更精確的時間尺度上制定相應(yīng)的管控方案。

3.3 東北地區(qū)氮流失管理建議

目前，我國的農(nóng)業(yè)污染管理依賴市、區(qū)級行政機構(gòu)執(zhí)行，以地級市作為氮流失風(fēng)險特征分類單元，有利于提出更符合現(xiàn)行管理機制的建議。一級控制區(qū)地級市為長春、盤錦、遼源、沈陽、四平和錦州，關(guān)鍵因子為施氮強度和農(nóng)田流失，主要作物為玉米(除盤錦以外)，需要重點管控其氮肥輸入，其次減少氮肥流失，如測土施肥、有機肥替代、保護性耕作、秸稈覆蓋等措施以減少農(nóng)田氮流失[41]；對于水田，可以施行淺施灌溉和間歇灌溉，通過減少稻田與排水溝之間的水位差異來減少滲漏量[42]。二級控制區(qū)地級市為通化、丹東、鞍山等8個市，玉米、水稻種植比例都較高，施氮強度及與產(chǎn)流相關(guān)的因子較高，需要在適當減少其農(nóng)業(yè)化肥流失的基礎(chǔ)上，選用如植被緩沖帶等攔截措施，加強對該區(qū)域土壤氮素隨徑流流失的攔截[43]，提高豐水年氮流失管控預(yù)警。三級控制區(qū)氮指數(shù)水平較低，極高風(fēng)險區(qū)呈零星分布，包括七臺河、雙鴨山和雞西等11個地級市，特別是黑龍江三江平原各市水田較為集中，需要重視其水田氮流失管理。四級控制區(qū)氮指數(shù)水平最低，各地級市氮流失風(fēng)險暫時無需特殊管控。

4 結(jié)論

東北地區(qū)作為我國重要的糧食生產(chǎn)基地，需要全面、系統(tǒng)的研究以優(yōu)化其農(nóng)田管理和非點源污染防控方案。通過構(gòu)建氮指數(shù)評估東北3省氮流失風(fēng)險，并分析其分布特征。通過相關(guān)分析確定了影響該地區(qū)氮流失風(fēng)險分布的關(guān)鍵因子。根據(jù)氮流失特征對地級市進行聚類劃分，著重探討作物類型和降雨年型對氮流失風(fēng)險的影響。東北3省極高、高風(fēng)險區(qū)域以耕地為主，集中在吉林省中部、遼寧省中部和南部以及黑龍江省中部和東部。在縣級區(qū)劃上，降水年型對降雨徑流因子和氮流失風(fēng)險的分布有顯著影響，特別是豐水年對氮流失風(fēng)險的影響較大。由于較高的施氮強度和農(nóng)田徑流損失，長春、盤錦、遼源、沈陽、四平和錦州氮指數(shù)較高，是東北3省需要首批重點管控的地級市。