不同地形梯度上的植被变化趋势及原因分析

马士彬 ，安裕伦 ，杨广斌 ，张勇荣

1. 贵州师范大学地理与环境科学学院，贵州 贵阳 550001；2. 贵州省山地资源与环境遥感应用重点实验室，贵州 贵阳 550001；3. 贵州师范大学喀斯特研究院，贵州 贵阳 550001；4. 喀斯特山地生态环境国家重点实验室培育基地，贵州 贵阳 550001；5. 六盘水师范学院旅游与历史文化学院，贵州 六盘水 553004

贵州省是亚热带喀斯特的连片分布区（Sweeting，1995），主要的生态环境问题之一是石漠化（Wang et al.，2002）。监测石漠化发生与演变需要客观分析植被的动态特征（熊康宁等，2002）。植被指数可以准确地反映植被的动态变化特征，如归一化植被指数（Normal differential vegetation index，NDVI）（Tucker et al.，1986），因此可以通过卫星遥感图像计算植被指数对陆地生态系统进行长时间、大尺度的观测。

卫星遥感监测显示，2000年以来贵州省NDVI值整体呈上升趋势，但存在空间差异，显著增加的区域主要分布在六盘水、毕节西部和安顺北部（Tian et al.，2017）。NDVI呈良性变化的原因目前存在争议，从植被对气候变化响应的角度分析认为，贵州省植被动态变化与气温呈正相关关系，并且气温的影响大于降水（Yang et al.，2004；Hou et al.，2015）；但从人类活动对植被变化影响的角度认为，贵州省植被动态变化与气候变化无显著相关性，主要是人为因素导致的（Cai et al.，2014；Tong et al.，2016）。显著减少的区域主要分布在贵阳、六盘水、遵义、凯里等城市周边和贵州省北部、东北和东南部等区域（Tian et al.，2017；Cai et al.，2014）。引起城市周边NDVI降低的原因主要是城市扩张，而贵州省北部、东北和东南部等非喀斯特区域的NDVI降低的原因目前还缺乏相关研究。

正确地分析植被变化的趋势和驱动力是掌握区域生态系统变化规律的基础（Pan et al.，2018）。目前大部分研究都是在线性趋势分析基础上，得到二维空间植被变化特征。但是，由于地形严格控制着人类活动的空间范围（马士彬等，2008；马士彬等，2012；郜红娟等，2015；周忠发等，2019）和自然地理要素的空间分布，如果在地形破碎的区域进行植被变化分析时仅考虑二维的空间特征，那么垂直结构上的特征有被掩盖的危险。Tong et al.（2016）以广西和贵州为研究区，在分析2001-2013年植被变化趋势的基础上评价了未来植被变化规律的可能性；研究得出，植被退化主要发生在低地形位上，中地形位上植被具有很好的恢复前景，高地形位上植被变化趋势的不确定性风险较大，应是今后生态工程实施重点关注的区域。上述研究虽然分析了植被变化趋势的垂直特征，但没有系统分析其形成原因。因此，本研究采用Mann-kendall趋势分析法和地形位指数，分析贵州省植被变化趋势的空间分布特征，并探讨不同地形梯度上植被变化的直接驱动力。

1 数据与方法

贵州省位于中国西南部，与云南、四川、湖南、广西、重庆接壤，国土面积1.76×105km2，其中喀斯特地貌占72%，海拔在1000-3000 m之间；根据贵州省农业地貌区划，将研究区划分为东部丘陵山地（Hill-mountain Eastern，HME）、中部丘原山地（Hilly plain-mountain Center，HPMC）和西部高原山地（Plateau-mountain Western，PMW）3个地貌区（图 1a）。贵州属于亚热带湿润季风气候，年平均气温14-16 ℃，年平均降水量1000-1300 mm。

1.1 数据源

1.1.1 MODIS-NDVI时序数据

为了分析 2000-2015年植被的变化趋势，从LPDAAC网站（https://e4ftl01.cr.usgs.gov/）收集了MODIS 16 d的NDVI数据，空间分辨率250 m。采用最大值方法将16 d NDVI合成月NDVI，利用月NDVI计算年平均NDVI和季节NDVI。

1.1.2 2000-2015年土地利用/覆被数据

土地利用/覆被数据来源于全国生态环境10年变化（2000-2010年）遥感调查评估（欧阳志云等，2014）和中国科学院国家资源环境科学数据中心（RESDC）（http://www.resdc.cn/）。土地利用/覆被分为9大类：亚热带常绿阔叶林、亚热带落叶阔叶林、亚热带针叶林、亚热带混合林、灌木、草地、耕地、水体和建设用地，图1b为研究区土地利用/覆被图。

1.1.3 气象数据

2000-2015年日气温和降水数据来源于中国气象局（http://data.cma.cn/site/index.html）。收集了51个气象站，包括贵州33个，四川2个，重庆5个，湖南3个，广西4个，云南4个，气象站点分布见图1b。由日值数据合成月值，然后分别计算年值和季节值，其中12月、1月、2月为冬季，3月、4月、5月为春季，6月、7月、8月为夏季，9月、10月、11月为秋季。

1.2 方法

1.2.1 预处理

图1 研究区数字高程模型（DEM）和地貌分区图（a），土地利用/覆被类型和气象站点分布图（b）Fig. 1 Distribution of landform regions overlaying a digital elevation model (DEM)(a); the spatial distribution of landuse/cover types and climate stations are shown(b)

MODIS NDVI时间序列中的异常值是由云污染和设备故障等因素造成的，为了减少这种噪声影响时间序列的实际趋势，预处理是必要的（Cai et al.，2015）。运用 TimeSat 3.1软件，利用Savitzky-Golay滤波法去除低质量观察点。图2展示了原始数据 2000-2015年月时间序列和重构后的月时间序列。

1.2.2 Mann-kendall变化检测及突变分析

应用 Mann-kendall趋势变化检测模型，从MOD13Q1 NDVI每个栅格的时间序列中计算出Mann-Kendall单调趋势显著性的 Z值统计量（Mann，1945；Kendall，1970）。

由正态分布表可知，当|z|≥2.23和|z|≥1.64表示时间序列在 0.01和 0.05水平上显著变化。根据 z值将植被变化趋势分为3个类型：

（1）Z≥2.23，极显著增加；

（2）1.64≤Z＜2.23，显著增加；

（3）-1.64＜Z＜1.64，无显著变化；

（4）-2.32＜Z≤-1.64，显著降低；

（5）Z≤-2.23，极显著降低。

1.2.3 Landtrendr扰动类型分析

LandTrendr（Landsat-based Detection of Trends in Disturbance and Recovery）算法可以同步检测出变化趋势和扰动事件（Kennedy et al.，2010）。基本原理是采用任意时间分割技术分割光谱轨迹，用直线线段来模拟时间轨迹的重要特征，重要特征分割点的时间和光谱值为生成变化图提供所需的基本信息。本文利用 LandTrendr算法检测植被变化的扰动类型，分为扰动后修复、扰动后不修复和未扰动。

1.2.4 地形位指数

地形是地表空间分异的重要影响因子之一，因此地表景观结构及其变化往往在地形梯度上表现出一定的分布规律。为了定量分析植被变化趋势的地形分异规律，本研究使用地形位指数来度量地形梯度。地形位指数可用下式进行计算：

1.2.5 分布指数

为了消除地形位梯度分段和植被斑块面积差异的影响，使用分布指数来描述植被变化趋势在地形位梯度上的分布特征。分布指数为：

图2 MODIS NDVI数据预处理效果图Fig. 2 NDVI derived from MOD13Q1

式中，Pie为第 i种植被变化趋势在地形位 e上的分布指数；Sie为第 i种植被变化趋势在地形位e上的面积；Si为研究区内第i种植被变化趋势的面积；Se为研究区内地形位e的总面积；S为整个研究区的总面积。分布指数是一个标准化、无量纲的指数，其取值范围为[0, S]。若Pie＞1，表示地形位e是该种植被变化趋势分布的优势地形位，Pie越大，优势越明显。反之若Pie＜1，表示地形位e不是该种植被变化趋势分布的优势地形位（喻红等，2001）。

2 结果与分析

2.1 NDVI变化趋势的空间分布特征

由图3可知，贵州省NDVI年平均值变化趋势的空间分布存在显著差异。呈现上升趋势的栅格占总面积的 48.2%，其中显著上升的栅格占 10.2%（Z≥1.64，P≤0.05），极显著上升的栅格占 3.7%（Z≥2.32，P≤0.01）。这些栅格主要分布在贵州省西部的威宁、毕节，西南部的盘县、普安、兴仁、贞丰，北部和东北部的沿河、印江、铜仁东部等地区，这些地区的主要植被类型为灌丛、草地和农田。呈现下降趋势的栅格占研究区域总面积的51.8%，其中呈现极显著下降趋势的栅格占 3.8%（Z≥2.32，P≤0.01）。这些栅格主要位于贵阳、六盘水、遵义、凯里、铜仁等城市周边，以及贵州省北部、东北和东南部等区域，这些区域的植被类型主要包括常绿阔叶林、落叶阔叶林和针叶林。

图3 研究区NDVI变化趋势空间分布图Fig. 3 Spatial patterns of NDVI change trend

2.2 植被变化趋势的地形效应

根据贵州省农业地貌区划，将研究区划分为 3个地貌区（图1a）。本文以3个地貌区为单元计算地形位指数，并基于等距划分将地形位指数划分为12个区段，依次从1-12进行编号，其中1-4为低地形位，5-8为中地形位，9-12为高地形位。3个地貌区域的地形位指数均呈正态分布，地形位指数集中分布于4-7区段（图4a-图4c）。

图4a反映出，贵州省东部山地丘陵区内，NDVI呈显著增加趋势的栅格主要分布在低地形位上，随着地形位升高，NDVI呈增加趋势的栅格的分布优势快速下降。NDVI呈显著降低趋势的栅格主要分布在地形位1区段和7-12区段上。由图4d可以看出，在东部山地丘陵区内，地形位1区段上土地利用类型以建设用地为主，2-4区段上是耕地和植被的混合分布区，5-12区段是森林的优势分布区。

图4b可以看出，中部丘陵山地地区NDVI呈显著增加趋势的栅格的优势分布区位于 1和 5-7区段，但是5-7区段的分布指数仅略高于1；呈显著降低趋势的栅格主要分布在 1-5区段和 9-12区段上，其中1区段上，NDVI呈降低趋势的栅格的分布优势大于呈增加趋势的栅格分布优势。图4e显示，1-3区段为建设用地和耕地的分布区，4-7为林地和耕地的优势分布区，8以上集中分布森林。

图4c显示，西部高原山地地区NDVI变化趋势的地形分布结构与中部丘陵山地区域基本相同，显著增加的栅格以5-7区段为优势分布区，显著降低的栅格以1-4区段为优势分布区。土地利用结构也是 1-4区段是建设用地的优势分布区，5-7是林地和耕地的优势分布区，其中包括灌木林、草地和稀疏林；10以上区段是森林的优势分布区。

3 讨论

3.1 植被变化空间特征的一致性分析

本研究得到的贵州省NDVI变化趋势的空间分布特征与前人的研究结果基本一致。NDVI呈显著增加趋势的栅格主要分布在贵州省威宁、毕节、盘县、兴义、贞丰、沿河、务川、德江和铜仁东部等区域；显著降低的栅格主要分布在赤水、习水、施秉、黔东南等地区及贵州省主要的城区周边，如贵阳市、六盘水市、遵义市、毕节市和凯里市。为了对NDVI的变化趋势进行一致性分析，计算了2000-2010年NDVI年平均值的变化趋势比例。在0.1显著性水平（Z=1.23）上，3种变化趋势类型（显著下降、不显著变化和显著增加）占比分别为6.9%、55.3%和 37.8%，与Cai et al.（2014）关于贵州省NDVI呈显著降低趋势的栅格比例和空间分布的研究结果基本一致。

图4 3个地貌区内NDVI变化趋势和土地利用/覆被类型的地形梯度分布特征Fig. 4 Distributions characteristic of NDVI changes trends and landuse/cover in three landform regions

本研究分析得出，贵州省NDVI变化趋势的垂直空间分布特征为：低地形位和高地形位是 NDVI呈显著降低趋势的栅格的优势分布区，中地形位是NDVI呈显著增加趋势的栅格的优势分布区。这一结果与Tong et al.（2016）以贵州和广西为研究区得到的研究结果基本一致。

3.2 人类活动对植被变化趋势地形分异的影响

人类活动对植被的扰动主要是通过建设开发、开垦和砍伐造成的，建设开发和开垦使土地利用类型从有林地或灌木林地转变为建设用地或耕地，造成生态系统结构发生变化，在扰动前后NDVI时间序列存在变化趋势的突变点（Verbesselt et al.，2010）。根据 Mann-Kendall趋势分析，贵州省 3个地貌区的低地形位上，扰动后无修复和扰动后修复的栅格占低地形位NDVI呈显著降低的栅格总数的 34.4%和 14.8%，NDVI平均降低速率分别为-0.012 a-1和-0.0091 a-1。从土地利用变化角度分析，低地形位上NDVI呈显著降低的栅格中，2000-2015年转变为建设用地、开垦和林地退化的面积比例达到 33.5%（表 1）。因此可以看出，人为活动造成的植被扰动是低地形位成为NDVI呈显著降低的栅格的优势分布区的主要原因。

贵州省中地形位上 NDVI呈显著增加的栅格中，扰动后修复和未扰动的栅格比例分别为30.7%和 43.2%，NDVI变化速率分别为 0.0163 a-1和0.0062 a-1。同时，这部分栅格中存在NDVI快速升高的突变性变化趋势，说明中地形位成为NDVI显著增加的栅格的优势分布区是由人为活动促使形成的。从土地利用变化角度分析，2000-2015年中地形位上NDVI显著增加的栅格内生态修复比例达到 30.7%（表 1）。根据贵州省土地利用的空间分布可知，中地形位是旱地、灌木林和草地的优势区域（郜红娟等，2015），也是石漠化治理、退耕还林还草和封山育林工程的主要实施区域（Yang et al.，2014），这些工程的实施都可以使NDVI快速上升。尽管中地形位是NDVI显著增加栅格的优势分布区，但是分布优势不强，这主要是由于在中地位上还存在土地开垦和植被退化两种土地利用转变类型，二者比例与生态修复比例基本相等。理论上中地形位上NDVI呈显著增加和显著降低的栅格的分布优势应该相等，但是土地开垦属于破坏原有植被覆盖后的重新耕种，扰动类型可以表现为扰动后修复，NDVI呈先降低后增加的变化趋势，这导致 NDVI显著增加栅格的分布优势略高于显著降低。虽然本研究反映出生态修复工程的实施促使中地位上NDVI值显著上升，但必须密切关注生态保护和修复工程的具体实施效果，防止边治理边破坏的双轨变化趋势的出现。

表1 3个地形梯度上植被扰动及土地利用变化表Table 1 Vegetation disturbance and the landuse change in the low, middle and high terrain gradient

3.3 自然因素对高地形位植被显著降低的影响

贵州省东部和中部地貌区的高地形位上是NDVI呈显著降低趋势的栅格的优势分布区。从NDVI时间序列轨迹看，86.2%属于未扰动类型，NDVI平均变化速率为-0.008 a-1。这部分栅格内，2000-2015年土地利用类型未发生变化，因此可以得出，人为因素不是高地形位上NDVI显著减少的直接原因。

图5 高地形位夏季（a）和春季（b）NDVI值与气温的相关关系Fig. 5 Percentage of pixel of correlation between NDVI value and seasonal average temperature

由图5可知，东部和中部的高地形位区内，夏季NDVI值与夏季温度呈负相关关系的栅格比例分别为97%和96.5%，春季二者关系与夏季基本相同。贵州省 2000-2015年夏季气温呈显著增高趋势，气温趋势斜率为0.0952 ℃·a-1，R2=0.49；夏季降水趋势斜率为-0.992 mm·a-1，R2=0.0022，夏季气候变化呈暖干趋势。土壤水分是植被生长直接有效的水源，植被对土壤水分的缺乏具有快速、强烈的反应（Rodriguez-Itur et al.，2001）。当气温上升，降水减少或不变时，土壤水分减少，植被长势减弱，NDVI呈下降趋势。但由于植被对气候变化的响应是一个缓慢的过程（De Jong et al.，2012），这种变化可能被人为活动促使的植被快速增长趋势所掩盖，因此，贵州省NDVI的二维线性趋势表现为上升趋势。最新的关于全球植被变化趋势的研究表明，20世纪90年代以来北半球所有纬度带上，由气候变干导致植被变褐的趋势有所增加，但这种增加趋势被变绿趋势所掩盖（Pan et al.，2018），本文通过植被变化趋势的地形梯度分异进一步证实上述研究的结论。同时表明，随着地形梯度的加大，人为活动减少或消失，植被变化由人为和自然双重因素控制转向单一自然因素控制，可以更加明显的反映出植被对气候变化的响应。

4 结论

（1）2000-2015年贵州省NDVI变化趋势的空间分布呈现为：显著降低的栅格主要分布在北部、东北和东南部的非喀斯特区域，显著增加的栅格主要分布在贵州省的西部、南部和东北部等地区。

（2）贵州省 NDVI变化趋势的地形分异规律显著。低地形位区是NDVI显著降低的优势区，中地形位区是NDVI显著增加的优势区域。东部和中部地貌区的高地形位区是NDVI显著降低的优势分布区。

（3）植被变化趋势的地形梯度分异特征能够反映出植被变化的直接驱动力，人为扰动是低地形位上NDVI显著降低的直接原因，石漠化治理工程和生态修复措施的实施是中地形位上NDVI显著升高的主要原因；高地形位区NDVI显著降低与春、夏季气候呈暖干的变化趋势有关。

（4）生态环境保护过程中，针对不同地形梯度上的植被应该采取不同保护措施：低地形位区应该重点关注人为活动对植被的扰动；中地形位区应重点关注生态工程的治理成效和不合理的土地利用方式，防止植被恢复与退化的同时发生；高地形位上应通过实地监测密切关注植被对全球气候变暖的响应。