河南省豫西黄土丘陵区2015—2020年水土保持功能时空变化

杜 军， 赵胜朝， 徐建昭， 郭 雷

（1.河南省科学院地理研究所，郑州 450052； 2.河南省水土保持监测总站，郑州 450003）

习近平总书记在黄河流域生态保护和高质量发展座谈会上发表重要讲话中指出，黄河流域是我国重要的生态屏障和经济地带，是打赢脱贫攻坚战的重要区域，在我国经济社会发展和生态安全方面具有十分重要的地位，保护黄河是事关中华民族伟大复兴的千秋大计. 豫西黄土丘陵区是河南省沿黄区域生态环境比较脆弱的区域之一，也是河南省水土流失最严重的地区之一［1-2］，区域水土保持功能的改善与否，对于黄河流域生态保护与高质量发展具有重要的影响［3-5］. 经过近些年国土绿化、封山育林、坡改梯、淤地坝等水保措施［6-8］，区域生态环境有了较大改善，但对于区域水土保持功能空间变化以及近五年的时间变化缺乏相关研究，因此本研究在前人研究的基础上，结合豫西黄土丘陵的实际情况，基于2015和2020年的遥感影像获取研究区土地利用动态变化数据，采用RUSLE 模型对豫西黄土丘陵地区2015—2020 年间区域水土保持功能的变化特征进行评估，为豫西黄土丘陵区土壤侵蚀的治理提供一定的参考依据.

1 研究区概况与数据来源

1.1 研究范围

豫西黄土丘陵区［9-10］（110°21′～113°30′、33°32′～35°17′）位于河南省西北部，行政区划上包括郑州市、洛阳市、焦作市、三门峡市和济源市的19 个县（市），总面积2.71 万km2，约占河南省总面积的16.23%（图1）. 属于黄土高原向黄淮海平原过渡的交错地带，水土流失较为严重，是河南省生态环境比较脆弱的区域之一. 区域内地形起伏较大、地貌复杂、土地利用类型多样，最低海拔85 m，最高海拔2392 m，地貌类型主要包括山地（中山和低山）、丘陵和平原等，土地利用类型主要以耕地、林地为主. 该区属于温带半湿润半干旱气候，多年平均气温8～15 ℃，多年平均降水约为574.0～634.4 mm.

图1 豫西黄土丘陵区概况Fig.1 General situation of Loess Hilly Area in western Henan

1.2 数据来源

本研究中利用2015年、2020年的遥感影像数据，对影像进行融合、校正等预处理，利用监督分类法并结合人工目视判读对2期利用进行解译［11］，将研究区分为林地、草地、耕地、水域、建设用地和未利用地6种类型，解译结果随机采样并进行野外核查，总体精度大于95%.

2 研究方法

利用RUSLE模型［12-14］的土壤侵蚀评价方法，土壤侵蚀模型RUSLE考虑了降雨、土壤可蚀性、作物管理、坡度坡长和水土保持措施五大因子，方程式如下：

式中：A为年平均土壤流失量，t/（km2·a）；R为降雨和径流侵蚀因子，MJ·mm/（km2·h·a）；K为土壤可蚀性因子，此因子是指在标准状态下（小区长22.1 m，坡度9%，休闲地且顺坡耕作）每年单位降雨侵蚀指标所产生的侵蚀数量，单位为t·h/（MJ·mm）；L，S为地形因子，其中L为坡长因子，S为坡度因子；C为作物管理因子；P为治理措施因子.

3 土壤侵蚀评价因子提取方法

3.1 降雨侵蚀力R因子

降水时雨滴击溅、分离土壤颗粒以及径流冲刷和转运是引起土壤侵蚀的主要动力因素. 美国在使用通用土壤流失方程预报土壤流失量时，一直用EI30作为经典算法来计算降雨侵蚀力指标，这种经典算法最初是由Wischmeier 为确定降雨侵蚀力指标［15］，对降雨量、降雨动能、最大15 min 雨强、最大30 min 雨强等变量与土壤流失量之间进行决定系数分析发现的，EI30 可以解释土壤流失量变化的70.2%～89.2%. 经典算法一般都需要降雨过程资料，必须用断点降雨资料计算降雨动能，工作量大且资料的获取难度大，因此采用经典算法计算降雨侵蚀力十分不易. 后来有学者逐步使用了简易算法，即采用年、月或日降雨资料进行降雨侵蚀力的计算. 简易算法资料获取更为容易，同时资料的整理及计算相对简单. 谢云等［16］对我国降雨情况进行分析得出日降雨一般能够替代次降雨计算降雨侵蚀力. 章文波和付金生［17］通过全国71 个气象站点的降雨数据分析得出，在利用年、月、日降雨资料等5 种代表性雨量资料估算降雨侵蚀力的结果中，日雨量计算多年平均侵蚀力的精度最高，并建立了利用日雨量估算降雨侵蚀力的简易算法模型，提出了模型的参数估计方法. 模型在全国的决定系数平均为0.718，估算多年平均侵蚀力的相对误差平均为4.2%.公式如下：

式中：Mi为第i个半月时段的侵蚀力值；α，β为模型参数；K为该半月时段内的天数；Dj为半月时段内第j天的日雨量，要求日雨量≥12 mm，否则以0计算，12 mm与中国侵蚀性降雨标准对应. 参数α和β反映了区域降雨特征，如下式：

式中：Pd12为日雨量≥12 mm的日平均雨量；Py12为日雨量≥12 mm的年平均雨量.

由于不是每次降雨都会产生侵蚀，少于一定数量的降雨经常不会引起土壤侵蚀，因此在统计降雨资料时，存在一个侵蚀性降雨标准. 美国侵蚀性降雨标准定为次降雨量大于6.5 mm，或一次降雨时间长于12.7 min 或15 min. 中国在黄土高原的多年研究中，建立了相应的侵蚀性降雨雨量标准为日降雨量大于12 mm，也是本文采用的标准.

本研究利用章文波等的日降雨量估算方法计算降雨侵蚀力R因子. 降雨数据来自河南省2015—2020年的逐日降雨资料，并整理其中的侵蚀性降雨，按照以上公式进行估算.

3.2 土壤可蚀性K因子

土壤可蚀性是土壤侵蚀预报中的一个重要方面，它是评价土壤是否易受侵蚀营力破坏的性能，也是土壤对侵蚀营力分离和搬运作用的敏感性的体现. 由于自然界存在的土壤类型多种多样，在确定其可蚀性值大小时，不论人力物力还是时间上都很难现实. 因此，国内外的很多研究都致力于土壤可蚀性的估算.Wischmeier 和Smith［18］根据实测的土壤可蚀性因子与土壤性质的相关性，建立了土壤可蚀性K值的诺谟方程，参数包括土壤质地、土壤有机质、土壤结构和土壤通透性等. 但是这种方法需要较多参数，特别是土壤结构级别和土壤渗透级别很难准确地获得，因此较为烦琐. Williams 和Sharply［19］在侵蚀-生产力评价模型EPIC（Erosion productivity impact calculator）中，把土壤可蚀性因子K值的计算公式简化成只与土壤砂粒含量、粉粒含量、黏粒含量和土壤有机碳含量有关的公式，建立了一种相对简便的土壤可蚀性因子K值的计算方法. 我国学者也对土壤可蚀性K值进行了大量研究，涉及研究区土壤K值的修正诺谟公式、修正EPIC模型以及只考虑几何平均粒径的计算方法等.

土壤数据（质地、有机质等）来源于黑河计划数据管理中心（http：//westdc.westgis.ac.cn）的“基于世界土壤数据库的中国土壤数据集（v1.1）”，从中截取了河南省地区数据.

式中：SANd为砂粒含量，%；SILi为粉粒含量，%；CLAi为黏粒含量，%；C为有机质含量，%；SN1=1-SAN/100.

3.3 地形因子LS因子

地形是自然地理要素中最基本的要素，对土壤侵蚀强度有重大影响. 各地形因子中，坡度、坡长应用最为广泛，因为从动力学角度分析，所有的地形因子的作用最终都取决于地形坡度. 因此在土壤侵蚀的研究中，使用坡度、坡长相乘来计算地形因子LS值. 坡度和坡长因子是降雨侵蚀作用中的加速因子，反映了地形对土壤侵蚀的影响，决定着径流的运动状态和方向.

数字高程模型（Digital Elevation Model，DEM）是描述地表每个单元的高程属性和空间位置分布的有序集合. DEM数据的采集方式主要有地面测量、摄影测量、激光扫描、地形图手扶跟踪数字化、立体遥感、GPS、地形图屏幕数字化、干涉雷达等. 其常用的表达形式有规则格网和不规则三角网两种. 本文选取空间分辨率为30 m的DEM数据.

3.3.1 坡度因子S坡度是影响土壤侵蚀的另一个重要参数. 在土壤侵蚀量与坡度的关系上，一般建立的是幂函数或二次多项式关系. McCool等［20］通过分析大量的野外径流试验小区资料和人工降雨资料，分别在0～9%和9%～25%坡度范围内建立了土壤侵蚀量与坡度的线性关系. Liu等［21］在McCool等的基础上，结合中国黄土高原绥德、安塞和天水三个水土保持试验站的试验资料，建立了25%以上的陡坡坡度和土壤侵蚀量的关系式. 本文采用Liu等和McCool等结合的方法，计算坡度因子.

式中：θ为像元坡度，%.

3.3.2 坡长因子L坡长L因子是反映地形地貌特征对土壤侵蚀的影响，是水土保持、土壤侵蚀等研究中的基本因子. 坡长能够直接影响地面径流的速度，进而影响对地面土壤的侵蚀力. USLE模型中，坡长L因子是指其他条件（降雨、土壤、耕作、水土保持措施、坡度等）相同的情况下，某一给定坡长的坡面上，土壤流失量与标准径流小区（坡长22.13 m，坡度9%）典型坡面土壤流失量的比值介于0～1之间，反映地形地貌特征对土壤侵蚀的影响. Wischmeier和Smith［18］对坡长指数的改进结果在USLE模型中得到了广泛应用.

式中：L为坡长因子；l为像元坡长；m为坡长指数.m取值如下式：

传统坡长计算一般是在野外或地形图上进行，不仅费时费力，而且不能获取研究区任意点的坡长. 地理信息系统、数字高程模型等技术的发展，使得坡长计算的效率得到大幅提高，本文即是通过ArcGIS 软件和DEM进行坡长提取. 由于GIS专业处理软件中没有直接求坡长的功能，本文参考利用汤国安和杨昕［22］在水文分析功能和空间分析应用时提出的一种快速近似计算坡长的方法，在假定坡面水流方向与山脊线（分水岭）垂直基础上，计算每个点沿垂直方向到山脊线的实际汇流路径长度为坡长的近似值.

计算步骤如下：

1）对DEM数据进行填挖，再对填挖后的DEM数据提取负地形，获得由原来沟道转化成的山脊线.

2）求负地形的水流方向并计算流水累积量.

3）提取负地形的河网. 使用条件函数工具con（accum>100，1），即像元流入量超过100 的所有像元均将成为河流网络的一部分.

4）对山脊线进行欧氏距离计算，得到近似坡长l.

5）使用条件函数计算m值.

6）采用公式计算出坡长因子L.

3.4 作物管理因子C

本研究区域的植被覆盖度估算的数据源选用Landsat-8卫星遥感影像，遥感植被指数选用归一化差异植被指数（Normal Differential Vegetation Index，NDVI），方法选用像元二分模型. 植被是土壤侵蚀的主要抑制因素，C因子是这种抑制作用的定量化表达指标之一.C因子指一定植被覆盖条件耕作农地上的土壤流失量与同等条件下实施清耕的连续休闲地的土壤流失量之比，其值大小变化于0～1 之间，C值越大，说明对应的土地利用类型的土壤侵蚀越严重，其计算公式采用蔡崇法等［23］的研究结果.

3.5 治理措施因子P

USLE模型中，治理措施因子P定义为采用特定措施土地上的土壤流失量与顺坡种植的土壤流失量的比值，介于0～1之间. 水土保持措施通过采取一定的耕作方式或改变小地形、拦蓄地表径流、增加降雨入渗以充分利用光热和水土资源，改善农业生产条件达到保水保土保肥、减少土壤侵蚀和提高农业生产的目的，主要包括耕作措施和工程措施两大类. 据调查，流域内的农地主要为退耕还林后保留的基本农田，大多已完成梯田改造；林地大部分是退耕还林工程造林，通常的侵蚀控制措施有梯田和淤地坝，未采取水土保持措施的土地利用类型P因子值为1. 具体参照刘宝元等的《土壤侵蚀预报模型》［24］并根据划分的土地利用类型结合坡度信息对P因子进行修正，确定修正后的P因子值，获得P因子图层.

4 结果与分析

4.1 土壤侵蚀强度等级变化

根据土壤侵蚀分类分级标准（SL 190—2007），豫西黄土丘陵地区的土壤侵蚀等级包括微度、轻度、中度、强度、极强度、剧烈等6 个级别［25-26］. 从表1 可以看出，2015 年到2020 年土壤侵蚀等级面积在豫西黄土丘陵地区发生变化，具体表现为：微度侵蚀中10 137.42 km2未变化，0.9 km2转变为轻度侵蚀；轻度侵蚀中6 566.85 km2未变化，100.08 km2转变为微度侵蚀，1.34 km2转变为中度侵蚀，0.31 km2转变为强度侵蚀；中度侵蚀中3 580.29 km2未变化，7.22 km2转变为微度侵蚀，69.59 km2转变为轻度侵蚀，0.9 km2转变为强度侵蚀，1.06 km2转变为极强度侵蚀；强度侵蚀中2 469.92 km2未变化，0.69 km2转变为微度侵蚀，35.33 km2转变为轻度侵蚀，17.39 km2转变为中度侵蚀，0.69 km2转变为极强度侵蚀，0.46 km2转变为剧烈侵蚀；极强度侵蚀中2 657.36 km2未变化，0.23 km2转变为微度侵蚀，0.04 km2转变为轻度侵蚀，41.28 km2转变为中度侵蚀，4.59 km2转变为强度侵蚀，0.76 km2转变为剧烈侵蚀；剧烈侵蚀中1 313.67 km2未变化，0.02 km2转变为微度侵蚀，9.49 km2转变为强度侵蚀，2.51 km2转变为极强度侵蚀.

表1 2015—2020年土壤侵蚀强度等级转移矩阵Tab.1 Transfer matrix of soil erosion intensity grade from 2015 to 2020 单位：km2

4.2 土壤侵蚀时间变化

从表2可以看出，豫西黄土丘陵地区2015年和2020年从微度到剧烈的土壤侵蚀面积和占比的变化均呈整体减少趋势，局部均呈先减少再增加再减少的双峰波动变化，峰值在强度和极强，土壤侵蚀面积和占比都在微度最大、剧烈最小. 从2015年到2020年，豫西黄土丘陵地区的土壤侵蚀面积和占比在不同侵蚀强度等级发生不同变化. 土壤侵蚀面积和占比分别为：微度分别增加105.62 km2、0.4%；轻度分别增加3.29 km2、0.02%；中度分别减少19.17 km2、0.07%；强度分别减少39.37 km2、0.14%；极强分别减少39.55 km2、0.15%；剧烈分别减少10.81 km2、0.04%. 由此可见，随着时间推移，豫西黄土丘陵地区的强土壤侵蚀不断减少，逐渐趋于轻微土壤侵蚀.

表2 2015年和2020年豫西黄土丘陵地区各土壤侵蚀强度的面积Tab.2 Area of soil erosion intensity in loess hilly area of western Henan from 2015 to 2020

4.3 土壤侵蚀空间变化

由图2 可以看出，豫西黄土丘陵地区的大部分地区发生轻、微、中程度的土壤侵蚀，少部分地区发生强度以上土壤侵蚀. 多数轻微程度土壤侵蚀出现在研究区的东北部、中部以及零散分布在其他地区，而强度以上土壤侵蚀出现在研究区的北部和南部边缘地带以及东南部的部分地区，从行政区上看，中度侵蚀在卢氏县、灵宝市最广泛，强度侵蚀在卢氏县、灵宝市、嵩县最广泛，极强度侵蚀在卢氏县、嵩县最广泛，剧烈侵蚀在嵩县最广泛.

图2 2020年豫西黄土丘陵地区的土壤侵蚀强度的空间分布图Fig.2 Spatial distribution of soil erosion intensity in loess hilly area of western Henan in 2020

5 结论

本文利用RUSLE 模型对豫西黄土丘陵地区2015—2020年间区域水土保持功能的变化特征进行研究，得到如下结论：

1）与2015年相比，区域的水土保持功能得到了显著了提升，区域平均侵蚀模数总体属于中度侵蚀量级，对全区域平均土壤侵蚀量下降率约为1.03%. 从不同侵蚀量级看，中度及以上面积都有所下降，总面积下降了约112 km2；微度面积增加明显，增加了约106 km2.

2）2020年，豫西区域土壤侵蚀量级总体上以中度及以下为主，占比为76.09%，其中微度为37.93%，轻度24.69%，中度13.47%；强度及以上占比为23.91%，其中强度占9.2%，极强占9.85%，剧烈占4.87%. 土壤侵蚀防治重点区域主要分布于区域南部洛宁、汝阳、嵩县，西部的灵宝，北部的济源、渑池等.