多情景分析的農業面源污染關鍵源區識別軟件開發及應用

覃苑，胡海棠，淮賀舉，李存軍*，張巧玲，楊鐵利，王佳宇

1.遼寧科技大學

2.北京市農林科學院信息技術研究中心

農業面源污染（agricultural non-point source pollution，ANSP）指由于農業生產中的氮磷營養元素、農藥及其他溶解態和固態污染物流失導致的污染。近年來，農業面源污染對水體污染的貢獻越來越受到重視。2021年3月生態環境部和農業農村部聯合發布了《農業面源污染治理與監督指導實施方案（試行）》，明確了“十四五”至2035年農業面源污染防治的總體要求、工作目標，提出主要任務包括深入推進農業面源污染防治、完善農業面源污染防治政策機制和加強農業面源污染治理監督管理。ANSP的污染負荷估算、關鍵源區源頭識別和治理決策是監督管理的重要內容。

ANSP在空間上分布分散且隱蔽，其形成受到區域的氣候、地形、土地利用類型和植被覆蓋等多方面因素的影響，導致ANSP在監測、模擬和控制等方面面臨較大挑戰。GIS技術在空間數據管理、集成、分析方面有優勢，已被廣泛應用于ANSP的數據管理、空間分析、面源污染模型集成、信息系統開發和可視化等方面[1]。國內外采用年度農業面源污染模型（annualized agricultural nonpoint source pollution model，AnnAGNPS），水土評價模型（soil and water assessment tool，SWAT），輸出系數模型（export coefficient model，ECM）等模型和平均濃度法等方法進行ANSP模擬與識別研究。如夏昊[2]結合AnnAGNPS模型，開發了基于GIS的面源信息系統，其擁有空間分析和可視化表達等功能，用于管理ANSP中不同來源的空間和屬性數據；喬衛芳等[3]基于研究區域數字高程模型（digital elevation model，DEM）數據，將GIS與SWAT模型結合構建了ANSP基礎信息庫，分析了研究區ANSP負荷的分布特征，并模擬了區域內主要土地利用類型單位面積ANSP的負荷量。但AnnAGNPS和SWAT等模型所需參數多，計算復雜，在大區域的模擬上略顯乏力；生態系統服務和權衡的綜合評估（InVEST）模型中的營養物傳輸率模型（NDR）所需參數少，機理清晰，適用于大區域的ANSP模擬研究[4]。

在進行ANSP治理時，關鍵源區識別是面源污染管理的重要組成部分。研究表明，少數區域輸出的污染物通常是整個流域污染負荷的主要貢獻者，對水體環境質量有著至關重要的作用[5]。《2015年中國環境狀況公報》顯示，海河、巢湖等流域的水質狀況存在明顯的空間異質性，表現為支流水質在重度污染的情況下，其所屬的干流僅為輕度污染[6]。目前主要是針對面源污染負荷空間分布的異質性特征進行ANSP關鍵源區識別[7]，但是單純考慮污染負荷易忽略流域面積與產水量對污染負荷貢獻的影響，導致關鍵源區的識別準確度下降，進而影響治理效果。應結合潛在氮、磷徑流濃度和入河負荷的空間分布進行ANSP防控的關鍵源區識別與分析[8]。目前農業面源污染關鍵源區識別涉及水文分析、農業面源污染模型、柵格矢量運算和統計，對于面源污染非專業人士和基層工作者存在跨平臺操作、數據類型多、資料匯總繁瑣、計算量大等問題，亟需集成式的軟件工具平臺。

ANSP包括農田種植業面源污染、畜禽養殖、農村生活污染等，其中農田種植業面源污染是農業面源污染的重要組成部分，污染源相對分散且影響因素多。筆者針對在農田種植業氮、磷化肥施用形成的面源污染關鍵源區識別過程中情景單一、數據量大、涉及環節多，存在計算復雜、操作繁瑣等問題，將GIS與InVEST模型結合，研發農業面源污染關鍵源區識別與治理模擬的集成計算機軟件系統，在多種情景下對關鍵源區進行識別，并模擬治理后的效果，以期為農田種植業面源污染治理提供有效的輔助手段。

1 ANSP關鍵源區識別軟件系統設計

基于GIS二次開發，集成水文網絡拓撲構建和面源污染關鍵源區識別技術，利用DEM數據提取河流節點、河段與子流域等水文地貌要素，整合InVEST產水量模型和營養物傳輸率模型，結合基于入河污染負荷、潛在徑流污染物濃度和負荷與產水量比值3種情景的面源污染關鍵源區識別方法，以決策支持為重點，建立農業面污染關鍵源區識別與治理模擬軟件系統。

1.1 系統目標

通過研究ANSP的污染負荷情況與成因，運用GIS技術，研制具有農業面源污染關鍵源區識別與治理模擬的輔助決策軟件系統。對研究區域內的水文網絡進行拓撲構建、分析和專題信息可視化表達，快速高效得到研究區內拓撲關系、子流域營養元素分布規律和子流域潛在營養元素徑流濃度等信息。在此基礎上，根據不同的方法分析面源污染關鍵源區空間位置，并模擬關鍵源區治理后的效果，為ANSP治理提供依據。

1.2 功能設計

以靈活的應用集成接口構建GIS軟件系統，對研究區內的水網數據、DEM數據和土壤數據等多源異構的資源環境數據進行集成與調度，根據研究區已知的地理數據和資源環境數據，識別面源污染關鍵源區。系統按照實際功能分為6個模塊單元：1)拓撲關系分析單元，用于基于目標區域的DEM數據構建水文網絡拓撲關系；2)產水量分析單元，用于利用InVEST模型的產水量子模型，根據水文網絡拓撲關系，確定子流域產水量；3)入河負荷估算單元，根據地表營養物負荷、傳輸率、次表層營養物傳輸率和各柵格單元的入河負荷，計算子流域營養物入河總負荷；4)潛在徑流濃度分析單元，用于根據子流域產水量，確定目標區域內各子流域的潛在氮、磷徑流濃度；5)關鍵源區識別單元，基于氮、磷入河負荷，潛在氮、磷徑流濃度，子流域氮、磷負荷與產水量比值3種情景，識別和確定目標區域內的面源污染關鍵源區；6)關鍵源區治理模擬單元，用于模擬對識別出的關鍵源區進行治理后的結果，確定治理效果。

ANSP關鍵源區識別與治理模擬統一建模語言如圖1所示。由圖1可知，1種關鍵源區識別結果對應1個關鍵源區識別模型，模擬治理結果依賴于關鍵源區識別結果；1個關鍵源區識別模型對應若干農田氮、磷入河負荷估算模型，子流域負荷與產水量比值估算模型與潛在氮、磷徑流濃度估算模型；農田氮、磷入河負荷估算模型包含1個氮、磷流失負荷分配結果，1組地表營養物傳輸率和1組次地表營養物傳輸率；子流域負荷與產水量比值估算模型包含1個產水量模型，依賴于農田氮、磷入河負荷估算模型結果；潛在氮、磷徑流濃度估算模型包含1個產水量模型，依賴于水網拓撲關系和農田氮、磷入河負荷估算模型結果。ANSP關鍵源區識別與治理模擬軟件結構如圖2所示。

圖1 ANSP關鍵源區識別與治理模擬統一建模語言Fig.1 Unified modelling language for the identification and control simulation of ANSP key source areas

圖2 ANSP關鍵源區識別與治理模擬軟件結構Fig.2 Software structure diagram for identification and simulation of the treatment of ANSP key source areas

2 ANSP關鍵源區識別主要方法與模型

2.1 水文網絡拓撲構建

獲取研究區的DEM數據，基于DEM數據建立水文網絡拓撲關系。利用計算機編程提取河段以及子流域，徑流節點是基于幾何網絡的方法提取的，分為源頭節點、匯流節點、流域出口節點[9]。具體流程見圖3。

圖3 水文地貌要素提取Fig.3 Extraction of hydrological and geomorphic elements

流域水文網絡拓撲結構的建立以流域匯流關系建立的基礎——河段為紐帶，對河網的徑流節點、河段、子流域進行矢量要素的空間分析，構建拓撲關系。河段屬性表里記錄對應河段的起始節點和終止節點，以及對應的子流域編碼；子流域屬性表里記錄與其相對應的河段編碼。定義節點的編碼為Point-ID，節點類型是NodeType，上下游節點分別為FROM_NODE、TO_NODE；河段編碼為 River_ID，長度為River_len，匯入和匯出河段分別為Riv_in、Riv_out；子流域編碼為Basin_ID，每個子流域的匯入和匯出流域編碼分別為Basin_in、Basin_out。

2.2 流域氮、磷入河負荷計算方法

流域氮、磷入河負荷的計算分為2步，首先計算出流域氮、磷流失負荷，再將其與InVEST模型中的NDR模型進行結合，計算得到流域氮、磷入河負荷。

2.2.1 氮、磷流失負荷模型

流域氮、磷流失負荷計算公式如下：

2.2.2 氮、磷入河負荷模型

通過模擬區域內營養物的來源和轉移過程，計算氮、磷入河系數，在考慮氮、磷流失負荷的情況下估算氮、磷入河負荷。首先按照式（2）將氮、磷流失負荷分配到表層和次表層：

最后，匯總計算氮、磷入河負荷，計算公式如下：

2.3 流域產水量計算方法

利用InVEST模型的產水量子模型計算流域產水量，計算公式如下：

為了反映流域多年平均產水量，降低年間的氣候差異和降水量差異造成不確定性的影響，潛在的蒸發散失量與降水量采用2001—2015年的平均值。

2.4 ANSP關鍵源區的識別

由于面源污染易受到地形、植被、氣候、施肥和土地利用類型等因素影響[12-15]，不同識別方法通常會呈現出不同程度的空間差異性。本研究采用3種不同情景識別關鍵源區。

2.4.1 基于入河負荷

該情景由分區統計獲得子流域氮、磷入河負荷，根據每個子流域氮、磷入河負荷大小進行關鍵源區的識別[16-19]，計算公式如下：

2.4.2 基于潛在徑流濃度

該情景通過計算河段氮、磷入河通量與徑流量的比值，得到河段潛在徑流濃度，根據潛在徑流濃度的貢獻率大小進行關鍵源區的識別[8]。

2.4.3 基于負荷與產水量比值

將子流域產水量與氮、磷入河負荷同時納入考慮，計算子流域產水量與氮、磷入河負荷的比值，根據該比值進行關鍵源區的識別[20-22]，計算公式如下：

2.4.4 關鍵源區識別

根據式（8）～式（11）計算得到入河負荷、潛在徑流濃度和負荷與產水量比值，對子流域或河段對應子流域進行降序排列。按順序累加子流域入河負荷，每累加1次就通過子流域累計入河負荷與流域總負荷之比（Per）來判斷是否滿足設定閾值，若大于等于設定閾值則停止累加并將這些區域標記為面源污染關鍵源區，流程如圖4所示。

圖4 關鍵源區識別流程Fig.4 Identification process for key source area

2.5 關鍵源區治理模擬

針對3種情境識別出的關鍵源區，設定相同的預期治理幅度，可直觀地觀察到在不同情景下對關鍵源區治理的效果。將3種情景下識別出的關鍵源區，按照與河流源頭的接近程度進行排序，根據設定治理幅度依次降低這些區域的氮或磷負荷，再從上游—下游的順序依次更新流域所有受影響的河段的氮或磷潛在徑流濃度。

3 農業面源污染關鍵源區識別軟件應用

3.1 研究區概況及特征

海河流域（112°E～120°E，30°N～43°N）地處我國華北地區，包括北京、天津、河南北部、山東北部、山西東部、內蒙古和遼寧小部分區域與河北省大部分區域。流域地勢為西北高、東南低，包含灤河、北三河、大清河、永定河、子牙河、漳衛南運河、黑龍港運東河、徒駭馬頰河與冀東沿海八大水系[23]。海河流域是我國三大糧食生產基地之一，流域土地利用類型主要為耕地，占流域總面積的48.9%。ANSP是海河流域主要污染源之一，也是海河水質污染的主要原因之一[24]。

3.2 總氮關鍵源區識別

研發的農業面源污染關鍵源區識別軟件主界面見圖5。以TN為例，使用該軟件對海河流域的農業面源污染關鍵源區進行分析。流域ANSP的面積與TN、TP的關系通常是呈對數上升趨勢，往往不到10%的面積就貢獻了30%以上的氮、磷污染負荷[25]，據此將軟件識別幅度設定為30%，分別基于入河污染負荷、潛在徑流污染物濃度、負荷與產水量比值3種情景進行TN面污染關鍵源區識別，結果如圖6、圖7所示。

圖5 農業面源污染關鍵源區識別軟件用戶界面Fig.5 User interface of the key source area identification software for ANSP

由圖6、圖7可知，TN入河負荷較高的區域主要分布在灤河、北三河、璋衛河水系與黑龍港運東河下游，子牙河、大清河水系中部以及徒駭馬頰河全域，呈現自西北向東南方向增長的分布特征。潛在TN徑流濃度升高的區域主要分布中部—東南部地區的黑龍港運東河、徒駭馬頰河與大清河水系。TN負荷與產水量比值高的區域集中分布在中部—南部地區的徒駭馬頰河上游、大清河水系及黑龍港運東河。基于入河負荷情景下識別出的關鍵源區較為分散，除了永定河與子牙河水系外皆有分布。基于潛在徑流濃度情景下的關鍵源區分布較為集中，多為中部—南部分布，其余在東南部也有分布。基于流域負荷產水比值情景下識別出的關鍵源區分布非常集中，以中部—南部分布為主，識別出的區域與基于潛在徑流濃度情景下的區域重合度較高。

圖6 3種情景下的TN污染區域分布情況Fig.6 Regional distribution of total nitrogen pollution under three scenarios

圖7 3種情景下識別的關鍵源區Fig.7 Key source areas identified by three scenarios

3.3 治理模擬結果與分析

在識別出關鍵源區的基礎上，對關鍵源區進行模擬治理。為了更好地反映治理效果，將識別出的關鍵源區TN、TP入河負荷清零，并在考慮水網的拓撲關系下，對河段TN、TP潛在徑流濃度進行重新統計。模擬結果如圖8所示。由圖8可知，基于入河負荷情景下，模擬治理的效果欠佳，流域中部水系的潛在徑流濃度下降不明顯；基于潛在徑流濃度、負荷與產水量比值2種情景下，治理效果相近，治理區域的潛在徑流濃度均有較為明顯的下降，但對東南部和南部的部分水系治理效果不明顯。

圖8 原始徑流濃度與3種情景下治理效果對比Fig.8 Comparison chart of original runoff concentration and that of three methods after treatment

計算了3種情景識別出的面污染關鍵源區模擬治理后TN、TP負荷的下降率，結果見表1。

由表1可知，在研究區的數據與條件下，基于入河負荷情景識別出的面污染關鍵源區的TN治理效果弱于基于潛在徑流濃度和基于負荷產水量比值情景的治理效果。由于海河流域西北地區入河負荷高，但林地和草地面積較大，攔截和稀釋了大量的氮、磷污染物，而中南部地區耕地面積占比高，且采取了一年兩熟的種植模式，化肥的施用量大且耗水嚴重，導致地表徑流量偏低，所以3種情景下識別出的關鍵源區有著較顯著的空間差異性。由于基于潛在徑流濃度情景反映的是河段與其上游輸入污染物的匯集和徑流過程共同作用的結果，而子流域負荷與產水量比值反映子流域產水量與氮、磷污染負荷比值，二者更能反映出ANSP對水質的影響程度。

表1 3種情景下關鍵源區模擬治理效果Table 1 Simulation of governance effects in key source areas under three scenarios

對海河流域的治理模擬結果表明，僅針對面源污染入河負荷識別出的氮和磷污染關鍵源區在治理的結果上并不理想，基于潛在徑流濃度和基于負荷與產水量比值情景下的治理效果優于前者，且結果接近。這是因為面源污染受到地形、植被和土地利用類型等因素影響較大，單一情景下識別效果并不一定能準確反映面源污染對水質的影響程度。考慮到不同地區地理環境的復雜性，單一識別方法具有一定的局限性，在實際使用中，ANSP防控關鍵源區的識別應將入河負荷、潛在徑流濃度和負荷產水量比值的空間分布情況結合分析，以期得到更好的治理效果，支持農業面源污染治理科學決策。

4 結論

(1)針對農業面源污染關鍵源區識別過程中情景單一、數據量大、涉及環節多，以及模型復雜、操作繁瑣和治理效果不明確等問題，設計開發了多情景分析的農業面源污染關鍵源區識別軟件。軟件具備拓撲關系分析、產水量分析、入河負荷估算、徑流濃度分析、關鍵源區識別和關鍵源區模擬治理等功能，可模擬3種常見情景下農業面源污染分布的特點，為農業面源污染治理決策提供有效的輔助手段。

(2)以海河流域為例，進行了農業面源污染關鍵源區的識別和治理模擬應用，關鍵源區識別結果顯示，流域TN負荷呈現自西北向東南方向增長的分布特征，中部—東南部地區的潛在TN徑流濃度較高，中部—南部地區的TN負荷與產水量比值明顯高于其他地區。模擬治理的結果表明，基于入河負荷情景對氮或磷的面源污染關鍵源區治理效果較不理想，而基于潛在徑流濃度、負荷與產水量比值情景的治理效果較好，且后者略優于前者。