余汛 王莉 景海涛 李新元 范晨雨

摘要：归一化的植被指数（NDVI）的时空变化可以反映生态环境的演变，以MODND1M中国500 m NDVI数据和气象数据为数据源，采用趋势分析和相关分析方法，研究了太行山地区的NDVI、气温、降水量的时空变化特征及其相关性。研究结果表明：（1）16年来，太行山地区的平均NDVI呈增加趋势，增长率为0.079/10年，其中NDVI增加区域占太行山地区的96.84%;（2）在空间尺度上，NDVI与降水量呈正相关性的像元数占像元总数的30.7%，NDVI与气温呈正相关的像元数占总像元数的26.8%;（3）在时间尺度上，NDVI与降水量的相关系数为0.43，NDVI与气温的相关系数为0.22;（4）研究区域降水量与气温对植被生长存在滞后效应。

关键词：MODIS NDVI;响应;ArcGIS;太行山地区;气候因子;克里金插值;相关性;时空变化

中图分类号：Q948.112;S127

文献标志码：A

文章编号：1002-1302（2020）16-0260-08

植被是地表物质循环与能量流动的主要载体[1]，与地质、地貌、水文、土壤、动物和微生物共同构成自然地理环境，且最能反映其他要素性质[2]。气候对植被的空间分布、生产力、物候[3-5]等具有重要影响。因此，深入研究植被的动态变化及与气候间的相关性，对揭示区域生态环境演变以及应对气候变化具有重要的现实意义。

归一化的植被指数（NDVI）能够比较准确地反映植被代谢强度、植被绿化程度、光合作用强度、季节变化和年际变化特征，因此该指数可用于植被动态监测、植被分类、作物长势检测、物候监测等研究[6-7]。不少研究显示，近年来全球植被呈变好的趋势，尤其在北半球中高纬度地区变好趋势显著[8-12]。朴世龙等利用气象卫星NOAA-AVHRR（美国国家海洋和大气管理局-甚高分辨率扫描辐射计）数据探讨了植被活动对全球变化的响应方式，发现我国20世纪80、90年代NDVI有明显的增加趋势[10]。刘爽等利用中分辨率成像光谱仪（moderate-resolution imaging spectroradiometer，MODIS）数据发现，2000—2010年中国植被绿度正在全面提高，荒漠化在减少[13]。

山地作为陆地生态系统中最为活跃的地理单元，是环境变化的驱动机和放大器，山地植被的变化对全球变化的响应更为敏感[14]。太行山脉作为京津和华北大平原的天然屏障，不仅能抵御西北寒潮的袭击，还能接受东南暖湿气流。近年来，由于太行山森林资源曾遭到严重破坏，水土流失严重，容易爆发山洪。相关部门为改善太行山的植被覆盖率实施了太行山绿化等一系列的植被恢复计划。目前，太行山地区植被变化的相关研究已有一部分成果，如李薇等关于太行山不同坡度NDVI变化趋势差异分析研究等[15]。已有的研究结果没有深入地分析太行山地区的植被覆盖率变化信息及其对气候因子的响应。目前，用于植被动态变化检测的遥感数据源日益丰富，常用的NDVI数据源主要有NOAA-AVHRR、陆地卫星（Landsat TM/ETM +）、SPOT、MODIS等，在植被动态检测中各有千秋，其中MODIS数据具有较好的空间分辨率和时间分辨率，研究者可以从中获取更详尽的关于研究区域土地覆盖和生态系统变化的信息[16]，因此，本试验基于MODIS NDVI数据，结合气温、降水等数据，对太行山地区2000—2015年这16年来的植被覆盖率的时空变化情况进行研究，并研究植被对气候因子的响应机制。

1 研究区与数据源

1.1 研究区概况

太行山脉位于山西省与华北平原之间，纵跨北京、河北、山西、河南4省（市），山脉北起北京市西山，向南延伸至河南与山西交界地区的王屋山，西接山西高原，东临华北平原，呈东北—西南走向，绵延约400 km。它是我国地形第二阶梯的东缘，也是黄土高原的东部界线，位于34°71′～40°34′N、110°60′～115°62′E，总面积12.7万km2，平均海拔为1 500～2 000 m，最高峰为小五台山，高3 099 m。

1.2 数据源与数据处理

本研究采用的NDVI数据来源于中国科学院计算机网络信息中心国际科学数据镜像网站（http：//www.gscloud.cn）气象卫星MODND1M中国500 m NDVI月合成产品，数据是由MODND1D计算得到，计算方法为取月内每天最大值。坐标系EPSG：（WGS84），空间分辨率500 m，时间分辨率为每个月，时间为2000—2015年。

本研究所使用的MODIS NDVI数据已经经过水、云、气溶胶等处理，减弱了其他因素对数据的影响，再根据太行山地区的矢量边界对NDVI图像进行裁剪，得到研究区的影像数据。

研究区的边界矢量图来自中国科学院山地表生过程与生态调控重点实验室[17]。

气象数据为中国气象数据共享网的（http：//www.nmic.cn/）太行山区内18个站点的2000—2015年日降水量及日平均气温数据，并对气象数据进行异常值剔除等预处理。

1.3 趋势分析法

基于像元尺度的趋势分析法能模拟研究区中每个栅格单元的变化趋势，从而反映植被NDVI变化的方向和速率[19]。计算公式为：

式中：Slop是NDVI像元線性回归方程的斜率;i代表年数，取值范围为1～16;n是时间跨度，在本研究中n的数值为16年。当Slop>0时，NDVI呈增加趋势;当Slop=0时，NDVI基本稳定，无明显变化;当Slop<0时，NDVI呈减少趋势。利用F检验显著性分析，公式为：

式中：U=∑ni=1（yi-y）2，为误差平方和;Q=∑ni=1（yi-y）2为回归平方和;yi为第i年的植被覆盖率观测值;y为覆盖率的16年来的平均值。将得到的结果划为3个等级：极显著相关（P<0.01）、显著相关（0.010.05）。

1.4 克里金插值方法

克里金插值是基于内蕴假设，应用变异函数来研究在空间上随机分布并且相关的变量的方法。其方法是根据待估测点周围的若干已知信息，以变异函数作为工具，确定待估测点周围已知点的参数对待估测点的加权值，然后对待估测点做出最优的、无偏差的估测。克里金插值法不仅能够预测表面，还能对预测的确定性提供评判标准[18]。其公式为：

式中：Z（x0）代表未知样点的值;Z（xi）代表未知样点周围的已知样点的值;N表示已知样点的数量;λi是第i个样本点的权重。在地统计中，一共有6种克里金插值法，分别是普通克里金、简单克里金、泛克里金、指示克里金、概率克里金、析取克里金法，其中普通克里金法是所有克里金插值方法中最基本、最重要、应用最广的插值方法[18]。

1.5 相关性分析

相关分析是测量2个变量之间相关程度的有效方法。本研究先对气温、降水量数据进行插值分析，得到了与NDVI数据分辨率相同的栅格数据集，再利用ArcGIS进行空间分析，基于像元尺度的NDVI与降水量、气温间的相关分析的结果能说明NDVI与降水量、气温间的相关程度以及分析结果的空间分布情况，相关系数的计算公式如下：

式中：rxy为相关系数;xi与yi分别代表x和y在第i时期的值;x和y代表x和y的平均值;n为样本容量。相关系数值为-1～1，正值代表变量间呈正相关关系，负值代表变量间呈负相关关系，相关系数的绝对值越接近于1，变量间相关性越强。

2 结果与分析

2.1 太行山地区植被NDVI时空变化分析

基于研究区2000—2015年各年年均NDVI数据和趋势斜率法，得到太行山区年均NDVI变化程度分布图，结果如图2所示，把研究区的NDVI变化情况划分为中度退化、轻度退化、基本不变、轻度改善、中度改善、明显改善6个等级。从表1可以看出，太行山地区植被退化地区的面积占总面积的0.163%，改善的面积占总面积的96.84%。

基于ArcGIS提取的各年年均NDVI值做一元线性回归分析，由图3可知，16年来太行山地区植被的NDVI值波动比较剧烈，NDVI值都在大于0.40，最高值（0.57）出现在2014年，最低值（0.41）出现在2005年，平均值为0.499。对研究区16年间的年均植被NDVI值进行线性趋势分析，得到该区域年均植被NDVI值的变化趋势方程，可以看出2000—2015年研究区NDVI呈微弱增长趋势，增长率为0.079/10年（P=0.000 078 6<0.01），通过显著性检验，这表明16年太行山地区的植被覆盖呈增长趋势。

基于ArcGIS提取各季季均NDVI值做一元线性回归分析，结果如图4所示，添加趋势线得到回归方程，并在SPSS中对16年来各季季均植被NDVI的变化趋势做显著性检验，从表2中可以看出，太行山地区植被覆盖变化特征具有明显的季节差异，夏季的植被生长状态最好，秋季的植被生长状态次之，冬季的植被生长状态最差，这符合植被生长的自然规律。从图4可以看出，16年来的太行山区域的植被都表现为不同程度的增长趋势，其中冬季（0.095/10年）>夏季（0.086/10年）>春季（0.068/10年）>秋季（0.014/10年），每个季节的增加速率都通过了显著性检验，表明太行山地区每个季节的植被生长率都有改善的趋势。

2.2 研究区气温和降水的变化特征

在ArcGIS图层属性中统计太行山地区2000—2015年各年的平均气温数据，绘制年际温度变化趋势，结果如图5所示，添加趋势线得到年均气温的变化趋势方程，然后在SPSS中进行显著性检验，研究区的年均气温的变化趋势为0.325 ℃/10年（P=0.304），没有通过0.05水平显著性检验，说明16年来太行山地区气温升高趋势不显著，最高气温出现在2014年，为11.34 ℃，最低气温出现在2012年，为9.34 ℃。

在ArcGIS图层属性中统计2000—2015年各年的季均气温数据，绘制出太行山16年来各季季均气温的年际变化趋势，结果如图6所示，将四季气温的变化趋势方程添加进去，并进行显著性检验，其中仅夏季气温（P=0.041<0.05）通过显著性检验，说明夏季气温上升趋势明显，而春季、秋季、冬季的氣温上升趋势不明显。春季最高气温出现在2010年，为17.02 ℃，最低气温出现在2011年，为10.66 ℃;夏季最高气温出现在2010年，为 22.96 ℃，最低值出现在2004年，为20.94 ℃;秋季的最高气温出现在2006年，为11.26 ℃，最低气温出现在2000年，为9.11 ℃;冬季的最高气温出现在2007年，为 -2.15 ℃，最低气温出现在2005年，为-6.09 ℃。整的来看，16年来太行山地区四季气温均有变暖趋势。

在ArcGIS图层属性中统计太行山地区2000—2015年各年年均降水量数据，绘制年际变化趋势，添加线性趋势线得到降水量的变化趋势方程，并进行显著性检验（图7）。太行山地区年均降水量的变化趋势为16.16 mm/10年（P>0.05），说明16年来太行山地区降水量无明显增加趋势。

在ArcGIS图层属性中统计太行山地区2000—2015年各年的季均降水量，绘制太行山区的16年间的各季季均降水量的年际变化趋势，添加各季的线性趋势线得到各级的降水量变化趋势方程，并在SPSS中进行显著性检验（图8）。太行山地区16年来的4个季节的降水量呈现不同程度的上升或下降趋势，春季、秋季的降水量呈上升趋势，变化趋势为春季（11.543 mm/10年）>秋季（11.185 mm/10年），夏季和冬季的降水量呈微弱下降趋势，变化趋势分别为冬季（-3.8 mm/10年）>夏季（-0.319/10年），其中夏季（P=-0.003<0.05）降水变化趋势通过显著性检验，说明夏季降水量减少趋势明显。

2.3 研究区NDVI变化对气温和降水量变化的响应

因为MODIS-NDVI数据是从2000年开始记录的，所以把年平均气温、年累积降水量数据集中的2000—2015年数据提取出来。在ArcGIS中分别进行NDVI与降水量、气温的空间相关性分析，并在SPSS软件分别做NDVI与气温、降水量的时间相关性分析。基于太行山地区2000—2015年年均NDVI、气温、降水量的栅格数据，在ArcGIS中借助栅格計算器工具，分别得到NDVI与气温、降水量的相关系数空间分布图，再利用t检验对相关系数进行显著性检验，分别得到NDVI与气温、降水量的t检验空间分布图。

从图9可以看出，16年来太行山区绝大部分区域的植被NDVI与气温的相关性较弱，只有1.09%的区域通过了0.01水平显著性检验，4.27%的区域通过了0.05水平显著性检验。呈显著负相关及极显著负相关的区域分布在襄垣、长治、太古境内，说明这些区域的植被与气温呈负相关关系，呈显著、极显著正相关区域分布在怀来境内，说明该地区的植被与气温呈正相关。整体来看，相关系数为正的区域所占比例为26.8%，主要分布在怀来、蔚县、盐湖境内;相关系数为负的区域所占比例为73.2%，主要分布在侯马东北方向至灵丘境内、以及延庆境内。

从图10可以看出，只有3.84%的区域通过了0.05水平显著性检验，0.94%的区域通过了0.01水平的显著性检验。呈极显著、显著正相关的区域集中分布在太谷、安阳境内，说明这些区域植被与降水量呈较好的正相关关系;呈显著、极显著负相关的区域主要分布在怀来境内。相关系数自西南向东北逐渐递减，相关系数为正的区域主要分布在太行山地区的西南部，所占比例为30.7%，相关系数为负的区域主要分布在太行山区的东北部，所占比例为69.3%。

将太行山地区2000—2015年年均植被NDVI分别和同时期的年均气温、年降水量在Excel中做双轴折线图分析，从图11可以看出，太行山地区植被NDVI的变化与年均气温的变化并不完全一致，2000—2001、2004—2005年这几年间的植被NDVI随气温的升高而表现出下降的趋势;2007—2008、2009—2011年期间植被NDVI随着气温的降低表现出上升的趋势，这表明太行山地区的年均气温过高或过低都会抑制植被的生长，整体来看，16年来太行山地区的年均NDVI与年均气温呈正相关关系，相关系数为0.217。

从图11可以看出，太行山地区年均NDVI值和年均降水量在大多年份的变化趋势表现一致。2006—2007、2009—2010年期间，随着降水量的增加，植被的NDVI反而减少;其他年份植被NDVI的变化与降水量变化呈正相关性，即植被NDVI随着降水量的增加而增加。整体看来，16年来太行山地区年均NDVI与降水量呈正相关性，相关系数为0.43。

将太行山地区的年均植被NDVI与年均降水量和年均气温变化趋势进行对比，得出太行山地区植被NDVI的变化是降水量和气温变化共同作用的结果。

为了揭示太行山地区各个季节植被的NDVI值与气候变化的关系，在SPSS中将2000—2015年各年的季平均植被NDVI分别与同时段各年季平均气温、季均降水量做相关性分析。因为气温、降水量变化对植被的生长会存在一定的滞后效应，所以还要对各季季均植被NDVI和其对应的前1季的季均降水量、季均气温做相关性分析（表3）。其中，春季、夏季的季均植被NDVI与气温呈负相关关系，且夏季负相关性较高，相关系数为-0.201，说明夏季气温变高，会影响植被的生长。秋季、冬季的季均植被NDVI与气温呈正相关关系，表明秋季、冬季气温变高，会促进植被的生长。4个季节的季均植被NDVI与前一季的气温均呈正相关关系，但是都没通过显著性检验，说明气温与植被NDVI无显著相关性。由结果可得，春夏秋冬四季的气温高，会促进下一季的植被生长。

将太行山地区2000—2015年各年季均植被NDVI与同季、前1季的季均降水量在SPSS中做相关性分析并进行显著性检验，发现春、夏、秋3个季节的季均植被NDVI与降水量呈正相关关系。其中，春季降水与植被NDVI值的相关性最大，为0.443，表明春季降水量会促进植被的生长;夏季与秋季的季均植被NDVI与季均降水量相关性虽然相对较小，但也能说明夏季与秋季的季均降水量能促进植被生长;冬季的季均植被NDVI与降水量呈现负相关性，表明冬季较多的降水会抑制植被的生长。春季季均植被NDVI与冬季降水量呈负相关性，表明冬季降水量会抑制春季植被的生长，夏季的季均植被NDVI与春季的降水量呈正相关性，且通过了0.05水平显著性检验，表明降水量增加会促进夏季的植被生长。秋季和冬季的季均植被NDVI与前1季的降水量呈正相关关系，其中秋季的季均植被NDVI与夏季降水量呈弱正相关性，表明夏季和秋季的较多的降水量会促进后一季的植被生长。综上所述，太行山地区的降水量和气温对植被生长存在滞后效应。

3 结论与讨论

本试验基于太行山地区2000—2015年的MODIS NDVI数据，结合气温、降水量数据，对该地区2000—2015年16年间的植被覆盖率时空变化进行研究，并对影响其变化的气候与植被覆盖情况进行相关性分析，得到如下结论：

（1）2000—2015年研究区NDVI呈微弱增长趋势，增长率为0.079/10年，植被增加的区域占总区域的96.84%，植被减少及基本不变的区域为3.16%，表明16年太行山地区的植被覆盖呈增长趋势。

（2）总的来看，16年来太行山地区的年均NDVI与年均气温呈正相关关系，相关系数为0.217。太行山地区年均NDVI与降水量呈正相关性，相关系数为0.43。将太行山地区的年均植被NDVI与年均降水量、年均气温变化趋势进行对比，可知太行山地区植被NDVI是降水量和气温变化共同作用的结果。

（3）通过相关性分析发现，夏季、春季的季均植被NDVI值与气温呈负相关关系，表明夏季、春季气温较高，会影响植被的生长。秋季和冬季的季均植被NDVI与气温呈正相关关系，表明秋季、冬季气温较高会促进植被的生长。春季、夏季与秋季的季均植被NDVI与季均降水量呈正相关关系，表明春季、夏季与秋季的季均降水能促进植被生长;冬季季均植被NDVI与季均降水量呈负相关关系，表明冬季的季均降水量抑制植被的生长。

（4）4个季节的季均植被NDVI与前1季的气温均呈正相关关系，表明春夏秋冬四季的气温较高，会促进下一季的植被生长。冬季的季均降水量与春季的季均NDVI呈负相关关系，表明冬季过高的降水量会抑制春季的植被生长，夏季、秋季、冬季的季均植被NDVI与前1季的降水量呈正相关关系，其中秋季的季均植被NDVI与夏季降水量呈弱正相关关系。表明春季、夏季、秋季较多的降水会促进后一季的植被生长。

（5）太行山地区的降水量和气温对植被生长存在滞后效应。植被覆盖变化是人为因素和自然因素共同作用的结果。本研究主要结合气温、降水量探讨太行山地区植被变化趋势以及驱动因素。本研究对影响植被变化的自然因素，只从气温和降水量2个方面进行研究，对湿度、下行短波辐、坡度等其他的自然因素没有研究，而且没有考虑人为因素对植被覆盖变化的影响。因此，在后续植被变化的研究中，对其驱动因素的选取上应该不断完善。

参考文献：

[1]赵 杰，杜自强，张 红，等. 中国季节性昼夜增温的不对称性及其对植被活动的影响[J]. 生态学报，2018，38（11）：3909-3919.

[2]杨 菊，熊康宁，肖时珍. 赤水丹霞地貌区植被覆盖度时空变化分析[J]. 人民长江，2018，49（18）：38-44.

[3]於 琍，李克让，陶 波，等. 植被地理分布对气候变化的适应性研究[J]. 地理科学进展，2010，29（11）：1326-1332.

[4]张景华，李英年. 青海气候变化趋势及对植被生产力影响的研究[J]. 干旱区资源与环境，2008，22（2）：97-102.

[5]候 静，杜灵通，刘 可，等. 1982—2012年北半球荒漠草原过渡带植被物候特征及其与气候因子的关系[J]. 气候变化研究进展，2017，13（5）：473-482.

[6]Xiao X M，Boles S，Liu J Y，et al. Characterization of forest types in Northeastern China，using multi-temporal SPOT-4 VEGETATION sensor data[J]. Remote Sensing of Environment，2002，82（2）：335-348.

[7]赵赛帅. 基于环境卫星NDVI时间序列的盐沼植被分类与提取[D]. 南京：南京大学，2015.

[8]Myneni R B，Tucker C J，Asrar G，et al. Interannual variations in satellite-sensed vegetation index data from 1981 to 1991[J]. Journal of Geophysical Research Atmospheres，1998，103（D6）：6145-6160.

[9]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-1563.

[10]朴世龙，方精云. 1982—1999年我国陆地植被活动对气候变化响应的季节差异[J]. 地理学报，2003，58（1）：119-125.

[11]朴世龙，方精云. 最近18年来中国植被覆盖的动态变化[J]. 第四纪研究，2001，21（4）：294-302.

[12]方精云，朴世龙，贺金生，等. 近20年来中国植被活动在增强[J]. 中国科学（C辑），2003，33（6）：554-565，578-579.

[13]刘 爽，宫 鹏. 2000—2010年中国地表植被绿度变化[J]. 科学通报，2012，57（16）：1423-1434.

[14]常 纯，王心源，杨瑞霞，等. 基于DEM-NDVI的高山植被带定量刻划[J]. 地理研究，2015，34（11）：2113-2123.

[15]李 薇，谈明洪. 太行山区不同坡度NDVI变化趋势差异分析[J]. 中国生态农业学报，2017，25（4）：509-519.

[16]李晓光. 基于MODIS-NDVI的内蒙古植被覆盖变化及其驱动因子分析[D]. 呼和浩特：内蒙古大学，2014.

[17]刘宪锋，任志远，林志慧，等. 2000—2011年三江源区植被覆盖时空变化特征[J]. 地理學报，2013，68（7）：897-908.

[18]牟乃夏，刘文宝，王海银，等. ArcGIS 10地理信息系统教程——从初学到精通[M]. 北京：测绘出版社，2012（12）：43-43.

[19]王丽霞，余东洋. 渭河流域NDVI与气候因子时空变化及相关性研究[J]. 水土保持研究，2019，26（2）：249-254.