0 引 言
区域土壤侵蚀状况直接对人类健康和经济社会发展产生重要影响[1]，定期开展水土流失动态监测可为水土流失状况评价、水土保持治理工程空间布局等提供重要支撑[2]。土壤侵蚀模型是开展水土流失动态监测的重要技术手段[3-4]。迄今为止，美国通用土壤流失方程（Universal Soil Loss Equation，USLE）及其改进模型（Revised Universal Soil Loss Equation，RUSLE）仍在土壤侵蚀领域广泛应用[5-6]；该模型综合土壤侵蚀多个影响因子及其相互作用关系，同时各因子具有一定物理意义且因子获取方法及模型计算过程较为简单[7-8]。Liu等[9-10]利用黄土高原丘陵沟壑区径流小区实测资料构建了中国土壤流失方程（China Soil Loss Equation，CSLE），该模型依据中国土壤侵蚀与水土保持工作特点，充分考虑了生物措施、工程措施和耕作措施对土壤侵蚀及水土流失过程的影响，与USLE相比更能反映中国土壤侵蚀的实际情况[11]。
2018年水利部印发《区域水土流失动态监测技术规定（试行）》，中国水土流失动态监测预报工作从定量与定性相结合的“三因子”综合判别法发展到以定量为主的“CSLE模型法”[12]。2018年结果表明中国仍有水土流失面积273.69万km2，约占国土总面积的28.54%[13]。CSLE模型广泛应用对中国水土流失动态监测应用工作和土壤侵蚀基础研究具有十分显著的促进作用，但也存在诸如模型计算精度偏低、模型因子区域适用性较差、地形及土壤等基础数据现势性不足等亟需解决的问题[14-15]。作为影响土壤侵蚀最重要因素之一：植被覆盖因子（B），可直接反映植被覆盖变化状况或植被措施调控土壤侵蚀的能力，也是CSLE模型中的关键因子[16-19]。植被覆盖因子指一定条件下有植被覆盖或实施田间管理的土地的土壤流失总量与同等条件下实施清耕休闲地的土壤流失总量的比值，是介于0和1之间的一个无量纲数[20]。植被覆盖度大小是影响植被覆盖因子的关键因素；此外，植被覆盖因子大小还受到植被不同生长期内侵蚀性降雨多寡的影响，对于植被覆盖因子年值的研究还应充分考虑降雨侵蚀力年内季节分布的影响。标准小区法、次因子法、模型反推法、与植被覆盖度关系式法、直接赋值或手册查询法等是较为常用的植被覆盖因子获取方法。前3种方法依赖野外实测数据，逻辑性强、计算精度高；关系式法虽难以直接反映植被减蚀作用机理，但其经验关系具有相对较强的通用性和适用性；直接赋值或手册查询法较为依赖大量丰富的外业数据观测成果，通用性相对较弱。目前中国年度水土流失动态监测项目采用《区域水土流失动态监测技术规定（试行）》推荐的B因子计算方法，利用监测年前3年内每年24个半月MODIS数据计算得到B因子值。植被覆盖因子与土地利用方式关系密切，合理的土地利用方式可显著降低植被覆盖因子值而减少土壤侵蚀量[21-22]。张雪花等[23]基于人工模拟降雨试验确定了东北黑土农业区不同土地利用植被因子的取值范围，发现该值普遍小于USLE方程推荐的最优值。蔡崇法等[24]基于野外实测资料构建了三峡库区植被因子与植被覆盖度间经验模型，为三峡库区土壤侵蚀预报提供了基础支撑。中国地域辽阔，各区域气候、地形地貌、土地利用方式、植被类型与分布各不相同，造成CSLE模型中现有植被覆盖因子计算方法区域适用性和可操作性亟需提高。目前中国学者在USLE或RUSLE模型植被覆盖因子研究方面开展大量工作，且主要集中在西北黄土高原区和东北黑土区[11-25]；针对CSLE模型B因子，特别是南方红壤区B因子研究成果仍较为鲜见。此外，由于植被生长期划分、降雨侵蚀力空间分布等差异而导致植被覆盖因子研究结果地区间差异性较大[16]，应加强典型区域B值的系统性研究，为提高年度水土流失动态监测计算精度提供可靠的基础参数。
林下水土流失严重、崩岗侵蚀剧烈、水土流失呈零星斑点状分布等是南方红壤区水土流失特点[26-27]。近10年该区域侵蚀面积总体呈递减趋势，其中江西、湖北等6个省（区）呈下降趋势，广东、广西等4个省（区）呈增加趋势[28]。小流域综合治理（如农业生态工程、“大封禁、小治理”等）、坡地水土保持工程（等高植物篱等）、林下水土流失治理（林下补种草灌、针阔混交等）、侵蚀劣地和崩岗治理等是红壤区常见的水土流失防治模式[29]。与狗牙根、阔叶雀稗草等全园覆盖措施相比，百喜草全园覆盖措施蓄水保土效果最为显著，是常见水土保持植物措施[30]；特别是百喜草植物篱对于防治坡地水土流失效果十分显著，且对土壤理化性质改善具有积极作用[31]。也有研究表明，尽管植物措施可显著拦蓄径流、减少泥沙，但红壤丘陵区梯田果园应以浅根系豆科牧草为主，且不适当的牧草种植模式对果树产量有负面影响[32]。胡建民等[33]基于野外标准径流小区3种不同耕作方式5 a连续观测数据研究发现，横坡间作水土保持防治效果最优，减流减沙率在75%～80%；顺坡间作减流减沙率为60%～65%。黄俊等[34]基于野外天然降雨试验研究发现，自然撂荒地次降雨产流产沙量均高于灌草地和人工林地。总体而言，植物、工程措施的有机结合、高效配置是红壤坡地水土保持防治切实有效途径[35-38]。
本研究基于南方红壤区自然坡面野外系列降雨试验，分析次降雨土壤侵蚀模数与土壤侵蚀各影响因子间定量关系，构建植被覆盖因子与植被覆盖度间定量关系，以期为南方红壤区水土流失动态监测、土壤侵蚀研究等提供有益启示。
1 材料与方法
1.1 试验地点、材料及设备
试验在广东省五华县水土保持试验推广站进行（115°37'E，24°05'N），该站点属珠江流域韩江上游典型的红壤丘陵区，亚热带季风气候，年降水量为1 300～1 900 mm；该区域土壤侵蚀以面蚀、沟蚀和崩岗侵蚀等为主，地表植物主要有马尾松（Pinus masso
niana Lamb.）、桉树（Eucalyptus gran dis Hill）、柚树（Citrus maxima (Burm.) Merr.）、木荷（Schima kwangtungensis Hung T. Chang）、绢毛蔷薇（Rosa sericea Lindl.）、芒萁（Dicranopteris dichoto ma (Thunb. ) Bernh.）、蔗鸪草（Eriachne pallescensR. Br.）等。研究区土壤田间持水量及饱和含水率分别为26.7%～36.7%和42.5%～44.4%，有机质质量分数为0.78%～1.69%；砂粒（0.05～1 mm）、粗粉粒（0.01～<0.05 mm）、中细粉粒（0.001～<0.01 mm）、黏粒（<0.001 mm）质量分数分别为12.7%～14.7%、10.8%～22.9%、26.3～15.8%、28.6%～36.1%。
1.2 试验材料及设备
本试验使用的人工模拟降雨装置由西安理工大学生产的下喷式人工模拟降雨器，主要由降雨喷头及供水管道、恒压供水系统、控制系统、铝合金支架等构成（如图1所示）。该装置有效降雨面积为2 m×3 m，降雨高度3.0 m，降雨强度范围0.5～2.5 mm/min，雨滴直径大小为0.6～2.2 mm（滤纸色斑法测定），雨滴终速度为2.85～7.14 m/s，次降雨降雨动能为19.89～30.06 J/(m2·mm)。该装置降雨均匀度测试结果表明，20场不同降雨强度下降雨均匀度均值为81.5%～93.6%。试验过程中降雨器周围设有遮风帘，以减少自然风对降雨试验的不利影响。
图1 人工模拟降雨装置示意图Fig.1 Schematic diagram of artificial rainfall simulation device
1.3 观测指标及方法
受雨坡面为自然灌草地或清耕休闲地（裸地）坡面，在“降雨区”投影下方使用PVC板材形成封闭的人工径流小区（1.5 m×2.5 m），坡度使用罗盘实测获取，植被覆盖度采用照相法获取[39]。各场次降雨试验观测指标均参照《径流小区和小流域水土保持监测手册》[40]收集。各场次降雨径流深由次降雨径流体积除以径流小区投影面积计算得到。次降雨土壤侵蚀模数等于次降雨产沙量除以降雨面积与次降雨历时乘积值。次降雨量由电子雨量计自动记录（记录间隔1 min）。
试验开始前采用烘干法测定1次雨前土壤含水率，测定土层深度为0～15 cm。每个小区仅开展人工模拟降雨1次。经统计本试验各小区雨前表层土壤含水率在10.5%～28.4%之间波动，平均值为17.3%±5.1%。
本试验各场次降雨历时、降雨强度和次降雨量统计结果如表1所示，可以看出本试验各场次降雨总体呈现短历时、大雨强、小雨量的基本特征。按照秦伟等[41]在红壤区次降雨雨型研究结果，本试验人工模拟降雨属于“A雨型”降雨，也是造成土壤侵蚀的主要降雨类型，与试验地点常见的“短历时大雨强小雨量”降雨雨型基本一致。
表1 本试验各场次降雨参数统计结果Table 1 Statistical results of each individual rainfall event in this study统计量Statistics降雨历时Rainfall duration/min降雨强度Rainfall intensity/(mm·min-1)次降雨量Individual rainfall amount/mm 最大值Maximum6.010.5412.37 最小值Minimum34.283.1779.82 平均值Average value16.95±6.831.58±0.6234.44±18.12 中位数Median16.251.4530.71
1.4 植被覆盖因子计算方法
本文植被覆盖因子（B）值计算方法基于降雨试验实测数据，借助CSLE方程反推获取，具体如下：
设某一清耕休闲地（裸地）1次人工模拟降雨土壤流失量为
同一区域某一植被覆盖度为VC的灌草地坡面1次人工模拟降雨土壤流失量为
植被覆盖度为VC的灌草地坡面植被覆盖因子Bs由如下方程计算：
式中A为土壤侵蚀模数，t/(hm2·a)；R为降雨侵蚀力因子，MJ·mm/(hm2·h·a)；K为土壤可蚀性因子，t·hm2·h /(hm2·MJ·mm)；L为坡长因子，无量纲；S为坡度因子，无量纲；B为植被覆盖与生物措施因子，无量纲；E为工程措施因子，无量纲；T为耕作措施因子，无量纲；下标c和s分别表示清耕休闲地（裸地）、植被覆盖度为VC的灌草地坡面对应的变量值。A、R、K、L、S因子可通过实测或计算获取，E、T因子可根据实际情况取值获取；清耕休闲地（裸地）植被覆盖因子Bc=1。
各场次降雨侵蚀计算采用章文波等[42]提出PI10指标，即为次降雨量（Pr，mm）和最大10 min降雨强度（I10，mm/h）的乘积。土壤可蚀性采用基于侵蚀/生产力的土壤侵蚀预报模型（Erosion Productivity Impact Calculator，EPIC）计算获得[43]。坡长和坡度因子计算公式如下[5,44]：
式中G为小区坡度值（本研究中人工模拟降雨试验小区坡度在10°～15°之间，自然降雨试验小区坡度为20°和30°），Len为小区实际坡长值，m；t为坡长因子指数为(sinG/0.0896)/(3sinG0.8+0.56)。
此外，本研究中所有径流小区均无水土保持耕作和工程措施，因此工程措施因子和耕作措施因子均为1，即为E=1，T=1。
1.5 数据收集与处理
2016年7―8月和2018年5―6月共计完成110余场次人工模拟降雨试验，经过分析筛选出98组包含完整次降雨产流产沙过程曲线的有效观测数据（包含对照处理清耕休闲地9场次）。此外，本研究还使用了《全国水土流失动态监测与公告项目》广东省五华县乌陂河综合观测站2016年、2017年和2019年标准人工径流小区（20 m×5 m）天然降雨观测数据81场次（包含对照处理清耕休闲地10场次）、29场次（包含对照处理清耕休闲地5场次）和33场次（不含对照处理清耕休闲地观测数据），数据来源于《全国水土流失动态监测与公告项目》广东省五华县乌陂河综合观测站2016年、2017年和2019年年度整编数据成果，未公开发表。
本研究中的人工模拟降雨试验与标准径流小区天然降雨试验的土壤理化参数完全一致、降雨雨型基本相同、植被类型高度一致，因此人工模拟降雨试验数据与标准径流小区观测数据具有较强的可比性，可一同用于开展植被覆盖因子研究。
使用Microsoft Excel 2013?软件包进行数据处理，使用Origin Pro 8.5?软件包绘图。
1.6 模型模拟精度评价
采用拟合方程决定系数（R2）、纳什统计系数（Nash-Sutcliffe Efficiency Coefficient，NS）、均方根误差（Root Mean Square Error，RMSE）、平均绝对百分误差（Mean Absolute Percentage Error，MAPE）、平均绝对误差（Mean Absolute Error，MAE）对模型模拟进行定量评价，计算公式如下：
（6）
（9）
式中YSIMi为数据序列第i个模拟值；YOBSi为数据序列第i个观测值；YOBSmean为数据序列观测值的算术平均值；n为数据样本总数。
2 结果与分析
2.1 次降雨土壤侵蚀模数
图2为次降雨土壤侵蚀模数（Individual Rainfall Event Soil Erosion Modulus，ISEM）与土壤侵蚀各影响因子散点图。当VC≤0.5区间内，ISEM随VC增加而急剧降低，而后降低速率逐渐放缓趋于平稳。ISEM与雨前表层土壤含水率呈对数变化关系。ISEM随地表坡度（G）呈先增加后减小的变化趋势，其中临界坡度值约为10.75°。当G≤6°或G≥15°时，ISEM相对较小均不超过0.5 t/(hm2·h)。ISEM与次降雨量、次降雨侵蚀力和次降雨径流深间呈现极显著正相关线性关系（R2=0.63～0.81）。
图2 次降雨土壤侵蚀模数与各因子散点图Fig.2 Scatters between Individual Rainfall Event Soil Erosion Modulus (ISEM) and each factor
2.2 B因子计算模型构建
基于2016年和2018年98场次人工模拟降雨试验、2016年81场天然降雨试验、2017年29场次天然降雨试验观测数据，根据式（3）计算各场次降雨植被覆盖因子（B）值，进而采用线性函数、指数函数、自然对数函数、幂函数，以及一阶、二阶和三阶指数衰减模型7种函数关系对B和VC进行了回归拟合，获得不同植被覆盖度（VC）和因子B之间定量关系，如表2所示。各拟合方程均达到统计学显著水平（P<0.01）。其中“模型Ⅵ”决定系数最大为0.947；“模型Ⅵ”和“模型Ⅶ”纳什系数最大且均为0.876、均方根误差最小且均为0.021、平均绝对误差最小约为0.02；“模型Ⅱ”平均绝对百分误差最小为33.1%。总体而言，“模型Ⅵ”和“模型Ⅶ”为B和VC回归拟合的最优模型。但“模型Ⅶ”相对更为复杂，因此本研究选择“模型Ⅵ”作为定量刻画B与VC间定量关系的数学模型。基于模型VI的拟合曲线如图3所示，随VC增加B呈先迅速降低、后逐渐放缓的变化过程。
表2 植被覆盖因子与植被覆盖度不同拟合方程Table 2 Fitted equations between vegetation cover factor and vegetation cover模型序号Model No.函数关系Function relationship拟合方程Fitting equationsR2NSRMSEMAPE/%MAE Ⅰ线性函数B=-0.327VC+0.1940.6230.3950.038108.20.025 Ⅱ指数函数B=0.364e-5.308Vc0.8130.7700.02433.10.018 Ⅲ自然对数函数B=-0.116lnVC-0.0540.8570.8340.02364.00.019 Ⅳ幂函数B=0.010VC-1.5030.7110.5620.09043.10.030 Ⅴ一阶指数衰减模型B=0.043+0.783e-13.97Vc0.8550.7550.03288.80.024 Ⅵ二阶指数衰减模型B=-0.008+0.316e-4.35Vc+0.691e-39.87Vc0.9470.8760.02136.90.016 Ⅶ三阶指数衰减模型B=-0.008+0.156e-4.35Vc+0.160e-4.35Vc+0.691e-39.87Vc0.9460.8760.02136.90.016
对于特定区域而言，在土壤侵蚀其他影响因素不变的情况下，当植被覆盖度从0增加到最大值，植物措施的水土保持效应也逐渐增强，对应的B因子也将从最大值降低到最小值，而该区域土壤侵蚀量也从理论最大值降低到理论最小值，这是植被措施调控土壤侵蚀作用的具体表现。这一变化过程中客观存在某一特定点VC值（为临界植被覆盖度，设为VC_p），当VC>VC_p，B值衰减速率逐渐放缓并趋于稳定某一值，植被措施调控土壤侵蚀的作用也达到最大化；此时，在其他影响因素恒定情况下，VC增加基本不会引起土壤侵蚀量改变。
图3 植被覆盖因子与植被覆盖度散点图Fig.3 Scatter of vegetation coverage vs vegetation cover and management factor
基于B和VC间二阶指数衰减模型无法直接写出f(VC)=0对应VC计算表达值，难以直接获取VC_p值。基于模型Ⅲ（自然对数函数），令f(VC)=0可直接得到VC=e0.054/0.116=62.59%，故而可近似认为VC_p为62.59%。
对模型Ⅵ（二阶指数衰减函数）求一阶导数，得到：
显然，任意植被覆盖度VC计算得到的B'均表示曲线B=f(VC)上对应点的斜率值，即为B值随植被覆盖度增加而下降的速率。将VC_p=62.59%代入式（4）得到B'=-9‰。由于二阶指数衰减函数是单调递减函数，随着植被覆盖度增加，其递减速率逐渐降低。
因此，当VC>62.59%时，△VC引起B因子变化量△B的绝对值<9‰。在VC∈[0.626 1.00]范围内，9‰≤B≤13.1‰，且B均值为5.9‰±4.3‰；可认为在VC∈[0.626 1.00]范围内B因子并未跟随VC持续增加而显著降低，即植被措施调控土壤侵蚀能力并未显著增加。因此，为减少计算工作量可简化认为当VC>62.59%时，B≈0。则B因子计算模型可以简化为
（11）
2.3 B因子计算模型验证
本研究成果为点尺度B因子研究结论，由于植被措施水土保持效应具有典型的时空特征，B因子计算模型应从点、面2个尺度进行验证。因此，本研究首先采用经典方法使用点尺度野外观测数据对B因子计算式（11）进行评价，以确定式（11）计算结果的精度和合理性；在点尺度验证结果精度和合理性满足要求的基础上，从面尺度对式（11）的适用性和精度进行评价。点、面2个尺度对式（11）的验证结果可相互支撑佐证。
点尺度采用乌陂河综合观测站2019年度33场次自然降雨试验观测数据对式（11）进行验证；面尺度是以五华县为例，使用式（11）计算得到的B因子（记为Bf）与广东省五华县2018和2019年“区域水土流失动态监测成果”B因子（记为Bp）进行对比统计分析。其中Bp计算使用前3年的24个半月MODIS数据，具体计算方法见《区域水土流失动态监测技术规定（试行）》。为保证计算结果的可比性，本文同样使用前3年24个半月MODIS数据计算植被覆盖度值，然后基于式（11）计算Bf值。图4为B因子计算式（11）在点尺度验证结果。可以看出，使用式（11）计算得到各次降雨土壤流失量模拟值与实际观测值基本吻合，二者间线性拟合方程达到了极显著水平（P<0.01）；模拟值和观测值相对误差全部分布在±0.40以内，其中相对误差绝对值小于0.20的模拟值占总数的54.55%，小于0.30的占90%；33组观测数据模拟值与观测值均方根误差为0.012 t/hm2，相对误差绝对值的平均值为0.186。总体而言，B因子计算式（11）在点尺度的计算精度良好，其中计算误差较大的几组样本主要是其次降雨量和次降雨侵蚀力偏小所致。
图4 植被覆盖因子（B）计算模型点尺度验证结果Fig.4 Verification results of vegetation cover and management factor calculation model at point scale
如前所述，B因子计算模型点尺度验证结果较为理想，土壤流失模拟值与观测值十分接近，计算精度良好。因此，从面尺度对B因子计算模型进一步验证。图5和表3分别为五华县2018和2019年Bf与Bp相对误差绝对值栅格图和统计结果。总体而言Bf与Bp十分接近，2018和2019年相对误差绝对值小于0.1的分别超过80%和70%。这表明式（11）实际可行的，其计算精度也较为理想，可使用植被覆盖度作为单一变量计算B值。需要注意的是2018和2019年Bf与Bp相对误差绝对值超过50%的也分别占到样本总量的17%和29%左右，这些误差较大的像元主要位于植被和非植被分界地带，而本研究仅将植被区域纳入计算，加上植被措施的空间尺度效应等综合导致位于这些边界地带的像元计算结果误差偏大。
图5 植被覆盖因子相对误差绝对值栅格图Fig. 5 Raster map of absolute relative error of vegetation cover and management factor
表3 植被覆盖因子相对误差绝对值统计结果Table 3 Statistical results of absolute relative error of vegetation cover and management factor相对误差绝对值范围Range of absolute relative error占总像元个数比例Proportion to total number of pixels/% 2018年2019年 0～0.181.4270.87 >0.1～0.20.130.03 >0.2～0.30.030.06 >0.3～0.40.030.13 >0.4～0.51.400.30 >0.5～0.60.090.99 >0.6～0.70.293.45 >0.7～0.80.638.92 >0.8～0.95.2713.68 >0.9～1.08.721.38 >1.01.990.19
本文基于点尺度降雨试验、基于中国土壤流失方程获取不同VC下B值，采用回归分析法构建了B值与VC间二阶指数衰减模型，并从点尺度和面尺度2个层次验证了该模型计算精度。尽管该模型在点尺度验证上表现出良好的结果，但面尺度验证中仍存在20%～30%样本值相对误差仍超过50%。尽管有诸多研究表明人工模拟降雨试验一定程度可代替自然降雨开展水土保持研究，但诸多雨滴速度、降雨动能、雨滴大小等指标仍不能完成与自然降雨一致，或者上述指标未能与自然降雨实现一致的比例关系，这是本研究结论中相对误差原因之一，需要后续更为丰富的野外观测数据对模型进行完善补充。本研究面尺度验证中仅考虑了植被区域B值计算，但由于遥感影像分辨率的客观问题，影像中植被/非植被交接区域像元的计算结果仍存在较大误差，也是误差的主要来源。此外，本研究成果直接“由点到面”的验证应用仍存在空间尺度转换的问题，需要不同空间尺度下更为丰富的观测数据继续开展深入研究。
3 结 论
本文分析了次降雨土壤侵蚀模数与影响水土流失各关键因子定量关系，可为红壤区水土流失防治、土壤侵蚀模型构建提供参考。
1）次降雨土壤土壤侵蚀模数随植被覆盖度增加而逐渐降低；次降雨土壤侵蚀模数随雨前土壤含水率增加而逐渐增大，二者呈正相关对数关系；次降雨土壤侵蚀模数与坡度间呈二次抛物线关系，次降雨土壤侵蚀模数随坡度增加呈先增加后减小变化趋势；次降雨土壤侵蚀模数随次降雨量、次降雨侵蚀力和次降雨径流深呈显著正相关线性变化关系。
2）基于系列次降雨实测资料，采用中国土壤流失方程反推法获取了不同植被覆盖度下植被覆盖因子；采用线性函数、指数函数、对数函数、幂函数和指数衰减函数对植被覆盖因子和植被覆盖度定量关系进行了拟合，基于5个统计量对各拟合方程进行了对比分析，优选确定二阶指数衰减函数为刻画植被覆盖因子和植被覆盖度间定量关系的数学模型，该拟合方程纳什系数、均方根误差和决定系数分别为0.876、0.021和0.947。该模型点尺度验证结果良好，计算模拟值与实际观测值十分吻合，二者线性拟合方程达到了极显著水平（P<0.01）；面尺度验证结果较为理想，70%～80%模型计算数据相对误差绝对值不超过0.1。
限于试验数据有限，研究中用于建立和验证植被覆盖因子计算模型的样本数据仍然相对较少，后续工作还需更为丰富的观测资料、继续完善相关研究成果。
致谢：感谢广东省梅州市五华县水土保持试验推广站对本研究的支持。
0 引 言
区域土壤侵蚀状况直接对人类健康和经济社会发展产生重要影响[1]，定期开展水土流失动态监测可为水土流失状况评价、水土保持治理工程空间布局等提供重要支撑[2]。土壤侵蚀模型是开展水土流失动态监测的重要技术手段[3-4]。迄今为止，美国通用土壤流失方程（Universal Soil Loss Equation，USLE）及其改进模型（Revised Universal Soil Loss Equation，RUSLE）仍在土壤侵蚀领域广泛应用[5-6]；该模型综合土壤侵蚀多个影响因子及其相互作用关系，同时各因子具有一定物理意义且因子获取方法及模型计算过程较为简单[7-8]。Liu等[9-10]利用黄土高原丘陵沟壑区径流小区实测资料构建了中国土壤流失方程（China Soil Loss Equation，CSLE），该模型依据中国土壤侵蚀与水土保持工作特点，充分考虑了生物措施、工程措施和耕作措施对土壤侵蚀及水土流失过程的影响，与USLE相比更能反映中国土壤侵蚀的实际情况[11]。
2018年水利部印发《区域水土流失动态监测技术规定（试行）》，中国水土流失动态监测预报工作从定量与定性相结合的“三因子”综合判别法发展到以定量为主的“CSLE模型法”[12]。2018年结果表明中国仍有水土流失面积273.69万km2，约占国土总面积的28.54%[13]。CSLE模型广泛应用对中国水土流失动态监测应用工作和土壤侵蚀基础研究具有十分显著的促进作用，但也存在诸如模型计算精度偏低、模型因子区域适用性较差、地形及土壤等基础数据现势性不足等亟需解决的问题[14-15]。作为影响土壤侵蚀最重要因素之一：植被覆盖因子（B），可直接反映植被覆盖变化状况或植被措施调控土壤侵蚀的能力，也是CSLE模型中的关键因子[16-19]。植被覆盖因子指一定条件下有植被覆盖或实施田间管理的土地的土壤流失总量与同等条件下实施清耕休闲地的土壤流失总量的比值，是介于0和1之间的一个无量纲数[20]。植被覆盖度大小是影响植被覆盖因子的关键因素；此外，植被覆盖因子大小还受到植被不同生长期内侵蚀性降雨多寡的影响，对于植被覆盖因子年值的研究还应充分考虑降雨侵蚀力年内季节分布的影响。标准小区法、次因子法、模型反推法、与植被覆盖度关系式法、直接赋值或手册查询法等是较为常用的植被覆盖因子获取方法。前3种方法依赖野外实测数据，逻辑性强、计算精度高；关系式法虽难以直接反映植被减蚀作用机理，但其经验关系具有相对较强的通用性和适用性；直接赋值或手册查询法较为依赖大量丰富的外业数据观测成果，通用性相对较弱。目前中国年度水土流失动态监测项目采用《区域水土流失动态监测技术规定（试行）》推荐的B因子计算方法，利用监测年前3年内每年24个半月MODIS数据计算得到B因子值。植被覆盖因子与土地利用方式关系密切，合理的土地利用方式可显著降低植被覆盖因子值而减少土壤侵蚀量[21-22]。张雪花等[23]基于人工模拟降雨试验确定了东北黑土农业区不同土地利用植被因子的取值范围，发现该值普遍小于USLE方程推荐的最优值。蔡崇法等[24]基于野外实测资料构建了三峡库区植被因子与植被覆盖度间经验模型，为三峡库区土壤侵蚀预报提供了基础支撑。中国地域辽阔，各区域气候、地形地貌、土地利用方式、植被类型与分布各不相同，造成CSLE模型中现有植被覆盖因子计算方法区域适用性和可操作性亟需提高。目前中国学者在USLE或RUSLE模型植被覆盖因子研究方面开展大量工作，且主要集中在西北黄土高原区和东北黑土区[11-25]；针对CSLE模型B因子，特别是南方红壤区B因子研究成果仍较为鲜见。此外，由于植被生长期划分、降雨侵蚀力空间分布等差异而导致植被覆盖因子研究结果地区间差异性较大[16]，应加强典型区域B值的系统性研究，为提高年度水土流失动态监测计算精度提供可靠的基础参数。
林下水土流失严重、崩岗侵蚀剧烈、水土流失呈零星斑点状分布等是南方红壤区水土流失特点[26-27]。近10年该区域侵蚀面积总体呈递减趋势，其中江西、湖北等6个省（区）呈下降趋势，广东、广西等4个省（区）呈增加趋势[28]。小流域综合治理（如农业生态工程、“大封禁、小治理”等）、坡地水土保持工程（等高植物篱等）、林下水土流失治理（林下补种草灌、针阔混交等）、侵蚀劣地和崩岗治理等是红壤区常见的水土流失防治模式[29]。与狗牙根、阔叶雀稗草等全园覆盖措施相比，百喜草全园覆盖措施蓄水保土效果最为显著，是常见水土保持植物措施[30]；特别是百喜草植物篱对于防治坡地水土流失效果十分显著，且对土壤理化性质改善具有积极作用[31]。也有研究表明，尽管植物措施可显著拦蓄径流、减少泥沙，但红壤丘陵区梯田果园应以浅根系豆科牧草为主，且不适当的牧草种植模式对果树产量有负面影响[32]。胡建民等[33]基于野外标准径流小区3种不同耕作方式5 a连续观测数据研究发现，横坡间作水土保持防治效果最优，减流减沙率在75%～80%；顺坡间作减流减沙率为60%～65%。黄俊等[34]基于野外天然降雨试验研究发现，自然撂荒地次降雨产流产沙量均高于灌草地和人工林地。总体而言，植物、工程措施的有机结合、高效配置是红壤坡地水土保持防治切实有效途径[35-38]。
本研究基于南方红壤区自然坡面野外系列降雨试验，分析次降雨土壤侵蚀模数与土壤侵蚀各影响因子间定量关系，构建植被覆盖因子与植被覆盖度间定量关系，以期为南方红壤区水土流失动态监测、土壤侵蚀研究等提供有益启示。
1 材料与方法
1.1 试验地点、材料及设备
试验在广东省五华县水土保持试验推广站进行（115°37'E，24°05'N），该站点属珠江流域韩江上游典型的红壤丘陵区，亚热带季风气候，年降水量为1 300～1 900 mm；该区域土壤侵蚀以面蚀、沟蚀和崩岗侵蚀等为主，地表植物主要有马尾松（Pinus massonianaLamb.）、桉树（Eucalyptus gran disHill）、柚树（Citrus maxima(Burm.) Merr.）、木荷（Schima kwangtungensis HungT. Chang）、绢毛蔷薇（Rosa sericeaLindl.）、芒萁（Dicranopteris dichoto ma(Thunb. ) Bernh.）、蔗鸪草（Eriachne pallescensR. Br.）等。研究区土壤田间持水量及饱和含水率分别为26.7%～36.7%和42.5%～44.4%，有机质质量分数为0.78%～1.69%；砂粒（0.05～1 mm）、粗粉粒（0.01～<0.05 mm）、中细粉粒（0.001～<0.01 mm）、黏粒（<0.001 mm）质量分数分别为12.7%～14.7%、10.8%～22.9%、26.3～15.8%、28.6%～36.1%。
1.2 试验材料及设备
本试验使用的人工模拟降雨装置由西安理工大学生产的下喷式人工模拟降雨器，主要由降雨喷头及供水管道、恒压供水系统、控制系统、铝合金支架等构成（如图1所示）。该装置有效降雨面积为2 m×3 m，降雨高度3.0 m，降雨强度范围0.5～2.5 mm/min，雨滴直径大小为0.6～2.2 mm（滤纸色斑法测定），雨滴终速度为2.85～7.14 m/s，次降雨降雨动能为19.89～30.06 J/(m2·mm)。该装置降雨均匀度测试结果表明，20场不同降雨强度下降雨均匀度均值为81.5%～93.6%。试验过程中降雨器周围设有遮风帘，以减少自然风对降雨试验的不利影响。
图1 人工模拟降雨装置示意图Fig.1 Schematic diagram of artificial rainfall simulation device
1.3 观测指标及方法
受雨坡面为自然灌草地或清耕休闲地（裸地）坡面，在“降雨区”投影下方使用PVC板材形成封闭的人工径流小区（1.5 m×2.5 m），坡度使用罗盘实测获取，植被覆盖度采用照相法获取[39]。各场次降雨试验观测指标均参照《径流小区和小流域水土保持监测手册》[40]收集。各场次降雨径流深由次降雨径流体积除以径流小区投影面积计算得到。次降雨土壤侵蚀模数等于次降雨产沙量除以降雨面积与次降雨历时乘积值。次降雨量由电子雨量计自动记录（记录间隔1 min）。
试验开始前采用烘干法测定1次雨前土壤含水率，测定土层深度为0～15 cm。每个小区仅开展人工模拟降雨1次。经统计本试验各小区雨前表层土壤含水率在10.5%～28.4%之间波动，平均值为17.3%±5.1%。
本试验各场次降雨历时、降雨强度和次降雨量统计结果如表1所示，可以看出本试验各场次降雨总体呈现短历时、大雨强、小雨量的基本特征。按照秦伟等[41]在红壤区次降雨雨型研究结果，本试验人工模拟降雨属于“A雨型”降雨，也是造成土壤侵蚀的主要降雨类型，与试验地点常见的“短历时大雨强小雨量”降雨雨型基本一致。
表1 本试验各场次降雨参数统计结果Table 1 Statistical results of each individual rainfall event in this study统计量Statistics降雨历时Rainfall duration/min降雨强度Rainfall intensity/(mm·min-1)次降雨量Individual rainfall amount/mm 最大值Maximum6.010.5412.37 最小值Minimum34.283.1779.82 平均值Average value16.95±6.831.58±0.6234.44±18.12 中位数Median16.251.4530.71
1.4 植被覆盖因子计算方法
本文植被覆盖因子（B）值计算方法基于降雨试验实测数据，借助CSLE方程反推获取，具体如下：
设某一清耕休闲地（裸地）1次人工模拟降雨土壤流失量为
同一区域某一植被覆盖度为VC的灌草地坡面1次人工模拟降雨土壤流失量为
植被覆盖度为VC的灌草地坡面植被覆盖因子Bs由如下方程计算：
式中A为土壤侵蚀模数，t/(hm2·a)；R为降雨侵蚀力因子，MJ·mm/(hm2·h·a)；K为土壤可蚀性因子，t·hm2·h /(hm2·MJ·mm)；L为坡长因子，无量纲；S为坡度因子，无量纲；B为植被覆盖与生物措施因子，无量纲；E为工程措施因子，无量纲；T为耕作措施因子，无量纲；下标c和s分别表示清耕休闲地（裸地）、植被覆盖度为VC的灌草地坡面对应的变量值。A、R、K、L、S因子可通过实测或计算获取，E、T因子可根据实际情况取值获取；清耕休闲地（裸地）植被覆盖因子Bc=1。
各场次降雨侵蚀计算采用章文波等[42]提出PI10指标，即为次降雨量（Pr，mm）和最大10 min降雨强度（I10，mm/h）的乘积。土壤可蚀性采用基于侵蚀/生产力的土壤侵蚀预报模型（Erosion Productivity Impact Calculator，EPIC）计算获得[43]。坡长和坡度因子计算公式如下[5,44]：
式中G为小区坡度值（本研究中人工模拟降雨试验小区坡度在10°～15°之间，自然降雨试验小区坡度为20°和30°），Len为小区实际坡长值，m；t为坡长因子指数为(sinG/0.0896)/(3sinG0.8+0.56)。
此外，本研究中所有径流小区均无水土保持耕作和工程措施，因此工程措施因子和耕作措施因子均为1，即为E=1，T=1。
1.5 数据收集与处理
2016年7―8月和2018年5―6月共计完成110余场次人工模拟降雨试验，经过分析筛选出98组包含完整次降雨产流产沙过程曲线的有效观测数据（包含对照处理清耕休闲地9场次）。此外，本研究还使用了《全国水土流失动态监测与公告项目》广东省五华县乌陂河综合观测站2016年、2017年和2019年标准人工径流小区（20 m×5 m）天然降雨观测数据81场次（包含对照处理清耕休闲地10场次）、29场次（包含对照处理清耕休闲地5场次）和33场次（不含对照处理清耕休闲地观测数据），数据来源于《全国水土流失动态监测与公告项目》广东省五华县乌陂河综合观测站2016年、2017年和2019年年度整编数据成果，未公开发表。
本研究中的人工模拟降雨试验与标准径流小区天然降雨试验的土壤理化参数完全一致、降雨雨型基本相同、植被类型高度一致，因此人工模拟降雨试验数据与标准径流小区观测数据具有较强的可比性，可一同用于开展植被覆盖因子研究。
使用Microsoft Excel 2013?软件包进行数据处理，使用Origin Pro 8.5?软件包绘图。
1.6 模型模拟精度评价
采用拟合方程决定系数（R2）、纳什统计系数（Nash-Sutcliffe Efficiency Coefficient，NS）、均方根误差（Root Mean Square Error，RMSE）、平均绝对百分误差（Mean Absolute Percentage Error，MAPE）、平均绝对误差（Mean Absolute Error，MAE）对模型模拟进行定量评价，计算公式如下：
（6）
（9）
式中YSIMi为数据序列第i个模拟值；YOBSi为数据序列第i个观测值；YOBSmean为数据序列观测值的算术平均值；n为数据样本总数。
2 结果与分析
2.1 次降雨土壤侵蚀模数
图2为次降雨土壤侵蚀模数（Individual Rainfall Event Soil Erosion Modulus，ISEM）与土壤侵蚀各影响因子散点图。当VC≤0.5区间内，ISEM随VC增加而急剧降低，而后降低速率逐渐放缓趋于平稳。ISEM与雨前表层土壤含水率呈对数变化关系。ISEM随地表坡度（G）呈先增加后减小的变化趋势，其中临界坡度值约为10.75°。当G≤6°或G≥15°时，ISEM相对较小均不超过0.5 t/(hm2·h)。ISEM与次降雨量、次降雨侵蚀力和次降雨径流深间呈现极显著正相关线性关系（R2=0.63～0.81）。
图2 次降雨土壤侵蚀模数与各因子散点图Fig.2 Scatters between Individual Rainfall Event Soil Erosion Modulus (ISEM) and each factor
2.2B因子计算模型构建
基于2016年和2018年98场次人工模拟降雨试验、2016年81场天然降雨试验、2017年29场次天然降雨试验观测数据，根据式（3）计算各场次降雨植被覆盖因子（B）值，进而采用线性函数、指数函数、自然对数函数、幂函数，以及一阶、二阶和三阶指数衰减模型7种函数关系对B和VC进行了回归拟合，获得不同植被覆盖度（VC）和因子B之间定量关系，如表2所示。各拟合方程均达到统计学显著水平（P<0.01）。其中“模型Ⅵ”决定系数最大为0.947；“模型Ⅵ”和“模型Ⅶ”纳什系数最大且均为0.876、均方根误差最小且均为0.021、平均绝对误差最小约为0.02；“模型Ⅱ”平均绝对百分误差最小为33.1%。总体而言，“模型Ⅵ”和“模型Ⅶ”为B和VC回归拟合的最优模型。但“模型Ⅶ”相对更为复杂，因此本研究选择“模型Ⅵ”作为定量刻画B与VC间定量关系的数学模型。基于模型VI的拟合曲线如图3所示，随VC增加B呈先迅速降低、后逐渐放缓的变化过程。
表2 植被覆盖因子与植被覆盖度不同拟合方程Table 2 Fitted equations between vegetation cover factor and vegetation cover模型序号Model No.函数关系Function relationship拟合方程Fitting equationsR2NSRMSEMAPE/%MAE Ⅰ线性函数B=-0.327VC+0.1940.6230.3950.038108.20.025 Ⅱ指数函数B=0.364e-5.308Vc0.8130.7700.02433.10.018 Ⅲ自然对数函数B=-0.116lnVC-0.0540.8570.8340.02364.00.019 Ⅳ幂函数B=0.010VC-1.5030.7110.5620.09043.10.030 Ⅴ一阶指数衰减模型B=0.043+0.783e-13.97Vc0.8550.7550.03288.80.024 Ⅵ二阶指数衰减模型B=-0.008+0.316e-4.35Vc+0.691e-39.87Vc0.9470.8760.02136.90.016 Ⅶ三阶指数衰减模型B=-0.008+0.156e-4.35Vc+0.160e-4.35Vc+0.691e-39.87Vc0.9460.8760.02136.90.016
对于特定区域而言，在土壤侵蚀其他影响因素不变的情况下，当植被覆盖度从0增加到最大值，植物措施的水土保持效应也逐渐增强，对应的B因子也将从最大值降低到最小值，而该区域土壤侵蚀量也从理论最大值降低到理论最小值，这是植被措施调控土壤侵蚀作用的具体表现。这一变化过程中客观存在某一特定点VC值（为临界植被覆盖度，设为VC_p），当VC>VC_p，B值衰减速率逐渐放缓并趋于稳定某一值，植被措施调控土壤侵蚀的作用也达到最大化；此时，在其他影响因素恒定情况下，VC增加基本不会引起土壤侵蚀量改变。
图3 植被覆盖因子与植被覆盖度散点图Fig.3 Scatter of vegetation coverage vs vegetation cover and management factor
基于B和VC间二阶指数衰减模型无法直接写出f(VC)=0对应VC计算表达值，难以直接获取VC_p值。基于模型Ⅲ（自然对数函数），令f(VC)=0可直接得到VC=e0.054/0.116=62.59%，故而可近似认为VC_p为62.59%。
对模型Ⅵ（二阶指数衰减函数）求一阶导数，得到：
显然，任意植被覆盖度VC计算得到的B'均表示曲线B=f(VC)上对应点的斜率值，即为B值随植被覆盖度增加而下降的速率。将VC_p=62.59%代入式（4）得到B'=-9‰。由于二阶指数衰减函数是单调递减函数，随着植被覆盖度增加，其递减速率逐渐降低。
因此，当VC>62.59%时，△VC引起B因子变化量△B的绝对值<9‰。在VC∈[0.626 1.00]范围内，9‰≤B≤13.1‰，且B均值为5.9‰±4.3‰；可认为在VC∈[0.626 1.00]范围内B因子并未跟随VC持续增加而显著降低，即植被措施调控土壤侵蚀能力并未显著增加。因此，为减少计算工作量可简化认为当VC>62.59%时，B≈0。则B因子计算模型可以简化为
（11）
2.3B因子计算模型验证
本研究成果为点尺度B因子研究结论，由于植被措施水土保持效应具有典型的时空特征，B因子计算模型应从点、面2个尺度进行验证。因此，本研究首先采用经典方法使用点尺度野外观测数据对B因子计算式（11）进行评价，以确定式（11）计算结果的精度和合理性；在点尺度验证结果精度和合理性满足要求的基础上，从面尺度对式（11）的适用性和精度进行评价。点、面2个尺度对式（11）的验证结果可相互支撑佐证。
点尺度采用乌陂河综合观测站2019年度33场次自然降雨试验观测数据对式（11）进行验证；面尺度是以五华县为例，使用式（11）计算得到的B因子（记为Bf）与广东省五华县2018和2019年“区域水土流失动态监测成果”B因子（记为Bp）进行对比统计分析。其中Bp计算使用前3年的24个半月MODIS数据，具体计算方法见《区域水土流失动态监测技术规定（试行）》。为保证计算结果的可比性，本文同样使用前3年24个半月MODIS数据计算植被覆盖度值，然后基于式（11）计算Bf值。图4为B因子计算式（11）在点尺度验证结果。可以看出，使用式（11）计算得到各次降雨土壤流失量模拟值与实际观测值基本吻合，二者间线性拟合方程达到了极显著水平（P<0.01）；模拟值和观测值相对误差全部分布在±0.40以内，其中相对误差绝对值小于0.20的模拟值占总数的54.55%，小于0.30的占90%；33组观测数据模拟值与观测值均方根误差为0.012 t/hm2，相对误差绝对值的平均值为0.186。总体而言，B因子计算式（11）在点尺度的计算精度良好，其中计算误差较大的几组样本主要是其次降雨量和次降雨侵蚀力偏小所致。
图4 植被覆盖因子（B）计算模型点尺度验证结果Fig.4 Verification results of vegetation cover and management factor calculation model at point scale
如前所述，B因子计算模型点尺度验证结果较为理想，土壤流失模拟值与观测值十分接近，计算精度良好。因此，从面尺度对B因子计算模型进一步验证。图5和表3分别为五华县2018和2019年Bf与Bp相对误差绝对值栅格图和统计结果。总体而言Bf与Bp十分接近，2018和2019年相对误差绝对值小于0.1的分别超过80%和70%。这表明式（11）实际可行的，其计算精度也较为理想，可使用植被覆盖度作为单一变量计算B值。需要注意的是2018和2019年Bf与Bp相对误差绝对值超过50%的也分别占到样本总量的17%和29%左右，这些误差较大的像元主要位于植被和非植被分界地带，而本研究仅将植被区域纳入计算，加上植被措施的空间尺度效应等综合导致位于这些边界地带的像元计算结果误差偏大。
图5 植被覆盖因子相对误差绝对值栅格图Fig. 5 Raster map of absolute relative error of vegetation cover and management factor
表3 植被覆盖因子相对误差绝对值统计结果Table 3 Statistical results of absolute relative error of vegetation cover and management factor相对误差绝对值范围Range of absolute relative error占总像元个数比例Proportion to total number of pixels/% 2018年2019年 0～0.181.4270.87 >0.1～0.20.130.03 >0.2～0.30.030.06 >0.3～0.40.030.13 >0.4～0.51.400.30 >0.5～0.60.090.99 >0.6～0.70.293.45 >0.7～0.80.638.92 >0.8～0.95.2713.68 >0.9～1.08.721.38 >1.01.990.19
本文基于点尺度降雨试验、基于中国土壤流失方程获取不同VC下B值，采用回归分析法构建了B值与VC间二阶指数衰减模型，并从点尺度和面尺度2个层次验证了该模型计算精度。尽管该模型在点尺度验证上表现出良好的结果，但面尺度验证中仍存在20%～30%样本值相对误差仍超过50%。尽管有诸多研究表明人工模拟降雨试验一定程度可代替自然降雨开展水土保持研究，但诸多雨滴速度、降雨动能、雨滴大小等指标仍不能完成与自然降雨一致，或者上述指标未能与自然降雨实现一致的比例关系，这是本研究结论中相对误差原因之一，需要后续更为丰富的野外观测数据对模型进行完善补充。本研究面尺度验证中仅考虑了植被区域B值计算，但由于遥感影像分辨率的客观问题，影像中植被/非植被交接区域像元的计算结果仍存在较大误差，也是误差的主要来源。此外，本研究成果直接“由点到面”的验证应用仍存在空间尺度转换的问题，需要不同空间尺度下更为丰富的观测数据继续开展深入研究。
3 结 论
本文分析了次降雨土壤侵蚀模数与影响水土流失各关键因子定量关系，可为红壤区水土流失防治、土壤侵蚀模型构建提供参考。
1）次降雨土壤土壤侵蚀模数随植被覆盖度增加而逐渐降低；次降雨土壤侵蚀模数随雨前土壤含水率增加而逐渐增大，二者呈正相关对数关系；次降雨土壤侵蚀模数与坡度间呈二次抛物线关系，次降雨土壤侵蚀模数随坡度增加呈先增加后减小变化趋势；次降雨土壤侵蚀模数随次降雨量、次降雨侵蚀力和次降雨径流深呈显著正相关线性变化关系。
2）基于系列次降雨实测资料，采用中国土壤流失方程反推法获取了不同植被覆盖度下植被覆盖因子；采用线性函数、指数函数、对数函数、幂函数和指数衰减函数对植被覆盖因子和植被覆盖度定量关系进行了拟合，基于5个统计量对各拟合方程进行了对比分析，优选确定二阶指数衰减函数为刻画植被覆盖因子和植被覆盖度间定量关系的数学模型，该拟合方程纳什系数、均方根误差和决定系数分别为0.876、0.021和0.947。该模型点尺度验证结果良好，计算模拟值与实际观测值十分吻合，二者线性拟合方程达到了极显著水平（P<0.01）；面尺度验证结果较为理想，70%～80%模型计算数据相对误差绝对值不超过0.1。
限于试验数据有限，研究中用于建立和验证植被覆盖因子计算模型的样本数据仍然相对较少，后续工作还需更为丰富的观测资料、继续完善相关研究成果。
致谢：感谢广东省梅州市五华县水土保持试验推广站对本研究的支持。
文章来源：南方园艺 网址: http://nfyy.400nongye.com/lunwen/itemid-12542.shtml
上一篇： 南方某镇污水处理厂三期工程设计
下一篇： 农业经济论文_黑龙江省主要农产品价格监测分析——以2021年第1季度为例
