Introduction to treeclim
R Workshop, AmeriDendro 2016
Christian Zang (edited by Andy Bunn)
26-March-2016
Introduction
Basic Usage
treeclim is an R package for the numerical calibration of proxy-climate relationships, developed with a focus on tree-ring data.
The workhorse function of treeclim is dcc, shorthand for “dendroclimatic calibration”. In its most simplest case, a growth-climate analysis can be performed like this:
dc1 <- dcc(muc_spruce, muc_clim)## Running for timespan 1950 - 2007...For this first example, we use data that comes with treeclim: a chronology of Picea abies growing near Munich (Bavaria, Germany) together with corresponding monthly climate data (temperature and precipitation). The default configuration of dccis to perform a bootstrapped response function analysis (Guiot 1991), using monthly climate variables from previous June to current September. The main difference here to other software (e.g., Dendroclim2002, Biondi and Waikul 2004) is that treeclim uses stationary bootstrapping (Politis and Romano 1994) by default, to account for potentially autocorrelated tree-ring data.
The resulting coefficients can be viewed by calling:
dc1## Coefficients (significance flags correspond to p < 0.05):
## id varname month coef significant ci_lower ci_upper
## prec.prev.jun 1 prec jun 0.014 FALSE -0.195 0.241
## prec.prev.jul 2 prec jul -0.032 FALSE -0.250 0.166
## prec.prev.aug 3 prec aug -0.037 FALSE -0.255 0.183
## prec.prev.sep 4 prec sep 0.050 FALSE -0.151 0.223
## prec.prev.oct 5 prec oct 0.086 FALSE -0.124 0.297
## prec.prev.nov 6 prec nov -0.032 FALSE -0.231 0.144
## prec.prev.dec 7 prec dec 0.143 FALSE -0.084 0.341
## prec.curr.jan 8 prec JAN 0.090 FALSE -0.138 0.345
## prec.curr.feb 9 prec FEB 0.113 FALSE -0.077 0.271
## prec.curr.mar 10 prec MAR 0.036 FALSE -0.226 0.325
## prec.curr.apr 11 prec APR -0.130 FALSE -0.347 0.067
## prec.curr.may 12 prec MAY 0.070 FALSE -0.145 0.257
## prec.curr.jun 13 prec JUN 0.042 FALSE -0.186 0.248
## prec.curr.jul 14 prec JUL 0.122 FALSE -0.088 0.379
## prec.curr.aug 15 prec AUG -0.008 FALSE -0.239 0.219
## prec.curr.sep 16 prec SEP 0.076 FALSE -0.177 0.291
## temp.prev.jun 17 temp jun -0.081 FALSE -0.272 0.125
## temp.prev.jul 18 temp jul -0.022 FALSE -0.217 0.165
## temp.prev.aug 19 temp aug -0.112 FALSE -0.323 0.098
## temp.prev.sep 20 temp sep -0.028 FALSE -0.245 0.191
## temp.prev.oct 21 temp oct 0.212 FALSE -0.070 0.423
## temp.prev.nov 22 temp nov -0.109 FALSE -0.306 0.140
## temp.prev.dec 23 temp dec 0.163 FALSE -0.062 0.419
## temp.curr.jan 24 temp JAN -0.070 FALSE -0.280 0.160
## temp.curr.feb 25 temp FEB 0.134 FALSE -0.125 0.363
## temp.curr.mar 26 temp MAR 0.017 FALSE -0.182 0.255
## temp.curr.apr 27 temp APR -0.154 FALSE -0.349 0.010
## temp.curr.may 28 temp MAY 0.027 FALSE -0.172 0.226
## temp.curr.jun 29 temp JUN -0.028 FALSE -0.343 0.297
## temp.curr.jul 30 temp JUL 0.044 FALSE -0.193 0.315
## temp.curr.aug 31 temp AUG -0.086 FALSE -0.320 0.236
## temp.curr.sep 32 temp SEP 0.004 FALSE -0.216 0.255We see a table listing the monthly climatic predictors that have been generated by treeclim with the value of their corresponding coefficients, significance flags and confidence intervals. In this example, none of the climatic predictors turned out to be significant at a level of p < 0.05. (One important notice: the object returned by dcc is actually more complex than the table we get when print it. To access the coefficients, e.g. for writing your own plotting functions, you could use coef(dc1).)
A basic plot is obtained by issueing:
plot(dc1)
We will see later how this plot can be styled to meet your personal preferences.
This was a very quick and basic introduction to the main functionality of treeclim. Everything else beyond that is refining the specification of your analysis. The next section details on the conceptual framework of treeclim and how it reflects in its design and usage.
Caution! Please make sure to have a recent version of treeclim (version >= 1.0.15) installed when you follow along the examples in this document.
Citation
I enjoy writing treeclim and helping users, and Franco Biondi (of Dendroclim2002 fame) enjoys mentoring my work and collaborating closely to set the direction for this software package. Additional to improving citation metrics, getting proper credits for our work through correct citation also helps us to get an idea about who uses treeclim and for what. So, please make sure to cite treeclim correctly when you use it in your analyses:
print(citation("treeclim"), style = "text")## Zang C and Biondi F (2015). "treeclim: an R package for the
## numerical calibration of proxy-climate relationships."
## _Ecography_, *38*(4), pp. 431-436. ISSN 1600-0587, <URL:
## http://doi.org/10.1111/ecog.01335>.The same applies to R in general and dplR, when you use it together with treeclim (which you definetively should, of course).
The Five Elements of Growth-Climate-Relationships
Everyone looking at tree-growth as a function of climate will have to consider at least five elements of the analysis to be carried out:
- Data: which tree-ring and climate data shall be used
- Design: which climatic variables are considered interesting a priori, i.e. based on your assumptions about what is limiting tree growth for the site you are studying
- Scope: which subset of tree-ring and climate data shall be considered for the analysis
- Method: which method to use, i.e. response functions, correlation functions, or seasonal correlations
- Dynamics: shall the analysis be carried out in a static fashion, yielding one set of coefficients, or shall it be done in a dynamic fashion, allowing to look at time-varying response of tree-growth to climate.
The design of treeclim allows to build an analytic strategy considering these five elements. In the following, we will look at each element in some detail.
Data
treeclim expects a tree-ring chronology as input, ideally as produced by the chron function in dplR. For climate data, two general formats are accepted:
- Year-Month format: here, the climate data is formatted as a
data.framewith columnsyear,month, and one column for each climatic variable like so:
## year month temp prec
## 1 1949 1 -0.4 44.9
## 2 1949 2 -0.7 19.2
## 3 1949 3 1.3 33.9
## 4 1949 4 9.8 58.9
## 5 1949 5 11.1 250.5
## 6 1949 6 13.1 70.7It is important to note that treeclim can only work with data for which all years are complete, which means you might have to truncate your climate data if this is not the case.
- 13-Column format: here, for each climate variable a
data.frameor numerical matrix is supplied with year as first column and monthly observations in columns 2 to 13 like so:
## Year Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
## 1 1948 114 52 30 29 21 82 129 45 7 19 18 16
## 2 1949 40 22 44 46 101 24 12 35 35 6 21 54
## 3 1950 43 68 12 70 59 45 76 103 59 38 149 29
## 4 1951 54 45 28 14 44 137 55 34 25 33 60 44
## 5 1952 60 48 132 54 47 61 7 51 112 77 144 73
## 6 1953 20 22 3 38 45 124 116 29 29 26 8 26To use more than one variable with this data format, you can supply them as a list to dcc like so: dcc(chrono = rt_spruce, climate = list(rt_temp, rt_prec))
For both climate data formats, there is no limit on the number of variables that can be used for the analyses. Please note that no missing values are allowed for both formats, so if your climate data contains missing values, consider interpolation or other imputation methods.
For the 13-Column format, there are no obvious variable names for the climatic variables. If you want to have names, you can supply them using the var_names argument like so: dcc(chrono = rt_spruce, climate = list(rt_temp, rt_prec), var_names = c("tmp", "pre")). This also works for renaming variables coming from the Year-Month format.
Design
Being able to design a customized matrix of predictors is probably the most interesting part about dcc. The default design matrix will select all months from previous June to current September from all variables supplied to dcc. Since this is most of the times a good starting point, but not necessarily the ideal design matrix, you can easily specify a custom design matrix using a set of simple commands.
Modifiers
In treeclim, there are three modifier functions for building the design matrix.
.range: for selecting all months individually.mean: for building the mean over the specified months.sum: for summing up over the specified months
All modifiers have the same structure of MODIFIER(.months, .variables), i.e. they take as arguments a number of months, and the variables for which the modifier shall be applied. Both arguments can be omitted; when .months is omitted, the default set of months (previous June to current September) is used, when .variables is omitted, the modifier will be applied to all available climate variables.
Here are two examples:
dc2 <- dcc(muc_spruce, muc_clim, .range(-6:9, c("temp", "prec")))## Running for timespan 1950 - 2007...rownames(dc2$coef)## [1] "prec.prev.jun" "prec.prev.jul" "prec.prev.aug" "prec.prev.sep"
## [5] "prec.prev.oct" "prec.prev.nov" "prec.prev.dec" "prec.curr.jan"
## [9] "prec.curr.feb" "prec.curr.mar" "prec.curr.apr" "prec.curr.may"
## [13] "prec.curr.jun" "prec.curr.jul" "prec.curr.aug" "prec.curr.sep"
## [17] "temp.prev.jun" "temp.prev.jul" "temp.prev.aug" "temp.prev.sep"
## [21] "temp.prev.oct" "temp.prev.nov" "temp.prev.dec" "temp.curr.jan"
## [25] "temp.curr.feb" "temp.curr.mar" "temp.curr.apr" "temp.curr.may"
## [29] "temp.curr.jun" "temp.curr.jul" "temp.curr.aug" "temp.curr.sep"This is identical with the default configuration and can be written more concise as dcc(muc_spruce, muc_clim, -6:9). Another example:
dc3 <- dcc(muc_spruce, muc_clim, .sum(5:9, c("temp", "prec")))## Running for timespan 1949 - 2007...dc3## Coefficients (significance flags correspond to p < 0.05):
## id varname month coef
## prec.sum.curr.may.cu...ul.curr.aug.curr.sep 1 prec.sum MAY...SEP 0.129
## temp.sum.curr.may.cu...ul.curr.aug.curr.sep 2 temp.sum MAY...SEP -0.129
## significant ci_lower ci_upper
## prec.sum.curr.may.cu...ul.curr.aug.curr.sep FALSE -0.098 0.308
## temp.sum.curr.may.cu...ul.curr.aug.curr.sep FALSE -0.310 0.098In this example, we did use all individual months, but as a sum of temperature and precipitation for the given period of time, resulting in much smaller design matrix with only two predictor variables. Since we only have the variables temp and prec in muc_clim, we could have written this more concise as dcc(muc_spruce, muc_clim, .sum(5:9)).
Exclusion
It can be desireable to exclude some time spans from analysis, for example dormant seaons. For this task, treeclim has the function exclude_from (or exfr for short) to exclude specific months from being integrated into the predictors build using the modifier functions.
dc4 <- dcc(muc_spruce, muc_clim, .range(exfr(-6:9, -10:2)))## Running for timespan 1950 - 2007...rownames(dc4$coef)## [1] "prec.prev.jun" "prec.prev.jul" "prec.prev.aug" "prec.prev.sep"
## [5] "prec.curr.mar" "prec.curr.apr" "prec.curr.may" "prec.curr.jun"
## [9] "prec.curr.jul" "prec.curr.aug" "prec.curr.sep" "temp.prev.jun"
## [13] "temp.prev.jul" "temp.prev.aug" "temp.prev.sep" "temp.curr.mar"
## [17] "temp.curr.apr" "temp.curr.may" "temp.curr.jun" "temp.curr.jul"
## [21] "temp.curr.aug" "temp.curr.sep"In this example, we excluded autumn and winter months October (previous year) to February (current year) from the analysis. The remaining months from .range(-6:9) were selected as individual months for building the design matrix, because we used the .range modifier. Often, as here, it is not necessary to explicitly use .range, and dcc(..., exfr(-6:9, -10:2)) would have done the same thing.
Modifier Chaining
The last component for building design matrix is introduced by the ability to chain modifiers. This is done using the + operator:
dc5 <- dcc(muc_spruce, muc_clim, .sum(-4:-9, "prec") +
.sum(4:9, "prec") + .range(exfr(-6:9, -10:3), "temp"))## Running for timespan 1950 - 2007...dc5## Coefficients (significance flags correspond to p < 0.05):
## id varname month coef
## prec.sum.prev.apr.pr...ul.prev.aug.prev.sep 1 prec.sum apr...sep 0.037
## prec.sum.curr.apr.cu...ul.curr.aug.curr.sep 2 prec.sum APR...SEP 0.179
## temp.prev.jun 3 temp jun -0.121
## temp.prev.jul 4 temp jul -0.062
## temp.prev.aug 5 temp aug -0.072
## temp.prev.sep 6 temp sep -0.017
## temp.curr.apr 7 temp APR -0.135
## temp.curr.may 8 temp MAY 0.099
## temp.curr.jun 9 temp JUN -0.005
## temp.curr.jul 10 temp JUL 0.004
## temp.curr.aug 11 temp AUG -0.047
## temp.curr.sep 12 temp SEP 0.033
## significant ci_lower ci_upper
## prec.sum.prev.apr.pr...ul.prev.aug.prev.sep FALSE -0.138 0.208
## prec.sum.curr.apr.cu...ul.curr.aug.curr.sep FALSE -0.028 0.400
## temp.prev.jun FALSE -0.291 0.087
## temp.prev.jul FALSE -0.239 0.111
## temp.prev.aug FALSE -0.281 0.127
## temp.prev.sep FALSE -0.243 0.223
## temp.curr.apr FALSE -0.315 0.037
## temp.curr.may FALSE -0.174 0.343
## temp.curr.jun FALSE -0.288 0.297
## temp.curr.jul FALSE -0.275 0.256
## temp.curr.aug FALSE -0.264 0.294
## temp.curr.sep FALSE -0.186 0.250Here, we are using two precipitation sums (corresponding to the previous and current growing seasons), and a range of month-wise temperatures excluding a dormant season. We can chain as many modifiers as we want. Using modifiers, exclusion and modifier chaining, we can create as complex (or simple) design matrices as needed.
Redundancy Check
While treeclim will check carefully that no two predictors are exactly the same, it can by no means check if the specification of the design matrix makes any sense at all.
Scope
By default, treeclim will use the maximum timespan available for both tree-ring and climate data. This can be changed using the timespan argument: dcc(..., timespan = c(1965, 2000)) will run the analysis using the timespan from 1965 to 2000. When the number of years for the analysis becomes smaller than the number of predictors in the design matrix, treeclim will raise an error.
Method
With dcc, one may choose between response function analysis (the default) and correlation function analysis using the method argument: dcc(..., method = "correlation") will switch to the correlation method. A third method is available: seasonal correlations. The rationale behind this approach developed by Dave Meko (Meko et al. 2011) is somewhat different, which is why it deserves its own function seascorr.
Dynamics
The temporal stability of climate-growth-relations can be analysed using moving response and correlation functions. By default, dcc computes static analyses, but by setting the movingargument to TRUE, the analyses will be repeated in shifting time windows.
dc6 <- dcc(chrono = muc_spruce, climate = muc_clim,
selection = 4:10, moving = TRUE, sb = FALSE)
dc6## Coefficients (significance flags correspond to p < 0.05):
## 1949-1973 1950-1974 1951-1975 1952-1976 1953-1977 1954-1978
## prec.curr.apr "-0.2" "-0.132" "-0.269" "-0.071" "-0.138" "-0.229"
## prec.curr.may "-0.042" "-0.176" "-0.031" "-0.048" "0.029" "0.06"
## prec.curr.jun "0.064" "0.105" "0.14" "0.319" "0.358" "0.305"
## prec.curr.jul "0.151" "0.158" "0.164" "0.033" "0.08" "0.215"
## prec.curr.aug "0.141" "0.106" "-0.065" "-0.085" "-0.088" "-0.161"
## prec.curr.sep "0.102" "0.045" "-0.084" "-0.13" "-0.143" "-0.251"
## prec.curr.oct "-0.073" "-0.048" "0.021" "0.161" "0.147" "0.08"
## temp.curr.apr "-0.197" "-0.184" "-0.213" "-0.191" "-0.197" "-0.198"
## temp.curr.may "0.137" "0.145" "0.071" "0" "-0.015" "0.011"
## temp.curr.jun "0.144" "0.177" "0.229" "0.033" "0.031" "0.027"
## temp.curr.jul "0.15" "0.16" "0.056" "0.006" "-0.04" "-0.008"
## temp.curr.aug "-0.106" "-0.122" "-0.126" "-0.014" "-0.052" "-0.132"
## temp.curr.sep "0.165" "0.12" "-0.016" "0.032" "0.099" "0.116"
## temp.curr.oct "-0.018" "0.059" "0.218" "0.175" "0.165" "0.147"
## 1955-1979 1956-1980 1957-1981 1958-1982 1959-1983 1960-1984
## prec.curr.apr "-0.167" "-0.092" "-0.141" "-0.136" "-0.098" "-0.125"
## prec.curr.may "0.005" "-0.018" "-0.039" "-0.063" "-0.077" "-0.091"
## prec.curr.jun "0.217" "0.21" "0.15" "0.072" "0.04" "0.029"
## prec.curr.jul "0.278" "0.23" "0.359" "0.239" "0.165" "0.144"
## prec.curr.aug "-0.188" "-0.169" "-0.26" "-0.21" "-0.151" "-0.072"
## prec.curr.sep "-0.236" "-0.241" "-0.097" "-0.131" "-0.041" "0.123"
## prec.curr.oct "0.059" "0.096" "0.213" "0.216" "0.181" "0.166"
## temp.curr.apr "-0.112" "-0.087" "-0.054" "-0.089" "-0.021" "-0.122"
## temp.curr.may "0.047" "0.042" "0.047" "0.132" "0.056" "-0.017"
## temp.curr.jun "-0.134" "-0.14" "-0.096" "-0.045" "-0.016" "-0.019"
## temp.curr.jul "-0.227" "-0.251" "-0.22" "-0.18" "-0.027" "-0.018"
## temp.curr.aug "-0.116" "-0.045" "-0.045" "0.007" "0.088" "0.111"
## temp.curr.sep "0.143" "0.137" "0.213" "0.287" "0.289" "0.281"
## temp.curr.oct "0.212" "0.221" "0.242" "0.3" "0.307" "0.396"
## 1961-1985 1962-1986 1963-1987 1964-1988 1965-1989 1966-1990
## prec.curr.apr "-0.117" "-0.096" "-0.132" "-0.185" "-0.116" "-0.082"
## prec.curr.may "-0.033" "-0.025" "0.011" "-0.031" "-0.064" "-0.018"
## prec.curr.jun "0.006" "0.003" "-0.027" "-0.024" "-0.052" "0.001"
## prec.curr.jul "0.118" "0.076" "0.092" "0.158" "0.167" "0.18"
## prec.curr.aug "0.028" "0.001" "-0.074" "-0.07" "-0.121" "-0.134"
## prec.curr.sep "0.085" "0.099" "0.078" "0.007" "0.041" "0.082"
## prec.curr.oct "0.16" "0.15" "0.09" "0.022" "0.071" "0.103"
## temp.curr.apr "-0.16" "-0.124" "-0.087" "-0.009" "0.02" "-0.023"
## temp.curr.may "0.021" "-0.059" "-0.1" "-0.067" "0" "0.034"
## temp.curr.jun "-0.114" "-0.118" "-0.134" "-0.113" "-0.115" "-0.079"
## temp.curr.jul "-0.021" "-0.02" "-0.038" "-0.015" "0.055" "0.064"
## temp.curr.aug "0.134" "0.093" "0.153" "0.179" "0.118" "0.164"
## temp.curr.sep "0.247" "0.295" "0.296" "0.276" "0.319" "0.297"
## temp.curr.oct "0.379" "0.396*" "0.368" "0.392" "0.335" "0.382"
## 1967-1991 1968-1992 1969-1993 1970-1994 1971-1995 1972-1996
## prec.curr.apr "-0.094" "-0.043" "-0.066" "-0.017" "0.055" "0.005"
## prec.curr.may "-0.05" "0.114" "0.09" "0.105" "0.117" "0.131"
## prec.curr.jun "-0.008" "-0.023" "-0.057" "-0.058" "-0.096" "-0.113"
## prec.curr.jul "0.184" "0.161" "0.031" "0.071" "0.065" "0.056"
## prec.curr.aug "-0.106" "-0.104" "-0.122" "-0.128" "-0.087" "-0.09"
## prec.curr.sep "0.092" "0.118" "0.131" "0.139" "0.164" "0.168"
## prec.curr.oct "0.105" "-0.027" "0.024" "0.02" "0.071" "0.049"
## temp.curr.apr "-0.022" "-0.106" "-0.161" "-0.136" "-0.183" "-0.159"
## temp.curr.may "0.058" "-0.01" "-0.057" "-0.093" "-0.107" "-0.091"
## temp.curr.jun "-0.115" "-0.187" "-0.202" "-0.169" "-0.131" "-0.088"
## temp.curr.jul "0.073" "-0.001" "0.022" "0.035" "0.005" "-0.021"
## temp.curr.aug "0.169" "0.028" "-0.024" "0.006" "0.008" "-0.001"
## temp.curr.sep "0.269" "0.183" "0.196" "0.189" "0.191" "0.168"
## temp.curr.oct "0.384" "0.427*" "0.448*" "0.437*" "0.406" "0.39"
## 1973-1997 1974-1998 1975-1999 1976-2000 1977-2001 1978-2002
## prec.curr.apr "0.017" "0.002" "-0.039" "-0.031" "0.004" "0.044"
## prec.curr.may "0.073" "0.073" "0.105" "0.118" "0.124" "0.124"
## prec.curr.jun "-0.111" "-0.121" "-0.122" "-0.093" "-0.204" "-0.199"
## prec.curr.jul "0.077" "0.071" "0.068" "0.039" "0.047" "0.085"
## prec.curr.aug "-0.15" "-0.149" "-0.177" "0.002" "-0.053" "0"
## prec.curr.sep "0.015" "0.014" "0.041" "0.125" "0.182" "0.22"
## prec.curr.oct "0.083" "0.053" "0.041" "0.052" "-0.003" "-0.021"
## temp.curr.apr "-0.151" "-0.162" "-0.157" "-0.142" "-0.158" "-0.229"
## temp.curr.may "-0.144" "-0.142" "-0.119" "-0.078" "-0.071" "-0.079"
## temp.curr.jun "-0.079" "-0.116" "-0.152" "-0.16" "-0.027" "0.07"
## temp.curr.jul "-0.061" "-0.06" "-0.074" "-0.146" "-0.087" "-0.065"
## temp.curr.aug "-0.025" "-0.03" "0.027" "0.03" "-0.129" "-0.134"
## temp.curr.sep "0.176" "0.172" "0.113" "0.205" "0.19" "0.116"
## temp.curr.oct "0.323" "0.282" "0.167" "0.121" "0.095" "0.078"
## 1979-2003 1980-2004 1981-2005 1982-2006 1983-2007
## prec.curr.apr "0.096" "0.1" "0.068" "0.054" "0.154"
## prec.curr.may "0.106" "0.12" "0.111" "0.107" "0.17"
## prec.curr.jun "-0.221" "-0.115" "-0.093" "-0.006" "-0.004"
## prec.curr.jul "0.094" "0.083" "0.073" "0.137" "0.114"
## prec.curr.aug "0.096" "0.139" "0.123" "0.142" "0.096"
## prec.curr.sep "0.203" "0.246" "0.266" "0.248" "0.291"
## prec.curr.oct "0.005" "0.002" "0.029" "-0.12" "-0.088"
## temp.curr.apr "-0.19" "-0.275" "-0.32" "-0.29" "-0.249"
## temp.curr.may "-0.076" "0.018" "0.002" "0.021" "0.045"
## temp.curr.jun "-0.147" "-0.107" "-0.122" "-0.139" "-0.182"
## temp.curr.jul "-0.083" "-0.057" "-0.074" "-0.148" "-0.187"
## temp.curr.aug "-0.266" "-0.353" "-0.311" "-0.148" "-0.13"
## temp.curr.sep "0.109" "0.088" "0.085" "-0.046" "-0.06"
## temp.curr.oct "0.116" "0.011" "0.012" "-0.086" "-0.087"plot(dc6)
The stars in the plot indicate windows with significant correlations for the given variable. Like in the non-dynamic case, you are free to create a customized design matrix using the treeclim modifiers. Further options comprise the window size, the window offset (the gap between consecutive window starts), and whether to start with windowing from the most recent end (default), or from the least recent end.
If method is set to correlation for computing a moving correlation function, we can use g_test to test the temporal fluctuations of the correlation coefficients for significance by comparing them to the magnitude of fluctuations that would occur by chance alone (using a multivariate extension of the approach by Gershunov et al. 2001):
dc7 <- dcc(chrono = muc_spruce, climate = muc_clim,
selection = .mean(3:5) + .mean(6:8) + .mean(9:10), method = "cor",
moving = TRUE, win_size = 30, sb = FALSE)
plot(dc7)
g_test(dc7, sb = FALSE)## Performing test without bootstrapping individual correlations.## varname month p
## 1 temp.mean SEP...OCT 0.104
## 2 temp.mean JUN...AUG 0.028
## 3 temp.mean MAR...MAY 0.228
## 4 prec.mean SEP...OCT 0.368
## 5 prec.mean JUN...AUG 0.390
## 6 prec.mean MAR...MAY 0.082In this example, only the fluctuations in the correlation coefficient for mean temperature from current June to current August are significant at a level of p < 0.05. This means, that the probability of an even larger fluctuation by chance is less than 0.05. For more details, look into ?g_test.
Seasonal Correlations
treeclim provides similar functionality to Dave Meko’s seascorr MATLAB program: partial correlations for a primary and a secondary climatic variable for seasons of different lengths. The tree-ring and climate input data is required in the same format as for dcc.
sc1 <- seascorr(muc_spruce, muc_clim)## Running for timespan 1949 - 2007...sc1## Results for a season length of 1 month:
## month type coef significant
## 1 aug primary -0.201353 FALSE
## 2 sep primary -0.046407 FALSE
## 3 oct primary 0.289829 TRUE
## 4 nov primary -0.229822 FALSE
## 5 dec primary 0.213082 FALSE
## 6 JAN primary -0.058630 FALSE
## 7 FEB primary 0.194547 FALSE
## 8 MAR primary 0.096562 FALSE
## 9 APR primary -0.216712 FALSE
## 10 MAY primary -0.002871 FALSE
## 11 JUN primary -0.081163 FALSE
## 12 JUL primary -0.040511 FALSE
## 13 AUG primary -0.122883 FALSE
## 14 SEP primary -0.007321 FALSE
## 15 aug secondary -0.121383 FALSE
## 16 sep secondary 0.143349 FALSE
## 17 oct secondary 0.099696 FALSE
## 18 nov secondary -0.110715 FALSE
## 19 dec secondary 0.160683 FALSE
## 20 JAN secondary 0.122742 FALSE
## 21 FEB secondary 0.141762 FALSE
## 22 MAR secondary 0.060445 FALSE
## 23 APR secondary -0.142786 FALSE
## 24 MAY secondary 0.074795 FALSE
## 25 JUN secondary 0.033870 FALSE
## 26 JUL secondary 0.158503 FALSE
## 27 AUG secondary 0.019846 FALSE
## 28 SEP secondary 0.220631 FALSE
##
## Results for a season length of 3 months:
## month type coef significant
## 1 aug primary -0.26670 FALSE
## 2 sep primary -0.21116 FALSE
## 3 oct primary 0.04095 FALSE
## 4 nov primary 0.01417 FALSE
## 5 dec primary 0.18663 FALSE
## 6 JAN primary -0.02033 FALSE
## 7 FEB primary 0.17070 FALSE
## 8 MAR primary 0.11878 FALSE
## 9 APR primary 0.09937 FALSE
## 10 MAY primary -0.03333 FALSE
## 11 JUN primary -0.14034 FALSE
## 12 JUL primary -0.05457 FALSE
## 13 AUG primary -0.10413 FALSE
## 14 SEP primary -0.08277 FALSE
## 15 aug secondary -0.08814 FALSE
## 16 sep secondary -0.02398 FALSE
## 17 oct secondary 0.08810 FALSE
## 18 nov secondary 0.05480 FALSE
## 19 dec secondary 0.07216 FALSE
## 20 JAN secondary 0.11419 FALSE
## 21 FEB secondary 0.30820 TRUE
## 22 MAR secondary 0.21256 FALSE
## 23 APR secondary 0.08474 FALSE
## 24 MAY secondary 0.04763 FALSE
## 25 JUN secondary 0.01954 FALSE
## 26 JUL secondary 0.17789 FALSE
## 27 AUG secondary 0.15294 FALSE
## 28 SEP secondary 0.23795 FALSE
##
## Results for a season length of 6 months:
## month type coef significant
## 1 aug primary -0.207292 FALSE
## 2 sep primary -0.218888 FALSE
## 3 oct primary -0.095932 FALSE
## 4 nov primary -0.185217 FALSE
## 5 dec primary -0.032638 FALSE
## 6 JAN primary 0.005526 FALSE
## 7 FEB primary 0.149260 FALSE
## 8 MAR primary 0.172287 FALSE
## 9 APR primary 0.057038 FALSE
## 10 MAY primary 0.100142 FALSE
## 11 JUN primary 0.030095 FALSE
## 12 JUL primary 0.043238 FALSE
## 13 AUG primary -0.086953 FALSE
## 14 SEP primary -0.134976 FALSE
## 15 aug secondary 0.022943 FALSE
## 16 sep secondary 0.036155 FALSE
## 17 oct secondary 0.067985 FALSE
## 18 nov secondary 0.026315 FALSE
## 19 dec secondary 0.106765 FALSE
## 20 JAN secondary 0.125854 FALSE
## 21 FEB secondary 0.226750 FALSE
## 22 MAR secondary 0.206780 FALSE
## 23 APR secondary 0.163925 FALSE
## 24 MAY secondary 0.229434 FALSE
## 25 JUN secondary 0.167334 FALSE
## 26 JUL secondary 0.208102 FALSE
## 27 AUG secondary 0.153896 FALSE
## 28 SEP secondary 0.199985 FALSEThis is the default configuration for seascorr: considering three different season lengths (1, 3, and 6 months), and current September as month where tree-ring growth is assumed to have finished. There is a plot method available as well:
plot(sc1)
By default, the primary variable is the first one found in the climate data supplied to seascorr. So with year, month, temp, prec in muc_clim, precipitation is chosen as primary variable. You can change this by either supplying names or positions, so to swap primary and secondary variables from the last example we can issue:
sc2 <- seascorr(muc_spruce, muc_clim, primary = "prec", secondary = "temp")## Running for timespan 1949 - 2007...plot(sc2)
Reconstruction Skills
treeclim offers some basic functionality for evaluating the reconstruction skills of tree-ring chronologies via the skills function. skills works on existing output from either dcc or seascorr and tests the prediction skills of the tree-ring data for a given climate target. For this example, we use a different built-in data set, a chronology of pines from Visdalen in Norway (consult ?norw015 for details) and associated climate data:
dc8 <- dcc(norw015, list(norway_prec, norway_temp), 3:9, var_names = c("prec", "temp"))## Running for timespan 1901 - 1983...plot(dc8)
Let us assume we are interested in reconstructing summer temperatures, and want to make use of the significant correlation with July temperature. We call skills on the dc8 object, and specify .range(7, "temp") as target to get July temperature as response variable (note that for skills, like in many dendroclimatic applications, climate is modelled subject to tree-ring data):
sk8 <- skills(dc8, .range(7, "temp"))
sk8## Call:
## skills(object = dc8, target = .range(7, "temp"))
##
## Using climate target temp.curr.jul and calibration model with 50 percent (= 42 years) of data starting at recent end as calibration period:
##
## r p-value
## 0.543895 0.009901
##
## Coefficients:
## intercept slope
## 8.572 3.238
##
## Verification statistics:
## Red. of Error (RE) Coeff. of Efficiency (CE)
## 0.032 0.004
##
## Durbin-Watson Test (DW): 1.628 (p = 0.09402)
##
## Model for whole period:
##
## r p-value
## 0.445744 0.009901
##
## Coefficients:
## intercept slope
## 9.157 2.968By default, skills will take the maximum available timespan and divide it into equally sized calibration and verification periods. Reported are the parameter estimates and goodness of fit for the calibration and full models, plus some established verification statistics. Here, both RE and CE are positive, so the reconstruction has some skills at least. We can also look at the reconstruction:
plot(sk8)
The specification of the calibration period in skills is quite flexible, you might want to look in to ?skills for more details. Furthermore, you can chose between an OLS regression model (default) or an RMA regression.
Plotting
treeclim uses ggplot2 (Wickham 2009) to build its plots. ggplot2-plots are themeable and easily be customized. This means, it only takes a few lines of code to adapt your treeclim output to your aesthetical preferences.
plot(dc1)
This is the default. If we feel the default style does not blend in nicely with our gothic inspired interior design, we can always change some things (we need to include ggplot2 for this):
library(ggplot2)
plot(dc1) + theme_dark() + scale_color_manual(values = c("#FFFFFF", "#F62525")) +
ylab("Response Coefficients") + ggtitle("Happy Worlddendro 2016!")
You may want to consult the ggplot2 documentation for further options (there are many).
Conclusion & Acknowledgements
treeclim tries to provide a flexible framework for the exploration of tree-climate-relationships. treeclim does not impose complexity on the user, and keeps simple things simple. When needed, however, functionality is available to adapt to various research questions. Most aspects of what treeclim can do has been touched in this vignette. There are few more options that are potentially interesting, such as the ability to use different bootstrapping schemes with dcc.
Areas where treeclim is currently lacking include better support of southern hemisphere data by providing access to climate data from the year before the previous one. Also, testing reconstructing skills is only basic.
Some of the features currently implemented derive from requests by the user community (seascorr being the most prominent example), and I am grateful for these and further suggestions.
Citations
Biondi, F., and K. Waikul. 2004. DENDROCLIM2002: A C++ program for statistical calibration of climate signals in tree-ring chronologies. Computers in Geosciences 30:303–311.
Gershunov, A., N. Schneider, and T. Barnett. 2001. Low-frequency modulation of the ENSO-Indian Monsoon rainfall relationship: Signal or noise? Journal of Climate 14:2486–2492.
Guiot, J. 1991. The bootstrapped response function. Tree Ring Bulletin 51:39–41.F
Meko, D. M., R. Touchan, and K. J. Anchukaitis. 2011. Seascorr: A MATLAB program for identifying the seasonal climate signal in an annual tree-ring time series. Computers & Geosciences 37:1234–1241.
Politis, D. N., and J. P. Romano. 1994. The stationary bootstrap. Journal of the American Statistical Association 89:1303–1313.
H. Wickham. ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag New York, 2009.