全 文 ：第 35 卷第 22 期
2015年 11月
生 态 学 报
ACTA ECOLOGICA SINICA
Vol.35,No.22
Nov., 2015
http: / / www.ecologica.cn
基金项目:国家自然科学基金资助项目(41171067)
收稿日期:2014鄄03鄄30; 摇 摇 网络出版日期:2015鄄04鄄20
*通讯作者 Corresponding author.E鄄mail: jiangy@ bnu.edu.cn
DOI: 10.5846 / stxb201403300597
赵守栋,江源,焦亮,王明昌,张凌楠,李文卿.ARSTAN程序和 R语言 dplR扩展包进行树轮年表分析的比较研究.生态学报,2015,35(22):7494鄄7502.
Zhao S D, Jiang Y, Jiao L, Wang M C, Zhang L N, Li W Q.A Comparative analysis of ARSTAN and the dplr package of R language in analyses of tree鄄
ring chronologies.Acta Ecologica Sinica,2015,35(22):7494鄄7502.
ARSTAN程序和 R 语言 dplR 扩展包进行树轮年表分
析的比较研究
赵守栋,江摇 源*,焦摇 亮,王明昌,张凌楠,李文卿
北京师范大学资源学院, 北京摇 100875
摘要:在树轮年代学领域,ARSTAN是去趋势处理和建立年表方面应用最为广泛的程序,而新兴的 R 语言 dplR 扩展包实现了
ARSTAN的主要功能,且具有源代码公开、扩展性强等优点,是传统程序的良好补充。 使用贺兰山青海云杉(Picea Crassifolia)树
轮宽度数据,分析了 ARSTAN和 dplR进行树轮年代学分析所得结果的差异。 结果显示,两种程序计算平均敏感度和一阶自相
关系数的平均误差为 0.005—0.008,但具有确定的转换关系;两种程序如果使用同种方法去趋势,拟合曲线的参数相近,建立标
准年表的平均误差为 0.002;拟合自回归模型时差异较大,其中时域上表现为差值年表起始 30a 内差异显著,在频域上表现为
dplR的差值年表保留的低频信息较少;年表统计量计算和公共区间分析中,不同程序计算样本总体代表性和信噪比的差异较
大。 分析表明,两程序在拟合生长趋势和自回归模型时存在算法上的较大差异,同时年表统计量和公共区间各指标的算法也不
尽一致,但存在较为确定的转换关系。 对开展不同来源数据的整合分析提出了建议,应明确不同研究中树轮数据的处理过程,
在条件允许时使用同一程序或算法重新处理数据,确保结果的可比性。
关键词:R语言;dplR;去趋势;公共区间分析;树轮年代学
A comparative analysis of ARSTAN and the dplr package of R language in
analyses of tree鄄ring chronologies
ZHAO Shoudong, JIANG Yuan*, JIAO Liang, WANG Mingchang, ZHANG Lingnan, LI Wenqing
College of Resources Science and Technology, Beijing Normal University, Beijing 100875, China
Abstract: Dendrochronology plays an important role in estimating past climatic conditions and predicting future climate
change. Detrending and chronology development are the fundamental steps of the study of dendrochronology. ARSTAN is the
most popular program used to accomplish this step, and it has played an important role in the development of
dendrochronology. However, ARSTAN uses the Fortran programming language, so users find it difficult to understand and
revise the algorithm of the source program to meet their needs. An emerging package of the R language named dplR provides
similar functions to ARSTAN. R and dplR爷 s source code is fully open to the public; thus, it has numerous users. When
scholars from different domains communicate and share the methods and results of dendrochronology, it can help them
improve those chronologies. In addition, R and dplR have become a good supplement to traditional analysis software. This
paper compares the different dendrochronological analysis algorithms and results provided by ARSTAN and dplR with tree
ring width data from Picea crassifolia on Helan Mountain, Ningxia Hui Autonomous Region, China. The results show that
the two programs calculated exactly the same means and standard deviations. The mean error of the mean sensitivities (MS)
and first鄄order autocorrelations (AC) were 0.005 and 0.008, respectively, but they had a clear conversion relationship.
When using the same method for detrending with both types of software, the parameters of fitting curves were generally
http: / / www.ecologica.cn
equal, and the corresponding standard chronologies developed by the two programs had a mean error of only 0. 002.
However, the residual chronologies were very different. In the time domain, a significant difference was observed in the
residual chronologies in the first 20—30 years. In the frequency domain, the residual chronologies created using ARSTAN
showed more low frequency information than that created using dplR. For example, the former showed periods of 32 years
with higher power than those of dplR. In the common interval analysis, ARSTAN gave a higher expressed population signal
(EPS) and signal鄄to鄄noise ratio (SNR) of chronologies than dplR. EPS error was 0.4% and SNR error was 30%—40%. By
comparing the algorithms of the two programs, we found that ARSTAN and dplR have different initial value setting rules and
nonlinear fitting methods to choose the best fitting model during detrending. When fitting an autoregression model, ARSTAN
used a pooled algorithm to find the integral growing pattern and used the same fitting order for different sequences. However,
dplR directly used different optimal fitting models for different sequences. In addition, the two programs used different, but
similar, formulas for calculating MS, AC, EPS, and SNR. Although the absolute value of the results was different,
calculation results of the same program using different data were comparable. In conclusion, this paper offers two suggestions
for the meta鄄analysis of tree ring data from different sources. First, if the source data are available, researchers should
choose a single program for statistical calculation, detrending, and common interval analysis based on their needs. Second,
if the source data are not available, information related to the chronologies is sufficient; researchers should use only a single
program to calculate EPS and SNR chronological statistics to ensure that the results will be comparable.
Key Words: R language; dplR; detrending; common interval analysis; dendrochronology
摇 摇 树轮年代学研究树木年轮的宽度、密度等性质构成的时间序列及其与环境因子的响应关系[1],在反演过去和预
测未来气候变化中扮演着重要作用[2鄄5]。 树轮年代学中,去除树木自身的生长趋势并建立树轮年表是开展研究的基
础步骤之一。 在去趋势处理和建立年表方面,Cook等编写的 ARSTAN程序是目前最权威的分析工具,在树轮年代
学的发展过程中发挥了巨大作用。 然而 ARSTAN基于 Fortran语言编写,同时缺乏系统而详细的说明文档,使用者
不容易理清其内部的算法,实现新的分析方法也比较困难,这在一定程度上限制了研究的深入和学科间的交流。
R语言是当前发展较快的统计分析工具之一,在众多研究领域得到了越来越广泛的应用。 Bunn 等使用
R语言编写了 dplR扩展包,实现了树轮数据读取、转换、分析、绘图等功能[6]。 相比于经典工具而言,R 语言
和 dplR等相关扩展包的源代码完全公开,且容易理解,研究者可以根据自身需要修改或编写程序;同时使用
者众多,网络社区发达,有利于不同领域的相互交流。 因此,R 语言及 dplR 扩展包是 ARSTAN 等传统程序的
良好补充,对于树轮年代学的进一步发展具有重要意义。
目前在国外研究中,dplR被广泛应用于去趋势处理和建立年表[7鄄9],在交叉定年[10]、小波分析[3]等方面
也发挥了作用;在其基础上又发展出交互式去趋势[11]、响应函数分析[12]、扰动事件重建[13]等相关扩展包,逐
渐形成了一套较为完整的树轮数据分析软件体系。 而在国内树轮年代学领域,使用 dplR开展的研究还较少,
使用的功能也相对单一[14]。 作为一个新兴的树轮年表分析工具,dplR处理结果与传统程序一致与否,将直接
影响到不同来源研究数据的可比性,因此有必要对 ARSTAN和 dplR在算法和结果上进行比较分析。
本文使用贺兰山青海云杉树轮宽度数据,分别借助 ARSTAN和 dplR进行原始序列统计量计算、去趋势处
理、年表建立、公共区间分析等操作,通过相关分析、小波分析等手段比较了两者所得结果的差异,并从算法的
角度探讨了差异的来源,为开展树轮年代学数据的整合分析提出了需要注意的问题及相应建议。
1摇 数据与方法
1.1摇 树轮宽度资料与气候资料
本文中使用的树轮资料为贺兰山青海云杉(Picea Crassifolia)树轮宽度数据[15],采样点位于贺兰山东坡
林线附近,共包含来自 25棵样树 50根样芯的树轮宽度测量结果,数据精度为 0.01 mm。 使用 COFECHA程序
检验交叉定年质量,结果显示序列平均长度为 83.9 a,平均敏感度为 0.393,序列间相关系数为 0.825。
5947摇 22期 摇 摇 摇 赵守栋摇 等:ARSTAN程序和 R语言 dplR扩展包进行树轮年表分析的比较研究 摇
http: / / www.ecologica.cn
1.2摇 ARSTAN处理
ARSTAN程序目前分为 Windows XP 版和 OS X 版,早期的 DOS 版已不再更新。 本文主要使用 Windows
XP 版 ARSTAN进行去趋势处理并建立年表,程序下载自美国哥伦比亚大学树木年轮实验室(www. ldeo.
columbia.edu / tree鄄ring鄄laboratory),版本号为 44。 选择的去趋势方法为修正的负指数函数,即:
f( i) = a 伊 exp( - b 伊 i) + d (1)
式中,f( i)为生长趋势估计值,i为年份的序数。 设定约束条件为 d>0,拟合失败则拟合斜率非正的线性函数。
公共区间设定为所有树芯均含有的共同时间段,即 1966—2009 年。 ARSTAN 最终建立标准年表( standard
chronology, STD)、差值年表 ( residual chronology, RES)和自回归年表 ( autoregressive standard chronology,
ARS),并给出年表的基本统计量和公共区间分析结果。
1.3摇 dplR处理
R语言主程序和 dplR扩展包均下载自综合 R档案网络(cran.r鄄project.org),其中 R语言版本号为 3.1.1,dplR版
本号为 1.6.1。 使用“rwl.stats冶函数计算原始树轮宽度序列和树轮宽度指数年表的基本统计量。 使用“detrend.series冶
函数对原始树轮宽度序列逐一进行去趋势操作,确保各样芯使用的去趋势方法与 ARSTAN中相同,使得结果具有可
比性。 使用“rwi.stats冶函数对树轮宽度指数序列进行公共区间分析,时间段同样设定为 1966—2009年,这也是程序
默认的公共区间。 使用“chron冶函数采用双权重平均法建立年表,最终得到 STD和 RES年表。
1.4摇 ARSTAN与 dplR结果对比
比较 ARSTAN和 dplR所得原始轮宽序列统计量的差异,使用 Spearman相关系数考察误差与序列长度的
关系。 对比 ARSTAN和 dplR拟合生长趋势所得曲线的参数及拟合结果,使用 Spearman相关系数考察两程序
所得 STD、RES年表的相关性,以 ARSTAN所得结果为准,计算 dplR所得年表的标准误;针对差异相对较大的
RES年表,使用 Morlet小波变换分析两程序所得 RES 年表在频域上的差异。 最后比较不同程序在年表统计
量和公共区间分析结果上的差异。 相关分析中,显著性水平设置为 0.05。
2摇 结果
2.1摇 原始轮宽序列统计量
针对 50个原始轮宽序列,ARSTAN 和 dplR 计算所得均值和标准差完全相等;而 dplR 所得平均敏感度
(mean sensitivity, MS)和一阶自相关系数( first鄄order autocorrelation, AC)均略小于 ARSTAN 所得结果,其中
MS的平均误差为 0.005,标准误为 0.0002,AC的平均误差为 0.008,标准误为 0.0003,两者的误差与序列长度
存在显著负相关关系(P<0.05)。 进一步计算发现,ARSTAN和 dplR所得 MS和 AC满足:
MSr = MSa 伊
n - 2
n - 1
(2)
ACr = ACa 伊
n - 1
n
(2)
式中,MSr为 dplR所得 MS,MSa为 ARSTAN所得 MS,ACr为 dplR 所得 AC, ACa为 ARSTAN 所得 AC,n 为序列
长度。 随着序列长度的增加,ARSTAN 与 dplR 所得结果的差异将逐渐减小,当序列长度在 100 年以上时,相
对误差可以降至 1%以下。
2.2摇 建立树轮宽度指数年表
ARSTAN和 dplR拟合所得负指数函数的参数 a、b完全相同,参数 d的误差在 0.0003以内,而线性函数参
数完全相同。 去趋势所得各样芯的两个标准化序列间存在差异的年份不超过 10 年,且误差均在 0.003 以内。
在此基础上,两程序建立的 STD年表相关系数为 0.999(P<0.05),对应年份宽度指数误差在 0.05 以内,平均
误差为 0.005,标准误为 0.0010,相位变化情况基本相同(图 1)。
在建立差值年表时,ARSTAN对各标准化序列均拟合为一阶自回归模型,所得差值序列与原序列长度相
等;而 dplR则拟合成阶数各异的自回归模型,包括 29 个一阶模型和 17 个二阶或更高阶模型,序列初始损失
6947 摇 生摇 态摇 学摇 报摇 摇 摇 35卷摇
http: / / www.ecologica.cn
的年份数与模型阶数相等,另有 4个样芯未拟合成自回归模型。 各样芯的两个差值序列间的平均误差为 0.
034,标准误为 0.0011。 两程序建立的 RES 年表间相关系数为 0.987(P<0.05),在 1914 年前(即序列前 30
年),两者平均误差为 0.048,标准误为 0.0096;而 1914 年后,两者平均误差为 0.028,标准误为 0.0021。 两个
RES年表仅在起始的 10年内存在相位不一致的情况,例如 1887—1889年(图 1)。
对两个 RES年表分别进行 Morlet小波变换(图 2)。 结果显示,两个 RES年表含有的高频信息基本相同,
图 1摇 ARSTAN和 dplR所得树轮宽度指数年表
Fig.1摇 Tree鄄ring chronologies produced by ARSTAN and dplR
(a)为标准年表;(b)为差值年表;(c)为不同程序建立的年表之差
图 2摇 差值年表小波功率谱
Fig.2摇 Wavelet power spectrum of the residual chronologies
7947摇 22期 摇 摇 摇 赵守栋摇 等:ARSTAN程序和 R语言 dplR扩展包进行树轮年表分析的比较研究 摇
http: / / www.ecologica.cn
均表现出 1940年前 2—11a尺度的震荡能量较小,周期变化不明显,在 1940 年后 2—11a 尺度的能量增强,周
期变化显著;而在低频信息方面,ARSTAN建立的 RES年表在 1940 年后稳定存在着 32a 左右的周期,而 dplR
建立的 RES年表 32a左右的周期震荡能量始终较小。
2.3摇 年表统计量和公共区间分析结果
年表统计量及公共区间分析结果见表 1。 在年表统计量方面,ARSTAN 的 STD、RES 年表和 dplR 的 STD
年表的序列起止年相同,而 dplR的 RES年表与这 3个年表相比少 1年(1884 年)。 两个 STD年表的均值、标
准差均相差 0.001;MS相差 0.003,相对误差约为 1%;AC 相差 0.005,相对误差约为 1.5%。 而两个 RES 年各
项统计量的误差均大于 STD年表,其中 dplR的 RES年表 MS比 ARSTAN低 4%,而 AC的绝对值低 36%。
在公共区间分析方面,两个 STD年表的所有样芯平均相关系数(all鄄series Rbar, Rtot)、树内平均相关系
数(within鄄trees Rbar, Rwt)、树间平均相关系数(among鄄trees Rbar, Rbt)和年表有效信号(effective chronology
signal, Reff)均完全相等,而 ARSTAN的 RES年表在这些指标上普遍略高于 dplR 的 RES 年表。 dplR 的 STD
和 RES年表的样本总体代表性(expressed population signal, EPS)比 ARSTAN 计算结果低 0.004,相对误差约
为 0.4%,而信噪比(signal鄄to鄄noise ratio, SNR)结果低 30左右,相对误差高达 30%—40%。
表 1摇 年表统计量及公共区间分析结果
Table 1摇 Chronologies statistics and common interval analysis
程序 Program ARSTAN dplR
年表类型 Chronology STD RES STD RES
年表统计量
Chronologies statistics
序列长度 Time span 1884—2010 1884—2010 1884—2010 1885—2010
Mean 0.966 0.986 0.967 0.993
SD 0.308 0.283 0.309 0.280
MS 0.298 0.348 0.295 0.333
AC 0.312 -0.137 0.317 -0.087
公共区间分析
Common interval analysis
公共区间
Common interval 1966—2009 1966—2009 1966—2009 1966—2009
样树 /样芯 Trees / Cores 25 / 50 25 / 50 25 / 50 25 / 50
Rtot 0.674 0.728 0.674 0.725
Rwt 0.808 0.850 0.808 0.843
Rbt 0.671 0.725 0.671 0.722
Reff 0.743 0.784 0.743 0.784
EPS 0.990 0.993 0.986 0.989
SNR 103.477 133.745 72.125 90.525
摇 摇 Mean: 平均指数 mean indices;SD: 标准差 standard deviation;MS: 平均敏感度 mean sensitivity;AC: 一阶自相关系数 first鄄order autocorrelation;
Rtot: 所有样芯平均相关系数 all鄄series Rbar;Rwt: 树内平均相关系数 within鄄trees Rbar;Rbt: 树间平均相关系数 among鄄trees Rbar;Reff: 年表有效
信号 effective chronology signal;EPS: 样本总体代表性 expressed population signal;SNR: 信噪比 signal鄄to鄄noise ratio
3摇 讨论
3.1摇 ARSTAN与 dplR计算序列统计量的算法差异
根据原始轮宽序列统计量比较结果,ARSTAN 与 dplR 在 MS 和 AC 的计算上存在着较为确定的换算关
系,这很可能是由于两者使用了不同的计算公式。 dplR计算 MS的公式为[16]:
MS = 2
n - 1移
n
i = 2
| w i - w i -1
w i + w i -1
(4)
式中,n为序列长度,w i为第 i年的树轮宽度。 根据结果推测,ARSTAN在计算 MS 时,将公式中的 n-1 换成了
n-2。
dplR调用 R中的“acf冶函数计算 AC,其计算公式可表达为[1,17]:
8947 摇 生摇 态摇 学摇 报摇 摇 摇 35卷摇
http: / / www.ecologica.cn
AC =
移
n
i = 2
(w i - —w )
移
n
i = 1
(w i - —w ) 2
=
移
n
i = 2
(w i - —w )
(n - 1) s2
(5)
式中,s为序列的标准差。 ARSTAN在计算 AC时,可能将 s处理成了总体标准差:
s =
移
n
i = 1
(w i - —w ) 2
n
(6)
尽管计算方法不同,由于两指标相对 n的取值具有单调性,因此同一个程序的计算结果具有可比性,而不同程
序的计算结果需要经过换算后才可相互比较。
3.2摇 ARSTAN与 dplR建立年表的算法差异
ARSTAN和 dplR均具有多种去趋势方法,两者的对应关系如表 2所示,其中在拟合负指数函数的程序实
现方面,两者具有较大差异。 在算法上,ARSTAN使用 Fritts等专为负指数函数拟合问题设计的算法[18],该方
法充分考虑到负指数函数的曲线特征,首先估计曲线的切线斜率 b,通过调节 b 的值进行迭代。 而 dplR 首先
估计 y轴截距 a和渐近线 d的初始值,b的初始值则根据经验直接设为 0.01,再调用 R内置的“nls冶函数进行
曲线参数估计。 如果不设定约束条件,当树轮宽度序列的趋势较符合负指数函数曲线时,两种算法均能得到
精确的拟合结果;而当树轮宽度序列的趋势近似于直线时,d 为负值,这时 dplR 的初始值设定严重偏离最优
解,导致迭代过程报错,转而使用线性函数进行拟合。 如果设定约束条件 d>0,ARSTAN首先在整个定义域范
围内寻找各参数的最佳估计值,然后将参数带入约束条件中进行判断,如果不符合条件则拒绝拟合负指数函
数,转而拟合线性函数;而 dplR则是在“nls冶函数中缩小参数的取值范围,最终拟合得到参数的局部最优解。
当拟合方法出现分歧时,去趋势所得序列将出现较大差异,而本文在操作中参考 ARSTAN 的运行结果,使用
dplR对每根样芯单独去趋势,确保两者使用的去趋势方法相同。
尽管拟合生长趋势的方法相同,ARSTAN和 dplR建立的标准年表仍存在一定差异,这与缺失值的处理方
式有关。 为防止运算中可能出现的错误,dplR将序列中的缺失值(一般赋值为 0)均替换成 0.001,因此 dplR
所得标准化序列不存在树轮宽度指数为 0 的情况;而 ARSTAN并未对缺失值进行特殊处理,因此标准化序列
中仍存在 0值。
表 2摇 ARSTAN和 dplR去趋势方法的对应关系
Table 2摇 Corresponding relationships between detrending methods of ARSTAN and dplR
去趋势方法 Detrending methods ARSTAN选项Detrending options
dplR中函数与参数设置
Functions and parameters in dplR
均值水平线
Horizontal line (arithmetic mean) 6 detrend(method
= “Mean冶, …)
修正的负指数函数
Modified negative exponential curve
摇 无约束条件
摇 Unconstrained 3 detrend(method
= “ModNegExp冶, …)
摇 有约束条件,拟合失败则采用斜率非正的直线
摇 Constrained, or a line with a negative or zero slope 2
detrend(method = “ModNegExp冶,
constrain.modnegexp = “always冶, …)
摇 有约束条件,拟合失败则采用任意斜率的直线
摇 Constrained, or a line with any slope 1
detrend(method = “ModNegExp冶,
constrain.modnegexp = “always冶,pos.slope = TRUE, …)
区域曲线标准化
Regional Curve Standardization
-2 rcs( ...)
步长为序列长度 n%的样条
Spline (n% of the series length)
-n detrend(method = “Spline冶, nyrs = n, …)
摇
9947摇 22期 摇 摇 摇 赵守栋摇 等:ARSTAN程序和 R语言 dplR扩展包进行树轮年表分析的比较研究 摇
http: / / www.ecologica.cn
建立差值年表,需要去除树轮宽度指数标准化序列中的自回归成分。 ARSTAN首先用所有标准化序列拟
合出一个 p阶合并自回归模型,将拟合值作为树木群体所共有的持续性生长量。 再对各标准化序列分别拟合
限定为 p阶的自回归模型,将拟合值作为整体生长模式下个体的持续性生长量,而将残差作为差值序列,并使
用双权重平均法建立 RES年表。 如果将群体的持续性生长量加回到 RES年表上,则得到 ARS年表[19]。 dplR
未考虑整体的生长模式,而是直接对各标准化序列分别拟合自回归模型,不同序列的模型可以具有不同的阶
数,再将残差作为差值序列,建立 RES年表。 针对本文所用数据,dplR 的计算过程对个体特有的持续性生长
量剔除得更为彻底,低频信息保留得更少,这种差异在频域上更为明显(图 2)。
拟合 p阶自回归模型会造成序列前 p年树轮宽度指数的损失。 针对这一情况,ARSTAN使用序列均值将
标准化序列向前延长 p年[19],从而使自回归模型计算出的拟合值序列与原序列等长,避免了序列长度上的损
失。 而 dplR未对初始年份的损失进行特殊处理,因此两个 RES 年表存在的差异相对 STD 年表更大,尤其是
样本量较小的最初若干年。 根据本文结果,起始值处理方式造成的差异会在约 30a 后消失,同时考虑到参与
“树轮鄄气候冶关系分析的年表区间通常远小于年表总长度,如果起始年份不参与后续的分析,则起始年份的不
同处理方式对最终结果的影响较小[19鄄20]。
3.3摇 ARSTAN与 dplR进行公共区间分析的算法差异
公共区间分析是指选择包含特定时间范围的树轮宽度序列,通过计算年表信号,评价最终所得年表对森
林总体生长的代表性,其中年表信号是衡量序列中包含的共同变化的统计量[21]。 DOS 版 ARSTAN 在导入测
量数据时会提供一个推荐的公共区间范围,其选择标准是公共区间内包含尽可能多的测量值,即满足区间内
完整样芯数与区间长度乘积取得最大值;而 Windows XP 版 ARSTAN不再提供这一结果。 dplR中对数据集进
行去趋势使用“detrend冶函数,其默认公共区间是所有样芯共有的公共区间,该区间长度一般小于 DOS 版
ARSTAN提供的公共区间长度。 dplR中的“common.interval冶函数可以给出各种标准对应的时间范围。
在进行公共区间分析时,最初使用方差分析对年表信号进行估计[1],后来的研究发现使用相关矩阵也可
以估计年表信号的大小[22]。 以序列间的平均相关系数(mean series intercorrelation, Rbar)作为年表信号的衡
量指标,可以通过计算 EPS和 SNR来评价年表的质量[17]:
EPS = t
伊 Rbar
t 伊 Rbar + (1 - Rbar)
(7)
SNR= t
伊Rbar
(1-Rbar)
(8)
式中,t为样树数量。 对于某树轮年表,如果每棵样树均取一根样芯,则 Rbar 的取值为 Rbt;如果每棵样树取
多根样芯,则 Rbar的取值应当为 Rbt与 Rwt加权平均结果,即 Reff,如果仍以 Rbt作为年表信号则会低估样本
中的公共信息[21]。
根据运算结果推测,ARSTAN可能使用样芯数量 m和 Rtot计算 EPS 和 SNR。 这种算法将所有样芯等同
看待,没有剔除同一样树的不同样芯包含着的重复信息,因此可能高估了年表质量。 EPS 除了用于评价年表
整体质量外,还用于计算子样本信号强度(subsample signal strength, SSS) [21],而 EPS 和 SSS 均是衡量年表可
靠长度的重要指标。 由于 ARSTAN 计算得出的 EPS 偏高,可能会令 EPS 或 SSS 高于阈值(一般设为 0.8 或
0.85)的最小样本量偏低,造成对年表可靠长度的高估。
目前的树轮年代学研究中,一棵样树往往取两个甚至多个样芯,因此 dplR采用 Reff计算 EPS 和 SNR,计
算结果更符合现有理论。 使用基于 Reff 的公式重新评估本文中 ARSTAN 建立的年表,则 STD 年表 EPS 为
0.986,SNR为 72.276,而 RES年表 EPS为 0.989,SNR为 90.741,这些结果均与 dplR所得结果接近。
3.4摇 对比不同来源结果时应当注意的问题
随着树木年轮基础数据日渐丰富,对前期样点尺度的研究成果进行整合,开展区域、大陆甚至全球尺度上
森林对于全球气候变化的响应分析已成为研究热点之一[23鄄24]。 根据本文所得结果,在树木年轮研究中如果
使用不同程序处理原始数据,所得年表及其质量评估需要考虑算法差异可能造成的影响。 在标准年表建立方
0057 摇 生摇 态摇 学摇 报摇 摇 摇 35卷摇
http: / / www.ecologica.cn
面,即使 ARSTAN与 dplR拟合生长趋势的参数设定相同,拟合结果可能差异很大。 而对于差值年表,在频域
上,ARSTAN建立的 RES年表比 dplR具有更多的低频信息;在时域上,起始年份处理方式的不同,造成年表在
最初的 20—30年内存在较明显的差异。 不同程序在 MS、AC等指标上可能存在计算方法差异;ARSTAN使用
Rtot计算 EPS和 SNR等指标,相对于文献中基于 Reff的算法而言,可能高估了年表中包含的共同信息[25鄄26]。
综合以上讨论,对开展不同来源树轮数据的整合分析提出如下建议:
(1)如果能够在国际树轮数据库等共享平台上获取原始数据,建议使用同一程序完成年表建立和公共区
间分析等处理步骤,并注意检查各样芯去趋势方法的具体参数。
(2)如果仅能获得文献资料,建议记录其处理过程的详细信息,包括软件的版本信息、去趋势的具体方法
等;当信息足够时,建议根据同一标准设定公共区间,并重新计算 MS、EPS等指标,确保结果的可比性。
4摇 结论
在树轮年代学研究中,新兴的 R语言 dplR 扩展包是对 ARSTAN 等传统程序的良好补充,其平台开放性
可以促进研究者对程序算法的认识和跨学科的相互交流。 本文对比了 ARSTAN程序和 dplR 扩展包,分别使
用两种程序对贺兰山青海云杉树轮宽度数据进行原始序列统计量估计,完成了去趋势处理并建立了标准年表
和差值年表,通过公共区间分析评价了年表的质量。 分析表明,ARSTAN与 dplR在拟合生长趋势和自回归模
型时存在算法上的较大差异,同时平均敏感度、一阶自相关系数、样本总体代表性和信噪比等指标的算法也不
一致,但存在较为确定的转换关系。 整合分析不同来源的树轮资料时,应明确其数据处理过程,在条件允许的
情况下使用同一程序或算法进行去趋势处理和公共区间分析。
参考文献(References):
[ 1 ]摇 Fritts H C. Tree Rings and Climate. London: Academic Press, 1976:10鄄11,254鄄260.
[ 2 ] 摇 彭剑峰, 勾晓华, 陈发虎, 刘普幸, 张永, 方克艳. 阿尼玛卿山地不同海拔青海云杉(Picea crassifolia)树轮生长特性及其对气候的响应.
生态学报, 2007, 27(8): 3268鄄3276.
[ 3 ] 摇 Pederson G T, Gray S T, Woodhouse C A, Betancourt J L, Fagre D B, Littell J S, Watson E, Luckman B H, Graumlich L J. The unusual nature
of recent snowpack declines in the North American Cordillera. Science, 2011, 333(6040): 332鄄335.
[ 4 ] 摇 Liu H Y, Park Williams A, Allen C D, Guo D L, Wu X C, Anenkhonov O A, Liang E Y, Sandanov D V, Yin Y, Qi Z H, Badmaeva N K. Rapid
warming accelerates tree growth decline in semi鄄arid forests of Inner Asia. Global Change Biology, 2013, 19(8): 2500鄄2510.
[ 5 ] 摇 芦晓明, 梁尔源. 灌木年轮学研究进展. 生态学报, 2013, 33(5): 1367鄄1374.
[ 6 ] 摇 Bunn A G. A dendrochronology program library in R (dplR). Dendrochronologia, 2008, 26(2): 115鄄124.
[ 7 ] 摇 Bader M K F, Leuzinger S, Keel S G, Siegwolf R T W, Hagedorn F, Schleppi P, K觟rner C. Central European hardwood trees in a high鄄CO2
future: synthesis of an 8鄄year forest canopy CO2 enrichment project. Journal of Ecology, 2013, 101(6): 1509鄄1519.
[ 8 ] 摇 Boden S, Kahle H P, Wilpert K V, Spiecker H. Resilience of Norway spruce (Picea abies (L.) Karst) growth to changing climatic conditions in
Southwest Germany. Forest Ecology and Management, 2014, 315: 12鄄21.
[ 9 ] 摇 Rodr侏guez鄄Gonz佗lez P M, Campelo F, Albuquerque A, Rivaes R, Ferreira T, Pereira J S. Sensitivity of black alder (Alnus glutinosa [L.]Gaertn.)
growth to hydrological changes in wetland forests at the rear edge of the species distribution. Plant Ecology, 2014, 215(2): 233鄄245.
[10] 摇 Martin A R, Caspersen J P, Fuller M M, Jones T A, Thomas S C. Temporal dynamics and causes of postharvest mortality in a selection鄄managed
tolerant hardwood forest. Forest Ecology and Management, 2014, 314: 183鄄192.
[11] 摇 Campelo F, Garc侏a鄄Gonz佗lez I, Nabais C. detrendeR – A Graphical User Interface to process and visualize tree鄄ring data using R.
Dendrochronologia, 2012, 30(1): 57鄄60.
[12] 摇 Zang C, Biondi F. Dendroclimatic calibration in R: The bootRes package for response and correlation function analysis. Dendrochronologia, 2013,
31(1): 68鄄74.
[13] 摇 Altman J, Fibich P, Dolezal J, Aakala T. TRADER: A package for tree ring analysis of disturbance events in R. Dendrochronologia, 2014, 32
(2): 107鄄112.
[14] 摇 王蔚蔚, 张军辉, 戴冠华, 王秀秀, 韩士杰, 张寒松, 王云. 利用树木年轮宽度资料重建长白山地区过去 240年秋季气温的变化. 生态学
杂志, 2012, 31(4): 787鄄793.
1057摇 22期 摇 摇 摇 赵守栋摇 等:ARSTAN程序和 R语言 dplR扩展包进行树轮年表分析的比较研究 摇
http: / / www.ecologica.cn
[15]摇 赵守栋, 何新, 张冰琦, 刘琦, 王辉, 江源. 贺兰山东坡高山林线青海云杉(Picea crassifolia)径向生长对气候因子的响应. 北京师范大学
学报: 自然科学版, 2013, 49(5): 501鄄505.
[16] 摇 Biondi F, Qeadan F. Inequality in paleorecords. Ecology, 2008, 89(4): 1056鄄1067.
[17] 摇 Cook E R, Pederson N. Uncertainty, emergence, and statistics in dendrochronology / / Hughes M K, Swetnam T W, Diaz H F, eds.
Dendroclimatology: Progress and Prospects. Netherlands: Springer, 2011: 77鄄112.
[18] 摇 Fritts H C, Mosimann J E, Bottorff C P. A revised computer program for standardizing tree鄄ring series. Tree鄄Ring Bulletin, 1969, 29(1鄄2): 15鄄20.
[19] 摇 Cook E R. A Time Series Analysis Approach to Tree Ring Standardization[D]. Tucson, Arizona, USA: University of Arizona, 1985:63鄄73.
[20] 摇 Cook E R, Shiyatov S, Mazepa V. Estimation of the mean chronology / / Cook E R, Kairiukstis L A. Methods of Dendrochronology: Applications in
the Environmental Sciences. Dordrecht: Kluwer Academic Publishers, 1990: 123鄄132.
[21] 摇 Briffa K R, Jones P D. Basic chronology statistics and assessment / / Cook E R, Kairiukstis L A. Methods of Dendrochronology: Applications in the
Environmental Sciences. Dordrecht: Kluwer Academic Publishers, 1990: 137鄄152.
[22] 摇 Wigley T M, Briffa K R, Jones P D. On the average value of correlated time series, with applications in dendroclimatology and hydrometeorology.
Journal of Climate and Applied Meteorology, 1984, 23(2): 201鄄213.
[23] 摇 郑景云, 邵雪梅, 郝志新, 葛全胜. 过去 2000年中国气候变化研究. 地理研究, 2010, 29(9): 1561鄄1570.
[24] 摇 PAGES 2k Consortium. Continental鄄scale temperature variability during the past two millennia: Supplementary Information. Nature Geoscience,
2013, 6(5): 339鄄346.
[25] 摇 Yang B, He M H, Melvin T M, Zhao Y, Briffa K R. Climate control on tree growth at the upper and lower Treelines: A case study in the Qilian
Mountains, Tibetan Plateau. PLoS ONE, 2013, 8(7): e69065.
[26] 摇 Ahmed M, Palmer J, Khan N, Wahab M, Fenwick P, Esper J, Cook E. The dendroclimatic potential of conifers from northern Pakistan.
Dendrochronologia, 2011, 29(2): 77鄄88.
2057 摇 生摇 态摇 学摇 报摇 摇 摇 35卷摇