1982—2020年黄土高原水蚀风蚀交错带土壤风蚀的时空变化及其主控因素

doi:10.7522/j.issn.1000-694X.2026.00025

1982—2020年黄土高原水蚀风蚀交错带土壤风蚀的时空变化及其主控因素

严云籍^,¹, 白飞², 张加琼^,¹^,³, 安鑫¹, 范晓强¹, 王椰¹^,³

1.西北农林科技大学水土保持科学与工程学院/水土保持与荒漠化整治全国重点实验室，陕西杨凌 712100

2.榆林市水土保持生态工程建设中心，陕西榆林 719000

3.中国科学院水利部水土保持研究所，陕西杨凌 712100

Spatiotemporal variations in wind erosion intensity and their main impact factors in the wind-water erosion crisscross region of the Loess Plateau from 1982 to 2020

Yan Yunji^,¹, Bai Fei², Zhang Jiaqiong^,¹^,³, An Xin¹, Fan Xiaoqiang¹, Wang Ye¹^,³

1.College of Soil and Water Conservation Science and Engineering / State Key Laboratory of Soil and Water Conservation and Desertification Control，Northwest A&F University，Yangling 712100，Shaanxi，China

2.Yulin Engineering Construction Center of Soil and Water Conservation Ecological，Yulin 719000，Shaanxi，China

3.Institute of Soil and Water Conservation，Chinese Academy of Sciences and Ministry of Water Resources，Yangling 712100，Shaanxi，China

通讯作者: 张加琼（E-mail： jqzhang@nwsuaf.edu.cn）

收稿日期: 2025-12-31 修回日期: 2026-02-11

基金资助:

国家自然科学基金联合基金项目. U2243210
国家自然科学基金面上项目. 42577396

Received: 2025-12-31 Revised: 2026-02-11

作者简介 About authors

严云籍（2000—），女，云南丘北人，硕士研究生，主要从事土壤侵蚀研究E-mail：2019011610@nwsuaf.edu.cn , E-mail：2019011610@nwsuaf.edu.cn

摘要

在退耕还林（草）等生态工程持续实施背景下，黄土高原水蚀风蚀交错带土壤风蚀强度与格局均明显改变。然而，目前缺乏该区域土壤风蚀强度和格局的系统评估，难以支撑变化环境下的风蚀精准防控。本文采用RWEQ估算1982—2020年的逐年风蚀模数，解析其时空演变及主控因素。结果表明：1982年以来区域年均风蚀模数可分为波动下降（1982—1989）、相对稳定（1990—2000）、快速下降（2001—2008）和低值稳定（2009—2020）4个阶段，阶段四的年均风蚀模数（681.8 t·km^-2·a^-1）较阶段一降低73.3%。各阶段春季均为风蚀高发期，风蚀量平均占全年的53.5%。近40年来，风蚀强度显著减小区域面积占比超过80%，尤其是中部地区风蚀强度显著减小的面积达该地区总面积的98.1%。此外，聚集分布于西北部和北部地区的中度及以上和强烈及以上风蚀区明显减少，二者在全区域的面积占比分别从阶段一的22.4%和12.9%降至阶段四的6.7%和2.8%。近40年来，土壤风蚀强度空间分异主要受土壤性质、气候变化与植被覆盖共同影响，退耕还林（草）等生态工程的实施使区域植被覆盖度平均提高16.5%，其对风蚀强度降低的贡献从29.6%升高到62.1%，植被建设的风蚀防控效应显著。

关键词： 水蚀风蚀交错带 ; RWEQ ; 风蚀强度 ; 时空变化 ; 影响因素

Abstract

The temporal and spatial patterns of soil erosion intensity in the wind-water erosion crisscross region of the Loess Plateau have obviously changed with the continuous implementation of ecological restoration projects represented by Grain-for-Green Program. However， understanding of the spatiotemporal variation in wind erosion intensity remains limited， hindering the development of precise and targeted control strategies under changing environmental conditions. This study applied the Revised Wind Erosion Equation （RWEQ） to estimate annual wind erosion modulus from 1982 to 2020， to clarify spatiotemporal variation in wind erosion intensity， and to determine its main impact factors. The results showed that annual wind erosion modulus exhibited four distinct stages since 1982， which were fluctuating decline （1982-1989）， relative stability （1990-2000）， rapid decline （2001-2008）， and low-value stability （2009-2020） stages. The mean modulus in the fourth stage （681.8 t·km^-2·a^-1） reduced by 73.3% compared to that of the first stage. Spring was consistently the high-risk season of wind erosion in all phases， accounting for 53.5% of annual total erosion modulus. In the past four decades， areas with significant reduction of wind erosion intensity exceeded 80% of the study region， particularly in the central area， where 98.1% of this area showed significant reduction in wind erosion intensity. The area of above-moderate and above-severe erosion concentrated in the northwest and north decreased from 22.4% and 12.9% in the first stage to 6.7% and 2.8% in the fourth stage， respectively. The spatial variation in wind erosion intensity was jointly influenced by soil properties， climate， and vegetation cover. However， the impact of vegetation cover attributed to ecosystem restoration projects such as the Grain-for-Green Program increased from 29.6% to 62.1% with mean vegetation cover increasing by 16.5% in the past four decades， which indicated substantial environmental improvements induced by ecological restoration.

Keywords： wind-water erosion crisscross region ; Revised Wind Erosion Equation （RWEQ） ; wind erosion intensity ; spatiotemporal variation ; impact factors

PDF (4479KB) 元数据多维度评价相关文章导出 EndNote| Ris| Bibtex 收藏本文

本文引用格式

严云籍, 白飞, 张加琼, 安鑫, 范晓强, 王椰. 1982—2020年黄土高原水蚀风蚀交错带土壤风蚀的时空变化及其主控因素. 中国沙漠[J], 2026, 46(2): 336-347 doi:10.7522/j.issn.1000-694X.2026.00025

Yan Yunji, Bai Fei, Zhang Jiaqiong, An Xin, Fan Xiaoqiang, Wang Ye. Spatiotemporal variations in wind erosion intensity and their main impact factors in the wind-water erosion crisscross region of the Loess Plateau from 1982 to 2020. Journal of Desert Research[J], 2026, 46(2): 336-347 doi:10.7522/j.issn.1000-694X.2026.00025

0 引言

黄土高原水蚀风蚀交错带是风力与水力在时空尺度上交替并叠加作用的典型区域，生态环境脆弱、土壤侵蚀过程复杂^［1-2］。该区地表物质组成复杂且质地疏松，植被覆盖度较低，对气候变化和人类活动扰动高度敏感。随着退耕还林（草）、三北防护林等生态工程的持续实施，该区域生态环境明显改善，侵蚀环境和土壤侵蚀空间分布格局明显改变^［3-5］。典型区域的植被盖度从20世纪70—80年代的约38%增加至2020年的52%，多年平均侵蚀模数从3 499 t·km^-2·a^-1降低至526 t·km^-2·a^-1［6］。然而，水蚀风蚀交错带仍是黄土高原土壤侵蚀最强烈的地区^［7］，是黄河下游河床粗泥沙的主要来源地^［8］，更是未来黄土高原生态环境建设的重点^［9］。因此，阐明黄土高原水蚀风蚀交错带土壤侵蚀强度的时空变化特征，对精准制定科学的水土流失治理措施体系具有重要意义。

目前，黄土高原水蚀风蚀交错带土壤侵蚀时空变化特征主要聚焦水蚀过程^［10-12］，缺乏对风蚀时空变化的系统刻画^［13］。针对该区域内毛乌素沙地（2000—2022年）、宁夏（1988—2018年）、无定河流域（2000—2018年）、北方农牧交错带（2001—2012年）等典型地区的研究显示，风蚀模数整体呈波动下降趋势，但在2003、2010年前后出现增大，且不同研究区的增大幅度差异明显^［14-17］。如宁夏草地禁牧（2003年）后平均风蚀模数由5 590 t·km^-2·a^-1（1988—2003年）降低至1 450 t·km^-2·a^-1（2003—2018年），波动上升期（2003—2010年）的平均风蚀模数（1 783 t·km^-2·a^-1）仍明显低于禁牧前。而在无定河流域及毛乌素沙地边缘，风蚀模数2010年增大到2 287 t·km^-2·a^-1，比2005—2009年提高了74%。在空间分布上，西北和北部风蚀较强，东南或南部较弱；且在退耕还林（草）等生态工程作用下^{［11，18-19］}，不同地区风蚀强度均呈现从局部减弱到总体持续减弱的变化，强度从强烈及以上区域连续分布逐渐转变为中度及以下为主^［20-21］，但西北和北部风蚀模数年际波动明显（540~2 650 t·km^-2·a^-1）且区域差异大（640~3 020 t·km^-2·a^-1）^［22-23］。这与风蚀严重区脆弱的环境对气候、土壤、植被等变化敏感密切相关^［24-25］。可见，前人多聚焦黄土高原水蚀风蚀交错带的典型、局部区域，且更关注2000年以后的风蚀变化，未明晰区域土壤风蚀时空变化特征，不利于识别风蚀强度及环境条件巨变的区域并明晰主控因素，难以为区域风蚀精准防控提供可靠的科学依据。因此，亟须聚焦于黄土高原水蚀风蚀交错带，系统揭示土壤风蚀强度的长时序演变特征、空间分异格局及其影响因素。

据此，本文针对变化环境背景下黄土高原水蚀风蚀交错带土壤风蚀时空变化特征认知不足的问题，融合多源气象、遥感与地理数据，选用修正土壤风蚀方程（RWEQ，Revised Wind Erosion Equation）估算黄土高原水蚀风蚀交错带1982—2020年土壤风蚀模数，明晰近40年来风蚀强度的时间演变特征及空间格局变化，识别风蚀变化的主控因子及其交互作用，为黄土高原水蚀风蚀交错带风蚀防控与生态治理提供科学依据。

1 材料与方法

1.1　研究区概况

黄土高原水蚀风蚀交错带（35°01′—40°44′N，100°52′—114°31′E）为东北-西南走向的带状区域，面积约为19.3万 km²，是典型的农牧交错带和黄河粗泥沙的主要来源区^［26-27］（图1）。该区属温带半干旱大陆性季风气候，温度、降雨等年际和年内变异大。该区域年均气温7~9 ℃，≥10 ℃年积温2 200~3 500 ℃，无霜期140~170 d。水蚀和风蚀全年交错发生，多年平均降水量250~450 mm，年际变幅可达2~4倍；70%~80%的降水集中于6—9月，且多以短历时强降雨的形式出现。风季主要在10月到次年5月，年均风速约2.2 m·s^-1，春季多大风沙暴，大风日数15~87 d，沙暴日数多在4 d以上^［28］。地貌类型以黄土丘陵沟壑、土石丘陵、风沙丘陵和山地为主^［29］，沟道密度3~6 km·km^-2，地表破碎度25%~50%^［26］。土壤类型主要为质地疏松、抗蚀性差的黄绵土、沙黄土、灰钙土和风沙土等。该区植被以草原为主，灌丛与人工林主要呈斑块状分布。随着退耕还林（草）、三北防护林等生态工程的持续实施，该区域植被覆盖度明显提升（提升15%~30%），生态环境明显改善^［30］。

图1

新窗口打开| 下载原图ZIP| 生成PPT

图1 黄土高原水蚀风蚀交错带地理位置及用于模拟的气象站点分布

注：基于自然资源部标准地图服务网站审图号GS（2023）2763号标准地图制作，底图边界无修改

Fig.1 Location of the wind-water erosion crisscross region of the Loess Plateau and meteorological stations used for simulation

1.2　数据来源与处理

研究时段内（1982—2020年），用于计算地表粗糙度因子的DEM数据来源于地球观测系统数据和信息系统的全球数字高程模型（https：//search.earthdata.nasa.gov），空间分辨率为30 m。用于计算天气因子的风速、降水、气温等气象数据来源于中国气象科学数据共享服务网（https：//data.cma.cn/metadata）的逐日站点观测资料。共包括研究区内的68个及其周边的125个站点（图1）。降水、温度和风速均基于站点气象数据插值。降水和温度插值时，考虑其与高程的相关性，采用考虑DEM（重采样到1 km）的协同克里金（Co-kriging）插值^［31］；风速插值时，鉴于其空间连续性较弱的特征，采用普通克里金（Ordinary Kriging）插值以降低过拟合误差^［32］。插值后的降水、气温、风速数据空间分辨率均为1 km。用于计算雪盖因子和土壤湿度因子的逐月雪深数据与地表潜在蒸散发数据来源于国家青藏高原科学数据中心（https：//data.tpdc.ac.cn/），空间分辨率分别为0.25°和1 km。国家青藏高原科学数据中心发布的中国区域250 m分辨率的NDVI数据（2000—2020年）未能完整覆盖本文的研究时段，本文将其与国家地球系统科学数据中心（https：//www.geodata.cn/）发布的中国8 km分辨率逐月NDVI数据集（1982—2020年）融合，构建连续时间序列的NDVI数据，并降尺度至1 km。用于计算土壤可蚀性因子与土壤结皮因子的土壤数据来源于世界土壤数据库（HWSD v2.0，https：//gaez.fao.org/pages/hwsd），空间分辨率为1 km。这些空间分辨率不一致的数据均使用最邻近采样统一到1 km空间分辨率。

1.3　研究方法

1.3.1　风蚀模数估算

采用修正土壤风蚀方程（RWEQ）估算研究区域的土壤风蚀模数。RWEQ基于气候条件、土壤可蚀性、地表粗糙度、植被覆盖、结皮盖度等要素估算风蚀模数^［33］：

S L = \frac{2 x}{s^{2}} Q_{m a x} e^{- {(x / s)}^{2}}

（1）

s = 150.71 {(W F \times E F \times S C F \times K' \times C O G)}^{- 0.3711}

（2）

Q_{m a x} = 109.8 (W F \times E F \times S C F \times K' \times C O G)

（3）

式中：SL为地块风蚀模数（t·km^-2·a^-1）；Q_max为最大输沙能力（kg·m^-1）；x为距上风向不可蚀边界的地块长度（取50 m）；s为关键地块长度（m）。

天气因子

W F = \frac{\sum_{i = 1}^{n} U_{i} (U_{i} - U_{t})^{2} \times N_{d} \times ρ}{N \times g} \cdot S W \cdot S D

（4）

U_{i} = U_{10} {(\frac{2}{10})}^{\frac{1}{7}}

（5）

式中：WF为气候因子（kg·m^-1）；U_i 和U₁₀分别为近地面2 m和10 m高度处风速（m·s^-1）；U_t为2 m高度的临界风速，本文取值5 m·s^-1；N为风速观测次数，通常为500次，由于气象站点日均风速数据无法满足该要求，采用郭中领等^［34］的方法进行降尺度以提高风速数据的时间分辨率；N_d 为实验天数，本文中风速估算在半月尺度进行，取值为15天；ρ为空气密度（kg·m^-3）；g为重力加速度，取值9.8 m·s^-2；SW为土壤湿度因子（无量纲）；SD为积雪盖度因子（无量纲）。

土壤可蚀性与结皮因子

E F = \frac{29.09 + 0.31 S a + 0.17 S i + 0.33 S a / C l - 2.59 O M - 0.95 C a C O_{3}}{100}

（6）

S C F = \frac{1}{1 + 0.0066 {(C l)}^{2} + 0.021 {(O M)}^{2}}

（7）

式中：EF为土壤可蚀性因子（无量纲）；SCF为土壤结皮因子（无量纲）；Cl、Si和Sa分别为黏粒、粉粒和砂粒含量（%）；OM为有机质含量（%）。

植被覆盖度因子。RWEQ依据作物的直立残茬、平铺残茬、作物覆盖估算该因子^［33］。为适应区域尺度的计算，选用作物覆盖估算方法计算植被覆盖度因子^［35］：

C O G = e^{- 5.614 C C^{0.7366}}

（8）

C C = \frac{N D V I - N D V I_{s o i l}}{N D V I_{v e g} - N D V I_{s o i l}}

（9）

式中：CC为植被覆盖度（%）；NDVI_soil和NDVI_veg分别为裸地和植被覆盖区域的NDVI值，在本文中分别选择累积的5%和95%为置信区间所对应的NDVI值代替^［36］。

地表粗糙度因子 $K'$ 。RWEQ使用随机粗糙度和垄作形成的定向粗糙度表示该因子。然而，在非农田、地形存在明显变化等条件下，原模型中该因子的适用性不足^［37］。因此，基于坡度（α，°）的变化计算该因子^［38］：

K' = c o s α

（10）

1.3.2　基于逐像元一元线性回归的NDVI数据空间分辨率的统一

由于覆盖本文全部时段的NDVI数据空间分辨率较低（8 km），基于从国家青藏高原科学数据中心发布的中国区域250 m分辨率的NDVI数据，选取两个数据集重叠的时段（2000—2020），采用逐像元一元线性回归分析构建重采样算法。其中，2000—2015年的数据用于构建一元线性回归模型算法，得到1982—2020年的GIMMS-MODIS NDVI数据集，2016—2020年的数据用于验证。结果表明GIMMS-MODIS NDVI数据集与MODIS数据的相关性较高（R²=0.82），与前人的研究结果一致^［39-40］。逐像元构建一元线性回归模型如下：

{\hat{M}}_{i} = a + b G_{i} + ε_{i}

（11）

b = \frac{\sum_{i = 1}^{n} (M_{i} - \bar{M}) (G_{i} - \bar{G})}{\sum_{i = 1}^{n} (G_{i} {- \bar{G})}^{2}}

（12）

a = \bar{M} - a \bar{G}

（13）

式中：a，b分别代表不同数据集中的数据经过线性拟合后的截距和斜率；ε_i 代表可能引起的误差；G_i 为第i个月GIMMS NDVI；M_i 为第i个月MODIS NDVI； $\bar{G}$ 为GIMMS NDVI月平均值； $\bar{M}$ 为MODIS NDVI月平均值； ${\hat{M}}_{i}$ 为构建的第i个月GIMMS-MODIS NDVI。

1.3.3　土壤风蚀模数时空变化特征

风蚀模数的时间变化特征通过其在研究期内的变化趋势及其变化的显著性表征。首先，采用Theil-Sen Median趋势分析阐明研究时段内风蚀模数的变化趋势^［41］。然后采用Mann-Kendall显著性检验方法对Theil-Sen Median趋势分析结果进行非参数显著性检验^［42］。

S l o p e = M e d i a n (\frac{x_{j} - x_{i}}{j - i}), \forall j > i

（14）

Z_{M K} = \{\begin{matrix} \frac{S - 1}{\sqrt[]{V a r (S)}}, & S > 0 \\ 0, & S = 0 \\ \frac{S - 1}{\sqrt[]{V a r (S)}}, & S < 0 \end{matrix}

（15）

式中：Median为中值；x_i 和x_j 分别为第i年和第j年的风蚀模数；Slope为趋势斜率，Slope>0、Slope<0和Slope=0分别表示风蚀模数呈上升、下降趋势和无明显变化；Z_MK为标准化的趋势检验统计量；S为趋势检验统计量；Var为方差；在α=0.05时，|Z_MK|>1.96表明风蚀模数具有显著变化（P<0.05）。依据Theil-Sen Median趋势分析和Mann-Kendall显著性检验结果将风蚀模数趋势变化划分为5个等级：显著增加（Slope>0，|Z_MK|>1.96），不显著增加（Slope>0，|Z_MK|≤1.96），无变化（Slope=0，|Z_MK|≤1.96），显著减少（Slope<0，|Z_MK|>1.96），不显著减少（Slope<0，|Z_MK|≤1.96）。

依据《土壤侵蚀强度分级标准》将研究区土壤风蚀强度划分为6个等级：微度（<200 t·km^-2·a^-1）、轻度（200~2 500 t·km^-2·a^-1）、中度（2 500~5 000 t·km^-2·a^-1）、强烈（5 000~8 000 t·km^-2·a^-1）、极强烈（8 000~15 000 t·km^-2·a^-1）和剧烈（>15 000 t·km^-2·a^-1）^［43］。结合风蚀强度的变化程度，将其划分为稳定不变、轻度改善、中度改善、大幅改善、轻微恶化、中度恶化、大幅恶化7个等级。稳定不变表示风蚀强度无变化，轻度至大幅改善表示风蚀强度等级分别降低1、2~3和4~5个等级；反之轻度恶化到大幅恶化分别表示风蚀强度等级增强1、2~3和4~5个等级。依据不同时期风蚀强度的变化特征，采用风蚀强度转移矩阵表示侵蚀强度变化区域及其面积：

E = (\begin{matrix} E_{11} & \dots & E_{1 i} \\ ⋮ & ⋮ & ⋮ \\ E_{j 1} & \dots & E_{i j} \end{matrix})

（16）

式中：E为转移矩阵；E_ij 为t时期i级土壤侵蚀强度转变为t+T时期j级土壤侵蚀强度的面积。

1.3.4　风蚀模数变化的主控因素解析

采用1.3.3中的时间变化特征分析方法刻画主控因素的变化趋势，采用地理探测器的因子探测器与交互探测器模块量化各因子及其交互作用对风蚀模数的影响大小。基于RWEQ中的主要因子，将风蚀模数变化的影响因子分为侵蚀力（年均风速X₁、年降水量X₂）、植被覆盖度X₃、土壤性质（黏粒含量X₄、有机碳含量X₅）、土地利用类型X₆、地形（用坡度X₇表示）相关因子。对于主要变化因素表征上述因子对风蚀模数变化解释力的核心指标q值通过以下公式计算：

q = 1 - \frac{\sum_{h = 1}^{L} N_{h} σ_{h}^{2}}{N σ^{2}}

（17）

式中：q值取值为［0，1］，越接近于1，表明影响程度越强；h为因子的分类序号；L为土壤风蚀模数或影响因子的分层总数；N_h 和N分别为进行风蚀模数估算的每个分区单元数和全区单元数； $σ_{h}^{2}$ 和 $σ^{2}$ 分别为每个分区和全区土壤风蚀模数的方差。

2 结果与分析

2.1　土壤风蚀的时间变化特征

近40年来，黄土高原水蚀风蚀交错带的土壤风蚀模数在剧烈波动（变异系数高达82.4%）中显著减小（P<0.05），多年平均模数为1 253 t·km^-2·a^-1，具有明显的阶段性和季节变化。年风蚀模数5 a滑动平均曲线显示，1982—2020年，风蚀特征变化可划分为4个阶段（图2A）：风蚀模数大且剧烈波动下降（1982—1989年）、相对稳定（1990—2000年）、快速下降（2001—2008年）和稳定在较低水平（2009—2020年）。1982—1989年年均风蚀模数（525.1~4 764.9 t·km^-2·a^-1）及其年际变化幅度（变异系数为59.1%）均为近40年来的最大值。1990—2000年，年风蚀模数较稳定（平均为1 291.1 t·km^-2·a^-1），比前一时期平均减小49.4%。2001—2008年风蚀模数快速下降，平均风蚀模数（757.2 t·km^-2·a^-1）分别较阶段一和阶段二减小70.3%和41.4%。2009—2020年风蚀模数维持在较低水平（平均为681.8 t·km^-2·a^-1），除2010年外，年际波动幅度较小（变异系数35.8%）。该阶段的平均风蚀模数分别较阶段一、阶段二和阶段三减小了73.3%、47.2%和10.0%。可见，1982年以来，黄土高原水蚀风蚀交错带经历了两次风蚀模数明显减小并维持相对稳定的变化，表明风蚀环境总体持续改善。

图2

新窗口打开| 下载原图ZIP| 生成PPT

图2 黄土高原水蚀风蚀交错带1982—2020年土壤风蚀模数年际与季节变化

Fig.2 Temporal variations of wind erosion modulus in the wind-water erosion crisscross region of the Loess Plateau during 1982-2020

近40年来风蚀模数的季节变化显示，春季（3—5月）为风蚀高发期，贡献了全年12.8%~84.1%（平均53.5%）的风蚀量。此外，冬季风蚀贡献较大，其变化趋势与春季相似但对全年风蚀贡献比例明显低于春季，占年风蚀总量的4.5%~71.0%，平均为33.1%；而夏、秋两季的风蚀贡献平均仅占全年的1.2%和12.3%（图2B）。可见，春季与冬季，尤其是春季风蚀模数的变化决定了年风蚀模数变化特征，表明春季风蚀防控是该区域风蚀防控的重点。

2.2　土壤风蚀的空间变化特征

2.2.1　土壤风蚀强度的空间分布格局

1982—2020年，研究区土壤风蚀强度长期以轻度及以下为主，总体持续减弱，由强烈及以上等级向中度及以下等级转移，高强度风蚀区范围持续收缩。区域风蚀强度总体西南部较低、东北部及北部较高，微度和轻度风蚀面积占区域总面积的53.4%~98.0%；强烈及以上风蚀主要集中于西北部和北部区域，该区域中度及以上和强烈及以上风蚀区在全区域的占比分别由阶段一的22.4%和12.9%下降至阶段四的6.7%和2.8%，2005年达到近40年来的最低值，中度及以上和强烈及以上风蚀区仅为全区域面积的4.6%和0.4%。为进一步明确风蚀强度的空间变化，结合地理区划、地形地貌、土壤质地等，将区域划分为东部、中部和西部地区。东部地区风蚀强度明显高于中、西部，且变化最剧烈。该地区风蚀强度呈西北高东南低的分布，中度及以上强度的风蚀区连片分布于西北部。1982年，东部地区轻度及以下风蚀区占该区面积的43.5%，强烈及以上风蚀区占36.6%。到2000年，该地区风蚀强度明显减弱，中度及以上风蚀区由连续分布逐渐变为斑块散布，但强烈及以上风蚀区相对面积占比仍然较高（14.4%~21.2%）。到2005年，该地区风蚀强度明显减弱，轻度及以下风蚀区相对面积高达91.2%，强烈及以上风蚀面积降至0.9%。2005年后，东部地区风蚀强度呈现波动特征，从2005年到2010、2015、2020年，中度及以上风蚀区的区内面积占比分别增加了6.3%、0.4%和1.1%（图3）。与东部地区类似，中部地区风蚀强度西高东低且以微度和轻度为主，风蚀强度较高的区域长期集中连片分布于其西北部。1982年，中部地区强烈及以上强度风蚀区的区内面积占比为19.4%，到2000年和2005年逐渐缩减到11.0%和0.1%，2005—2010年出现扩张，达到12.7%，而后持续缩减至2020年的0.1%。西部地区是研究区域中风蚀强度最小的地区，轻度及以下风蚀区的区内面积占比为88.3%~99.4%，强烈及以上风蚀区零星分布（区内面积占比为0.1%~3.4%）。可见，研究区风蚀强度总体呈现东北高西南低的空间格局。

图3

新窗口打开| 下载原图ZIP| 生成PPT

图3 黄土高原水蚀风蚀交错带1982—2020年风蚀强度空间变化特征（以5年间隔为主的典型年份为例）

Fig.3 Spatial variation characteristics in wind erosion intensity in the wind-water erosion crisscross region of the Loess Plateau during 1982-2020 （A case of typical year with interval of 5 years）

2.2.2　不同时期土壤风蚀强度的转移特征

近40年来，黄土高原水蚀风蚀交错带土壤风蚀强度总体呈明显改善，改善区由中东部向西部扩展，局部恶化区集中于西北部。从阶段一（1982—1989年）到阶段二（1990—2000年），风蚀强度明显改善，西部、中部和东部地区轻度改善及以上区域的区内面积占比分别为18.9%、40.0%和49.3%，稳定不变的区域分别为81.1%、60.0%和50.7%，中度改善与大幅改善的区域主要呈片状集中分布在东部地区的西北地带，区内面积占比为12.2%和0.01%。阶段二到阶段三（2001—2008年），轻度改善及以上区域与稳定不变区域在全区域面积占比为22.5%和76.8%，轻度恶化及以上区域零星分布于东部和西部地区，其区内面积占比为0.5%和1.5%。阶段三到阶段四（2009—2020年），总体形成了中部轻度改善和东西两边轻度恶化的空间分异格局。风蚀强度改善区域主要分布于中部、西部地区的交界区以及东部地区的东部，其在全区的面积占比（11.6%）较前一阶段减少了10.9%；同时，分布于东部地区和西部地区西侧的轻度恶化及以上区域面积分别为该地区面积的5.8%和6.7%，较上一阶段增加了5.3%和5.2%（图4）。可见，风蚀强度改善主要发生在前两个阶段，后期改善范围缩小并伴随东部和西部地区的局部波动，但未改变整体改善趋势。

图4

新窗口打开| 下载原图ZIP| 生成PPT

图4 黄土高原水蚀风蚀交错带1982—2020年土壤风蚀强度转移空间分布

Fig.4 Spatial distribution of soil erosion intensity transition in the wind-water erosion crisscross region of the Loess Plateau during 1982-2020

2.2.3　土壤风蚀空间变化趋势

Theil-Sen Median趋势分析与Mann-Kendall显著性检验结果显示，1982—2020年黄土高原水蚀风蚀交错带风蚀强度总体显著减小（P<0.05）。显著减少区域连续分布，范围最广（约16万 km²）。不显著减少（P≥0.05）区域次之（约3万 km²），主要分布于西部地区和东部地区的西部（图5A）。风蚀模数增大（包括显著与不显著增大）及无显著变化区零星散布于西部地区边缘以及东部地区的西部。东部、中部和西部地区的风蚀模数均以显著减少（P<0.05）为主。中部地区风蚀模数显著减少（P<0.05）占该地区总面积的98.1%，为研究区内最高值，且风蚀模数增大的面积<0.1%，表明该地区风蚀减弱最明显。西部和东部地区虽同样以显著减少（P<0.05）为主（分别占各地区面积的80.3%和75.9%），但不显著减少（P≥0.05）区域占比相对较高，分别达到16.4%和21.6%，风蚀模数增大的面积分别为该地区面积的3.0%和2.5%，表明部分区域风蚀减弱幅度相对有限，环境条件相对稳定（图5B）。

图5

新窗口打开| 下载原图ZIP| 生成PPT

图5 基于Theil-Sen Median趋势分析与Mann-Kendall显著性检验的黄土高原水蚀风蚀交错带1982—2020年风蚀模数的空间变化趋势

Fig.5 Spatial variation trends of wind erosion modulus in the wind-water erosion crisscross region on the Loess Plateau during 1982-2020 based of Theil-Sen median trend and Mann-Kendall significance analysis

2.3　不同时期土壤风蚀模数变化的主要影响因素

单因子及分阶段交互探测结果表明，近40年来，区域风蚀模数主要受土壤性质、气候和植被覆盖共同影响。土壤性质对风蚀模数空间分异的重要解释表现为黏粒含量及其与年均风速、年降水量、植被覆盖度、土地利用类型、坡度的交互效应在各时期均较高。在4个阶段，土壤黏粒含量解释力为0.13~0.37，其与上述各因子交互作用的q值依次为0.24~0.46，此外，其与土壤有机质含量的交互作用解释力为0.25~0.34。气候因子对区域风蚀模数变化亦有重要影响，年均风速和年降水量的单因子解释力分别为0.08~0.28和0.02~0.19，二者的交互作用的解释力更高（0.26~0.31），且二者与植被覆盖度交互作用的解释力（平均分别为0.27和0.18）仅次于二者与土壤黏粒含量交互作用的解释力（平均为0.44和0.39）。可见，植被对研究区风蚀模数空间分异具有重要作用（图6）。

图6

新窗口打开| 下载原图ZIP| 生成PPT

图6 基于地理探测器的土壤风蚀关键影响因子分析

注：X₁为年均风速；X₂为年降水量；X₃为植被覆盖度；X₄为土壤黏粒含量；X₅为土壤有机碳含量；X₆为土地利用类型；X₇为坡度

Fig.6 Main impact factors for wind erosion analysis based on the geographical detection results

分阶段比较植被及其与其他因子交互作用对风蚀模数变化的解释力发现，从阶段一到阶段二，区域年均风速及其与年降水量的交互作用解释力均呈降低趋势的背景下，植被覆盖度与两气候因子的交互作用同步下降，二者的q值分别从0.25和0.22降低到0.24和0.18。植被覆盖度与土地利用类型、坡度交互作用的解释力均从0.11降低到0.10。而从阶段三到阶段四，在年均风速q值较高（0.19~0.28）、年降水量q值明显降低（从0.10减小到0.06）且二者的交互作用q值稳定在较高值（0.31）的条件下，植被覆盖度与年均风速交互作用的q值从0.27增大到0.32，其与年降水量交互作用的q值从0.18降低到0.14；其与土地利用类型、坡度的交互作用的解释力分别从0.08和0.09增大到0.10和0.11。可见，2000年以来植被盖度的风蚀防控作用明显增强，反映了区域退耕还林（草）等植被建设措施的成效。

3 讨论

3.1　基于RWEQ模拟的水蚀风蚀交错带风蚀模数可靠性

对比本文与前人估算的黄土高原水蚀风蚀交错带土壤风蚀模数发现，本文通过RWEQ估算的区域风蚀模数总体与现有结果基本一致。前人采用Cs-137示踪、Be-7示踪、定位观测、模型模拟等方法获得的结果显示，不考虑风积的条件下，黄土高原水蚀风蚀交错带不同区域的风蚀模数为141~4 230 t·km^-2·a^-1，本文在与前人研究区域相同或类似区域的风蚀模数为295.2~3 414.6 t·km^-2·a^-1（表1），尤其是在陕西靖边、陕西神木、水蚀风蚀交错带丘陵沟壑区等区域结果一致性较高。然而，本文估算的风蚀模数范围（434.1~1 429.2 t·km^-2·a^-1）与前人在山西大同、陕西定边等地基于Cs-137示踪方法估算的风蚀模数有明显差别。这可能与该部分风蚀模数采用风蚀和水蚀对风水复合侵蚀贡献的总比例进行估算相关。结合前人的多项结果，黄土高原水蚀风蚀交错带土壤风蚀约占总侵蚀的10%~40%，普遍在25%左右^［44-46］，且不同研究中区域、时段、环境与气候条件等的变化也可能造成明显差异。可见，采用不同方法在同一研究区或相似环境下的结果大多可与本文相互印证，表明本文基于RWEQ估算的区域风蚀模数具有较高可信度。

表1 本文与前人量化的风蚀模数结果的差异

Table 1 Comparison of wind erosion modulus between this study and previous results

区域/研究时期	研究方法	土壤总侵蚀模数/（t·km^-2·a^-1）	土壤风蚀模数/（t·km^-2·a^-1）	本文相应区域的风蚀模数/（t·km^-2·a^-1）	数据来源
山西大同/2008年	Cs-137示踪	17 500	4 375^*	1 052.7	姜洪涛^[47]
宁夏彭阳/2007年		2 500~10 000	625~2 500^*	434.1	马远远等^[48]
陕西绥德/2005年		12 000	3 000^*	282.5	李勉等^[49]
陕西靖边/2020年左右		4 849~9 463	1420~3 138^*	1 210.6	Zou等 ^[44]
陕西神木/2014年		8 060~13 370	519~4 230^*	1 429.2	Zhang等^[50]
陕西定边/2014年		1 513~8 314	378~2 079^*	3 414.6	脱登峰^[51]
陕西神木/2015—2016年	Be-7示踪	—	141~771^*	1 362.2	Zhang等^[45]
建设三北防护林的黄土高原水蚀风蚀交错带丘陵沟壑区/2000—2018年^§	RWEQ	—	600~1 700	295.2~1 792.5	张雄一等^[21]
黄土高原水蚀风蚀交错带/1980、1990、2000、2010年和2020年	RWEQ	—	246	1 253.0	黎恩丹^[13]

注：*风蚀模数根据区域风蚀和水蚀的比例估算，已有结果表明黄土高原水蚀风蚀交错区风蚀在总侵蚀中的贡献约为10%~40%，多为25%^［44-46］。因此，按照风蚀平均贡献为25%估算。§黄土高原建设三北防护林工程的丘陵沟壑区略大于本文的区域，因此，本文中采用的风蚀模数为整个水蚀风蚀交错带2000—2018年的年均风蚀模数范围。

新窗口打开| 下载CSV

3.2　气候变化和植被建设对区域风蚀强度时空变化的作用

虽然土壤性质一直是影响区域土壤风蚀的重要因素，但近40年水蚀风蚀交错带风蚀强度时空变化与该区域的气候变化与植被建设密切相关^{［16，52］}。现有结果表明以植被建设为主的人为活动对风蚀模数减小的贡献为23%~91%^［52-53］。在水蚀风蚀交错带，气候因子和植被覆盖度在不同时段的变化趋势不相同，导致二者对风蚀强度阶段性变化的影响程度也具有明显差异^［25］。对比2000年前后近20年，水蚀风蚀交错带的年均风速不显著下降（P≥0.05），年降水量由下降趋势转为上升趋势，但均未达到显著水平（P≥0.05），植被覆盖度显著上升（P<0.05），阶段四的平均植被覆盖度较阶段一（42.0%）提高16.5%（图7）。对比不同时期气候和以退耕还林（草）等生态工程为主的植被建设对不同阶段风蚀模数变化的贡献发现，从阶段二到阶段四，气候变化的贡献从70.4%降低到37.9%，而植被建设的贡献从29.6%增大到了62.1%（表2）。可见，2000年后，在风速未明显变化，降雨不显著增大的条件下，植被盖度增大在区域风蚀防控中发挥了巨大作用。

图7

新窗口打开| 下载原图ZIP| 生成PPT

图7 1982—2020年年均风速、年降水量、植被覆盖度年均值变化趋势

Fig.7 Variation trends in annual mean wind speed， precipitation， and vegetation cover from 1982 to 2020

表2 气候变化和植被建设对黄土高原水蚀风蚀交错带土壤风蚀模数的影响

Table 2 Effects of climate change and vegetation restoration on wind erosion modulus in the wind-water erosion crisscross region of the Loess Plateau

时段	风蚀模数/（t·km^-2·a^-1）	气候变化贡献/%	植被建设贡献/%
1982—1989（阶段一，基准期）	2 553.04	—	—
1990—2000（阶段二）	1 291.07	70.4	29.6
2001—2008（阶段三）	757.19	53.0	47.0
2009—2020（阶段四）	681.84	37.9	62.1

新窗口打开| 下载CSV

4 结论

1982—2020年，黄土高原水蚀风蚀交错带年均风蚀模数显著下降（P<0.05），多年平均模数为1 253 t·km^-2·a^-1，在退耕还林（草）前后均有一个明显下降而后趋于稳定的变化。1982—1989年剧烈波动下降，而后在1990—2000年保持相对稳定，稳定期的平均风蚀模数为1 291.1 t·km^-2·a^-1，较前一阶段平均减小49.4%；2001—2008年风蚀模数快速下降，而后在2009—2020年维持较低水平且相对稳定（681.8 t·km^-2·a^-1），较阶段一降低73.3%。空间上，风蚀强度总体呈现东部中部高西部低、西北部高东南部低的空间格局，强烈风蚀区由集中连片分布变为零星分布。区域风蚀强度总体以轻度和微度为主，风蚀强度明显由高强度向轻度及以下转移，风蚀模数显著减少区占研究区总面积的80%以上。中度及以上、强烈及以上风蚀区面积占比由1982—1989年的22.4%和12.9%降至2009—2020年的6.7%和2.8%，风蚀环境持续改善。区域风蚀强度变化主要受土壤性质、气候和植被的共同作用，从阶段一到阶段四，退耕还林（草）为主的人为活动的贡献从29.6%增大到62.1%，充分体现了区域生态工程建设的良好风蚀防控效应。

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

唐克丽

黄土高原水蚀风蚀交错区治理的重要性与紧迫性

［J］.中国水土保持，2000（11）：11-12.