土壤侵蚀破坏地表植被，造成土壤养分流失，降低土地生产力，导致土地资源破坏、流域水质污染、水库泥沙堆积和洪涝灾害，危及生态安全，是全球性的生态环境问题之一^[1-3]。USLE和RUSLE模型的大力发展^[4]，广泛应用于土壤侵蚀研究，多在坡面和小流域尺度上进行。当前国内学者对于流域土壤侵蚀时空变化特征以本地小流域的研究较多，且大多局限于短时间尺度的预测，对于大型流域长时间尺度的研究较少。蒸水流域曾经是湖南省水土流失严重的地区之一，流域内土壤母质风化程度高，防洪、排涝设施少，洪水来临导致河水漫滩，造成较大财产损失^[5]。为了遏制水土流失，当地相关部门在蒸水流域开展了水土保持和水环境综合治理工作，如蒸水南堤草桥-蒸水桥防洪工程建设等。要了解流域土壤侵蚀历史及现状，需要在蒸水流域土壤侵蚀状况进行定量评价与分析。笔者基于土壤侵蚀的基本理论，利用GIS和RS技术，将CSLE模型与长时序地面观测数据结合，定量评价蒸水流域1995、2000、2005、2010和2015年的土壤侵蚀情况，了解流域内土壤侵蚀环境背景和土壤侵蚀差异，掌握其水土流失整体现状，研究结果可望为水土流失防治、水土保持规划和技术管理提供科学参考，以便更好的利用水土资源。

蒸水属洞庭湖水系湘江的一级支流，流域介于E 111°53′~112°37′、N 26°52′~27°10′之间，流域面积3480km²，河长194km，形状呈“乙”字型，河流坡降0.54‰。北、东、南3面环山，山地丘陵为2720km²，约占流域总面积的71%；中、下游地区为平原、低地，面积为750km²，约占总面积的20%^[6]；其余为水域。流域内及流域周边有双峰站、邵东站等9个气象站，神山头站为水文站(图 1)。

图 1(Fig. 1)

图 1 蒸水流域高程图 Fig. 1 Elevation map of Zhengshui River Basin

研究采用中国土壤流失方程CSLE模型来估算蒸水流域的土壤侵蚀量，CSLE是Liu等^[7]基于我国土壤侵蚀现状，以USLE和RUSLE为原型，改进提出的适用于我国的预测坡耕地年土壤流失的经验模型。该模型结构简单、参数要求低，也适用于沟壑纵横，山高坡陡地区。模型为

式中：A为土壤侵蚀模数，表示多年平均土壤流失速率，t/(hm²·a)；R为降雨侵蚀力因子，MJ·mm/(hm²·h·a)；K为土壤可蚀性因子，t·hm²·h/(hm²·MJ·mm)；L为坡长因子；S为坡度因子；B为植被覆盖与生物措施因子；E为工程措施因子；T为耕作措施因子。参考水利部办公厅文件办水保(2018)189号文件《关于印发区域水土流失动态监测技术规定(试行)的通知》及水利部监测中心颁布的技术规定^[8]：R因子采用殷水清等^[9]的利用冷暖季日雨量资料计算半月降雨侵蚀力公式，累加半月降雨侵蚀力得到年降雨侵蚀力，经过普通克里金空间插值后，形成蒸水流域降雨侵蚀力因子分布图。K因子数据由刘宝元教授提供，计算方法采用的国务院第1次全国水利普查计算K值方法^[10]。LS因子采用Liu等^[7]的CSLE算法。

参照《区域水土流失动态监测技术规定(试行)》附录7的要求，并结合前人对B因子的研究和对华中区林下盖度的野外调查成果^[11-13]，得到不同用地类型的B因子值，其中林地和草地的B因子利用TM多光谱影像(1995—2010年利用Landsat4-5 TM卫星影像，2015年利用Landsat 8 OLI_TIRS卫星影像)中红波段和近红外波段，计算归一化植被指数(normalized difference vegetation index，NDVI)，基于NDVI数据计算得到30m空间分辨率的植被覆盖度(fractional vegetation cover，FVC)。

灌木林地和草地的B值采用计算公式：

果园、其他园地、有林地和其他林地的B值采用计算公式：

式中：GD为乔木林的林下盖度。根据前人对华中区林下盖度的研究和野外调查^[12-13]，分别拟定有林地、疏林地和其他林地的林下盖度值，代入到上式计算。耕地、城镇建设用地、农村建设用地、其他建设用地、水域、裸土地和裸岩的B值分别直接赋值为1.000、0.010、0.025、0.100、0、1.000和0^[8]，形成蒸水流域植被覆盖与生物措施因子图。

蒸水流域内工程措施主要为土坎水平梯田，参照《区域水土流失动态监测技术规定(试行)》，将水田赋E因子值为0.084，其余用地类型赋值为1.000。蒸水流域所在区域在全国轮作区内均处于2湖平原丘陵水田中三熟二熟区，将耕地T因子赋值为0.312，其他用地类型赋值为1。形成蒸水流域工程措施因子图和耕作措施因子图。

由于土壤侵蚀变化受到了自然和人类活动的综合作用影响，其中自然因素主要考虑于降水变化，人类活动主要考虑土地利用变化和水土保持措施实施的变化。为了区分降水和人类活动对蒸水流域土壤侵蚀的不同影响，利用控制变量方法，分别模拟只考虑降水变化和只考虑土地利用变化情景下2015年的土壤侵蚀状况^[14]。具体步骤如下：1)不考虑其他因素变化，只考虑降水变化情景下的2015年土壤侵蚀模数：A_2015R=R₂₀₁₅KLSB₁₉₉₅E₁₉₉₅T₁₉₉₅；2)不考虑其他因素变化，只考虑土地利用变化情景下的2015年土壤侵蚀模数：A₂₀₁₅=R₁₉₉₅KLSB₂₀₁₅E₂₀₁₅T₂₀₁₅；3)同时考虑降水变化和土地利用变化情景下的2015年土壤侵蚀模数：A_2015RBET=R₂₀₁₅KLSB₂₀₁₅E₂₀₁₅T₂₀₁₅；4)将计算的1995年实际土壤侵蚀模数与模拟情景下的2015年土壤侵蚀模数进行对比。

基于CSLE的因子算法，利用1995—2015年的5个年度数据，获得了蒸水流域30m分辨率R、K、LS、B、E、T因子图。

统计蒸水流域周边9个气象站每年侵蚀性降雨的出现时间，绘制柱状图(图 2)。1995—2015年各年的平均侵蚀性降雨时间分别为32.6、44.6、44.1、47.9和48.3d；非侵蚀性降雨时间分别为332.4、321.4、320.9、317.1和316.7d。1995年侵蚀性降雨占全年比例最小，而2010和2015年侵蚀性降雨占全年比例最大。

图 2(Fig. 2)

图 2 蒸水流域各气象站5个年度各级侵蚀性降雨全年出现次数 Fig. 2 Frequency of different erosive rainfall of 5 periods in each climate station of Zhengshui River Basin

如表 1与图 3所示，2010年R值最大，均值为5934MJ·mm/(hm²·h·a)，地区间年降雨侵蚀力分布差异最大，标准差为1280MJ·mm/(hm²·h·a)；1995年R值最小，均值4275MJ·mm/(hm²·h·a)，标准差为215MJ·mm/(hm²·h·a)；2005年R值为5177MJ·mm/(hm²·h·a)，降雨侵蚀力空间分布最为均匀，标准差为199MJ·mm/(hm²·h·a)。

图 3(Fig. 3)

图 3 蒸水流域5个年度降雨侵蚀力因子图 Fig. 3 Rainfall erosivity factor map of 5 periods in Zhengshui River Basin

表 1(Tab. 1)

表 1 蒸水流域降雨侵蚀力统计特征值 Tab. 1 Statistics of rainfall erosivity in Zhengshui River Basin   MJ·mm/(hm2·h·a) 年份 Year 最大值 Max 最小值 Min 平均值 Mean 标准差 Std Dev 1995 4763 3835 4275 215 2000 6908 4128 5526 589 2005 5602 4183 5177 199 2010 9186 4285 5934 1280 2015 6169 5095 5528 225

表 1 蒸水流域降雨侵蚀力统计特征值 Tab. 1 Statistics of rainfall erosivity in Zhengshui River Basin

蒸水流域表层土壤可蚀性因子K值主要介于0~0.006748t·hm²·h/(hm²·MJ·mm)之间，均值为0.004333t·hm²·h/(hm²·MJ·mm)。高值区主要分布在流域南部、东南部以及中部区域，低值区主要分布在流域西北部地区。造成这种空间异质性的原因主要取决于土壤类型的异质性；高值区的主要土壤类型为水稻土，低值区的主要土壤类型为石灰土(图 4)。

图 4(Fig. 4)

图 4 蒸水流域土壤可蚀性因子K值图 Fig. 4 Soil erodibility factor K map in Zhengshui River Basin

利用DEM数据提取得到蒸水流域坡度分布图和地形因子图(图 5和6)蒸水流域LS值均值为4.45。流域东北部、北部、西南山区LS值较高，中部、南部与东南部LS值较低，呈马蹄形分布。

图 5(Fig. 5)

图 5 蒸水流域坡度分布图 Fig. 5 Slope map of Zhengshui River Basin

图 6(Fig. 6)

图 6 蒸水流域坡度坡长因子图 Fig. 6 Slope and slope length factor LS map of Zhengshui River Basin

1995、2000、2005、2010和2015年蒸水流域的年均B因子分别为0.4511、0.4507、0.4496、0.4485和0.4456。B因子逐渐减小，年际差异主要来源于土地利用变化情况。西南、北部区域B值较低，东南部区域B值较高(图 7)。

图 7(Fig. 7)

图 7 蒸水流域5个年度植被覆盖与生物措施因子图 Fig. 7 Vegetation cover and biological practice factor B map during 5 periods in Zhengshui River Basin

1995年蒸水流域内以耕地和林地为主，其面积分别占总面积的42.79%和54.16%，呈现以农、林地为主的土地利用结构特征。此后随着退耕和建设占用，耕地面积逐年下降。到了2015年，耕地和林地为主要用地类型的格局没有变化，分别占流域面积的42.22%和53.89%。相比于1995年，水田面积占比减少0.37%，旱地面积占比减少0.19%，林地面积比例减少0.27%，而城镇用地、农村居民点等建设用地面积比例增加0.80%，这种变化是流域内城市化发展的结果。

1995—2015年蒸水流域的年均E因子分别为0.7109、0.7115、0.7119、0.7121和0.7143，西南、北部区域E值较高，中部与东南部区域E值较低。E因子的空间异质性来源于梯田，而西南及北部地区为非水田区域，工程措施因子都为1(图 8)。1995—2015年蒸水流域的年均T因子分别为0.7056、0.7061、0.7068、0.7071和0.7095。西南、北部区域T值较高，中部与东南部区域T值较低。T因子的空间异质性来源于2湖平原丘陵轮作区的耕地，西南及北部地区用地类型为非耕地，耕作措施因子都为1(图 9)。

图 8(Fig. 8)

图 8 蒸水流域2015年工程措施因子E值图 Fig. 8 Engineering practice factor E map of Zhengshui River Basin in 2015

图 9(Fig. 9)

图 9 蒸水流域2015年耕作措施因子T值图 Fig. 9 Tillage practice factor T map of Zhengshui River Basin in 2015

运用CSLE模型计算得到蒸水流域1995、2000、2005、2010和2015年的土壤侵蚀模数(图 10)。各年的平均土壤侵蚀模数分别为412、520、479、530和528t/(km²·a)，1995与2005年属于微度侵蚀等级，2000、2010和2015年属于轻度侵蚀等级。根据水利部行业标准SL190—2007《土壤侵蚀分类分级标准》规定，蒸水流域所在的南方红壤丘陵区土壤侵蚀模数所推荐的容许土壤流失量为500t/(km²·a)，这5年的平均侵蚀模数为494t/(km²·a)，低于容许土壤流失量，可以判定蒸水流域水土流失情况较轻。

图 10(Fig. 10)

图 10 蒸水流域5个年度土壤侵蚀模数空间分布 Fig. 10 Spatial distribution of soil erosion modulus during 5 periods in Zhengshui River Basin

利用ArcGIS软件统计功能，得到各土壤侵蚀强度等级的侵蚀面积与比例。按照水利侵蚀强度分级对土壤侵蚀模数进行划分，蒸水流域1995—2015年发生的主要侵蚀类型为微度侵蚀，均占到全流域面积的70%以上；轻度侵蚀区域均占全流域面积的20%以上；中度侵蚀区域占全流域面积的2%以内；强烈及以上等级侵蚀区域占全流域面积的1%以内(表 2)。将土壤侵蚀强度空间分布图(图 10)与流域坡度分布图(图 5)比较，发现土壤侵蚀强度在空间上的分布与坡度区间在空间上的分布相似，坡度较大的区域，土壤侵蚀强度较高。

表 2(Tab. 2)

表 2 蒸水流域5个年度各土壤侵蚀强度的面积 Tab. 2 Area of each soil erosion intensity during 5 periods in Zhengshui River Basin   km2 土壤侵蚀强度 Soil erosion intensity 1995 2000 2005 2010 2015 微度Slight 2675 2542 2574 2500 2559 轻度Light 735 835 834 880 826 中度Moderate 50 69 48 66 60 强烈Intense 10 19 13 21 19 极强烈Severe 7 12 8 11 11 剧烈Violent 3 3 4 3 6

表 2 蒸水流域5个年度各土壤侵蚀强度的面积 Tab. 2 Area of each soil erosion intensity during 5 periods in Zhengshui River Basin

2015年与1995年相比，时间跨度较长，R值与B、E和T值差异较大，A值差异也较大。通过模拟考虑降水变化和考虑土地利用变化情景下2015年的土壤侵蚀状况，结果发现：1)只考虑降雨变化情景下的2015年单位面积土壤侵蚀量均值为526t/(km²·a)；2)只考虑土地利用变化情景下的2015年单位面积土壤侵蚀量均值为412t/(km²·a)；3)同时考虑降水和土地利用变化情景下的2015年单位面积土壤侵蚀量均值为528t/(km²·a)。通过对比可以看出：从1995年到2015年，降水变化加剧了蒸水流域的土壤侵蚀，土地利用变化对土壤侵蚀变化不明显，其中降水变化使流域侵蚀量增加27.67%，土地利用变化使流域侵蚀量增加0.18%，降水和土地利用变化共同作用使侵蚀量增加28.16%。

分区统计不同用地类型的土壤侵蚀情况(表 3)。流域内面积大小占比靠前的土地利用类型分别为有林地(35.94%)、水田(31.24%)、疏林地(15.20%)、旱地(11.00%)，以上4类用地占到研究区总面积的93.38%。土壤侵蚀强度从大到小分别为旱地 > 有林地 > 疏林地 > 水田，旱地土壤侵蚀模数大于有林地土壤侵蚀模数，主要原因是旱地较有林地的植被覆盖度小，其B值远大于有林地的B值；有林地土壤侵蚀模数大于疏林地土壤侵蚀模数，主要原因是有林地的平均海拔高、坡度大，其LS值远大于疏林地的LS值；疏林地土壤侵蚀模数略大于水田土壤侵蚀模数，主要原因是水田在工程措施与耕作措施的共同影响下，保土效果略优于植被覆盖与生物措施对疏林地的影响。

表 3(Tab. 3)

表 3 蒸水流域5个年度主要用地类型与土壤侵蚀模数 Tab. 3 Soil erosion modulus under different land use type during 5 periods in Zhengshui River Basin 主要用地类型 Main land use type 面积 Area/km2 比例 Percentage/% 年均土壤侵蚀模数Annual soil erosion modulus/(t·km-2·a-1) 1995 2000 2005 2010 2015 有林地Forest land 1251 35.94 542 645 642 667 670 水田Paddle field 1088 31.24 128 220 150 229 205 疏林地Sparse woodland 529 15.20 314 328 349 337 282 旱地Rain fed cropland 383 11.00 980 1298 1189 1316 1249

表 3 蒸水流域5个年度主要用地类型与土壤侵蚀模数 Tab. 3 Soil erosion modulus under different land use type during 5 periods in Zhengshui River Basin

根据模拟估算结果表明，土壤侵蚀现象在蒸水流域并不常见。蒸水流域涉及的行政区面积前3的区县为衡阳县、邵东县和衡南县，面积占比分别为64.98%、16.47%和11.84%，其余各县区总面积比例不足7%(表 4)。南方红壤丘陵区所推荐的容许土壤流失量为500t/(km²·a)，衡阳县仅在2010年达到轻度侵蚀等级，平均侵蚀模数为462t/(km²·a)；衡南县在2000年、2005年和2015年达到轻度侵蚀等级，平均侵蚀模数为517t/(km²·a)；邵东县在5个年度均为轻度侵蚀，平均侵蚀模数为591t/(km²·a)。流域内土壤侵蚀模数>500t/(km²·a)的地区集中分布在邵东县东部，衡阳县西部、北部、东北部，衡南县西北部，南岳区，以上地区是未来水土保持的重点地区。

表 4(Tab. 4)

表 4 蒸水流域5个年度各县区面积与土壤侵蚀模数 Tab. 4 Soil erosion modulus and area of each county or district under different region during 5 periods in Zhengshui River Basin 行政区 Region 面积 Area/km2 面积比例 Percentage/% 年均土壤侵蚀模数Annual soil erosion modulus/(t·km-2·a-1) 1995 2000 2005 2010 2015 衡南县Hengnan county 411.92 11.84 449 631 513 459 535 衡山县Hengshan county 14.12 0.41 517 670 652 970 640 衡阳县Hengyang county 2 260.83 64.98 373 492 429 518 497 南岳区Nanyue district 11.81 0.34 1406 1799 1739 1988 1593 祁东县Qidong county 15.89 0.46 544 650 788 919 838 邵东县Shaodong county 573.07 16.47 553 552 639 571 639 石鼓区Shigu district 61.97 1.78 240 399 267 351 276 双峰县Shuangfeng county 54.71 1.57 506 593 530 762 555 蒸湘区Zhengxiang district 74.89 2.15 190 288 233 287 371 全流域Whole basin 3 479.26 100.00 412 520 479 530 528

表 4 蒸水流域5个年度各县区面积与土壤侵蚀模数 Tab. 4 Soil erosion modulus and area of each county or district under different region during 5 periods in Zhengshui River Basin

神头山站是蒸水流域的水文站，其控制面积为2857km²，约占全流域面积的82.10%，利用ArcGIS分区统计工具来计算神头山水文站控制面积的侵蚀模数，并将计算侵蚀量与神头山站观测的输沙量进行对比(表 5)。

表 5(Tab. 5)

表 5 神山头站控制范围内外各年统计特征值 Tab. 5 Statistical characteristic values of each year within and outside the control area of Shenshantou Station 年份 Year 控制范围侵蚀模数 Soil erosion modulus under monitoring range/(t·km-2·a-1) 实测输沙量 Observed sediment/(104t) 输沙模数 Sediment modulus/(t·km-2) 计算侵蚀量 Calculated soil erosion/(104t) 泥沙输移比 Sediment delivery ratio 1995 427 26.4 92.4 121.99 0.22 2000 528 26.4 92.4 150.85 0.18 2005 489 12.4 43.4 139.71 0.09 2010 552 11.4 39.9 157.71 0.07 2015 536 10.2 35.7 153.14 0.07

表 5 神山头站控制范围内外各年统计特征值 Tab. 5 Statistical characteristic values of each year within and outside the control area of Shenshantou Station

结果表明1995—2015年神山头水文站控制范围侵蚀量呈波动变化，输沙量逐渐减少，由1995年的26.4万t减少到2015年的10.2万t，各年侵蚀量远远大于该年输沙量，泥沙输移比逐年下降，由0.22下降到了0.07。神山头站控制范围内输沙量远小于侵蚀量，实测输沙量逐年降低。究其原因：一是因为本研究的计算侵蚀量是基于坡面模型的每个单元侵蚀量的总和，没有考虑每个单元侵蚀量在输移过程中的沉积^[15]，神头山水文站的控制面积占全流域的80%以上，侵蚀土壤需要通过长距离运输，输移过程中有大量侵蚀土壤发生沉积，绝大部分侵蚀土壤因就地沉积或沿途沉积等因素最终未能抵达神头山站；二是由于流域内林草覆盖度达到50%以上，植被覆盖除减少坡面侵蚀且拦截侵蚀土壤的作用显著；三是流域除梯田外还有一些其他水土保持工程措施，在本次评价中未考虑，他们也能拦截侵蚀土壤^[16]。

1) 蒸水流域土壤侵蚀在空间上呈现由西北至东南逐渐减弱的特征。5个年度的侵蚀模数分别为412、520、479、530和528t/(km²·a)，呈现波动上升的趋势。5个年度的平均侵蚀模数为494t/(km²·a)，低于容许土壤流失量，属于微度侵蚀等级，全流域平均水土流失情况较轻，但存在局部侵蚀严重的情况。

2) 侵蚀性降雨发生频次较高，侵蚀性降雨变化是蒸水流域近20年侵蚀时空变化最主要的驱动因素，降水变化使流域侵蚀量增加了26.67%。

3) 蒸水流域主要用地类型的土壤侵蚀情况为旱地＞有林地＞疏林地＞水田。

4) 神山头水文站控制范围内侵蚀计算结果远大于实测输沙量，输沙量逐年降低；流域内高覆盖度的林草，水土保持工程参与拦截侵蚀土壤过程，使绝大部分侵蚀土壤就地或沿途发生沉积，最终未能到达水文站泥沙出口，使泥沙输移比很小，而且逐年下降。

5) 流域内土壤侵蚀模数>500t/(km²·a)的地区集中分布在邵东县东部，衡阳县西部、北部、东北部，衡南县西北部，南岳区，以上地区是未来水土保持的重点地区。

本次研究原基于研究区土壤粒径组成以及土壤有机碳含量，采用Williams的EPIC模型法^[17]对蒸水流域的土壤可蚀性K值进行计算，但直接利用EPIC公式估算K值与实测值存在较大差异，最大可达10倍，与其他学者研究结果类似^[18]，因此直接使用了现有的土壤可蚀性研究成果评价了蒸水流域的土壤侵蚀。