1990—2015年内蒙古高原土壤风蚀时空差异特征

吴晓光 姚云峰 迟文峰 闫旭东 张宏飞 金 良 白文科 秦富仓

(1．内蒙古农业大学 沙漠治理学院，呼和浩特 010011；2．内蒙古自治区土地整治中心，呼和浩特 010020；3．内蒙古财经大学 资源与环境经济学院，呼和浩特 010017；4．西华师范大学 生态研究院，四川 南充 637002)

随着联合国千年生态系统评估计划(MA)在全球范围的广泛开展，生态评估已成为当今国际生态学领域研究的前沿课题之一[1-2]。中国北方生态环境脆弱，气候变化与不合理的人类活动共同作用，导致了以土壤风蚀为主要过程的水土流失与土地荒漠化问题[3-5]。因此，防风固沙服务功能的提升是其生态安全屏障的重要基础[6]。而植被覆盖是定量评价干旱与半干旱区生态质量的重要指标[7-8]，也是评价和预测土壤风蚀变化的基础[9]。有研究表明不同土地利用/覆盖类型条件下植被覆盖状况的季节差异明显，林地与草地比耕地具有更好的水土保持效果，春季与冬季表现最为明显[10]。为此，有研究指出调整地表植被覆盖结构及措施是防治土壤风蚀的有效途径[11-12]。耕地留茬免耕可显著提高地表的春秋冬季抗风蚀能力[13]。但这些研究大多集中在样地或小区域尺度，而在大尺度上(例如内蒙古高原地区)土壤风蚀的季节变化研究仍相对薄弱。

内蒙古高原作为祖国北方重要的生态屏障，近年来荒漠化严重。内蒙古属于典型的农牧交错带，作物以一年一熟为主，耕地每年近7个月处于裸露状态，春秋季地表干燥，大风天气频发加剧土壤风蚀过程，致使扬沙、沙尘暴等灾害频发[13-14]。内蒙古东西跨度较大，土地利用/覆盖类型由西向东总体呈荒漠戈壁、荒漠草原、典型草原过渡到大兴安岭林区，经向分布植被覆盖差异明显，季相变化对土壤风蚀影响巨大。区域内严重的土壤风蚀不仅破坏了地表，导致水土流失、土地退化、荒漠化等生态系统问题加剧，也严重威胁着国家生态环境安全[7,15]。基于此，本研究在遥感-地面同步观测试验、遥感定量反演与验证、模型系统参数化的基础上，构建1990—2015年内蒙古高原土壤风蚀数据集，辅以植被覆盖、土地利用/覆盖类型、气候、基础地理背景等数据，定量刻画土壤风蚀时间差异以及相应的季节变化特征与植被覆盖之间相互关系，以期界定有效抑制土壤风蚀的植被覆盖度调控阈值，为内蒙古高原土壤风蚀防治提供科学依据。

1 材料与方法

1.1 研究区概况

内蒙古高原土地面积约124.7×104km2，东至大兴安岭和苏克斜鲁山，西至马鬃山，南界祁连山麓和长城，北接蒙古国与俄罗斯东部，地理坐标为92°68′～125°83′E，36°16′～52°99′N，平均海拔1 000 m；研究区属干旱与半干旱气候，降水较少且年际变率大，东部以林地为主，中部以草地和农田生态系统为主，西部以沙漠和戈壁为主(图1)。受地形、植被、土壤、气候等自然原因和过牧、矿产开发等人为因素作用，该地区成为中国生态环境退化和水土流失最严重的地区之一。

图1 内蒙古高原地理位置及土地覆被示意图Fig.1 Location and land cover of Inner Mongolia Plateau

1.2 数据来源与处理

研究基础数据主要包括1990、2000、2005、2010和2015年土地利用/覆盖数据，数字高程(DEM)、归一化植被指数(NDVI)、气象、植被类型数据、土壤类型与理化性状、雪深数据。其中，土地利用/覆盖数据源于中国科学院“全国土地利用1∶10万数据库”，数据综合精度达91.2%[16]；DEM从USGS(http:∥glovis.usgs.gov/)获取的30米分辨率数字高程数据产品；NDVI数据则从NASA(https:∥ecocast.arc.nasa.gov)获取1990—2015年GIMMS3g和MODIS-NDVI数据集。采用通用“像元二分模型”获取年、月植被覆盖度数据集[17-18]；气象数据源于中国气象科学数据共享服务网(http:∥cdc.cma.gov.cn/home.do)，为了保证气象资料的完整性及持续性，构建稳定的气候因子数据，共选取研究区及中国周边气象台站(图1)的1990—2015年逐日平均气温、降水量、平均风速、风向、沙尘暴和日照时数等要素，数据经过质量检查和控制处理，消除奇异值等非气候因素造成的影响，并对缺失数据采用临近站点数值进行插补，将整理后的数据通过ANUSPLIN软件插值，获取长时间序列空间数据集；植被类型数据通过中国科学院资源环境科学数据中心获取(http:∥www.resdc.cn/)；土壤数据(1∶100 万土壤图)源自中国西部环境与生态科学数据中心(http:∥westdc.westgis.ac.cn)，包含土壤类型、理化性状、有机质含量、砂粒、土层厚度及粘粒含量等信息[19]；雪深源自中国西部环境与生态科学数据中心(http:∥westdc.westgis.ac.cn)，将获取的长时间序列数据集计算后形成雪盖因子数据。本研究的评价单元将不同精度的数据同化到1 km网格。

1.3 土壤风蚀模型

本研究采用修正土壤风蚀方程(RWEQ)评估土壤风蚀量，在充分考虑气候条件、植被状况、地表粗糙度、土壤可蚀性、土壤结皮等要素情况下，RWEQ模型可较好定量评估土壤侵蚀模数[17,20-21]，模型基本形式为：

(1)

S=150.71·(WF·EF·SCF·K′·COG)-0.371 1

(2)

Qmax=109.8[WF·EF·SCF·K′·COG]

(3)

式中：SL为土壤风蚀量，g/m2；z为下风向最大风蚀出现距离，m；S为关键地块长度，m；Qmax为风沙最大转运容量，kg/m；EF(无量纲)为土壤可蚀性因子；WF为气候因子，kg/m；K′(无量纲)为地表粗糙度因子；SCF(无量纲)为土壤结皮因子；COG(无量纲)为植被因子，主要涉及不同季节植被生长、枯萎、农作物及其他植被(草地)残茬等状况。

1)气象因子WF即风速、温度、降水、雪盖因子等多气象因子对土壤风蚀综合作用，公式如下：

(4)

(5)

(6)

(7)

SD=1-P(Hsnow>25.4 mm)

(8)

式中：WF为气象因子，kg/m；U2为2 m处风速，m/s；Ut为2 m处的起沙临界风速，m/s(假定为5 m/s)；ρ为空气密度，kg/m3；EL为海拔高度，km；T为绝对温度，K；N和Nd分别代表风速和试验的天数，d；SW为土壤湿度因子(无量纲)；ETp为潜在相对蒸发量，mm；R为降雨量，mm；I为灌溉量，mm；Rd为降雨或灌溉天数，d；Nd2为各月风速>5 m/s天数，d；SR为太阳辐射总量，cal/cm2；DT为平均温度，℃；SD为雪盖因子；Hsnow为雪覆盖深度，mm。

气象因子数据通过中国气象科学数据共享服务网获取，要素包括日均降水、气温、日照数和风速(http:∥cdc.cma.gov.cn/home.do)等，采用ANUSPLIN软件对长时间序列进行插值处理；雪盖因子利用中国西部环境与生态科学数据中心(http:∥westdc.westgis.ac.cn)下载的中国雪深长时间序列数据集来计算。

2)土壤可蚀性因子(EF)是土壤理化、机械组成等因素作用下引起的土壤风蚀量，公式如下：

(9)

式中：EF(无量纲)为土壤可蚀性因子；Sa为土壤粗砂含量，%；Si为土壤粉砂含量，%；CL为土壤粘粒含量，%；OM为有机质含量%；CaCO3为碳酸钙含量，%。土壤类型、有机质含量、砂粒和粘粒含量数据来源于中国西部环境与生态科学数据中心(http:∥westdc.westgis.ac.cn)提供的1∶100万土壤图及所附土壤属性；CaCO3数据源来自于地球系统科学数据共享服务网(http:∥www.geodata.cn)的 1∶400 万土壤碳酸钙含量数据。

3)土壤结皮因子SCF表征土壤表层的坚硬结皮有效抵抗风蚀能力的大小，公式如下：

(10)

式中：CL为土壤粘粒含量，%；OM为有机质含量，%。

汪队长的眼睛，睁得像铜铃。一听男兵要吻他的女兵李晓英，顿时“啊”了一声，高声骂道：“你混蛋，战场上调戏女兵是要掉脑袋的，你不知耻还来说情，你这个连长也当到头了。”

4)植被因子COG表示一定植被覆盖条件下对土壤风蚀的抑制程度，公式如下：

COG=e-0.048 3(SC)

(11)

式中：SC为植被覆盖度，%。

植被因子主要由土地利用/覆盖数据和NDVI数据综合处理获取，NDVI数据基于美国国家航空航天局网站(NASA: https:∥ecocast.arc.nasa.gov)获取。

5)地表粗糙度因子K′是由地形引起的土地表面粗糙程度对土壤风蚀的影响，公式如下：

K′=e(1.86Kr-2.41Kr0.934-0.124Crr)

(12)

(13)

式中：Crr为随机粗糙度因子，cm；Kr为地形粗糙度因子，cm；L为地势起伏参数；ΔH为离L范围内的海拔高程差，m。

1.4 模型修正

为了提高RWEQ模型模拟精度，采用中国科学院“全国1∶10万土地利用数据库”获取的1990、2000、2005、2010和2015年数据集对模型进行修正[16,22]。将耕地、林地、草地等多种土地利用类型处理成 1 km×1 km 栅格面积比例数据[23](图2(a))，在计算土壤风蚀影响因子时计算到模型参数图层中，并以加权方式获取RWEQ模型参数因子[21]。

(14)

在详细分析内蒙古高原土地利用/覆盖类型及地貌特点、气候差异的基础上，野外调查和收集不同区域地表粗糙度、不同季节植被覆盖度等信息，以提升RWEQ模型相应参数精度。采用野外直尺测定试验方法[24]测定不同区域植被覆盖类型的地表粗糙度(图2(b))，将野外观测获取的信息值与基于遥感的土地利用/覆盖类型相结合(式12)；针对不同季节植被覆盖差异特征，对生长期植被、枯萎期植被、植被直立残茬等状况进行野外采样，获取四季不同区域植被覆盖度信息，综合土地利用/覆盖数据、野外采样数据与植被覆盖度遥感数据加权估算不同季节地表植被覆盖信息(图2(c)和图2(d))，利用回归方程对不同时期遥感数据进行校正处理。

2 结果与分析

2.1 RWEQ模型精度验证

土壤风蚀模数监测方法诸多，成熟方法包括风洞试验、土壤剖面分析、剖面粒度分析、插钎法、137Cs 测定、集沙仪测定法等[25-26]，本研究根据文献收集与实地地面测定的中国北方不同地区风蚀模数数据的结果进行验证[27]。本研究基于RWEQ模型，获取了1990—2015年内蒙古高原土壤风蚀模数数据集(图3)。

收集和获取已有研究者在内蒙古高原区域以137Cs 法推算土壤风蚀模数的研究成果[25-27]，并对本研究反演成果进行验证(表1)。结果显示：基于RWEQ模型反演结果与采样结果拟合性较好(R2=0.85,P<0.01)，说明该模型在内蒙古高原区域尺度上在理论上可行。

2.2 土壤风蚀季节差异分析

从内蒙古高原多年平均的春、夏、秋、冬土壤风蚀量统计上看，春季土壤风蚀量最大，冬季次之，夏季土壤风蚀量最小；春季、夏季、秋季、冬季土壤风蚀量比例为5.42∶1.00∶1.36∶3.00(图4)。

总体上看，土壤风蚀春季与冬季变化表现最为突出。土壤风蚀量表现为春季>冬季>秋季>夏季；2001年春季与冬季土壤风蚀达到峰值，其中春季与冬季土壤风蚀量分别为32.34×108和 55.40×108t，第二次峰值出现在2010年，相对于2001年风蚀量明显减少；土壤风蚀量时序变化过程中，波谷均出现2008年，究其原因主要因为风场强度的降低(风蚀力低)以及气候的干燥度下降，尤其是2004年后生态建设致使植被恢复明显(图5)。而2010年不同季节土壤风蚀表现明显，主要由风场强度大幅升高引起的。内蒙古高原春、冬季植被覆盖率低，农田植被覆盖受季节影响最大，春和冬季农田基本为裸露状态，极易造成土壤风蚀；草地覆盖在非生长季的春、冬季作用极为重要，内蒙古高原以草地生态系统类型为主，在人为过牧(山羊破坏最为突出)的干扰下会致使春和冬季的大面积草地形成裸露状态，进而增加土壤风蚀。

图3 研究区1990—2015年土壤风蚀平均模数Fig.3 Average soil wind erosion modulus from 1990 to 2015

表1 实测结果与RWEQ模型反演结果对比Table 1 Comparison of measured results and RWEQ model inversion results

图4 内蒙古高原1990—2015年春夏秋冬土壤风蚀多年平均值统计Fig.4 Average soil wind erosion of different seasons in Inner Mongolia Plateau from 1990 to 2015

2.3 不同土地利用/覆盖类型植被覆盖季节差异特征

由表2可知，15年期间，10种土地利用/覆盖类型下不同季节(冬季除外)平均植被覆盖的相关系数、变化量及变化速率都呈显著增加趋势，夏季与秋季土地利用/覆盖类型增加速率均大于春季；耕地尤其是水田植被覆盖增加速率夏季与秋季最大，旱地次之，进一步揭示了内蒙古高原植被生长季覆盖变化差异明显；中、低覆盖度草地和荒漠植被在不同季节的变化方式和大小呈相似性，草地和荒漠植被生长季与非生长季植被覆盖度差异较为显著。

图5 内蒙古高原年和季节土壤风蚀量变化特征Fig.5 Characteristics of soil wind erosion amount of different seasons in Inner Mongolia Plateau

在典型区(图1)，春季至次年春季不同土地利用/覆盖类型转化过程中对土壤风蚀影响研究结果表明(表3)：低覆盖度向高/中覆盖度草地类型转化有效抑制土壤风蚀(-61.42%)，尤其是春季植被覆盖度的增加有效抑制土壤风蚀；荒漠周边区域，中覆盖度向低覆盖度草地类型转换过程中，土壤风蚀量呈增加(51.39%)态势，草地转化为沙地、旱地土壤风蚀量增加。

2.4 植被覆盖对土壤风蚀响应阈值

基于土地利用/覆盖类型对土壤风蚀的影响分析，揭示了植被覆盖变化对土壤风蚀起到重要的作用。研究结果发现，1990—2000年不同程度植被退化对土壤风蚀造成的影响差异较大，植被重度退化区域土壤风蚀表现出侵蚀模数增加(1.59 t/(hm2·年))，轻度退化致使土壤风蚀模数增加0.40 t/(hm2·年)；2000—2010年内蒙古高原植被覆盖总体呈增加趋势，其中植被轻微改善对土壤风蚀的作用不明显(减少0.10 t/(hm2·年)，中度改善区域土壤风蚀模数减少1.07 t/(hm2·年)，明显改善区域土壤风蚀模数减少1.33 t/(hm2·年)。基于植被覆盖多年均值与土壤风蚀多年均值建立回归关系，发现植被覆盖度达20%～25%时可迅速降低土壤风蚀模数，植被覆盖度超过75%时，随植被覆盖度的增加土壤风蚀模数变化幅度较小，植被覆盖达到75%时有效抑制土壤风蚀(图6)。秋季耕地和草地收割过程中，适当留茬从而增加植被残茬，对抑制土壤风蚀具有积极作用。

表2 2000—2015年不同土地利用/覆盖类型植被覆盖变化统计Table 2 Statistical analysis of vegetation cover changes of different land use/cover types from 2000 to 2015

表3 典型区不同季节土地利用/覆盖结构变化对土壤风蚀影响Table 3 Effect of land use/cover changes on soil wind erosion in different seasons of typical areas %

图6 内蒙古高原植被覆盖与土壤风蚀相关性阈值分析Fig.6 Relevance threshold analysis of Inner Mongolia Plateau vegetation cover and soil erosion

3 结论与讨论

3.1 结论

1)不同季节土壤风蚀差异较大，春季、夏季、秋季、冬季土壤风蚀量比例为5.42∶1.00∶1.36∶3.00。内蒙古高原在2001年春季与冬季土壤风蚀达到峰值，其中春季与冬季土壤风蚀量分别为32.34×108和55.40×108t；第二次峰值出现在2010年，相对于2001年风蚀量明显减少。

2)基于土地利用/覆盖变化遥感监测信息[16,22]，1990—2000年，“粗放”经营管理下，耕地转变成未利用地达380.24 km2；其中，转变成沙地与盐碱地的面积分别为31.29与318.44 km2。1998年实施“退耕还林还草”政策，植被覆盖恢复明显，耕地转换成林地的面积为567.07 km2，转换成草地的面积达2 737.42 km2，研究区中部与西部区域草地与林地面积增加明显，为土壤风蚀的减缓起到了关键性的作用，退耕还林还草有效的固化了土壤，起到防风固沙的积极作用。

3)2000—2015年，不同土地利用/覆盖类型的植被覆盖度都呈显著增加趋势，夏季与秋季土地利用/覆盖类型增加速率均大于春季。春季至次年春季表现为低覆盖度向高覆盖度草地类型转化能够有效抑制土壤风蚀(-61.42%)；荒漠周边区中覆盖度向低覆盖度草地类型转换过程中，土壤风蚀量增加(51.39%)，草地转化为沙地、旱地土壤风蚀量增加。

4)植被覆盖达到75%时可有效抑制土壤风蚀，低覆盖度植被对降低土壤风蚀作用敏感性高于高覆盖度植被，当植被覆盖度达20%～25%时可迅速降低土壤风蚀模数。

3.2 讨论

区域土壤风蚀、退耕还林还草对土壤风蚀的影响是一个复杂而且系统的研究，涉及到多背景、多学科的交叉融合[28]。因此鉴于作者知识背景和时间限制，仍然在几个方面需要改进和深入研究：

土壤风蚀模型应用的不确定性分析。RWEQ模型主要是应用在农田生态系统，且参与计算的因子较多，精度要求高等[28-29]，一方面，区域内的气象站点较少，部分采用周围的气象站点数据，另一方面气象数据采用空间插值法实现区域覆盖，插值的精度具有不确定性。基于此在小尺度研究土壤风蚀过程，应大量补充气象数据。

土壤风蚀是一种无法消除的过程，只能是控制，而不可能使区域内土壤风蚀几乎为零[29]。在土壤风蚀驱动因子贡献中，在同一自然条件下，控制人类活动降低土壤风蚀是最有效的解决办法。本研究从理论上界定了自然因素与人为因素加速和改善土壤风蚀作用过程的贡献程度，并未对人类活动贡献程度的显著性进行研究，在进一步研究的过程中，将继续对该部分进行深入研究，界定人类活动贡献显著性，为自然资源可持续发展和生态修复等工作提供科学的参考和借鉴依据。

本研究中并未对退耕为草地、灌木林地、牧草地等进行分类更深层次的研究，内蒙古高原地理单元特殊，水资源条件有限，充分考虑气候、水资源等条件，退耕后地类的选择对区域生态效应更为显著、更为适宜，有待进一步研究。