基于GIS的太原市土壤侵蝕定量研究

鹿晨昱，張 琳，2，薛 冰，張子龍，逯承鵬，李勇進

（1.西北師范大學 地理與環境科學學院，甘肅 蘭州730070；2.中國科學院 沈陽應用生態研究所，遼寧 沈陽110016；3.蘭州大學 中國西部循環經濟研究中心，甘肅 蘭州730000）

土壤侵蝕是土地退化的根本原因，也是導致生態環境惡化的重要因素。對土壤侵蝕進行定量研究，能夠為區域土地的持續利用提供科學決策依據。太原市地處山西省腹部，地形以山地、丘陵為主，中南部為汾河沖積扇平原，地勢平坦，土壤肥沃，灌溉方便，農業發達。

太原市是山西省農業發達地區，隨著區域人口的增加和社會經濟的發展，人類加大了對土地資源開發的力度，引起土地利用景觀格局發生變化，對土地資源過度或無科學指導的利用，導致生態環境惡化，產生了諸如土地退化、水土流失等嚴重問題。本研究選取太原市區作為研究區，在行政區劃上包括杏花嶺、小店、迎澤、尖草坪、萬柏林、晉源6個區，基于GIS與RS技術，進行土壤侵蝕定量評價及分析，研究太原市土壤侵蝕分布格局，以期能為區域土地的持續利用提供科學決策依據。

1 材料與方法

1.1 研究方法

土壤侵蝕計算模型的選擇正確與否直接關系到計算結果的科學性，到目前為止最為廣泛應用的計算模型是通用土壤流失方程（USLE），但其是在美國國情的基礎上研究得到的。后來科學家在全面考慮各影響因子的基礎上，提出了修正的土壤流失方程（RUSLE），該模型是建立在土壤侵蝕理論及大量實地觀察數據統計分析的基礎上，形式簡單、使用方便、易于推廣接收。因此，本文采用RUSLE模型進行研究。其表達式為:

式中:A——年土壤侵蝕量；R——降雨侵蝕力因子；K——土壤可蝕性因子；L——坡長因子；S——坡度因子；C——植被覆蓋因子；P——水土保持措施因子。

1.2 技術路線

基于RS和GIS空間分析平臺，采用RUSLE模型，在對區域土壤侵蝕影響因子進行分析的基礎上，對太原市2010年土壤侵蝕進行定量研究，并根據土壤侵蝕強度進行分級。

（1）降雨侵蝕力因子R。主要依據現有的降雨資料，計算得到太原市所轄各區的降雨侵蝕力。（2）植被因子C。基于Landsat－TM影像計算的植被覆蓋度（NDVI），參考 DeJong［1］在歐洲土壤侵蝕研究中提出的利用NDVI計算C值的方法，計算太原市的C值。（3）地形因子LS。選擇基于SRTM－DEM的坡度坡長因子LS的計算方法［2］，在ArcGIS平臺下計算，獲取太原市的地形因子值。（4）土壤可蝕性因子K。主要參考中國科學院水土保持研究所關于太原土壤特性的研究成果，在此基礎上進行提取。（5）水土保持因子P。由于土地利用/土地覆被能夠基本反映植被覆蓋率的高低與水土保持措施的差異，因此依據土地利用類型賦值的方法來加以確定。

2 土壤侵蝕影響因子分析

2.1 降雨因子R

2.1.1 降雨侵蝕力 降雨是引起土壤侵蝕最直接的驅動因素。降雨時雨滴所攜帶的動能對土壤顆粒進行沖擊，產生擊濺侵蝕，而后還產生徑流，沖刷地表，導致土壤侵蝕的發生。降雨侵蝕力就是反映降雨對土壤剝離及搬運的動力指標。降雨侵蝕力與降雨量、降雨動能、降雨強度、降雨歷時、雨型等因素有關，是降雨物理性質的函數［3］。

2.1.2 降雨侵蝕力計算 降雨侵蝕力反映了降雨對土壤侵蝕的影響。長期以來，國內外學者關于降雨侵蝕力做了大量的研究工作，總體來說有經典算法和簡易算法兩種計算方法。采用的數據源為6個雨量站點的月降雨量數據，由于雨量站有限，因此采用章文波、付金生［4］提出的利用逐月雨量估算降雨侵蝕力的簡易算法。計算公式為:

式中:Pij——第i年j月的降雨量（mm）；N——年數（a）；R——降雨侵蝕力﹝MJ·mm/（hm2·h·a）﹞；α，β——模型參數。α＝0.183 3，β＝1.995 7。

降雨侵蝕力是反映降雨對土壤流失的影響，降雨侵蝕力越大的地方，土壤流失相對越嚴重。根據公式（3）計算太原市所轄各區的年降雨侵蝕力，可以得出，晉源區的降雨侵蝕力最大﹝289.097MJ·mm/（hm2·h·a）﹞，其次是萬柏林區﹝267.742MJ·mm/（hm2·h·a）﹞，杏花嶺區﹝233.665MJ·mm/（hm2·h·a）﹞、迎澤區﹝225.253MJ·mm/（hm2·h·a）﹞和尖草坪區﹝223.031MJ·mm/（hm2·h·a）﹞降雨侵蝕力相仿，小店區﹝178.634 MJ·mm/（hm2·h·a）﹞的降雨侵蝕力最小。

2.2 地形因子LS及其計算

地形因子是影響土壤侵蝕的基本自然驅動要素，在RUSLE及其它許多土壤侵蝕模型中，坡度和坡長是度量地形因子對土壤侵蝕影響的主要指標，通常用LS來表示。在RUSLE中，坡度因子被定義為某一坡度的單位面積土壤流失量與坡度為5.130，其它條件都一致的坡面產生的單位面積土壤流失量之比率；坡長因子被定義為某一坡面的單位面積土壤流失量與坡長為22.13m，其它條件都一致的坡長產生的土壤流失量之比值；LS因子則是指某一坡面的土壤流失量與坡度為5.130，坡長為22.13m，其它條件都一致的坡面所產生的土壤流失量之比值［5－6］。

傳統坡面尺度上的土壤侵蝕評價，其坡度和坡長因子可以靠野外實測來計算，但是在區域尺度上，這種方法是不現實的，其地形指標要通過DEM來提取。因此，選擇能夠較為真實地表達研究區實際地形的DEM，是提取地形因子的基礎［7］。

關于利用DEM來提取LS因子，國內外學者進行了大量的研究，其中 Hickey［8］和 VanRemortel［9］提出了在ArcGIS下計算LS值的算法，而國內學者的研究則主要是對LS算法進行修正，以使其適用于我國不同的區域。

本文在ArcGIS軟件柵格模塊的支持下，利用柵格計算器，實現LS值的自動生成。算法為:

式中:S——坡度因子；θ——由ArcGIS提取的坡度值；L——坡長因子；λ——由 DEM 提取的坡長；α——坡長指數；β——參數。計算得到太原市坡度坡長因子。分析可知，坡度越大，徑流能量越大，對坡面的沖刷能力越強，土壤侵蝕就越嚴重；坡長與土壤侵蝕的關系雖較為復雜，但也表現出坡長增加，徑流量增大，侵蝕作用加強的趨勢。太原處于盆地地帶，西、北、東三面環山，中、南部為河谷平原，形成從中部延伸到南部的汾河流域，汾河作為具有較強侵蝕能力的徑流，對周圍土壤進行長期侵蝕。在地形背景的影響下，中部和南部的坡長、坡度是土壤侵蝕的最主要影響因素。由計算結果可知，市區中部和南部的坡度坡長因子較高，表明該區域的土壤侵蝕相對較為嚴重。

2.3 土壤可蝕性因子K

土壤是指陸地表層具有一定的肥力并且能夠生長植物的疏松表層，它是土壤侵蝕發生的主體，是被侵蝕的對象。土壤受侵蝕的影響主要是由其各種理化特性決定的［10］。土壤可蝕性K定義為標準小區單位降雨侵蝕力引起的年土壤流失率。本研究利用Wiiliams［11］提出的土壤可蝕性因子K與土壤性質之間的關系式確定土壤因子K值:

式中:Sa——砂粒含量（%）；Si——粉砂含量（%）；Cl——黏粒含量（%）；C——有機碳含量（%）；Sn＝（1－Sa）/100。

參照《太原市土壤志》，根據研究區各種土壤的機械組成、粒級含量和有機質含量，利用公式計算出不同土類的可蝕性因子K值（表1）。在此基礎上，根據中國科學院水土保持研究所關于太原市土壤特性的研究成果，在土壤類型圖上進行賦值，得到土壤可蝕性因子。

表1 太原市土壤性質及K值

由表1可知，太原市土壤可蝕性K值分布于0.35～0.42（t·hm2·h）/（hm2·MJ·mm）。K值越大，土壤抗蝕性越差，反之亦然。太原市區K值最大約為0.42（t·hm2·h）/（hm2·MJ·mm），主要分布在太原市中南部地區。該地帶處于河谷平原，土壤類型以廣泛分布在河漫灘地的潮土為主，而潮土結構疏松，遇水易分解，抗蝕性差，因此K值較高。北部及西部的土壤可蝕性K值次之，約為0.38（t·hm2·h）/（hm2·MJ·mm）。其余區域的可蝕性K值相對最低，約為0.35（t·hm2·h）/（hm2·MJ·mm）。

2.4 植被覆蓋與管理因子C及其計算

植被是土壤侵蝕的主要抑制因素。應用廣泛的通用土壤流失方程中的植被覆蓋與管理因子C（在1965年版的USLE方程中被稱為種植管理因子）是這種抑制作用的定量化表達指標之一。C因子指一定植被覆蓋條件下，耕作農地上的土壤流失量與同等條件下實施清耕的連續休閑地的土壤流失量之比，為一無量綱參數，其值大小變化介于0～1［12］。C值越大，意味著與其對應的地域的土壤侵蝕越嚴重。

目前已有的40多種植被指數中，歸一化植被指數NDVI是水土流失預測和研究中應用最廣泛的植被指數［6］。一般來講，通用土壤侵蝕模型計算C值的方法是通過建立植被指數與C值的關系來實現的，其計算公式為:

式中:N——NDVI數值；α，β——決定NDVI與C關系曲線的參數，較合理的設置是α＝2，β＝1。

計算得到太原市植被因子。總體上看，西南部的植被覆蓋與管理因子C值相對較高，表明該區域的土壤侵蝕程度相對較為嚴重；東北部的植被覆蓋與管理因子C值相對較低，說明該區域的土壤侵蝕程度相對較輕。

2.5 水土保持因子P

水土保持措施因子P是指特定水保措施下的土壤流失量與相應未實施水保措施的順坡耕作地塊的土壤流失量之比。P值范圍在0～1之間，0代表根本不發生侵蝕的地區，而1值代表未采取任何控制措施的地區。

根據太原市水土保持試驗站資料以及縣志資料，并結合當地土地利用及農事活動情況，通過賦值法確定太原市不同土地利用類型的P值。其中按土地利用類型劃分，林地是指生長喬木、竹類、灌木、沿海紅樹林的土地，此類區域不采取控制措施，因此將P值賦為1；建設用地是指建有建筑物、構筑物的土地，此類區域的土地利用率最大，因此將P值賦為0；研究區的水體是由貫穿太原市區的汾河及部分人工水庫組成的，該類區域P值也賦為0；耕地中，由于坡耕地是山區泥沙來源的策源地，因此將坡耕地的P值賦為0.5，梯田及其它耕地的P值賦為0.4；將未利用地（包括荒草地、鹽堿地、裸巖石礫地）P值賦為0.7。通過對太原市的Landsat TM遙感影像進行非監督分類與監督分類，得到太原市土地利用類型圖，經驗證分類的精度可達90%以上。之后，將相應的P值賦給不同的土地利用類型，得到水土保持措施因子。由結果可知，太原市由南向北P值逐漸升高，南部的P值相對較低，意味著土壤侵蝕較輕或不發生土壤侵蝕。其中東部和西部的邊緣地帶多為以梯田為代表的耕地，因此這些地帶土壤侵蝕程度較中部和北部要輕。總體上看，采取水土保持控制措施的面積約占總面積的2/3。

3 太原市土壤侵蝕評價

3.1 基于RUSLE的土壤侵蝕評價

在ArcGIS軟件中利用地圖代數的計算方法運行RUSLE模型進行太原市土壤侵蝕評價，并參考水利部頒布的《土壤侵蝕分類分級標準》，結合研究區土壤侵蝕量的實際分布特征，確定具有區域特色的土壤侵蝕強度分級標準。最終將太原市土壤侵蝕等級分為微度侵蝕、輕度侵蝕、中度侵蝕、強度侵蝕、極強度侵蝕5個等級（表2），并得到2010年太原市土壤侵蝕模數圖（圖1）。

圖1 太原市土壤侵蝕模數圖

表2 太原市土壤侵蝕程度分級

3.2 土壤侵蝕空間分布規律

由表2及圖1可知，太原市不同等級的土壤侵蝕具有一定的空間分布規律:（1）極強度侵蝕分布在市區的邊緣地帶，約占侵蝕面積的2%，該地帶的土地類型為未利用地，呈零星狀。其主要影響因素是地形因素，太原市發源于河谷平原，由于地形導致邊緣的坡面較大，在一定程度上土壤流失量與坡面呈正相關，故侵蝕程度最大；其次是水土保持因素，該區域土地質量較差，多為荒廢地，相關部門對此地帶的管理程度較為欠缺，因此加重了土壤侵蝕程度。（2）強度侵蝕主要集中分布于東北部山地和西部零散地區，約占侵蝕面積的21%，此地帶的土地類型為未利用地。首先，太原市的地形是北高南低呈簸箕形，因此地形因素為主要影響因素；其次，該區域土質的植被覆蓋極少，屬于荒廢地，荒漠化嚴重，儲水能力極低，故侵蝕強度較大。（3）中度侵蝕主要分布在中部的建設用地和水體周邊（以汾河流域、晉陽湖附近區域為主），以及呈零星分布狀的坡耕地，約占侵蝕面積的22%。首先，晉陽湖位于晉源區，該地區的降雨侵蝕力相對最高，居6個區之首，因此降雨因素對土壤侵蝕的影響較大；其次，太原市地形是北高南低，這種地形特點也是一個重要影響因素。（4）輕度侵蝕主要集中分布在東部和南部的以梯田為代表的耕地地帶（坡耕地除外），約占侵蝕面積的25%。此地區土質為潮土，由于其結構疏松，遇水易分解，抗蝕性差，因此土壤可蝕性對土壤侵蝕有一定影響；但該區域地形相對平緩，土壤有機質含量較高，此種情況下不易發生土壤侵蝕，故該區域屬于輕度侵蝕。（5）微度侵蝕主要分布在中南部的林地地區，約占侵蝕面積的30%，在侵蝕分類等級中所占侵蝕面積最大、侵蝕程度最低。該區域地勢低，土壤肥沃，植被覆蓋大，各侵蝕因子的影響均較小，因此屬于微度侵蝕。

4 結論

（1）太原市土壤侵蝕以微度和輕度所占比例最大，主要分布在中南部的林地，以及東部和南部的以梯田為代表的耕地地帶（坡耕地除外），約占市區總面積的1/2以上，土壤侵蝕較嚴重的地帶主要分布在北部以及邊緣地帶。總體上看，太原市大部分地區的土壤侵蝕程度相對較低、土地利用相對較充分，但亦有一部分地區的土壤侵蝕較為嚴重，盡管所占面積較小，但不可忽視。

（2）總體上看，荒草地、裸巖石礫地等未利用地的侵蝕強度最大，這類區域的植被覆蓋少，荒漠化嚴重，水土保持措施不到位，故侵蝕強度最為強烈；其次是建設用地、水體和耕地，侵蝕強度次之，處于中間層次；林地的侵蝕強度最小。

［1］De Jong S M.Applications of Reflective Remote Sensing for Land Degradation Studies in a Mediterranean Environment［M］.Utrecht:Universiteit Utrecht，1994:6－8.

［2］汪凌.美國航天飛機雷達地形測繪使命簡介［J］.測繪通報，2000，13（12）:38－40.

［3］劉寶元，謝云，張科利.土壤侵蝕預報模型［M］.北京:中國科學技術出版社，2001:12－18.

［4］章文波，付金生.不同類型雨量資料估算降雨侵蝕力［J］.資源科學，2003，25（1）:35－41.

［5］李晶，任志遠.基于GIS的陜北黃土高原土地生態系統水土保持價值評價［J］.中國農業科學，2007，40（12）:2796－2803.

［6］劉志紅.基于遙感與GIS的全國水土流失評價［D］.陜西楊凌:中國科學院 水利部 水土保持研究所，2007.

［7］楊勤科，趙牡丹，劉詠梅，等.DEM與區域土壤侵蝕地形因子研究［J］.地理信息世界，2009，7（1）:25－31.

［8］Hickey R，Smith A，Jankowski P.Slope Length Calculations from a DEM within ARC/INFO GRID［J］.Computers，Environment and Urban Systems，1994，18（5）:365－380.

［9］Van Remortel R，Hamilton M，Hickey R.Estimating the LS factor for RUSLE through iterative slope length processing of digital elevation data［J］.Cartography，2002，30（1）:27－35.

［10］張黎明.我國南方不同類型土壤可蝕性K值及相關因子研究［D］.海南 海口:華南熱帶農業大學，2005.

［11］Williams J R，Renard K G，Dyke P T.EPIC:A new method for assessing erosion′s effect on soil productivity［J］.Journal of Soil and Water Conservation，1983，38（2）:381－383.

［12］張巖，袁建平，劉寶元.土壤侵蝕預報模型中的植被覆蓋與管理因子研究進展［J］.應用生態學報，2002，13（8）:1033－1036.