基于NDVI的1982—2012年黄河流域多时间尺度植被覆盖变化

颜 明，贺 莉，王随继，郑明国，孙莉英，许炯心

(中国科学院地理科学与资源研究所陆地水循环及地表过程重点实验室,100101,北京)

陆地生态系统中的植被不仅在全球物质与能量循环中起着重要作用,而且影响着全球碳的收支平衡,稳定全球气候等方面具有重要作用[1],对于流域来说,同样具有重要意义,不仅维持流域的生态平衡,而且对于流域的水文及地貌过程起着制约和调节作用。NDVI(normalized difference vegetation index,归一化差分植被指数)是当前常用的表征植被覆盖状况的指标,它与植被的一系列指数密切相关,如植被覆盖度、生物量、净初级生产量、叶面积指数等。近30年来,国内外学者利用NDVI数据在不同的空间尺度上对植被年际变化进行了深入研究[2-5]。研究发现:全球植被的覆盖在增加,以北半球中高纬地区最为明显,在我国的青藏高原及北方等地区植被覆盖也明显增加[6-7]。在月和季节尺度上关于NDVI的变化也有一些研究开展[8]。这些研究结论大多是针对单一时间尺度,没有分析不同时间尺度间所表现出来的植被变化是否具有一致性或存在抵消关系；因此,需要更为细致的研究,发掘NDVI在月、季和年3个时间尺度之间的趋势性变化是否存在一致性和抵消关系,以确定所计算的更大的时间尺度趋势性是否能代表更小时间尺度的趋势性。

黄河流域作为我国的第2大流域,水沙问题十分突出,黄河上游水源涵养和中游的水土流失是黄河长期以来备受科研工作者和广大人民关注的现实问题,植被一方面能缓冲径流的汇聚过程,同时能保护地表阻止降雨—径流的冲刷能量减少土壤流失。自新中国成立以来,黄河流域的植被建设工作从未停止,尤其是黄河中游为防治水土流失开展的一系列水土保持工作；因此有必要检测和评估这些措施对黄河流域的植被覆盖的影响[9-11],发现过去工作中的不足,为提出新的治理措施奠定基础。而且对黄河流域的植被覆盖在不同时空尺度上开展研究,有助于揭示黄河流域植被对气候变化和人类活动的响应[12],为针对水土流失和水资源紧缺等问题制订规划时提供植被覆盖相关的资料,为流域生态文明建设提供空间信息和理论支撑。

笔者以黄河流域花园口水文站以上为研究区域,利用反映植被覆盖状况的NDVI数据分析黄河流域1982—2012年不同时间尺度NDVI的趋势性,主要分析以下3个内容:1)不同时间尺度下NDVI的空间特征,包括年和不同季节的空间分布情况,并对比它们之间的异同；2)分析年、季和月3个时间尺度NDVI的趋势性以及不同时间尺度趋势性的空间特征；3)通过相关矩阵和空间对应关系来探讨不同时间尺度植被覆盖趋势性变化之间的相互关系,分析不同时间尺度间植被覆盖的趋势性存在的一致性和差异。

1 研究区概况

黄河流域西起巴颜喀拉山,东临渤海,南至秦岭,北抵阴山,大部分区域位于中国西北部,由一系列平均海拔在4 000 m以上高山组成[13]；北部为内蒙古高原,海拔在1 300～2 200 m之间,高原北侧是阴山山脉,西侧是贺兰山,中间由一些列断陷盆地和鄂尔多斯高原构成；东部地区为黄土高原,海拔在1 000～2 000 m之间,地貌上主要是黄土覆盖的塬梁峁,水土流失严重。黄河花园口以下是被大堤束窄的河道,由河道、河漫滩和大堤构成,覆盖面积较小,仅2.3万km2,且以耕地为主,对于全流域的植被覆盖变化影响较小；因此,本研究以花园口水文站以上作为研究区位于 E 95°52′～113°57′,N 32°8′～41°50′(图 1)。

黄河是我国第2大河,是贯穿北方众多省份的复杂水系,以年均不足600亿m3水资源量承担本流域和下游引黄灌区的农业灌溉、工业生产和城乡居民生活用水[14]。流域大部分地区属于干旱与半干旱地区,在气候和人类活动的影响下,生态环境脆弱,干旱是流域的基本特征。研究区大部分为干旱半干旱和高寒区,多年平均降水量不足450 mm,平均气温为2.68℃,气温的地区差异显著,东南高,西北低,年际变化较大,蒸发力强。这些地貌及气候特点决定了黄河流域的植被种类及覆盖程度的区域差异,河源受低温控制,植被以高寒植物为主,主要植被类型为高寒灌丛、高寒草甸、高寒草原、沼泽及水生植被[15]；黄土高原为温带大陆性季风气候,植被以森林和草原为主,从西部的草原向东逐渐转变为森林草原。

2 数据来源和处理

本研究使用了表征 NDVI的2个数据集:GIMMS和MODIS。GIMMS数据来源于美国航空航天局全球监测与模型研究组,该数据集为半月最大值合成数据[16-17],分辨率是8 km,该数据集经过了校准,几何纠正,消除了火山气溶胶及其它因素的影响,时间段为1982年1月—2006年12月；MODIS数据下载于 EOSDataGateway(NASA-EOS 2006)[18],数据来自于载有“中分辨率成像光谱仪(MODIS)”的Terra和Aqua两颗卫星向地面发送的数据。MODIS-NDVI数据是AVHRR-NDVI的后续产品[19],为地面植被覆盖监测研究提供更长时间序列的资料,这里选用的是1 km月值数据,这个月植被指数产品用无云、经大气订正好的16 d合成的反射率值加权平均得到,时间段为2000年1月—2012年12月。

图1 研究区概图Fig.1 Sketch map of the studied area

由于采用的MODIS产品的分辨率为1 km,与GIMMS的8 km分辨率不一致,为将像元大小统一,将MODIS产品在ArcGIS平台上重采样为8 km。GIMMS数据截止时间为2006年12月,为将数据通过MODIS产品进行延续,将MODIS产品与GIMMS数据进行了拟合,它们具有相同的时间段为2000—2006年,首先将1—12月逐月取6年逐像元的多年平均值,然后基于像元进行拟合,得到每个月的拟合公式,经验证,采用线性拟合是最优的拟合方式,公式为

式中:NDVG为GIMMS的NDVI值,NDVM为MODIS的NDVI值；a为系数,b为常数。

从相关系数平方来看,所有拟合后的R2值均高于0.81,说明MODIS产品与GIMMS数据在空间上超过80%的像元能很好地对应,它们在空间上具有较好的一致性。取得逐月的拟合值后,根据月的拟合值按年和季分别取算术平均,将年和季的NDVI数据延续到2012年。

3 研究方法

3.1 趋势性检验

为分析黄河流域不同时间尺度NDVI的趋势性,采用了非参数检验方法Mann-Kendall趋势性检验法[20-21],Mann-Kendall检验的统计S,利用下式计算:

式中:xj和xk分别为第j时间(年、季、月)和第k时间(年、季、月)的NDVI值；n为系列的记录长度(个数),sgn(xj-xk)为表征函数。

随机序列Si(i=1,2,…,n)近似地服从正态分布,则Si的均值和方差为:

利用下式计算统计检验值

以上为计算Mann-Kendall趋势值的一般过程,为了计算黄河流域每个像元上不同时间尺度的趋势值,基于年、季和月3个时间尺度每一年的栅格数据,在ArcGIS平台上利用Raster Calculator进行不同时间尺度的趋势值计算。

利用ZS的值进行趋势统计的显著性检验。ZS值为正时表明植被覆盖有上升趋势,反之为下降趋势。在双边趋势检验中,如果则接受零假设(无变化趋势),反之则拒绝零假设。Z1-α/2可以从标准正太分布函数获得,为便于比较,这里取90%,95%和99%3个置信水平,它们对应的ZS的绝对值分别为大于等于1.28、1.64和2.32。

3.2 不同时间尺度间的趋势性关系

为观察不同时间尺度趋势的一致性,建立不同时间尺度趋势的相关矩阵,通过不同时间尺度间的相关系数大小来判断两时间尺度间的相似性。为了进一步从像元数量来分析不同时间尺度间的一致性关系,统计了具有相同趋势、相反趋势和无趋势的个数,将不同时间尺度上每个像元能通过90%置信水平的具有增加趋势的像元值记为1,减少趋势的记为-1,不具有趋势的像元值记为2。在ArcGIS中利用Raster Calculator对不同时间进行相加处理,得到2个时间尺度每个像元的相加值(表1)。以此统计两两时间尺度间具有增加趋势像元数、不具有趋势的像元数以及趋势性不一致的像元数占据研究区的比例。

表1 两时间尺度间的趋势性关系表Tab.1 Trend relationships between two different time scales

4 结果与分析

4.1 NDVI的空间分布特征

根据1982—2012年的多年平均NDVI空间分布图来看,总体分布情况是东南高,西北低,西南青藏高原河源区的植被覆盖较差(图2)。东部和南部为山区,NDVI值都较高,达到0.5以上,西部及北部为沙漠和裸岩山区,植被覆盖较少,NDVI值不足0.1。不同季节多年平均NDVI空间分布的总体格局与年平均分布基本相同,随着季节的变化,高值区和低值区的变化都较小,干旱-半干旱过渡区(风沙黄土过渡区)NDVI的变化较大,冬春2冷季的空间分布比较相似,夏秋2暖季比较相近。冬季NDVI值以0.1～0.3为主,大于0.5的高值和小于0.1的低值区所占的区域极小,春季时东部高值区的值有一定的提高,面积也稍有扩大；夏季时低值区的面积明显缩小,只存在毛乌素沙地和库布齐沙漠的西部及阴山西段和黄河河源区北侧的极小区域,处于中间值0.2～0.3的区域明显增加,主要在沙漠的东侧和南部得到了拓展。秋季随着植被的凋零,NDVI值又开始向低值区发展,空间分布特征介于夏季和冬季之间。各季节内不同月份间的NDVI空间分布的相似性都较高(图3),在年内,随着冷季向暖季迁移,气温升高和降雨转多,NDVI的空间转变基本表现为:高值区增加、低值区减少,中值区域被压缩,当气温转向降低和降雨减少时,流域NDVI的空间分布变化相反。

图2 多年平均NDVI值和一年趋势值(Z)的空间分布图Fig.2 Spatial distribution of average annual NDVI value and the trend of average annual NDVI value(Z)

图3 多年平均四季NDVI值空间分布图Fig.3 Spatial distribution of seasonal average NDVI value

4.2 不同时间尺度NDVI的Mann-Kendall趋势性

4.2.1 年的趋势性 黄河流域大部分区域近几十年NDVI变化表现出较强的增加趋势,其中能通过90%显著性水平具有增加趋势的区域面积为68.04%,在95%显著性水平下的像元数也达到了54%,仅有3.5%的区域表现出下降趋势,不具有明显变化的面积略多于1/4(表2)。不具有变化趋势的区域多处于沙漠区的中心地带、裸岩山区和较高的山区,也就是植被覆盖较差的和植被覆盖良好的地区变化相对较小(图2)。表现出增加趋势的主要位于沙漠向黄土转变的过渡区,集中在毛乌苏沙地和库布齐沙漠向黄土过渡的黄河中游多沙粗沙区,这里的NDVI值表现出了明显的增加趋势,其他具有增加趋势的区域处于山地中的谷地和盆地内,具有减少趋势的在研究区内没有出现比较集中的区域。

4.2.2 不同季节的趋势性 4个季节的Mann-Kendall趋势性可以分为2类:具有一定比例增加趋势的秋冬2季和不具有明显变化的春夏2季。从表2中可以看到:冬季具有增加趋势的区域在研究区内的面积比例最大,为59.82%；其次是秋季,占研究区面积的 39.76%；春季所占的面积最小,仅4.65%；夏季具有增加趋势的面积也不超过10%。4个季节中NDVI具有减少趋势性的面积都较小,春季最高,占全流域面积16.63%的像元表现出了减小趋势,冬季和秋季都仅有1%左右的像元有减少趋势。4个季节NDVI随时间不具有明显趋势性变化的像元在全流域中都占有相当比例:夏季所占范围最高,达到83.78%；春季的比例也较高,所占百分数为78.72%；也就是说,夏季和春季黄河流域大部分区域都不具有明显趋势性变化,冬季范围最小,但也占到了38.74%。

表2 一年和季节不同类型趋势(Z)占研究区面积比例Tab.2 Area percentage of annual and different seasonal NDVI trend(Z)in different significant level %

4.2.3 不同月份的趋势性 从不同月份趋势值所占的面积比例(表3)来看,在冷暖转换比较明显的3月和10月具有增加趋势,所占的面积最多,分别为65.12%和70.75%,3月份NDVI表现出增加的区域分布在黄土高原西侧及其以西的区域,也就是毛乌素沙地周边及青藏高原上,10月的NDVI值表现出增加趋势的区域分布在黄土高原及鄂尔多斯高原上(图4)。与它们相邻的2月和9月也有相当比例的区域表现出了增加趋势,而处于温度较高时间段的6月最少,仅10.98%的像元数表现出了增加趋势。所有的月份NDVI值随时间表现出减少趋势的像元数都很少,仅8月具有减少趋势的像元数超过了10%,其他月份均不足10%。不具有变化趋势的面积数在大部分月份所占的比例较高,尤其是1月、4月和6月,都超过了70%,较低的是3月、8月、9月和10月,最少的10月仅27.84%。

4.3 不同时间尺度趋势性的关系

建立月、季和年3个时间尺度的相关系数矩阵来分析不同时间尺度间的关系(表4),相关矩阵显示年的趋势性与季节和月份的趋势性的相关系数都较低,季节与其所包含的月份之间具有一定的相关性,其中,冬季与12月、1月和2月之间都表现为中度相关,说明冬季的植被覆盖增加是3个月的共同体现,而其他3个季节仅与其中的部分月份表现出中度相关,余下月份与相应季节为弱相关,反映出春夏秋3季的植被覆盖变化是由部分月份体现的。进一步分析不同时间尺度对应像元趋势的一致性时,由于每个时间尺度上具有减少趋势的像元所占的比例都较少；因此,统计时将之忽略,为进一步概化统计值,分别统计了2个时间尺度上在90%显著性水平下具有增加趋势的像元占据整个研究区的比例,都不具有趋势性的像元比例以及趋势性不一致的比例。从植被覆盖增加趋势的像元占据研究区的比例来看,年与其他时间尺度共同拥有大约35%具有增加趋势的像元数,其中10月最高,接近50%,最低的是夏季,仅16.17%。冬季和秋季与它们包含的月份之间共同拥有增加趋势的像元比例较高,大部分都超过38%,而春季和夏季与它们包含的月份之间共同拥有增加趋势的像元数较低,都不足20%(表4)。

表3 不同月份不同类型趋势占研究区面积比例Tab.3 Area percentage of different monthly NDVI trend in different significant level %

图4 3月和10月NDVI的 Mann-Kendall值(Z)空间分布图Fig.4 NDVI spatial distribution of Mann-Kendall value of NDVI(Z)in March and October

表4 不同时间尺度间趋势性的相关系数及具有共同增加趋势的像元占据研究区的比例Tab.4 Matrix of correlation coefficients between different temporal scales and percentage of pixels with common trend between different temporal scales__

从统计出的表5看,年与季节和月份间趋势性不一致的像元数在3类统计中所占比例最高,平均值达到47.77%,接近于研究区的半数,春季和夏季与它们包含的月份间趋势不一致的比例相对较高,而冬季和秋季与它们所包含月份之间趋势不一致的比例较少,在20% ～35%之间。从统计表5看,共同不具有趋势性的像元数所占的比例总体上不高,均值仅20.64%,最高值为秋季与1月之间的比例达到59.45%,其他时间之间的比例都远低于这个值,多居于10% ～30%之间。

表5 不同时间尺度间共同拥有趋势性不一致的像元比例及不同时间尺度间都不具有趋势性的像元比例Tab.5 Percentages of pixels in inconsistent trend among different temporal scales(Left)and percentages of pixels in non-trend among different temporal scales(Right)

5 结论

1)黄河流域3个时间尺度(年、季、月)NDVI的空间分布特征为:东部及南部山区NDVI值较高,西北部沙漠及土石山区的NDVI值较低,黄河源区的NDVI值较低。冬春2季和夏秋2季的NDVI空间分布分别具有比较高的相似性。

2)NDVI趋势性在年尺度上反应最好,超过60%的像元表现出了增加趋势,季节上,仅冬春部分区域表现出良好的增加趋势,而夏秋具有增加趋势的像元数都较少。冷暖转换月份3月和10月NDVI增加的趋势性明显高于其他月份。

3)相关性分析显示年与其他任一时间尺度间的NDVI趋势性的相关性不强,不同季节与其所包含月份间的NDVI趋势性在相关性上存在差异,冬季与包含月份间的都具有中度相关,而春夏秋3季与它们所包含的部分月份间存在中度相关。不同趋势的像元数统计显示,年和冬季与其他时间尺度间共同具有增加趋势的像元数较多,而夏秋2季与其他时间尺度间占有的比例较少。

综上所述,黄土高原植被覆盖在不同时间尺度上的变化存在明显差异,在利用年的趋势值来说明植被变化时应该谨慎对待,需要更为清晰地说明更小时间尺度的植被变化是否能够体现更高时间尺度的植被变化。

6 参考文献

[1] CRAMER W P,LEEMANS R.Assessing impacts of climate change on vegetation using climate classification systems[M]∥SOLOMON A M,SHUGART H H.Vegetation Dynamics and Global Change.London:Chapman and Hall,1993:190.

[2] MYNENI R B,KEELING C D,TUCKER C J,et al.Increased plant growth in the northern high latitudes from 1981 to 1999[J].Nature,1997,386(6626):698.

[3] MYNENI R B,TUCKER C J,ASAR G,et al.Interannual variations in satellite-sensed vegetation index data from 1981 to 1991[J].Journal of Geophysical Research Atmospheres,1998,103(D16):6145.

[4] ZHOU L,TUCKER C J,KAUFMANN R K,et al.Variations in northern vegetation activity inferred from satellite data of vegetation index during 1981- 1999[J].Journal of Geophysical Research Atmospheres,2002,106(17):20069.

[5] NEMANI R R,KEELING C D,HASHIMOTO H,et al.Climate-driven increases in global terrestrial net primary production from 1982 to 1999[J].Science,2003,300(5625):1560.

[6] 方精云,朴世龙,贺金生,等.近20年来中国植被活动在增强[J].中国科学(C辑:生命科学),2003,33(6):554.FANG Jingyun,PIAO Shilong,HE Jinsheng,et al.Activity of vegetation in China has been augmented in recent 20 years[J].Science in China Series C:Biological Sciences,2003,33(6):161.

[7] 陈效逑,喻蓉.1982—1999年我国东部暖温带植被生长季节的时空变化[J].地理学报,2007,62(1):41.CHEN Xiaoqiu,YU Rong.Spatial and temporal variations of the vegetation growing season in warm-temperate eastern china during 1982 to 1999[J].Acta Geographica Sinica,62(1):41.

[8] ZHANG Yili,GAO Jungang,LIU Linshan,et al.NDVI-based vegetation changes and their responses to climate change from 1982 to 2011:A case study in the Koshi River Basin in the middle Himalayas[J].Global and Planetary Change,2013,108:139.

[9] 刘志红,郭伟玲,杨勤科,等.近20年黄土高原不同地貌类型区植被覆盖变化及原因分析[J].中国水土保持科学,2011,9(1):16.LIU Zhihong,GUO Weiling,YANG Qinke,et al.Vegetation cover changes and their relationship with rainfall in different physiognomy type areas of Loess Plateau[J].Science of Soil and Water Conservation,2011,9(1):16.

[10] 钟莉娜,赵文武.基于NDVI的黄土高原植被覆盖变化特征分析[J].中国水土保持科学,2013,11(5):57.ZHONG Lina,ZHAO Wenwu.Detecting the dynamic changes of vegetation coverage in the Loess Plateau of China using NDVI data[J].Science of Soil and Water Conservation,2013,11(5):57.

[11] 何远梅,姚文俊,张岩,等.黄土高原区植被恢复的空间差异性分析[J].中国水土保持科学,2015,13(2):63.HE Yuanmei,YAO Wenjun,ZHANG Yan,et al.Spatial variability of vegetation restoration on the Loess Plateau based on MODIS/NDVI[J].Science of Soil and Water Conservation,2015,13(2):63.

[12] 信忠保,许炯心,郑伟.气候变化和人类活动对黄土高原植被覆盖变化的影响[J].中国科学(D辑:地球科学),2007,37(11):1504.XIN Zhongbao,XU Jiongxin,ZHENG Wei.Impact of climate change and human activities in vegetation cover change on the Loess Plateau[J].Science in China Series D:Earth Sciences,2007,37(11):1504.

[13] CUO Lan,ZHANG Yongxin,GAO Yanhong,et al.The impacts of climate change and land cover/use transition on the hydrology in the upper Yellow River Basin,China[J].Journal of Hydrology,2013,502:37.

[14] CAI Ximing,ROSEGRANT M W.Optional water development strategies for the Yellow River Basin:Balancing agricultural and ecological water demands[J].Water Resources Research,2004,40(8):1.

[15] 武吉华,张绅.植物地理学[M].3版.北京:高等教育出版社,1995:1.WU Jihua,ZHANG Shen.Geobotany[M].3rd Edition.Beijing:Higher Education Press,1995:1.

[16] PINZON J.Using HHT to successfully uncouple seasonal and interannual components in remotely sensed data[C]∥SCI 2002 Conference Proceedings,2002:14.

[17] TUCKER C J,PINZON J E,BROWN M E,et al.An extended AVHRR 8-km NDVI dataset compatible with MODIS and SPOT vegetation NDVI data[J].International Journal of Remote Sensing,2005,26(20):4485.

[18] NESTLER M S,PFISTER R.The NASA earth observing system data gateway[C]∥AGU Fall Meeting Abstracts,2002.

[19] HUETE A,DIDAN K,MIURA T,et al.Overview of the radiometric and biophysical performance of the MODIS vegetation indices[J].Remote Sensing of Environment,83(1 and 2),2002:195.

[20] MANN H B.Nonparametric tests against trend[J].Econometrica,1945,13:245.

[21] KENDALL M G.Rank correlation methods[M].2nd edition: Oxford, England: HafnerPublishing Co.,1955:85.