1850年以来川西高原北部植被气候生产潜力时空变化特征

王春学, 秦宁生, 周 斌, 庞轶舒, 罗 玉, 茅海祥

(1.中国气象局 成都高原气象研究所/高原与盆地暴雨旱涝灾害四川省重点实验室,成都 610072; 2.四川省气候中心, 成都 610072; 3.铜仁市气象局, 贵州 铜仁 554300)

植被气候生产潜力是指在其他条件均适宜的情况下，自然生长发育的植被由气候资源所决定的在单位时间单位面积上由光合作用产生的生物学产量或经济产量[1]。生产潜力是表征陆地生态过程的关键参数和评价陆地生态系统可持续发展的一个重要指标[2]。对生产潜力的研究不但可以揭示生产力与气候要素之间的内在联系，还可以用来预测生产力的未来发展趋势[3]。我国人多地少，资源相对贫乏，对生产潜力的研究显得更加迫切和重要，可以为合理开发利用自然资源和应对全球变化提供科学依据[4]。

Miami模型是经典的、使用最广泛的生产潜力估算模型[5]。Miami模型主要针对自然森林草原植被，利用五大洲53个地点的气象资料，获得了气象因子模拟自然植被生物量的关系表达式，之后中国学者对模型进行了验证，认为该模型适合中国自然植被生产潜力的研究应用[6]。树木年轮(树轮)数据拥有长时间序列和高空间分辨率的特点，利用树轮资料重建生产潜力有其独特的优势，能够获取长时间尺度的数据[7]。许多研究表明，树轮资料结合传统的研究方法能够较好地表征生产潜力的逐年变化，对生产潜力的估算研究起到了推动作用[8-12]。

川西高原位于青藏高原东部，地广人稀，具有丰富的植被，是生态脆弱区和气候变化敏感区，其气候及生态特征对青藏高原乃至中国西南地区都有重要影响[13-15]。陈卓奇等[16]利用遥感数据反演了青藏高原生产潜力，发现降水量小于450 mm的区域生产潜力变化的主导因子为降水量，降水量大于450 mm的区域，主导因子变为气温，并且随着气温的升高生产潜力显著提高。刘刚等[17]的研究也表明青藏高原生产潜力与温度有很好的正相关关系。朴世龙等[18]发现青藏高原生产潜力自东南向西北递减，并且1982—1999年有逐渐上升的变化趋势。李宗善等[19]重建了1788—2010年川西米亚罗地区的生产潜力，发现夏季至秋季温度是限制森林生长的最为重要的气候要素。

目前川西高原的生产潜力估算时间尺度主要集中在近几十年，更长时间尺度的估算较少，而且个别研究也只是针对单个站点展开的，非常缺乏对川西高原大范围、长时间生产潜力的科学认识。近些年川西高原树轮研究工作取得了比较丰富的成果，重建了川西高原多个地点的气候变化序列，揭示了川西高原部分地区的历史气候变化特征，同时也积累了大量的树轮资料，这为利用树轮资料建立大范围、长时间生产潜力奠定了良好的理论和资料基础[20-24]。本文将利用川西高原北部阿坝藏族羌族自治州(阿坝州)多条树轮资料重建1850年以来阿坝州生产潜力，并分析其时空变化特征，提高对川西高原历史时期生态系统的科学认识，为自然资源的合理开发利用和可持续发展提供科学支撑。

1 资料和方法

1.1 研究区概况及气象资料来源

阿坝州位于四川省西北部，地貌以高原和高山峡谷为主，岷江和大渡河贯穿全境，是黄河上游的重要水源地，野生动植物资源极其丰富。阿坝州是中国五大牧区之一的川西北牧区的重要组成部分，境内的若尔盖生态区是我国第一大高原沼泽湿地，也是气候变化的敏感区[25]。阿坝州受西风南支急流、东南季风以及西南季风等综合影响，形成了冬寒夏凉、降水适中的山地气候[24]。

本文使用的气象资料为阿坝州13个国家气象站(图1)1961—2010年逐月气温和降水量资料。分析发现阿坝州年平均气温8.5℃，1月最冷，平均气温-1.6℃；7月最热，平均气温17.3℃；年平均降水量663.5 mm，12月最少，平均2.7 mm；6月最多，平均119.5 mm。

图1 阿坝州气象站和树轮采样点分布

1.2 树轮资料来源及树轮年表建立

本文所用的树轮样本由中国科学院地理科学与资源研究所、国家气候中心和中国气象局成都高原气象研究所联合采集。采样点位于阿坝州境内(图1)，区域内森林呈带状分布于河谷两侧，受人类活动影响较少。本文所用4个采样点的优势树种均为冷杉[Abiesfabri(Mast.) Craib]，除二道海采样点(EDH)的海拔高度接近3 300 m外，其余3个采样点海拔高度均超过3 500 m，其中热务沟采样点(RWG)最高(3 680 m)。树轮采样对象为树龄较长的活树，每株至少采集2个树芯，其中除九寨沟采样点(JZG)只有11株29个样芯外，其余3个采用点的样本数均在25棵树、50个样芯以上，热基沟采样点(RJG)样芯最多(68个)，4个样本点共采集了88棵树的213条树芯，复本量能够满足树轮气候研究的要求(表1)。

表1 树轮采样点概况

将采集到的样芯按照国际通用的基本程序进行干燥(自然晾干)、固定、磨光处理和交叉定年，并利用COFECHA程序对定年结果进行检查和验证。为了去除原始树轮宽度中所包含的非气候信息，利用ARSTAN程序，采用步长为样本长度67%的样条函数拟合生长趋势，然后用双权重平均法进行合成，得到每个采样点的标准化年表(STD)。表2给出了各采样点标准年化表统计特征及公共区间分析结果，可以看到4个标准化年表的统计特征比较相似，其中平均敏感度略偏低，但是一阶自相关系数和总体代表性都比较高，具有一定的表征过去气候变化的潜力。从SSS>0.85的起始年份来看，EDH最晚(1849年)，RJG最早(1514年)，为了保证分析结果的可靠性和一致性，后文统一选取1850—2010年的树轮宽度资料进行分析。

表2 各采样点标准年表统计特征及公共区间(1800-2010年)分析结果

1.3 气候生产力的估算模型及遥感资料验证

本文采用Miami模型[3]估算阿坝州气候生产潜力，相关计算公式为：

MT=30000/[1+e(1.315-0.119T)]

(1)

式中:MT为由年平均温度决定的气候生产潜力[kg/(hm2·a)]；T为年平均气温(℃)。

MR=30000[1-e(-0.000664R)]

(2)

式中:MR为由年平均降水量决定的气候生产潜力[kg/(hm2·a)]；R为年降水量(mm)。

为了验证模型估算结果，本文使用了中国科学院地理科学与资源研究所的中国植被净初级生产力(Net Primary Productivity，NPP)遥感观测资料，提取得到了九寨沟、若尔盖、壤塘、马尔康和理县5个点的2001—2010年逐年资料。图2给出了模拟和遥感观测的对比情况，可以看到MT与遥感观测NPP的对应关系比较好，从量值上看MT为NPP的1～1.5倍，二者散点图(图2A)呈显著的正相关(R2=0.78)，从年平均变化(图2B)来看，2001—2010年二者的变化趋势和年际波动情况基本一致。MR与NPP的对应关系则不好，可以看到(图2C)二者为弱的负相关，年际变化(图2D)虽然都为增多趋势，但是具体年际波动的对应关系较差。通过与遥感资料的对比验证，发现利用温度估算的气候生产潜力比利用降水估算的气候生产潜力更据研究价值和代表性。

1.4 经验正交函数分解

气候统计诊断中应用最为普遍的方法是把原变量场分解为正交函数的线性组合，构成为数很少的互不相关的典型模态，代替原始变量场，每个典型模态都含有尽量多的原始场信息。经验正交函数(Empirical Orthography Function，EOF)分解就是这样一种方法，EOF的优点是没有固定函数，能在有限区域对不规则分布的站点进行分解，很容易将变量场的信息集中到几个模态上，而且分离出的空间结构具有一定的物理意义[26]。

1.5 MTM-SVD

MTM-SVD方法是由Mann等[27]提出的一种多变量频域分解技术。这是一种将谱分析的多锥度方法(Multi-Taper Method，MTM)和变量场的奇异值分解(Singular Value Decomposition，SVD)方法结合在一起的气候周期信号检测技术，详细内容参阅相关文献[28-30]。

图2 Miami模型估算的MT，MR与遥感NPP的散点图(A，C)和逐年对比图(B，D)

2 结果与分析

2.1 树轮和气候生产潜力相关分析

利用Miami模型分别建立1961—2010年13个气象站MT和MR，并分别统计其与树轮宽度序列的相关特征。由于树木生长情况往往与前期气候状况密切相关，所以本文计算了当年树轮序列与12个月累积生产力的相关系数，累积值从上年1月开始逐月滑动至当年1月，找出每个气象站与树轮序列相关系数的最大值，然后统计13个气象站中超过一定显著性水平的台站数，并计算其平均值。

图3A为MT与树轮序列的相关统计情况，可以看到MT与树轮宽度都为正相关，其中上年10月—当年9月的MT与树轮相关系数最大(0.57)，同时13个站均超过了0.01显著性水平，其余时段的相关系数在0.52上下波动，显著台站数在12个左右。

图3B为MR与树轮序列的相关统计情况，可以看到MR与树轮宽度都为负相关，其中上年8月—当年7月的MR与树轮相关系数最大(-0.41)，但是只有1个台站达到了0.05的显著性水平，其余时段的相关系数在-0.35上下波动，显著台站有7个左右。

由相关分析可知，树木生长越旺盛时，MT往往越大，反之则越小，即MT可以代表阿坝州的气候生产潜力。另外相关分析发现树轮与上年10月—当年9月的MT相关性最好，所以后文主要针对该时段的MT展开分析。

图3 MT和MR与树轮序列相关分析

表3为上年10月—当年9月MT与树轮序列的具体相关系数情况，可以看到13个气象站与4个树轮序列的相关系数都较高，除个别序列间的相关系数达到0.05显著性水平外，绝大多数相关系数都通过了α=0.01的显著性检验。

表3 上年10月-当年9月MT与树轮序列相关系数

2.2 建立回归方程

从表3中不难发现，除了理县站与JZG的相关系数最高外，其余站均与RJG的相关系数最高，但是如果只用RJG单一序列来重建所有站的MT，那么就无法体现MT空间变化的差异，所以为了更多地提取时空特征信息，在保证回归方程具有一定显著性水平的前提下，尽量多地选取树轮序列参与重建。经过反复试验，发现选择相关系数最高的2个树轮序列参与回归方程重建时，可以满足上述要求，即MT=B0+B1X1+B2X2,式中：B0，B1，B2为回归系数；X1，X2分别为选择的两个树轮年表指数。

为了检验重建方程是否稳定可靠，采用国际惯用的严格显著性检验方法(逐一剔除法，leave one out)对方程进行检验，表4给出了各重建方程的回归系数及相应检验统计量，可以看到所有方程的相关系数、F检验值、乘积平均数(t)和误差缩减值(RE)均达到了0.01的显著性水平，说明重建方程具有一定的稳定性和可靠性。符号检验有3个达到0.01显著性水平，5个达到0.05显著性水平，5个接近0.05的显著性水平，但同时一阶差符号检验都没有通过0.05的显著性检验，表明重建方程对低频变化的重建要好于高频变化的重建。

表4 回归方程系数及逐一剔除检验统计量

2.3 重建结果与分析

2.3.1 EOF模态分析 图4A给出了1850—2010年阿坝州MT平均值分布图，可以看到阿坝州南部的小金—若尔盖—汶川一线MT超过了15 000 kg/(hm2·a)，九寨沟—松潘—马尔康一线为11 000～14 000 kg/(hm2·a)，西北部若尔盖—阿坝—壤塘一线只有7 000～10 000 kg/(hm2·a)，即阿坝州MT分布不均，由东南向西北递减，南部大值区MT可达北部低值区的1倍以上。从1850—2010年阿坝州年平均MT变化来看(图4B)，1850—1950年MT主要表现为年际波动，没有明显的变化趋势，1950—1980年出现波动减少的变化趋势，1980—2010年出现显著增多的变化趋势。

对重建的阿坝州13个站1850—2010年MT进行EOF分析，第一模态通过了显著性检验，解释方差贡献率达92%，从其空间型(图5)可以看到整个区域都为正值，即全区一致变化型，其中西部和北部为大值带，对应时间系数与其年平均曲线变化基本一致。进一步利用滑动t检验方法对其时间系数进行突变检验，图6给出了20 a滑动的统计量曲线，可以看到1850年以来统计量有2处超过0.01的显著性，一处为正值(1967年前后)，另一处为负值(1990年前后)，即过去161 a中阿坝州MT出现过两次突变，20世纪60—80年代经历了一次显著的由多到少的转变，20世纪90年代开始又经历了一次明显的由少到多的变化。

图4 1850-2010年MT平均值分布及其逐年变化曲线

图5 1850-2010年阿坝州MT的EOF分析第一模态空间型及其对应的时间系数

注：虚线为0.01显著性水平。

2.3.2 周期分析 MTM-SVD方法可以对场空间进行周期分析，图7给出了1850—2010年阿坝州MT场的LFV谱分析结果，可以看到年际尺度上的2，2.6 a周期通过了0.01的显著性检验，4～5 a周期通过了0.05的显著性检验，年代际尺度在12，17，80 a都出现峰值，但是没有达到一定的显著性水平。从对应的10 a滑动窗口分析可以发现(图8)，准2 a周期有3个显著时段，分别是1860—1880年、1910—1940年和1970—1990年；4～5 a周期也有3个显著时段，主要出现在1890—1940年、1950—1980年和1990—2010年。

注：虚线为蒙特卡洛置信度。

注：阴影表示通过置信度检验。

3 结论与讨论

通过与遥感资料的对比，发现在川西高原地区利用Miami模型估算的气候生产潜力具有可靠的研究价值和代表性，尤其是与温度的关系要好于降水，说明川西高原植被气候生产潜力对温度的变化更为敏感。同时发现该区域上年10月—当年9月的MT与树轮宽度指数有很好的对应关系，建立的回归方程通过了显著性检验，重建的MT稳定可靠。利用树轮宽度资料重建川西高原历史时期气候生产潜力是一种弥补观测资料年限短的有效方法，有助于增强对高原地区植被生态系统的科学认识，对自然资源的开发利用有重要参考价值，而川西高原森林资源丰富，在获得更多树轮样本资料的情况下，可以重建时间尺度更长、时空精度更高的植被气候生产潜力场。

分析表明，阿坝州气候生产潜力分布不均，由东南向西北递减，1850年以来总体具有明显的年际和年代际变化特征，尤其20世纪90年代以来，气候生产潜力出现快速增大的变化趋势，虽然1999年以后出现了全球变暖停滞现象，但是川西高原MT仍然保持着较快的增大速度，这与高寒地区对全球变暖的响应较低海拔地区偏慢是一致的[31]。另外1850年以来阿坝州MT出现过两次突变，20世纪60—80年代经历了一次显著的由多到少的转变，20世纪90年代开始又经历了一次明显的由少到多的变化，随着全球变暖减缓21世纪是否会再次出现较大的突变也值得密切关注。本研究增强了对阿坝州MT空间分布、变化趋势、周期和突变等规律的科学认识，这对高原农业、牧业、林业的开发利用和发展规划都有重要意义。

由于阿坝州历史时期气候生产潜力的相关研究比较少，而且缺乏相关的公开数据，所以本文没有利用其他资料对建立的1850—2010年气候生产潜力序列进行对比验证，后期还需要不断搜集相关研究成果，开展对比验证分析。另外本文只是探索了利用树轮资料建立阿坝州气候生产潜力的可能性，并对其分布特征和变化规律等事实进行了分析，没有探讨其变化机理，尤其是为什么气温决定的气候生产潜力更具代表性等问题，值得进一步深入研究。

Miami模型在我国被广泛应用于气候生产潜力的估算，尤其是在西北干旱区[32-33]、内蒙古草原[34-35]以及我国南方地区[36]已有大量的研究成果，但是在川西高原的应用还比较有限，而本研究则证明了其温度模型在高原地区的适用性。利用Miami模型可以方便快捷的给出植被气候生产潜力的时空动态演变，但是该模型受限于气象站点的分布，而本文使用了13个国家气象站，站点分布比较均匀，但是站点密度还是相对稀疏，如果能够利用更多的气象加密站点资料，无疑会获得更加准确的结果。