黄河源植被时空变化及其对地形和气候的响应

覃巧婷,陈建军,2*,杨艳萍,赵晓宇,周国清,2,尤号田,2,韩小文,2 (.桂林理工大学测绘地理信息学院,广西 桂林54004；2.桂林理工大学广西空间信息与测绘重点实验室,广西 桂林 54004)

黄河源区位于青藏高原腹地,是黄河流域的主要产流区、水源涵养区,也是我国重要的生态屏障[1-2].黄河源区植被以高寒草甸、高寒草原和高寒灌丛为主,其草地生态系统是长期自然演化的结果,生态环境十分脆弱,容易因外界的干扰和破坏而发生变化,恢复难度极大且恢复过程缓慢[3].近几十年来,在气候变化和人类活动的影响下,黄河源区的生态环境发生了显著的变化,主要表现为高寒草地严重退化[4]、土地沙漠化[5]、生物多样性和数量锐减[6]、水土流失加剧[7]等,引起了众多学者和决策部门的广泛关注[8-10].作为青藏高原陆地生态系统的重要组成部分,黄河源区植被对维护黄河流域生态环境安全及可持续发展具有非常重要的作用[11].因此,研究黄河源区植被的变化及其驱动机制具有重要意义.

利用遥感卫星数据监测植被的生长有助于更好地了解植被信息.当前用于度量植被变化的植被指数有很多,例如归一化植被指数(NDVI)、增强型植被指数(EVI)、土壤调节植被指数(SAVI)、差值植被指数(DVI)等,但NDVI是应用最广泛的植被指数[12].植被变化监测的方法主要有两种,一种是站点长期监测,一种是遥感监测.站点监测精度较高,但空间代表性有限;遥感监测能够较精确地反映地面上植被的变化情况,是当前研究区域植被变化的主要技术手段[13].近年来,许多学者利用 NDVI时间序列从月尺度[14]、年尺度[15]、像元尺度和区域尺度[16]等不同方面对植被变化进行动态监测,认为全球气候变化导致植被发生了显著的变化.影响植被变化的气候因子主要有气温、降水和辐射等,而这些气候变化主要是通过改变植被生长所需要的能量和水分,影响碳的积累过程、水循环过程和土壤有机碳的分解转化过程,从而影响植被的生长和分布格局[17].目前关于植被与气候变化的研究相对较多,大多数学者采用气温、降水来研究植被对气候变化的敏感性[18-20],指出与人类活动、地形因素等相比,气候变化对植被的生长起决定性作用.以往的研究大多讨论了整体植被与气候变量之间的相关性,但气温、降水对植被的影响存在空间异质性,且不同植被类型对气候变化的响应特征存在差异[21].同时,相关研究表明植被在垂直方向上的变化是不同的,高海拔地区气候变化对植被生态系统的影响更大[22].黄河源区位于青藏高原东北部,其平均海拔为4000m左右,地形起伏较大,地理条件复杂[23].而目前关于不同地形条件下气候驱动影响的植被分布及演化的研究较少,因此,综合考虑地形因子和气候因子对不同植被类型的影响,对于更好地了解黄河源区植被动态具有重要意义.

本文基于 2000～2019年的时序 MODIS NDVI数据,采用一元线性回归及趋势分析等方法,分析黄河源区植被的时空变化特征;同时,结合地形数据和气候数据,利用相关分析和统计分析等方法,探讨植被对地形和气候因子的响应机制,以期为黄河源区生态环境保护提供理论依据.

1 材料与方法

1.1 研究区概况

黄河源区位于青藏高原东北部,东经 95.5°～103.5°,北纬 31.5°～36.5°,涉及青海、四川、甘肃 3 省的6个州、18个县,面积约为12.2万km2[24-25],占黄河流域总面积的 16%[26].黄河源区海拔约在 3000m以上,地势西高东低,地貌复杂,分布有高山、盆地、峡谷、草原、沙漠和众多的湖泊、沼泽、冰川及多年冻土等.该区属于典型的高原大陆性高寒气候,冷热两季交替,干湿季节分明,日照时间长,辐射强烈[27].黄河源区多年平均气温为-4.0～5℃,自西北向东南逐渐降低,海拔越高气温越低;多年平均降水量为220～780mm,降水年内分配不均,年际变化大,表现为冬干秋旱、夏秋降水集中[28].区内高寒植被分布广泛,主要有高寒草甸、高寒草原、高寒灌丛和高寒沼泽等,其中高寒草甸和高寒草原约占总面积的 80%.黄河源区主要以畜牧业为主,产业结构单一,生产力不发达[29].

1.2 数据来源与处理

1.2.1 NDVI数据 本研究使用的 MODIS NDVI数据是中分辨率成像光谱仪植被指数产品MOD13Q1,下载自美国国家航空航天局(NASA)(https://lads web.modaps.eosdis.nasa.gov/search/),空间分辨率为 250m,时间分辨率为 16d.采用 MODIS重投影工具 (MRT)对原始数据进行拼接、重采样、格式转换等处理,并利用黄河源区矢量边界对处理后的数据进行裁剪.为了消除云层、大气和太阳高度角等的影响,将每年 NDVI数据采用最大值合成法进行合成.

1.2.2 气象数据 采用气温和降水作为气象资料,气象数据来自欧洲中期天气预报中心(ECMWF)第三代再分析资料 ERA-Interim (https://apps.ecmwf.int/datasets/data/interim-full-daily/levtype=sfc/),空间分辨率为0.5°.ERA-Interim提供了自1979年以来的再分析资料,并实时更新.相比第二代产品 ERA-40数据同化模型,ERA-Interim采用了四维变分分析(4D-Var),并结合改进的湿度分析,卫星数据误差校正等技术,实现了再分析资料质量的提升[30].采用克里金插值法对气象数据进行空间插值,将气象数据重采样为空间分辨率为250m的栅格数据,实现气象数据与NDVI数据空间分辨率的统一.

1.2.3 辅助数据 DEM 数据来源于地理空间数据云(http://www.gscloud.cn/),空间分辨率为 30m,地理坐标为GCS_WGS_1984,对DEM数据进行拼接、裁剪、重投影的处理,获得黄河源区的高程数据,并提取坡度数据,研究黄河源区植被随高程、坡度的变化特征.植被类型数据来源于中国科学院资源环境科学数据中心(http://www.resdc.cn/),其植被类型主要包括针叶林、阔叶林、灌丛、荒漠、草原、草甸、沼泽、高山植被、栽培植被等,本文仅选取主要的3种植被类型(高寒草甸、高寒草原、高寒灌丛)进行分析,其空间分布如图1所示.

图1 黄河源区主要植被类型空间分布Fig.1 Spatial distribution of main vegetation types in the source region of Yellow River

1.3 研究方法

1.3.1 最大值合成法 采用最大值合成法获取黄河源区年NDVI的最大值,从而去除残云、残雾、大气等的影响.计算方法如下所示:

式中:i为年序号,取值范围为2000～2019;NDVIi为第i年的最大 NDVI值;NDVI1、NDVI2、NDVI3、NDVI4分别表示每年7月份上半月和下半月、8月份上半月和下半月的NDVI数据.

1.3.2 NDVI变化趋势 采用线性回归趋势分析方法计算植被 NDVI的变化趋势.多年回归方程中的趋势斜率代表年际变化,通常利用最小二乘法求解.计算公式如下:

式中:θslope为 NDVI的回归斜率;n为研究时间段累计年数,本文中n=20;i代表年序号;NDVIi为第i年平均 NDVI值.若θslope＞0,表明 NDVI在 n 年间呈增加趋势;θslope＜0,表明 NDVI在 n年间为下降趋势;θslope=0,表明NDVI在n年间无明显变化[31].

1.3.3 NDVI与气温、降水的偏相关分析 偏相关分析是探索植被生长与单一气候因子之间的联系,同时剔除其他气候因子的影响,被广泛用于检测影响植被生长的主要气候驱动因子[32].计算公式如下:

式中:x为 NDVI;y/z为气温/降水,单位为℃/mm;rxy,z为将降水/气温固定后气温/降水与 NDVI的偏相关系数.若 r＞0,表明气温/降水与 NDVI呈正相关;r＜0表明气温/降水与 NDVI呈负相关;r=0表明气温/降水与NDVI不存在线性相关性.rxy、rxz、ryz分别代表NDVI与气温、NDVI与降水、气温与降水的相关系数.其中,相关系数的计算公式如下:

式中:x为 NDVI;y为气温/降水;rxy为气温/降水与NDVI之间的相关系数;xi与yi分别为NDVI、气温/降水第 i年的值;¯x和¯y分别表示 NDVI、气温/降水两变量n年的平均值;n为样本数,本文中n=20.

1.3.4 NDVI与气温、降水的复相关分析 复相关分析是研究一个变量与多个变量之间的相关程度.植被是受到多种因素影响的,因此运用复相关分析计算出NDVI与气温、降水的复相关系数,以此分析NDVI受两者共同作用的影响有多大.复相关系数计算公式如下所示:

式中:x为NDVI;y为气温;z为降水;rx,yz为NDVI与气温、降水的复相关系数.

复相关系数的显著性检验采用 F检验法,其统计量计算公式如下所示.

式中:n为样本数,本文中n=20;k为自变量个数.

1.3.5 NDVI与高程、坡度的关系 为了探究NDVI随海拔梯度的分布和变化特征,将黄河源区分为 37个 100m(2564～6295)的高程带,统计每一高程带的平均 NDVI值和所占的像元比例,其中整个研究区共有 156056034个像元.同理,将黄河源区分为30 个 2°(0～60°)和 1 个 30°(60～90°)的坡度带,统计不同坡度带所占的像元比例和 NDVI值,分析 NDVI随坡度的变化特征.

2 结果与分析

2.1 植被NDVI的时空变化特征

2.1.1 植被NDVI的时间变化特征 如图 2所示,2000～2019年黄河源区植被NDVI呈波浪式的增加趋势,线性增长率为 0.0014/a.在 2003年植被 NDVI出现最小值,在2018年植被NDVI出现最大值.黄河源区各植被类型NDVI的年际变化与植被整体基本一致.其中高寒草原的增加速率最大,为 0.0027/a;高寒草甸次之,为 0.0013/a;高寒灌丛的增加速率最小,为0.0008/a.表明近20a来研究区高寒草原的改善趋势更为明显.在数值上,3种植被类型NDVI值从大到小依次为高寒灌丛、高寒草甸、高寒草原.其中高寒草甸与植被整体的变化趋势最为接近.

图2 黄河源区植被NDVI的年际变化Fig.2 Annual variation of vegetation NDVI in the source region of the Yellow River

2.1.2 植被 NDVI的空间分布格局 将黄河源区植被 NDVI划分为＜0.2、0.2～0.4、0.4～0.6、0.6～0.8、＞0.8 5个等级.如图3所示,2000～2019年黄河源区平均NDVI为0.65,表明多年植被NDVI处于中高水平,但不同地区差异较大,呈现出由东南向西北递减的变化趋势.其中中高植被覆盖区主要分布在东南部的红原县、若尔盖县、玛曲县、久治县等区域;低植被覆盖区主要分布在西北部的玛多县、曲麻莱县的部分区域.

图3 黄河源区植被NDVI的空间分布Fig.3 Spatial distribution of vegetation NDVI in the source region of the Yellow River

2.1.3 植被NDVI的变化趋势 对黄河源区2000～2019年植被NDVI做一元线性回归分析,根据回归显著性检验结果(P＜0.01),将黄河源区植被 NDVI变化情况划分为显著减少、无显著变化和显著增加3类.

如图4所示,总体上植被NDVI显著增加的面积大于植被 NDVI显著减少的面积.分析其原因,一方面是近几十年来实施的生态保护与建设工程取得的成效,例如围封、退牧还草、湿地保护等;另一方面是降水增多、气温升高等使得植被生长的气候环境变好.其中植被 NDVI显著增加的面积占源区总面积的 20.42%,主要呈块状分布在黄河源区北部的兴海县、同德县和泽库县的部分区域以及扎陵湖和鄂陵湖周边地区;植被 NDVI显著减少的面积占研究区总面积的 3.14%,主要零星分布在东南部曲麻莱县和南部达日县和久治县等部分区域;植被NDVI基本不变的面积占研究区总面积的 76.44%,分布于黄河源区大部分区域.

2.2 NDVI随地形因子的变化特征

2.2.1 NDVI随海拔的变化特征 黄河源区地势西高东低,最高海拔达到6295m,高海拔地区主要分布在黄河源区西北部,低海拔地区主要分布在东南部.如图 5所示,黄河源区的高程主要分布在3400～4700m之间,其所占比例达到整个研究区的92.26%.植被 NDVI随着高程带的不同而存在一定的差异,在高程小于3500m范围内,植被NDVI随着高程的增加而增加;高程在 3500～4100m 之间,植被 NDVI处于一个较平缓的趋势,之后随着高程的增加迅速下降;高程在 5300m 以上时,植被 NDVI极小,接近于0.

图5 黄河源区植被NDVI随海拔的变化关系Fig.5 The relationship of vegetation NDVI with altitude in the source region of the Yellow River

2.2.2 NDVI随坡度的变化特征 由图6可知,研究区内坡度小于20°的区域所占比例较大,占研究区总面积的 72.42%,主要分布在黄河源区的西北部和东南部.植被 NDVI随着坡度的增加呈现出单峰型曲线,在 24～26°坡度带达到峰值.当坡度小于 24°时,植被 NDVI随着坡度的增加而逐渐增大;当坡度大于26°时,植被NDVI随着坡度的增大而逐渐减小.

图6 黄河源区植被NDVI随坡度的变化关系Fig.6 The relationship of vegetation NDVI with slope in the source region of the Yellow River

2.3 气候因子对植被NDVI的影响

2.3.1 植被NDVI与气温、降水的时间相关性 植被NDVI与气温、降水时间序列的相关性采用皮尔逊系数来表示,如表 1所示,黄河源区植被 NDVI与气温、降水的相关性均显著,相关系数分别为0.5109和0.5048,均通过了0.05的显著性检验,表明气温、降水对植被的生长起促进作用.各植被类型与气温、降水均表现出正相关关系,其中气温与高寒草甸的相关性最好,且通过了0.05的显著性检验;降水与高寒灌丛的相关性最好,高寒草甸次之;高寒草原与气温、降水的相关性相对较差,其与降水的相关系数通过了 0.05的显著性检验,与气温的相关系数不具有统计上的显著性.

表1 黄河源区植被NDVI与气温、降水的相关系数Table 1 The correlation coefficient between vegetation NDVI and temperature and precipitation in the source region of the Yellow River

2.3.2 植被NDVI与气温、降水的空间相关性 由于不同地区特殊的自然条件,植被对气温、降水的响应程度各不相同.由图 7可知,黄河源区植被 NDVI与气温、降水在空间上呈现出以正相关为主的趋势.NDVI与气温、降水的最大偏相关系数分别为0.922和0.907.植被NDVI与气温呈现正相关的区域主要分布在黄河源区西北部,表明在这些区域气温是影响植被变化的重要因素;与气温呈现负相关的区域主要分布在黄河源区东南部.植被 NDVI与降水呈现正相关的区域主要分布在扎陵湖、鄂陵湖西北部及龙羊峡水库周边,与降水呈现负相关的区域主要分布在研究区中部和东南部部分地区.

图7 黄河源区植被NDVI与气温和降水的偏相关系数Fig.7 Partial correlation coefficient between vegetation NDVI and air temperature and precipitation in the source region of the Yellow River

除了单一的气候因子对植被的影响之外,植被的变化还受到多种因素共同作用的影响.利用复相关公式对气温、降水进行计算并用F检验法进行显著性检验,分析气温、降水共同作用下植被的变化.由图8可知,黄河源区植被NDVI与气温、降水的复相关系数介于0～0.93之间.整体上,植被NDVI与气候因子复相关性较强的区域主要集中在黄河源区西北部的曲麻莱县、玛多县、达日县;复相关性较弱的区域分布在龙羊峡水库周围及若尔盖县等海拔较低的区域.

图8 黄河源区植被NDVI与气温、降水的复相关系数Fig.8 The multiple correlation coefficient between vegetation NDVI and temperature and precipitation in the source region of the Yellow River

3 讨论

3.1 植被NDVI动态变化

本文基于MODIS NDVI数据,采用最大值合成法、趋势分析法研究黄河源区植被的时空变化特征.研究结果表明,2000～2019年黄河源区植被NDVI呈增长趋势,与刘启兴[24]、Yang[33]等的研究结论相一致,表明近年来黄河源区生态环境在不断改善.植被的变化是气候因子和人类活动共同作用的结果[34],对于生态系统较为脆弱的敏感区域,气候的突变或人类活动的加剧,都会使植被的变化出现较大的波动[35].黄河源区地处高海拔地区,是我国生态环境脆弱的地区之一.近年来由于人类活动加剧和全球气候变暖,区域内植被变化在不同年份间出现了较大的波动性,整体上植被呈增加趋势.但植被的改善主要依赖于局部地区的增长,而不是整个区域.其中植被 NDVI增加的区域主要分布在研究区北部,例如泽库县、兴海县和玛多县等;植被NDVI减少的区域主要分布在研究区南部和东南部,例如红原县、若尔盖县、达日县和甘德县等.除了气候因素,植被变化与人类活动也有着不可分割的联系.黄河源区在过去很长一段时间一直处于严重超载放牧的状态,过度放牧是引起该区草地生态系统退化的主要原因

[36].过度放牧会破坏黄河源区脆弱的草地生态系统,降低其保水性能,从而导致草场退化,地表裸露,出现沙砾地.近年来黄河源区植被的好转可能与多种生态保护与建设工程有关,例如三江源自然保护区、退耕还草、生态移民等.

3.2 植被与气候因子的相关性

气候变化对草地生态系统的影响主要是通过降水和温度的变化,改变土壤湿度和温度,同时影响光合作用和植物呼吸作用,进而调节草地的生长和生态系统生产力[37-38].在全球变暖的大背景下,对植被变化与气候因子的研究呈现多样性,尤其是在年尺度上.通过对黄河源区过去20a气温和降水的变化曲线进行分析(图 9),发现近 20a间研究区气温和降水整体上均呈增加趋势,线性增长率分别为 0.0344℃/a和 0.4949mm/a,其中气温增长速率与冯晓莉[39]、刘光生[40]等的研究结果一致,降水增长率低于两位作者的研究结果,可能是研究区域范围以及研究时间段不一致导致的结果.植被 NDVI随气温和降水上下波动,与降水相比,NDVI与气温的变化曲线走势较为一致,表明植被对气温的响应更强烈.

图9 黄河源区植被NDVI与气温和降水的变化曲线Fig.9 Variation curve of NDVI, temperature and precipitation in the source region of the Yellow River

经过计算可知,植被NDVI与气温、降水均呈线性正相关关系,与气温、降水的相关系数均通过了0.05的显著性检验(表 1).不同植被类型对气候变化的响应存在差异,其中高寒草甸受到气温的影响较大,高寒灌丛对降水的敏感性强于气温.高寒草甸主要分布在较高海拔地区,气温较低,降水量较大,植被的生长主要受气温的影响,所以高寒草甸对气温的敏感性更强.高寒灌丛主要分布在较低海拔地区,气温高,蒸发量大,可能由于水分条件的限制,植被的生长受到影响,所以其生长更依赖于降水[41].结合相关分析,发现近20a来黄河源区植被NDVI与气温、降水均呈现出以正相关为主的趋势,但与降水相比,气温对黄河源区植被变化的影响更为显著,这与相关研究结果一致[42-43].在一定程度上,温度的升高可以增强植物的光合强度,延长植物的生长期[44-45],从而促进植被的生长.因此,黄河源区植被覆盖的提高可能是气候变暖的结果.

3.3 植被与地形因子的相关性

相关研究表明,海拔高度对植被覆盖有一定的影响,植被覆盖与温度的关系随着海拔的变化具有不同的特征[46-47].在本文的研究结果中,植被 NDVI随着海拔高度的增加呈现出单峰型曲线,峰值出现在 2700～2800m 和 3500～3600m.植被 NDVI较大的区域主要集中在3100～4800m的中高海拔地区.在高程＜3600m 的范围内,植被 NDVI与高程呈正相关,NDVI随着海拔高度的增加而增加.海拔在3600～5300m的范围内,植被NDVI随着高程的增加呈现出减小的趋势;当海拔＞5300m时,植被NDVI趋于稳定,且数值较低,这可能反映了气候变化条件下高寒植被对环境的响应.分析其原因,可能是在低海拔地区,地势相对平缓,水热条件优越,人类活动较为频繁,严重影响了植被的恢复和生长,使植被覆盖降低.而高海拔地区人口稀少,人类活动相对较少,从而黄河源区植被呈现出增长的趋势,并达到一个峰值.而后,随着海拔的进一步升高,地势陡峭,降水量增加,气温下降,光照强烈,风速大,蒸散迅速,植被的生长受到限制,因此植被呈减少的趋势.

坡度是影响地表物质流动的重要因子,直接影响物质交流与能量转换的方式和程度,改变土壤的基本属性,从而影响植被的分布态势[48].本文的研究结果表明,随着坡度的增大,植被 NDVI呈现出先增大后减小的变化趋势,且在 24～26°坡度带 NDVI达到最大值(图6).相对于坡度较大的区域,坡度平缓更有利于植被的生长,但黄河源区大部分坡度相对较缓的地区受到的人为干扰相对较多,所以水热条件较好的缓坡地区其植被覆盖反而较低.坡度影响地表径流和排水状况,所以坡度越高的区域自然降水就越容易流失,土壤矿质营养等稀薄,植被生长差,植被覆盖也较低.

3.4 本研究的局限性

植被的动态演变受到气候因素、地形因素、土壤因素、人类活动等多个方面因素的共同影响,作用机制复杂.本文由于缺少相关资料,只考虑了地形因子和气候因子对黄河源区植被变化的影响,尚未考虑到其他的因子,若能把影响植被变化的其他因子共同进行研究,则更能深刻地揭示各种因素对植被生长的作用机理.

此外,目前NDVI时间序列数据主要有NOAA/AVHRR NDVI、MODIS NDVI、SPOT/VGT NDVI等,而研究植被变化及其影响机制的方法也有许多种,例如趋势分析法、相关分析法、残差分析法、Hust指数法、统计分析法等.本研究使用的数据为MODIS NDVI,采用的方法为趋势分析法、偏相关分析法、复相关分析法,研究结果因使用的数据、方法不同可能存在一些差异.

4 结论

4.1 黄河源区植被 NDVI整体呈增加趋势,线性增长率为0.0014/a,波动性较大.多年平均NDVI为0.65,但空间差异显著,呈现出由东南向西北递减的空间分布格局.NDVI显著增加的区域主要呈块状分布在黄河源区北部及扎陵湖、鄂陵湖周边地区,NDVI显著减少的区域主要零星分布在研究区东南部及南部的部分区域.

4.2 高海拔地区主要分布在黄河源区西北部,低海拔地区主要分布在东南部.随着高程的增加,植被NDVI呈现出先增加后减少的趋势,但在 3500～4100m区间 NDVI变化不显著.随着坡度的增大,植被 NDVI呈现出单峰型曲线,并在 24～26°坡度带达到最大值.

4.3 2000～2019年气温和降水整体上均呈上升趋势.NDVI与气温、降水的相关性均通过了0.05的显著性检验,表明气温、降水对植被的生长起促进作用.与降水相比,NDVI对气温的敏感性更强.