陈方磊 王计平 程复 谢海燕
摘要：采用日降雨量、DEM、土壤类型、泥沙含量及多期NDVI等数据，基于修正通用土壤流失方程（RUSLE）和地理探测器，研究了国家生态退耕还林还草工程实施前后近41年延河甘谷驿流域土壤侵蚀动态与驱动因子。结果表明，1980—2020年研究区土壤侵蚀强度总体呈波动变化趋势，1980年、1990年、2000年、2010年和2020年平均侵蚀模数分别为6 746.30、5 740.28、6 389.56、5 450.46、5 480.56 t/（km2·年）。1980—2000年研究区整体侵蚀强度逐渐增强，强烈及以上等级侵蚀面积占比逐渐增加，表现为“增蚀升级”的特点；
2000年后研究区内土壤侵蚀强度开始降低，强烈及以上等级的侵蚀面积减少，总体表现为“减蚀降级”的特点。研究区土壤侵蚀强度随着坡度的升高而加剧，同时发现海拔1 000～ 1 200 m和1 200～1 400 m是研究区内侵蚀发生的主要高程带。2020年土地利用类型因子解释力最为突出，表明退耕还林还草工程实施效果显著，大面积的耕地向林草地转换是使得研究区2000年后土壤侵蚀强度降低的最主要原因。土壤侵蚀各影响因子的协同作用明显强于单一因子的影响。
关键词：土壤侵蚀；
修正通用土壤流失方程（RUSLE）；
地理探测器；
驱动因子；
延河甘谷驿流域
中图分类号：S157.1 文献标识码：A
文章编号：0439-8114（2024）06-0027-08
DOI：10.14088/j.cnki.issn0439-8114.2024.06.005 开放科学（资源服务）标识码（OSID）：
Evaluation of soil erosion and analysis of driving factors in the Ganguyi Watershed of Yanhe River from 1980 to 2020
CHEN Fang-lei1， WANG Ji-ping2，3， CHENG Fu4， XIE Hai-yan1
（1.College of Resources and Environment， Xinjiang Agricultural University， Urumqi 830052， China；
2.Institute of Ecological Protection and Restoration， Chinese Academy of Forestry， Beijing 100091， China；
3.Research Center of Saline and Alkali Land， National State Forestry and Grassland Administration， Beijing 100091， China；
4.Water and Soil Conservation Monitoring Center， Ministry of Water Resources of the Peoples Republic of China， Beijing 100053， China）
Abstract：
The daily rainfall data， DEM data， soil type data， sediment content data and multi period NDVI data were used to study the soil erosion dynamics and driving factors in the Ganguyi Watershed of the Yanhe River in the past 41 years before and after the implementation of the national ecological rehabilitation project of returning farmland to forest and grassland based on the Revised General Soil Loss Equation （RUSLE） and geographic detectors. The results showed that， from 1980 to 2020， the overall soil erosion intensity in the study area showed a fluctuating trend， with an average erosion modulus of 6 746.30 t/（km2·a）， 5 740.28 t/（km2·a）， 6 389.56 t/（km2·a）， 5 450.46 t/（km2·a） and 5 480.56 t/（km2·a） in 1980， 1990， 2000， 2010 and 2020， respectively. From 1980 to 2000， the overall erosion intensity in the study area gradually increased， and the proportion of erosion areas at the strong level and above gradually increased， which was characterized by “erosion increase and upgrading”. After 2000， the intensity of soil erosion in the study area began to decrease， and the area of erosion at the strong level and above decreased， which was characterized by “erosion reduction and degradation”. The intensity of soil erosion in the study area increased with the increase of slope. At the same time， it was found that 1 000～1 200 m and 1 200～1 400 m were the main elevation zones for erosion occurrence in the study area. The explanatory power of land use type factors was most prominent in 2020， indicating that the implementation of the project of returning farmland to forests and grasslands had a significant effect. The conversion of large areas of farmland to forests and grasslands was the main reason for the decrease in soil erosion intensity in the research area after 2000. The synergistic effect of various influencing factors on soil erosion was significantly stronger than that of a single factor.
Key words：
soil erosion；
Revised Universal Soil Loss Equation （RUSLE）；
geographic detector；
driving factor；
Ganguyi Watershed of Yanhe River
土壤侵蚀是指结构不稳定的土壤因流水或降雨而被冲刷、剥蚀、搬运，因重力作用而失稳移动，因风力作用而悬浮移动、沉积，或因冻融作用而形成泥质流体的现象［1］。中国是世界上土壤侵蚀状况最为严重的国家，目前受侵蚀的土壤总面积达26万km2，西北黄土高原地区尤为严重［2］。黄土高原不仅是中国水土保持工程治理与生态修复的重点地区，同时也是土壤侵蚀机理、治理理论、方法及防治技术研究与实践关注的热点区［3］。因此，开展黄土高原小流域范围内水土流失动态变化分析变得尤为重要。
为了开展土壤侵蚀的定量评估等相关工作，国内外水土保持研究相关学者提出了众多土壤侵蚀相关模型，这些模型可分为物理过程性模型、产量模型和经验模型。物理过程性模型主要有WEPP模型等［4］；
产量模型有EPIC［5］和PI［6］；
原理简单和适用性强是经验模型的优势所在，主要有土壤流失方程（USLE）［4］、修正通用土壤流失方程（RUSLE）［7］和中国土壤流失方程（CSLE）等［8］。黄土高原土壤侵蚀的相关研究表明，RUSLE可以更为准确地模拟黄土高原土壤侵蚀状况［9-11］。
延河流域是黄土高原水土流失最为严重的区域［12］，整个流域受水土流失影响，土壤肥力下降，水库淤积，河床上升，旱灾频发，农业生产条件落后，生态环境脆弱，严重制约经济社会高质量发展［13］。以小流域为单元的水土流失综合治理是中国目前防治水土流失的主要方式［14］。本研究以延河流域甘谷驿水文站控制区为研究区，运用RUSLE模型，分析该地区1980—2020年长时间研究序列的土壤侵蚀时空变化特征，阐明研究区近41年土壤侵蚀变化特征及驱动因素，以期为该区未来制定土壤侵蚀应对策略提供理论依据。
1 数据与方法
1.1 研究区概况
延河流域甘谷驿水文站控制区（36°22′N—37°19′、108°39′—109°48′E，图1）位于黄河流域黄土高原中部，控制区总面积为5 872 km2。延河流域地表破碎，黄土侵蚀极为剧烈，沟间地以长梁、斜梁、梁峁和峁为主，同时也有少量残塬存在，部分地区有石质丘陵与裸露基岩，整个流域土壤稳定性较差，因此长期遭受强烈剥蚀［15］。
1.2 数据来源
数据包括逐日降雨量数据、土地利用类型数据、土壤类型数据、DEM数据、归一化植被指数（Normalized difference vegetation index，NDVI）数据、泥沙含量数据和人口数据。降雨量数据下载于国家气象科学数据中心-中国气象数据网（http：//data.cma.cn） ，选用“中国地面气候资料日数据集（V3.0）”中延河流域甘谷驿水文站控制区内及周边共13个气象站点（大宁、靖边、清涧、横山、甘泉、志丹、延长、延川、延安、宜川、安塞、子长、吴起）1980—2020年逐日降雨量数据；
水文站输沙模数数据来源于水利部黄河水利委员会编制的黄河泥沙公报（www.yrcc.gov.cn）。土壤类型和质地数据来源于世界土壤数据库中国土壤数据集（Harmonized world soil database version 1.1，HWSD）。地形数据为30 m分辨率的DEM数据，来源于地理空间数据云（http：//www.gscloud.cn/）。归一化植被指数数据来自Google Earth Engine平台，由于研究内容时序跨度大，目前没有单一数据产品可以满足。本研究采用250 m分辨率MODIS（Moderate-resolution imaging spectroradiometer）-NDVI产品和8 km分辨率GIMMS（Global inventory modelling and mapping studies）-NDVI产品作为数据源，通过合成处理得到1981—2020年（1980年数据缺失，1980—1981年未发生重大变化）NDVI数据。1980年、1990年、2000年、2010年和2020年5期土地利用数据均来源于中国科学院地理科学与资源研究所/地理国情监测云平台（http：//www.dsac.cn）。人口数据来自LandScan数据库。
1.3 土壤侵蚀模数计算
计算土壤侵蚀模数选用修正后的通用土壤流失方程（RUSLE），表达式如下。
A=R×L×S×K×C×P （1）
式中，A为土壤侵蚀模数；
R为降雨侵蚀力因子；
L、S分别为坡长、坡度因子；
K为土壤可蚀性因子；
C为植被覆盖与管理因子；
P为水土保持措施因子。
1.3.1 降雨侵蚀力因子 采用章文波等［16］提出的基于逐日降雨量的降雨侵蚀力计算方法，其表达式如下。
[R半月=αk=1n（Pk）β] （2）
[R年=i=124R半月i] （3）
[β=0.836 3+18.144/Pd12+24.455/Py12] （4）
[α=21.586β-7.189 1] （5）
相关学者根据黄河流域径流实验站降雨和径流观测资料的分析结果拟定侵蚀性降雨量标准定为12 mm［17］。式中，[R半月]为半月降雨侵蚀力；
[R年]为年降雨侵蚀力；
k为某半月内侵蚀性降雨日数；
[Pk]为半月内第k天的侵蚀性降雨日雨量；
[Pd12]为一年内侵蚀性降雨日雨量均值；
[Py12]为侵蚀性降雨年总量的多年均值；
[α]、[β]为该模型的两个参数。
由于降雨具有偶然性，为了减少单个年份降雨存在的偶然性，故将1980—2020年分为5个时段计算多年平均降雨侵蚀力，第1时段为1980—1985年（时段Ⅰ），第2时段为1986—1995年（时段Ⅱ），第3时段为1996—2005年（时段Ⅲ），第4时段为2006—2015年（时段Ⅳ），第5时段为2016—2020年（时段Ⅴ）。最后采用反距离权重法［18］进行空间插值，得到研究区各个时段降雨侵蚀力的空间分布。
1.3.2 土壤可蚀性因子 以世界土壤数据库中国土壤数据集为基础，运用EPIC模型计算研究区内各类型土壤的K，计算式如下。
[K=0.2+0.3exp-0.025 6SAN1-SIL100×SILCLA+SIL0.3×1-0.25CC+exp3.72-2.95C×1-0.7SN1SN1+exp-5.51+22.9SN1] （6）
式中，SAN为沙粒含量；
SIL为粉粒含量；
CLA为黏粒含量；
C为有机碳含量；
SN1=1-SAN/100。
1.3.3 坡长、坡度因子 坡长、坡度因子是评估土壤侵蚀的重要参数，本研究采用符素华等［19］提出的坡长、坡度因子计算方法，计算式如下。
[S=10.8sin θ+0.03 θ≤5°16.8sin θ-0.05 5°<θ≤10°21.9sin θ-0.96 θ>10°] （7）
[L=λ/22.1m， m0.2 θ≤1° 0.3 1°<θ≤3°0.4 3°<θ≤5°0.5 θ>5°] （8）
式中，λ表示坡长；
m为坡长系数；
[θ]为坡度。
1.3.4 植被覆盖与管理因子 为满足1980—2020年长时间序列的研究需求，需选取2种数据源相同时序（2002—2006年）NDVI月最大合成数据，将MODIS-NDVI（250 m）平均聚合至8 km空间分辨率，得出逐月空间分布系数，并将MODIS-NDVI（8 km）与GIMMS-NDVI（8 km）进行回归统计，得到二者的回归方程；
利用所得的回归拟合方程对GIMMS-NDVI（8 km）数据进行修正，随后将修正后的结果与MODIS-NDVI（250 m）逐月数据的空间分布系数相乘，最终得到空间分辨率为250 m的GIMMS-NDVI数据［20］。由于1980年GIMMS-NDVI数据缺失，该地区历史相关资料表明研究区在1980—1981年并未发生大规模土地利用类型转变，同时未发生重大气候灾害，故用1981年NDVI数据作为1980年缺失数据代入公式计算植被覆盖与管理因子。数据降尺度前后对比如图2所示。
土壤侵蚀强度在不同的植被覆盖类型和管理方式下呈现出差异。植被覆盖与管理因子在土壤侵蚀计算模型中反映植被覆盖和管理措施对土壤侵蚀的影响，取值范围为0～1。本研究参考相关研究［21］成果对研究区C进行赋值和计算，耕地赋值为0.44，水域和生产建设用地的植被覆盖度较低，土壤侵蚀强度较小，故将其C赋值为0，未利用地赋值为1。根据江忠善等［22］提出的方法，草地和林地的植被覆盖与管理因子计算式如下。
[CG=1 fvc≤5%exp-0.041 8fvc-5 fvc>5%] （9）
[CF=1 fvc≤5%exp-0.008 5fvc-51.5 fvc>5%] （10）
式中，CG为草地植被覆盖与管理因子；
CF为林地植被覆盖与管理因子，fvc表示植被覆盖度。根据研究区土地利用数据和基于NDVI数据计算的植被覆盖度，得到研究区5个时期的植被覆盖与管理因子。
1.3.5 水土保持措施因子 土地利用类型可以反映出水土保持措施因子，本研究参考文献［23］，采用经验赋值法得到研究区不同时期的水土保持措施因子P分布，林地、草地和耕地分别赋值0.05、0.16和0.31，水域和未利用地均赋值1.0，生产建设用地赋值0。
1.4 地理探测器
地理探测器是探测空间分异性并揭示其背后驱动因子的统计学方法［24］。采用该工具中的单因子探测器来定量评估延河流域甘谷驿水文站控制区土壤强度变化的影响因子，定量分析各因子对土壤侵蚀强度的影响程度。在地理探测器中，运用因子探测器来探测因变量的空间分异性，用q表示自变量对因变量的解释度，表达式如下。
[q=1-?=1LN? σ2?Nσ2=1-SSWSST] （11）
[SSW=?=1LN?σ2?] （12）
[SST=Nσ2] （13）
式中，L为自变量或因变量的分层或分区；
Nh和N分别为层h和全区的单元数；
[σ2?]和[σ2]分别为层h和全区Y值的方差；
SSW和SST分别为层内方差之和与全区总方差；
q的值域为［0，1］［25］。
交互作用探测器通过对比单因子q和双因子q来确定2个因子间的交互作用方向和方式，交互作用的判定依据见表1。
2 结果与分析
2.1 流域土壤侵蚀强度及其时空变化总特征
按照上述各因子方法逐个计算1980年、1990年、2000年、2010年和2020年共5期延河流域甘谷驿水文站控制区单位面积土壤侵蚀模数，根据《土壤侵蚀分类分级标准》对研究区土壤侵蚀强度计算结果进行分级，利用GIS生成延河流域甘谷驿水文站控制区土壤侵蚀强度分布（图3）。查阅黄河水资源和泥沙公报中延河流域甘谷驿水文站年输沙模数数据，对比发现本研究土壤侵蚀计算结果变化趋势与延河甘谷驿水文站实测输沙模数数据变化趋势一致。同时发现1987—2010年年均输沙模数为5 510 t/（km2·年），计算出的1980—2010年年均土壤侵蚀模数约为6 082 t/（km2·年），计算结果比输沙模数高出10.4%，呈现出这样结果的主要原因是由于RULSE是基于坡面来计算坡面侵蚀的模型，在工程措施以外同时存在泥沙沉淀和重力侵蚀，因此通过模型计算所得的侵蚀模数与水文站观测到的输沙模数有一定的出入。
由土壤侵蚀强度分级可以看出，与1980年相比，1990年土壤侵蚀强度空间分布特征变化较小，中度及以下的侵蚀主要分布在河谷周围和研究区南部植被覆盖度较高的地带；
2000年与前2期结果不同，研究区中部地带土壤侵蚀强度明显增大；
2010年较2000年侵蚀情况变化最为明显，沿河道向外两侧侵蚀强度显著降低；
2020年流域西北部高海拔地区侵蚀强度较高，研究区中部和南部区域侵蚀强度较低，基本上以中度、轻度和微度侵蚀为主。
土壤侵蚀模数计算结果表明，1980—2020年延河流域甘谷驿水文站控制区年均土壤侵蚀强度呈波动变化的发展趋势，1980年平均土壤侵蚀模数为6 746.30 t/（km2·年）；
1990年平均侵蚀模数减小到5 740.28 t/（km2·年），比1980年减少14.9%；
2000年平均侵蚀模数增加到 6 389.56 t/（km2·年），比1990年增加11.3%；
2000年后研究区平均侵蚀模数开始大幅降低，2010年和2020年平均侵蚀模数分别为5 450.46、5 480.56 t/（km2·年），相较于2000年分别减少14.7%和14.2%。
表2为研究区1980—2020年各侵蚀强度等级的面积及占比情况。通过分析5期土壤侵蚀强度的面积分布可知，1980年、1990年、2000年、2010年和2020年微度和轻度侵蚀面积占比较大，分别为48.42%、50.58%、46.31%、56.02%和62.85%，强烈及以上等级侵蚀面积占比分别为40.06%、37.36%、42.36%、31.99%和27.90%，均小于各时期微度和轻度侵蚀所占面积。
从5期土壤侵蚀强度面积变化来看，1980—1990年，土壤侵蚀强度只有剧烈侵蚀的面积减小，面积占比也降低3.29个百分点，可见剧烈侵蚀面积减少是该段时间内流域土壤侵蚀量减小的主要原因；
1990—2000年，强烈及以上等级侵蚀面积占比累计增加了5.00个百分点，微度、轻度、中度侵蚀面积占比分别减少了1.84、2.43、0.73个百分点，可见该时间段内土壤侵蚀量的增加与强烈及以上等级侵蚀面积增加有关；
2000—2010年，强烈及以上等级侵蚀面积占比共计减少了10.37个百分点，可见该时间段内土壤侵蚀强度的大幅减弱与强烈及以上等级侵蚀面积的减少有关；
2010—2020年，剧烈侵蚀面积增加0.64个百分点，土壤侵蚀量的细微增加可能与剧烈侵蚀面积的增加有关。
2.2 流域土壤侵蚀强度的空间分布变化
为进一步分析延河流域甘谷驿水文站控制区退耕还林还草措施实施前后土壤侵蚀的空间分布变化，利用GIS对侵蚀结果进行叠加分析，得出土壤侵蚀强度空间转化，结果如图4所示。由图4可知，1980—2000年，土壤侵蚀强度空间分布变化明显，整体上呈侵蚀强度升级的趋势，流域中部地区坪桥镇、建华镇、化子坪镇、真武洞镇尤为显著；
2000—2020年土壤侵蚀强度减弱面积明显大于强度增大的面积，土壤侵蚀强度减弱的面积占变化总面积的67.5%，土壤侵蚀量在这个时段下降明显，侵蚀强度减弱主要发生在化子坪镇、招安镇、真武洞镇和建华镇等流域中部地区的乡镇。
通过分析1980—2020年土壤侵蚀强度空间分布变化结果，得到1980—2000年、2000—2020年2个时段研究区的土壤侵蚀强度等级转移矩阵，见表3和表4。1980—2000年，流域侵蚀强度整体呈增大趋势，土壤侵蚀强度增大的区域面积占比为16.1%，说明在这期间研究区土壤侵蚀强度加剧。2000—2020年，47.7%的区域土壤侵蚀强度发生变化，其中高侵蚀等级向低侵蚀等级转换的面积占比为35.5%，这期间土壤侵蚀强度降低明显。
2.3 流域土壤侵蚀强度的地形分布变化
根据研究区的海拔范围，分为<1 000 m、1 000～1 200 m、1 200～1 400 m、1 400～1 600 m和1 600～ 1 800 m共5个海拔段，基于GIS进行叠加分析和统计，得出研究区1980—2020年不同海拔段的土壤侵蚀参数，结果见表5。可以看出，研究区内土壤侵蚀强度与海拔间密切相关，随着海拔的升高，土壤侵蚀强度先明显增大后减小。在<1 000 m高程带，1980—2020年土壤侵蚀面积占比较低且强度变化不明显，微度侵蚀面积占比最高。从侵蚀面积占比来看，1 000～1 200 m和1 200～1 400 m是研究区内侵蚀发生的主要高程带，在这2个高程带内，2000年之前各侵蚀强度所占比例变化并不明显，2000年后微度和轻度侵蚀等级面积占比明显增大，强烈及以上等级侵蚀面积减少，表明导致2000年后研究区土壤侵蚀状况向好发展的主要改善区域在该高程带内。研究区1 400～1 600 m高程带内各等级土壤侵蚀强度面积占比变化不大，强烈和极强烈侵蚀面积占比呈减小趋势。1 600～1 800 m高程带内各时期土壤侵蚀面积占比几乎不发生变化，可见该高程带并非近年来改善水土保持措施实施的重点区域。
坡度是影响坡面土壤侵蚀模型计算结果中的重要影响因子［19］。利用流域高程数据提取坡度信息，将坡度按照≤5°、5°～8°、8°～15°、15°～25°和>25°进行分级，随后将坡度分级结果与5期土壤侵蚀强度结果带入ArcGIS软件进行叠加分析，得到5期不同坡度的土壤侵蚀分布状况，结果见表6。总体来看，研究区土壤侵蚀强度随着坡度的增大而加剧。坡度
≤5°的侵蚀强度较弱，不发生强烈及以上等级的侵蚀，15°～25°和>25°是研究区内土壤侵蚀强度最高的地区。1980—2000年，8°～15°坡度带中度侵蚀面积占比增加，这也是导致在此期间土壤侵蚀强度整体上升的原因之一。2000年以后15°～25°和>25°坡度带内强烈、极强烈和剧烈侵蚀面积明显减少，这与当地在该时段所实施的生态恢复相关工程有密切关联。
2.4 流域土壤侵蚀驱动因子分析
为进一步正确认识在黄土高原地区开展的“退耕还林”工程背景下，工程后研究区内土壤侵蚀的自然与社会因素对土壤侵蚀的影响，选取研究区平均降雨量（X1）、植被覆盖度（X2）、海拔（X3）、坡度（X4）、土地利用类型（X5）、土壤类型（X6）和人口密度（X7）7个因子作为代入运算的自变量，以2020年土壤侵蚀强度作为因变量Y代入地理探测器中进行运算。地理探测器中要求输入的自变量为类型数据，本研究采用王劲峰等［26］提出的数据离散化方法，使用自然断点法将时段内降雨量、高程和人口密度数据离散化为6类，土壤类型数据和土地利用按各自类别分为6类，植被覆盖度数据分为≤0.3、0.3～0.4、0.4～0.5、0.5～0.6、0.6～0.7、0.7～0.8、0.8～0.9、0.9～1.0共8类，坡度分为≤5°、5°～10°、10°～15°、15°～20°、20°～25°、>25°共6类，采用GIS软件的渔网功能，将研究区划分成1 km×1 km格网，共提取5 870个采样点代入地理探测器运行，影响土壤侵蚀结果的7个因子q由高到低依次为土地利用类型（0.380 0）>坡度（0.145 0）>海拔（0.026 1）>植被覆盖度（0.026 0）>土壤类型（0.024 0）>平均降雨量（0.011 0）>人口密度（0.009 0），结果显示所有因子均通过显著性检验（P<0.000 1）。
因子探测器的应用结果表明，退耕还林还草工程实施后，不同因子对土壤侵蚀强度的解释力存在差异，土地利用类型因子的解释力最强，q为0.380 0，是影响研究区土壤侵蚀空间分布的主导因子。延河流域地理位置和环境特殊，流域内土壤稳定性较差，生态环境脆弱。自1999年开始的退耕还林还草工程实施以来，黄土高原土地利用类型发生根本性改变［27，28］，根据计算分析得到研究区2000—2020年耕地向林地和草地转移面积分别为204.75 km2和502.68 km2，植被覆盖度显著增加，生态环境得到明显改善，根据探测结果，研究区内坡度因子也是影响土壤侵蚀的另一主导因子。上述结果表明，自然和人为因素共同影响着延河流域甘谷驿水文站控制区土壤侵蚀空间分布格局。
利用交互探测器研究了延河流域甘谷驿水文站控制区土壤侵蚀和各驱动因素间的交互作用，发现各因素间并不存在相互独立的作用，而是表现为非线性增强和双因子增强2种交互作用，结果如表7所示。平均降雨量与土壤类型、平均降雨量与人口密度、海拔与坡度、海拔与土壤类型、海拔与人口密度、坡度与土壤类型间的交互作用以双因子增强的形式影响土壤侵蚀分布格局变化，土壤类型与人口密度交互探测结果为单因子非线性增强，其他因子交互探测结果均为非线性增强。土地利用类型与植被覆盖度、海拔、坡度、土壤类型的交互影响力均在0.40以上，其中土地利用类型与坡度的交互影响力为0.600 0，表明这2个因子组合对研究区土壤侵蚀分布格局影响最大。
表7 2020年土壤侵蚀各因子交互探测结果
[q1 q2 q1+q2 q1∩q2 结果 X1=0.011 0 X2=0.026 0 0.037 0 0.045 0 非线性增强 X1=0.011 0 X3=0.026 1 0.037 1 0.043 0 非线性增强 X1=0.011 0 X4=0.145 0 0.156 0 0.166 0 非线性增强 X1=0.011 0 X5=0.380 0 0.391 0 0.398 0 非线性增强 X1=0.011 0 X6=0.024 0 0.035 0 0.034 0 双因子增强 X1=0.011 0 X7=0.009 0 0.020 0 0.016 0 双因子增强 X2=0.026 0 X3=0.026 1 0.052 1 0.056 0 非线性增强 X2=0.026 0 X4=0.145 0 0.171 0 0.198 0 非线性增强 X2=0.026 0 X5=0.380 0 0.406 0 0.427 0 非线性增强 X2=0.026 0 X6=0.024 0 0.050 0 0.060 0 非线性增强 X2=0.026 0 X7=0.009 0 0.035 0 0.047 0 非线性增强 X3=0.026 1 X4=0.145 0 0.171 1 0.167 0 双因子增强 X3=0.026 1 X5=0.380 0 0.406 1 0.410 0 非线性增强 X3=0.026 1 X6=0.024 0 0.050 1 0.043 0 双因子增强 X3=0.026 1 X7=0.009 0 0.035 1 0.034 0 双因子增强 X4=0.145 0 X5=0.380 0 0.525 0 0.600 0 非线性增强 X4=0.145 0 X6=0.024 0 0.169 0 0.166 0 双因子增强 X4=0.145 0 X7=0.009 0 0.154 0 0.161 0 非线性增强 X5=0.380 0 X6=0.024 0 0.404 0 0.420 0 非线性增强 X5=0.380 0 X7=0.009 0 0.389 0 0.393 0 非线性增强 X6=0.024 0 X7=0.009 0 0.033 0 0.013 2 单因子非线性增强 ]
3 小结
1）利用RUSLE模型对延河流域甘谷驿水文站控制区1980—2020年土壤侵蚀强度进行计算，结果表明，1980年平均土壤侵蚀模数为6 746.30 t/（km2·年）；
1990年平均侵蚀模数减小到5 740.28 t/（km2·年），比1980年减少14.9%；
2000年平均侵蚀模数增加到 6 389.56 t/（km2·年），比1990年增加11.3%；
2000年后研究区平均侵蚀模数开始大幅降低，2010年和2020年平均侵蚀模数相较于2000年分别减少14.7%和14.2%。
2）延河流域甘谷驿水文站控制区1980—2000年土壤侵蚀强度逐渐增强，强烈及以上等级侵蚀面积占比逐渐增加，表现为“增蚀升级”的特点。2000年后研究区内土壤侵蚀强度开始降低，强烈及以上等级的侵蚀面积减少，轻度和微度侵蚀面积增大，总体表现为“减蚀降级”的特点。结合DEM数据分析得出研究区土壤侵蚀强度随着坡度的升高而加剧，坡度≤5°的侵蚀强度较弱，不发生强烈及以上等级的侵蚀，15°～25°和>25°是研究区发生土壤侵蚀最为严重的地区，同一坡度区间内各时段土壤侵蚀强度变化特征同区域整体变化特征一致。同时发现1 000～1 200 m和1 200～1 400 m是研究区内侵蚀发生的主要高程带。
3）通过地理探测器对研究区退耕还林还草工程实施背景下土壤侵蚀强度的影响因素分析发现，土地利用类型因子解释力较为突出，表明退耕还林还草工程实施后水土流失治理效果显著，大面积的耕地向林草地转换是研究区2000年后土壤侵蚀强度降低的最主要原因。交互探测结果表明，各影响因子的协同作用明显强于单一因子的影响。
本研究以RUSLE模型计算结果作为基础，对研究区1980—2020年土壤侵蚀的时空变化规律进行探索，使用地理探测器对研究区7个影响因子进行探测分析，通过对研究结果分析可知，研究中采用的方法合理，预期研究任务基本完成。相较以往研究，本研究着重突出近41年长时间序列数据的变化，可以较好地降低研究时段内数据突变的影响，同时能更加全面地掌握研究区土壤侵蚀的变化特点和归因，使研究结果更有价值。本研究仍存在不足之处，在计算结果验证中，如果可以获取研究区水土保持实地监测数据，模型计算结果会更具说服力。针对地理探测器的使用，如果在今后研究中可以加入GDP和种植作物类型数据等社会经济因子，研究结果会更有价值。
参考文献：
［1］ 洪宇辰，姬鑫慧，戴晓爱.基于GIS/RS和USLE模型的土壤侵蚀研究及主要侵蚀因子识别［J］.测绘， 2017， 40（6）：
277-283.
［2］ CHEN R Y，YAN D C，WEN A B， et al. The regional difference in engineering-control and tillage factors of Chinese Soil Loss Equation［J］. Journal of mountain science， 2021， 18（3）：
658-670.
［3］ 火 红.基于地形因子的延安市土壤侵蚀及其地表驱动因素的时空分异特征研究［D］.西安：长安大学， 2021.
［4］ 陈朝良，赵广举，穆兴民，等.基于RUSLE模型的湟水流域土壤侵蚀时空变化［J］.水土保持学报， 2021， 35（4）：
73-79.
［5］ 方兴义.基于EPIC模型的农业旱灾风险模糊评估方法［J］.灾害学，2020， 35（3）：
55-58.
［6］ 段兴武，谢 云，张玉平，等. PI模型在东北松嫩黑土区土壤生产力评价中的应用［J］.中国农学通报， 2010， 26（8）：
179-188.
［7］ RENARD K G，FERREIRA V A. RUSLE model description and database sensitivity［J］. Journal of environmental quality，1993，22（3）：458-466.
［8］ 周正朝，上官周平.土壤侵蚀模型研究综述［J］.中国水土保持科学， 2004，2（1）：
52-56.
［9］ 郭 达，宋小宁，董 震，等.基于RUSLE与GIS的黄土高原水土流失评价研究——以宁夏中卫地区为例［J］.泥沙研究， 2020， 45（5）：
55-60.
［10］ 王泽宇，陈旭阳，马彩诗，等.陕北榆林市退耕还林前后土壤侵蚀及生态服务价值变化［J］.西北林学院学报， 2021， 36（3）：
59-67.
［11］ 张 素，熊东红，吴 汉，等.基于RUSLE模型的孙水河流域土壤侵蚀空间分异特征［J］.水土保持学报， 2021， 35（5）：
24-30.
［12］ 钟 雪.延河流域社会水文耦合规律研究［D］.西安：西北大学， 2020.
［13］ 张永峰，王明华.延安市水土保持高质量发展先行区示范创建［J］.中国水土保持， 2022（9）：
37-40.
［14］ 王福岭，吴 鹏，郭子琦.以小流域为单元的水土流失强度评估［J］.中国水土保持， 2023（1）：
44-47.
［15］ 孙 虎.陕西延河流域地貌组合类型的模糊聚类划分［J］.陕西师范大学学报（自然科学版）， 1996，24（4）：
83-88.
［16］ 章文波，谢 云，刘宝元.利用日雨量计算降雨侵蚀力的方法研究［J］.地理科学， 2002，22（6）：
705-711.
［17］ 谢 云，刘宝元，章文波.侵蚀性降雨标准研究［J］.水土保持学报，2000，14（4）：
6-11.
［18］ 邹 强.流域降雨空间插值方法比较［J］.节水灌溉， 2009（7）：
12-14.
［19］ 符素华， 刘宝元， 周贵云， 等. 坡长坡度因子计算工具［J］. 中国水土保持科学， 2015， 13（5）：
105-110.
［20］ 林锦阔.河西地区土壤侵蚀时空分异及其驱动因素［D］.兰州：兰州大学， 2020.
［21］ 秦 伟，朱清科，张 岩.基于GIS和RUSLE的黄土高原小流域土壤侵蚀评估［J］.农业工程学报， 2009， 25（8）：157-163，4.
［22］ 江忠善，郑粉莉，武 敏.中国坡面水蚀预报模型研究［J］.泥沙研究， 2005（4）：
1-6.
［23］ 闫 瑞，张晓萍，李够霞，等.基于RUSLE的北洛河上游流域侵蚀产沙模拟研究［J］.水土保持学报， 2017， 31（4）：
32-37.
［24］ 宋丽洁，戴昭鑫，刘 新. 2000—2020年东营市人口分布的时空特征演变规律及影响因素分析［J］.北京测绘， 2022， 36（3）：
260-265.
［25］ 姜旭海，韩 玲，白宗璠，等.内蒙古自治区沙漠化敏感性时空演变格局和趋势分析［J］.生态学报， 2023， 43（1）：
364-378.
［26］ 王劲峰，徐成东.地理探测器：原理与展望［J］.地理学报，2017，72（1）：116-134.
［27］ 龚直文，姚顺波，党晶晶.黄土高原退耕前后土地利用/覆被变化及驱动力分析——以志丹县为例［J］.林业经济，2017，39（2）：54-58， 66.
［28］ 杨 灿，魏天兴，李亦然，等.黄土高原水蚀风蚀交错区退耕还林工程前后NDVI时空变化特征［J］.北京林业大学学报， 2021， 43（6）：
83-91.