Trees do not always act their age: size-deterministic tree ring standardization for long-term trend estimation in shade-tolerant trees

With increasing awareness of the consequences of climate change for global ecosystems, the focus and application of tree ring research have shifted to reconstruction of long-term climate-related trends in tree growth. Contemporary methods for estimating and removing biological growth trends from tree ring series (standardization) are ill-adapted to shade-tolerant species, leading to biases in the resultant chronologies. Further, many methods, including regional curve standardization (RCS), encounter significant limitations for species in which accurate age estimation is difficult. In this study we present and test two tree ring standardization models that integrate tree size in the year of ring formation into the estimation of the biological growth trend. The first method, dubbed size-deterministic standardization (SDS), uses tree diameter as the sole predictor of the growth trend. The second method includes the combined (COMB) effects of age and diameter. We show that both the SDS and COMB methods reproduce long-term trends in simulated tree ring data better than conventional methods; this result is consistent across multiple species. Further, when applied to real tree ring data, the SDS and COMB models reproduce long-term, time-related trends as reliably as traditional RCS and more reliably than other common standardization methods (i.e. C-method, basal area increments, conservative detrending). We recommend the inclusion of tree size in the year of ring formation in future tree ring standardization models, particularly when dealing with shade-tolerant species, as it does not compromise model accuracy and allows for the inclusion of unaged trees.

1 Introduction

Tree rings have long-served as a record of environmental change in forest ecosystems. Early dendrochronological studies used tree ring chronologies from climate-sensitive species to elucidate the dynamics of growth–climate relationships and reconstruct climate anomalies from periods before the existence of instrumental records. However, with increasing awareness of the consequences of climate change for global ecosystems, the focus and application of tree ring research have shifted to reconstruction of low-frequency climate-related trends in tree growth (Gedalof and Berg, 2010; Boisvenue and Running, 2006; Jacoby and D'Arrigo, 1997). As it stands, previous optimism regarding the benefits of carbon fertilization for forest growth (Battipaglia et al., 2013; Norby et al., 2005) has been quelled by a lack of consistent evidence in real forests. While many studies have noted increases in long-term growth rates over time in temperate forests (Gedalof and Berg, 2010; Huang et al., 2007; Martinelli, 2004), others suggest no change (Giguère-Croteau et al., 2019; Camarero et al., 2015; Granda et al., 2014; Silva et al., 2010; Peñuelas et al., 2011). Further, in boreal and drought-prone species, growth decline (Chen et al., 2018; Dietrich et al., 2016; Girardin et al., 2011; Silva and Anand, 2013) and increased mortality (Herguido et al., 2016; Liang et al., 2016) in response to climate stress have been prevalent. Central to all of these studies is the assumption that long-term growth trends can be accurately and unbiasedly estimated from tree ring data.

As it stands, accurate estimation of long-term growth trends in forests may be limited by poorly adapted tree ring standardization (age trend removal) methods (Briffa et al., 1996) and inappropriate sampling methods (Nehrbass-Ahles et al., 2014; Brienen et al., 2012). Early standardization methods (i.e. conservative detrending) were designed to maintain high-frequency variation in tree ring series and discard long-term, low-frequency variation. It is accepted that these methods are inappropriate for estimating long-term, climate-related growth trends (Briffa et al., 1992); however, they are still used in situations where contemporary standardization methods are not applicable due to restrictive data requirements (e.g. Villalba et al., 2012; Gedalof and Berg, 2010; Wang et al., 2006).

Modern standardization methods are designed to estimate biological age- and size-related effects on tree growth independent of time-related variance, thus theoretically maintaining long-term trends in the final chronologies. Among these, the conversion of tree ring widths to basal area increments (BAIs), and the closely related C-method (Biondi and Qeadan, 2008), as well as the use of regional curve standardization (RCS; Briffa et al., 1992), and its many variants (see Helama et al., 2017), have become commonplace (Peters et al., 2015). Traditional RCS relies on the assumption that the species-specific biological growth trend of local trees can be estimated, and thus removed, from a sufficiently large sample of trees using tree age alone. Alternatively, the BAI method assumes that the biological growth trend is sufficiently related to basal area accrued in a given year and, as such, chronologies presented as BAI (instead of raw ring width) contain minimal biological effects. In practice, it is unlikely that this strict relationship accounts for all the variation in ring width that is related to biological size and age effects. As such, some studies have proposed explicit models of BAI that attempt to include variables related to tree age and size or environmental conditions (i.e. tree density, soil fertility, etc.; see, e.g. Linares et al., 2009; Nock et al., 2011). Similarly, the C-method (CM) assumes that tree-wise BAI (tree ring area) distributed over a growing surface in time is constant and as such, annual deviations from this trend can represent the standardized chronology (free from biological trend, Biondi and Qeadan, 2008). Both BAI and CM are best suited to open-growth, shade-intolerant trees, where the strict relationship between annual growth and expected BAI is not impeded by early competition for light.

However, due to the difficulties in separating climate-related trends that vary on long timescales from those related to biological tree growth and/or succession-related environmental change, none of these methods are likely to produce accurate estimates of external forcing when trees from only a single age or size class are sampled (Brienen et al., 2012; Briffa and Melvin, 2011). It follows that studies which only sample even-aged stands or dominant trees are likely to produce biased estimates of long-term growth. While increased awareness of sample biases has led to better prescriptions for study design (see Nehrbass-Ahles et al., 2014; Brienen et al., 2012), systematic tests of the ability of these models to accurately reproduce long-term trends are still limited (e.g. Sullivan et al., 2016; Peters et al., 2015; Esper et al., 2003).

Despite these limitations, RCS remains the standard method for estimating long-term growth trends in tree ring data (Helama et al., 2017). However, the standard RCS approach encounters large limitations for many species in which accurate age estimation is difficult. Additionally, we suggest the inherent assumption of RCS that biological growth trends are sufficiently determined by tree age may not be appropriate in all species. More specifically, this assumption is problematic for shade-tolerant trees. Shade-tolerant species exhibit relatively low low-light mortality and thus can persist in forest understories for variable amounts of time before release from overstory light suppression. In these cases, traditional age-deterministic models exhibit high variance, and thus low precision, in the period following tree establishment and leading up to the age when most trees have been released from suppression (Fig. 1). This period of ill-fit means that trees that are released relatively early (or late) from light suppression will exhibit inflated (or deflated) growth relative to the chronology. As a result, the final chronology will show less agreement than would be expected in a shade-intolerant species. Even more problematically, if trees are sampled according to minimum size thresholds, the youngest trees in the chronology are likely to be early-release trees, leading to an artificial inflation of modern growth rates in the final chronology. While modifications to traditional RCS that address variance in contemporaneous growth rates and regional environmental conditions have been prevalent in shade-intolerant species (see Helama et al., 2017) there has been little to no focus on the improvement of standardization techniques specific to shade-tolerant tree species.

Figure 1 (a) In shade-tolerant species, young trees are stochastically released from low-light suppression in the understory. (b) Since release from suppression is not strictly related to tree age, widely used communal age trend models (RCS) poorly model tree growth in the period following establishment and leading up to the age when most trees have been released from suppression. (c) Poor model fit in this period implies that the biological growth trend is not entirely removed from individual series and leads to high residual variance when standardized tree ring series are aligned according to calendar year.

Alternatively, in the field of forest growth and yield modelling, size-deterministic, rather than age-deterministic, predictive growth models are ubiquitous. It is well-understood that tree size regulates the capacity for resource acquisition, namely, light (Canham et al., 2004), water and nutrients (Homann et al., 2000), resource allocation (Lehnebach et al., 2018), and metabolic costs (West et al., 2001). As such, the notion of radial growth being deterministic according to size rather than age is logical from both a physiological and ecological perspective. Tree size in a given year is dependent on its previous size and annual growth, so shade-tolerant trees that have yet to be released from overstory light suppression remain small as they grow older. This relaxes the period of “ill-fit” that would be observed in an age-based model. Accordingly, we propose that a size-deterministic model for tree ring standardization may be more appropriate than traditional RCS for shade-tolerant tree species. The application of size-deterministic models has been limited, with few examples of tree size in a given year being incorporated into BAI models (e.g. Marqués et al., 2016; Camarero et al., 2015; Nock et al., 2011; Martínez-Vilalta et al., 2008) and even fewer of uniquely size-based tree ring models (e.g. Bontemps and Esper, 2011; Gavin et al., 2008). Further, there have been no systematic evaluations of the ability of size-based models to accurately estimate long-term trends in tree ring series.

We present two tree ring standardization models that integrate tree size in the year of ring formation into estimation of the biological growth trend. The first model uses tree diameter as the sole predictor of the communal growth trend while the second includes the combined effects of both age and diameter. It follows that the objective of this study is to determine the efficacy of both models in estimating long-term growth trends in their resultant tree ring chronologies. First, we use modelled tree ring data from shade-tolerant and shade-intolerant species to make the inappropriateness of age-based models explicit for shade-tolerant trees. Further, we investigate the performance of size-based models relative to contemporary standardization methods in the presence of size thresholds in tree sampling. Last, we apply the developed models to tree ring data from shade-tolerant temperate species to evaluate model performance relative to contemporary methods, on the basis of accurate reconstruction of known long-term time-related trends in the series.

2 Methods

2.1 Model formulation

Traditional RCS makes two assumptions about tree growth: first, that trees of the same species in a given region exhibit a common growth trend as they age, and, second, that growth of an individual tree in a given year is thus a product of its age and common climatic or environmental forcing in that year (Esper et al., 2003; Briffa et al., 1992). We present a variant of the RCS method that uses tree size, measured by diameter at breast height (DBH), in the year of ring formation as the primary determinant of the common biological growth trend. As with RCS we assume that the relationship between expected growth and tree size is non-linear and can be approximated for a region from a sufficiently large sample of trees from the species in question. Further, we assume that using a sample of trees from a range of size and age classes ensures estimation of the common trend is not confounded by underlying low-frequency climate or environmental forcing in the chronology (Brienen et al., 2012). The size-based regional curve model, hereafter referred to as the size-deterministic standardization (SDS) model, takes the following form:

$$\begin{array}{}\text{(1)}& E\left({\mathrm{RW}}_{y,i}\right)=B_o+f_{\mathrm{1}}\left({\mathrm{DBH}}_{y,i}\right)+e_{yi},\end{array}$$

where E(RW_yi) represents the expected ring width of a given tree (i) in year (y), and f₁ represents a non-linear function relating DBH of a given tree (i) in year (y) to E(RW_yi). As in RCS, the communal non-linear relationship is estimated communally for all local trees of interest. In our study we estimate f₁ with a penalized thin plate regression spline in a generalized additive model (GAM); however, this relationship could be estimated by a number of different spline-fitting or non-linear regression techniques (i.e. ffcsaps function in dplR, Bunn et al., 2018, time-varying splines, Melvin et al., 2007). Under this paradigm the model residuals (e_yi) represent individual standardized ring width indices and, by extension, individual tree response to climatic or environmental forcing. Annual model residuals subject to a robust mean thus represent the final standardized chronology. This approach differs slightly from traditional RCS, whereby standardized ring width indices are occasionally produced by division of raw measurements by the expected value. Calculation of standardized ring width indices by subtraction from the expected value, as in the case of residuals, is now commonly used, as it tends to reduce bias in the resultant chronology (Helama et al., 2004) and eases in the formulation of more complex tree ring standardization models. However, unlike division methods, the subtraction method does not provide any stabilization of variance in the resulting residuals; as such, it may be necessary to use a stabilization procedure (i.e. log transformation, power transformation) on raw ring width data beforehand.

Tree size in a given year can be estimated by outside-in or inside-out techniques. If the pith of a tree is present in the core (or reasonably close to) DBH_y is a simple summation of all previous ring widths since the year of origin, multiplied by two. Alternatively, if the pith is missed, DBH_y can be calculated via subtraction of more modern ring widths (multiplied by two) from the inside-bark diameter. In this case, inside-bark diameter is calculated as measured DBH minus bark thickness (multiplied by two), where bark thickness can be directly measured or estimated using species-specific allometric equations (e.g. Stayton and Hoffman, 1970).

Similar to the model formulation for SDS, RCS models were estimated with GAMs of the following form:

$$\begin{array}{}\text{(2)}& E\left({\mathrm{RW}}_{yi}\right)=B_o+f_{\mathrm{1}}\left({\mathrm{Age}}_{yi}\right)+e_{yi},\end{array}$$

where Age_yi is the age of an individual tree in a given year and the resultant standardized tree ring indices are derived from model residuals (e_yi).

In addition, a more complex model that integrated independent size and age effects was also evaluated for comparison. This model, hereafter referred to as the combined model (COMB), took the following form:

$$\begin{array}{}\text{(3)}& E\left({\mathrm{RW}}_{yi}\right)=B_o+f_{\mathrm{1}}\left({\mathrm{Age}}_{y{i}_{\mathrm{a}}}\right)+f_{\mathrm{2}}\left({\mathrm{DBH}}_{yi}\right)+e_{iy}.\end{array}$$

In a large variety of long-lived tree species, accurate age estimation (pith sampling) is difficult or impossible; rendering traditional RCS or combined models inappropriate for all trees sampled. To address this issue, the above model can incorporate unaged trees. Here f₁ represents the non-linear function relating age to expected ring width for the subset of all trees that are aged (i_a). In this model, ring widths from unaged trees are assigned arbitrary ages that do not contribute to the linear approximation of the smooth term for age (i.e. $f_{\mathrm{1}}\left({\mathrm{Age}}_{y{i}_{\mathrm{a}}}\right)$) but these trees still contribute to the smooth term for size f₂(DBH_yi). Syntax for missing data in GAMs follows the protocol provided in mgcv (Wood, 2011). In this study all GAMs were fit using the mgcv package (Wood, 2011) in the R statistical program (v.3.5.0).

In addition to the models presented above we investigated three additional standardization methods: conservative detrending (CD), CM, and BAI. Conservative detrending describes functions (i.e. negative exponentials, straight lines) or flexible splines fit to individual tree ring series (see Cook and Kairiukstis, 2013). In this study we use spline-fitting techniques rather than modified negative exponentials as they are more appropriate for shade-tolerant tree species. As above, the individual standardized tree ring width indices are derived from model residuals. The C-method estimates tree-specific expected ring widths by assuming constant annual basal area increment (tree ring area) over the life span of the tree (see Biondi and Qeadan, 2008). Annual deviations from expected values thus represent standardized ring width indices. For consistency, the standard CM approach in dplR (Bunn et al., 2018) was modified in order to calculate indices via subtraction (residuals) instead of division. Tree ring widths were converted to BAI using the dplR package in R (Bunn et al., 2018). R code for worked examples of all standardization procedures used in this study is available (https://doi.org/10.5281/zenodo.3565813).

2.2 Simulated tree ring data

We simulated tree ring data using a well-established gap-phase model. The SORTIE-ND model was chosen over other similar gap-phase models as it better emulates understory light conditions and low-light mortality, both of which are central to the notion of age being an inappropriate determinant of growth in shade-tolerant species. In SORTIE, annual radial tree growth is calculated as an asymptotic function of light availability and previous tree diameter. As such, the underlying growth trend in SORTIE-simulated data should be well-approximated by a flexible curve estimated on the basis of tree size (SDS). As such, we use this analysis solely to elucidate the problematic nature of age-based standardization methods for shade-tolerant species and not to confirm the efficacy of size-based standardization methods.

For simplicity, a stand 100 % dominated by sugar maple (Acer saccharum) was simulated, as sugar maple is a model shade-tolerant species that grows in self-replacing stands. All living trees (>5 cm DBH; n=3657) in the final year of the model run were used for further analysis. Additionally, to elucidate our claim that age-deterministic growth estimation is more problematic in shade-tolerant species, we completed a similar SORTIE simulation for the shade-intolerant species white pine (Pinus strobus). Again, the stand was 100 % white pine, standard model parameters were used, and the simulation was run for 1000 years. All living trees (>5 cm DBH; n=7362) in the final year of the model run were used for further analysis. Additional details regarding model parameters for the SORTIE simulations are provided in the Supplement (Sect. S1).

To simulate a low-frequency climate-related growth trend, a logistic trend was added to raw tree ring width of individual trees produced by both SORTIE simulations. The logistic trend simulated an initial rapid increase in growth and subsequent levelling off that aimed to represent a period of carbon fertilization and eventual acclimation. The logistic model was applied to the last 100 years of growth and took the following form, where RWt_yi represents ring widths with the simulated long-term trend and RWr_yi are raw ring widths:

$$\begin{array}{}\text{(4)}& {\mathrm{RWt}}_{yi}={\mathrm{RWr}}_{yi}\left({\displaystyle \frac{k}{\mathrm{1}+a{e}^{-ry}}}+\mathrm{1}\right).\end{array}$$

The logistic trend parameters (r=0.12, k=0.629, a=20) were chosen such that increases in individual tree growth averaged approximately 5 % per decade. Additionally, we tested the standardization models in their ability to detect simulated negative trends in tree growth, as previous studies have noted a failure of contemporary methods to accurately reproduce declining growth trends (Peters et al., 2015). The simulated negative logistic trend took the form of Eq. (4) with parameters (r=0.12, $k=-\mathrm{0.421}$, a=20) chosen such that decreases in growth averaged 5 % per decade. For completeness, we also simulated positive and negative linear trends. Results of those analyses are provided in the Supplement (Sect. S3).

A total of 60 trees were randomly selected, without replacement, from the simulated tree populations and subject to each of the six standardization methods (SDS, RCS, COMB, CD, BAI, CM). Model residuals (in the case of RCS, SDS, COMB, CD and CM) or transformed (BAI) tree ring widths were compiled into an annual mean chronology using Tukey's biweight robust mean. The resultant chronologies were then tested for significant correlation with the imposed trends using Spearman's rank correlation coefficient. This process was bootstrap resampled (with replacement) 100 times, in order to produce confidence intervals for the resultant mean chronologies and their respective correlation coefficients.

To examine the effect of minimum size sampling thresholds on the accuracy of long-term trend reconstruction by each of the standardization methods, we completed the same analysis on trees from the simulated populations that exceeded certain size thresholds. The thresholds employed were 10 cm DBH, which represented a practical minimum size threshold for sampling, and 30 and 50 cm DBH which represented thresholds for mature and dominant trees, respectively. The CD method was only applied when size thresholds exceeded 10 cm DBH due to the troublesome nature of fitting splines to excessively short time series. The mean Spearman's rho for all detrending methods and sampling thresholds were compared using two-way ANOVA and post hoc tests.

2.3 Real tree ring data

Additionally, we evaluated the performance of the six standardization methods in real tree ring data from shade-tolerant species. We collected tree ring data from seven mature sugar maple-dominated stands in Ontario, Canada (Table 1). Further, tree ring data sets from the shade-tolerant species red spruce (Picea rubens) were obtained from the DendroEcological Network database (https://www.uvm.edu/femc/dendro, last access: 6 December 2019; Table 1). Red spruce was chosen as it had sufficient replication across studies in the database. Descriptions of the sampling strategies and data-processing methods for all sites considered are provided in either the Supplement (Sect. S3) or in their respective references (i.e. Kosiba et al., 2013, 2017). Data were considered suitable for this study if age and DBH estimates were provided and if a minimum 10 trees per site and species were sampled and accurately aged. All cores in which pith offset was estimated to be greater than 10 years were considered unaged. To simplify comparisons of the resultant chronologies, unaged trees were not included in the models.

Prior to model application, a time-deterministic thin plate regression spline was applied to all raw ring widths from each site. This ensured there was no underlying time trend present in the data. Since trees of multiple ages/sizes were sampled in each study we assume the removed time trend is therefore independent of biological trends in the series. For each site, residuals from the regression spline were centred according to the site-wise mean and standard deviation of raw ring widths prior to analysis.

Table 1 Location, sample size, chronology length, and source of tree ring data sets used in this study.

Again, increasing and decreasing logistic trends (Eq. 4), as well as linear trends (Supplement, Sect. S3), were added to the (recentred) tree ring residuals. Trend parameters were chosen such that the increase (or decrease) in tree growth averaged 5 % per decade over the last 50 years of growth (r=0.12, k=0.276, a=20, positive trend; r=0.12, $k=-\mathrm{0.226}$, a=20, negative trend). For each site all trees were subject to each of the six standardization methods (SDS, RCS, COMB, CD, BAI, and CM). Model residuals (in the case of RCS, SDS, COMB, CD, and CM) or transformed (BAI) tree ring widths were compiled into an annual mean chronology using Tukey's biweight robust mean. The resultant chronologies were then tested for significant correlation with the imposed trends using Spearman's rank correlation coefficient. In both species (sugar maple and red spruce) one-way ANOVA and Tukey post hoc comparisons were used to test for significant differences in model performance, as estimated by chronology correlation with the imposed trend.

3 Results

3.1 Comparisons of methods in simulated data

In order to evaluate the efficacy of each standardization method we calculated correlations between chronologies produced by each method and a variety of imposed trends in simulated sugar maple and white pine tree ring data. Bootstrapped confidence intervals for chronologies from each of the standardization methods are provided in Fig. 2a and b for sugar maple and red pine, respectively. Distributions of the respective Spearman's rank correlation coefficients between the chronologies and the imposed trends are provided in Fig. 3a for sugar maple and Fig. 3b for white pine.

Figure 2 The 95 % confidence intervals for standardized chronologies produced by each standardization method (legend on the right side) applied to SORTIE-simulated (a) sugar maple and (b) white pine tree ring data. Confidence intervals obtained via bootstrap resampling (rep = 100) of 60 trees (>10 cm DBH) from the SORTIE-simulated populations. Dotted lines indicate the standardized positive (left side) or negative (right side) logistic trend that was added to the raw tree ring data.

Figure 3 Spearman's rho correlation between chronologies produced by each of the five standardization methods and the imposed positive (left column) or negative (right column) logistic trend in SORTIE-simulated (a) sugar maple and (b) white pine tree ring data. Correlation distribution created by bootstrap resampling 60 trees (rep = 100) from SORTIE-simulated tree populations. The horizontal axis denotes minimum tree size (DBH) thresholds for sampling from the population. Horizontal lines indicate threshold for significant Spearman's rho (a=0.05) for correlation between chronologies and the imposed trend.

3.2 Simulated sugar maple tree ring data

In the simulated sugar maple data, two-way ANOVA suggested a significant effect of both the standardization model (p<0.001) and minimum size sampling threshold (p<0.001) on average correlation with the positive logistic trend. Alternatively, for the negative logistic trend there was a significant effect of the standardization model (p<0.001) but not of the size sampling threshold. For both positive and negative logistic trends SDS ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.974}\pm \mathrm{0.037}$, $\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.954}\pm \mathrm{0.068}$, respectively) and COMB ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.965}\pm \mathrm{0.039}$, $\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.894}\pm \mathrm{0.123}$, respectively) models produced chronologies with significantly higher correlations than all other models (p<0.001 for all) but that were not significantly different from each other (p=0.998, p=1.000, respectively). For the positive imposed trend, BAI ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.864}\pm \mathrm{0.236}$) and RCS ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.900}\pm \mathrm{0.162}$) produced chronologies with correlations significantly higher than CD ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=-\mathrm{0.503}\pm \mathrm{0.329}$) and CM ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.746}\pm \mathrm{0.306}$) (with p<0.001 for all examples) but that were not significantly different from each other (p=0.996). Notably, correlations exhibited by BAI chronologies were dependent on size sampling thresholds with BAI chronologies performing best when size thresholds exceeded 50 cm DBH (Fig. 3a). At this threshold BAI chronologies produced significantly higher correlations than when all trees were sampled (p=0.003) and when trees >10 cm DBH were sampled (p<0.001). The CD method produced chronologies that exhibited the lowest average correlation with the imposed positive trend of all models (p<0.001 for all).

Alternatively, when considering negative imposed trends, BAI ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.745}\pm \mathrm{0.426}$) chronologies performed significantly worse than RCS ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.706}\pm \mathrm{0.281}$, p<0.001) but still better than CD ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=-\mathrm{0.609}\pm \mathrm{0.291}$) and CM ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.666}\pm \mathrm{0.364}$) (p<0.001 for both). Again, CD chronologies exhibited significantly lower correlations than all other models (p<0.001 for all). Notably, RCS chronologies produced at the 50 cm DBH sampling threshold exhibited significantly lower correlations than all other sampling thresholds (p<0.001; Fig. 3a). All other models exhibited similar correlation distributions across the various size thresholds used for sampling.

3.3 Simulated white pine tree ring data

In simulated white pine data, two-way ANOVA suggested a significant effect of both standardization model (p<0.001) and minimum size sampling threshold (p<0.001) on average correlations for both the positive and negative logistic trend analyses. For the positive trend, chronologies produced by SDS ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.977}\pm \mathrm{0.026}$), RCS ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.932}\pm \mathrm{0.091}$), COMB ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.956}\pm \mathrm{0.052}$), and CM ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.953}\pm \mathrm{0.045}$) produced high correlations across all sampling thresholds, with SDS performing significantly better than CM (p=0.006) and RCS (p=0.001). All four models produced significantly higher correlations than those produced by BAI ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.899}\pm \mathrm{0.222}$) or CD ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.767}\pm \mathrm{0.126}$) chronologies, with CD producing the lowest correlations of all models. Contrasts suggested that the significant effect of minimum size threshold was driven by significant differences in correlations from BAI chronologies across sample thresholds, whereby BAI chronologies exhibited significantly lower correlations when no minimum size thresholds (i.e. all trees sampled) were employed (p<0.001 in all cases; Fig. 3b).

When examining negative imposed trends, SDS ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.942}\pm \mathrm{0.090}$) and COMB ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.904}\pm \mathrm{0.097}$) models produced chronologies with significantly higher correlations than all the other models but that are not significantly different from each other (p=0.594). BAI ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.750}\pm \mathrm{0.390}$) and RCS ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.772}\pm \mathrm{0.245}$) produced chronologies with correlations significantly higher than CD ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=-\mathrm{0.505}\pm \mathrm{0.316}$) and CM ($\stackrel{\mathrm{‾}}{r_{\mathrm{s}}}=\mathrm{0.623}\pm \mathrm{0.362}$) (p<0.001 for all) but that were not significantly different from each other (p=1.00). CD chronologies exhibited significantly lower correlations than all other models (p<0.001 for all). Contrasts suggested that the significant effect of minimum size threshold was driven by significant difference in correlations of chronologies produced by BAI and CM among sampling thresholds. As evident in Fig. 3b, BAI chronologies performed significantly better when sampling thresholds exceeded 50 cm DBH, and CM chronologies performed best when sampling thresholds exceeded 30 cm DBH.

3.4 Comparisons of methods in real tree ring data

Standardization methods were evaluated on the basis of correlations between their resultant chronologies and known time-related trends in tree ring series from shade-tolerant species.

Confidence intervals surrounding chronologies produced from each of the standardization methods applied to the tree ring series from six sugar maple stands are provided in Fig. 4a for both positive and negative logistic trends. The corresponding distributions of Spearman's rank correlation coefficients are provided in Fig. 5a with significant differences (p<0.05) being denoted by letters. Chronologies and corresponding correlation coefficients for the identical analysis performed on 12 red spruce stands are provided in Figs. 4b and 5b.

Figure 4 Standardized chronologies produced by each standardization method (legend on the right side) applied to tree ring series from (a) sugar maple (n=6) and (b) red spruce (n=12) stands. Solid lines represent the resultant model-wise mean chronologies across all stands considered, while ribbons represent respective 95 % confidence intervals. Dotted lines indicate the standardized positive (left side) or negative (right side) logistic trend that was added to the raw tree ring data.

Figure 5 Standardized chronologies produced by each standardization method applied to tree ring series from (a) sugar maple (n=6) and (b) red spruce (n=12) stands. Solid lines represent the resultant model-wise mean chronologies across all stands considered while ribbons represent respective 95 % confidence intervals. Dotted lines indicate the standardized positive (left side) or negative (right side) logistic trend that was added to the raw tree ring data.

Regardless of trend direction, RCS, COMB, and SDS chronologies exhibited comparable and consistent results across both species (Fig. 5). In general, chronologies produced by all three methods exhibited conservative but reliable estimations of the imposed trends (Fig. 4). SDS produced chronologies with correlations as high or higher (Fig. 5b, negative trend) than traditional RCS chronologies. Notably, the BAI and CM methods produced strong positive correlations between chronologies and the imposed trend only when the imposed trend was increasing (Figs. 4, 5) but both consistently failed to reproduce negative trends (Fig. 4). Finally, across both species, CD chronologies exhibited low correlations with the imposed trend regardless of direction (Figs. 4, 5).

4 Discussion

4.1 Size-deterministic vs. age-deterministic models for long-term trend reconstruction

Using simulated tree ring data from the shade-tolerant species sugar maple, we have shown that standardization models that include tree size in the year of ring formation (SDS, COMB) produced chronologies that retain long-term and low-frequency variation better than those produced by models that only included age as a predictor (RCS). Alternatively, in the shade-intolerant white pine species, chronologies produced by the RCS and COMB models showed no significant difference in their estimation of long-term trends, though SDS chronologies slightly outperformed RCS chronologies. As discussed previously, the finding that size-based standardization models perform well in simulated tree ring data is not surprising given that the SORTIE model calculates annual tree growth as function of tree size. Thus, the underlying growth trend would be well-approximated by a flexible curve estimated on the basis of tree size. As such, we use these results solely to elucidate the problematic nature of age-based standardization methods for shade-tolerant species. SORTIE's use of diameter, rather than age, as a determinant of tree growth is not arbitrary; it is well-established that tree metabolic processes are directly related to size (West et al., 2001). Additionally, there is little evidence for a unique effect of age on tree growth that is independent of size (Munné-Bosch, 2007). With the exception of dendrochronological models, the vast majority of individual tree growth and process models are indeed size-based. It follows that the ubiquitous use of age or calendar year in tree ring standardization methods (RCS, signal-free standardization, CD, Hugershoff curves) is a practice born out of convenience rather than physiological consideration. As such, we agree with previous accounts that this assumption may be especially problematic in shade-tolerant trees where age and size may not be perfectly correlated (Peters et al., 2015; Bontemps and Esper, 2011).

Unfortunately, all systematic comparisons of tree ring standardization methods in real tree ring data (e.g. Sullivan et al., 2016) are limited by their inability to validate long-term trends estimated by chronologies. In this study we evaluate standardization methods on their ability to reconstruct artificial trends in tree ring data. We show that SDS and COMB models are as reliable as the traditional RCS method in accurately detecting long-term trends in shade-tolerant species. Further, SDS appears to provide more reliable reconstructions when the underlying trend is negative. To our knowledge, only one other study has evaluated size-deterministic models on the basis of long-term trend reconstruction in chronologies. Bontemps and Esper (2011) compared RCS and SDS chronologies in common beech (Fagus sylvatica L.) and conclude that both exhibit similar variations, with the magnitude of difference varying between 3 % and 7 %. However, other studies have examined the influence tree size in explicit models of BAI. In tropical tree species of varying shade-tolerance, Nock et al. (2011) note that linear mixed models of BAI that included tree diameter had more support than those that included age. This result is corroborated by analyses of mixed models of BAI in Mediterranean pine species that suggest that the effect of DBH on BAI is more important than the effect of tree age (Marqués et al., 2016). In line with discussion above, Nock et al. (2011) attribute this finding to size being a more important determinant of light capture as it relates to tree height and crown size (King et al., 2005).

The resultant chronologies are indeed more likely to be influenced by the sample of the underlying tree population than by choice of standardization model. Tree age can be difficult or impossible to accurately estimate for some trees. In contrast, annual tree size can be reliability estimated from DBH and tree ring measurements more ubiquitously. We note that in this study only 66 % of sugar maple trees could be accurately aged. Since unaged trees are likely to be the oldest trees in the chronology, it follows that RCS chronologies may exhibit poor sample replication (especially in early years) and may be significantly shorter than those typically produced by SDS or COMB models. This has obvious implications for data quality and suitability. Considerably problematic is the “segment length curse”, whereby almost all standardization methods are ill-equipped to estimate long-term trends on timescales greater than or equal to the length of the chronology itself (Cook et al., 1995). Excessively short RCS chronologies are therefore limited in their application. A large advantage of SDS and COMB models is that they can incorporate otherwise inadmissible tree ring data.

This study does not explicitly test the efficacy of COMB models relative to SDS in the presence of unaged trees. We have also not provided evidence to suggest that the added complexity of COMB models relative to SDS is beneficial to accurate reconstruction of trends in the resultant chronologies. Given the merit the of size-deterministic models presented here, we suggest that future research explores the implications of the trade-off between model information and complexity in the presence of unaged trees.

4.2 BAI, CM, and CD methods for long-term trend reconstruction

The finding that CD did not produce accurate long-term trends in simulated tree ring data is consistent with our expectations (Peters et al., 2015; Briffa et al., 1992). We maintain CD should be avoided if the goal is long-term reconstruction from tree ring data. More interestingly, we have shown that CM and BAI, although designed for shade-intolerant open growth trees, do not reliably reconstruct negative long-term trends in simulated white pine tree ring data. Further, our analysis suggests BAI is less reliable when small and young trees are sampled. This result is corroborated in our study by a failure of both methods to reconstruct negative trends in real sugar maple and red spruce tree ring data. Further, this finding is in line with Peters et al. (2015), who note low reliability of BAI and that BAI is likely to produce erroneous trends when the underlying trend is of low signal, as would be the case for young and small trees that have low BAI rates and low climate sensitivity.

Both BAI and CM impart a strict relationship between tree size and growth. It has been suggested that this relationship may not account for the entire biological growth trend, leading to the maintenance of erroneous long-term trends in the resultant chronologies (Peters et al., 2015). Erroneous increasing trends are indeed noted in both sugar maple (Fig. 4a) and red spruce (Fig. 4b) chronologies produced by BAI and CM in our study. Accordingly, we caution future studies in their interpretation of BAI and CM trends in low-signal tree ring series.

Other studies have explicitly modelled size and/or age effects on BAI using a mixed-effect modelling approach (e.g. Marqués et al., 2016; Camarero et al., 2015; Nock et al., 2011; Martínez-Vilalta et al., 2008). We suggest this approach may better account for species- and site-specific factors that influence expected growth rates, leading to more accurate estimates of long-term trends in the resultant chronology. While our findings regarding the importance of inclusion of size in tree ring standardization models are presented in the context of raw tree ring width models, they are also directly relevant to explicit models of BAI.

A more thorough discussion of the limitations of the CD, BAI, and CM methods, as relevant to reconstruction of long-term trends, is beyond the scope of this study. The interested reader is directed to Peters et al. (2015).

4.3 Other considerations and future research

It is important to note that the goal of this study was not to explicitly test the effect of sample biases (i.e. modern sample bias, selection bias, etc.) on trend reconstruction but instead to assess reliability across different underlying sampling distributions. Accordingly, our results do not suggest that any of the discussed standardization methods are immune to sample biases (i.e. big tree selection bias, slow grower survivorship bias), as our study is not designed to detect and isolate the effects of contemporaneous differences in growth among trees that produce these biases. There is now substantial evidence to suggest that the long-standing practice of sampling only dominant trees or trees exceeding a minimum size threshold within a stand leads to considerable bias in the resultant chronology (Nehrbass-Ahles et al., 2014; Brienen et al., 2012; Briffa and Melvin, 2011). This bias is consistent across standardization methods (Duchesne et al., 2019; Nehrbass-Ahles et al., 2014). We maintain that in cases of long-term trend reconstruction, stands should be sampled according to the underlying stand age and size distribution, either through use of fixed-plots or random tree selection, regardless of the standardization procedure used.

Given the underlying physiological justification of the models presented here, we have no reason to suggest they are not broadly applicable to species of all shade-tolerance levels. We recommend future studies investigate the applicability of SDS and COMB models to both tree ring width and BAI data in a wider range of species. That said, shade-tolerant and broadleaf species and their applicable standardization procedures are underrepresented in dendrochronological studies (Zhao et al., 2019). Further, the applicability of enhanced tree ring standardization models (including traditional RCS and BAI) to global tree ring data sets are limited by widely unavailable metadata (i.e. tree age and DBH) in tree ring databases. Accordingly, we recommend more stringent requirements on the inclusion of applicable metadata in global databases in order to accommodate more complicated standardization models. We advocate for continued refinement of tree ring standardization procedures that are relevant to the ecological questions they aim to address.