密云水庫流域土壤侵蝕模數時空動態監測

郭 晴，任永康，陳 濤*

（1. 中國地質大學（武漢） 地球物理與空間信息學院，湖北 武漢 430074）

土壤侵蝕對土壤分級、農業生產、水文系統和環境的不利影響，長期以來被認為是影響人類可持續性發展的嚴重問題[1]。截至2019 年，我國土壤侵蝕總面積已達到271.08 萬km2，成為世界上土壤侵蝕最嚴峻的國家之一。其中，水力侵蝕面積為113.47 萬km2，占總面積的41.86%；風力侵蝕面積為157.61 萬km2，占總面積的58.14%[2]。土壤侵蝕會引起土壤肥力衰退，進而威脅農業生產安全。徑流泥沙搬移的污染源將影響侵蝕鄰近區生態和經濟的發展，且侵蝕泥沙的搬運使得土壤碳、氮、磷等元素的含量與組分產生變化，從而影響全球生源要素的循環，甚至成為全球氣候變化的驅動要素之一[3]。

土壤侵蝕研究通常采用實地測算和構建土壤侵蝕模型來開展。由于實地觀測難以直接監測空間大范圍土壤侵蝕情況，且易受主觀因素和監測技術的限制，因此不能滿足土壤侵蝕治理工作對時效性、快速性與準確性的需求[4]；而土壤侵蝕模型不僅可以解決實地觀測中的矛盾，還保障了土壤侵蝕預測和評價的準確性。根據研究方法，土壤侵蝕模型可分為物理模型和經驗統計模型[5]，物理模型中最具代表性的是水蝕預報模型[6]，雖對土壤侵蝕的描述十分全面，但所需數據量大且難以在大區域中實現[7]；經驗統計模型中認可度最高、應用范圍最廣的是Kenneth G R[8]等提出的土壤流失方程（USLE）和修正后的通用土壤流失方程（RUSLE）[6]。隨著土壤侵蝕研究的不斷深入，許多學者結合研究區的真實狀況對模型中各因子進行了修正，并與GIS 和遙感技術進行集成，這樣不僅具有了RUSLE 模型結構簡單、所需數據量少的特點[9]，而且具有了GIS 和遙感技術快速獲取、有效管理并實時分析時空大數據的優勢，為土壤侵蝕提供了科學、快速、準確的監測和評估。例如，劉寶元[10]等根據我國國情，結合GIS 和遙感技術，輔以遙感圖像，修正了RUSLE 中的各個因子；楊佳佳[11]等基于遙感和GIS 技術，利用修正后的RUSLE模型對黑龍江典型黑土區土壤侵蝕進行了監測；王丹媛[12]等利用遙感數據修改了桂南沿海諸河流域USLE 模型中的因子，并評價了土壤侵蝕強度。近年來，眾多學者將RUSLE模型與GIS相結合以研究大尺度土壤侵蝕和水土保持格局與規律[13-14]，如陳峰[15]等利用RUSLE 模型研究了滇南山區土壤侵蝕的時空演變情況；王志杰[16]等運用RUSLE模型探究10 a間貴陽市土壤侵蝕空間演變規律。

作為北京市最大的飲用水源供應地，密云水庫擔負著北京生產生活用水的重要任務。然而，密云水庫流域生態環境較脆弱，土壤侵蝕不利于流域內的水循環，進而將威脅北京市民的飲用水安全[17]。因此，研究密云水庫的土壤侵蝕情況，為改進水庫周邊水土保持工作提供科學性的指導，對于保護其生態環境具有重要意義。目前，針對密云水庫流域土壤侵蝕的研究主要為單一年份，缺乏時空動態變化的分析，在一定程度上不利于決策和管理，阻礙了流域的可持續發展。本文基于遙感和GIS技術，利用RUSLE模型對密云水庫流域2001—2020年土壤侵蝕情況進行了動態監測，分析了時空動態變化特征及其影響因素，為密云水庫流域的水土保持提供了科學合理的依據。

1 研究區概況

密云水庫流域位于北京市城北 13 km2（41°14′~41°05′N、116°07′~117°30′E），面積為180 km2， 包括北京市的密云、懷柔、延慶以及河北省的灤平、豐寧、赤城等10個區縣（圖1）；坐落在潮河以及白河中游偏下，流域面積為15 788 km2。水庫流域年均氣溫為9～10 ℃，年均降水量為488.9 mm，主要集中在汛期。流域內土壤主要包括褐土、棕壤、草甸土和粟鈣土，植被以闊葉混交雜木林、油松、刺槐和落葉松為主。

圖1 密云水庫概要圖

2 數據來源與研究方法

2.1 數據來源與預處理

本文采用的遙感數據包括2001 年、2005 年、2009 年 、 2011 年 、 2015 年 和 2020 年 6 期 共 24 景Landsat8/TM 遙感影像。由于部分年份遙感影像存在云量較多、噪聲明顯和數據條帶丟失等問題，因此在等間隔的基礎上從研究年份中選擇數據精度較高的6 期數據作為研究對象；由于研究區遙感影像數據具有多時相、多光譜、多傳感器等特點，因此利用ENVI 5.3對原始影像進行輻射定標、大氣校正等預處理，形成標準化的多源遙感數據集。本文采用的地形數據為ASTER GDEM 30M產品；降水量數據包括2001—2020年流域內豐寧、赤城、沽源、崇禮、承德、懷來、灤平、懷柔、密云、興隆、延慶和張家口地區的年均降水量；土地利用數據包括2001年、2005年、2009年、2011 年、2015 年和 2020 年 6 期 1 km 分辨率的中國土地利用數據（包含耕地、林地、草地、水域、建設用地和未利用地6個一級地類）；土壤數據為土壤科學數據庫中的中國第二次普查數據，包括土壤質地砂礫、粉粒和有機碳含量等數據；各區縣生產總值來源于當地政府官方網站數據；人口數據來源于中國統計信息網（表1）。

表1 數據列表

2.2 研究方法

基于GIS和遙感技術，本文利用RUSLE模型計算得到土壤侵蝕模數，并利用轉換矩陣分析了2001—2020年土壤侵蝕強度的空間動態變化情況。RUSLE模型的計算公式為：

式中， A 為土壤侵蝕模數，單位為t/(hm2·a)；R 為降雨侵蝕力因子， 單位為(MJ·mm)/(hm2·h·a)； K 為土壤可蝕性因子 ，單位為(t·hm2·h)/(hm2·MJ·mm) ； LS 為坡度坡長因子，無量綱；C 為植被覆蓋與管理因子，無量綱；P 為水土保持措施因子，無量綱。

1）降雨侵蝕力因子表示降雨引起土壤侵蝕的潛在能力。本文基于監測區域代表站點連續30 min的最大降雨強度和每次降雨的降雨量，利用基于日降雨量的擬合模型[18-19]計算降雨侵蝕因子。為了得到流域內降雨侵蝕力因子的連續數據，在各測站降雨侵蝕力值和空間分布的基礎上，利用GIS 的樣條插值實現。其計算公式為：

式中，Pf為汛期雨量；I30B為侵蝕性降雨最大30 min雨強的年代表值，單位為mm；i 為侵蝕性降雨次數；Pt為汛期月份的月降雨總量。

2）土壤可蝕性因子是指土壤對侵蝕的敏感性，與降雨、土壤自身的理化特性等相關。采用Wisch?meier W H[20]等提出的公式，根據土壤質地、土壤中的有機質含量等主要土壤性質數值指標估算土壤可蝕性因子。其計算公式為：

式中， M 為優勢粒徑組成的乘積；a 為有機質含量百分比；b 為土壤結構等級；c 為土壤滲透等級。

3）坡度坡長因子。坡長因子是指任意坡長的單位面積土壤流失量與標準坡長單位面積土壤流失量之比；坡度因子是指任意坡度下的單位面積土壤流失量與標準坡度下單位面積土壤流失量之比。由于坡度從20%增加到40%和60%時，坡長指數無變化[6]，因此當坡度<21%時，采用式（5），當坡度＞21%時，采用式（6）[21]。

4）植被覆蓋與管理因子是指在土壤、坡度和降雨條件相同的情況下，具有特定植被覆蓋的土地的土壤流失量與裸地的土壤流失量的比值。本文采用蔡崇法[22]等提出的利用植被覆蓋度與定量估算管理因子的回歸方程，即

式中， f 為植被覆蓋度；C 為定量估算管理因子

植被覆蓋度計算采用李苗苗[23]等提出的改進的像元二分模型，則有：

式中， NDVI 為歸一化植被指數值，無量綱；NDVImax、NDVImin為遙感圖像中給定置信度置信區間內的最大值與最小值，本文采用的置信區間為[5%，95%]。

5）水土保持措施因子是指采取一定的水保措施后的土壤流失量與標準狀況下土壤流失量的比值。水土保持措施因子目前主要通過布設天然小區試驗得到。由于難以獲得水土保持措施因子，在考慮財力物力人力等綜合因素后，結合前人的研究結果[24-25]，本文將水土保持措施因子設定為1。

2.3 轉移矩陣

土壤侵蝕轉移矩陣反映了研究區在某一時段期初和期末不同侵蝕等級面積之間的轉化信息，包括期初各類面積的轉出量和期末各類面積的轉入量[26]。本文計算得到2001 年和2020 年不同侵蝕等級土壤侵蝕模數的轉移矩陣，并探究了侵蝕等級的變化特征。

2.4 相關性分析

相關性分析能反映因變量與自變量之間相互關系的密切程度。本文將土壤侵蝕模數分別與年均降雨量與平均植被覆蓋度進行統計回歸分析[27]，以探究年均降雨量、平均植被覆蓋度與土壤侵蝕模數的相關性。

3 研究結果分析

3.1 土壤侵蝕模數的時間變化特征

密云水庫流域年均土壤侵蝕模數呈“快速減小、小幅增長、又迅速回落”的趨勢，20 a間研究區土壤侵蝕狀況整體得到明顯控制，2009—2015年土壤侵蝕呈現惡化趨勢。根據國家水利部發布的SL190-2007《土壤侵蝕分類分級標準》[28]，本文將研究區土壤侵蝕分為微度、輕度、中度、強烈、極強烈和劇烈侵蝕6個等級，并統計分析不同時期密云水庫流域土壤侵蝕強度的面積比例。20 a間研究區土壤侵蝕以輕度和微度侵蝕為主，二者占總面積的70%以上，且呈波浪式變化；其次為中度侵蝕；除強烈侵蝕和極強烈侵蝕面積以外，其他侵蝕強度面積呈總體下降趨勢。2001—2020年各侵蝕強度等級的年變化率（降幅）均在5%~10%，但2001—2020年強烈和極強烈侵蝕增幅明顯（圖2）。

圖2 2001—2020年不同等級平均土壤侵蝕模數折線圖

由面積轉移矩陣可知（表2）， 2001—2020 年研究區土壤侵蝕以微度侵蝕少量轉出、其他侵蝕等級不同程度轉入為主要特征；侵蝕強度等級的轉出方向以向低等級侵蝕轉化為主要特征，微度侵蝕轉入了14.69 km2，占總微度侵蝕轉出面積的67.35%，輕度侵蝕轉入微度侵蝕的面積為2.91 km2，占輕度侵蝕轉出面積的56.7%，中度侵蝕轉入微度和輕度侵蝕的面積為0.43 km2，占總中度侵蝕轉出面積的67.2%，劇烈侵蝕幾乎全部轉向低等級侵蝕；但極強烈侵蝕大部分面積轉入為更高等級，強烈侵蝕轉向低等級侵蝕的面積占比為44.7%，因此流域內20 a 間土壤侵蝕呈“整體好轉、局部惡化”的趨勢。

表2 2001—2020年研究區土壤侵蝕強度面積轉移矩陣/km2

3.2 土壤侵蝕模數的空間變化特征

3.2.1 研究區不同行政區內土壤侵蝕變化分析

本文疊置分析了密云水庫流域研究期土壤侵蝕強度分布圖與研究區行政區劃圖，獲得流域內研究期間各縣區土壤侵蝕強度的面積比例（表3），可以看出，流域內土壤侵蝕最嚴重的區域為密云區和灤平縣，均值都在1 500 t/(km2·a)以上；密云、崇禮和灤平縣的年均侵蝕模數降幅最大，分別為53.9%、63.5%和42.5%；崇禮縣水土保持情況最好，平均土壤侵蝕模數接近 700 t/(km2·a)，且降幅較大。

表3 密云水庫流域2001—2020年各縣區侵蝕強度面積/km2

3.2.2 區域內不同海拔下土壤侵蝕變化分析

本文疊置分析了土壤侵蝕分布圖與數字高程模型分級圖，得到密云水庫流域2001—2020年不同海拔下土壤侵蝕變化情況（表4、5），可以看出，隨著海拔的升高，研究區土壤侵蝕模數整體呈“高—低—高”的特點，海拔500 m 以下以及海拔1 500 m 以上的地區侵蝕相對嚴重。在低海拔地帶，由于地形平坦有利于開發，導致人口分布密集、城市化程度較高，進而對地表植被造成嚴重破壞；而在1 500~3 000 m高海拔地區，土地類型以高山裸土區和荒草地為主，地勢陡峻，高強度的降雨易產生滑坡、泥石流等自然災害。

表4 密云水庫流域2001—2009年不同海拔下土壤侵蝕變化分析

表5 密云水庫流域2011—2020年不同海拔下土壤侵蝕變化分析

3.2.3 區域內不同土地利用類型下土壤侵蝕變化分析

本文統計分析了不同土地利用類型的侵蝕狀況，結果如表6所示。各種土地利用類型中，耕地的平均土壤侵蝕強度普遍較高，而城鄉居民用地和水域的平均土壤侵蝕強度較低。2001年、2009年耕地和林地的侵蝕強度接近于未利用地，處于較高水平；2015年耕地、林地和草地的侵蝕強度均大于1 500 t/（km2·a），林地的侵蝕強度最大，未利用地的侵蝕強度未達到1 500 t/（km2·a）；2020年耕地、林地和草地的侵蝕強度均有所回落，與逐步治理耕地以及大力推行退耕還林政策密不可分。

表6 密云水庫流域2001—2020年不同土地利用類型下土壤侵蝕變化/（t/（km2·a））

3.3 影響因素分析

3.3.1 自然因素

本文對密云水庫流域研究期間的平均土壤侵蝕模數與平均植被覆蓋變化、平均降雨侵蝕模數進行回歸分析，結果表明，平均土壤侵蝕模數與平均植被覆蓋變化存在明顯的負相關性，相關系數為-0.77（圖3）；而與平均降雨侵蝕模數呈正相關性，相關系數為0.71（圖4）。由于平均土壤侵蝕模數與植被覆蓋度的決定系數大于其與降雨侵蝕模數的決定系數，因此平均土壤侵蝕模數與植被覆蓋度的相關程度更高。

圖3 植被覆蓋度與土壤侵蝕模數關系散點圖

圖4 平均降雨侵蝕模數和土壤侵蝕模數關系散點圖

3.3.2 人為因素

分析研究區各縣市人口數據（圖5）和土壤侵蝕強度的面積比例統計表可知，研究期內密云區、興隆縣、灤平縣和懷柔區的平均侵蝕面積均在1 000 km2以上，且人口逐年遞增，其中密云區和懷柔區的人口基數大、增速也較快，人口的迅速增長促使農民大量開荒種糧，導致土地墾荒率越來越高，加劇了水土流失；后期（2015—2020年）崇禮縣、沽源縣和延慶區等地區人口出現零增長甚至負增長的情況，人口對環境的壓力下降了，水土流失也有所緩解。

圖5 研究區各縣區常駐人口變化圖

分析各區縣GDP 變化情況（圖6）和土壤侵蝕強度的面積比例統計表可知，土壤侵蝕與社會經濟發展具有相關性，研究區初期和中期（2001—2011年）各行政區的GDP增速較大，如懷柔區、密云區和灤平縣，其土壤侵蝕面積相對較大，原因在于其地理位置離北京主城區較近，經濟社會活動較頻繁，經濟的發展伴隨著開發區建設等土建工程的急劇增加，使得大量草地、林地遭到了破環，加劇了水土流失；后期各行政區GDP增速逐漸放緩，政府踐行綠色發展理念，加大對水土保持工作的投入，使得土壤侵蝕面積較前期有所回落。

圖6 研究區各縣區GDP變化圖

4 結 語

本文利用RUSLE 模型對密云水庫流域2001—2020 年的土壤侵蝕情況進行了定量研究。結果表明，2001—2020年密云水庫流域土壤侵蝕模數大致經歷了先下降后回升又大幅下降的變化過程，最小值出現在2009 年，2020 年比2001 年的土壤侵蝕面積有所減少，說明整體情況在好轉；流域內侵蝕等級主要由輕度轉為微度，但強烈侵蝕和劇烈侵蝕主要集中在灤平縣和密云區；隨著海拔的上升，流域內侵蝕模數先升高后降低又升高；土壤侵蝕模數與土地利用類型相關性程度較高；植被覆蓋度升高，土壤侵蝕模數下降，降雨侵蝕加劇，土壤侵蝕模數降低，但前者的影響程度較大；人口與社會經濟的發展也會對土壤侵蝕模數造成影響。