2015, Vol.35
1 引言
高分辨率的树轮气候代用资料以其对气候变化响应敏感、 精度高、 连续性强、 分布广泛及样本易得等优势，在古气候研究中正被广泛采用[1, 2, 3, 4, 5, 6, 7, 8]。树轮宽度的年际间变化不仅可以重建环境变化的历史[9]，还可以反映气候变化对生态系统的影响。但年轮宽度与气候因子之间有着复杂的关系，这种关系受气候因子间的相互制衡和因物种而异的树木生长节律的共同影响。在利用树木年轮开展历史时期气候变化的研究中，剔除树木年轮与年龄相关的生长趋势是准确获取气候信号的先决条件，又被称为年轮宽度序列标准化[10, 11]，其中以线性或负指数函数和样条函数法应用最为广泛，当前我国的树轮气候学者多采用这些传统的方法[12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25]。
然而，随着树轮气候学研究的需求和发展，已有的方法不断受到质疑和挑战。譬如负指数和样条函数的结合弥补了单一线性函数对生长趋势拟合的偏离，但样条函数步长的选取较依赖于研究者的经验，存在一定的随机性和主观性，同时会剔除较多的低频气候信息。Melvin和Briffa[26]建立了“signal-free”方法用于解决10年至百年尺度的异常变化对传统的树木生长拟合曲线的影响，但也同样存在低频变化的缺失以及使用不同拟合方法导致的不确定性问题[27]； 徐岩和邵雪梅[28]提出总体生长趋势曲线拟合法，运用树木径向生长模型及大样本量树木生长总体趋势的估计估算区域的总体生长曲线，仍然要求能够在同一气候区域获得同一树种且足够多的包含髓心的样本的情况下才可使用； Briffa等[29, 30]和Esper等[31]开展了关于树轮低频信息提取方面的研究，提出了区域曲线标准化(Regional Curve Standardization，简称RCS)方法，较好的保留了气候低频信息，在Mann等[32]利用树轮重建的北半球过去近千年来的温度距平序列的工作中得到了应用和检验，保留了百年以上尺度的气温变化信息，中世纪暖期的争议得到了解释。但Yang等[33, 34, 35, 36]通过对都兰地区长达2000年的祁连圆柏树轮数据的分析进一步指出了RCS方法应用中的缺陷，认为即使同一地区的同一树种在同一生命周期生长速率仍然有较大差异，传统的RCS方法会高估(或低估)与年龄相关的生长速率，源于对RCS方法改进的思路提出了特征值分析，对选定的祁连圆柏的特征向量进行归类，分析该树种与年龄相关的生长速率的共性及差异，尝试建立一套可以区分年龄相关的信号及气候相关的信号的方法； 并在此研究的基础上，结合Signal-Free Regional Curve Standardization(SF-RCS)方法，去除生长趋势，重建了青藏高原东北部过去3500年的降水变化[7]。
随着树种和环境的不同，越来越多的研究者普遍认为确定出每一棵树木的生理年龄及准确的生长量订正非常困难，靠一种方法或一种模式解决所有树种、 各种生态环境中树木生长量订正的问题不太可能。经验模态分解(Empirical Mode Decomposition，简称EMD)是1998年Huang等[37]提出的一种有效处理时间序列信号的新方法，能够比较有效地对信号进行多分辨率分析； 从本质上讲是将一个信号进行平稳化处理，按照自然变化把不同特征尺度(频率)或层次的波动或趋势逐级从原信号中分解出来，得到一系列具有不同特征尺度的本征模函数(Intrinsic Mode Function，简称IMF)分量，研究者可根据序列的自然特征，分析不同特征尺度序列表征的物理意义。该方法被认为是近年来线性和稳态谱分析的一个重大突破。目前EMD方法已成功应用于信号处理、 图像处理、 湍流、 地震及大气科学等非线性科学领域[38, 39, 40, 41]。在树轮学研究中，Fang等[42]利用EMD对所有树轮序列进行分解，分析可靠的高频和低频信号分别重建了高频和低频年表。本文进一步运用EMD方法对选用的祁连山青海云杉样芯宽度序列进行逐一分解，分析所得各序列的趋势项对该树种年龄相关的生长趋势的代表性，通过与传统去趋势方法的比较，综合分析该方法在去趋势中的适用性和可能存在的问题，尝试建立一套可参考的树轮宽度资料生长量订正的方法，供读者多一个新的选择。
2 资料与方法 2.1 研究区概况及数据处理祁连山位于我国西北干旱区、 东部湿润区和青藏高原区的过渡地带，受大陆性气候和高原气候的综合影响，属温带半干旱区。作为重要的地理气候分界线，在我国自然地理上占有非常重要的地位。本研究所用青海云杉(Picea crassifolia Kom.)树轮资料采自甘肃肃南县青龙乡孔岗木(Kgm)林场，位于祁连山中部地区，标号为Kgm03和Kgm04。采样时多选择土壤贫瘠、 岩石裸露，生长条件相对较差的小生境。为保证复本量，一般要求一个点能采集到15株以上，树上取两根样芯，以便于更准确的进行实验室交叉定年工作，表 1为具体的采样信息表。因Kgm04样点样芯株数低于15，又与Kgm03样点较为接近，故后期分析中将两个样点合并为一个样点。
|
表 1 孔岗木样点基本信息 Table 1 Information of Kgm sampling site |
所有样本在实验室晾干、 固定、 打磨等预处理之后，进行交叉定年。以0.01mm精度的轮宽测量仪进行轮宽量测，之后用COFECHA程序进行检验和校准，以保证所有样本定年和量测的准确无误。
2.2 青海云杉生长量订正上述方法得到的逐年树木生长量受生物遗传特性影响遵循一定的规律，已有研究表明[43, 44]，单株云杉林木初期生长缓慢，20年后生长迅速增加，径向生长的峰值一般在50年左右，当到达树木生长最旺盛的时期时，平均的年轮宽度常为整个生长期内的极大值，后逐步减小并呈规律性下降趋势，最后达到一个相对稳定的水平。但气候环境提供并限制了树木生长的大背景及长尺度上的轮宽交替变化； 且随着立地微环境及郁闭度的变化，在同一生命周期年龄相关的生长速率可能存在较大的分异[45]，与年龄相关的生长趋势也会呈现多种变化形式及微差异。
下文将采用经验模态分解(EMD)对每个云杉样芯的宽度测量序列进行分解及趋势项分析。EMD分解过程如下: 利用波动上、 下包络的平均值去确定“瞬时平衡位置”，进而分解得到IMF分量，每次提取都将高频信号分离，剩下较低频信号，最后得到趋势分量或者定常值，最终把原序列分解成有限个IMF分量，第一个IMF分量代表原序列中最高频的成分。将原始序列减去第1个IMF分量依次进行上述平稳化过程就可以得到第二个IMF分量，直到不能分解为止，最后一个分量代表原序列的均值或趋势，不同层次的IMF分量都可能对应某一物理背景，而原始序列就可以由这些IMF分量均值之和表示。
研究中首先利用EMD运算程序对祁连山孔岗木(Kgm03和Kgm04)两个云杉样点58个样芯(见表 1)的宽度测量序列进行逐一分析。其中本研究中最长的样芯Kgm0317a(1723~2000年)和Kgm0317c(1790~2000年)分解所得分量见下图(图 1和2)，由图 1和2可以看出，EMD分解将两个样芯序列由高频至低频分别分解得到了6个本征模态分量，每个IMF分量对应于一个窄波段信号，各分量分别与树木生长的限制性因子(生物效应、 气候效应和微环境)相联系，反映了树轮序列记录的不同时间尺度信号的自然波动。其中IMF1分量总是代表序列最高频的变化成分(树龄为200年左右的云杉径向生长的年际变化)，与年际的气候变化信号较为接近，可能主要反映了树木生长与气候相关的年际变化特征， 是树轮气候研究中最可靠的信号； IMF 2至IMF5依次显示了两个样芯年代际尺度至多年代际尺度的周期性振荡及振幅随时间的变化过程； 每个IMF分量的变化规律符合自然信号的非线性变化特点，波动具有相对稳定的准周期，但各尺度信号对原始数据总体特征的影响程度有很大差异，各分量振幅的平方可反映该分量在原始数据中的信号强度和能量大小，与原序列的相关系数反映各分量对原序列变化特征的影响程度； 而IMF6为两个样芯的最低频趋势项(图 1f和图 2f蓝线)，表示的是信号中存在的线性项或缓慢变化的周期大于信号资料长度的成分，与传统理论研究中认为的年龄相关的生长速率变化曲线较为接近，是否主要反映了树木生长的生物效应，下文将通过与传统拟合趋势线进行对比，分析EMD趋势项在祁连山青海云杉宽度序列生长量订正中的适用性； 以上各分量之和可还原为原始变量序列(图 1g，1h和图 2g，2h黑线)。
|
图 1 Kgm0317a宽度序列EMD分解各分量 Fig. 1 IMF components of the Kgm0317a based on EMD |
|
图 2 Kgm0317c宽度序列EMD分解各分量 Fig. 2 IMF components of the Kgm0317c based on EMD |
由IMF6(图 1f蓝线)可以看出，本研究中主要分析的趋势项变化特征如下:1723~1760年生长速率较低(一般云杉前20年生长缓慢)，1760~1800年缓慢增加(原则上该阶段呈指数式增长，但基于取样方向、 郁闭度等微环境影响下呈现出的生长速率变化表现了一定的独特性)，1800~1890年逐渐加快并达到生长的峰值，1900年后增长速率减缓，1900~2000年缓慢减速至平稳状态(慢-快-慢，渐趋于平稳)，该曲线基本与树木自幼龄期开始的生长波动较为接近，在幼龄期的低生长大多认为可能和周围大树遮挡有一定关系，1800~1890年可能存在一次采伐后的生长释放[45]，后期的生长则与树木的生物学遗传效应较为吻合。一般的线性函数拟合的生长趋势曲线(图 1g红线)在整个生长过程中呈持续上升趋势，可能对树木幼龄期生长速率微小变化信息的拟合有一定偏差，通常通过变换样条函数进行补充以达到更好的效果； 进一步对比了样条函数的拟合曲线(图 1h绿线)，步长选取序列长度的2/3，结果显示，样条函数的拟合曲线与EMD分解的趋势项较为接近。但样条函数的步长常常依赖于研究者的经验和对研究样芯的情况分析而定，具有一定的主观性，拟合结果随研究者选取步长的改变得到的曲线差异较大，因此，结果有一定的不确定性。
采自同一株树的两个样芯往往钻取方向不同，获得的长度不同，其包含的年轮数量也不同，图 2所示的Kgm0317c样芯展示了与Kgm0317a不同的情况。该序列起始年为1790年，而图 1所示的Kgm0317a起始年为1723年，表明反方向获取的样本远未达到髓心，即树木幼龄期的生长轮在采样中由于采集困难而缺失，因此，该序列的趋势项不能完全显示树木自幼龄期的生物趋势，参照图 1，Kgm0317a的趋势项自1790年开始由缓慢生长转为生长增速，1900年后生长速率减缓且渐趋于平稳。与图 1不同的是，获得的Kgm0317c样芯可能丢失了幼龄期的低生长阶段年轮，其趋势项的变化过程基本分为3个阶段，第一阶段1790~1900，生长速率缓慢增加，与Kgm0317a的同阶段生长速率变化特征较为一致； 第二阶段1900~1950趋于平稳，第三阶段1950年之后趋于水平； 同时对该样芯进行了线性函数和样条函数拟合，所得到的线性拟合趋势线(图 2g红线)整个时段接近于水平，与EMD得到的趋势线在第一阶段(1790~1900年)呈现出一定差异； 而样条函数拟合趋势线(图 2h绿线)仍然选取 2/3 步长，得到的拟合线在第一阶段与EMD趋势线大致接近，而在后两个阶段与前两者均有不同，与EMD趋势线的主要差别表现在第三阶段(1950~2000年)，样条曲线呈现为上升趋势，与实际需要拟合的树木随年龄变化的生长效应理论相悖，可能会产生虚假的过度拟合，导致原本树轮资料记录的上升的气候信号丢失。结合与Kgm0317a得到的一些共性，综合分析3种拟合方法得到的生长曲线，我们初步认为EMD的趋势项对Kgm0317c的年龄相关的生长趋势的表达可能更接近于真实情况。
基于上文的结果，我们进一步利用EMD方法对本研究中Kgm03和Kgm04样点的所有样芯(见表 1)宽度序列进行了逐一分析，得到了各序列的时间变化分量和最低频趋势项，初步检验该方法对青海云杉树轮宽度资料年龄相关的生长趋势提取的代表性、 适用性和可能限制因素。通过对所有样芯的对比归类，可以发现两种方法对绝大部分样芯的去趋势维持了相似的结果，但对个别样芯在一些时段的差异也值得重视。我们挑选出传统拟合曲线和EMD趋势线及相应得到的树轮生长量订正后的指数序列有较大差异的4种情况进行了较为详细的讨论(图 3)，分析其可能的原因和两种方法的利弊，主要如下: 1)生长效应低估(髓心端放大效应): 样芯Kgm0316b(图 3a黑线)序列接近髓心部分生长速率较低(幼龄期可能受到环境抑制)，这类样芯在生长条件完全满足的阶段常出现突然的生长释放，这种快速的转折在线性拟合时因为函数的斜率较大，对髓心端的生长效应低估而导致指数序列在髓心端异常放大，产生奇异值(图 3c红线)； EMD因为自适应的特性，在这种情况有一定优势，得到的序列在总体方差上更加平稳。2)生长效应高估(树皮端过度拟合): Kgm0317a为采到的最长样芯，1723~2000年树木由幼龄期至成熟期经历了先加速至顶峰，之后下降渐趋平稳，EMD分解的趋势线基本反映了上述变化(图 3b蓝线)； 而传统线性拟合得到的趋势线(图 3b红线)表现为从拟合期开始向一个方向持续变快或变慢的特点，对该序列的拟合在1775~1810年、 1850~1900年和1975~2000年3个时段均与EMD趋势线呈现出一定差异。传统方法得到的指数值(图 3d红线)在1850~1900年时段略高于EMD方法得到的指数值(图 3d蓝线)，而1775~1810年和1975~2000年时段低于EMD得到的指数值，初步分析认为第三时段1975~2000年(树皮端)基本为树木年龄相关的生长速率进入衰减阶段，较少可能仍呈现为持续上升的状态，因此EMD方法实现的生长量订正可能更可靠，更接近于真实情况(图 3d蓝线)，线性拟合得到的指数序列引起的树皮端的降低可能是生长效应被高估的一个表现(图 3d红线)。早期和中期两个时段的差异所反映的两种方法对生长量订正的优劣，需要与更多长序列样芯对比才能客观评价。3)线性或负指数拟合失效: 树木生长受到严重抑制时线性拟合可能完全失效，拟合函数系数出现负值，被迫选用非线性拟合(图 3e红线)，非线性拟合至目前为止采用较少； 作为对传统线性拟合方法的补充，EMD也是一个新的选择(图 3e蓝线)。两种方法得到的指数值主要在1910~1925年时段有一定差异，非线性拟合的指数值此时段为整个序列的峰值区，高于EMD得到的指数值，其他时段基本一致。而1910~1925年时段的峰值很难找到相对应的该时段树木生长受气候影响快速生长的研究结论，因此，分析认为可能因为非线性曲线对树龄决定的生长速率的急剧变化偏离产生了虚假高值。4)线性拟合方向单一: 线性拟合由于从拟合期开始主要向一个方向持续变快或变慢的特点，如对Kgm0411a(图 3f黑线)的拟合呈持续下降(图 3f红线)，受斜率影响导致末端上翘，出现虚假极大值(图 3h红线)； EMD趋势项末端趋于平缓(图 3f蓝线)，得到的树轮指数(图 3h蓝线)避免了最后5~10年的异常上升。
|
图 3 EMD方法与传统方法在拟合生长趋势时的比较原始宽度(a，b，e和f黑线)和线性拟合(a，b，e和f红线)及EMD分解趋势线(a，b，e和f蓝线)；传统拟合的标准化序列(c，d，g和h红线)及EMD方法得到的标准化序列(c，d，g和h蓝线) Fig. 3 Comparison of fitting growth trends between EMD and traditional methods. The original ring-width measurement series(black of a，b，e，and f)； the growth trends fitted by traditional(red of a，b，e，and f)； the growth trends fitted by EMD(blue of a，b，e，and f). The de-trended ring-width index series based on the traditional method(red of c，d，g，and h) and the de-trended ring-width index series based on the EMD method(blue of c，d，g，and h) |
将祁连山孔岗木Kgm03和Kgm04两样点树轮资料合并(58芯)，从58芯中剔除部分较短样本，选用质量较好的44根样芯宽度序列，参照常规的树轮曲线标准化订正思路(原始曲线除以趋势项)用EMD方法计算了树轮宽度年表(图 4a)，最终年表序列长度为1723~2000年，并与传统方法计算的树轮宽度年表(图 4b)进行了对比和偏差分析(图 4c)。图 4中竖直红线表示样本信号强度EPS大于0.85时对应的年份为1845年。年表在EPS<0.85部分由于样本量较少，故仅用作参考，最后确定可靠年表为1845~2000年。
|
图 4 不同方法得到的树轮宽度指数序列 细线为年值，平滑曲线为5年滑动平均值 (a)EMD方法得到的序列(the series developed using EMD method)；(b)传统方法得到的序列(the series developed using traditional method)；(c)两条序列之差(the difference between two series，a-b) Fig. 4 Ring-width indices derived by two methods. The thin line represents annual values and the smooth line represents the 5-year running means |
由EMD方法(图 4a)和传统方法(图 4b)对比结果可以看出，用不同方法得到的序列在低频和高频变化上大部分时段表现出了较好的一致性。两条序列一致的低值期主要有1850~1880年和1925~1940年； 一致的高值期主要有1825~1850年、 1880~1925年和1940~1960年。两种方法对幼龄期树木生长趋势的拟合常有一定差别，因此得到的年表在早期表现了一定的差异，1723~1750年传统方法得到的年表值略高于EMD方法得到的年表值，但该时段仅有一芯，代表性较差； 而1750~1845年间传统年表值低于EMD年表值，传统年表在该时段呈现出较大的方差，但因该时段EPS<0.85，样芯数仍较少，传统年表所指示的低值是否反映了真实的气候信息，而非生长趋势拟合效果引起的虚假信息尚具有一定的不确定性，该差异所反映的两种方法的优劣也有待进一步研究确定。
图 4c更直观的显示了两条序列的绝对差值，仅1750~1840年差异较大，1850~1960年两者的差值序列接近零值，说明该时段两种方法对树轮生长量的订正极为相似； 而在树皮端(1960~2000年)两者的差值略微升高，说明两者对树皮端的拟合略有不同。综上所述，两种方法对生长量的订正差异主要体现在前后两端，因此，应用中可重点分析拟合线对前后两端生长效应的拟合效果而决定最后需用的去趋势方法。
3.2 年表与气象资料的相关分析为了进一步显示EMD方法建立年表的可靠性，我们进一步分别计算了EMD方法建立年表和传统方法建立年表与样点最近点甘肃肃南气象站(38°50′N， 99°37′E； 海拔2311m； 1957~2000年)器测记录的月平均气温、 月降水量之间的相关系数(r)，结果见表 2和表 3(表中黑体表示r通过了99％信度检验)。
|
|
表 2 EMD年表及传统年表与肃南气象站月均温的相关结果对比 Table 2 Correlation comparison of monthly mean temperature of Sunan Meteorological Station with EMD index and traditional index |
|
|
表 3 EMD年表及传统年表与肃南气象站月降水的相关结果对比 Table 3 Correlation comparison of monthly precipitation of Sunan Meteorological Station with EMD index and traditional index |
表 2和表 3表明，两种方式建立的年表与最近点气象站气温和降水的相关时段及相关系数极为一致，均与前年7月、 8月和当年6月气温显著负相关； 与前年8~9月降水量呈显著正相关，随着气温变化，10~12月树木基本停止生长，过量的降水将不再与年轮序列呈明显正相关变化，相关系数由正转负； 当年1月的降水可能将以冰雪的形式存储于地表及土壤层，4月初随着气温达到较为适宜的状态，树木年轮快速增长，年轮指数与当年4~6月降水量呈较为连续的正相关。与前年8月至当年7月整年的降水量相关系数显著提高，通过了99％的信度检验。两个指数序列与气候因子的相关结果完全一致，说明两种方法得到的树轮指数在高频变化上的高度一致性，这与实际情况比较吻合，表明EMD方法适用于树轮宽度序列的生长量订正，可以较好的去除生长趋势，尤其是对于传统线性函数在特殊样芯前后几十年拟合效果较差时，EMD方法提供了一个较好的补充。
综上两个年表与肃南气象站气温和降水的分析结果表明，降水是该区云杉生长的主要气候限制性因子，云杉宽度指数序列主要反映了区域降水量的历史变化信息。下文进一步利用EMD对各样芯分解所得的不同尺度的IMF原始分量的聚类集合序列分析各分量所反映的对应尺度的区域降水变化特征。
3.3 祁连山云杉序列多尺度频谱分析选用功率谱对EMD方法和传统去趋势方法建立的两种指数年表的周期变化进行对比分析，结果如下图 5所示， 两条序列的频谱分布特征较为一致，其中高频2~3年的变化周期，可能主要受QBO(准2年震荡)以及QTO(准3年震荡)影响所致，即受赤道地区平流层的26~30个月为周期的东向风和西向风之间的不断转换形成的准2年及准3年震荡[46, 47]； 5~8年的周期可能主要与ENSO海气耦合周期相联系[48]； 另外还有33年左右的共同周期。不同之处是，EMD方法获取的年表显示了通过信度检验的16年周期，该信号与太阳活动的11年准周期[49]较为接近，可能主要反映了太阳活动对该区域年降水量的影响。
|
图 5 EMD方法(黑线)和传统方法(灰色线)计算出的 树轮指数序列的功率谱对比 Fig. 5 The spectral analysis of EMD index (black line) and Traditional index(gray line) |
将EMD方法对从两个样点中选出的质量较好的44根样芯分解所得的所有IMF分量分别进行功率谱分析，并按小于10年，10~30年，30~60年及大于60年尺度分别进行分类集合运算(即高频IMF分量平均，低频IMF分量平均)，获取高频至低频不同尺度的原始测量序列的IMF分量的分频集合序列，考虑到可靠序列长度仅为150年左右，故大于60年尺度的集合序列未做讨论(周期重复较少)，图 6a，6b和6c分别展示了44根样芯前3个尺度范围上原始IMF分量的集成序列。
|
图 6 EMD分解得到的各尺度分量集合序列 Fig. 6 The integrated components of ring-width series in different time-scales |
图 6a为小于10年的年代际变化集合序列，这是树轮气候研究中被认为最为可靠的信息，将去掉了低频趋势的该序列与肃南气象站的降水进行了进一步的相关分析，结果如表 4(黑体表示相关系数通过了99％信度检验)。提取的高频集合序列与降水量的相关时段与上文两个指数年表与气候因子的相关时段完全一致，但相关系数显著提高，尤其是与前年8月至当年7月，相关系数提高到0.734。而高频上两序列的相关更能客观的反映树轮序列与气候因子相关意义的真实性，避免了因低频趋势相关的假性相关。也为进一步深入研究树木年轮与气候在不同尺度的相关关系提供了一个思路。
|
|
表 4 高频集合序列1和肃南气象站月降水及月气温的相关结果 Table 4 Correlation coefficients of the first integrated component with monthly precipitation and monthly mean temperature of Sunan Meteorological Station |
集合序列1如图 6a，较好的显示了1845~2000年间降水的高频变化特征: 最大的波动主要发生在20世纪10~20年代(1907~1920年)，不仅振幅超过5个(标准差=0.161)且变率较快，是百年内一次较大的气候波动； 1980~2000年振幅显著减小。
EMD指数年表的功率谱显示16年的周期通过了信度检验，而传统方法建立的指数年表由于低频信息保留较少，16年的周期未达到信度检验水平(见图 5)。集合序列2为所有样本10~30年的年代际尺度变化分量(图 6b)，该序列与EMD树轮宽度年表 16年左右通过信度检验的主周期较为接近。1925~1933年的树木低生长与甘肃境内较大范围的干旱背景紧密相关，与文献记录的该时期的极端干旱事件极为吻合[3, 4, 5]。
30~60年尺度的多年代际集合序列3如图 6c，1845~1880年期间振幅变化较小，随后振幅快速变化，高低生长期交替出现: 高值期有1880~1924年、 1937~1957年和1974~1993年； 低值期有1924~1937年和1957~1974年及20世纪末。
4 结论本文基于EMD方法对青海云杉树轮宽度生长量的订正研究，分析了EMD方法在分离树木生长趋势及多尺度气候信息提取中的可行性及相对优势。初步分析表明，EMD分解得到的趋势项可以较好的反映祁连山青海云杉与年龄相关的生长趋势，适用于云杉生长量的订正，可以获取能较好反映研究区气候变化的可靠树轮年表。但该方法仍然过分依赖数学计算，缺乏更多生物学的耦合分析，且对不同树种及较长树龄的树轮样芯生长量订正的效果仍有待探讨，相对于传统方法仍然得不到绝对的优势和提高，有待进一步的研究。但是本研究旨在提供一种可选择的新的去趋势思路，为读者在研究中提供可靠的选择依据。初步得出以下结果:
(1)基于EMD方法的生长量订正与传统方法的差异大致可归为4种情况:EMD方法对于传统线性或负指数函数可能引起的样芯髓心端放大效应在一定程度上有所减缓； 长树龄样芯树皮端生理效应减弱，而线性函数由于斜率的持续性会导致末端生长效应被高估(或被低估)年表值偏低(或偏高)，EMD方法在这方面有一定改善； 传统线性及负指数函数拟合失效时，本文介绍的EMD方法的生长量订正可作为一个较为适用的补充选项； 同时EMD方法还可以较好的避免传统函数因方向单一导致的生长效应高估或低估出现奇异值。
(2)器测记录期内EMD方法和传统去趋势得到的两个序列与气候因子的相关分析表明，两者对气候因子的响应结果极为一致，表明EMD方法可以得到能反映气候信息的可靠年表； EMD分解得到的所有原始序列年代际尺度的高频集合分量与气候因子的相关系数显著提高，对分析树轮序列和气候因子之间不同尺度的响应关系提供了一个新思路。
(3)EMD年表的主周期对应的原始序列的3个集合分量分别显示了1845~2000年间区域降水的不同尺度变化特征。小于10年尺度的降水高频序列波动期主要发生在20世纪10~20年代； 10~30年尺度上1925~1933年的显著低值期与甘肃境内较大范围的干旱紧密相关，与高频序列一致的显示出20世纪末降水振幅趋于平稳的特征； 30~60年尺度的多年代际集合分量较好的显示了1845~2000年间对应的3次干湿交替变化特征。
| 1 | Cook E R, Anchukaitis K J, Buckley B M et al. Asian monsoon failure and mega-drought during the last millennium. Science, 2010,328 (5977):486~489 |
| 2 | Jones P D, Briffa K R, Osborn T J et al. High-resolution palaeoclimatology of the last millennium:A review of current status and future prospects. The Holocene, 2009,19 (1):3~49 |
| 3 | Liu Xiaohong, Qin Dahe, Shao Xuemei et al. Temperature variations recovered from tree-rings in the middle Qilian Mountain over the last millennium. Science in China(Series D), 2005,48 (4):521~529 |
| 4 | Zhang Y, Tian Q H, Gou X H et al. Annual precipitation reconstruction since A.D. 775 based on tree rings from the Qilian Mountains, Northwestern China. International Journal of Climatology, 2010,,31 (3):371~381 |
| 5 | Tian Q H, Gou X H, Zhang Y et al. May-June Temperature reconstruction over the past 300 years based on tree rings on the Qilian Mountains of northeastern Tibetan. IAWA Journal, 2009,30 (4):421~434 |
| 6 | Tian Qinhua, Zhou Xiuji, Gou Xiaohua et al. Analysis of reconstructed annual precipitation from tree-rings for the past 500 years in the middle Qilian Mountain. Science China:Earth Sciences, 2012, doi:10.1007/s11430~012-4375~6 |
| 7 | Yang B, Qin C, Wang J L et al. A3500-year tree-ring record of annual precipitation on the northeastern Tibetan Plateau. Proceedings of the National Academy of Sciences of the United States of America, 2014,111 (8):2903~2908 |
| 8 | Yang B, Qin C, Huang K et al. Spatial and temporal patterns of variations in tree growth over the northeastern Tibetan Plateau during the period AD 1450~2001. The Holocene, 2010,20 (8):1235~1245 |
| 9 | Fritts H C. Tree Rings and Climate. London:Academic Press, 1976.246~311 |
| 10 | Cook E R, Briffa K R, Shiyatov S et al. Tree-ring standardization and growth-trend estimation. In:Cook E R, Kairiukstis L A eds. Methods of Dendrochronology. Dordrecht/Boston/London:Kluwer Academic Publishers, 1990. 104~123 |
| 11 | Cook E R. A time series analysis approach to tree-ring standardization. Tucson:The Dissertation for the Doctoral Degree of University of Arizona, 1985.1~171 |
| 12 | Liu Y, Linderholm H W, Song Huiming et al. Temperature variations recorded in Pinus tabulaeformistree tings from the southern and northern slopes of the central Qinling Mountains, Central China. Boreas, 2008, doi:10.1111/j.1502~3885.2008.00065.x |
| 13 | 刘 禹, 蔡秋芳, 宋慧明. 关于青藏高原2485年温度的季节和空间代表性问题. 第四纪研究, 2013,33 (1):108~114 Liu Yu, Cai Qiufang, Song Huiming. Seasonal and spatial representativeness in the northern Tibetan Plateau. Quaternary Sciences, 2013,33 (1):108~114 |
| 14 | Shao X M, Xu Y, Yin Z Y et al. Climatic implications of a 3585-year tree-ring width chronology from the northeastern Qinghai-Tibetan Plateau. Quaternary Science Reviews, 2010,29 (17~18):2111~2122 |
| 15 | 邵雪梅, 王树芝, 徐 岩等. 柴达木盆地东北部3500年树轮定年年表的初步建立. 第四纪研究, 2007,27 (4):477~485 Shao Xuemei, Wang Shuzhi, Xu Yan et al. A 3500-year master tree-ring dating chronology from the northeastern part of the Qaidam Basin. Quaternary Sciences, 2007,27 (4):477~485 |
| 16 | Gou Xiaohua, Deng Yang, Chen Fahu et al. Tree ring based streamflow reconstruction for the upper Yellow River over the past 1234 years. Chinese Science Bulletin, 2010,55 (36):4179~4186 |
| 17 | Fan Zexin, Bräuning Achim, Tian Qinhua et al. Tree ring recorded May-August temperature variations since A.D. 1585 in the Gaoligong Mountains, southeastern Tibetan Plateau. Palaeogeography, Palaeoclimatology, Palaeoecology, 2010,296 (1~2):94~102 |
| 18 | 刘昶智, 勾晓华, 方克艳等. 甘肃南部公元1824年以来降水重建. 第四纪研究, 2013,33 (3):518~525 Liu Changzhi Gou Xiaohua, Fang Keyan et al. Precipitation reconstruction in southern Gansu Province since A.D. 1824. Quaternary Sciences, 2013,33 (3):518~525 |
| 19 | 蔡秋芳, 刘 禹. 山西五鹿山油松树轮宽度年表的建立及过去百余年5~6月平均气温变化. 第四纪研究, 2013,33 (3):511~517 Cai Qiufang, Liu Yu. The development of a tree-ring width chronology and the May-June mean temperature variability in Wulu Mountain, Shanxi Province of north-Central China. Quaternary Sciences, 2013,33 (3):511~517 |
| 20 | 李宗善, 刘国华, 傅伯杰等. 利用树木年轮宽度资料重建川西米亚罗地区过去200年夏季温度的变化. 第四纪研究, 2011,,31 (3):522~534 Li Zongshan, Liu Guohua, Fu Bojie et al. Tree-ring based summer temperature reconstruction over the past 200 years in Miyaluo of western Sichuan, China. Quaternary Sciences, 2011,31 (3):522~534 |
| 21 | 陈 峰, 袁玉江, 魏文寿等. 福建沙县马尾松树轮宽度与夏季亚洲-太平洋涛动指数的关系. 第四纪研究, 2011,,31 (1):96~103 Chen Feng, Yuan Yujiang, Wei Wenshou et al. Correlations between the summer Asia Pacific Oscillation Index and the tree-ring width of Pinus massiniana from Sha County, Fujian Province. Quaternary Sciences, 2011,,31 (1):96~103 |
| 22 | Liang E Y, Liu X, Yuan Y et al. The 1920s drought recorded by tree rings and historical documents in the semi-arid and arid areas of Northern China. Climatic Change, 2006, doi:10.1007/s10584~006-9082-x |
| 23 | Zhang Q B, Guo D C, Tan D Y et al. A 2326-year tree-ring record of climate variability on the northeastern Qinghai-Tibetan Plateau. Geophysical Research Letters, 2003,30 (14):1739~1742 |
| 24 | 时兴合, 秦宁生, 刘洪滨等. 青海治多公元1374年以来树轮记录的青藏高原地面加热场强度变化. 第四纪研究, 2013,33 (1):115~125 Shi Xinhe, Qin Ningsheng, Liu Hongbin et al. Surface heating field intensity change over Tibetan Plateau recorded Qinghai Zhiduo tree rings since 1374A.D.Quaternary Sciences, 2013,33 (1):115~125 |
| 25 | 任军莉, 刘 禹, 宋慧明等. 甘肃临夏地区过去195年最高温度历史重建——基于紫果云杉树轮宽度资料. 第四纪研究, 2014,34 (6):1270~1279 Ren Junli, Liu Yu, Song Huiming et al. The historical reconstruction of the maximum temperature over the past 195 years, Linxia region, Gansu Province——Based on the data from Picea purpurea Mast. Quaternary Sciences, 2014,34 (6):1270~1279 |
| 26 | Melvin T M, Briffa K R. A "Signal-Free" approach to Dendroclimatic Standardisation. Dendrochronologia. 2008, 26 :71~86, doi:10.1016/j.dendro.2007.12.001 |
| 27 | Fang K, Wilmking M, Davi N et al.An ensemble weighting approach for dendroclimatology:Drought reconstructions for the northeastern Tibetan Plateau PloS One, 2014,9 (1):e86689, doi:10.1371/journal.pone.0086689 |
| 28 | 徐 岩, 邵雪梅. 柴达木盆地东缘祁连圆柏轮宽序列标准化的方法研究. 地理学报, 2006,61 (9):919~928 Xu Yan, Shao Xuemei. Standardization of Qilian Juniper ring-width series in the eastern margin of the Qaidam Basin. Acta Geographica Sinica, 2006,61 (9):919~928 |
| 29 | Briffa K R, Jones P D, Schweingruber F H. Tree-ring density reconstructions of summer temperature patterns across western North America since 1600. Journal of Climate, 1992,5 (7):735~754 |
| 30 | Briffa K R, Jones P D, Schweingruber F H et al. Tree-ring Variables as Proxy-climate Indicators:Problems with Low Frequency Signals. Berlin:Springer-Verlag, 1996.9~41 |
| 31 | Esper J, Cook E R, Schweingruber F H. Low-frequency signals in long tree-ring chronologies for reconstructing past temperature variability. Science, 2002,295 (5563):2250~2253 |
| 32 | Mann M E, Bradley R S, Hughes M K. North hemisphere temperature during the last millennium:Inferences, uncertainties and limitations. Geophysical Research Letters, 1999,26 (6):759~762 |
| 33 | Yang B, Sonechkin D M, Datsenko N M et al. Eigen analysis of tree-ring records:Part 1, A limited representativeness of regional curve. Theoretical and Applied Climatology, 2011,106 (3~4):489~497 |
| 34 | Yang B, Sonechkin D M, Datsenko N M et al. Eigen analysis of tree-ring records:Part 3, Taking heteroscedasticity and sampling effects into consideration. Theoretical and Applied Climatology, 2012,107 (3~4):519~530 |
| 35 | Yang B, Sonechkin D M, Datsenko N M et al. Eigen analysis of tree-ring records:Part 2, Posing the eigen problem. Theoretical and Applied Climatology, 2012,107 (1~2):131~141 |
| 36 | 刘晶晶, 杨 保. 西藏浪卡子地区树轮宽度特征值年表的建立. 第四纪研究, 2014,34 (6):1280~1287 Liu Jingjing, Yang Bao. Development of a tree-ring-width chronology using the Eigen analysis method in Langkazi regions of Tibet. Quaternary Sciences, 2014,34 (6):1280~1287 |
| 37 | Huang N E, Shen Z, Long R S et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary series analysis. Proceedings of the Royal SocietyA: Mathematical, Physical & Engineering Sciences, 1998,454 (1971):903~995, doi:10.1098/rspa.1998.0193 |
| 38 | 杨 周, 林振山, 俞鸣同. 东亚冬季风演变和亚洲内陆干旱化信号的多尺度分析. 第四纪研究, 2011,,31 (1):73~80 Yang Zhou, Lin Zhenshan, Yu Mingtong. Multi-scale analysis of East Asian winter monsoon evolution and Asian inland drying force. Quaternary Sciences, 2011,,31 (1):73~80 |
| 39 | 郭志永, 翟秋敏. 北大西洋IODP U1313站深海沉积中环境敏感粒度组分研究. 第四纪研究, 2012,32 (6):1107~1120 Guo Zhiyong, Zhai Qiumin. A study of environmentally sensitive grain-size groups of deep-sea sediment at IODP site U1313, North Atlantic. Quaternary Sciences, 2012,32 (6):1107~1120 |
| 40 | 冯 雯, 李可军, 翟 倩等. 中国北方地区(33~41°N, 108~115°E)过去500年夏季降水量的周期成分分析. 第四纪研究, 2013,33 (1):126~135 Feng Wen, Li Kejun, Zhai Qian et al. Multi-scale analysis of the precipitation index in north-Central China(33~41°N, 108~115°E)over the past 500 years. Quaternary Sciences, 2013,33 (1):126~135 |
| 41 | 林振山, 汪曙光. 近四百年北半球气温变化的分析:EMD方法的应用. 热带气象学报, 2004,20 (2):90~96 Lin Zhenshan, Wang Suguang. EMD analysis of northern hemisphere temperature variability during last 4 centuries. Journal of Tropical Meteorology, 2004,20 (2):90~96 |
| 42 | Fang Keyan, Frank David, Gou Xiaohua et al. Precipitation variations over the past four centuries in the Dieshan Mountain area as inferred from tree rings:An introduction to an HHT-based method. Global and Planetary Change, 2013, 107 :109~118, doi:10.1016/j.gloplacha.2013.04.010 |
| 43 | 《青海森林》编委会. 青海森林. 北京: 中国林业出版社, 1993. 233~238 Editorial Committee of "Forest in Qinghai". Forest in Qinghai. Beijing:China Forestry Press, 1993. 233~238 |
| 44 | 刘兴聪. 东祁连山西段云杉生长规律的调查研究. 甘肃农大学报, 1983,1 (3):60~69 Liu Xingcong. On investigation and studies on the growth laws of Picea crassifolia in the western of the east Qilianshan. Journal of Gansu Agricultural University, 1983,1 (3):60~69 |
| 45 | 何志斌, 赵文智, 张立杰等. 祁连山青海云杉林采伐干扰与恢复过程. 林业科学, 2009,45 (2):12~16 He Zhibin, Zhao Wenzhi, Zhang Lijie et al. Logging disturbance and the restoration process of Picea crassifolia forest in Qilian Mountain. Scientia Silvae Sinicae, 2009,45 (2):12~16 |
| 46 | Gordon A H, Byron-Scott R A, Bye J A. A note on QBO-SO interaction, the quasi-triennial oscillation and the sunspot cycle. Journal of Atmospheric Sciences, 1982,39 (9):2083~2087 |
| 47 | Naujokat B. An update of the observed Quasi-Biennial Oscillation of stratospheric winds over the tropics. Journal of Atmospheric Sciences, 1986,43 (17):1873~1877 |
| 48 | Allan R. ENSO and climatic variability in the past 150 years. In:Diaz H, Markgraf V eds. El Niño and the Southern Oscillation:Multiscale Variability and Global and Regional Impacts. Cambridge:Cambridge University Press, 2000. 3~35 |
| 49 | Prabhakaran Nayar S R, Radhika V N, Revathy K et al. Wavelet analysis of solar wind and geomagnetic parameters. Solar Physics, 2002, 208 (2):359~373 |
Abstract
The relationship between climate factors and tree growth are very complex. Traditional standardization methods still face challenges in retaining the climate signal while removing age-related growth trends and non-climatic disturbances. In this paper, we introduce an Empirical Mode Decomposition(EMD)based standardization method to develop the ring-width chronology, which has the advantages of excluding non-climatic disturbances in individual tree-ring series. First, tree-ring measurements were decomposed by the EMD to produce a number of Intrinsic Mode Functions(IMFs) with different physically-meaningful frequency band. Second, the lowest frequency IMF band, mainly representing age-dependent growth trend, was used to calculate tree-ring index divided by the original tree-ring width records. Based on this method, a new tree-ring chronology was developed from 44 Spruce trees(Picea crassifolia Kom.)at Kong-gang-mu(Kgm)Forest of Sunan County of Qinghai Province in the Qilian Mountains, and the elevations ranged from 2620m to 2770m, which covered the period 1723A.D.to 2000A.D.However, the chronology was more reliable before 1845A.D. when the EPS was greater than 0.85.
Compared two ring-width chronologies from EMD method and traditional method, we can found that the EMD method has the significant advantages in extracting the age-dependent growth trend, especially for improving fitting on the early "slow-growth" and the late "flat-growth" periods. It is noted that two chronologies exhibited the coincident significant correlation with meteorological factors of Sunan Meteorological Station(38°50'N, 99°37'E; 2311m a.s.l.)in observation period(1957~2000), showed the highest correlation coefficient(r=0.664, p <0.01) with the total precipitation of the year from the previous August to current July, which indicated that the EMD method successfully preserved climatic signals. Although the EMD method appears to offer both practical and conceptual advantages, we emphasize that it is premature at this stage to conclude that the EMD-based method is better than the traditional methods. However, EMD method provides an additional option to mitigate bias in low-frequency variations longer than mean segment lengths of tree rings, which is inevitable with traditional methods.
The spectral analysis showed three main wavebands of this chronology: less than 10-yrs period, 10~30yrs, 30~60yrs, which would be integrated to three components of IMFs of all original ring width series derived from EMD method. The 10-yrs scale of IMF component was relatively stable in the time-frequency domain, representing a major portion of the high-frequency precipitation variations, with the maximum amplitude around the 1910s~1920s. The 10~30yrs scale of IMF component displayed the obvious lowest amplitude during the 1925~1933, corresponding to an anomalously dry period at that time of Gansu Province. And the 30~60yrs scale of IMF component showed that three dry-wet fluctuation. The integration components of the tree-ring width series with different main periods revealed the temporal distribution and detailed fluctuations of the major periods by their corresponding amplitudes in the time-frequency domain, which provide more information of frequency-domains of reconstructed climate series at multiple time-scales.