基于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.