三江源区多尺度水文干旱特征及植被的响应

肖祖香，朱 双，罗显刚，娄连惠，李 江，罗顺根，操 丽

(1.中国地质大学(武汉)地理与信息工程学院，湖北 武汉 430074；2.国家地理信息系统工程技术研究中心，湖北 武汉 430074； 3.中国地质环境监测院，北京 100081；4.湖北省自然资源厅信息中心，湖北 武汉 430071； 5.武汉中地云申科技有限公司，湖北 武汉 430074)

干旱事件是一种高频的自然灾害，具有最广泛的社会、经济和生态影响[1-2]。研究地区干旱演变特性对了解区域气候变化，指导未来水资源规划和利用具有重要意义[3-4]。Mckee等[5]提出的标准化降水指数(standardized precipitation index，SPI)被广泛应用于评估不同时间尺度的气象干旱事件，Sobral等[6]认为SPI在干旱时空特征研究中有非常好的作用。Wang等[7]采用归一化植被指数(normalized difference vegetation index, NDVI)来研究植被与干旱的关系。增强植被指数(enhanced vegetation index，EVI)是对NDVI的改进指数，减弱了大气和土壤噪音对植被的影响，反应植被情况更具稳定性。自20世纪90年代以来，三江源区气候呈暖湿现象，长江流域湿地面积呈整体上升趋势，黄河流域和澜沧江流域呈下降趋势[8]。三江源区的植被变化、气候变化和生态保护是重点研究方向，Wang等[9]研究发现草地退化比较明显的地区为黄河源区。目前大部分研究是针对2010年前后的地区干旱与植被演变，对三江源区的多尺度干旱演变与植被变化研究还很少[10]。

本文采用SPI和遥感植被指数(NDVI、EVI)量化1961—2018年三江源区的生长季植被和水文干旱的演化特性，基于Mann-Kendall(M-K)趋势检验法分析1961—2018年三江源区多尺度干旱的时空特征和变化趋势以及2000—2018年植被时空演变趋势，并探讨了三江源区的植被状态与干旱事件的联系，可为后期探究三江源区的植被-干旱的响应关系提供参考。

1 研究区域与数据

1.1 研究区域概况

三江源区位于中国青海省南部，地理位置在31°39′～36°12′N、89°45′～102°23′E之间，是长江、黄河和澜沧江的发源地，也是青藏高原的腹地(图1)，面积约30万km2，占青海省总面积的43%，地势高耸、气候寒冷，年平均气温在-5. 4～4.2 ℃之间，年降水量由东南部的 770 mm 向西北逐渐递减至 260 mm，气候条件具有明显的区域分异特征。三江源区中黄河、长江和澜沧江年平均径流量分别为201亿m3、124亿m3和150亿m3[11]，年蒸发量为730～1 700 mm，超过80%的降水量发生在5—10月[12]。由图1(b)可知2015年该地区的土地利用类型主要为林地和草地。

图1 三江源区的气象站点以及土地利用类型Fig.1 Geographic location, meteorological stations and land-use type of the TRSR region

1.2 数据来源与处理

三江源区数字高程模型(DEM)来自地理空间数据云网站(http://www.gscloud.cn/)，在ArcGIS中以三江源区边界来裁剪和镶嵌得到研究区的DEM。青藏高原土地利用覆盖数据来自中国科学院资源与环境科学数据中心(http://www.resdc.cn/)，通过裁剪和重分类得到三江源区的土地利用类型。

遥感卫星数据MOD13Q1来自美国国家航空航天局(NASA)EOS / MODIS数据产品(https://modis.gsfc.nasa.gov/)，空间分辨率为250 m×250 m，时间分辨率为16 d。MOD13Q1植被指数产品(NDVI和EVI)用水、云、重气溶胶等处理，以确保数据质量，这两种产品有效地表征了全球植被状态和过程。为了消除异常值的影响，采用最大合成法(MVC)合成月度NDVI数据，将月度NDVI数据合成为最大年NDVI影像[13]。

1961—2018年月降水量数据来自中国气象科学数据共享服务网(http://cdc.cma.gov.cn/)，其中1961—2016年数据来自中国0.5°×0.5°网格数据集(V2.0)产品，三江源区共有150个格点数据集。数据经过交叉验证，错误分析，质量状况良好；2017—2018年数据来自中国地面气候资料月值数据集。由研究区20个气象站点的地面月降水量数据，采取最优Kriging插值方法重采样2017—2018年150个格点的月降水量数据，得到1961—2018年的三江源区月降水数据。本文依据气象站点位置以及区域分布特点选取了40个格点数据进行研究。

2 研究方法

2.1 SPI

SPI是常用的多时间尺度气象干旱指数[5]，它可以计算多个时间尺度[14]。根据GB/T 20481—2006《气象干旱等级》的定义，SPI是表征某时段降水量出现概率的指标，适用于月尺度以上的干旱监测与评估。该指数假设降水量服从Γ分布，考虑了降水服从偏态分布的实际，然后进行正态标准化处理，使得同一个干旱指数可以反映不同时间尺度和不同类型的水资源状况。一般当SPI值小于0且连续观察到该值小于等于-1时，定义为发生干旱事件[15]。SPI1 通常用来反映短期地表水分异常变化，SPI3可反映干旱的季节变化。表1为中国的干旱分类标准[16]。SPI的详细介绍和计算机程序可在国家干旱减灾中心网站(https://drought.unl.edu/)上获得。

表1 基于SPI指数的干旱等级划分

2.2 植被指数

NDVI是目前应用最广泛的植被指数，可以反映植被的覆盖度、植被的基本生长状态等。EVI是增强型植被指数，相对于NDVI在大气噪声、土壤背景、饱和度等问题上有所改进[17]。本文根据三江源区植被生长的实际情况，确定干旱对植被生长状况的影响主要表现在植被生长季(5—9月)，提取每年5—9月的NDVI和EVI数据进行植被时空演变分析。

2.3 M-K趋势检验法

M-K非参数检验法常用于评估降水、径流、气温等数据时间序列的趋势变化，它的优点在于样本不需要遵循某一特定的分布，而且很少受到异常值的干扰，且计算简便[18]。该方法以标准正态分布变量Z表示长时间序列SPI和植被指数的变化趋势，Z>0表示序列有上升趋势，Z<0表示序列有下降趋势。考虑到正态分布的对称性，Z的绝对值在大于1.28、1.64、2.32时，分别表示通过了显著性水平α为0.1、0.05和0.01的显著性检验。本文利用MATLAB软件，编写M-K趋势分析的循环计算程序，得到所有格点的Z值。

2.4 斯皮尔曼等级相关系数

斯皮尔曼等级相关(Spearman rank correlation)在统计学中用于评估两个变量之间的相关性[19]，变量间的相关性可以使用单调函数来描述[20]。当变量x增大时，变量y也趋向于增大，则相关系数ρ>0；当x增大时，y趋向于减小，则ρ<0；当x增大时y没有任何趋向性，则ρ=0。斯皮尔曼方法估计的相关系数更合理且误差更小，较之其他如Kendall方法更具优势[21]。

3 结果与分析

3.1 干旱时空特征

利用三江源区的植被生长季(5—9月)的SPI1来表示生长季月度干旱，春季(3—5月)、夏季(6—8月)、秋季(9—11月)和冬季(12—2月)的SPI3极小值来表示四季干旱。

3.1.1 SPI时间演变

图2为1961—2018年三江源区的长时间序列的干旱特征SPI值(5—9月)，SPI1和SPI3线性拟合斜率分别为0.004和0.005，说明在此期间三江源区逐渐变得湿润。月度严重干旱事件(SPI1)发生在1963年、1966年、1969年、1984年、1997年和2017年，季节性中度干旱事件(SPI3)发生在1963年、1966年、1984年、1992年、2000年和2017年。

图2 1961—2018年SPI1和SPI3时间序列Fig.2 Time series of SPI1 and SPI3 from 1961 to 2018

3.1.2 SPI趋势空间分布

图3为1961—2018年SPI在月度与季度的M-K趋势检验值Z的空间分布。5月、9月、春、夏、秋和冬季平均趋势特征值Z分别为1.64、0.63、2.28、1.58、0.65和1.03，其中SPI上升趋势占比分别为92.68%、70.00%、97.50%、90.00%、57.50%和82.50%。可知三江源区呈月度湿润化趋势，5月在长江西北部有显著湿润化趋势，但9月在长江与黄河南部出现不显著的干旱化趋势(Z<0)。三江源区的季节性SPI有上升趋势，不同流域差异大。春季平均Z=2.28且通过了α=0.05的显著性检验，表示1961—2018年三江源区春季有显著的湿润化趋势。

图3 1961—2018年SPI的M-K趋势检验的空间分布Fig.3 Spatial distribution of M-K trend test value Z of SPI from 1961 to 2018

3.2 植被时空特征

图4是每年植被生长季的NDVI和EVI像元最大合成值(年最大值(5—9月)、5月和9月的值(共计19年))，得到三江源区的植被时间变化特征。由图4可知NDVI值大于EVI值。植被指数在5月的突变点发生在2001年和2017年，年际的突变点发生在2011年，9月的突变点发生在2002年和2013年。2010年以来，年际植被指数表现出植被先增后减再增现象。通过相关性分析得到两个指数呈高度正相关性(相关系数大于0.9，通过α=0.01显著性检验)。

图4 植被指数的时间变化Fig.4 Time evolution of vegetation index (NDVI/EVI)

图5为2000—2018年期间植被指数在5月、9月和年际的M-K趋势检验值Z的空间分布。5月、9月、年际的NDVI平均Z值分别为1.54、0.82和1.36，NDVI显著上升趋势的占比分别为50%、72.5%和55%；EVI的平均Z值分别为1.44、0.48和1.01，显著上升趋势的占比分别为52.5%、70%和65%。5月和年际的植被指数显示，黄河以北地区有显著的植被增加趋势；9月植被指数和年际植被指数显示，长江中部与澜沧江植被有不显著的下降趋势。2000—2018年三江源区植被变化趋势为北部强于南部、西部强于东部,该地区大部分区域植被状况有好转。

图5 2000—2018年植被指数M-K趋势检验值Z的空间分布Fig.5 Spatial distribution of trend test value Z of vegetation index from 2000 to 2018

3.3 相关性分析

表2为研究区不同SPI与不同植被指数的斯皮尔曼等级相关系数ρ及其显著性水α(表中NDVI5、EVI5、NDVI9和EVI9分别为5月和9月的植被指数，NDVImax和EVImax为5—9月植被指数的最大值。

表2 不同SPI与不同植被指数的斯皮尔曼相关系数

由表2可知5—9月SPI3的平均值与NDVI5、EVI5和NDVImax均有显著正相关性, NDVImax的相关系数最大(ρ=0.584)，其次是EVI5和NDVI5；整体来看5—9月SPI3的平均值与植被指数呈正相关关系，且NDVI的相关系数明显大于EVI。由表2可知5—9月SPI3的极小值与NDVI5、EVI5显著正相关，NDVI5的相关系数最大(ρ=0.525)，其次是EVI5；整体来看5—9月SPI3的极小值与植被指数呈正相关关系，且NDVI的相关系数略大于EVI。综上可知，三江源区的SPI3与植被指数有显著的正相关性且与NDVI的相关性更高。

4 讨 论

对于长时间的植被]检测，通常可以使用卫星遥感数据源。NDVI与植被变化之间有很好的关系，EVI则有空间相对较高、大面积的近实时数据、分辨率高的优点，已经在植被状况监测方面得到广泛地使用[22]。Ji等[23]研究表明SPI3与NDVI的相关性最高，Haro-Monteagudo等[24]提到3个月时间尺度的干旱指标与植被变化的相关性最高，这些结论与本文结果一致，表明植被变化与干旱事件之间有一定的联系和互相影响[25]。本文发现三江源区西北部和东北部大部分地区的植被呈显著增加趋势，这与该地区的干旱变化趋势一致。植被增加区域主要分布在长江源西北部，其中澜沧江源北部和黄河源西部某些区域有植被退化迹象，这与钱拴等[26]的研究结论一致。本文研究发现21世纪以来，三江源区整体呈湿润化趋势，但同时监测到2017年存在月度和季度干旱事件，以及黄河源西南部地区也在一定的植被度减少趋势。这也推动着我们今后探究不同时间尺度干旱特征变量以及其他气象因素对植被的影响，探讨更有效的方法来量化和解释该地区植被-干旱的响应关系。

5 结 论

a.1961—2018年不同时间尺度的 SPI均体现三江源区在逐渐变得湿润，时间尺度越大，干湿交替周期越长。21世纪以来，三江源区干旱事件出现在2002年、2006年、2012年、2015年和2017年，秋冬季节发生干旱的频率更高，其次是春季。

b.长江源西北地区有显著湿润化趋势，长江源与黄河源南部在秋冬季节有不显著的干旱化威胁，春季有显著的湿润化趋势，整体呈湿润化趋势。

c.三江源区的植被在21世纪以来有所增加，在2002年和2013年增加量出现减少现象；NDVI与EVI在三江源区有明显增加的趋势，空间上变化趋势是北部强于南部、西部强于东部。

d.三江源区的SPI3与植被指数均呈正相关关系: 5—9月SPI3的平均值与NDVImax有最显著的正相关性，5—9月SPI3的极小值与NDVI5、EVI5有显著正相关性。