吳漢卿，萬 煒，單艷軍，陳延華，李 強(qiáng)，李存軍，胡海棠，張寶貴

基于磷指數(shù)模型的海河流域農(nóng)田磷流失環(huán)境風(fēng)險評價

吳漢卿1，萬 煒1，單艷軍1，陳延華2，李 強(qiáng)1，李存軍3，胡海棠3，張寶貴1※

（1. 中國農(nóng)業(yè)大學(xué)土地科學(xué)與技術(shù)學(xué)院，北京 100193；2. 北京市農(nóng)林科學(xué)院植物營養(yǎng)與資源研究所，北京 100097；3. 北京農(nóng)業(yè)信息技術(shù)研究中心，北京 100097）

農(nóng)田面源磷流失是導(dǎo)致水體富營養(yǎng)化的主要原因，識別農(nóng)田磷流失的關(guān)鍵源區(qū)、影響因子是農(nóng)田面源污染防治的重要環(huán)節(jié)。該研究以海河流域為研究區(qū)，采用磷指數(shù)模型，選取土壤有效磷含量、磷肥施用量作為源因子，以土壤侵蝕模數(shù)、年徑流深、農(nóng)田和水體間歸一化距離指數(shù)作為遷移因子，結(jié)合GIS技術(shù)評估海河流域農(nóng)田磷流失風(fēng)險，并利用結(jié)構(gòu)方程模型研究農(nóng)田磷流失風(fēng)險指數(shù)與各影響因子間關(guān)系。結(jié)果表明：1）海河流域農(nóng)田土壤有效磷、磷肥施用量、土壤侵蝕模數(shù)、年徑流深及歸一化距離指數(shù)處于中-低、中-高、極低、中-高和高級別風(fēng)險等級的區(qū)域面積占比最高，分別占農(nóng)田總面積的66.5%、61.1%、99.0%、54.2%和64.8%；2）影響農(nóng)田磷流失的關(guān)鍵因子為遷移因子，其中關(guān)鍵的遷移因子為年徑流深及歸一化距離指數(shù)；3）源因子與遷移因子間呈極顯著負(fù)相關(guān)（<0.01），土壤性質(zhì)（包括土壤質(zhì)地及有機(jī)碳含量等）與源因子呈極顯著負(fù)相關(guān)（<0.01），與遷移因子呈極顯著正相關(guān)（<0.01）；4）海河流域農(nóng)田磷流失關(guān)鍵源區(qū)位于黃河北岸的山東省和河北省東南部的平原農(nóng)耕區(qū)、海河流域西北部的山區(qū)地帶。該研究結(jié)果對流域尺度上農(nóng)田磷流失研究的方法創(chuàng)新有所裨益，可為海河流域農(nóng)田面源污染防治提供科學(xué)參考。

農(nóng)田；模型；面源污染；磷流失；磷指數(shù)模型；關(guān)鍵源區(qū)；結(jié)構(gòu)方程模型；海河流域

0 引 言

磷（P）是作物生長的主要限制性養(yǎng)分，農(nóng)業(yè)生產(chǎn)過程中施用大量磷肥[1]。農(nóng)田土壤磷素在降水沖刷及土壤侵蝕作用下，以徑流、壤中流、吸附于泥沙表面等方式進(jìn)入受納水體[2]。據(jù)估計，全球每年有3×106~4×106t P2O5從土壤遷移到水體中[3]，其中約19.5 kg/hm2來自中國農(nóng)田[4]。農(nóng)業(yè)生產(chǎn)引發(fā)的磷流失是水體富營養(yǎng)化的最直接原因，已成農(nóng)田面源污染防控研究的熱點問題[5]。

農(nóng)田面源污染因受土壤類型、地形特征、氣候條件和人類活動等諸多因素影響具有高時空異質(zhì)性，農(nóng)田面源污染防控難度較大[3,6]。在流域尺度上，識別出高污染風(fēng)險的關(guān)鍵源區(qū)（Critical Source Areas，CSAs）是農(nóng)田面源污染防控的前提[7-8]。磷指數(shù)模型（Phosphorus Index Model，PIM）是識別關(guān)鍵源區(qū)及評估污染風(fēng)險的經(jīng)典模型之一[9-10]。該模型作為磷素綜合管理工具，已被國內(nèi)外學(xué)者廣泛應(yīng)用[10-12]。磷指數(shù)模型最初根據(jù)源因子與遷移因子特征用以評估田間尺度上的磷污染風(fēng)險等級，后來Gburek和Sharpley等將農(nóng)田到河網(wǎng)的水文連通性、距離等因子也考慮在內(nèi)，對流域尺度上的磷流失潛在風(fēng)險指數(shù)進(jìn)行了定性化評估[9-10]。Buchanan等開發(fā)了基于徑流行程時間的磷指數(shù)模型，以改進(jìn)磷指數(shù)表征排水對磷流失潛在風(fēng)險的影響[13]。在中國，張淑榮等最先采取磷指數(shù)模型，選取土壤侵蝕模數(shù)、地表徑流、潛在污染源距水體的距離、土壤磷含量、磷肥施用量和施用方法等因子，識別出了小流域磷流失高風(fēng)險區(qū)[14]。其后，學(xué)者們通過將土壤磷吸附指數(shù)（Phosphorus Sorption Index，PSI）、磷飽和度（Degree of Phosphorus Saturation，DPS）、溶解態(tài)磷、顆粒態(tài)磷及地形指數(shù)因子等指標(biāo)納入評價體系，對磷指數(shù)模型進(jìn)行修正，對流域尺度上磷流失潛在風(fēng)險進(jìn)行評價并識別出關(guān)鍵源區(qū)[15-16]。綜上，修正的磷指數(shù)模型通過定性評估流域尺度上農(nóng)田磷流失潛在風(fēng)險指數(shù)，從而識別關(guān)鍵源區(qū)，但對如何識別關(guān)鍵影響因子尚不清楚。

2011年“重點流域‘十二五’農(nóng)業(yè)面源污染防治研究報告”指出，在8個重點流域中，海河流域的總磷總負(fù)荷量及排放強(qiáng)度分別位居第一和第二位，總磷總負(fù)荷量占8個重點流域總量的43.1%，排放強(qiáng)度是8個重點流域平均值的2.1倍。2018年《中國生態(tài)環(huán)境狀況公報》指出海河流域為中度污染，地表水中劣Ⅴ類占20.0%。因此，海河流域農(nóng)田磷流失環(huán)境風(fēng)險評價對于磷污染防治具有重要意義。盡管磷指數(shù)模型在小流域尺度上農(nóng)田磷流失關(guān)鍵源區(qū)識別方面得到了廣泛應(yīng)用和發(fā)展，而海河流域空間范圍廣，空間分異規(guī)律明顯，農(nóng)田面源污染的驅(qū)動機(jī)制難以較好地把握，因此在海河流域農(nóng)田面源污染防治中的應(yīng)用極少，關(guān)于農(nóng)田磷流失關(guān)鍵源區(qū)、關(guān)鍵影響因子識別及磷流失潛在風(fēng)險指數(shù)與各影響因子間關(guān)系研究亦是鮮見報道。

因此，本研究以海河流域為研究對象，以1 km柵格尺度為最小地域單元（傳統(tǒng)方法主要以行政區(qū)為最小地域單元進(jìn)行空間分析），利用磷指數(shù)模型并結(jié)合GIS技術(shù)，通過分析影響海河流域農(nóng)田磷流失的影響因子及潛在風(fēng)險指數(shù)的空間分布，識別出海河流域農(nóng)田磷流失的關(guān)鍵源區(qū)，在此基礎(chǔ)上采用結(jié)構(gòu)方程模型方法定量地探討海河流域農(nóng)田磷流失潛在風(fēng)險與各影響因子間關(guān)系，進(jìn)而識別出關(guān)鍵影響因子。本研究結(jié)果希望對流域尺度上農(nóng)田磷流失研究方法創(chuàng)新有所裨益，并為海河流域農(nóng)田面源污染防治提供科學(xué)參考。

1 研究區(qū)概況及數(shù)據(jù)源

1.1 研究區(qū)概況

海河流域地處華北（112°～120°E、35°～43°N），共轄5?。ê颖笔?、山西省、山東省、河南省、遼寧?。?直轄市（北京市、天津市）及1自治區(qū)（內(nèi)蒙古自治區(qū)），流域總面積達(dá)3.18×105km2。該區(qū)域整體地勢為西北高、東南低，西北部為燕山-太行山山脈，東南部則為黃河及海河形成的沖積平原（圖1a）。境內(nèi)主要河流為海河與灤河。氣候類型以暖溫帶大陸性季風(fēng)氣候為主，年均溫為1.5～14 ℃，多年平均降水量為539 mm。土壤類型主要有華北山地棕壤與褐土，海河平原潮土、黃壚土、鹽堿土及風(fēng)沙土，內(nèi)蒙古高原栗鈣土與綿土。

土地利用類型以耕地為主，其中耕地面積占區(qū)域總面積的48.9%（圖1b）。該區(qū)典型的農(nóng)作物種植方式為冬小麥-夏玉米輪作。該區(qū)域是中國三大糧食生產(chǎn)基地之一和京、津等大中城市“菜籃子”主要生產(chǎn)基地，人口密度較高，人為作用強(qiáng)度較大，農(nóng)田面源污染嚴(yán)重，引發(fā)水體富營養(yǎng)化等環(huán)境問題[17]。

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

1）縣域主要作物的糧食種植面積及磷肥施用量數(shù)據(jù)：從研究區(qū)各市統(tǒng)計年鑒（2019年）提取出化肥折純量，再根據(jù)2019年中國統(tǒng)計年鑒中氮肥、磷肥及鉀肥施用量分別為2065.43、728.88和590.28萬t，參考全國平均比例折算出研究區(qū)磷肥施用量數(shù)據(jù)。

2）土地利用、DEM及土壤質(zhì)地數(shù)據(jù)：由中國科學(xué)院資源環(huán)境數(shù)據(jù)云平臺（http://www.resdc.cn/）下載得到。土地利用數(shù)據(jù)為2015年的柵格數(shù)據(jù)產(chǎn)品。研究區(qū)數(shù)字高程模型（DEM）利用ArcGIS水文分析模塊提取研究區(qū)河網(wǎng)空間分布[18]，并用于提取出研究區(qū)內(nèi)的坡度及坡長因子。土壤質(zhì)地數(shù)據(jù)用于提取研究區(qū)土壤的砂粒、粉粒和黏粒含量。上述數(shù)據(jù)的空間分辨率均為1 km。

3）土壤全磷含量和土壤有機(jī)質(zhì)含量數(shù)據(jù)：空間分辨率為1 km，從國家地球系統(tǒng)科學(xué)數(shù)據(jù)中心：土壤分中心（http://soil.geodata.cn/）下載得到。本研究在前人研究基礎(chǔ)上[19-20]，通過將有機(jī)質(zhì)含量除以1.724[19]和土壤全磷含量乘以轉(zhuǎn)化系數(shù)（6%）[20]的方式分別得到土壤有機(jī)碳及有效磷含量。

4）1980—2015年逐月降水量數(shù)據(jù)：由中國氣象數(shù)據(jù)網(wǎng)（http://data.cma.cn/）下載得到，先通過ANUSPLIN氣象插值軟件[18]經(jīng)薄板樣條函數(shù)處理得到研究區(qū)1 km空間分辨率的研究區(qū)逐月降水量數(shù)據(jù)，并處理得到多年平均逐月降水量。

5）2019年逐月的MODIS-NDVI數(shù)據(jù)：通過美國國家航空航天局（National Aeronautics and Space Administration，NASA，https://ladsweb.modaps.eosdis.nasa.gov/search/）下載得到。通過MRT軟件對其進(jìn)行鑲嵌、投影轉(zhuǎn)化、重采樣及裁剪等處理，再通過最大值合成法[21]得到研究區(qū)NDVI數(shù)據(jù)，用以計算作物生育期的植被覆蓋度。

2 研究方法

2.1 磷指數(shù)模型

磷指數(shù)模型是一種綜合考慮農(nóng)田磷流失的主要影響因子，進(jìn)而評價流域內(nèi)不同區(qū)域發(fā)生磷流失潛在風(fēng)險高低的經(jīng)典模型[9-11,22]。本研究基于磷指數(shù)模型來評估海河流域農(nóng)田磷流失風(fēng)險，參照國內(nèi)外相關(guān)研究成果[9-11,22]，將影響磷流失的因子劃分為源因子與遷移因子。在本研究中，選取了土壤有效磷含量及磷肥施用量作為源因子，這是因土壤有效磷能反映土壤本身磷含量及磷流失潛在能力，而磷肥施用量作為主要的外源磷投入，對土壤磷素含量具有重要影響[11]；選取土壤侵蝕模數(shù)、年徑流深及農(nóng)田與河網(wǎng)間歸一化距離指數(shù)作為遷移因子，這是由于此三者為農(nóng)田磷流失至受納水體過程中最主要的影響因素，分別反映了土壤侵蝕、地表徑流及農(nóng)田至受納水體的距離[11-12]。

采用ArcGIS中的自然斷點法對源因子、遷移因子進(jìn)行等級劃分[12]，分為5級（極低、低、中、高、極高），并對各等級進(jìn)行風(fēng)險賦值，因子風(fēng)險值及權(quán)重值確定，均參照前人研究結(jié)果進(jìn)行確定[9-11]，具體數(shù)值見表1，進(jìn)而采用磷指數(shù)模型計算出海河流域農(nóng)田磷流失風(fēng)險指數(shù)來反映研究區(qū)農(nóng)田磷流失風(fēng)險等級情況，再根據(jù)風(fēng)險等級識別出海河流域農(nóng)田磷流失關(guān)鍵源區(qū)[9-11]。

表1 農(nóng)田磷流失風(fēng)險等級劃分

注：RV代表風(fēng)險值。

Note: RV represents Risk Value.

2.2 源因子計算

2.2.1 土壤有效磷含量

表層土壤有效磷含量作為土壤磷流失的物質(zhì)基礎(chǔ)，含量越高則磷流失風(fēng)險越大[14]。本研究利用ArcGIS生成海河流域1 km柵格土壤有效磷含量空間分布。

2.2.2 磷肥施用量

2015年農(nóng)業(yè)農(nóng)村部339個國家級基層肥料零售網(wǎng)點農(nóng)戶調(diào)查資料顯示，海河流域磷肥施用類型主要以化肥為主，施肥方式以基肥形式為主。因此，本研究利用ArcGIS將統(tǒng)計年鑒中的縣域磷肥施用量數(shù)據(jù)按照土地利用類型對研究區(qū)內(nèi)耕地進(jìn)行空間插值，精確至1 km空間分辨率的耕地柵格中，從而生成海河流域磷肥施用量空間分布圖。

2.3 遷移因子計算

2.3.1 土壤侵蝕模數(shù)

本研究選用修正通用土壤流失方程（Revised Universal Soil Loss Equation, RUSLE）來估算海河流域土壤侵蝕模數(shù)，作為農(nóng)田土壤磷流失風(fēng)險的重要遷移因子之一[9-11]。利用如下公式進(jìn)行計算，并獲取海河流域土壤侵蝕模數(shù)（）分布

式中為土壤侵蝕模數(shù)，t/(km2·a)；為降雨侵蝕力因子，MJ·mm/(hm2·h·a)；為土壤可蝕性因子，t·hm2·h/(MJ·mm·hm2)；為地形因子，無量綱；為植被覆蓋與管理因子，無量綱；P為水土保持措施因子，無量綱。

1）降雨侵蝕力因子

降雨侵蝕力因子（）是用來評估降雨潛在引起的土壤侵蝕能力的動力學(xué)因子。本研究基于海河流域1980–2015年的平均逐月降雨量數(shù)據(jù)，采用Wischmeier的公式計算[23]，并獲取海河流域分布

式中為多年平均降雨量，mm；p為第月多年平均降雨量，mm。

2）土壤可蝕性因子

土壤可蝕性因子（）定義為單位降雨侵蝕力造成的單位面積侵蝕土壤數(shù)量，反映了不同土壤的侵蝕速率，是影響土壤侵蝕的內(nèi)在因子，與土壤物理性質(zhì)如質(zhì)地和有機(jī)質(zhì)含量等緊密有關(guān)。本研究利用研究區(qū)1 km柵格土壤砂粒、粉粒、黏粒及有機(jī)碳數(shù)據(jù)，采用Williams等建立的EPIC（Erosion Productivity Impact Calculator）模型計算[24]，并獲得海河流域分布

式中SAN、SIL和CLA分別為0~20 cm土層砂粒、粉粒和黏粒含量，%；和SOM分別為土壤有機(jī)碳和有機(jī)質(zhì)含量，%；S為粉粒和黏粒含量總和，%。

3）地形因子

地形因子（）亦稱坡長坡度因子，主要反映坡長因子（）與坡度因子（）對土壤侵蝕的影響。定義為其他條件一致情況下，特定坡長的坡地土壤流失量與標(biāo)準(zhǔn)小區(qū)坡長的坡地土壤流失量比值；定義為其他條件一致情況下，特定坡度的坡地土壤流失量與標(biāo)準(zhǔn)小區(qū)坡度的坡地土壤流失量比值。本研究利用ArcGIS水文分析功能從數(shù)字高程模型中提取出坡長與坡度，再根據(jù)下列公式依次計算出、與[25]，并獲得海河流域地形因子（）分布

式中和分別代表坡長因子和坡度因子，無量綱；是由DEM提取出的水平投影坡長，m；22.13是RUSLE采用的標(biāo)準(zhǔn)小區(qū)坡長，m；是由DEM提取出的坡度，(°)；是坡度坡長指數(shù)，無量綱。

4）植被覆蓋與管理因子

植被覆蓋與管理因子（）是指在一定條件下，有作物或植物覆蓋及實施管理情況下坡地土壤流失量與耕種過后連續(xù)休閑地土壤流失量之比。本研究先利用ENVI軟件獲取歸一化植被指數(shù)（NDVI），再計算植被覆蓋度（），最后利用估算[26]，獲得海河流域分布

式中為植被覆蓋度，%；NDVI為歸一化植被指數(shù)；NDVImax和NDVImin分別為研究區(qū)內(nèi)NDVI的最大值和最小值。

5）水土保持措施因子

水土保持措施因子（P）是指采取特定水土保持措施下的土壤流失量與相應(yīng)未采取措施順坡耕作下的土壤流失量之比。本研究采用被廣泛應(yīng)用于大尺度流域P值計算的經(jīng)驗?zāi)Ｐ蚚27-28]，獲得海河流域P分布

式中是由DEM提取出的坡度，%。

2.3.2 年徑流深

年徑流深（R）定義為每年將徑流量均勻地鋪在整個流域面積上所相當(dāng)?shù)乃畬由疃?，是一種度量地表徑流的經(jīng)典方法。本研究利用海河流域近35 a（1980–2015年）的多年平均逐月降雨量數(shù)據(jù)計算出多年平均降雨量，以下列公式計算R，并獲得海河流域R分布

式中R為年徑流深，mm/a；為年徑流系數(shù)。參考《海河流域水利手冊》[29]和劉春蓁等研究結(jié)果[30]，本研究中取值為0.12。

2.3.3 歸一化距離指數(shù)

農(nóng)田與受納水體的距離是影響磷流失的重要因子，距離越遠(yuǎn)，在遷移過程中被稀釋及截留的可能性越大，即潛在的磷流失風(fēng)險愈低，反映了農(nóng)田與河網(wǎng)間的距離及潛在的磷流失風(fēng)險[31]。本研究在前人[32]的研究基礎(chǔ)上，將研究區(qū)內(nèi)農(nóng)田與河網(wǎng)間距離進(jìn)行歸一化處理，首次提出了歸一化距離指數(shù)（Normalized Differential Distance Index，NDDI）這一概念，得到農(nóng)田與河網(wǎng)間距離指數(shù)。NDDI指數(shù)與傳統(tǒng)的主觀設(shè)定農(nóng)田與水體間距離相比，具有客觀定量評價的優(yōu)勢，使研究結(jié)果更加精準(zhǔn)與可信。

式中為各像元中農(nóng)田到河網(wǎng)間距離，max和min分別為農(nóng)田與河網(wǎng)間距離的最大值和最小值。

2.4 磷流失風(fēng)險指數(shù)計算

本研究在前人提出的磷指數(shù)模型評價指標(biāo)體系基礎(chǔ)上[11-12]，根據(jù)研究區(qū)源因子及遷移因子的分級標(biāo)準(zhǔn)與權(quán)重，計算出磷流失風(fēng)險指數(shù)。

式中PI表示磷流失風(fēng)險指數(shù)；S和W分別為源因子評價指標(biāo)對應(yīng)的風(fēng)險值及權(quán)重，T和W分別為遷移因子評價指標(biāo)對應(yīng)的風(fēng)險值及權(quán)重。PI的計算通過對源因子指標(biāo)與遷移因子指標(biāo)相乘，以此反映遷移因子對源因子的限制作用。

2.5 統(tǒng)計分析

結(jié)構(gòu)方程模型（Structural Equation Modeling，SEM）是一種驗證性的多元統(tǒng)計方法，能反映潛變量與外顯變量之間的關(guān)系，同時探討不同因子對同一過程直接或間接的影響，判別不同因子的相對重要性。因此，本研究利用Amos 21.0構(gòu)建SEM，分析土壤性質(zhì)（包括土壤質(zhì)地及有機(jī)碳含量等）與磷流失風(fēng)險指數(shù)、源因子（包括土壤有效磷含量、磷肥施用量等）及遷移因子（包括土壤侵蝕模數(shù)、年徑流深及耕地與水體間距離因子等）間的關(guān)系。在構(gòu)建SEM之前，采用SPSS 21.0對PI及各因子進(jìn)行相關(guān)分析。使用卡方值/自由度（CHI/DF）、配適度指數(shù)（Goodness-of-Fit Index，GFI）、調(diào)整后配適度指數(shù)（Adjusted Goodness-of-Fit Index，AGFI）、赤池信息量準(zhǔn)則（Akaike Information Criteria，AIC）和漸進(jìn)誤差均方根（Root Mean Square Error of Approximation，RMSEA）對SEM的擬合度進(jìn)行評價。

3 結(jié)果與分析

3.1 農(nóng)田磷流失源因子空間分布

海河流域內(nèi)耕地土壤有效磷質(zhì)量分?jǐn)?shù)范圍為0～193 mg/kg，均值為45 mg/kg。高值區(qū)主要分布在燕山-太行山山脈地區(qū)及渤海灣灘涂地區(qū)（圖2a），占農(nóng)田總面積的33.5%，前者可能是由于山區(qū)土壤有機(jī)質(zhì)含量高，在分解過程中會產(chǎn)生有機(jī)酸、腐殖酸等物質(zhì)，減少土壤對磷的吸附，導(dǎo)致土壤中磷含量高[33]，后者可能是由于渤海灣灘涂地區(qū)沉積物中磷含量較高，且黏粒、粉粒含量較高，對磷吸附能力較強(qiáng)[34]；中值區(qū)主要分布在海河及黃河沖積平原，且與耕地分布基本吻合，這可能是人為施肥及土壤磷素累積所致；低值區(qū)主要分布在燕山-太行山山麓地區(qū)，這可能是由于地形起伏，坡度較大，土壤侵蝕嚴(yán)重，導(dǎo)致含量低。整體上，研究區(qū)內(nèi)土壤有效磷含量中等的區(qū)域占農(nóng)田總面積為66.5%。

海河流域內(nèi)磷肥施用量范圍為0～115 kg/hm2，均值為52 kg/hm2?？臻g分布與土壤有效磷空間分布相反（圖 2b）。高值區(qū)主要分布在海河與黃河沖擊平原區(qū)域的農(nóng)耕區(qū)，其中山東省與河北省東南部磷肥施用量顯著高于其它省份，這是因為黃河及海河的沖積平原地勢平緩，適宜耕作，因此耕作強(qiáng)度大，化肥施用量高；低值區(qū)主要分布在燕山-太行山山脈地區(qū)，這主要是由于山區(qū)耕地較少，且土壤有效磷含量較高，因此磷肥施用量較低。整體上，海河流域內(nèi)磷肥施用量較高的區(qū)域占農(nóng)田總面積的61.1%。

3.2 農(nóng)田磷流失遷移因子空間分布

3.2.1 土壤侵蝕模數(shù)

海河流域降雨侵蝕力因子（）范圍為54～278 MJ·mm/(hm2·h·a)。東西向差異顯著，表現(xiàn)出東高西低的規(guī)律。高值區(qū)主要分布于河北省東北部區(qū)域，低值區(qū)主要分布于研究區(qū)西北部（圖3a）。海河流域內(nèi)值較高的區(qū)域占農(nóng)田總面積的51.9%。而土壤可蝕性因子（）范圍為0～0.5 t·hm2·h/(MJ·mm·hm2)（圖3b）。高值區(qū)主要分布于海河與黃河沖積平原及燕山-太行山山麓地區(qū)，海河流域內(nèi)值較高的區(qū)域占農(nóng)田總面積的67.0%。

地形因子（）范圍為0～362。高值區(qū)主要分布在燕山-太行山山脈地區(qū)，低值區(qū)主要分布在平原區(qū)域（圖 3c）。整體上，海河流域內(nèi)值較低的區(qū)域占農(nóng)田總面積的91.6%，絕大部分區(qū)域地勢平緩。而植被覆蓋與管理因子（）高值區(qū)主要分布于東南部及西北部（圖3d）。其中東南部區(qū)域值空間分布主要與建設(shè)用地的空間分布吻合，具有較強(qiáng)地人為擾動；西北部區(qū)域主要是農(nóng)牧交錯區(qū)，以草地為主，植被覆蓋度比平原耕地與林地低。整體上，海河流域內(nèi)值較低的區(qū)域占農(nóng)田總面積的77.9%。水土保持措施因子與地形相一致，高值區(qū)主要分布在燕山-太行山山脈地區(qū)（圖3e）。整體上，海河流域內(nèi)P值較低區(qū)域占農(nóng)田總面積的75.0%。

土壤侵蝕模數(shù)（）的分布狀況與地形因子（）及植被覆蓋與管理因子（）有關(guān)。和越大，值愈高。高值區(qū)主要分布在研究區(qū)的西北部（圖3f），其中山地區(qū)域由于地勢起伏大，易發(fā)生侵蝕；草地區(qū)域植被覆蓋度相對較低，亦易發(fā)生侵蝕。整體上，海河流域值較低的區(qū)域占農(nóng)田總面積的99.0%。

3.2.2 年徑流深和歸一化距離指數(shù)

海河流域年徑流深（R）的空間分布與降雨侵蝕力基本一致（圖4a和圖3a），R范圍為44.2～99.7 mm/a。R高值區(qū)主要分布于東部沿海地區(qū)及太行山山脈地區(qū)，前者由于沿海，降水量大，從東部沿海到內(nèi)陸呈降低趨勢，這與海陸地帶性規(guī)律相契合；后者主要是受山地丘陵的地形因素影響，在降水后因山地徑流的匯水作用致使徑流深度較多。整體上，海河流域R值較高的區(qū)域占農(nóng)田總面積的54.2%。而海河流域歸一化距離指數(shù)（NDDI）較低的區(qū)域占農(nóng)田總面積的64.8%，即絕大多數(shù)的耕地與河網(wǎng)間距離較近（圖4b）。

3.3 農(nóng)田磷流失關(guān)鍵源區(qū)及因子識別

源因子指數(shù)（Source Factor Index，SI）與磷肥施用量相一致，遷移因子指數(shù)（Transport Factor Index，TI）與歸一化距離指數(shù)（NDDI）相一致，磷指數(shù)（Phosphorus Index，PI）與TI相一致（圖5）。海河流域SI和TI空間分布圖（圖5a和5b）中農(nóng)田磷流失高風(fēng)險及以上區(qū)域分別占農(nóng)田總面積的40.2%和30.9%（表2）。海河流域農(nóng)田磷流失高風(fēng)險區(qū)域（即關(guān)鍵源區(qū)）占農(nóng)田總面積的20.5%（表2），主要分布于黃河北岸的山東省和河北省東南部的平原農(nóng)耕區(qū)、海河流域西北部的山區(qū)地帶（圖5c）。

PI、SI及TI與各源因子、遷移因子相關(guān)分析結(jié)果如表3所示，PI與年徑流深、降雨侵蝕力因子顯著正相關(guān)（<0.05），與有機(jī)碳極顯著正相關(guān)（<0.01）；SI與土壤有效磷、砂粒含量及土壤侵蝕模數(shù)、水土保持措施因子呈極顯著正相關(guān)（<0.01），與年徑流深、降雨侵蝕力因子、土壤可蝕性因子、粉粒含量、磷肥施用量極顯著負(fù)相關(guān)（<0.01），與黏粒含量顯著負(fù)相關(guān)（<0.05）；TI與磷肥施用量、年徑流深、歸一化距離指數(shù)、降雨侵蝕力因子、土壤可蝕性因子、黏粒含量呈極顯著正相關(guān)（<0.01），與粉粒含量呈顯著正相關(guān)（<0.05），與水土保持措施因子、砂粒含量和有機(jī)碳含量均呈極顯著負(fù)相關(guān)（<0.01）。

圖3 海河流域土壤侵蝕模數(shù)及其影響因子空間分布

圖5 海河流域農(nóng)田磷流失風(fēng)險指數(shù)空間分布

表2 海河流域農(nóng)田磷流失不同風(fēng)險等級面積占比

表3 農(nóng)田磷流失風(fēng)險指數(shù)與各影響因子間相關(guān)分析

注：*和**分別表示在0.05和0.01水平上顯著相關(guān)，ns表示顯著不相關(guān)（>0.05）。下同。

Note: * and ** indicate significant correlation at 0.05 and 0.01 levels, respectively. And ns represents no significant correlation (>0.05). Same as below.

結(jié)構(gòu)方程模型（SEM）配適度指標(biāo)（CHI/DF=0.304；GFI=1.000；AGFI=0.999；AIC=18.304；RMSEA=0.000）表明模型擬合效果較好（圖6）。SEM結(jié)果表明，PI與遷移因子顯著正相關(guān)（<0.05），源因子與遷移因子間呈極顯著負(fù)相關(guān)（<0.01），土壤性質(zhì)與源因子呈極顯著負(fù)相關(guān)（<0.01），與遷移因子呈極顯著正相關(guān)（<0.01）。

注：箭頭的寬度表示標(biāo)準(zhǔn)化路徑系數(shù)強(qiáng)度。實線表示正路徑系數(shù)，虛線表示負(fù)路徑系數(shù)。

4 討 論

評估農(nóng)田土壤磷流失潛在風(fēng)險對于制定有效的磷污染防控措施至關(guān)重要[7]。本研究識別出的農(nóng)田磷流失關(guān)鍵源區(qū)主要集中在海河流域東南部區(qū)域的平原農(nóng)耕區(qū)及山區(qū)河谷地帶，與孔佩儒等研究結(jié)果相吻合[35]。這是由于：1）東南部區(qū)域的平原農(nóng)耕區(qū)主要為海河與黃河沖積平原，該區(qū)域地勢平坦（圖3c），土壤肥沃，耕作強(qiáng)度及化肥投入量大（圖2b），加之農(nóng)田與河網(wǎng)間距離較近（圖 4b），在中高強(qiáng)度的降雨（徑流）作用下，土壤及肥料中磷素易隨地表徑流向受納水體遷移，導(dǎo)致農(nóng)田磷流失風(fēng)險高；2）海河流域西北部的山區(qū)地帶，由于其地勢起伏，坡度較大，土壤侵蝕較為嚴(yán)重，因此農(nóng)田磷流失風(fēng)險高（圖3c和圖3f）。

本研究中，海河流域磷指數(shù)（PI）與遷移因子顯著正相關(guān)（<0.05），而與源因子未達(dá)到顯著相關(guān)水平（圖6，>0.05），且PI與TI空間分布高度相關(guān)。這說明，海河流域農(nóng)田磷流失主要影響因子為遷移因子，這與Li等研究結(jié)果相一致[11]。從遷移因子角度來看，海河流域土壤侵蝕模數(shù)整體上處于極低水平，這與前人研究結(jié)果相一致[17]。土壤侵蝕模數(shù)值范圍較廣，這可能是由于自然條件、空間和時間尺度、生態(tài)系統(tǒng)相互作用的復(fù)雜性及土壤侵蝕地塊等因素導(dǎo)致土壤侵蝕模數(shù)的高度變異性[36]。而對于中高值區(qū)域仍需采取一定的水土保持措施，研究發(fā)現(xiàn)保護(hù)性耕作和等高種植能增強(qiáng)了土壤對高強(qiáng)度降水的適應(yīng)力，降低土壤侵蝕模數(shù)[37-38]。研究區(qū)內(nèi)年徑流深及農(nóng)田與河網(wǎng)間歸一化距離指數(shù)高值區(qū)所占比例較高，分別達(dá)31.0%和64.8%（圖4）。且TI和降雨侵蝕力因子、年徑流深及歸一化距離指數(shù)極顯著正相關(guān)（表 3，<0.01）。這說明徑流及農(nóng)田與河網(wǎng)間距離為主要遷移因子。相關(guān)研究表明，降水變化會影響磷投入與輸出間的平衡[39]，主要是由于高強(qiáng)度降水會影響年徑流深、降雨侵蝕力及磷的相對釋放。而海河流域不同區(qū)域年降水量變化較大，最大差異可達(dá)462.5 mm。本研究中降水量較少的區(qū)域（西北部）土壤有效磷含量高于高降水量區(qū)域（東南部），這可能是由于低降水區(qū)域徑流導(dǎo)致的磷損失較低[39]。有研究報道采取一定的水土保持措施及管理實踐，有助于減少水土流失及徑流作用[40]。從源因子角度來看，SI空間分布結(jié)果顯示40.2%的農(nóng)田潛在磷流失達(dá)到高風(fēng)險及以上級別（圖5a）。因此，源因子也是影響土壤磷流失的重要因子。SI空間分布與磷肥施用量一致，而處于海河與黃河沖擊平原區(qū)域的農(nóng)耕區(qū)（山東省及河北省東南部）磷肥施用量普遍過量，因此，在不影響糧食安全的前提下可以通過采取適量減肥、肥料替代及綜合養(yǎng)分管理等措施來減少農(nóng)田磷流失風(fēng)險[41]。

與此同時，本研究結(jié)果表明土壤性質(zhì)（土壤質(zhì)地及有機(jī)碳含量）與源因子、遷移因子間均呈極顯著相關(guān)（圖6，<0.01）。表3中SI與黏粒含量顯著負(fù)相關(guān)（<0.05），與粉粒含量極顯著負(fù)相關(guān)（<0.01），且SI與土壤有效磷含量、砂粒含量極顯著正相關(guān)（<0.01）。這可能是由于黏?；蚍哿：枯^高的土壤中磷淋洗量少，因此有效磷含量通常比砂粒含量高的土壤中磷含量高[34]。值得注意的是，表3中SI與磷肥施用量呈極顯著負(fù)相關(guān)，然而通過圖2b與圖5a結(jié)果對比分析，二者應(yīng)呈正相關(guān)關(guān)系，這一矛盾結(jié)果的出現(xiàn)可能是由于：1）SI與土壤有效磷呈極顯著正相關(guān)（表3，<0.01），而土壤有效磷與磷肥施用量間呈極顯著負(fù)相關(guān)（表3，<0.01）。農(nóng)業(yè)實際生產(chǎn)過程中，土壤有效磷含量較高的地區(qū)，磷肥施用量相對來說會較低。反之，土壤有效磷含量較低的地區(qū)，磷肥施用量相對來說較高。磷肥施用量一定程度上受制于土壤有效磷含量；2）海河流域范圍廣，空間異質(zhì)性較強(qiáng)，諸多地表要素所造成的不確定性因素也會對結(jié)果產(chǎn)生一定偏差。綜上兩點導(dǎo)致相關(guān)分析中出現(xiàn)這一結(jié)果矛盾現(xiàn)象。TI與土壤可蝕性因子、質(zhì)地（砂粒、黏粒和粉粒含量）、有機(jī)碳和水土保持措施因子顯著相關(guān)（表 3，<0.05），而這些指標(biāo)均與土壤性質(zhì)緊密相關(guān)。這說明土壤性質(zhì)對農(nóng)田磷流失具有重要影響[42]。大量研究報道，提升土壤質(zhì)量一定程度上也能有效減少水土流失[42-43]。土壤理化性質(zhì)會影響土壤中磷形態(tài)及分布[44]；較低的耕作頻率及土壤質(zhì)量增強(qiáng)提高了土壤的入滲能力和土壤結(jié)構(gòu)（增加團(tuán)聚體穩(wěn)定性），使土壤抗水土流失能力增強(qiáng)[45]。因此，通過一定的保護(hù)性耕作等水土保持措施提升土壤質(zhì)量，緩沖降水（徑流）影響及增加土壤粗糙度，減少在一定降雨（徑流）強(qiáng)度下農(nóng)田土壤磷流失對于農(nóng)田面源污染防治具有重要意義。

農(nóng)田磷流失的空間變異性通常歸因于多個空間變量的綜合作用，如土壤磷含量、地理位置、水文條件、降水因素及農(nóng)業(yè)耕作、管理措施等。一些變量間的協(xié)同作用會加劇磷流失，而另一些變量間的拮抗作用會減少磷流失。如山區(qū)地形和降水因子通常會導(dǎo)致高水土流失和磷流失風(fēng)險，而茂密植被又可能會減少土壤侵蝕，同時低施肥量也可以通過降低磷盈余來減少磷流失潛力，農(nóng)田與河網(wǎng)間距離越遠(yuǎn)越能通過提供更多的磷滯留機(jī)會來減少磷流失至受納水體。平原區(qū)域地勢起伏較小、坡度平緩，具有較低的侵蝕潛力，然而由于磷投入量較大，在中高強(qiáng)度的降雨（徑流）作用下，仍面臨高風(fēng)險級的磷流失[15]，這與本研究中結(jié)果相一致。本研究結(jié)果發(fā)現(xiàn)，源因子與遷移因子間呈極顯著負(fù)相關(guān)關(guān)系（圖6，<0.01），這可能是由于：1）磷在農(nóng)田中為作物必需養(yǎng)分，并非污染物；2）影響海河流域農(nóng)田磷流失的主導(dǎo)因子為遷移因子，而年徑流深及歸一化距離指數(shù)作為遷移因子中的關(guān)鍵因子。地表徑流強(qiáng)度越大、農(nóng)田與河網(wǎng)間距離越近，農(nóng)田土壤磷流失風(fēng)險越高，越易使土壤及肥料中的磷遷移到受納水體造成污染，進(jìn)而降低SI。此外，在TI與SI雙高區(qū)域如本研究中的關(guān)鍵源區(qū)（黃河北岸山東省區(qū)域、河北省東南部的平原農(nóng)耕區(qū)），農(nóng)田磷流失受源因子及遷移因子的綜合作用，因此需要整合該區(qū)域關(guān)鍵的源和遷移因子并進(jìn)行綜合考慮，制定出TI與SI雙減的最佳管理措施以減少農(nóng)田磷流失[46]。綜上，對于海河流域農(nóng)田磷流失關(guān)鍵源區(qū)中的平原農(nóng)耕區(qū)，可通過采取適量減肥、綜合養(yǎng)分管理、保護(hù)性耕作及在農(nóng)田與河網(wǎng)間布設(shè)緩沖區(qū)等措施來降低農(nóng)田磷流失風(fēng)險；對于山區(qū)地帶，可采取如等高種植、梯田種植等水土保持措施來減少土壤侵蝕，進(jìn)而降低農(nóng)田磷流失風(fēng)險。

本研究風(fēng)險評估采用經(jīng)典成熟的磷指數(shù)模型，并在參考國內(nèi)外研究的基礎(chǔ)上進(jìn)行了一定的修正，如改進(jìn)了磷指數(shù)模型中農(nóng)田到水體間距離因子，首次提出歸一化距離因子的概念，并將源因子、遷移因子與磷指數(shù)歸一化后進(jìn)行相關(guān)分析、建立結(jié)構(gòu)方程模型，研究了土壤性質(zhì)（如質(zhì)地及有機(jī)碳）、源因子、遷移因子及農(nóng)田磷流失風(fēng)險指數(shù)間關(guān)系，得出最終研究結(jié)果。由于磷指數(shù)模型中各因子權(quán)重的確定存在一定主觀性，權(quán)重賦值的差異會顯著影響PI數(shù)值大小，但對于定性化識別關(guān)鍵源區(qū)影響不大。研究過程中由于數(shù)據(jù)受限，在源因子中未考慮施肥時間及方式，這些可能會給研究結(jié)果帶來一定偏差，因此，對于各因子及磷流失風(fēng)險評價有待結(jié)合實地情況進(jìn)一步研究。本研究數(shù)據(jù)均來源于公共平臺，易獲取，這為流域尺度上直觀地識別農(nóng)田磷流失關(guān)鍵源區(qū)和因子，進(jìn)而采取有針對性的防控措施提供了快速實用的方法。

5 結(jié) 論

本研究利用磷指數(shù)模型評估了海河流域農(nóng)田磷流失環(huán)境風(fēng)險，識別出海河流域磷流失的關(guān)鍵源區(qū)、因子。研究結(jié)果對流域尺度上農(nóng)田磷流失研究方法創(chuàng)新有所裨益，為海河流域農(nóng)田面源污染防治提供科學(xué)參考。本研究得出結(jié)論如下：

1）海河流域農(nóng)田土壤有效磷、磷肥施用量、土壤侵蝕模數(shù)、年徑流深及歸一化距離指數(shù)的空間分布結(jié)果表明，大多數(shù)區(qū)域分別處于中-低、中-高、極低、中-高和高風(fēng)險級別，分別占農(nóng)田總面積的66.5%、61.1%、99.0%、54.2%和64.8%；

2）影響海河流域農(nóng)田磷流失的關(guān)鍵因子為遷移因子，其中關(guān)鍵遷移因子為年徑流深及歸一化距離指數(shù)；

3）源因子與遷移因子間呈極顯著負(fù)相關(guān)（<0.01），土壤性質(zhì)（包括土壤質(zhì)地及有機(jī)碳含量等）與源因子呈極顯著負(fù)相關(guān)（<0.01），與遷移因子間呈極顯著正相關(guān)（<0.01）；

4）海河流域農(nóng)田磷流失關(guān)鍵源區(qū)主要位于黃河北岸的山東省和河北省東南部的平原農(nóng)耕區(qū)、海河流域西北部的山區(qū)地帶，需要重點防控。

[1] Werner F, de la Haye T R, Spielvogel S, et al. Small-scale spatial distribution of phosphorus fractions in soils from silicate parent material with different degree of podzolization[J]. Geoderma, 2017, 302: 52-65.

[2] Granger S J, Harris P, Peukert S, et al. Phosphate stable oxygen isotope variability within a temperate agricultural soil[J]. Geoderma, 2017, 285: 64-75.

[3] Schoumans O F, Chardon W J, Bechmann M E, et al. Mitigation options to reduce phosphorus losses from the agricultural sector and improve surface water quality: A review[J]. Science of the Total Environment, 2014, 468/469: 1255-1266.

[4] Cao D, Cao W, Fang J, et al. Nitrogen and phosphorus losses from agricultural systems in China: A meta-analysis[J]. Marine Pollution Bulletin, 2014, 85(2): 727-732.

[5] Conley D J, Paerl H W, Howarth R W, et al. Controlling eutrophication: Nitrogen and phosphorus[J]. Science, 2009, 323(5917): 1014-1015.

[6] Sun B, Zhang L, Yang L, et al. Agricultural non-point source pollution in China: Causes and mitigation measures[J]. Ambio, 2012, 41(4): 370-379.

[7] Huang L, Wang H, Li Y, et al. Spatial distribution and risk assessment of phosphorus loss potential in urban–suburban soil of Lishui, China[J]. Catena, 2013, 100: 42-49.

[8] Wang F E, Sun Z L, Zheng S Y, et al. An integrated approach to identify critical source areas of agricultural nonpoint-source pollution at the watershed scale[J]. Journal of Environmental Quality, 2018, 47(4): 922-929.

[9] Gburek W J, Sharpley A N, Heathwaite L, et al. Phosphorus management at the watershed scale: A modification of the phosphorus index[J]. Journal of Environment Quality, 2000, 29(1): 130-144.

[10] Sharpley A N, Weld J L, Beegle D B, et al. Development of phosphorus indices for nutrient management planning strategies in the United States[J]. Journal of Soil and Water Conservation, 2003, 58(3): 137-152.

[11] Li Z, Zhang R, Liu C, et al. Phosphorus spatial distribution and pollution risk assessment in agricultural soil around the Danjiangkou reservoir, China[J]. Science of the Total Environment, 2020, 699: 134417.

[12] 劉潔，龐樹江，何楊洋，等. 基于小流域產(chǎn)流特征的磷流失關(guān)鍵源區(qū)識別[J]. 農(nóng)業(yè)工程學(xué)報，2017，33(20)：241-249. Liu Jie, Pang Shujiang, He Yangyang, et al. Critical area identification of phosphorus loss based on runoff characteristics in small watershed[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2017, 33(20): 241-249. (in Chinese with English abstract)

[13] Buchanan B P, Archibald J A, Easton Z M, et al. A phosphorus index that combines critical source areas and transport pathways using a travel time approach[J]. Journal of Hydrology, 2013, 486: 123-135.

[14] 張淑榮，陳利頂，傅伯杰，等. 農(nóng)業(yè)區(qū)非點源污染潛在危險性評價：以于橋水庫流域磷流失為例[J]. 第四紀(jì)研究，2003，23(3)：262-269. Zhang Shurong, Chen Liding, Fu Bojie, et al. The risk assessment of nonpoint pollution of phosphorus from agricultural lands: A case study of Yuqiao Reservoir Watershed[J]. Quaternary Sciences, 2003, 23(3): 262-269. (in Chinese with English abstract)

[15] Zhou H, Gao C. Assessing the risk of phosphorus loss and identifying critical source areas in the Chaohu Lake Watershed, China[J]. Environmental Management, 2011, 48(5): 1033-1043.

[16] 李文超，劉申，雷秋良，等. 高原農(nóng)業(yè)流域磷流失風(fēng)險評價及關(guān)鍵源區(qū)識別：以鳳羽河流域為例[J]. 農(nóng)業(yè)環(huán)境科學(xué)學(xué)報，2014，33(8)：1591-1600. Li Wenchao, Liu Shen, Lei Qiuliang, et al. Assessing risk of non-point source phosphorus loss and identifying critical source areas in a Chinese highland agricultural watershed[J]. Journal of Agro-Environment Science, 2014, 33(8): 1591-1600. (in Chinese with English abstract)

[17] 李曉松，吳炳方，王浩，等. 區(qū)域尺度海河流域水土流失風(fēng)險評估[J]. 遙感學(xué)報，2011，15(2)：372-387. Li Xiaosong, Wu Bingfang, Wang Hao, et al. Regional soil erosion risk assessment in Hai Basin[J]. Journal of Remote Sensing, 2011, 15(2): 372-387. (in Chinese with English abstract)

[18] 萬煒，李含微，王佳瑩，等. 基于空間平滑法的旱作區(qū)糧食產(chǎn)量時空變化與影響因素研究[J]. 農(nóng)業(yè)工程學(xué)報，2019，35(16)：284-296. Wan Wei, Li Hanwei, Wang Jiaying, et al. Spatio-temporal changes and influencing factors of grain yield based on spatial smoothing method in dryland farming regions[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2019, 35(16): 284-296. (in Chinese with English abstract)

[19] Pribyl D W. A critical review of the conventional SOC to SOM conversion factor[J]. Geoderma, 2010, 156(3/4): 75-83.

[20] 魯如坤. 土壤磷素水平和水體環(huán)境保護(hù)[J]. 磷肥與復(fù)肥，2003，18(1)：4-8. Lu Rukun. The phosphorus level of soil and environmental protection of water body[J]. Phosphate & Compound Fertilizer, 2003, 18(1): 4-8. (in Chinese with English abstract)

[21] 萬煒，鄧靜，王佳瑩，等. 基于潛力衰減模型的東北-華北平原旱作區(qū)耕地生產(chǎn)力評價[J]. 農(nóng)業(yè)工程學(xué)報，2020，36(5)：270－280. Wan Wei, Deng Jing, Wang Jiaying, et al. Evaluation of cultivated land productivity based on potential attenuation model in the dryland farming regions of Northeast and North China Plain[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2020, 36(5): 270－280. (in Chinese with English abstract)

[22] Sharpley A N, Mcdowell R W, Peter J A K. Phosphorus loss from land to water: Integrating agricultural and environmental management[J]. Plant and Soil, 2001, 237(2): 287-307.

[23] Wischmeier W H, Smith D D. Predicting Rainfall Erosion Losses. A guide to conservation planning[M]. Washington: United States Department of Agriculture, 1978.

[24] Williams J, Renard K, Dyke P. EPIC: A new method for assessing erosion’s effect on soil productivity[J]. Journal of Soil and Water Conservation, 1983, 38: 381-383.

[25] Phinzi K, Ngetar N S. The assessment of water-borne erosion at catchment level using GIS-based RUSLE and remote sensing: A review[J]. International Soil and Water Conservation Research, 2019, 7(1): 27-46.

[26] 蔡崇法，丁樹文，史志華，等. 應(yīng)用USLE模型與地理信息系統(tǒng)IDRISI預(yù)測小流域土壤侵蝕量的研究[J]. 水土保持學(xué)報，2000，14(2)：19-24. Cai Chongfa, Ding Shuwen, Shi Zhihua, et al. Study of applying USLE and geographical information system IDRISI to predict soil erosion in small watershed[J]. Journal of Soil and Water Conservation, 2000, 14(2): 19-24. (in Chinese with English abstract)

[27] Fu B J, Zhao W W, Chen L D, et al. Assessment of soil erosion at large watershed scale using RUSLE and GIS: A case study in the Loess Plateau of China[J]. Land Degradation & Development, 2005, 16(1): 73-85.

[28] Lufafa A, Tenywa M M, Isabirye M, et al. Prediction of soil erosion in a Lake Victoria basin catchment using a GIS-based Universal Soil Loss model[J]. Agricultural Systems, 2003, 76(3): 883-894.

[29] 任憲韶，戶作亮，曹寅白. 海河流域水利手冊[M]. 北京：中國水利水電出版社，2008.

[30] 劉春蓁，劉志雨，謝正輝. 近50年海河流域徑流的變化趨勢研究[J]. 應(yīng)用氣象學(xué)報，2004，15(4)：385-393. Liu Chunzhen, Liu Zhiyu, Xie Zhenghui. Study on trends in runoff for the Haihe River Basin in recent 50 years[J]. Journal of Applied Meteorological Science, 2004, 15(4): 385-393. (in Chinese with English abstract)

[31] 張世文，王朋朋，葉回春，等. 基于數(shù)字土壤系統(tǒng)的縣域土壤磷素流失風(fēng)險簡化評估[J]. 農(nóng)業(yè)工程學(xué)報，2012，28(11)：110-117. Zhang Shiwen, Wang Pengpeng, Ye Huichun, et al. Risk simplified assessment on phosphorus loss risk based on digital soil system at county scale[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2012, 28(11): 110-117. (in Chinese with English abstract)

[32] Bechmann M E, St?lnacke P, Kv?rn? S H. Testing the Norwegian phosphorus index at the field and subcatchment scale[J]. Agriculture, Ecosystems & Environment, 2007, 120(2/3/4): 117-128.

[33] Chacón N, Dezzeo N, Mu?oz B, et al. Implications of soil organic carbon and the biogeochemistry of iron and aluminum on soil phosphorus distribution in flooded forests of the lower Orinoco River, Venezuela[J]. Biogeochemistry, 2005, 73(3): 555-566.

[34] Momtaz H R, Ghasemi N, Sepehr E, et al. Spatial distribution of soil-available phosphorus and single-point phosphorus buffering index in the Khoy region, Iran[J]. Soil Use and Management, 2017, 33(1): 13-24.

[35] 孔佩儒，陳利頂，孫然好，等. 海河流域面源污染風(fēng)險格局識別與模擬優(yōu)化[J]. 生態(tài)學(xué)報，2018，38(12)：4445-4453. Kong Peiru, Chen Liding, Sun Ranhao, et al. Identification and optimized simulation of non-point source pollution risk patterns in the Haihe River Basin[J]. Acta Ecologica Sinica, 2018, 38(12): 4445-4453. (in Chinese with English abstract)

[36] Boix F C, Martínez M M, Calvo C A, et al. Concise review of interrill erosion studies in SE Spain (Alicante and Murcia): Erosion rates and progress of knowledge from the 1980s[J]. Land Degradation & Development, 2005, 16(6): 517-528.

[37] Liu R M, Zhang P P, Wang X J, et al. Assessment of effects of best management practices on agricultural non-point source pollution in Xiangxi River watershed[J]. Agricultural Water Management, 2013, 117: 9-18.

[38] Beniston J W, Shipitalo M J, Lal R, et al. Carbon and macronutrient losses during accelerated erosion under different tillage and residue management[J]. European Journal of Soil Science, 2015, 66(1): 218-225.

[39] Liu Z, Shao M, Wang Y. Spatial patterns of soil total nitrogen and soil total phosphorus across the entire Loess Plateau region of China[J]. Geoderma, 2013, 197/198: 67-78.

[40] Eekhout J P C, de Vente J. Assessing the effectiveness of Sustainable Land Management for large-scale climate change adaptation[J]. Science of the Total Environment, 2019, 654: 85-93.

[41] Heathwaite L, Sharpley A, Gburek W. A conceptual approach for integrating phosphorus and nitrogen management at watershed scales[J]. Journal of Environment Quality, 2000, 29(1): 158-166.

[42] Martínez-Mena M, Carrillo-López E, Boix-Fayos C, et al. Long-term effectiveness of sustainable land management practices to control runoff, soil erosion, and nutrient loss and the role of rainfall intensity in Mediterranean rainfed agroecosystems[J]. Catena, 2020, 187: 104352.

[43] Ouyang W, Huang H B, Hao F H, et al. Evaluating spatial interaction of soil property with non‐point source pollution at watershed scale: The phosphorus indicator in Northeast China[J]. Science of the Total Environment, 2012, 432: 412-421.

[44] Bechmann M, Stalnacke P, Kvoerno S, et al. Integrated tool for risk assessment in agricultural management of soil erosion and losses of phosphorus and nitrogen[J]. Science of the Total Environment, 2009, 407(2): 749-759.

[45] Haruna S I, Anderson S H, Nkongolo N V, et al. Soil hydraulic properties: Influence of tillage and cover crops[J]. Pedosphere, 2018, 28(3): 430-442.

[46] Rao N S, Easton Z M, Schneiderman E M, et al. Modeling watershed-scale effectiveness of agricultural best management practices to reduce phosphorus loading[J]. Journal of Environmental Management, 2009, 90(3): 1385-1395.

Environmental risk assessment of phosphorus loss from farmland based on phosphorus index model in the Haihe River Basin

Wu Hanqing1, Wan Wei1, Shan Yanjun1, Chen Yanhua2, Li Qiang1, Li Cunjun3, Hu Haitang3, Zhang Baogui1※

(1.,,100193,; 2.,,100097,; 3.,100097,)

The non-point phosphorus (P) loss from farmland is the main cause of eutrophication, and identifying the Critical Source Areas (CSAs) and Critical Influence Factors (CIFs) of the P loss from farmland is a pivotal part of controlling agricultural non-point source pollution. The objective of this study was to assess the environmental risk of P loss from farmland based on P index model in the Haihe River Basin. In this study, we took the Haihe River Basin as the research region, and applied the P index model. Thereinto, the soil available P content and fertilizer-P application rate were selected as the source factors. Synchronously, the soil erosion modulus, annual runoff, the Normalized Differential Distance Index (NDDI) between farmland and river network were chosen as the transport factors. Additionally, this research aims at identifying the CSAs combined with GIS technology, distinguishing the CIFs and exploring the relation between P index and influence factors such as soil properties (mainly including soil texture and organic carbon), source and transport factors via Structural Equation Modeling (SEM). The spatial distributions of soil available P content, fertilizer-P application rate, soil erosion modulus, annual runoff and NDDI in the Haihe River Basin indicated that the most areas were at the risk levels of low-medium, medium-high, very low, medium-high and high in the Haihe River Basin, respectively. Namely, these areas accounted for 66.5%, 61.1%, 99.0%, 54.2%, and 64.8% of the total farmland area in the Haihe River Basin, respectively. Moreover, we revealed that the CIF affecting P loss from farmland in the Haihe River Basin is the transport factors in which the key transport factors are the annual runoff and NDDI. Furthermore, SEM prescribed that there was an extremely significant negative correlation between the source factors (<0.01) and extremely significant positive correlation with transport factors (<0.01), and soil properties had highly significant negative correlations with source factors and transport factors (<0.01). Subsequently, the findings illustrated that the CSAs of P loss from farmland in the Haihe River Basin were mainly located in the plain farming areas of Shandong Province and the southeast of Hebei Province on the north bank of the Yellow River, and the mountain areas in the northwestern of the Haihe River Basin. Due to the plain farming areas have high cultivation intensity and fertilizer-P application rate, flat terrain, and close distance to river network, P of soil and fertilizer could easily migrate into the rivers following the surface runoff with the medium-high intensity. Thus, for the plain farming areas, some practices should be taken to decrease the risk level of P loss from farmland, such as appropriate reduction of the fertilizer-P application rate, comprehensive nutrient management, conservation tillage and deployment of buffer zone between farmland and river network and so on. Because mountain areas have severe soil erosion, therefore, the risk level of P loss from farmland is relatively high. For this, proper soil and water conservation measures should be taken to reduce the soil erosion and risk level of P loss from farmland, such as contour tillage and terrace cropping and so on. In conclusion, the results are beneficial to the method innovation of P loss from farmland at the watershed scale, and can provide a scientific reference for controlling agricultural non-point source pollution in the Haihe River Basin.

farmland; models; non-point source pollution; phosphorus loss; phosphorus index model; critical source area; structural equation modeling; Haihe River Basin

吳漢卿，萬煒，單艷軍，等. 基于磷指數(shù)模型的海河流域農(nóng)田磷流失環(huán)境風(fēng)險評價[J]. 農(nóng)業(yè)工程學(xué)報，2020，36(14)：17-27.doi：10.11975/j.issn.1002-6819.2020.14.003 http://www.tcsae.org

Wu Hanqing, Wan Wei, Shan Yanjun, et al. Environmental risk assessment of phosphorus loss from farmland based on phosphorus index model in the Haihe River Basin[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2020, 36(14): 17-27. (in Chinese with English abstract) doi：10.11975/j.issn.1002-6819.2020.14.003 http://www.tcsae.org

2020-06-22

2020-07-10

國家重點研發(fā)計劃項目（2016YFD0800906）

吳漢卿，博士生，主要從事農(nóng)田面源氮磷污染研究。Email：wuhanqing920705@163.com

張寶貴，教授，博士生導(dǎo)師，主要從事農(nóng)田面源污染和土壤重金屬污染對生態(tài)系統(tǒng)服務(wù)影響研究。Email：zhangbg@cau.edu.cn

10.11975/j.issn.1002-6819.2020.14.003

X52；K903

A

1002-6819(2020)-14-0017-11