基于SWAT模型的布爾哈通河流域面源污染的變化研究

石金昊，朱衛(wèi)紅，田 樂，金日

（延邊大學 地理與海洋科學學院/長白山濕地生態(tài)系統(tǒng)功能與生態(tài)安全重點實驗室，吉林 延吉133002）

0 引 言

【研究意義】據(jù)2010年《第一次全國污染源普査公報》數(shù)據(jù)顯示，面源污染已經(jīng)成為我國水體的主要污染源，是我國湖泊水質惡化的主要原因。面源污染具有隨機性、復雜的遷移轉化機理、時空變化幅度大等特點[1]。面源污染是流域不同土地利用類型對降雨徑流的綜合響應，其產(chǎn)生、輸出、遷移和轉化與土地利用類型及其不同組合與空間格局特征密切相關[2]，流域景觀的組成格局和空間格局變化也會影響面源污染物的產(chǎn)生、遷移和轉化。降雨對地表的沖刷通過陸地匯流將污染物帶入水體，土地利用和土地覆蓋方式影響了地表各類污染物的分布和負荷，用地類型的改變也影響了天然的水文過程。因此，運用景觀生態(tài)學原理探究景觀格局變化與面源污染之間的關系成為面源污染控制的關鍵。

【研究進展】隨著水文模型發(fā)展和應用的不斷成熟，其中SWAT 模型被越來越多地用來評估流域尺度下的面源污染，分析其時空分布進而識別關鍵污染區(qū)域和關鍵污染期，也用于分析和評價污染控制管理措施對水環(huán)境的影響[3]。Uriarte 等[4]研究了1977—2011年流域土地利用變化對水質的影響，發(fā)現(xiàn)河流水質惡化是因為城市化和牧場的發(fā)展。Kibena 等[5]發(fā)現(xiàn)城市化進程的加快和農(nóng)業(yè)用地的增加，使該地區(qū)的面源污染加劇。而李懷恩等[6]、宋林旭等[7]運用SWAT 模型模擬了多年徑流、泥沙、氮磷等污染物濃度，表明不同的土地利用條件下產(chǎn)出的氮磷污染物濃度有較大差異，且植被覆蓋的增加可以很好地降低面源污染負荷。Wang 等[8]、孫麗娜等[9]利用SWAT 模型和CLUE-S 模型相結合的方法，分析了長時期土地利用變化下面源污染的輸出特征，并結合2 種模型模擬未來不同土地利用情景下的面源污染特征?！厩腥朦c】上述研究通常集中在流域內的土地利用類型的構成如何解釋面源污染的變化，很少有學者研究土地利用模式（景觀格局）的空間配置，例如斑塊的大小、形狀，邊緣配置或景觀內土地利用類型的斑塊的空間連通性等。但土地利用類型指標只能粗略預測水環(huán)境的變化，大多數(shù)研究沒有區(qū)分不同的景觀格局與面源污染的關系，所以難以全面反映土地利用變化對面源污染的影響?！緮M解決的關鍵問題】因此，基于SWAT模型分析布爾哈通河流域1986、1996、2006、2016年的面源污染時空分布特征，并結合實測的水文水質數(shù)據(jù)進行率定及驗證?；谒鶚嫿ǖ牟紶柟ê恿饔騍WAT 模型，將重分類后的土地利用數(shù)據(jù)代入SWAT模型，分別模擬1986、1996、2006 和2016年4 期土地利用現(xiàn)狀下的面源污染時空分布，計算1986年和2016年布爾哈通河各子流域的景觀格局指數(shù)，分析子流域面源污染和景觀格局之間的關系，以期為該流域有效減少非點源污染提供科學依據(jù)。

1 研究區(qū)概況

布爾哈通河位于42°27′N—43°23′N，129°46′E—129°38′E 之間，是國際性河流圖們江的重要支流之一（圖1），源于安圖縣哈爾巴嶺山脈東南麓，橫穿延邊朝鮮族自治州中部平原盆地，最終匯入嘎呀河。河流全長約172 km，流域面積約7 064 km2，多年平均徑流量約63 198 萬m3[10]。全流域內有耕地11.3 萬hm2，耕地率16%，其中水田2.85 萬hm2，占耕地面積的25.2%，且多集中在布爾哈通河的河谷盆地，是吉林省著名的水稻產(chǎn)區(qū)。由于開發(fā)歷史較早，墾殖率高達20%～30%，河谷內的土地大多被開發(fā)成水田和建筑用地，河谷二側坡地開墾的旱田、果園非常普遍，原有的森林資源面積急劇減少，破碎化程度不斷加劇，水土流失面積占流域總面積的39%～49%[11]，加重流域面源污染風險。

圖1 布爾哈通河流域DEMFig.1 DEM of Burhatong River Basin

圖2 布爾哈通河流域1986—2016年4 期土地利用Fig.2 The land use map for the fourth period of 1986—2016 in the Burhatong River Basin

2 數(shù)據(jù)與方法

2.1 數(shù)據(jù)來源

①土地利用數(shù)據(jù)：本研究基于1986、1996、2006年和2016年布爾哈通河流域的Landsat系列衛(wèi)星遙感影像，先對遙感影像數(shù)據(jù)做處理工作，主要包括幾何校正、波段融合、波段合成、影像拼接與裁剪等。在ArcGIS10.1 和Cognition902 軟件的利用下，運用面向對象分類方法結合人工目視解譯方法進行土地利用信息提取工作。從TM 圖像中把土地利用類型分為7 類：農(nóng)田、林地、草地、建設用地（住宅、商業(yè)和工業(yè)用地）、濕地、水域、裸地（圖2）；②為了劃分研究所需要的流域，利用NASA 提供的ASTER DEM（分辨率30 m）；③土壤數(shù)據(jù)來源于寒區(qū)旱區(qū)科學數(shù)據(jù)中心，基于世界土壤數(shù)據(jù)庫中國土壤數(shù)據(jù)集（V1.1）；④氣象數(shù)據(jù)來自中科院水利水電科學研究院的SWAT 模型中國大氣同化驅動數(shù)據(jù)集[12]（CMADS，http://www.cmads.org/），時間尺度為2008—2017年，數(shù)據(jù)包含了SWAT 模型所需的降水、氣溫、風速、太陽輻射、相對濕度；⑤水文數(shù)據(jù)來自榆樹川站點月徑流數(shù)據(jù)，時間尺度2010—2013年。月水質監(jiān)測數(shù)據(jù)（總磷），時間尺度2011—2013年，該研究僅模擬流域總磷的時空分布。

2.2 分析方法

2.2.1 流域面源污染模擬

利用上述數(shù)據(jù)構建SWAT 模型后，采用SWAT-CUP 進行校準和不確定性分析，以月徑流量、月總磷量為目標對模型進行校準和驗證，為判斷模型的模擬效果，分別選取決定系數(shù)R2和納什效率系數(shù)NSE評價模型精度[13]。校準模型之后，將1986、1996、2006年的土地利用數(shù)據(jù)分別帶入模型，模擬4 期土地利用下的流域面源污染特征。

2.2.2 景觀格局分析

將1986年和2016年的土地利用數(shù)據(jù)由矢量數(shù)據(jù)轉化為柵格數(shù)據(jù)，然后按照SWAT 劃分出的子流域分別裁剪成44 份柵格數(shù)據(jù)，最后輸入Fragstats 4.2 軟件運行，計算景觀格局指數(shù)。本研究在斑塊類型水平上選取PLAND（Percentage of Landscape）和PD（Patch Densit）來表征景觀格局的組成特征，這2 個指數(shù)是確定子流域優(yōu)勢斑塊和分析子流域不同類型斑塊破碎度的重要依據(jù)。在景觀水平上選取表征斑塊形狀復雜程度的LSI指數(shù)、表征斑塊聚集程度和連通性的CONTAG指數(shù)及COHESION指數(shù)、表征子流域景觀異質性的SHDI指數(shù)，這幾個簡單且通用的景觀度量，能比較全面地反映子流域景觀空間特征。

2.2.3 統(tǒng)計分析

冗余分析（RDA）現(xiàn)在被廣泛用于確定環(huán)境因素和生態(tài)指標之間的關系，以44 個子流域為樣本，使用CANOCO5.0軟件分析面源污染負荷和景觀格局的關系，其中箭頭的長度代表該變量被排序圖解釋的程度，箭頭越長影響程度越高，當各景觀指標箭頭與水質指標箭頭之間夾角小于90°，二者關系為正相關，大于90°為負相關，等于90°則不存在相關性[14]。

3 結果與分析

3.1 流域景觀格局變化

林地作為布爾哈通河流域的主要土地利用類型，2016年的林地面積占流域面積的73%以上，并且從1986—2016年林地的面積有著小幅度的增加。由于開墾時間較早，農(nóng)業(yè)用地在流域內的面積占比僅次于林地，在1986年達到21.45%，但到2016年縮減到了20.57%。受城市化進程和經(jīng)濟發(fā)展的需要，建設用地從1986—2016年呈不斷擴張的趨勢，其面積占比增長了0.92%，其余土地利用類型，僅有小幅變動。從1986—2016年，面積占比最大的林地有5 503.45 hm2被開墾成了農(nóng)田（表1），而大量的濕地被侵占成了農(nóng)田。因為2002 退耕還林工程在吉林省全面啟動，到了2016年，共有6 718.72 hm2的農(nóng)田轉換為林地，其次有4 577.50 hm2的農(nóng)田變?yōu)榱私ㄔO用地。1986—2016年，共有1 182.28 hm2的建設用地轉換成了農(nóng)田，農(nóng)田成為建設用地的主要轉變方向。

表1 1986—2016年土地利用轉移矩陣Table 1 1986—2016 Land use transfer matrix

本研究主要選取了林地、農(nóng)田、建設用地這3 種面積占比前3 的用地類型來分析其斑塊破碎化程度的改變，PD是表征流域斑塊破碎化程度的重要指標，由表2 可知，1986—2016年，建設用地的PD平均值基本保持平衡，林地和農(nóng)田的PD平均值有著上升的趨勢，其中農(nóng)田的PD平均值上升幅度較大，破碎化程度較高。景觀水平上表征子流域斑塊聚集程度的CONTAG指數(shù)平均值呈下降趨勢，而斑塊連通性指數(shù)COHESION平均值基本維持不變。LSI和SHDI分別表征流域斑塊形狀復雜程度和景觀異質性的指數(shù)，1986—2016年，二者的平均值都呈上升趨勢。

表2 1986—2016 布爾哈通河景觀格局變化總體趨勢Table 2 General trends of land use and landscape pattern changes in Burhatong River from 1986 to 2016

注 n 表示景觀中斑塊類型的數(shù)量。

3.2 流域面源污染時空分布特征

3.2.1 模型校準及驗證

本研究通過SWAT-CUP 2012 進行參數(shù)的全局敏感性分析，最后篩選出敏感程度較大的23 個參數(shù)，利用SWAT-CUP 程序里的SUFI-2 參數(shù)估計最優(yōu)化方法對榆樹川站點2010—2013年的月徑流量和月總磷量進行校準及驗證。徑流在校準期的R2和NSE 分別達0.86、0.81，在驗證期的R2和NES 分別達0.72、0.72?？偭自谛势诘腞2和NSE 分別達0.62、0.57，在驗證期的R2和NSE 分別達0.71、0.69。一般情況下，R2＞0.6，NSE＞0.5 時，模型能有較好的模擬效果。綜合來看，本研究構建的布爾哈通河流域SWAT模型達到了模擬精度的要求，可以很好地模擬流域水文水質狀況。

圖3 1986—2016年4 期土地利用下子流域面源污染分布特征Fig.3 Distribution characteristics of non-point source pollution in each sub-basin under land use from 1986 to 2016

3.2.2 不同土地利用的面源污染時空變化特征

時間尺度上看，流域總磷流失負荷峰值主要發(fā)生在6—8月，降雨對總磷負荷的作用方式主要通過對徑流和輸沙量的影響產(chǎn)生，一般情況下降雨量和總磷流失量呈正相關?？偭棕摵膳c徑流量大致一致，表現(xiàn)為總磷負荷隨徑流量增加而增加。已有研究[15]表明，在夏季降雨量多，總磷負荷量也隨之增多，汛期面源污染負荷量通常較高?？臻g尺度上看（圖3），總磷的流失主要發(fā)生在流域北部和東部地區(qū)，集中在9、10、32、33、34 號子流域，在2016年，這5 個子流域的總磷流失負荷占總負荷的30%。1986、1996、2006、2016年的總磷負荷分別為74.04、73.78、82.50、128.31 t。1986—1996年，總磷流失量保持在相對穩(wěn)定的狀態(tài)，但從2006年以后，總磷流失負荷開始呈急劇上升趨勢，這一差異的原因可能是進入快速城市化發(fā)展階段后，城市的快速擴張和人類生活生產(chǎn)污水的排放直接加速了總磷流失。分析結果表明，建設用地和農(nóng)田都與總磷流失負荷正相關（圖4（b））。

3.2.3 景觀格局與面源污染的關系

通過RDA 排序分析探討1986年和2016年景觀格局與非點源污染過程之間的關系。由圖4 所知，林地與TP負相關，而建設用地和農(nóng)田與TP正相關。1986年，景觀指標對總磷污染的解釋率高達86.43%（表3），總磷流失負荷貢獻度較高的景觀指標主要是農(nóng)田、PD（林地）及PD（農(nóng)田）（圖4（a））。而在2016年，景觀指標對總磷污染的解釋率為77.60%（表3），總磷流失負荷貢獻度排前二的景觀指標變成了建設用地和PD（林地）（圖4（b））。1986—2016年，農(nóng)田的PD值呈上升趨勢，并且農(nóng)田的PD值與總磷負相關，表明農(nóng)田斑塊的破碎化在一定程度上減緩了面源污染的流失。林地在布爾哈通流域的占比達到了73%，盡管林地對改善面源污染有著積極作用，但林地的PD值與總磷正相關，因為林地斑塊的破碎化會造成水土流失，從而加大面源污染流失風險。CONTAG與TP負相關，這表明破碎化嚴重且景觀斑塊較為分散的區(qū)域面源污染負荷更高。COHESIO與TP負相關，表明水質惡化更可能與流域內分散分布的土地利用有關，所以提高流域內土地利用斑塊之間的連通性，減少斑塊破碎化對養(yǎng)分流失改善有積極作用。SHDI 和LSI都與總磷正相關，表明流域內景觀異質性的增加和形狀的復雜化，能夠導致污染負荷輸出風險水平的增加。

表3 冗余分析各排序軸方差解釋率Table 3 Total variance explained by the ordination axis

圖4 1986 和2016年污染物流失量與景觀指標的冗余分析Fig.4 Redundant analysis of pollutant loss and landscape pattern index in 1986 and 2016

4 討論

通過構建布爾哈通河流域SWAT 模型發(fā)現(xiàn)，流域降雨和徑流過程對面源污染流失有著重要影響，總磷流失量和徑流量正相關，流域豐水期多出現(xiàn)在6—9月，總磷產(chǎn)生量峰值也主要在6—9月，所以該時期是非點源污染的關鍵期。農(nóng)村和城市是非點源污染的主要來源，但在不同的歷史發(fā)展階段，二者對面源污染的影響也不同。基于4 期的流域土地利用數(shù)據(jù)，通過構建SWAT 模型分析因土地利用變化而導致的面源污染時空分布特征，發(fā)現(xiàn)受地形、農(nóng)業(yè)生產(chǎn)及城市發(fā)展影響的9、10、32、33 和34 號子流域總磷污染負荷較高，屬于流域污染的關鍵控制區(qū)。農(nóng)田和建設用地是面源污染的主要來源，而森林對面源污染改善有著積極影響，這與前人研究的結果基本一致[16]。具有污染輸出效應的建設用地和農(nóng)田一直被視為面源污染的“源”景觀，其中建設用地對水質的影響主要取決于流域內城市化水平的高低。大量的生活垃圾、工業(yè)污水的排放及農(nóng)業(yè)生產(chǎn)使用的化肥、農(nóng)藥等會經(jīng)過城市不透水表面進入河流。布爾哈通河中游到下游流域二岸分布城市、村屯、耕地等，來自城區(qū)的生活污水、農(nóng)業(yè)污水及工業(yè)廢水對流域面源污染有著負面影響。該流域傳統(tǒng)的農(nóng)業(yè)活動（即耕作或過度施肥）增加了土壤侵蝕的風險，并且過量的養(yǎng)分應用可能導致非點源污染經(jīng)地表徑流后流入到河流。此外，通過對流域的景觀格局的進一步分析，發(fā)現(xiàn)了農(nóng)田對面源污染貢獻度降低的可能原因。有研究表明，農(nóng)田斑塊的破碎化程度越低，斑塊廊道聯(lián)通性越好，農(nóng)田污染物容易形成面源污染[17]。盡管33 號子流域的建設用地面積最大，但33 號子流域的污染輸出負荷比32 號和34 號小。已有研究表明，流域“源”和“匯”景觀的空間位置不同，其對水環(huán)境的影響也不同[18]。33號子流域的建設用地和農(nóng)田基本位于上游區(qū)域，下游區(qū)域則是大面積的森林，這對該流域的面源污染起到了更好的凈化作用。

森林一直被認為與農(nóng)業(yè)面源污染有著直接的聯(lián)系，其植被覆蓋可以減少土壤侵蝕，具有削減暴雨徑流、減少水土流失、吸附污染物的作用，可以有效減少地表徑流攜帶營養(yǎng)鹽進入河流，對水質惡化在一定程度上有削減作用[19]，但往往局限于景觀組成格局的角度分析得出該結論。以具有污染物消減效應的森林為優(yōu)勢斑塊的流域景觀格局中，其對景觀破碎程度最為敏感[20]。林地景觀越破碎，斑塊越分散，越不利于抑制土壤侵蝕、徑流等過程，從而造成污染元素的流失。本研究也通過冗余分析發(fā)現(xiàn)林地景觀斑塊破碎化的增加，會增大流域面源污染產(chǎn)生的風險，而林地作為布爾哈通河流域的主要用地類型，其景觀空間格局變化引起的生態(tài)變化，更值得注意和警惕。

區(qū)域水質是多尺度環(huán)境因子的綜合反映，不單單受土地利用的影響，景觀格局的變化也被認為是影響區(qū)域水質的主要原因[21]，景觀格局與流域非點源污染研究是探究景觀格局與地表環(huán)境變化過程的典型模式。景觀形狀指數(shù)（LSI）、香濃多樣性指數(shù)（SHDI）與總磷污染呈正相關，主要因為這些指標反映了子流域中斑塊的多樣性，這些多樣性的增大又間接地體現(xiàn)了景觀的異質性與人類活動干擾的增強，因此也增加了面源污染的風險，其他相關研究結果也表明這些指標與水質的惡化有密切聯(lián)系[22]。斑塊聚集度指數(shù)（CONTAG）反映的是景觀的分離與散布程度，其高值反映了高度聚集。例如，當不同的土地利用最大限度地聚集時，CONTAG接近100，而當土地利用最大限度地分散時接近與0。斑塊連通性指數(shù)（COHESIO）表示斑塊內細胞的空間連通性和流域內土地利用的物理連通性，高蔓延度值和高連通性說明景觀中的某種優(yōu)勢斑塊類型高達集聚并且形成了良好的連接性，景觀的破碎化程度較低，所以面源污染風險也會隨之減小。有相關研究表明坡度主要通過影響產(chǎn)流量和產(chǎn)沙量進而影響磷素的流失總量[23]，而在布爾哈通河流域坡度較高區(qū)域，在雨水充沛且多暴雨出現(xiàn)的條件下，會加重水土流失，致使大量未被植物吸收利用的氮磷營養(yǎng)元素隨水土流失進入河流，是造成流域面源污染的重要原因。

5 結論

1）SWAT 模型能較好地模擬布爾哈通河流域總磷污染時空分布，反映研究區(qū)的總磷污染流失負荷。

2）布爾哈通河流域面源污染防治的關鍵控制區(qū)域在坡度較高的北部區(qū)域和農(nóng)田和城鎮(zhèn)集中的東部區(qū)域。

3）通過冗余分析發(fā)現(xiàn)林地與總磷負荷負相關，農(nóng)田和建設用地與總磷負荷正相關；林地斑塊破碎化越高，流域面源污染流失越嚴重，農(nóng)田斑塊的破碎化在一定程度上減緩了面源污染的擴散。流域斑塊之間的高度流通性和聚集性，對養(yǎng)分流失改善有積極作用。流域內景觀異質性的增加和形狀的復雜化，能夠導致污染負荷輸出風險的增加。對該區(qū)域進行土地利用規(guī)劃與景觀格局優(yōu)化時，應從景觀組成格局和空間格局二方面綜合考慮，才能更好地改善流域面源污染。