1999—2015年长江流域上游植被覆盖特征及其对气候和地形的响应
2019-09-17
(成都理工大学 国土资源部地学空间信息技术重点实验室，成都 610059)
1 研究背景
植被是生态系统的重要组成部分，在拦截降水、减缓径流、防治荒漠化和水土保持等方面发挥着重要作用，并且受到了广泛关注[1]。因此，利用有效的监测指标分析区域植被覆盖状况，对实现地区土地资源的科学管理和环境保护相关政策的制定具有重要的参考价值。归一化植被指数(NDVI)可以补偿光照条件、土壤颜色、坡度和观测方向差异造成的变化，被广泛认为是监测区域和全球植被变化的有效指标[2]。气候变化给植被生态系统带来了巨大的影响[3]。现有研究表明，植被覆盖与气候因子之间具有密切的内在联系，植被是全球气候变化的敏感指示器，尤其是温度和降水因子[4]。同时，许多研究也表明地形亦是影响植被覆盖变化的重要因素之一[5-6]。目前关于植被覆盖对气候响应的研究多集中于单一动态变化的角度，通常将整个研究区视为一个统一体，通过计算整个区域时间序列的平均值探讨植被与气候之间的关系[7-9]，而该方法却忽略了不同地区温度和降水在空间和时间上的异质性，无法从多角度判断不同地区植被对气候的动态响应及静态响应，并且没有根据区域差异性进行分区讨论。此外，许多植被对地形的响应研究多将地形因子如高程、坡度进行区间划分，通过统计不同区间范围内植被覆盖状况研究两者之间的关系，分级区间不同，得到的结果也不同[10-12]。因此，该方法不仅受分级方法、分级区间的主观影响较大，而且无法更加准确地显示出两者之间的连续变化关系，对结果造成影响。为此，本文特引入格网法通过散点图进行数据分析。目前散点图在许多领域得到应用[13]，但是把研究区划分成格网单元，将散点图运用于植被覆盖和气候与地形因子之间关系的研究还十分罕见。
本文主要解决以下问题：掌握1999—2015年长江流域上游的植被覆盖空间格局及变化特征；探讨植被覆盖对温度、降水变化的动态响应，并进行分区探讨；引入格网法探讨植被覆盖对温度、降水、高程和坡度在空间分布上的静态响应。
2 研究区概况
长江流域上游(24.46°N—35.75°N，90.54°E —111.46°E)位于中国大地形的第一、二阶梯，是整个长江流域的生态屏障，具有重要的水源涵养功能。该地区地形复杂，地势西高东低，西部为青藏高原，海拔在6 000 m以上，主要植被为高原草甸；东部主要为四川盆地的平原丘陵地带，海拔仅为几百米，主要植被类型为阔叶林、农作物等。年均气温在8～10 °C之间，年降水量约为900 mm。东西地形、植被、气候状况存在显著差异。目前该地区由于气候变化和人类活动的影响，植被覆盖发生了显著变化。大规模的生态恢复项目，例如“退耕还林”计划使得部分地区植被覆盖率提高。但同时，人口增长、经济发展和城市扩张等使得部分地区植被大量退化。因此，有必要全面分析长江流域上游植被覆盖特征，并了解植被对气候和地形的响应。
3 数据与方法
3.1 数据来源与处理
(1)遥感数据：研究区1999—2015年的SPOT NDVI数据集，空间分辨率为1 km，时间分辨率为10 d，采用最大值合成法(MVC)对其进行处理得到各年度NDVI最大值，MVC可以降低大气效应、减少云污染，提高数据的准确性[14]； 研究区2000年Landsat 5 TM及2015年Landsat 8 OLI影像，空间分辨率为30 m，经过波段合成、裁剪等操作得到；地形数据空间分辨率为90 m的数字高程模型(DEM)。
(2)气象数据: 1999—2015年年均温度、年降水量的栅格数据是基于全国2 400多个气象站点的日观测数据，通过整理、计算和空间插值处理得到。考虑到常用的插值方法由于地形起伏等因素使得数据精度较低，本文的数据使用ANUSPLIN(澳大利亚利用薄板平滑样条函数对气候数据进行空间插值开发的一个软件程序)进行插值，其已经被证明比其他方法更适合于气候的空间插值[15]。
遥感数据、气象数据和DEM来源于中国科学院资源环境科学数据中心(http://www.resdc.cn)；长江流域上游矢量边界数据来源于国家地球系统科学数据共享服务平台(http://www.geodata.cn)； landsat影像数据获取自地理空间数据云(http://www.gscloud.cn/)。此外，为提高研究结果的准确性，空间分辨率统一为1 km，投影坐标系采用Krasovsky_1940_Albers等面积割圆锥投影。
3.2 研究方法
3.2.1 趋势分析法
根据年度时间序列与NDVI的相关性，可以确定植被覆盖的年际变化率[16]。以像素为基础，得到植被覆盖的变化趋势，其表达式为
(1)
式中：Ks表示植被变化趋势的斜率；n为年份数， 本研究期为1999—2015年，n=17；NDVIi为第i年的NDVI。正(负)斜率表示植被覆盖增加(减少)，Ks的绝对值越大，表示变化越大。变化趋势的显著性采用相关系数显著性检验表进行检验。参照已有研究[2]，植被覆盖变化趋势分为5个类别：显著退化(Ks<0，p<0.01)；轻度退化(Ks<0，0.01
3.2.2 偏相关系数法
地理系统是一个多因素的复杂庞大系统，其中一个因素的变化会影响到其他因素的改变。在多因素体系中，研究一个因素与另一个因素的相关性时，不考虑其他因素的影响，这种相关性称为偏相关[17]。偏相关系数用来表示多个因素动态变化之间的相关性，本文基于像元利用偏相关系数法分析1999—2015年植被覆盖与年均温度和年降水量动态变化之间的偏相关性，即
(2)
简单相关系数计算公式为
(3)
3.2.3 格网法
研究区总体面积约98万km2，以栅格作为评价单元进行统计分析的方法数据量过大，且容易因数据误差导致分析结果不准确。为有效地避免上述问题，通过在研究区范围内构建格网评价单元，将研究区切分成数个相等的单位地块，分析植被覆盖对气候和地形因子的静态响应。结合研究区的实际情况，以面积为252.86 km2的格网单元地块作为基准单元，对于边界部分的处理，遵循质心在格网内即算作一个分析单元，质心不在格网内则不纳入统计的原则，最终整个研究区被划分为3 898个基准单元进行统计分析。在建立上述基准单元的基础上，基于ArcGIS 10.2进行区域统计，计算每个单元块内所有像元平均值，将空间数据转化为统计数据，用表格的形式记录下每一个分析单元的平均值。最后，用其中一个元素基准单元块内的平均值作为横坐标，另一元素基准单元块内的平均值作为纵坐标，绘制出分析变量与因变量之间关系的散点图，并以此为依据进行分析、讨论。
图1 1999—2015年研究区平均年度NDVI分布Fig.1 Average annual NDVI distribution in the study area from 1999 to 2015
4 结果与分析
4.1 植被覆盖特征分析
4.1.1 植被覆盖空间分布特征
基于像元计算得到1999—2015年长江流域上游平均年度NDVI的空间分布图(图1)。参照已有研究并结合地区的实际情况将植被覆盖分为5个类别[18-19]，如表1。研究区1999—2015年的平均NDVI为0.68，整体上表现为东部地区植被覆盖度高于西部地区。植被高覆盖区(NDVI>0.7)占总面积的66.14%，主要分布在中东部；植被低覆盖区(0.15≤NDVI<0.3)仅占4.54%，主要分布在西部青藏高原高寒草甸地带，平均NDVI值为0.26，草地植被覆盖低，这归因于草本植被生物量积累有限。此外，中东部的成都市、重庆市等人口密集，工业和经济发达的地区植被覆盖也较低。
表1 植被覆盖分类标准及数据统计Table 1 Classification standards and statistics of vegetation cover
图2 1999—2015年植被覆盖动态变化及呈贡区 Landsat影像Fig.2 Dynamic changes in vegetation cover from 1999 to 2015 and Landsat images of Chenggong
4.1.2 植被覆盖年际变化特征
基于像元的长江流域上游植被覆盖变化趋势空间分布如图2所示。由图2(a)可知，NDVI变化趋势具有空间差异性。1999—2015年该地区60.59%的植被覆盖状况有所改善，主要分布在东部地区，其中48.17%的地区显著增长，占植被覆盖变化的比例最大，12.42%为轻度增长。变化不明显的区域占36.44%，主要分布在东部的成都平原和西部的青藏高原地区。相反，仅2.97%的植被已发生退化，占很小的比例，呈星散状分布，其中一部分发生在成都、重庆、昆明等城市扩张的高度工业化地区，另一部分植被退化区主要分布在西北部。通过将NDVI计算结果与分辨率更高的Landsat影像对比分析发现，植被覆盖退化受人类活动影响较大。以研究区南部显著退化的区域(昆明市呈贡区)为例，通过选取2000年和2015年2期的Landsat影像(图2(b)、图2 (c))对比发现该地区植被减少主要由城市扩张造成。2011年根据国务院批准，对呈贡撤县设区，昆明市政府和主城区大学校园整体搬迁呈贡，并建立呈贡大学城片区(大学城占地43.15 km2)[20]。此后，新区经济发展导致流动人口增加，全区总人口相比2000年年平均增长率达到10.78%，人口的增加又使得建筑面积不断增加，植被不断减少。可见，政府的规划决策对植被覆盖变化产生了重要的推动作用。
4.2 植被覆盖对气候和地形的响应
4.2.1 植被覆盖对气候的响应
4.2.1.1 动态响应
焦珂伟等[4]指出：在气候变化影响下，植被总是在不断地去适应外界条件的改变使其自身活动更为有利，该过程是动态、非线性的，尤其在极端气候条件下，植被活动随气候因子的变化过程存在明显的阈值。当温度和降水等气候条件的改变超过植被的适应能力时，植被将会在结构功能等多方面受到抑制作用[21]。气候变化对植被活动的影响机制如表2所示。
表2 温度和降水变化对植被生长的影响机制[4]Table 2 Effects of temperature and precipitation changes on vegetation growth[4]
为探讨植被覆盖对气候变化的动态响应，本文计算了每个像元的NDVI与年均温度、年降水量之间的偏相关系数，如图3所示。年均温度与NDVI的偏相关系数范围在-0.49～0.65之间(图3(a))，年降水量与NDVI的偏相关系数在-0.40～0.43之间(图3(b))。通过统计，年NDVI与年均温度呈正相关的地区占61.40%，与年降水量呈正相关的地区占51.69%，但在不同区域相关性方向及程度不同，这与卓嘎等[22]利用MODIS NDVI数据集得到的不同地区植被覆盖与温度和降水的相关性不同的结果相符。
综上所述：植被覆盖对温度和降水变化的动态响应具有较强的空间异质性，不同地区表现出不同的动态响应关系。根据年均温度、年降水量与NDVI在空间上偏相关性方向的不同，对其进行分区(图3(c))讨论。植被覆盖对温度、降水变化的动态响应分区及统计如表3所示。
图3 1999—2015年年均温度、年降水量与NDVI变化 的偏相关系数及响应关系分区Fig.3 Partial correlation coefficient and response pattern of annual temperature, annual precipitation and NDVI variations from 1999 to 2015
分区温度降水面积比率/%分区1++15.79分区2+-45.61分区3-+35.90
注：“+”表示偏相关系数>0；“-”表示偏相关系数<0。
分区1：该区域占总面积的15.79%，分布在西部，NDVI与年均温度和年降水量的偏相关系数较高，且都为正相关。分析认为该地区为青藏高原草甸区，海拔高，气温降水较低，人类活动少，气候是影响该地区植被生长的主导因素。分区2：该区域占总面积的45.61%，分布在西南部，为青藏高原与四川盆地过渡地带，地貌主要为山区，植被以温带针叶林、阔叶林和草地为主。该地区植被覆盖变化与年均温度呈正相关，与年降水量呈负相关，但与温度的偏相关程度小于降水。根据已有成果[23-26]分析认为这是由于高山植被生长的可用水量取决于降水和大气蒸发需求的季节性因素，植物生长所需要的水分充足，若降水继续增加将会使植被光照条件下降，减弱其光合作用。此外，高山地区地形起伏度大，降水量增加会导致泥石流、滑坡等自然灾害，加重土壤侵蚀，土壤有机质减少，从而与植被覆盖变化表现出负相关。分区3:该区域占总面积的35.9%，主要分布在东部，NDVI与年均温度呈负相关，与年降水量呈正相关，但年NDVI与降水的偏相关程度低于NDVI与温度的偏相关程度，表明该地区NDVI对降水的变化相对不敏感。分析认为，该地区主要为四川盆地平原丘陵区，城镇化程度高，人类活动频繁，主要用地类型为城镇用地或灌溉农田，植被类型大多为农作物，在农作物生长季节主要用地下水灌溉，导致植被覆盖与降水量的相关性较弱。
图4 1999—2015年平均NDVI与年均温度、年均 降水量相关性散点图及年均温度、年均降水量空间分布Fig.4 Scatter plot of the correlation between annual average NDVI and temperature, annual average NDVI and precipitation, and spatial distribution of annual average temperature and precipitation from 1999 to 2015
4.2.1.2 静态响应
通过格网法，建立1999—2015年长江流域上游17 a的平均NDVI与年均温度及年均降水量相关性的散点图(图4(a)、图4(b))，从静态角度探讨植被覆盖对气候因子的响应。图4(a)显示了每个格网单元的NDVI与年均温度之间的关系，NDVI值随着温度的升高呈连续非线性变化，整体表现出二次函数关系(R2=0.713 5，p<0.001)，当温度升高到一定阈值时，NDVI随温度变化不再明显。观察散点图，取阈值为2.5 ℃，当温度≤2.5 ℃时，NDVI随温度升高而迅速增加，通过线性拟合(R2=0.56，p<0.001)发现NDVI以0.06/(1 ℃)的速率增加。通过图4(c)发现，在空间分布上，年均温度≤2.5 ℃的区域主要分布在研究区西部青藏高原高海拔地区，>2.5 ℃的区域主要分布在研究区东部四川盆地低海拔地区。据图4(b)，NDVI与降水的关系和NDVI与温度之间的关系类似，都为连续非线性关系(R2=0.700 8，p<0.001)，观察散点图取阈值为800 mm，当年降水量≤800 mm时，NDVI随年降水量升高而迅速增加，通过线性拟合(R2=0.69，p<0.001)发现NDVI以0.001/(1 mm)的速率增加。据图4(d)发现，以800 mm年降水量为界限的空间分布特点与温度一致。该现象表明植被覆盖对温度和降水的静态响应存在一定的阈值，根据该阈值分割得到的空间分布特征都表现为：西部响应明显，东部响应不明显。
4.2.2 植被覆盖对地形的响应
通过格网，建立了1999—2015年长江流域上游17 a的平均NDVI与高程及坡度相关性的散点图(见图5)。
图5 1999—2015年平均NDVI与高程、坡度相关性 散点图及高程、坡度空间分布Fig.5 Scatter plot of the correlation between average NDVI and elevation, NDVI and slope gradient, and spatial distribution of elevation and slope gradient from 1999 to 2015
由图5(a)探讨植被覆盖对高程的静态响应，多年平均NDVI值在R2=0.617 6的拟合度(p<0.01)显示不同海拔高度呈现连续非线性下降趋势，与NDVI对温度和降水的静态响应整体变化趋势相反，但同时也存在一定的响应阈值。取高程阈值为4 000 m，低于该值时，随着高程的增加，NDVI变化不明显，该区域多年平均NDVI >0.75；而高于该值时，NDVI急速下降。据图5(c)发现，海拔高度>4 000 m的区域主要分布在研究区西部青藏高原及附近的山区。综合表明高海拔地区会抑制植被的生长，这主要是因为高海拔地区气候寒冷，高原山地地区常年冰雪覆盖，冻融侵蚀严重，使土壤的物理性质遭到不同程度的破坏，影响坡面土体的稳定，加重土壤侵蚀[27-28]。
图5(b)表明NDVI与坡度之间的关系与高程的连续非线性变化趋势不同，呈不连续的非线性变化趋势。在坡度<5°时，NDVI<0.6的区域随着坡度的增加呈现明显的增加趋势，而当NDVI>0.6时，随着坡度的增加，NDVI无明显变化趋势。从坡度空间分布来看(图5(d))，坡度大的区域主要分布在研究区中部和西南部山区，而NDVI<0.6的区域主要分布在研究区西部青藏高原高海拔地区， NDVI>0.6的区域主要分布在研究区中部及东部低海拔地区。因此，结合由散点图得到的结果可知，在西部青藏高原高海拔高寒草甸区，NDVI对坡度变化的响应明显，在东部低海拔地区，NDVI对坡度变化无明显响应。
5 结论及讨论
(1)长江流域上游植被覆盖呈现自西向东逐渐增加的梯度分布特征。1999—2015年间，植被覆盖变化整体呈增长趋势，增长面积为60.59%。通过比较NDVI计算结果和Landsat影像发现，在东部及南部地区，植被覆盖显著退化主要由人类活动尤其是城镇扩张造成。
(2)植被覆盖对气候的响应具有空间异质性。从动态角度考虑，植被覆盖与温度和降水的变化在不同地区的偏相关性方向和程度不同，不同地区的响应关系不同；从静态角度考虑，植被覆盖对温度和降水的响应表现为连续非线性增加趋势，存在一定的阈值界限。
(3)植被覆盖对地形的静态响应为非线性。植被覆盖在不同海拔高度的变化趋势不同，其阈值约为4 000 m，当海拔超过该阈值时，随着海拔的升高，NDVI显著下降。植被覆盖对坡度的响应为不连续的非线性变化关系。整体来看，植被覆盖对温度、降水、高程、坡度的静态响应在空间分布上都表现出西部青藏高原地区响应明显，东部四川盆地地区无明显响应的一致特征。
(4)与将全部数据分成不同区间进行研究的传统方法相比，本研究中引入格网法进行空间统计，能较好地反映植被覆盖对气候和地形因子在空间分布上的静态响应。格网法是一种有效和可靠的工具，为未来植被覆盖方面的研究提供了一种新思路。
本研究中使用的遥感数据的空间分辨率为1 km，时间分辨率为1 a，这种相对较低的分辨率在分析中可能会忽略一些细节。此外，在选取响应要素方面，其它气候因素(日照时数、湿度等)和地表因素(岩性、土壤等)也是影响植被覆盖的重要因素，但由于数据获取原因，只选取了常规的因子进行分析。未来的研究需要探索高分辨率及其它更多因素对植被覆盖的响应。
