基于RUSLE模型的齐齐哈尔市土壤侵蚀时空分布

Spatiotemporal Distribution of Soil Erosion in Qiqihar City Based on the RUSLE Model

摘要: 为评估齐齐哈尔市近30年土壤侵蚀状况并明确其主导驱动因子，本文基于Google Earth Engine平台，融合多源遥感与地理数据，构建RUSLE模型，对1990~2020年间的土壤侵蚀模数进行定量估算，并引入地理探测器分析其影响机制。结果表明：土壤侵蚀以微度和轻度为主，空间上呈“中部低、两侧高”格局；坡度为主控因子，高坡地侵蚀显著增强。耕地因人为扰动普遍存在轻中度侵蚀，林地在局部退化区亦有高强度侵蚀。坡度与降雨、土地利用间存在显著交互作用。建议未来加强坡耕地治理与水土保持工程，提升土壤资源可持续利用水平。

Abstract: To assess soil erosion in Qiqihar City over the past 30 years and identify its key driving factors, this study utilized the Google Earth Engine platform to integrate multi-source remote sensing and geographic data. The RUSLE model was applied to quantify soil erosion modulus for the period 1990~2020, and the GeoDetector model was introduced to analyze influencing mechanisms. Results show that soil erosion was mainly slight and mild, with a spatial pattern of “low in the center and high on both sides.” Slope was the dominant factor, with significantly intensified erosion in steep areas. Due to frequent human disturbance, cultivated land exhibited widespread mild to moderate erosion, while forested areas also experienced severe erosion in some degraded and steep regions. Significant interaction effects were observed between slope, rainfall, and land use. It is recommended to strengthen slope farmland management and soil conservation practices to promote sustainable use of soil resources.

1. 引言

土壤侵蚀强度与空间分布是制定科学水土保持策略的基础，而修正通用土壤流失方程(RUSLE)作为国际公认的经典经验模型，由于其结构简单、参数物理意义明确，且易于与地理信息系统(GIS)集成，该方法在区域土壤侵蚀研究中得到了广泛应用。该模型通过耦合五个因子，分别是降雨侵蚀力(R)、土壤可蚀性(K)、坡长坡度(LS)、植被覆盖与管理(C)和水土保持措施(P) [1]，以此来定量预测年均土壤流失量，为区域生态保护与规划提供重要依据。

RUSLE模型被广泛应用于不同尺度的土壤侵蚀研究[2]。在区域尺度上，学者们针对西南山区、黄土高原、长江流域等典型区域，结合GIS空间分析技术，定量评估了土壤侵蚀的时空演变特征。在流域尺度上，浙东四明山、祖厉河上游等案例通过多因子动态监测，揭示了小流域土壤侵蚀与景观格局的互馈机制。此外，模型还被应用于矿区复垦、城市水土流失等特殊场景，实现对其跨领域的适应性。

目前，在土壤侵蚀研究领域，RUSLE模型的应用研究主要围绕三个关键方面：一是获取更为精确的模型参数；通过高分辨率DEM提取地形起伏度替代传统坡长坡度因子计算，结合多源遥感数据动态更新植被覆盖因子(C)，提升空间分辨率。二是多尺度的耦合分析；在流域尺度上耦合气象、土壤、土地利用等数据，揭示不同景观格局对侵蚀的调控作用；三是模型与其他技术的集成。将RUSLE与WEPP、EPIC等物理过程模型结合，弥补经验模型在陡坡侵蚀与泥沙输移模拟中的不足。未来研究仍需深化模型参数本地化修正，强化多源数据融合与多尺度耦合分析，以提升RUSLE模型在复杂环境下的预测能力，为全球土壤侵蚀防控提供科学支撑。

本研究旨在以齐齐哈尔市为研究对象，利用现有的GEE平台，结合RUSLE模型和GIS技术，对齐齐哈尔市近30年土壤侵蚀状况进行评估检测，探究齐齐哈尔市土壤侵蚀特征及其影响因素的时空变换规律，旨在为齐齐哈尔市开展水土保持工作以及推进侵蚀沟治理提供可靠依据。

2. 材料与方法

2.1. 研究区概况

齐齐哈尔市位于中国的东北地区，其地理位置大致介于46˚N~48˚N，122˚至126˚E之间。主要地貌类型为平原，处于温带大陆性季风气候带，四季气候特征鲜明，冬季持续时间长，表现为寒冷且干燥的气候状况；夏季相对短暂，呈现出温暖湿润的气候特点。7~9月份降雨集中，年均降水量介于400 mm~550 mm之间。齐齐哈尔市以草原和森林为主要自然植被类型，土壤类型丰富多样，主要为黑土、草甸土及沼泽土等。

2.2. 数据与方法

2.2.1. 数据来源

数据源信息如表1所示。

Table1. Data sources

表1. 数据来源

数据	数据来源
降雨量数据	CHIRPS数据集 (https://code.earthengine.google.com/)
DEM数据	ALOSDEM (https://search.asf.alaska.edu/)
土地覆盖类型数据	GLC_FCS30D数据集 (https://data.casearth.cn/)
Lansat5TM/Landsat8OLI数据	美国地质勘探局 (http://www.usgs.gov/)

2.2.2. 数据预处理

根据RUSLE模型，土壤侵蚀模数是由降雨侵蚀力、土壤可蚀性、坡度坡长、植被覆盖与管理以及水土保持措施五个影响因子累乘所得到的结果[3]，综合考虑了自然以及人为因素等对土壤侵蚀程度的影响，进而能较精准地对区域土壤侵蚀量进行估算。

在实际应用中，由于数据来源的多样性以及不同因子测量和获取方式的差异，每种侵蚀因子对应的栅格数据中单个栅格像元的空间大小并不一致，需统一重采样至同一空间分辨率为30 m进行计算。通过重采样操作，使得每个侵蚀因子的栅格数据具有相同的空间尺度和网格结构，从而能够准确地进行累乘计算，最终得到土壤侵蚀模数。这一过程为后续的土壤侵蚀空间分析、风险评估以及防治措施的制定提供了可靠的数据基础。

2.2.3. RUSLE模型

依据修正后的通用土壤侵蚀方程RUSLE，可计算得出土壤侵蚀模数，具体计算公式如下[4]：

$A=R\times K\times LS\times C\times P$

式中：A为年均土壤流失量[t/(hm²·a)]；R为降雨侵蚀力因子[(MJ·mm)/(hm²·h·a)]；K为土壤可蚀性因子[(t·hm²·h)/(MJ·mm·hm²)]；LS为坡长坡度因子，无量纲；C为植被覆盖因子，无量纲；P为水土保持措施因子，无量纲。

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

利用获取的年平均降雨相关数据资料，可依据下列式子计算出各年的降雨侵蚀力因子(R) [5]：

$R=0.053P_j^{1.6548}$

式中：R为降雨侵蚀力因子[(MJ·mm)/(hm²·h·a)]；$P_j$ 为年均降雨量(mm)。

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

土壤可蚀性因子计算方法较多，鉴于对数据的可获得性以及公式的简洁性的考虑，本研究采用Wishchmeier等人提出的算法进行计算[6]，见式：

$K_{EPIC}=\left\{0.2+0.3\exp \left[-0.0256SAN\left(1-\frac{SIL}{100}\right)\right]\right\}{\left(\frac{SIL}{CLA+SIL}\right)}^{0.3}$

$\left(1.0-\frac{0.25C}{C+exp\left(3.72-2.95C\right)}\right)\left(1.0-\frac{0.7SN}{SN+exp\left(22.9SN-5.51\right)}\right)$

式中：SAN为砂粒(0.050 mm~2.000 mm)含量(%)；SIL为粉粒(0.002 mm~0.050 mm)含量(%)；CLA为黏粒(<0.002 mm)含量(%)；C为有机碳含量(%)；SN = 1 − SAN/100。

计算后的结果应乘以0.1317，由美国制转换为国际制，单位应为(t·hm²·h)/(MJ·mm·hm²) [7]。

(3) 坡度坡长因子(LS)

采用符素华等人[8]开发的土壤侵蚀模型地形因子计算工具(2.0版本)中涉及的公式对坡长坡度因子进行计算。

坡度因子的计算，见下式：

$\{\begin{array}{l}S=10.8sin\theta \theta <5˚\\ S=10.8sin\theta -0.55˚\le \theta <10˚\\ S=21.91sin\theta -0.96\theta \ge 10˚\end{array}$

式中：S为坡度因子，无量纲；$\theta$ 为坡度(˚)。

坡长因子的计算，见下式：

$\{\begin{array}{l}m=0.2\theta <0.5˚\\ m=0.30.5˚\le \theta <1.5˚\\ m=0.41.5˚\le \theta <3˚\\ m=0.53˚\le \theta \end{array}$

式中m为坡长指数，无量纲；$\theta$ 为坡度(˚)。

坡长因子L随着汇流面积$A_{out}$ 增大而增大，并受到指数m的控制，而m又由坡度决定。

$\text{L=}\{\begin{array}{c}\begin{array}{cc}\frac{A_{out}^{m+1}-0^{m+1}}{A_x^m\cdot {22.13}^m}& A_{out}>0\end{array}\\ \begin{array}{cc}0& A_{out}\le 0\end{array}\end{array}$

式中：$A_{out}$ 为单位宽度上的流域累计流量(flowaccumulation)；m为坡长指数，无量纲；$A_x$ 为参考坡长，通常设为30 m，对应标准小区长度。

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

本研究采用蔡崇法等人[9]提出的公式先计算出植被覆盖度(FVC)进而估算C值，其值介于0~1之间，无量纲。

$FVC=\frac{NDVI-NDV{I}_{\min }}{NDV{I}_{\max }-NDV{I}_{\min }}$

式中：FVC为植被覆盖度；NDVI为归一化植被指数；NDVImax为归一化植被指数的最大值，即植被完全覆盖像元的NDVI值；NDVImin为归一化植被指数的最小值，即无植被覆盖像元的NDVI值。

$C=\{\begin{array}{l}\begin{array}{cc}1& FVC\le 0.095\end{array}\\ \begin{array}{cc}0.6508-0.3436lgFVC& 0.095<FVC\le 0.783\end{array}\\ \begin{array}{cc}0& FVC>0.783\end{array}\end{array}$

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

结合现有研究，针对研究区土地覆盖的具体情况，本文对水土保持措施因子P值分类赋值，具体见表2 [10]。

Table 2. P factor assignment

表2. P因子赋值

土地覆盖类型	林地	灌丛	草地	农田 < 5˚	农田 ≥ 5˚	裸地	冰雪	水体	建设用地	湿地	其他
P因子赋值	0.05	0.16	0.05	0.25	0.4	0.8	0	0	0	0	0.8

2.2.4. 地理探测器

地理探测器(GeoDetector)是由中国科学院王劲峰研究员团队提出的一种用于识别空间异质性及其驱动机制的统计分析工具。其核心理念基于假设：如果某一影响因子(自变量)对某现象(因变量)具有显著影响，那么二者的空间分布应具有一致性。该方法通过对空间数据进行分层和方差分解，定量测度各驱动因子的解释力，广泛应用于土壤侵蚀、生态环境、城市扩展等领域的驱动因素识别与解释。

本研究选取分异探测器与交互探测器对齐齐哈尔市1990~2020年间土壤侵蚀影响因子进行定量探测与交互分析，明确主要驱动因素及其耦合关系，为揭示区域土壤侵蚀形成机制提供依据。分异探测器用于衡量某一因子对因变量空间分异的解释力，采用q值进行度量，q值范围为0~1，值越大表示解释力越强。交互探测器则是用于判断两个因子联合作用是否具有增强或减弱效应，并区分线性增强、非线性增强、独立作用等类型。

3. 土壤侵蚀时空演变与坡度、土地利用类型的响应分析

3.1. 土壤侵蚀时空分布特征

根据SL190-2007《土壤侵蚀分类标准》[11]，将土壤侵蚀强度划分为6个等级分别为：微度侵蚀为<200 (t/km²·a)，轻度侵蚀为200~2500 (t/km²·a)，中度侵蚀为2500~5000 (t/km²·a)，强烈侵蚀为5000~8000 (t/km²·a)，极强烈侵蚀为8000~15,000 (t/km²·a)，剧烈侵蚀为>15,000 (t/km²·a)。

结合图1的对比分析与表3中不同侵蚀强度的面积统计可知，研究区土壤侵蚀状况于这30年期间展现出特定的时空演变规律。从时间维度看，1990~2020年间齐齐哈尔市土壤侵蚀总体呈现出轻度向中高强度等级转化的趋势。其中，轻度侵蚀面积在2000年达到最高值(27695.79 km²)，但自2010年起开始明显下降，至2020年已降至20314.10 km²，表明在一定程度上水土保持措施产生了积极效应。与此同时，中度及以上侵蚀等级面积自2010年起显著增加，尤其是强烈侵蚀面积在2020年跃升至467.05 km²，较2010年增长近6.9倍，极强烈和剧烈侵蚀也同步扩大，显示出部分地区生态系统稳定性下降，土壤侵蚀风险加剧。

在空间格局上，土壤侵蚀等级总体呈现出“中部低、两侧高”的分布态势。研究区中部为松嫩平原腹地，地势相对平坦，土壤类型较稳定，农业耕作集中但管理措施相对完善，因而微度和轻度侵蚀占据主导地位。相比之下，东部与西部丘陵区因受地形起伏大、坡度陡峭以及植被破碎等因素影响，更易发生中度及以上等级侵蚀，特别是在坡耕地密集、土地退化严重区域，出现了小范围的极强烈与剧烈侵蚀斑块。

Table 3. Area statistics by erosion intensity from 1990 to 2020

表3. 1990~2020年不同侵蚀强度的面积统计

侵蚀强度	1990年	2000年	2010年	2020年
微度侵蚀	11568.97	11918.98	11729.62	14612.79
轻度侵蚀	26985.43	27695.79	27059.58	20314.10
中度侵蚀	1141.97	232.61	1005.11	4338.54
强烈侵蚀	135.21	19.27	67.90	467.05
极强烈侵蚀	34.44	0.80	5.16	121.03
剧烈侵蚀	1.43	0.00	0.08	13.95

注：单位：km²。

Figure 1. Spatial distribution of soil erosion intensity from 1990 to 2020

图1. 1990~2020年土壤侵蚀强度分布图

3.2. 不同的坡度下土壤侵蚀强度的变化

在坡度与土壤侵蚀强度的关系分析中，不同坡度等级下的土壤侵蚀强度分布表现出显著差异。从图2可以看出，坡度越大，土壤侵蚀等级整体呈现由轻向重转变的趋势。

在坡度小于1˚的地区，土壤侵蚀以微度侵蚀为主，说明平坦区域的土壤侵蚀风险极低。随着坡度的增加，轻度和中度侵蚀的比例逐渐升高，在坡度为2˚~3˚和3˚~5˚的区域，轻度侵蚀分别占比80%和67%，说明低缓坡地虽然存在一定程度的侵蚀，但强度仍以轻度为主。随着坡度的增加，中度及以上等级的侵蚀显著增加。例如，在坡度为5˚~8˚区域，中度侵蚀与强烈侵蚀的面积占比分别达到了34%和6%；而在坡度大于15˚的区域，除微度侵蚀外，还出现了6%的剧烈侵蚀和16%的极强侵蚀，表明陡坡区域的土壤侵蚀问题更为严重。

总体来看，齐齐哈尔市1990~2020年间土壤侵蚀在空间上具有明显的坡度分异性。随着坡度的增加，土壤侵蚀等级由轻逐渐向中、重度发展，尤其在中高坡度区域(>5˚)，中度及以上等级侵蚀逐渐成为主要类型，反映了坡度是影响区域土壤侵蚀强度的重要驱动因子。

Figure 2. Proportion of annual average soil erosion intensity area under different slope classes

图2. 不同坡度等级下的年均土壤侵蚀强度面积占比

3.3. 不同的土地利用类型下土壤侵蚀强度的变化

为定量分析不同土地利用类型对土壤侵蚀的影响，本文以齐齐哈尔市耕地、林地和草地为代表，基于1990~2020年平均土壤侵蚀强度的面积占比(图3)开展对比分析。结果显示，不同土地利用类型的土壤侵蚀强度差异明显，且表现出一定的区域地理特征性。

其中，耕地土壤侵蚀以轻度为主，占比66%，微度侵蚀占29%，中度及以上侵蚀等级约占5%。作为齐齐哈尔市的主要农耕用地，耕地分布广泛，集中于中部平原区，地势较为平缓，但由于长期耕作、地表扰动频繁、秸秆覆盖不足等人为活动影响，使其在降雨季节更易发生表层侵蚀。林地的微度侵蚀占比达85%，轻度侵蚀仅为11%，其余等级侵蚀面积极低，显示出较好的生态防护效能。然而，相比草地，林地中出现了略高比例的中度及以上侵蚀现象，这在一定程度上与齐齐哈尔市林地的空间分布及地形条件有关。研究区内林地多分布于东部丘陵和山地地带，坡度较大，部分区域林下植被稀疏，加之近年来部分区域存在林地退化、采伐干扰等现象，使得局部地区土壤抗蚀能力下降，侵蚀等级升高。草地类型的土壤侵蚀则以轻度侵蚀为主(60%)，微度侵蚀占比为29%，中度及以上侵蚀等级略高于林地，但整体仍处于较低水平。这与草地主要分布于平坦地区、地表扰动较小、部分区域有天然草被覆盖等因素有关。

林地的水土保持能力优于草地，但在齐齐哈尔地区，由于林地所处地势陡峭、易发生坡面侵蚀，且部分区域存在人为干扰，导致其在局部区域的侵蚀强度反而略高于草地。这一现象表明：土地利用类型对土壤侵蚀的影响不仅与其本身属性有关，还与空间分布、地形条件和人为干预密切相关。

Figure 3. Proportion of annual average soil erosion intensity area under different land use types

图3. 不同土地利用类型下的年均土壤侵蚀强度面积占比

4. 土壤侵蚀多因子时空特征与驱动机制探究

4.1. RUSLE模型因子时空分布特征分析

4.1.1. 降雨侵蚀力因子(R)

Figure 4. Distribution map of rainfall erosivity R factor

图4. 降雨侵蚀力R因子分布图

从图4可以看出，在时间维度上，1990年最高值为349.501(MJ·mm)/(hm²·h·a)，2000年降至186.114 (MJ·mm)/(hm²·h·a)，2010年为291.758 (MJ·mm)/(hm²·h·a)，2020年又上升至484.074(MJ·mm)/(hm²·h·a)，这表明齐齐哈尔市的降雨侵蚀力在这40年间存在波动变化，并非单调递增或递减。进一步分析，1990~2000年R值整体呈下降趋势，而2010~2020年则呈上升趋势，这种变化可能与气候变化、大气环流等因素有关，为研究气候变化对区域土壤侵蚀的影响提供了科学参考。

从空间变化分布格局来看，不同年份的降雨侵蚀力因R值在各区域间的分布存在显著差异。具体而言，齐齐哈尔东部地区在各个年份均呈现出较高的R值，这表明该区域属于降雨侵蚀力的高值区域，相较于其他区域可能更易发生土壤侵蚀现象；反观中部地区，其R值则相对较为稳定[11]。

4.1.2. 土壤可蚀性因子(K)

如图5所示，齐齐哈尔市土壤可蚀性因子(K)在空间分布上呈现显著的不均衡性上，呈现出明显的东高西低、北高南低的格局。K值的整体变化范围在0.0157至0.0412 (t·hm²·h)/(MJ·mm·hm²)之间，反映了不同区域土壤对侵蚀的敏感性，其中东部和东北部丘陵地带明显K值普遍偏高，土壤颗粒细、结构松散，易于降雨或地表径流作用下发生侵蚀。而西部及中南部地区K值较低，多为平坦开阔的耕地或草地，土壤黏粒含量高，结构更稳定，抗侵蚀能力强。

从空间格局来看，K因子的分布状况与地形起伏、土壤质地特性及土地利用方式紧密相连。K值较高的区域往往兼具坡度较大、植被覆盖稀疏等特征，这些区域应当在未来水土保持与农业管理工作中被列为重点防护对象。因此，在规划区域水土流失防治策略的过程中，必须充分了解不同土壤类型在空间上的差异性，对K值较高的敏感区域实施更为精细化的分区管理措施，从而有效提升水土资源的可持续利用效率的同时，增强生态系统的平衡稳定性。

Figure 5. Distribution map of soil erodibility K factor

图5. 土壤可蚀性K因子分布图

4.1.3. 坡度坡长因子(LS)

地势平坦区域因坡度小、坡长短，水流速度慢、汇聚少，导致侵蚀力弱，坡度坡长因子对其土壤侵蚀影响较小。如图6蓝色区域所示，该类地形有利于水土保持与农业生产。而坡度坡长因子高值区域则表现为坡陡坡长，水流加速汇聚，侵蚀能力显著增强，成为水土流失的高风险区。

Figure 6. Distribution map of slope length and steepness LS factor

图6. 坡度坡长LS因子分布图

4.1.4. 植被覆盖与管理因子(C)

齐齐哈尔市植被覆盖区域及裸露或低覆盖区域的范围和分布在不同年份呈现出明显变化(图7)。1990年，市域内存在较大范围的低植被覆盖区，可能与大量裸露土地、农业用地处于非生长季节等因素有关。植被覆盖较好的区域零星分布，主要对应林地、草地及部分生长季的农田。

至2000年，绿色区域的面积和连片程度显著增加，表明高植被覆盖区有所扩展。这一变化可能归因于生态保护措施的推进、农业结构调整或气候条件的改善，提升了原本植被较差区域的生长环境。

2010年，低植被覆盖区再次扩大，可能受到干旱等不利气象条件影响，抑制了农田和草地植被生长，也可能因城市扩张、农业开垦等人类活动加剧了植被破坏。到2020年，植被覆盖度有所恢复与扩展，说明近年来实施的退耕还林还草、湿地保护等生态政策取得了积极成效，促进了植被恢复。

从空间分布看，植被覆盖在市域内呈不均衡格局。边缘地区由于地形、土壤条件适宜且人类干扰较少，植被覆盖相对较好；而市中心及人口密集区域，受城市建设和工业活动影响，植被覆盖度较低。河流、湖泊等水域周边植被状况则相对较好，受水资源充足的影响明显。

Figure 7. Distribution map of vegetation cover and management C factor

图7. 植被覆盖与管理C因子分布图

4.1.5. 水土保持措施因子(P)

Figure 8. Distribution map of conservation practice P factor

图8. 水土保持措施P因子分布图

P因子取值根据土地覆盖类型确定，范围为0~0.8，数值越小表明水土保持效果越显著。1990年齐齐哈尔市P因子高值区呈零散分布，集中于部分人类活动强度较高区域，如城乡结合部、工业集中区等。如图8所示，2000年高值区范围变化不大，仅在局部位置和形态上有所调整，反映出土地利用和生态环境的局部变化。2010年起高值区数量与分布范围明显增加，可能受城市扩张与农业结构调整影响。至2020年，高值区进一步扩展，分布更广，反映出区域发展、人类活动增强及土地利用优化等因素对P因子空间格局的影响。

4.2. 驱动力分析

4.2.1. 单因子探测分析

本研究选取坡度、坡向、高程、植被覆盖度、土地利用类型和年均降雨量六个因子，借助地理探测器模型对1990、2000、2010和2020年齐齐哈尔市土壤侵蚀模数的空间分异驱动机制进行定量评估。通过计算q值，揭示各因子对土壤侵蚀的解释力。

图9结果表明，自然地形因子(坡度、坡向、高程)对土壤侵蚀的控制作用显著高于生态和气候因子，其中坡度为最主要驱动因子，其q值始终高于0.45，在2000年达到最大值0.5296。高坡地形显著增加地表径流速度与侵蚀强度，是土壤侵蚀的核心影响因子。高程和土地利用类型的解释力次之，并呈上升趋势，反映出地貌起伏及人类活动强度变化对侵蚀格局的逐步强化。随着城市化发展，土地利用方式改变导致部分区域侵蚀风险升高。植被覆盖度作为生态因子，对土壤侵蚀的解释力保持在0.07~0.11之间，稳定发挥减缓水流、固土保水的生态功能。年均降雨量和坡向的q值较低(<0.10)，空间变异性弱，虽为直接或间接影响因素，但在本区域尺度上解释力有限。

总体而言，齐齐哈尔市土壤侵蚀格局以地形因子主导，生态与土地利用因素次之，气候因子作用相对较弱。未来水土保持应聚焦高风险坡地与不合理用地区域，优先实施退耕还林(草)、优化土地利用结构及增强植被恢复等措施，以实现土壤侵蚀的有效控制与资源可持续利用。

Figure 9. q-values of single-factor explanatory power

图9. 单因子解释力q值

4.2.2. 交互因子探测分析

为深入探究土壤侵蚀驱动机制的复合效应，全面揭示各驱动因子之间的协同作用对土壤侵蚀过程的影响，本研究运用地理探测器的交互探测模块，对齐齐哈尔市1990~2020年期间各主导因子之间的交互作用关系展开了系统且深入的分析。地理探测器的交互探测模块能够有效地评估不同因子之间相互作用对目标变量(本研究中为土壤侵蚀)的综合影响，相较于单因子探测，更能揭示因子间复杂的相互作用关系以及由此产生的非线性和空间异质性效应。该方法不仅提升了因子识别的解释能力，也为多因子协同机制提供了量化依据，特别适用于地理环境系统中复杂耦合过程的研究。在面向区域尺度的生态环境管理中，多因子交互分析为科学制定土壤保护与利用策略提供了理论支撑和决策基础。

研究结果表明，齐齐哈尔市土壤侵蚀驱动因子间普遍存在增强型或非线性增强型交互效应，多因子组合对土壤侵蚀的解释力明显优于单因子，体现出土壤侵蚀过程的多因子协同特征与复杂性。

由表4可以看出，“年均降雨量∩坡度”的交互作用在1990~2020年间始终具有较强解释力(平均 q  =  0.568)，并呈非线性增强关系。降雨与坡度共同影响地表径流强度，在陡坡地区更易引发强烈侵蚀。“土地利用类型∩坡度”交互效应更为显著(平均 q  =  0.644)，反映出人类土地利用方式在不同坡度条件下对侵蚀风险的放大作用。坡耕地的耕作方式直接影响侵蚀强度，合理管理(如等高耕作、梯田等)可有效降低风险。此外，“高程∩植被覆盖度”(平均q  =  0.304)与“植被覆盖度∩降雨量”(平均q  =  0.282)等组合亦表现出稳定增强效应，说明地形、生态与气候因子共同调控土壤侵蚀空间格局。

总体而言，多因子交互探测有效揭示了土壤侵蚀的非线性与空间异质性特征。建议未来水土保持策略重点关注坡度与降雨、土地利用等因子耦合区域，通过分区治理与生态修复相结合措施，实现对高风险区的精准防治与生态系统的持续优化。

Table 4. Interaction detection of driving factors

表4. 驱动因子交互探测

交互因子	1990年	2000年	2010年	2020年
高程∩植被覆盖度	0.2895^↑↑	0.2902^↑↑	0.3381^↑↑	0.3248^↑↑
高程∩土地利用类型	0.2023^↑	0.2370^↑	0.3124^↑	0.2774^↑
高程∩年均降雨量	0.1498^↑	0.1949^↑	0.2395^↑	0.1956^↑
高程∩坡度	0.5388^↑	0.5960^↑	0.5966^↑	0.5438^↑
植被覆盖度∩土地利用类型	0.1946^↑↑	0.2087^↑↑	0.1945^↑	0.2981^↑↑
植被覆盖度∩年均降雨量	0.1826^↑↑	0.1564^↑↑	0.1967^↑↑	0.1680^↑↑
植被覆盖度∩坡度	0.8053^↑↑	0.8091^↑↑	0.6989^↑↑	0.7868^↑↑
土地利用类型∩年均降雨量	0.1373^↑	0.1394^↑↑	0.2220^↑	0.2117^↑↑
土地利用类型∩坡度	0.6693^↑↑	0.7283^↑↑	0.7118^↑↑	0.7074^↑↑
年均降雨量∩坡度	0.5149^↑	0.5773^↑	0.5988^↑↑	0.5112^↑↑

注：^↑表示双因子增强；^↑↑表示非线性增强。

5. 结论

本文基于Google Earth Engine平台，集成多源遥感数据与地理信息技术，采用修正通用土壤流失方程(RUSLE)对齐齐哈尔市1990~2020年间的土壤侵蚀强度与空间分布格局进行了系统评估，并结合地理探测器模型探讨其主导驱动因素及交互机制。主要结论如下：

(1) 齐齐哈尔市土壤侵蚀以微度与轻度为主，整体侵蚀强度在1990~2020年间呈现“先降后升”的波动特征。空间上，侵蚀高值区主要分布于东部丘陵与西部坡耕地集中区域，中部平原区则整体侵蚀程度较低。

(2) 不同坡度等级对土壤侵蚀等级结构具有显著影响，随着坡度的增大，中度及以上强度侵蚀比例显著上升，>15˚区域成为剧烈侵蚀高发区，显示出坡度作为主控因子的作用。

(3) 土地利用类型显著影响土壤侵蚀特征。耕地因人为扰动频繁、表层裸露程度高，其轻度及中度侵蚀发生率远高于林地与草地。虽然林地具备较强的水土保持能力，但由于主要分布于东部山地丘陵区、坡度较大，且存在局部退化、采伐干扰等问题，导致部分林地区域中重度侵蚀比例高于草地。

(4) 地理探测器分析结果表明，坡度是最主要的驱动因子，其q值在各年份均为最高；高程、土地利用与植被覆盖度亦具一定解释力。交互探测结果进一步揭示：坡度与降雨、土地利用等因子之间存在显著的非线性增强关系，说明土壤侵蚀是多因子耦合影响下的复合过程。

综上，齐齐哈尔市土壤侵蚀受地形地貌、人类活动与气候因素共同驱动，呈现出典型的空间异质性与时序变异性。建议今后水土保持工作应重点关注中高坡耕地、林地退化区等重点敏感区域，优化土地利用结构，强化坡耕地生态治理，提升区域生态系统的稳定性与土壤资源的可持续利用效率。