黄土丘陵区切沟侵蚀发育过程中土壤侵蚀阻力变化

Variations in Soil Erosion Resistance during Gully Development in the Hilly and Gully Regions of the Loess Plateau

摘要

目的 切沟形态特征是量化切沟侵蚀发育的重要参数，其演变过程受制于水流剪切力与土壤侵蚀阻力之间的动态平衡。然而，切沟侵蚀发育过程中土壤侵蚀阻力在切沟不同部位如何随发育阶段演变并反馈于形态演化，还缺乏系统的量化研究。 方法 以黄土丘陵区典型坡面切沟为研究对象，结合野外调查与室内分析，系统比较不同发育阶段（前、中、后期）、沟道形态（宽型、窄型）及发育部位（沟头、沟壁、沟底）的土壤侵蚀阻力特征。 结果 随切沟发育，土壤临界剪切力（τ_c ）显著上升［发育前、中、后期分别为（2.19±1.39）、（3.11±1.92）、（3.62±2.77） Pa］，而细沟可蚀性（k_d ）显著下降［发育前、中、后期分别为（320.46±424.21）、（121.28±145.99）、（51.20±36.93） cm³/（N·s）］，表明发育后期土壤侵蚀阻力显著增强。随机森林模型揭示，此增强过程主要与土壤与根系属性演变有关，有机质含量是影响τ_c 与k_d 的共同关键因子；τ_c 的提升与根重密度密切相关，而k_d 的降低则主要得益于团聚体稳定性的增加。表明切沟发育过程中的土壤结构重建与根系发育，共同促进土壤侵蚀阻力的总体提高。此演变趋势在空间上呈现显著分异，不同地貌部位的土壤侵蚀阻力差异显著，表现为沟头τ_c 最高［（4.59±2.21） Pa］、k_d 最低［（44.81±25.79） cm³/（N·s）］，沟壁τ_c 最低［（1.34±0.61） Pa）］，k_d 最高［（380.06±406.89） cm³/（N·s）］，沟底居中，反映沟壁坍塌与沟底下切是主要的活跃侵蚀形式。沟道形态在此基础上进一步调节侵蚀阻力的时空格局，相比于宽型沟，窄型沟具有更高的k_d 和更低的τ_c，表现出更强的侵蚀敏感性。线性混合效应模型从发育阶段、地貌部位和沟道形态等宏观地貌因子视角量化土壤侵蚀阻力的变异来源，结果表明，三者共同解释k_d 的52.2%和τ_c 的58.2%变异，且发育阶段与地貌部位存在极显著交互作用（p<0.001），说明土壤侵蚀阻力随切沟发育在不同部位表现出差异化响应，并受沟道形态调节。综合来看，发育阶段决定土壤与根系属性演化的总体趋势，地貌部位和沟道形态通过塑造局地水动力与物源条件调节此演化过程，共同控制土壤侵蚀阻力的时空格局。 结论 研究揭示切沟发育与土壤侵蚀阻力耦合的定量规律，为黄土丘陵区侵蚀沟防控与生态修复提供科学依据。

Abstract

Objective Gully morphological characteristics are key parameters for quantifying gully erosion development， and its evolution is governed by the dynamic balance between hydraulic shear stress and soil erosion resistance. However， a systematic quantitative understanding is lacking regarding how soil erosion resistance evolves at different gully positions with developmental stages and feeds back to morphological evolution during gully erosion development. Methods This study investigated typical gullies on slopes in the loess hilly region. Field surveys and laboratory analyses were conducted to systematically compare soil erosion resistance characteristics across different developmental stages （early， middle， late）， gully morphologies （wide， narrow）， and positions （head， wall， floor）. Results As gully development progressed， critical shear stress （τ_c ） of soil significantly increased from （2.19±1.39） Pa in the early stage to （3.11±1.92） Pa in the middle stage and to （3.62±2.77） Pa in the late stage. Rill erodibility （k_d ） significantly declined from （320.46±424.21） cm³/（N·s） in the early stage to （121.28±145.99） cm³/（N·s） in the middle stage and to （51.20±36.93） cm³/（N·s） in the late stage， indicating that soil erosion resistance was significantly enhanced in the later stage of development. Random forest model revealed this enhancement was primarily related to changes in soil and root properties， with soil organic matter content being a key factor influencing both τ_c and k_d . The increase in τ_c was closely associated with root weight density， while the decrease in k_d was mainly attributed to increasing aggregate stability. These patterns suggested that soil structure reconstruction and root development during gully development jointly contributed to the overall increase in soil erosion resistance. This evolution trend exhibited significant spatial differentiation， with pronounced differences in soil erosion resistance across different geomorphological positions. Gully head showed the highest τ_c ［（4.59±2.21） Pa］ and lowest k_d ［（44.81±25.79） cm³/（N·s）］， gully wall had the lowest τ_c ［（1.34±0.61） Pa］ and highest k_d ［（380.06±406.89） cm³/（N·s）］， and gully floor displayed intermediate values. These values indicated that wall collapse and floor incision were the primary active erosion modes. Gully morphology further regulated the spatiotemporal pattern of erosion resistance. Compared to wide gullies， narrow gullies exhibited higher k_d and lower τ_c， suggesting greater sensitivity to erosion. Linear mixed-effects modelling quantified the sources of variation in soil erosion resistance from the perspective of macro-geomorphological factors such as developmental stage， geomorphological position， and gully morphology. The results showed that these three factors jointly explained 52.2% of k_d and 58.2% of τ_c variation， with a significantly pronounced interaction between developmental stage and geomorphological position （p<0.001）. These findings demonstrated that soil erosion resistance exhibited differentiated responses at different positions as gully development progressed and was conditioned by gully morphology. Overall， the developmental stage determined the overall trend of soil and root property evolution， whereas geomorphological position and gully morphology modulated this evolution by shaping local hydrodynamics and provenance conditions， thereby jointly controlling the spatiotemporal pattern of soil erosion resistance. Conclusion This study reveals the quantitative law of the coupling between gully development and soil erosion resistance， providing a scientific basis for gully erosion control and ecological restoration in the loess hilly region.

Graphical abstract

关键词

沟道侵蚀 / 临界剪切力 / 细沟可蚀性 / 沟道形态 / 发育阶段 / 黄土高原

Key words

gully erosion / critical shear stress / rill erodibility / gully morphology / developmental stage / Loess Plateau

土壤侵蚀是土地退化的主要过程和驱动力，其中切沟侵蚀作为最严重的侵蚀方式之一，严重威胁人类社会可持续发展，是全球关注的重大环境问题^［1-2］。黄土高原切沟发育广泛，活跃期切沟沟头溯源、沟壁扩张与沟底下切作用强烈，是流域侵蚀产沙的主要来源^［3］。已有研究^［4］表明，沟道侵蚀产沙可占集水区总沉积物质量的10%~95%，在严重水土流失区甚至可达80%以上^［3，5］。切沟形态特征是量化切沟侵蚀发育的重要参数，其演变过程受制于水流剪切力与土壤侵蚀阻力之间的动态平衡，在一定程度上反映侵蚀过程与发育阶段^［6］。其中，水流剪切力是侵蚀发生的驱动力，需克服由土壤性质（如体积质量、质地等）决定的临界剪切力，才能启动土壤颗粒的剥离与输移^［7］。已有研究^［7］表明，土壤侵蚀速率通常可以表示为外加剪切力超过临界剪切力后的函数；在相同水力条件下，侵蚀阻力越低的土体往往表现出更高的侵蚀敏感性和更快的侵蚀速率。因此，沟头溯源、沟底下切与沟岸扩张等侵蚀过程是该平衡的具体表现形式，共同决定切沟的形态格局与演化方向^［8-9］。深入揭示切沟侵蚀演变与其土壤侵蚀阻力之间的互馈机制，对完善切沟侵蚀动力模型具有重要的意义。

切沟侵蚀的发育过程具有高度复杂性。由于沟头、沟底和沟岸等关键部位在土体结构、理化性质及根系分布存在显著差异，其土壤侵蚀阻力呈现明显的空间分异，导致侵蚀过程在空间上表现出非均衡发育的特征^［1，4］。因此，切沟形态既是其土壤侵蚀阻力空间分异演变的直观反映，也为从形态特征推断其发育阶段和活跃程度提供重要依据。目前，研究普遍认识到地形、植被与土壤性质对切沟土壤侵蚀阻力的重要调控作用^{［4-5，10］}，并常通过沟道形态来经验性地推断其发育阶段与活跃性。然而，土壤侵蚀阻力作为连接环境属性与侵蚀过程的关键力学枢纽，其随切沟发育阶段的系统演变规律、空间分异格局及受沟道形态等宏观因子的调控作用，目前尚缺乏系统的定量揭示。本研究以黄土丘陵区典型坡面切沟为研究对象，结合野外调查与室内实验，旨在：1）定量刻画土壤侵蚀阻力随切沟发育的演变轨迹；2）揭示驱动该演变的关键土壤与根系属性；3）阐明地貌部位与沟道形态对其时空分异格局的调控作用。最终揭示切沟侵蚀发育演变与土壤侵蚀阻力之间的耦合机制，为切沟侵蚀的定量评估与精准防治提供科学依据。

1 材料与方法

1.1 研究区概况

本研究选取陕北黄土高原丘陵沟壑区的沙堰沟小流域（37°54′28″~37°55′02″N，110°08′59″~10°09′39″E）作为研究对象（图1）。该流域总面积69 hm²，海拔1 024~1 124 m。气候类型属半干旱大陆性季风气候，年平均气温9.4 ℃，年平均降水量451 mm。区域内沟壑密布、地表破碎、地貌形态典型，主要土壤类型为黄绵土，水土流失严重^［11］。

1.2 典型沟道选取与样品采集

本研究以黄土高原丘陵沟壑区不同发育阶段的典型坡面切沟为对象，于2025年4月采用无人机高精度测绘（DJI P4 RTK）获取其几何形态参数（沟长、沟宽、沟深）及集水区地貌信息（坡度、面积、坡长）。切沟是由间歇性股流冲刷形成，边壁陡峻，深度介于0.50~30 m，断面呈“V”形或“U”形，纵剖面与所在坡面基本一致，多呈跌水状，无法横跨耕作的永久性侵蚀沟^［1-2］。切沟发育阶段的划分参照THWAITES等^［12］提出的分类框架，综合考虑平面形态、剖面形态及主导侵蚀过程。在此基础上，结合李建军等^［13］基于形态特征（如长宽比等）的分类思路，采用3个关键形态比值量化切沟发育程度：沟长与集水区长度之比（L_r，表征沟头溯源）、沟宽与集水区宽度之比（W_r，表征沟道拓宽），以及沟头下切高度与集水区相对高差之比（D_r，表征沟道下切）。各比值临界值基于实地观测设定：L_r与W_r的均为0.5，对应溯源或侧蚀作用强烈时沟头或沟岸可侵蚀至分水岭的最终稳定状态；D_r为0.25，对应由于坡面纵坡降，沟道剖面由“V”形转为稳定“U”形、下切显著减弱的阶段（如集水区高差30~50 m、下切深度8~10 m情形）。依据上述阈值（L_r=0.50，W_r=0.50，D_r=0.25）组合定义切沟形态类型。理论上存在8类组合，但实地未发现“短窄深”（L_r≤0.5，W_r≤0.5，D_r≥0.25）与“短宽深”（L_r≤0.50，W_r≥0.50，D_r≥0.25）样本，故最终选取的6条切沟分别代表短宽浅、长宽浅、长宽深、短窄浅、长窄浅、长窄深6类（表1），可归为宽型与窄型沟的发育前、中、后期，确保样本在研究区形态谱系中具有代表性。实地调查表明，不同发育阶段切沟沟道的植被生长状况存在明显差异。发育早期切沟植被以草本为主；中期开始出现低矮灌木；后期沟道坡面植被覆盖度显著提高，形成较为稳定的灌草群落，局部甚至出现乔木。

为系统研究沟头溯源、沟壁拓宽和沟底下切3类主导侵蚀过程对土壤可蚀性的影响，本研究以切沟发育阶段（前期、中期、后期）为首要分类框架。在每一发育阶段内分别选取宽型和窄型典型切沟，在每条切沟的沟头、沟壁和沟底3处各布设3个空间重复采样点（图1），构成“发育阶段×沟道形态×发育部位×重复”的采样设计，共采集54个土壤样品。沟头点沿沟缘等距布设；沟壁与沟底点沿沟向自沟头至沟尾等间距设置，沟底点避开显著沉积区，沟壁点垂直沟壁采集。所有点位记录经纬度与高程。采样前清除植被、表层结皮与枯落物，每点使用1 540 cm³环刀采集原状土，用于测定土壤侵蚀阻力。

1.3 指标测定和数据处理

细沟可蚀性（k_d ）和临界剪切力（τ_c ）是表征土壤侵蚀阻力的2个关键力学参数，其测定对理解沟蚀起动机制尤为重要^［14-15］。本研究采用的射流侵蚀仪测定法，通过垂直射流模拟沟头溯源或沟壁冲刷的核心剥蚀机制，能够直接、标准化地获取土壤的侵蚀阻力。该方法在可控条件下剥离野外复杂因素的干扰，为机理研究提供可重复、可比较的基础参数，已在土壤侵蚀阻力表征及相关模型参数化中得到广泛应用^［16］。射流侵蚀仪装置见图2^［17］。该装置由供水单元（含三脚架、管路及溢流槽）与测定单元（含浸没槽、喷头、挡片及测量杆）组成。测定时，将原状土样固定于浸没槽中，在恒定水头射流持续冲刷下，记录不同时间对应的侵蚀深度，直至深度趋于稳定，获得平衡冲刷深度。基于冲刷过程的水动力-侵蚀数据序列，k_d 与τ_c 的计算通常采用3种数据分析方法：Blaisdell法^［18］、迭代法和冲刷深度法^［19］。鉴于迭代算法在数据拟合上优于Blaisdell法^{［16，19-20］}，本研究采用该方法进行分析。

使用环刀取样并测定土壤体积质量（BD）、总孔隙度（TP）；使用马尔文激光粒度仪测定粒度组成，并计算黏粒（Clay）、粉砂（Silt）与砂粒（Sand）含量；采用湿筛法测定>0.25 mm水稳性团聚体含量（WSA），并据此计算团聚体平均重量直径（MWD）；使用重铬酸钾外加热法测定土壤有机碳（SOC），并按需换算为土壤有机质（SOM）；利用烘干法获得根重密度（RWD）。

采用Microsoft Excel 2016软件进行数据整理与统计。采用单因素方差分析检验组间差异，显著性水平设置为0.05。采用随机森林量化土壤与根系属性对侵蚀阻力的贡献大小。采用线性混合效应模型分析发育阶段、沟道形态和发育部位的固定效应及其交互作用对土壤侵蚀阻力的影响，并通过计算方差解释比例量化各因素的独立与联合贡献。采用ArcGIS 10.8、Origin 2022和Microsoft PowerPoint 2016软件进行图形绘制。

2 结果与分析

2.1 不同发育阶段切沟的土壤理化性质和根系特征

由表2可知，不同切沟发育阶段下土壤理化性质及根系指标均存在显著差异（p<0.05）。土壤物理性质在切沟发育前期与中后期差异显著，而土壤化学性质与根系指标则在发育中期与发育前、后期差异显著。就物理性状而言，黏粒占比由发育前期的（9.53±2.44）%增至中期的（14.62±3.02）%，后期为（11.48±3.00）%；砂粒由（40.74±7.84）%降至中期的（32.57±5.38）%，后期回升至（35.75±8.67）%。土壤体积质量从发育前期的（1.30±0.08） g/cm³降低至发育后期的（1.20±0.11） g/cm³，且总孔隙度由（51.09±2.99）%增加至（54.83±4.00）%（表2），反映出发育后期切沟的土壤结构更为疏松。水稳性团聚体含量与团聚体平均重量直径在中期达到最大值［水稳性团聚体含量：前期（9.96±7.26）%、中期（23.93±11.06）%、后期（18.48±10.70）%；团聚体平均重量直径：前期（0.22±0.18） mm、中期（0.49±0.35） mm、后期（0.31±0.25） mm］，说明发育中后期团聚体结构相对前期更加完善。有机质质量分数和根重密度亦在中期达到峰值且发育后期高于前期［有机质：前期（5.70±2.80） g/kg，中期（8.24±5.09） g/kg，后期（7.11±3.56） g/kg；根重密度：前期（0.29±0.26） mg/cm³，中期（0.57±0.47） mg/cm³，后期（0.31±0.31） mg/cm³］，二者趋势一致，表明发育中后期植被生长相对旺盛，根系发育可能是土壤养分与有机质积累的重要驱动因子。

由表2可知，切沟沟头、沟壁与沟底在大多数土壤理化性质及根系指标上存在显著差异（p<0.05）。在土壤质地方面，沟壁的黏粒占比［（10.68±2.67）%］显著低于沟头［（12.39±3.67）%］与沟底［（12.56±3.92）%］，而粉粒与砂粒占比无显著差异，表明质地差异主要由黏粒占比变化驱动。沟壁土壤体积质量较高［（1.34±0.12） g/cm³］且总孔隙度较低［（49.39±4.31）%］，分别显著异于沟头［（1.24±0.09） g/cm³；（53.34±3.35）%］与沟底［（1.22±0.08） g/cm³；（54.00±3.01）%］，表明沟壁土壤更为紧密、孔隙连通性较差。沟壁水稳性团聚体含量也较低［（12.28±8.93）%］，表明沟壁的结构稳定性较弱。根重密度分布规律与表征土壤结构稳定性的水稳性团聚体含量和团聚体平均重量直径相似［沟头（0.65±0.38）mg/cm³、沟底（0.44±0.33） mg/cm³、沟壁（0.08±0.07） mg/cm³］，表明根系分布可能是影响切沟系统结构稳定的重要因素；根系稀疏的沟壁区域因此总体土壤结构稳定性较差，而沟头与沟底因根系较发达而表现出更强的稳定性与侵蚀阻力。

2.2 不同发育阶段切沟的土壤侵蚀阻力特征

由图3可知，不同发育阶段的沟道在k_d 和τ_c 上表现出显著差异。总体上，τ_c 随发育阶段呈上升趋势［发育前期（2.19±1.39） Pa；中期（3.11±1.92） Pa；后期（3.62±2.77） Pa）］，而k_d 则呈下降趋势［发育前期（320.46±424.21） cm³/（N·s）；中期（121.28±145.99） cm³/（N·s）；后期（51.20±36.93） cm³/（N·s）］。2项指标均存在较大组内变异，表明不同沟型之间可能存在较大的差异。具体而言，窄型沟的τ_c 在中期达到峰值后略有回落［前期（1.70±1.00） Pa；中期（2.46±1.61） Pa；后期（2.43±0.95） Pa］，而宽型沟τ_c 则随发育持续增加［前期（2.69±1.59） Pa；中期（3.76±2.08） Pa；后期（4.81±3.49） Pa］；相应地，窄型沟的k_d 从前期的（448.95±548.68） cm³/（N·s） 降至后期的（55.47±15.77） cm³/（N·s），降幅显著大于宽型沟［宽型沟由（191.96±210.23） cm³/（N·s）降至（46.93±51.07） cm³/（N·s）］。总体平均值显示宽型沟τ_c ［（2.99±2.58） Pa］显著高于窄型沟［（1.75±1.23） Pa］，而k_d ［（115.42±159.94） cm³/（N·s）］显著低于窄型沟［（237.5±359.74） cm³/（N·s）］，表明窄型沟总体上更易被侵蚀，侵蚀活跃度更高。

不同发育部位（沟头、沟壁、沟底）的τ_c 与k_d 变化规律与切沟整体变化趋势相对一致（图4）。平均τ_c 在沟头、沟壁、沟底分别为（4.59±2.21）、（1.34±0.61）、（2.98±1.88） Pa，平均k_d 分别为（44.81±25.79）、（380.06±406.89）、（68.07±59.25） cm³/（N·s）。具体而言，宽型沟沟头与沟底的τ_c ［沟头：（5.84±2.90） Pa；沟底：（3.93±2.42） Pa］显著高于窄型沟［沟头：（3.35±1.50） Pa；沟底：（2.03±1.00） Pa］，分别约为后者的1.74、1.93倍；窄型沟沟头和沟壁的k_d ［沟头：（60.29±30.72） cm³/（N·s）；沟壁：（495.88±359.74） cm³/（N·s）］显著高于宽型沟［沟头：（29.32±16.82） cm³/（N·s）；沟壁：（264.23±207.21） cm³/（N·s）］，分别约为其2.05、1.88倍。窄型沟各部位τ_c 通常在中期达峰后趋于稳定或略回落，而宽型沟τ_c 则持续上升；同时2类沟道各部位k_d 均随切沟发育呈显著下降趋势。总体来看，τ_c 表现为沟头>沟底>沟壁，k_d 表现为沟壁>沟底>沟头，表明沟壁拓宽和沟底下切可能是比沟头溯源更主要的侵蚀方式。

2.3 土壤与根系属性对侵蚀阻力的影响

随机森林模型分析结果表明，土壤与根系特征对τ_c 与k_d 具有较强的解释力，R²分别为58%和46%（图5）。所有因子均与k_d 显著相关；除黏粒和粉粒占比外，其余指标均与τ_c 显著相关。对于τ_c，前5个变量解释82.88%的相对重要性，其中有机质质量分数（%IncMSE =18.29）、根重密度（%IncMSE=12.21）和水稳性团聚体含量（%IncMSE=8.6）是影响其变化的关键预测因子，其相对重要性分别为28.47%、19.00%和13.39%，均呈极显著相关（p<0.001）。土壤体积质量、平均重量直径和总孔隙度等土壤理化性质也贡献较大，表明土壤结构特征对τ_c 影响较大。对于k_d，前5个变量解释77.26%的相对重要性，其中有机质质量分数（%IncMSE=14.30）、水稳性团聚体含量（%IncMSE=6.48）和总孔隙度（%IncMSE=5.30）是影响其变化的关键因子，其相对重要性分别为31.91%、14.47%和11.82%，且均达到显著相关水平（p<0.01）。根重密度、黏粒占比和土壤体积质量等根系指标与土壤理化性质也贡献较大，表明土壤根系、粒径和结构特征共同影响k_d。综合来看，根重密度和有机质质量分数等根系与养分指标对τ_c 重要性更高，而k_d 则更偏向受土壤团聚体和孔隙度等指标影响，表明土壤抗剪强度更依赖于根系固结作用，而可蚀性更易受土壤结构稳定性的影响。

2.4 土壤侵蚀阻力对切沟侵蚀演变的响应

线性混合效应模型表明，发育部位、发育阶段与沟道形态共同解释k_d 变异的52.2%（R²=0.522）和τ_c 变异的58.2%（R²=0.582）（图6）。其中，发育部位为主导性因子，独立解释k_d 变异的28.9%和τ_c 变异的37.0%；发育阶段与沟道形态的贡献次之，但对τ_c 而言，沟道形态的影响显著增强（独立解释13.0%）。凸显切沟演化中不同发育部位土壤侵蚀阻力变化的显著差异。边际效应分析显示，随着发育推进，k_d 由前期的320.5显著下降至后期的51.2，累计降幅84.0%（表3）；在窄型沟中该趋势更明显（由前期的449.0降至后期的55.5，降幅87.6%）。相反，τ_c 由发育前期的2.19显著上升至后期的3.62（增幅65.3%），在宽型沟中增幅更大（由发育前期的2.69至后期的4.81，增幅78.8%）。发育阶段与地貌部位存在显著交互作用（p<0.001），表明不同发育部位对发育过程的响应模式各异（表4）。沟头呈缓慢且渐进的k_d 下降，沟底在发育前期出现快速下降后保持稳定，而沟壁部位k_d 值从792.0急剧下降至71.7，降幅达91.0%，成为发育效应体现最为突出的地貌单元。各部位的τ_c 均随发育阶段稳步增加（增幅57.3%~85.1%）。综合来看，黄土丘陵区切沟发育演化与土壤侵蚀阻力显著相关：k_d 和τ_c 随发育阶段推进而系统性降低和增加，且此过程在不同地貌部位表现出差异化响应，并受到切沟形态的调节。

3 讨论

τ_c 和k_d 作为土壤侵蚀模型中的关键参数，分别表征侵蚀启动所需的临界水动力应力和土壤对侵蚀作用的敏感性与易蚀性^［21-22］。传统上，τ_c 和k_d 常通过Parker公式和WEPP模型等根据土壤物理属性（如粒径分布、密度、比重等）或者经验关系估算^［23］；而JET仪器通过受控喷射水流直接作用于土体表面，实测侵蚀速率与剪切应力间的关系，从而更真实地反映土壤属性与植被作用对土壤侵蚀阻力的影响^［16］。本研究中τ_c 和k_d 随切沟发育在整体发育阶段和各发育部位均分别呈现增加和降低的系统性变化，反映土壤侵蚀阻力在空间与时间上的系统性增强。从过程机制来看，发育中后期土壤黏粒与有机质累积、土壤体积质量下降、孔隙度与团聚体稳定性增加，以及根重密度的峰值出现，共同构成侵蚀阻力增强的物质与结构基础（表2）。其中，有机质为团聚体形成和孔隙结构优化提供黏结与胶结介质，根系则通过缠结、加筋和根系分泌物促团聚作用提高土体结构完整性和抗拉强度；二者协同作用使土体在同等水动力作用下更难侵蚀和剥离^{［10，24-25］}。随机森林模型进一步表明，有机质、根重密度及水稳性团聚体是调控τ_c 和k_d 的最关键因子，其相对重要性显著高于其他土壤理化指标，在统计层面印证上述以有机质（化学胶结）、根系（生物加筋与分泌）和水稳性团聚体（物理结构稳定性）为核心的生物-物理-化学过程的综合效应。因此，发育阶段对τ_c 和k_d 的影响，本质上是通过土壤结构重建和植被根系发育提高土壤侵蚀阻力。

土壤侵蚀阻力随发育的增强过程，在不同地貌部位表现出显著的差异化响应，体现局部侵蚀动力与物源属性的交互作用。其中，沟壁是响应最为突出的单元。沟壁在发育早期受边坡不稳定性及岸脚集中径流侵蚀影响，土体松散、结构破碎，植被难以定殖，表现出极低的侵蚀阻力；随着植被恢复、成壤作用及结构重组的推进，有机质累积、团聚体稳定性增强，根重密度显著提高，发育后期侵蚀阻力显著增强^［1］。沟底的侵蚀阻力则呈现前期快速上升、后期趋稳的特征，反映出初期切沟在强烈水流冲击下经历松散物质的迅速剥离与水力分选，随后床面粒级与形态趋于稳态，沉积与再分配作用减弱^［26］。沟头变化相对平缓，可能因坡面植被抑制溯源侵蚀，使径流冲刷作用有限，土壤颗粒分布较均匀、质地较好^［21］。总体上，地貌部位的差异通过塑造各部位特定的水动力作用方式和物源供给条件，进而影响土壤、根系属性的演化路径，最终形成“沟头侵蚀阻力最高、沟底次之、沟壁最低”的空间格局。此结果在趋势上与王亚娟等^［27］提出的“沟床>沟岸>沟头”可蚀性排序一致，但在沟壁数值特征上仍存在差异，可能与本研究采样位置取自垂直裸露面有关，表明结论对采样位置较为敏感，同时揭示沟壁作为潜在高风险侵蚀单元的重要性。

沟道形态是调节上述发育效应的另一关键因素。窄型沟因流量集中，单位面积水流速度和剪切力显著高于坡面流，长期维持较强的水力扰动，相对抑制土壤结构的稳定化^［21］，从而降低侵蚀触发的临界应力τ_c，使其在发育过程中保持较高的易蚀性。相比之下，k_d 主要受土壤理化性质（如质地、有机质含量等）支配^［27］，因此受形态控制较弱，但其下降过程在窄型沟中更为显著，可能与窄沟内更快速的物源选择性剥离有关。总体而言，沟道形态并非直接改变τ_c 和k_d，而是通过调节局地水动力场及其对土壤结构与植被根系恢复的制约，间接影响土壤侵蚀阻力的演化速率和空间格局，即初始土壤侵蚀阻力和沟道形态共同决定侵蚀热点部位与形态演化路径，而随发育推进，沟道形态的反馈又反过来通过水动力、土壤和植被条件来重塑侵蚀阻力的时空分布。

综上所述，黄土丘陵区切沟侵蚀发育可视为一个在发育阶段、地貌形态及土壤与根系属性等多因子耦合控制下演化的动力学过程。发育阶段、地貌部位和沟道形态通过塑造局地水动力与物源条件，控制土壤与根系属性的演化路径和空间分异；而土壤与根系属性则通过调节τ_c 和k_d，决定局部土壤侵蚀阻力水平，并通过影响侵蚀速率和物质重分配，在更长时间尺度上反过来作用于切沟的进一步发育与形态调整。从治理实践而言，应重点关注发育初期及窄型沟道的侵蚀控制，并对沟壁等高风险单元实施针对性固土与植被恢复措施，以顺应自然演替规律，阻断侵蚀恶性循环，促进土壤结构与生态功能的协同恢复。

4 结论

1） 切沟发育过程中土壤侵蚀阻力整体呈增强趋势。τ_c 随切沟发育显著上升（增幅65.3%），而k_d 显著下降（降幅84.0%），表明土壤侵蚀阻力随切沟发育而强化，系统趋于稳定。

2） 土壤侵蚀阻力变化由关键土壤与根系属性驱动。随机森林模型表明，有机质是影响τ_c 和k_d 的共同关键因子；τ_c 更依赖于根系生物量，而k_d则更受团聚体稳定性控制。

3） 地貌部位和沟道形态显著调节土壤侵蚀阻力的空间分异。侵蚀阻力表现为沟头最强、沟底次之、沟壁最弱的格局，且沟壁对发育阶段响应最为敏感，表明沟壁侧蚀与沟底下切是主要的活跃侵蚀方式。窄型沟整体土壤侵蚀阻力（τ_c 更低，k_d 更高）弱于宽型沟，且其τ_c 增幅与k_d 降幅更显著，表明初始形态通过改变水动力与物质输移条件影响土壤侵蚀阻力重建的速率与路径。

4） 切沟侵蚀演变受多因子协同控制。地貌部位、发育阶段和沟道形态共同解释k_d 变异的52.2%（R²=0.522）和τ_c 变异的58.2%（R²=0.582）。其中，发育阶段与地貌部位存在极显著交互作用（p<0.001），表明k_d 和τ_c 随切沟发育在不同地貌部位表现出差异化响应，并受到切沟形态的协同调节。

参考文献

[1] HAYAS A， GÓMEZ J A. Modelling of gully widening， a review. Implications for research on gully evolution and restoration［J］.Earth-Science Reviews，2024，255：e104836.

[2] 刘宝元，杨扬，陆绍娟.几个常用土壤侵蚀术语辨析及其生产实践意义［J］.中国水土保持科学，2018，16（1）：9-16.

[3] LIU B Y， YANG Y， LU S J. Discriminations on common soil erosion terms and their implications for soil and water conservation［J］.Science of Soil and Water Conservation，2018，16（1）：9-16.

[4] CHEN R， ZHOU Y， WANG Z T， et al. Towards accurate mapping of loess waterworn gully by integrating google earth imagery and DEM using deep learning［J］.International Soil and Water Conservation Research，2024，12（1）：13-28.

[5] CHEN Y X， JIAO J Y， YAN X Q， et al. Response of gully morphology and density to the spatial and rainy-season monthly variation of rainfall at the regional scale of the Chinese Loess Plateau［J］.Catena，2024，236：e107773.

[6] 张光辉.切沟侵蚀研究进展与展望［J］.水土保持学报，2020，34（5）：1-13.

[7] ZHANG G H. Advances and prospects for gully erosion researches［J］.Journal of Soil and Water Conservation，2020，34（5）：1-13.

[8] 宋爽，范昊明，牛天一.东北黑土区不同季节侵蚀沟形态发育比较分析［J］.水土保持学报，2022，36（1）：18-23.

[9] SONG S， FAN H M， NIU T Y. Comparative analysis of the development the geometric forms of gully in different seasons in black soil area of northeast China［J］.Journal of Soil and Water Conservation，2022，36（1）：18-23.

[10] YUAN X Y， YE F， FU W X， et al. Estimating the critical shear stress for incipient particle motion of a cohesive soil slope［J］.Scientific Reports，2022，12：e9736.

[11] THOMAS J T， IVERSON N R， BURKART M R， et al. Long-term growth of a valley-bottom gully，western Iowa［J］.Earth Surface Processes and Landforms，2004，29（8）：995-1009.

[12] POESEN J， NACHTERGAELE J， VERSTRAETEN G， et al. Gully erosion and environmental change： Importance and research needs［J］.Catena，2003，50（2/3/4）：91-133.

[13] ROBERTS M E， BURROWS R M， THWAITES R N， et al. Modelling classical gullies： A review［J］.Geomorphology，2022，407：e108216.

[14] 王义霞.黄土丘陵区水力侵蚀和植被恢复对土壤有机碳输移及分布的影响［D］.陕西 杨凌：中国科学院大学（中国科学院教育部水土保持与生态环境研究中心），2018.

[15] WANG Y X. Effects of water erosion and vegetation rehabilitation on the soil organic carbon transportation and distribution in small catchments of the hilly loess region［D］.Yangling， Shaanxi： University of Chinese Academy of Sciences （Research Center of Soil and Water Conservation and Environment， Chinese Academy of Science and Ministry of Education），2018.

[16] THWAITES R N， BROOKS A P， PIETSCH T J， et al. What type of gully is that the need for a classification of gullies［J］.Earth Surface Processes and Landforms，2022，47（1）：109-128.

[17] 李建军，焦菊英，陈一先，等.青藏高原侵蚀沟分类初探［J］.水土保持学报，2025，39（1）：1-9.

[18] LI J J， JIAO J Y， CHEN Y X， et al. Elementary study on gully classification on the Qinghai-Xizang Plateau［J］.Journal of Soil and Water Conservation，2025，39（1）：1-9.

[19] 卢超，马周加态，李佳辉，等.冻融条件下土壤侵蚀阻力影响因素［J］.水土保持学报，2023，37（6）：25-33.

[20] LU C， MA Z J T， LI J H， et al. The factors influencing soil erosion resistance under freeze-thaw conditions［J］.Journal of Soil and Water Conservation，2023，37（6）：25-33.

[21] LEE S， CHU M L， GUZMAN J A， et al. Modeling soil erodibility and critical shear stress parameters for soil loss estimation［J］.Soil and Tillage Research，2022，218：e105292.

[22] KARAMIGOLBAGHI M， GHANEEIZAD S M， ATKINSON J F， et al. Critical assessment of jet erosion test methodologies for cohesive soil and sediment［J］.Geomorphology，2017，295：529-536.

[23] 郭金金，方怒放，刘前进，等.小麦秸秆粉碎还田配施氮肥与多酚对土壤侵蚀阻力的影响［J］.水土保持学报，2021，35（4）：122-128.

[24] GUO J J， FANG N F， LIU Q J， et al. Effects of incorporated wheat straw returning combined application with nitrogen and polyphenols on soil erosion resistance［J］.Journal of Soil and Water Conservation，2021，35（4）：122-128.

[25] BLAISDELL F W， HEBAUS G G， ANDERSON C L. Ultimate dimensions of local scour［J］.Journal of the Hydraulics Division，1981，107（3）：327-337.

[26] DALY E R， FOX G A， MILLER R B， et al. A scour depth approach for deriving erodibility parameters from jet erosion tests［J］.Transactions of the ASABE，2013，56（6）：1343-1351.

[27] KHANAL A， FOX G A， AL-MADHHACHI A T. Variability of erodibility parameters from laboratory mini jet erosion tests［J］.Journal of Hydrologic Engineering，2016，21（10）：e04016030.

[28] KNAPEN A， POESEN J， GOVERS G， et al. Resistance of soils to concentrated flow erosion： A review［J］.Earth-Science Reviews，2007，80（1-2）：75-109.

[29] SINGH H V， THOMPSON A M. Effect of antecedent soil moisture content on soil critical shear stress in agricultural watersheds［J］.Geoderma，2016，262：165-173.

[30] XING H， HUANG Y H， CHEN X Y， et al. Comparative study of soil erodibility and critical shear stress between loess and purple soils［J］.Journal of Hydrology，2018，558：625-631.

[31] WANG C， CAI C F， DENG Y S. Tiny niche terrain induces gully headcut retreat［J］.Earth Surface Processes and Landforms，2024，49（7）：2260-2277.

[32] WU J， YANG G， MA Y X， et al. Effects of vegetation restoration on soil aggregate characteristics and soil erodibility at gully head in loess hilly and gully region［J］.Scientific Reports，2024，14（1）：e31149.

[33] BARTLEY R， POESEN J， WILKINSON S， et al. A review of the magnitude and response times for sediment yield reductions following the rehabilitation of gullied landscapes［J］.Earth Surface Processes and Landforms，2020，45（13）：3250-3279.

[34] 王亚娟，陈云明，孙亚荣，等.宁南山区侵蚀沟不同部位土壤理化性质及可蚀性研究［J］.水土保持学报，2023，37（2）：11-18.

[35] WANG Y J， CHEN Y M， SUN Y R， et al. Study on physical and chemical properties and erodibility of soil in different parts of erosion gully in the mountainous area of south Ningxi［J］.Journal of Soil and Water Conservation，2023，37（2）：11-18.

基金资助

国家自然科学基金项目(42177335)

陕西省自然科学基础研究计划项目(2024ZY-JCYJ-02-18)