近35年青藏高原植被带变化对气候变化及人类活动的响应

魏彦强，芦海燕，王金牛，孙 建，王旭峰

（1.中国科学院西北生态环境资源研究院 甘肃省遥感重点实验室，甘肃 兰州 730000；2.兰州大学 经济学院，甘肃 兰州 730000；3.兰州财经大学 会计学院，甘肃 兰州 730101；4.中国科学院成都生物研究所，四川 成都 610041；5.中国科学院地理科学与资源研究所，北京 100101）

近些年来，当代全球暖温化问题(contemporary global warming)，尤其是在高海拔、高纬度等气候变化敏感区的环境变化已引起了全球的广泛关注[1-2]。在高纬度和高海拔地区的积雪、冰川、冻土和植被等对气候变化的响应显得尤为明显，是较为典型和敏感的气候指示因子[3-4]。目前的研究结果表明，高海拔地区气候暖化较低海拔地区明显[5-7]，全面的暖化对寒冷的高海拔地区并不意味着都是负面的影响，例如增温使得严寒地区开始改善植物的生长环境，有利于植物向高海拔地区迁移和动物及人类的进入，是全球变暖中物种的避难所[8-10]。因此，对高海拔地区气候变化的研究在当代气候变化研究中具有代表性和典型性。

气候制约着植被的地理分布，植被是区域气候特征的反映和指示，二者之间存在密不可分的联系。植被带往往被认为是气候带的典型指示物。在大尺度的植被特征研究中，植被带的概念与因海拔、经度、纬度变化而产生的水热梯度差异相结合，以指示植被类型的划分[11-15]。由于植物的分布受控于气候模式，最低温度对限制物种的分布范围和向南北极高纬度地区及高海拔地区的极化扩展起着决定性的作用[16-17]。可以预测，全球全面暖化将改变物种分布的适宜环境和生存范围的极限，因而使得物种在暖化过程中得到重新分布[16-17]。当前在高海拔地区关于全球变暖引起的分布界限变化研究主要集中于以下几个方面：1)雪线/零物质平衡线 (snowline / equilibrium line altitude, ELA)随冰川退缩的变动[18-19]；2)树线/林线(timberline)和生态交错带(ecotone)随全球变暖向两极高纬度地区及高海拔地区的极化移动[20-22]；3)气候带(climate belt)或物种(species)向两极高纬度地区及高海拔地区的迁移[8, 16-17]；4)土地利用变化 (land-use and land-cover change, LUCC)[23-25]。Colwell等[16]和 Jump 等[17]研究认为，由于在垂直方向上气温递减率的存在以及在南北方向上温度随纬度的递减，在全球增温情形下，树线的分布范围将向高纬度高海拔地区扩展。植被的重分布及其他物种的迁移将随之发生[21, 24, 26]。与此同时，暖温化使得物种的适宜生存条件下限随之改变，人口的迁移及以游牧为主的牲畜也随之扩展其活动范围，表现为向两极及高海拔地区的迁入，人类活动增强[16, 27-28]。因此高山区和高纬度地区生态系统对全球暖化的响应往往表现为植被带的分布范围及其面积的变化，以及植被生长状况的变化。而目前在该领域长时间序列的监测和研究，归一化植被指数 (normalized difference vegetation index, NDVI)被认为是一个很好的监测指标[9, 29]。

青藏高原 (Qinghai-Tibetan Plateau, QTP)因其对气候变化的高敏感性和在气候变化中其气候因子的指示意义较同纬度低海拔地区更具代表性和典型性而倍受关注[30-31]。根据最近的研究结果，青藏高原地区在过去的几十年中经历了一个加速的暖化过程，其降水变化较为平缓，而潜在蒸散能力则有下降趋势[32-34]。根据最新对青藏高原冰川及冻土退化和环境效应的调查显示，从1950s开始青藏高原暖温化趋势明显，其冻土和山岳冰川遭受了大规模的退缩和面积缩减[35-37]。大量事实已表明，因高原的全面暖化，其雪线/零物质平衡线及树线的位置因冰川退缩而得到了抬升，植被覆盖发生了明显的变化[19, 22]。气候条件的变化已经引起了其生态环境和生态系统的变化[38-39]，而人为因素如人口增加、基础设施工程建设、城镇化发展、过度放牧等的干扰已对生态系统产生了明显的影响，并在一些区域超过了气候变化影响本身[27-28, 40]。另外，大量的气候模拟情景和预测结果表明，青藏高原在21世纪将经历持续的暖化过程，并且会进一步缓解其严寒而干冷的气候环境[5, 41-42]。目前已有很多关于青藏高原植被变化及其对气候变化响应的研究[38, 42-44]，但通过林线及生态交错带位置变化而研究植被带变动，以及除了考虑气候因素之外，综合考虑人为活动影响的研究还相对较少[27-28]。且在人为影响因素的解释上，指标的选取具有较强的任意性，概念化、定性化的描述较多而定量化的数值模拟和分析较为少见[27-28, 35]，往往与气候影响难以剥离，缺乏客观性评价[43-45]。近些年来，随着人口规模、基础设施建设及城镇化的发展，越来越多的研究表明，人类活动增强是一个不可忽视的重要影响方面，在区域环境变化中往往超过 气 候 变 化 影 响 本 身[24, 27-28, 44, 46]， 而 目 前 在 青 藏 高原地区对该方面的研究相对较少，相关研究较为薄弱。因此，本研究从NDVI的变化出发，分析植被带变化趋势和气候变化、人类活动之间的关系，进而剖析植被带变化对气候变化及人类活动的响应。

1 数据和方法

1.1 研究区域

本研究选择气候变化较为剧烈且垂直地带性较为典型的青藏高原地区作为研究区。该区域海拔普遍较高，且拥有相对较为宽阔的海拔范围、垂直带谱较为发育，其对气候变化响应较为敏感。依据Zhang等[47]界定的包括其南部的喜马拉雅山脉、北部的昆仑山及祁连山、西至帕米尔高原及喀喇昆仑山、东至成都平原的范围目前在学界普遍比较接受。本研究以此为参考，南部以中国国界线为界，其北部及东部选择海拔在3 500 m以上的地区为其边缘。本研究中青藏高原总面积约2.57 ×106km2，平均海拔 4 390 m 左右，其经度在 70°-105° E，纬度为 25°-40° N(图 1)。

1.2 数据

1.2.1 NDVI数据

由于植被地上生物量及净初级生产力与归一化植被指数(NDVI)之间有很高的相关性，借用遥感数据对地上植被状况的分析通常以NDVI来代替植被覆盖状况[9, 48-50]。本研究选择美国航空航天局(NASA)全球监测与模型研究组GIMMS (Global Inventory Monitoring and Modelling Studies)提 供的NOAA AVHRR-NDVI-3g V1.0 数据集 (https://ecocast.arc.nasa.gov/data/pub/gimms/3g.v1/)，其空间分辨率为8 km × 8 km，属半月合成数据，时间跨度为 1981 年7月1日-2015年12月31日，共计69个NetCDF-4半月合成图像文件。该数据集的NDVI已经过几何精纠正、辐射校正、大气校正等预处理，且已采用最大值合成法(MVC)以减少云、大气、太阳高度角的影响。

1.2.2 气候及其他数据

本研究的气象台站数据由中国气象局国家气候中心提供，时间范围为1961-2015年。其中位于高原上的站点91个，因站点变动，记录缺失等原因，其中有完整月平均气温和月降水记录的站点87个。由于高原上的站点分布较少，尤其是广大的中西部地区长时间序列的站点较少，本研究另外选取了环境数据分析中心(CEDA)的长序列再分析数据 CRUTS 4.01 (http://catalogue.ceda.ac.uk/intro)。该数据集包括了每天的平均温度、降水等常用指标，数据包括1900年1月1日到2016年12月31日的逐日气象要素，空间分辨率为 0.5° × 0.5°。人口、牲畜数量等经济统计数据均来自各个省市的统计年鉴。由于该类数据均以行政单元(县)来统计，但青藏高原是一个自然地理边界，在高原上的部分以其所占总县域面积的比例来进行分割。另外DEM 数据由 NASA提供，空间分辨率为90 m × 90 m。植被覆盖/土地利用类型数据为基于GLC2000的中国科学院1∶100万土地利用数据集(图1A)；生态区划数据采用Zheng[51]拟定的中国地理生态区域系统的框架方案，将青藏高原分为11个自然区(图1B，表1)，包含了3个温度带和4个干湿区。

图1 青藏高原地区根据GLC2000分类的植被覆盖/土地利用类型(A)，多年平均NDVI以及依照郑度等[51]拟定的生态气候分区(B)Figure 1 Vegetation types/land covers in the Tibetan Plateau according to the GLC2000 categories (A),the mean NDVI in recent years and eco-climatic zones in the TP according to [51] (B)

1.3 数据处理

首先对GIMMS-NDVI的半月合成数据通过国际通用的最大值合成法 MVC (maximum value composites)合成月NDVI，该方法可以进一步消除云、大气、太阳高度角等部分干扰。并利用12个月月数据合成4个季节NDVI数据和年数据。合成公式为：

式(1)中，NDVIi为半月NDVI数据或已经合成的月值数据，NDVIm为合成后的NDVI季节数据或年数据。考虑到青藏高原地区植被覆盖普遍较低，且大多为温带草原及高寒草甸，本研究取NDVI的阀值为0.05[9, 52]，即小于此值的地区被认为是无植被带覆盖。为清楚地了解NDVI在空间上的变化趋势，逐栅格的趋势分析能直观地反映其在空间上的变化趋势。一般地，可以通过计算植被绿度变化率 (greenness rate of change, GRC)来表征植被带变化趋势[10, 52]。而常用的则是采用容易实现的最小二乘法拟合年均NDVI的变化斜率[10]。但最小二乘法线性回归的缺点则是样本量较小时，极易受到异常观测点的影响。本研究采用了Theil-Sen线性趋势分析法。该方法是一种非参数的估计算法，能有效消除长时间序列分析中的异常值的影响[53-54]，是对最小二乘线性回归方法的改进。其计算公式为：

式(2)中：slope为植被绿度变化率，变量i，j为年的序号，即i，j= 1, 2, …, 35。n为研究的时间序列长度，n= 35；NDVIi为第i年的年最大 NDVI值。其中slope＞ 0说明NDVI在这35年中的变化趋势是增加的，反之则是减少的，其值的大小反映出植被带变化的剧烈程度。以上线性趋势分析均在PyCharm平台上利用Python 3.7.1计算完成。

2 结果

2.1 植被带的年际变化

时间序列的变化分析可以从植被带演化更替的角度来分析整个青藏高原整体的状况。为了综合分析植被带随时间的变化，本研究计算了整个青藏高原NDVI平均值的历年变化情况，并按季节平均值统计其在过去35年的变化。从统计结果来看，NDVI在以月为尺度的变化中呈现出倒“U”字形的周期性变化特点，即年内的变化在冬春两季较低，平均值约为0.16，而夏秋两季达到年内植被生长的最好水平，平均值夏季为0.29，秋季为0.26左右。年内的变化显现出地上生物量与四季的更替及气温的变化较为一致的特点。从总体变动趋势上看，最近35年的NDVI年际变化较小，且波动不大，总体上处于平稳的波动态势。

表1 青藏高原生态地理系统分区Table 1 Eco-climate zones in the Qinghai-Tibetan Plateau

由于青藏高原整体上全年气温普遍较低，冬季漫长而严寒，生长季(5-9月)相对短暂，因此植被带的发育主要集中在以夏季为主的生长季。从生长季的NDVI变化趋势来看，过去35年生长季的NDVI呈现出一定的波动性，但总体趋势平缓，最大值出现在1988年，为0.296；最小值出现在1995年，为0.267；多年生长季平均值为0.282左右，有轻微的上升趋势。以上结果表明，一方面，较低的NDVI反映了青藏高原总体上植被覆盖较差，整个柴达木盆地、藏北高原等地区植被分布稀疏，几近于裸地。由于生长季节较短，地上生物量普遍相对较小，使得区域平均值较低。另一方面，青藏高原内部植被类型多样、植被带较多、生态气候区种类较多，气候的垂直地带性分布较为明显，复杂的地形和气候类型状况使得区域间的植被发育状况差异明显。有些因增温而使得植被带生长条件改善的地区，很可能因一些地区因生长环境的恶化以致植被覆盖减少，但在整个高原求其平均值时而被相互抵消掉，使得总体上植被带变化趋势波动性不大。因此，以青藏高原为整体来研究其平均的植被带生长状况的变化显得意义甚微。为此，在空间上分区域和植被带类型的生长状况变化研究显得颇为必要。

2.2 植被带的空间变化

2.2.1 总体变化

逐栅格的趋势分析能直观地反映植被带在空间上的变化趋势，能较好地反映出每个栅格在时间序列上的变化差异，因此是在空间上研究植被带对气候及人为因素响应差异的有效途径。利用公式(2)可以计算出最近35年植被带的空间变化趋势(图2)。

图2 青藏高原1981-2015年间植被带生长状况变化情况Figure 2 Vegetation trends during 1981-2015 in the Tibetan Plateau

从整个区域来看，大部分地区的植被覆盖趋于好转。在高原西南部及中部、柴达木盆地周围的广大地区植被都出现了增长情况。部分地区的年增长率达到了0.002 5及以上。尤其是在H1C1、H1C2两个高原亚寒带高寒草甸草原区及H2C2藏南高山灌丛草原区的增长尤为明显。而处于恶化趋势的植被带分布于山南林芝、云南北部、川西北及青海湖附近的低海拔地区，高原的东部边缘是这一趋势较为集中的区域。另外H1B1果洛那曲高原亚寒带灌丛草甸区的植被覆盖退化也十分明显(图2)。

2.2.2 不同分区间的差异

由于植被带的分布与海拔、植物种类、气候类型等紧密相关，本研究依据海拔带、植被覆盖种类和生态气候区统计在各个类型区中的植被带生长状况变化情况(图3)。

从生态气候区的内部统计结果来看(图3A1)，在整个青藏高原的11个生态气候区中，植被带生长状况平均水平处于较明显改善的有6个区，其中有H1C1、H1C2和H1B1三个高原亚寒带气候区以及H2C2、H2D1、H2D2三个高原温带气候区。而植被带平均生长状况趋于退化的有4个区，均为植被生长条件相对较好的亚热带阔叶林区(H3)和2个温度较高的高原温带针叶林(H2AB1和H2C1)及高原温带荒漠区(H2D3)。从统计结果来看，位于高原亚寒带的H1D1昆仑山高寒荒漠区变化不明显。由于趋于改善的生态区集中于高原亚寒带和高原温带，而退化的区域处于平均温度等相对较高的亚热带及高原温带地区，因此推测，全球暖化对改善高海拔地区严寒条件下的植被带生长状况有积极作用，温度升高有可能是树线和林线在这一地区升高、植被趋于改善的主导因素。而在气候条件相对较好的亚热带、高原温带地区的植被退化极有可能是升温后蒸散发加剧或人为因素影响的结果。但具体原因还待进一步详细分析。

为了进一步的深入分析植被带生长状况在空间上的变化差异，以植被覆盖类型来统计各个类型区内的植被带变化情况显得十分必要。植被带的分布与植被覆盖/土地利用变化之间有着直接的联系，植被覆盖/土地利用变化的种类界定了植被带的基本类型[16-17, 55]。本研究以 GLC 2000 的用地分类系统为参考，对其6个大类的25个小类进行重新合并调整，将水体和冰川等进行归并，最后形成12个土地覆盖类型。以此类型为基础，分别对各个区内的植被带生长状况变化率进行了统计(图3B1)。

从统计的结果可以看出，在12个大的植被覆盖/土地利用类型中，有4个显著增长的类型区，即L31高覆盖度(＞ 50%)草地、L32中覆盖度(20%～50%)草地、L4水体及L33低覆盖度(5%～20%)草地。其中平均增长速度最快的为L33低覆盖度(5%～20%)草地。从植被带的平均变化速度来看，保持基本不变的为L6沙漠/戈壁，这与实际的情况相吻合。而在其他所剩的几个植被带减小的区域中，L21森林、L22高覆盖度灌丛及L23疏林地等的减少较为明显，是3种减少最为剧烈的类型区。从植被覆盖在空间上的分布来看，增长的区域L31、L32、L33基本上是位于高海拔地区的草地，由于水体中在分类时已包含了冰川和积雪，其覆盖度的增加推测应当是气温增加后雪线升高、植被带的上线苔原线随之升高的结果；而高海拔地区整体上植被的好转则反映了对气候转暖后植被带生存状况改善、植被长势趋好及覆盖范围扩大的结果。在减少的几个区域中，L21森林、L22高覆盖度灌丛及L23疏林地均分布于高原边缘谷地及南部海拔相对较低、温度、降水等生长环境相对较好的地区。

在以上几个分析结果中，无论是基于生态气候分区的分析还是从用地类型分区分析，均显示出植被带的空间变化和海拔之间的密切联系，在不同的海拔上其变化趋势表现出明显的差异性。由于本研究的植被带位置移动假设是以植被带的垂直带谱变化为依据，据此，按海拔分区来分析植被带的空间变化显得很有必要。本研究以500 m等间隔做海拔分区，分别划分出各个范围的海拔带，即 0-500 m、500-1 000 m 和 1 000-1 500 m等不同的分区，分别对每个区内植被带的变化率进行了统计(图3C1)。从计算的结果来看，在划分的15个海拔区中，有5个均表现出明显的增长趋势，除去最高海拔的一个分区，几个植被好转的地区均集中分布于4 500-6 500 m的高海拔地区，并且占整个高原绝大部分的面积(图3C2)。而4 000-4 500 m似乎是植被带变化的一个临界位置，该范围内没有统计意义上的明显变化，而在该海拔以下各个植被带均表现出了较为明显的退化趋势。从平均减少速率来看，1 000-2 500 m 之间的 3 个较低海拔区减少较快，处于明显的退化水平。

从以上结果可以看出，植被带生长状况的变化与生长季节、植被带类型、生态气候区和海拔等有关。以青藏高原为整体的时间序列分析并不能体现出植被带的变化差异，而在考虑到按生态气候区、植被覆盖类型和海拔后，其变化的差异性过程才得到了明显体现。以分区为视角的空间差异分析体现出青藏高原植被带在最近35年的变化中与生态气候区、植被覆盖类型及海拔之间的密切关系。

3 讨论

3.1 植被带变化对气候变化的响应

一般认为，植被类型、生长状况及其分布受气候因子的直接影响，在特定的气候模式制约下在地理分布上形成特定的植被带和生态系统类型[16, 55]。因而植被带类型及其生长状态的变化是对一定气候模式和生长环境变化的直接反应，气候条件的改变使得植被的种类和分布范围随之发生相应的变化以适应新的气候特点[26, 55-56]。在众多的气候因子中，光照、温度、降水、蒸散发等被认为是植被生长的最重要且起决定性作用的气候因子，而通常以日照、温度和降水的影响最为直接[57]。考虑到数据的可获得性，本研究在分析植被带与气候变化的相互关系时，选取了与之相对应的年平均温度和年总降水作为主要气候参考因子。

从台站的气象记录来看(图4)，整个青藏高原的年平均温度(87个站点的平均值)呈现出明显的增长趋势，在Theil-Sen法线性拟合的结果中，年温度增长率为 0.035 3 ℃·a-1(R2= 0.688,P＜ 0.001)。表明整个青藏高原最近几十年经历了较为快速的增温暖化过程，这一暖化过程与NDVI缓慢增长的趋势基本一致，可以认为是对气候暖化的正向响应。

与温度相比较，年总降水量在总体上呈现出略微的增长，但增长幅度较小，为每年0.831 mm(R2=0.041,P＜ 0.32)，线性增长趋势不太明显，呈现出较大的年际变化差异。因此，从整个高原总体的NDVI及温度和降水变化的趋势来看，NDVI的增长趋势与温度有较好的一致性，但与降水增加的关系不太明显。

从同化的数据来看(图5A、B)，气温的增长率在H1C2、H1C1、H2C2、H2AB1及H2DI等藏北高原、高原南部河谷及柴达木盆地等地区相对较高，而在H2C1、H1B1高原的东北部及果洛那曲亚寒带高原草甸区的增长率相对较小。而降水在整个高原的增加则主要集中在西藏和青海两省 (区)的中部地区，尤其以H1C1、H2C2和H2C2等区的增加为主。而其他生态气候区内则基本上呈现出减少趋势。以温度和降水的变化区分气候变化的类型，即暖湿、暖干、冷湿和冷干4种气候变化模式。其中以暖湿为主的地区主要集中于高原中部及南部的H1C1、H2C2、H2C2等高海拔地区。在其他几个暖干化为主的区域中NDVI变化出现了不同的变化态势，其中H2D1柴达木盆地周边地区有较为显著的增长，但在H2C1及H1B1两区中总体上有不同程度的减少。因此，总体上可以判断，气候的暖湿化有利于植被的生长，植被带在暖湿气候条件下出现了明显的恢复和扩展；而暖干化则在带来部分地区植被好转的同时，退化较为明显。值得一提的是，在气候条件较好的高原南部林芝地区及滇西北的NDVI减小趋势明显。虽然年均温度上升较快，但近些年降水减少明显，尤其是近些年的夏季干旱极有可能是这一地区植被退化的主要原因[58-59]。

以上分析可以总结为，最近35年整个青藏高原的暖温化对植被状况的改善具有明显的促进作用，尤其是在高原中部及西南部的H1C1、H1C2、H2C2三个高海拔区暖湿化促进作用明显，但在H2AB1区中尤其是川西北地区，无论是年均温还是降水都不明显。由此可以推断，暖化的过程使得在高海拔地区的植被带生长环境改善，有利于植被的生长和植被带向高海拔，高纬度地区扩展。但在H2AB1区中，尤其是川西北地区，无论是年均温还是降水都与NDVI变化的关系并不明显，则说明这些地区的植被带变化原因中，气候要素变化并不是其最主要原因，可能是人为因素影响等其他因素作用的结果。

3.2 植被带变化对人类活动的响应

图5 2010s年平均气温T(A)和平均降水量P(B)与1980s年平均气温T(A)和平均降水量P(B)的差值Figure 5 The mean temperature (A2) and mean precipitation increasing maps (B2) were generated from the differences between the two periods (Mean T/P2010-2016 - Mean T/P1981-1990)

由于植被覆盖状况受气候条件和人为干扰两个方面的影响，尤其是在生态环境比较脆弱的青藏高原高寒草甸及草原生态区，其生长季较短，植被地上生物量较少，破坏后恢复能力弱，很容易受气候及人为因素的影响而使植被带退化[27, 55, 60]。除了气候因素之外，人为因素的影响其实一直为人 们 所 关 注[24, 27-28, 43]。 由 于 青 藏 高 原 地 区 是 我 国 四大牧区之一，畜牧业在该区域中的比重历来较高[28]。近些年来，随着气候转暖及西部大开发等国家发展政策的实施，青藏铁路的开通等使得这一地区的迁移及流动人口逐渐增多，人口规模持续地增加[35]，旅游业等第三产业的发展使得该地区的城镇化发展速度加快[44, 46]。由于人口的增加及城镇的发展，伴随人口增长的是对奶类、肉类、毛皮制品等畜产品的需求增加，并且牧民的收入与牲畜数量直接相关，收入需求的刺激使得畜牧业超规模发展，人类活动对高寒草原草甸生态系统的影响在持续加强[28, 43-44]。近些年来在牲畜数量持续增长的同时，畜牧业的发展对高寒草甸草地生态系统的压力进一步加大，过度放牧、鼠虫害等使草地退化的报道越来越多[35, 43-44]。为了分析人类活动对青藏高原地区植被带变化的影响，本研究选取了人为因素影响中最重要的两个指标，即反映人为因素规模和城镇规模的人口规模，以及代表畜牧业规模的牲畜数量，并以二者与NDVI变化间的相互关系进行分析，进而探讨人类活动对植被带的影响。

从统计的结果来看，近些年来，青藏高原地区无论是人口规模还是牲畜规模都得到了较快的增长(图6)。人口数量从1971年的813万人很快增长为2010年的1 553.5万人，平均年增长率达到18.51万人·年-1，总人口增长了近两倍。而以羊单位换算的整个青藏高原牲畜的数量则由1971年的1.002 8 ×108羊单位增长到 2010 年的 1.343 6 × 108羊单位，平均增长速度为85.2万羊单位·年-1，总规模增长了1.34倍。由于青藏高原地区以畜牧业为主的第一产业占有较大的比重，从第一产业的发展及整个经济的增长规模来看，青藏高原地区近些年的经济增长较为明显，尤其是第二、三产业增长速度一度处于加速阶段(图7)。另外近些年来，以矿产开采及能源开发为主的第二产业以及以旅游和服务业为主的第三产业发展势头较快，其比重越来越高，并逐渐占据了经济发展的主体地位。另外人口的快速发展也进一步地推动了青藏高原城镇数量和规模的发展，近些年城镇化的发展速度较快[46]。

由于青藏高原地区的植被生长季节较短，高海拔高寒及缺氧的生态环境使得生态系统较为脆弱，极易受人为活动的影响和干扰。另外与其他低海拔气候较适宜地区相比，青藏高原地区无论是从植被的生物总量、活动状况、污染的自然降解水平还是破坏后的自我恢复能力等都比其他暖湿气候条件下的生态系统差，脆弱的生态环境使得这一地区对人为因素的干扰尤其敏感[43, 60]。

图6 青藏高原地区1971-2010年人口规模及以羊单位换算的牲畜数量的增长情况Figure 6 Population and sheep-unit livestock numbers for the Qinghai-Tibetan Plateau during 1971-2010

图7 1978-2010年青海、西藏地区三次产业的规模(以当年价格计算)及城镇化率Figure 7 The sizes of the three industries and urbanization rate of the Qinghai Province and Tibet

图8 青藏高原地区最近35年的人口(A)及牲畜(B)密度线性增长率Figure 8 Mean population (A) and mean livestock increase rates (B) over the past 40 years

为了从空间上探讨人为因素对植被带的影响，本研究分析了青藏高原人口的分布及近40年的人口增长速度(图8)，以及以羊单位换算的整个青藏高原的牲畜增长速度(图8)。受高海拔及高寒缺氧等恶劣的自然条件的限制，人口基本上分布于海拔高度相对较低、植被状况较好的低山河谷地区，并在高原边缘接近低海拔的地区有相对较为集中的分布。从生态气候分区上来看，H2C1、H1B1、H2AB1及H3等地区的人口密度较高。而在高原广阔的中西部地区人口分布则较为稀疏。从人口的年平均增长率也可以看出，高原边缘地区、藏南及中部地区的人口增长较快。牲畜的空间分布与人口的分布相似，主要分布于海拔高度相对较低、气候环境相对较好的低海拔河谷地区，H2C1、H1B1、H2AB1及H2C2等地区是分布相对较为集中的地区。与人口分布不同的是，牲畜的增长速度近些年在川西北及高原的西部高海拔地区增长较快，而在高原中部的H1C1、H1B1两个地区有减少的态势。由于这两个地区是三江源自然保护区，近些年来的退牧还草等生态移民及植被保护政策可能是其主要原因。

以上无论是人口规模的分布还是牲畜规模的分布，如果对比分析最近几十年的植被带生长状况在空间上的变化(图2)，则可以发现，在H2C1、H1B1、H2AB1等人口及牲畜密度较高的地区，虽然植被的生长环境相较高原北部及西部高海拔地区为好，但与其植被改善相反，这些人口及牲畜密度较高的地区其植被却普遍出现了统计意义上的明显下降态势。由于人为活动对植被带的影响除了人为的对矿山开采、森林的砍伐及城镇、基础设施建设对土地的侵占之外，最为直接的则是以发展畜牧业为主的牲畜放牧对植被带的影响[27, 43-44]。为了检验人为因素对植被带的压力和影响，本研究按县统计了最近35年来青藏高原人口及牲畜的密度及其与NDVI变化率之间的关系。由于密度的数量级差异较大，本研究在取密度的自然对数后与NDVI增长率之间进行了Theil-Sen线性回归分析 (图 9A，9B)。

从在青藏高原上的207个县区的统计分析结果来看，人口、牲畜密度对数和NDVI变化率之间均出现了较为一致的负相关性。其中NDVI变化率NDVItrend和人口密度自然对数ln(Population Density)之 间 的 回 归 方 程 为y= -0.000 2x+ 0.000 2(R2=0.166，P＜ 0.001)，NDVI变化率与牲畜密度自然对数ln(Livestock Density)之间的回归关系为y= -0.000 3x+0.009(R2= 0.142，P＜ 0.001)，二者的回归方程检验均较为显著。回归结果表明，随着人口、牲畜密度的增加，NDVI增长率出现了明显的递减，反映出随着人口、牲畜对植被带的压力加大后，植被带出现了退化现象。这与在NDVI变化生态气候分区统计中的结果基本一致。

图9 以及1981-2015年以县为统计单位人口密度(取自然对数)与NDVI增长率之间的回归关系(A)及牲畜密度(取自然对数)与NDVI增长率之间的回归关系(B)Figure 9 The correlations between the NDVI trend and ln(population density) (A)/ln(livestock density) (B)during 1981-2015 in each county

由于人口及牲畜分布较为集中的地区为H2C1、H2AB1及H3等海拔相对较低、气候环境相对较好的地区，但这些地区的植被生长状况却普遍出现了退化。在分析人口/牲畜密度与NDVI变化率之间的关系后可知，气候的暖化使这些地区严寒的气候条件缓解，变得更适宜于人类的游牧和定居，人口的迁入和畜牧业的快速发展却对这一地区的植被带带来了压力和干扰、产生了负面的影响，干扰加大从而使植被带总体上处于退化状态。在这些地区植被带生长状况恶化的原因中，人为因素的影响已经超过了平稳而缓慢的气候变化影响本身，已成为其退化的主要原因。

4 结论

本研究以对气候变化和人为影响较为敏感的青藏高原地区植被带变化趋势为视角，采用GIMMSNDVI来反映植被带的生长状况，并就其时空间变化与气候因子、人为影响之间的响应关系进行了分析，得出以下结论：

ⅰ)最近35年来整个青藏高原的植被带生长状况有缓慢的增长趋势，这与近些年来的气候全面暖化和降水的轻微增加有正向的相关性。

ⅱ)在高原中部及西南部等高海拔地区，暖湿化的气候变化趋势是其植被好转的直接原因。即随着气候暖化，高海拔地区严寒的植物生长环境得到改善，林线上升，植被带向高海拔、高纬度地区扩展，NDVI相应地增加。

ⅲ)在高原东北部及其东部边缘、东南部等海拔高度相对较低、气候环境相对较好的地区植被普遍出现恶化。与缓慢的气候暖化相比，人口及牲畜数量的增加等人为影响增强是其主要原因。即随着气候的暖化，一些地区因严寒的气候环境改善而开始适宜于游牧及定居，人口迁入，人为活动的影响加剧，当超过了植被带的承载能力时，植被带开始退化，表现为NDVI相应地减少。

致谢：感谢两位匿名审稿人对本文提出的修改建议。中国气象局国家气象信息中心(NMIC)及英国环境数据分析中心(CEDA)提供了气象资料，在此一并致谢。