2001-2020年黄土高原干旱时空动态及其对气候变化的响应

樊莉莉, 耿 斌, 王吉林, 陈 杰

(1.石家庄职业技术学院, 石家庄 050081; 2.黄河水利科学研究院, 郑州 450003; 3.武汉大学 水资源与水电工程科学国家重点实验室, 武汉 430072)

黄土高原作为中国重要的农业区之一，大部分地区属于半干旱和亚湿润气候，沿东南—西北方向的年平均降水量(1961—2009年)为143.6～811.8 mm，干旱事件的频繁发生成该地区21世纪面临最复杂的挑战之一[1-2]。但干旱事件的发生是一个复杂的过程，受气候变化、植被退化、人类活动等多种因素的影响，由于其发生频率高、持续时间长、影响范围广等原因，目前干旱已成为影响农业、生态最严重的气象灾害[3-4]。近年来，在全球变暖的影响下，黄土高原地区的温度将在整个世纪持续升高，并随后加速该地区的水分流失，使土地退化与恢复、干旱成为严重影响黄土高原生态环境的气候因子[5]。因此，及时、有效监测黄土高原干旱发展规律已成为亟待解决的问题。

目前已有诸多学者对黄土高原的干旱进行了研究[6-10]。在干旱变化趋势的研究中，有学者发现1960—2016年黄土高原干旱整体呈上升趋势，山西西部、宁夏北部和甘肃中东部地区的干旱加重最为显著，而西北部和西南部则表现为变湿趋势[6-9]。也有学者发现未来时期(2018—2100年)黄土高原干旱发生次数明显增多，未来干旱趋势显著加重的区域占黄土高原面积的51.62%～99.90%，80%的地区将发生从无旱到轻旱、从无旱到中旱的等级转换[10]。因黄土高原降水主要集中在夏、秋季，因此不同地区的干旱程度也因季节不同而存在显著差异。如有研究发现黄土高原东南部和西部在春、冬季干旱发生频率较其他地区频繁，而西北部干旱在夏、秋季发生频率较其他地区高，这主要与降水量、气温和蒸散发存在显著相关性[6,10-12]。也有研究发现黄河上游的干旱危害趋于严重，中游干旱危害过去40 a中已得到一定程度缓解[13]。虽然目前针对黄土高原干旱时空动态的研究很多，但很多是基于“点”尺度进行的，不同研究学者的结论不一致，且很少探讨未来干旱的持续情况和季节性气候对干旱的影响情况。而黄土高原干旱作为我国长期的一些生态环境监测任务，因此需要继续关注。温度植被干旱指数(TVDI)综合考虑了温度与植被双重特征，可以较为综合的反映一个地区干旱状况[14-18]。如刘一哲[19]、齐述华[20]、程伟[21]、候静[22]等分别对藏北地区、全国、内蒙古草地的干旱变化特征进行分析，发现综合了植被覆盖信息和陆地表面温度信息的TVDI指数能很好模拟土壤表层水分变化情况，在地表干旱监测方面具有较大应用潜力[19-21]。

因此，本文基于NASA提供2001—2020年黄土高原的MOD13A-NDVI和MOD11A2-LST数据，基于Python编程计算得到2001—2020年黄土高原季节性TVDI指数，并基于TVDI指数探究了不同植被分区不同季节干旱时空动态和季节性气温、降水对季节干旱的影响情况，以期为我国黄土高原生态工程制定和未来生态环境治理提供参考意见。

1 研究区与研究方法

1.1 研究区概况

黄土高原位于中国北方(33.7°—41.3°N，100.8°—114.6′E)，穿越黄河上中游，总面积达6.49×105km2，包括山西、内蒙古、陕西、河南、宁夏、甘肃、青海7省。区域处在沿海向内陆、平原向高原过渡地带，自南向北兼跨暖温带、中温带两个热量带，自东向西横贯半湿润、半干旱两个干湿区，受经纬度和地形的双重制约。降水量少而不均，从西北向东南逐渐增加，6—9月，该地区60%～70%的年降水量以高强度风暴的形式出现。年平均气温为3.6℃，西北部为14.3℃。该地区全年气温日较差较大，冬季干冷，夏季湿热，秋季气温迅速下降，春季气温迅速上升。据估计，该地区的年潜在蒸发量(865～1 274 mm)超过了降水水平[3](图1)。

注：Ⅰ为温带南部荒漠草原区，Ⅱ为温带南部典型草原区，Ⅲ为温带南部森林(草甸)草原区，Ⅳ为暖温带北部落叶栎林区，Ⅴ为暖温带南部落叶栎林亚区，下同。图1 黄土高原植被分区及气象站点空间分布

1.2 数据来源

(1) 遥感数据。NDVI作为植被绿色度和活动性的衡量标准，已被广泛用于在大空间尺度上检索植被生长监测方面。我们使用了由星载传感器MODIS观测得出的的NDVI数据计算TVDI指数。MODIS NDVI数据集是从MOD13A3产品版本6产品中提取的，空间分辨率为1 km，时间间隔为16 d。地表温度数据(LST)选取MOD11A2产品提供的白天地表温度数据，该数据空间分辨率为1 km，时间分辨率为8 d，数据处理方法与NDVI数据一致。为进一步去除异常值的影响，本文基于Timesat3.3插件提供的S-G滤波法对NDVI和LST数据进行去噪，然后采用最大值合成法合成为月值数据[23]。LST在不同的海拔高度上会因气温和大气湍流的影响而存在显著变化[24-25]。而黄土高原地形起伏明显，气温和降水的时空差异性大，因此本文对地表温度进行地形校正才能满足后续要求，地表温度订正公式见文献[26—27]。

(2) DEM数据。因海拔会影响地表温度，因此本文使用DEM对LST数据进行地形校正。我们使用了Google Earth Engine(GEE)平台提供的2000年的90 m的NASA DEM和2015年的30 m的AW3D作为LST地形校正的DEM数据。首先将其采样成1 000 m分辨率，2001—2010年的LST采用NASA DEM进行校正，2011—2020年的LST采用AW3D进行校正。

(3) 气象数据。本文的气候要素主要为气温和降水，来源于2001—2020年黄土高原及其周边地区的气象站(图1)，数据可从中国气象局数据网(http:∥data.cma.cn)获取。本研究利用ANUSPLINE以不同时段的DEM为协变量，将气温和降水数据从点尺度插值成为1 000 m分辨率的栅格数据。将12—2月定义为冬季，3—5月定义为春季，6—8月定义为夏季，9—11月定义为秋季[28]。

(4) 中国植被区划数据来源于1∶100万植被图。该植被区划数据是依据植被类型和地理分布特征进行划分区域，可以反映植被地理分布的规律性及其与环境的关系，提供区域或全球的植被地理图式[29]。

1.3 研究方法

1.3.1 温度植被干旱指数(TVDI) LST对裸露土壤或稀疏植被条件下的地表加热过程(热容量和热导率)有影响，因此对地表土壤含水量很敏感。随着植被覆盖度的增加，LST 对地表土壤含水量的变化响应缓慢。Price等[30]发现NDVI和LST的遥感像素值散点图通常为三角形或梯形。后来，Sandholt等[31]基于LST-NDVI典型的关系构建了TVDI指数用于反映地表干旱程度。理论上，在三角域中，与NDVI轴平行的三角形底边对应于研究区域内蒸散量最大的像素点，三角形的顶边(斜边)对应于研究区域内的零蒸散量像素点。在三角形顶部和底部之间，三角形的大部分对应于不同干旱条件下的蒸发量。随着NDVI沿x增加轴，最大LST减小，可以使用最小二乘法拟合负斜率，定义为干边缘。湿边缘由一组点组成，这些点为不同的地表植被覆盖形成一条水平线或倾斜线。在垂直方向上，NDVI值恒定的三角形域的底部到顶部，由于表层土壤的水分胁迫，LST逐渐增加，从湿边处的最小值到湿边处的最大值[32]。干旱边缘地表土壤含水量相应地从最大值减小到最小值，TVDI从0增加到1，表明地表从极端潮湿变为极端干旱。TVDI使用以下等式估算：

(1)

式中：Ts为观测到的地表温度;Tsmax和Tsmin分别为像素的最小和最大地表温度，它们在研究区域中具有相同的NDVI值，并被定义为湿边缘和干边缘。Tsmax和Tsmin通过散点图上下限的点组计算，他们的方程式如下：

Tsmax=a+b·NDVI

(2)

Tsmin=c+d·NDVI

(3)

式中：a，b，c，d为在整个研究区域内估计的干湿边缘参数，包括整个表面水分含量范围。

1.3.2 TVDI变化趋势分析 本文采用Sen斜率法分析近20 a黄土高原TVDI变化趋势。并基于Mann+Kendall对Sen趋势分析的显著性进行检验。Sen的计算公式为[33]：

(4)

式中:β为干旱变化趋势;median为取中值函数;当β<0时,TVDI呈增加趋势;β>0时,TVDI呈减小趋势。

1.3.3 TVDI未来持续状态分析 本文利用Hurst指数对黄土高原不同季节的TVDI进行未来变化趋势预测。Hurst指数是根据时间序列长期依赖性或持续性的一种有效方法，一般采用R/S分析方法计算Hurst指数，具体算法见参考文献[34]。其中，Hurst>0.5，表示序列TVDI的长期相关性为正持续性；Hurst=0.5，表明TVDI变化为随机序列，与过去没有关系；Hurst<0.5，表示序列TVDI长期相关性特征为反持续性[35]。

1.3.4 气温、降水对TVDI的影响分析 本文采用偏相关分析法分析季节性气温和降水对不同季节TVDI的影响情况[36]。气温、降水和TVDI之间的偏相关系数计算公式如下：

(5)

式中:rij为i和温度j之间的偏相关系数，控制变量为降水量，反之亦然；i为TVDI；j和n分别为季节性气温和降水量。

2 结果与分析

2.1 TVDI空间分布特征

参阅吴黎等[36]对研究区TVDI进行分级，TVDI≤0.46，地表湿润，无旱；0.46

使用均值法计算得到20 a间黄土高原TVDI的时间特征发现(表1，图2)：黄土高原TVDI平均每年的增长速率为0.000 2，TVDI集中在0.6～0.7，属中旱等级，即土壤表面干燥、植被叶片有缺水干黄现象。2005年TVDI有最大值(0.738 2)，2012年有最小值(0.591 8)。但不同季节TVDI的变化特征各不相同：春季TVDI的增加速率为四季最小(增长速率=0.000 2/a，R2=0.000 3)，变化速率与年均值相同，变旱趋势不明显，但从每年均值的变化来看，春季TVDI的波动性最大，且其波动走势与年均值走势一致，但波动性比年均值大；夏季TVDI的增加速率为0.001 2/a，R2为0.453,，变化速率大于春季和年均值；秋季TVDI的增加速率最大(增长速率=0.004 4/a，R2=0.283 6)，变化较为明显，且与年均值变化基本呈相反趋势；冬季TVDI在20 a间的变化呈现变湿趋势，TVDI呈减小趋势(减小速率=-0.004 3，R2=0.339 5)，减小速率明显，且与年均值变化趋势基本呈相反趋势。总体而言，2012年为TVDI变化的突变年份。

图2 TVDI时间曲线

表1 TVDI年际变化相关系数

从TVDI均值空间图发现，近20 a黄土高原干旱整体为北部和东南部高，西部和西南部低的空间分布特征，最低值集中在黄土高原西部的青海省的西宁、门源、贵德、甘肃的临夏、临洮及东部山西的太原地区，TVDI基本在0.57以下，即属于轻旱、无旱级别。最高值集中在黄土高原东南部河南和山西的长武等地区，TVDI值基本在0.76以上，属重旱级别。从不同年份来看，2005年的干旱程度最大，除黄土高原西部的部分地区，其他TVDI值均在0.76以上；其次为2007年，其空间分布特征与2005年相似，TVDI大范围表现为重旱；2019年的干旱程度仅次于2007年，尤其在内蒙古的银川等地区，干旱程度为重旱。TVDI最小的为2014年和2012年，TVDI值集中在0.7以下，干旱程度较低(图3)。

图3 2001-2020年黄土高原TVDI空间分布

2.2 TVDI变化趋势分析

基于Sen和M-K检验方法等到20 a间全年及不同季节的TVDI的变化趋势的空间分布、不同植被分区中不同变化趋势等级中所占面积百分比(图4)。其中，春季TVDI以增加趋势为主，增加速率集中在0～0.005/a，占黄土高原总面积的85.33%。TVDI呈不显著减少趋势的区域集中在青海西部、甘肃中部、陕西和山西的少部分地区，减少速率集中在0～0.005/a，占黄土高原总面积的14.67%；夏季TVDI变化速率除在山西少部分地区与春季不同，其他地区与春季的变化趋势相近；秋季98.21%的像元数表现为增加趋势，增长速率集中在0～0.005/a的像元占41.67%，56.54%的像元上增长速率大于0.005/a，呈不显著减小趋势的像元零星分布在陕西南部地区；20 a间冬季TVDI大部分像元上表现为减小趋势，即土壤湿度在增加，尤其在青海东部大部分地区、甘肃西南部少部分地区的减小速率大于0.005/a，呈增加趋势的像元零星分布在陕西南部少部分地区。从全年变化趋势空间分布来看，20 a间黄土高原TVDI以不显著增加为主，呈减少趋势的区域在青海中部、甘肃大部分地区、陕西中部、南部及山西有零星分布。

从不同植被分区下四季TVDI在变化速率分布区间数据发现(图4)：全年TVDI在分区Ⅰ和Ⅱ均以0～0.005/a变化速率增加，且像元百分比均小于20%，分区Ⅲ，Ⅳ和Ⅴ不显著减小和增加趋势的像元所占百分比各级基本为50%左右，不显著减小的像元数明显增多。从不同季节来看，春季的像元统计数与全年的较为相近；夏季呈减小趋势的像元数明显小于春季和全年；秋季TVDI整个以增加趋势为主，且不同植被分区的像元数的差异较小；冬季除Ⅴ区统计数与其他4区差别较大，其他4区相差不大，各个区呈增加趋势的像元占比均小于20%。

图4 2001-2020年黄土高原TVDI变化趋势及不同分区TVDI变化趋势

2.3 TVDI未来持续趋势分析

从黄土高原Hurst指数年变化、各季节变化持续性空间分布(图5)发现：黄土高原春季Hurst指数均值为0.41，表示将来TVDI变化特征与过去20 a黄土高原春季TVDI长期相关性特征为反持续性，即将来黄土高原春季TVDI将以减小趋势为主，在陕西境内的大部分像元上的这种持续性较其他地区更强；夏季Hurst集中在0.5以下的像元占比为61.76%，但在青海省西部、内蒙古的中部和北部、山西和陕西的少部分地区的Hurst集中在0.5以上，即这些区域的未来变化趋势与过去20 a的趋势相反，结合趋势分析发现青海省未来TVDI变大，其他几个地区未来TVDI均继续变小；秋季整个地区TVDI的Hurst指数为0.47，即秋季TVDI未来以变小趋势为主，呈增加趋势的主要集中在山西北部、内蒙古南部；冬季Hurst值为0.51，大部分像元变化趋势与过去20 a相反，即冬季TVDI未来呈增加趋势，在陕西中部地区的这种趋势更为明显。

图5 2001-2020年黄土高原Hurst指数空间分布

从全年来看，Hurst指数均值为0.45，即未来一段时间内的持续状态与过去20 a的变化呈相反趋势，空间分布上存在明显的空间异质性，存在明显的分区差异，在子植被分区Ⅰ，Ⅱ的Hurst集中在0.5以上，根据变化趋势结果表明在Ⅰ区域未来一段时间内，TVDI将呈减小趋势。Ⅲ，Ⅳ，Ⅴ三区的Hurst指数集中在0.5以下，即未来一段时间的变化趋势与过去20 a的变化趋势相反。

2.4 气温、降水对TVDI的影响情况

在控制气温的条件下，春季TVDI与降水以负相关性为主，占黄土高原总面积的85.41%，但只有0.21%的区域达到p<0.05的显著性检验。TVDI与降水呈正相关性区域集中在海拔较高的青海北部、陕西南部地区；夏季TVDI与降水以正相关性为主，青海西部、内蒙古大部分地区及山西北部少部分地区的相关系数集中在0.3以上，相关性较高。在海拔较低的陕西南部少部分地区，TVDI与降水呈极显著负相关关系。秋季TVDI与降水以负相关性为主，在陕西西部和甘肃西北部和东部的部分像元上相关系数达到-0.5以下，但在海拔较低的青海省西部、内蒙古地区的相关系数以正相关为主。冬季TVDI与降水主要以负相关为主，甘肃陇南、青海的青海湖区域少部分像元上表现为正相关关系。整年TVDI与降水以负相关性为主，仅在河南、陕西南部、山西东部部分像元上表现为正相关关系(图6)。

注：上为降水偏相关系数空间分布图，下为气温偏相关系数空间分布图。图6 TVDI与气温、降水偏相关系数

在控制降水的条件下春季气温与TVDI以负相关为主，即温度升高，TVDI减小，空间差异较小；夏季TVDI与气温以正相关关系为主，即温度升高会导致TVDI增大，尤其在陕西、河南、山西南部及甘肃中东部地区的正相性较大；秋季TVDI与气温在黄土高原大部分地区表现为正相关关系，只有在海拔较低的河南、陕西南部部分像元上呈负相关性；冬季TVDI与气温以负相关性为主，尤其在宁夏中部地区。就全年尺度来看，气温对TVDI主要为促进作用，即温度的升高会是TVDI增大，地表干旱程度加深，仅在海拔较高的内蒙古部分地区、甘肃陇南和青海青藏高原境内地区部分像元上表现为负相关关系。

由图7可知，就降水来看，不同植被分区TVDI与降水的相关系数与季节存在较为明显的关系，全年除子植被分区Ⅴ的相关系数为正，其他4区的TVDI与降水均为负相关关系，尤其是分区Ⅰ，Ⅱ的相关系数达到-0.2以下，负相关性明显，Ⅲ区相关系数为-0.14，Ⅳ区的相关系数均值为-0.008，相关性接近0，负相关关系不明显；春季TVDI与降水在5个区均呈负相关，除Ⅴ区TVDI与降水的负相关关系不明显，其他4区的相关系数集中在-0.2～-0.1；夏季TVDI与降水呈正相关关系，Ⅱ区的相关系数均值达到0.2以上，正相关关系明显，Ⅴ区TVDI与降水的关系不明显。秋季TVDI与降水均呈负相关关系，Ⅰ，Ⅳ区的负相关性较高，即TVDI随降水量的增多，明显下降，其他3区的相关性极低；冬季与全年的统计结果较为接近，但其相关性相较全年变低。

图7 黄土高原TVDI与季节性气候不同植被分区的平均相关系数分布

不同分区TVDI与气温的相关关系存在较大差异，全年TVDI与气温在5个区的相关系数均为正相关，尤其是Ⅳ区的相关系数均值达到0.24，相关性较高，其他4区的相关系数在0.1上下波动；春季TVDI与气温除在Ⅴ区的相关系数在0以上，其他4区TVDI与气温均为不显著的负相关关系；夏季在Ⅰ区的TVDI与气温的相关系数接近0，即气温的变化对TVDI的影响不大，Ⅱ，Ⅲ区的相关系数接近0.1，较Ⅳ，Ⅴ区的相关性(0.25以上)小很多；秋季除植被分区Ⅰ，Ⅴ区外，其他3区TVDI与气温的相关系数接近0.4，相关性极高；冬季TVDI与气温的相关性在不同植被分区间的关系为正负相关并存，但总体而言，TVDI与气温的相关性较秋季低很多。

3 讨论与结论

3.1 讨 论

地处半湿润、半干旱、干旱的黄土高原地区，因生态环境十分脆弱、干旱灾害严重、土壤侵蚀剧烈而成为世界上水土流失最严重的地区，且90%的黄河泥沙来自黄土高原，因此黄土高原的生态健康直接维系黄河的健康运行。因此，20世纪80年代以后，国家开展了一系列水土保持工程、退耕还林(草)工程、坡耕地整治工程控制水土流失，改善黄土高原生态环境建设。但随着植被覆盖度大幅度增加也使黄土高原部分地区土壤蒸发量加大，土壤干燥化程度加深[28,37-39]。本文通过研究2001—2020年黄土高原干旱时空分布特征、过去20 a呈现的趋势、未来一段时间可能持续的趋势及其对气温和降水变化的响应。通过研究近20 a黄土高原时空分布特征，发现黄土高原过去20 a整体为变旱趋势，干旱整体呈现出北部、东南高、西部、西南低的空间分布特征，最低值集中在黄土高原西部的青海省的西宁、门源、贵德、甘肃的临夏、临洮及东部山西的太原地区，属于轻旱、无旱级别，但在2005年、2007年、2019年大面积呈现为重旱，仅有海拔较低的青海没有出现重旱现象。最高值在黄土高原东南部河南和山西的长武等地区，属重旱级别。刘英等[39]通过TVDI指数对陕西省的干旱监测也发现2005年比相近年份的受旱面积大，也有报道显示2007年在内蒙古、宁夏、河南、陕西等地区受到严重的旱灾，且干旱范围较大，局地干旱非常严重，全国因旱受灾面积多于1991年以来的年均值，这一研究结果和报道均为本文的研究结果提供支持。

20 a间黄土高原TVDI以不显著增加趋势为主，即干旱程度在不显著增加，呈减小趋势的区域在青海中部、甘肃大部分地区、陕西中部、南部及山西有零星分布，未来一段时间内的干旱程度呈减弱趋势，但存在明显的分区差异。但春、夏两季干旱以增加趋势较秋季显著，冬季整体干旱程度减小。未来一段时间内黄土高原春季大部分区域的干旱程度有所缓解，夏季在大面积上表现为干旱缓解趋势，冬季未来一段时间内会在大范围变旱，尤其在陕西中部地区的这种趋势更为明显。降水量与TVDI呈负相关，仅在河南、陕西南部、山西东部部分像元上表现为正相关关系；气温对TVDI起促进作用，即温度的升高会使TVDI增大，地表干旱程度加深，仅在海拔较高的内蒙古部分地区、甘肃陇南和青海青藏高原境内地区部分像元上表现为负相关。究其原因可能是因为黄土高原位于东亚夏季风影响区域边缘，东南部降水受夏季风制约，东亚夏季风强度显著减弱，导致黄土高原东、南部干旱增加。

3.2 结 论

(1) 20 a间黄土高原TVDI整体呈现出北部、东南高，西部、西南低的空间分布特征，最低值集中在青海省的西宁、门源、贵德，甘肃的临夏、临洮及东部山西的太原地区，TVDI基本在0.57以下，即属于轻旱、无旱级别。最高值集中在黄土高原东南部河南和山西的长武等地区，TVDI值基本在0.76以上，属重旱级别。

(2) 20年间，黄土高原TVDI以不显著增加为主，在青海中部、甘肃大部分地区、陕西中部、南部及山西零星分布有呈减少趋势的像元。春季85.33%的面积上TVDI呈现出不显著的增加趋势，增加速率集中在0～0.005/a，夏季TVDI以不显著增加趋势为主；秋季98.21%区域为增加趋势，56.54%的像元上增长速率大于0.005/a，冬季TVDI大部分像元上表现为减小趋势。

(3) 黄土高原Hurst指数均值为0.45，即未来一段时间内黄土高原TVDI可能呈现出不显著的减少趋势，这说明黄土高原的干旱程度会有所缓解。春季未来一段时间内黄土高原大部分区域的TVDI将呈减小趋势，在陕西境内的大部分像元上的这种持续性较其他地区更强；夏季未来一段时间内青海省未来TVDI变大，其他几个地区未来TVDI均呈变小趋势；未来秋季TVDI将在一段时间内大部分像元数呈变小趋势，继续呈增加趋势的主要集中在山西北部、内蒙古南部的少部分地区；冬季未来一段时间冬季TVDI将呈增加趋势，尤其在陕西中部地区的这种趋势更为明显。

(4) 20年间，黄土高原TVDI与降水除在河南、陕西南部、山西东部少部分地区呈现出正相关，其他大部分地区表现为负相关，气温对TVDI主要为促进作用，但在海拔较高的内蒙古部分地区、甘肃陇南和青海青藏高原境内地区部分像元上表现为负相关关系。春季85.41%的面积上TVDI与降水之间的关系以负相关为主，相关系数主要集中在-0.5～0，与春季气温也呈负相关；夏季TVDI与降水以正相关为主，在海拔较低的陕西南部少部分地区与降水呈极显著负相关，与气温呈正相关；秋季TVDI与降水以负相关为主，在海拔较低的青海省西部、内蒙古地区的相关系数以正相关为主，在海拔较高的内蒙古部分地区、甘肃陇南和青海青藏高原境内地区部分像元上表现为负相关；冬季TVDI与降水主要以负相关为主，与气温以负相关为主。