北京市山区小流域治理前后土壤侵蚀强度及空间格局分析

赵猛，姚吉利，*，王建，胡信志，王家晖，张超，李彩林

1. 山东理工大学建筑工程学院，山东淄博 255049

2. 北京市水土保持工作总站，北京 100036

0 前言

北京山区作为首都重要的生态屏障，土壤侵蚀的防治问题不容忽视。北京山区坡度在10°以上的区域约占山区总面积的3/4，传统的坡度坡长因子计算中只分了三个区间，而对于坡度在 10°以上的区域并未充分考虑区间性。小流域梯田的水土保持措施因子也受到地形因子的影响，须得出适合于该小流域梯田的水土保持措施因子值。继 USLE[1]和RUSLE[2-3]提出后，国内外学者对流域坡面土壤侵蚀开展了大量探索性试验和深入研究，取得了一些显著性的研究成果: Kheir等[4]提出结合土壤可蚀性、地形特征、土地覆盖/利用、降雨侵蚀力与GIS技术开发区域侵蚀模型。唐磊等[5]应用景观生态学理论与 USLE分析岷江上游杂古脑流域 1990年和2015年的景观格局变化与土壤侵蚀变化。张养安等[6]基于GIS和RS技术，运用USLE评估杨凌区土壤侵蚀强度以及杨凌区土壤侵蚀的时空变异性。温莉等[7]以RS、GIS、e Cognition、SPSS 等软件为平台，将基于 RUSLE模型计算得到的年水土流失量数据与官方公布数据进行回归分析，建立重庆市土壤侵蚀模型。刘宝元等[8]结合北京山区土壤侵蚀发生的特点，建立了用于评价坡面土壤流失量、评价水土保持措施效益及其环境影响的北京土壤流失方程BJSLE。此外小流域水土保持效益的评价多依赖于实地考察取样分析，时间较长浪费人力物力，而基于遥感数据的快速效益评价较少[9]。

本文以中心村小流域水土保持措施本底调查为基础，通过利用高分辨率遥感影像，结合实地调查，选用像元二分模型[10]和K均值算法分别获取中心村小流域综合治理前后的植被覆盖、土地利用状况，结合这些基础数据，在BJSLE的基础上，改进了地形因子LS的计算方法，推算出中心村小流域梯田的水土保持工程措施因子P，定量评价中心村小流域治理前后(2015年、2017年)的土壤侵蚀状况，并分区进行对比分析与评价，得出水土保持措施效益情况。

1 研究区概况

中心村小流域位于北京市平谷区东北部(117°17'—117°21'E，40°12'—40°17'N，图 1)，隶属金海湖镇，流域面积 22.49 km2，项目区治理面积 10 km2；小流域属于燕山山脉，总体地势北高南低，东高西低，流域内以片麻岩、硅质白云岩、碎硝岩、碳酸盐岩类构成的单面山、馒头山为主，土壤以褐土为主，海拔在138 m—947 m之间。中心村小流域属暖温带半湿润大陆性季风区，多年平均降水量644 mm，降水年际变化大，时空分布极不均匀，降水多集中在 7、8、9三个月，汛期降水量占年降水的 83%左右。多年平均蒸发量1425 mm，平均气温11.6℃，极端最高气温40.2℃，极端最低气温－26.6℃。

2 研究方法

2.1 研究数据

图1 研究区位置及概况Figure 1 Location of the study region

本次研究使用的数据主要包括: 2015年7月和2017年7月的GF2号1米全色、4米多光谱影像数据一景，用于获取覆盖与管理因子和水土保持工程措施因子；中心村小流域综合治理状况基础数据、中心村小流域土壤质地资料、中心村小流域降雨资料，土壤质地资料结合全国二次普查结果获取土壤可蚀性因子，结合ArcGIS中的克里格插值法获取降雨侵蚀力因子。中心村小流域 DEM 数据通过北京市10米DEM数据裁剪获得，用于提取地形因子。

2.2 土壤侵蚀量的估算方法

本研究基于刘宝元[8]等提出的北京土壤流失方程BJSLE(Beijing soil loss equation)，在此基础上发展了适用于中心村小流域的土壤流失方程。其形式为:

式中，A为土壤侵蚀模数，t·(hm2·a)-1；K为土壤可蚀性因子，t·hm2·h·(hm2·MJ·mm)-1；L为坡长因子，无量纲(0—1)；S为坡度因子，无量纲(0—1)；C为覆盖与管理因子，无量纲(0—1)；R为降雨侵蚀力因子，MJ·mm·(hm2·h·a)-1；P为水土保持措施因子，无量纲(0—1)。

2.3 各因子提取方法

(1)降雨侵蚀力因子(R)

降雨是引起土壤侵蚀最直接的外部驱动因素[11-12]。降雨时雨滴所携带的动能对土壤颗粒进行冲击，产生击溅侵蚀，而后还产生径流，冲刷地表，再次导致土壤侵蚀的发生。R因子反映降雨引起土壤流失的潜在能力。R因子由小流域降雨量数据提取。本研究采用刘宝元等提出的北京土壤流失方程中建立的月降雨侵蚀力的估算公式[8]

式中，Pm为多年月平均降雨量，mm；r为回归方程决定系数。

(2)土壤可蚀性因子(K)

K因子是衡量土壤抗侵蚀性的指标，用于反映土壤对侵蚀的敏感性。影响K因子的因素是多方面的，一般说来，质地越粗或越细的土壤有较低的K值，而质地适中的反而有较高的K值[13]。结合土壤质地数据，K值估算采用 Williams等在EPIC 模型中的方法，从土壤质地数据提取，根据全国二次普查结果和中心村小流域的土壤质地状况得出中心村小流域主要土属土壤可蚀性因子K，见下表:

(3)地形因子(LS)

地形因子是影响土壤侵蚀的基本自然地理要素，是形成有侵蚀能力的径流的最主要因素。LS因子表示在其他条件不变的情况下，某给定坡长和坡度的坡面上土壤流失量与标准径流小区典型坡面上土壤流失量的比值，对土壤侵蚀起加强作用[14]。L因子从DEM数据中提取，其计算公式采用刘宝元[8]等提出的适用于北京山区的坡长因子修正公式:

表1 中心村小流域主要土属土壤可蚀性因子KTable 1 The main soil erosive factorKof the Zhongxin Village small watershed

式中，λ为像元坡长，m为坡长指数，随坡度发生变化，取值如下:m= 0 .2，θ＜ 0 .5°；m= 0 .3，0.5° ≤θ＜1.5°；m= 0 .4，1.5° ≤θ＜3°；m= 0 .5，3°≤θ。因子从DEM数据中提取，由于中心村小流域的坡度在 10°以上的区域所占比例为 75.59%，其计算公式采用刘滨涛[15]等修正的刘宝元算法，修正的坡度因子计算公式为:

(4)覆盖与管理因子(C)

中心村小流域的植被覆盖度估算的数据源选用GF2卫星遥感影像，遥感植被指数选用归一化差异植 被 指 数 (Normal Differential Vegetation Index，NDVI)，方法选用像元二分模型[10]。

植被是土壤侵蚀的主要抑制因素，C因子是这种抑制作用的定量化表达指标之一。C因子指一定植被覆盖条件耕作农地上的土壤流失量与同等条件下实施清耕的连续休闲地的土壤流失量之比，其值大小变化于0—1之间，C值越大，说明对应的土地利用类型的土壤侵蚀越严重[17]，其计算公式采用蔡崇法[18]等的研究结果

式中，VFC为植被覆盖度。

(5)水土保持工程措施因子(P)

P因子是采取水保措施后，土壤流失量与顺坡种植时的土壤流失量的比值，其值大小变化于 0—1之间[19]。在大尺度的土壤侵蚀研究中，通常采用结合实地考察给不同的土地利用类型赋予不同的P因子值。

中心村小流域在本次治理的水土保持工程措施中，梯田是北京山区最有效的水土保持措施之一，中心村小流域梯田整修共有两处: 中心村村南梯田10.86 hm2；中心村村西梯田，面积 6.59 hm2，共17.45 hm2，两处全部为石埂梯田。中心村村南梯田面净宽10—40m，田坎高度约为 0.6m—1.0m，坡度5°—8°。中心村村西梯田田面净宽3 m—12 m，田坎高度约为1 m—1.5 m，坡度10°—12°。由于该区域缺少实测资料，本研究将基于北京2004—2014年各坡面径流小区的实测资料，计算出适用于该区域梯田的水土保持工程措施因子。计算公式为

式中:Pm为措施小区的水土保持工程措施因子；Ac为对照小区的有效降雨年侵蚀量；Am为措施小区的有效降雨年侵蚀量。由于中心村小流域整修的梯田与北京山区的坡面径流小区的坡度坡长有所差异，需要根据公式将梯田的水土保持措施因子推算到该治理区域中，用下式

式中:Lm Sm为措施小区的坡长坡度因子；Lp Sp为治理区的坡长坡度因子。其中坡长坡度因子的计算方法同公式(3)、(4)。

其他水土保持措施中，村庄美化、挡土墙等措施也能够有效地减少小流域内的水土流失，参考前人研究[20-21]得出其他水土保持措施保土效益为: 村庄美化的年保土量为 18 t·(km2)-1；挡土墙的年保土量为 0.2 t·m-1。

中心村小流域其他不同土地利用类型的P值结合前人研究[22-23]整理见表2。

2.4 土壤侵蚀模数计算

基于 ArcGIS软件，利用土壤侵蚀因子计算值，运用北京土壤流失方程，对降雨侵蚀力R、土壤可蚀性因子K、地形因子LS、覆盖与管理因子C、水土保持措施因子P，进行图层栅格乘积运算。获取优于 10 m空间分辨率的土壤侵蚀模数计算值栅格图层。

3 结果与分析

3.1 中心村小流域综合治理前后土壤侵蚀强度分级

按照《北方土石山区水土流失综合治理技术标准》(SL665—2014)等技术标准，判断出每个土壤侵蚀地块的土壤侵蚀强度等级，北方土石山区容许土壤流失量定为 200 t·(km2·a)-1。统计微度(＜200)、轻度(200—1000)、中度(1000—2500)、强烈(2500—4000)、极强烈(4000—6000)和剧烈(＞6000)(分级单位t·(km2·a)-1)各级土壤侵蚀强度的面积，再计算其水土流失面积。中心村小流域治理前和治理后土壤侵蚀强度分级如图2。

表2 不同土地利用类型P值Table 2 P value of different land use types

中心村小流域治理前和治理后，对不同等级土壤侵蚀强度的侵蚀面积、侵蚀总量和侵蚀模数指标进行对比，见表3和图3。

治理前(2015年)，中心村小流域土壤侵蚀总量为 12565.36 t，平均侵蚀模数为 558.71 t·(km2·a)-1，属于轻度侵蚀。从侵蚀量上来看，轻度侵蚀和中度侵蚀所产生的侵蚀量较大，分别占土壤侵蚀总量的44.59%和 44.71%；其次是强烈侵蚀和微度侵蚀，分别占土壤侵蚀总量的5.40%和4.92%；极强烈侵蚀所产生的侵蚀量最小，占比为0.39%。从侵蚀面积上来看，轻度侵蚀面积最大，占小流域总面积的 45.47%；紧随其后的是微度侵蚀面积，占小流域总面积的35.88%；中度侵蚀面积、强烈侵蚀面积和极强烈侵蚀面积较小，分别占比为17.56%、1.02%、0.07%。从侵蚀分布上来看，微度侵蚀主要分布在坡度较小的河流流经处的平坦地带；轻度侵蚀主要分布在植被覆盖度较高的低海拔山坡，山谷和山口处；中度和强烈侵蚀则主要分布在植被覆盖度较低的高海拔地带和沟道分布较多的小流域西北部。

图2 中心村小流域治理前和治理后土壤侵蚀强度分级Figure 2 Grading of soil erosion intensity before and after treatment in the small watershed of the Zhongxin Village

表3 中心村小流域治理前和治理后不同等级土壤侵蚀强度对比Table 3 Comparison of soil erosion intensity of different grades before and after treatment in small watershed of Zhongxin Village

图3 不同等级土壤侵蚀强度面积比率和侵蚀量比率对比Figure 3 Comparison of soil erosion intensity area ratio and erosion ratio of different grades

治理后(2017年)，中心村小流域土壤侵蚀总量为 11167.85 t，平均侵蚀模数为 496.57 t·(km2·a)-1，属于轻度侵蚀，达到了治理目标。从侵蚀量上来看，轻度侵蚀所产生的侵蚀量最大，占土壤侵蚀总量的50.26%；中度侵蚀量较其他等级侵蚀量也比较大，占比为 41.12%，其次是微度侵蚀和强烈侵蚀，分别占土壤侵蚀总量的5.80%和2.72%；极强烈侵蚀所产生的侵蚀量最小，占比为0.10%。从侵蚀面积上来看，轻度侵蚀面积最大，占小流域总面积的 45.85%，主要分布在坡度较小的河流流经处的平坦地带；紧随其后的是微度侵蚀面积，占小流域总面积的 38.57%，主要分布在植被覆盖度较高的低海拔山坡、山谷和山口处；中度侵蚀面积、强烈侵蚀面积和极强烈侵蚀面积较小，分别占比为 14.94%、0.63%、0.01%，主要分布在植被较稀疏的高海拔地带。

3.2 中心村小流域综合治理前后土壤侵蚀强度变化及空间格局分析

通过对比治理前和治理后中心村小流域土壤侵蚀量(表3)，治理后的中心村小流域的土壤侵蚀总量比治理前减少了1397.51 t，综合治理前年平均土壤侵蚀模数为 558.71 t·(km2·a)-1，综合治理后年平均土壤侵蚀模数为 496.57 t·(km2·a)-1，年平均土壤侵蚀模数减小了 62.14 t·(km2·a)-1。从各级侵蚀量角度来看，综合治理后的微度侵蚀量和轻度侵蚀量较综合治理前有所增加，分别增加了 29.81 t、9.99 t；而中度侵蚀量、强烈侵蚀量和极强烈侵蚀量有所减少，其中中度侵蚀量减少幅度最大，减少了 1026.08 t。从各级侵蚀面积来看，治理后微度侵蚀面积和轻度侵蚀面积都有增加；中度侵蚀、强烈侵蚀、极强烈侵蚀面积均有减小。通过对综合治理前后各等级侵蚀面积结果进行转移矩阵(表4)分析，侵蚀等级正逐步由高向低过渡。由此可见通过综合治理后的中心村小流域土壤侵蚀状况正在逐步改善，侵蚀量减少。

中心村小流域综合治理前后的土壤侵蚀量的变化跟植被覆盖度有较大的关系，按照植被覆盖度分级标准将小流域的植被覆盖度情况划分为五个等级(表5)，各级覆盖面积对比情况见图4，由图4可以看出，低覆盖、中低覆盖、中覆盖所占得面积相对于治理前都有所减少，中高覆盖和高覆盖所占面积有所增加，其中中高覆盖面积增加显著，由治理前的4.31 km2增加到治理后的7.88 km2。从空间分布上来看，小流域坡中、坡下部和坡脚地区的植被覆盖度变化明显，有利于减少此地区的土壤流失量，植被覆盖度的变化也影响了覆盖与管理因子C的分布与大小。

在中心村小流域的土地利用类型中，林地所占面积最大，耕地、建筑用地等主要分布在流域河道两侧。通过中心村小流域治理前后土地利用结果对比可以看出，小流域经过综合治理，该流域建筑用地、交通运输用地和水体面积保持不变，林地面积有所增加，增加了 0.50 km2；但裸地面积和耕地面积分别减少了0.09 km2和0.34 km2。通过对两期土地利用结果进行矩阵分析(表 6)，可以看出: 综合治理前的耕地，有0.50 km2转化成了林地，有0.18 km2经过修筑和整治转化成了梯田；而裸地则主要转化为了耕地。由此可以看出，土地利用现状的转化也影响了水土保工程持措施因子 P，特别是耕地的转化与整治，很大程度上减少了水土的流失。

表4 中心村小流域综合治理前后侵蚀等级面积转移矩阵Table 4 Erosion level transfer matrix before and after comprehensive management of small watershed in Zhongxin Village

表5 植被覆盖度等级划分Table 5 Vegetation coverage classification

图4 中心村小流域治理前后植被覆盖度面积Figure 4 Vegetation coverage area before and after treatment of small watershed in Zhongxin Village

从空间分布格局来看，综合治理将小流域划分为生态修复区、生态治理区和生态保护区。通过对比分析(图5)生态修复区的侵蚀量由治理前的10207.55 t减小到治理后的 9260.91 t，年平均侵蚀模数由 689.70 t·(km2·a)-1减小到 625.74 t·(km2·a)-1，生态修复区主要位于小流域的山顶和坡上部，坡度一般都在25°以上，面积为14.80 km2，占小流域总面积的65.81%，该区域主要实施封山育林治理，因为该区域人类活动较少，没有开发建设及大规模的农业生产活动等人为干扰，通过封山育林能将人类活动趋近于 0，靠自然修复达到水土保持效果。生态治理区的侵蚀量由治理前的2046.00 t减小到治理后的1836.88 t，年平均侵蚀模数由277.61 t·(km2·a)-1减小到 249.24 t·(km2·a)-1，生态治理区主要位于坡中、坡下部和坡脚地区，坡度在 25°以下，面积为7.73 km2，占小流域总面积的32.77%，该区域人类活动频繁，生态环境干扰、破坏程度相对较大，此为小流域综合治理的重点区域，该区域治理的主要实施包括梯田整修、节水灌溉、挡土墙、村庄美化、田间生产道路治理，通过数据对比分析，该区域的治理措施成效显著，梯田整修能进一步的减少水土流失，提高农田的保土能力；节水灌溉能减少因灌溉带来的泥土冲刷；挡土墙能够有效的遏制裸露边坡坍塌，减少水土流失；村庄美化和田间生产道路的整修也有效的减少了水土流失。生态保护区的侵蚀量由治理前的76.61 t减小到治理后的 70.07 t，年平均侵蚀模数由239.40 t·(km2·a)-1减小到 218.96 t·(km2·a)-1，该区域主要位于河道周边地带，面积为 0.32 km2，占小流域总面积的 1.42%，河道周边实施防护措施，能减少河道冲刷带走的泥沙量。

4 讨论

坡度是地貌形态特征的主要要素，是影响坡面土壤侵蚀的重要因素[24]，北京山区不同坡度的侵蚀量存在一定的差异，但总体上表现出随着坡度的增加而增大的趋势，这是因为在相同的条件下，随着地面坡度增加，径流速度加快，径流在坡面滞留时间变短，使得入渗量减少，坡面径流量增加，加剧了土壤侵蚀过程[25]。这与本文研究结果类似，中心村小流域的坡度在 10°以上的区域所占比例为75.59%，坡度相差较大，由图2可知，中度及以上侵蚀大都分布在坡度较大的区域，并且随着坡度的增大，侵蚀强度也随之增大，再者生态修复区主要位于小流域的山顶和坡上部，坡度一般都在 25°以上，治理前后的侵蚀量均占总侵蚀量的 80%以上，更能说明研究区坡度的增加侵蚀过程也随之加剧。坡面土地利用方式及水土保持综合措施的合理布置是控制水土流失最重要的影响因素，通过合理的土地利用方式与适当的坡面水土保持措施可以有效地防治土壤侵蚀[25-26]。中心村小流域的治理措施成效显著，各项措施都发挥了应有的作用，梯田作为北京山区最有效的水土保持措施之一，在中心村小流域综合治理中发挥作用最大。本文研究结果从总体上来看，侵蚀量和侵蚀面积的变化并不是很大，究其原因，本研究是基于治理后一年来进行分析评价，各项措施完成时间较短，只体现短期效益，并没有发挥治理措施的长远效益，所以应当对治理区进行长期的观测与评价，才能体现综合治理措施最终的效益状况。

表6 中心村小流域综合治理前后土地利用转移矩阵Table 6 Land use transfer matrix before and after comprehensive management of Zhongxin Village small watershed

图5 中心村小流域三区划分治理前后平均侵蚀模数Figure 5 The average erosion modulus before and after the division of the three districts in the Zhongxin Village

5 结论

本文在北京土壤流失方程(BJSLE)的基础上，通过改进地形因子(LS)、梯田的水土保持工程措施因子(P)的计算方法，结合 GIS、RS，定量的评价了中心村小流域综合治理前后的土壤侵蚀状况，并从植被覆盖度、土地利用现状和三区划分的角度对综合治理前后土壤侵蚀时空变化格局进行了探讨，最后得出以下结论:

(1)中心村小流域综合治理前后年平均土壤侵蚀模数分别为 558.71 t·(km2·a)-1、496.57 t·(km2·a)-1，属于轻度侵蚀。小流域治理前后大部分的侵蚀量均来源于轻度侵蚀和中度侵蚀，分别达到总侵蚀量的89.30%和90.58%，所以应该进一步加强对该侵蚀区域的治理。

(2)中心村小流域综合治理前后的土壤侵蚀量有一定的变化，土壤的侵蚀级别正在从高级别向低级别转化，特别是强烈侵蚀和极强烈侵蚀的侵蚀量减少了治理前的一半以上。

(3)中心村小流域三区划分治理，因地制宜，不同治理区的侵蚀量均有较大变化，治理措施的效益达到了预期的效果。