Interactive comment on “ Forward modelling of tree-ring width and comparison with a global network of tree-ring chronologies ” by P .

The second issue is the focus of the paper. If the Discussion were much more clearly and closely aligned with the questions raised in the Abstract and Introduction, this would help. A tightened focus of the paper and deeper exploration of the results would provide the basis for a much more insightful discussion. This would make the paper a considerably more useful reference source for others wishing to explore the use of process-based models for simulation purposes, and/or the subsequent development of transfer functions for specific species/regions.


Introduction
Information derived from tree rings is one of the most powerful tools available for studying past climatic variability as well as identifying fundamental relationships between tree growth and climate (e.g.Fritts, 1976;Briffa et al., 2002a, b;IPCC, 2007).Tree-ring archives offer many advantages for climate reconstructions including their wide spatial distribution, annual resolution, calendar-exact dating, and high climate sensitivity (Hughes, 2002;Jones et al., 2009).However, relationships between tree growth and climate are complex and classical calibration methods (i.e.establishing a statistical relationship between instrumental and proxy data during their period of overlap) are limited for several reasons.Firstly, empirical methods generally treat climate as a function of the proxy rather than quantify how proxy signals are driven by climate variability.Secondly, statistical calibration is bound by the assumption that the relation between proxies and the target variable is stationary and independent of other climate variables -conditions that may not be adequately met (Raible et al., 2006).Thirdly, most tree-ring-based reconstructions are limited to geographic regions where growth variability is predominantly constrained by a single environmental variable.For example, temperature reconstructions are limited to regions where growth is Published by Copernicus Publications on behalf of the European Geosciences Union.source: https://doi.org/10.7892/boris.48741| downloaded: 26.12.2019P. Breitenmoser et al.: Forward modelling of tree-ring width primarily limited by cold temperatures.This leaves the vast areas away from high mountains and high-latitude environments poorly represented by high-quality temperature proxy data (Wahl and Frank, 2012).
In recent years, data assimilation approaches have found an increasingly important role in palaeoclimatic reconstruction efforts (Hughes and Ammann, 2009;Trouet et al., 2009;Goosse et al., 2010;Widmann et al., 2010;Franke et al., 2011;Bhend et al., 2012;Tingley et al., 2012;Steiger et al., 2014; for an introduction and overview see Hakim et al. (2013) and Brönnimann et al. (2013)).Combining information from climate models and proxy archives, also including associated errors and their covariance, provides a physically plausible representation of past climate consistent with available proxy data.An observation operator, which expresses the proxy as a function of the model data, is the formal link between model and proxy data.So-called forward proxy models (process-based or empirical) could potentially be used as an observation operator in data assimilation approaches (Hughes and Ammann, 2009).
The Vaganov-Shashkin (VS;Vaganov et al., 2006Vaganov et al., , 2011) ) process-based forward model has been demonstrated to skilfully simulate tree-ring width in North America and Siberia (Anchukaitis et al., 2006;Evans et al., 2006;Vaganov et al., 2011), China (Shi et al., 2008;Zhang et al., 2011), and Tunisia (Touchan et al., 2012).This forward-model of tree-ring growth was developed to quantify tree-ring formation as a function of climate and environmental variables (Vaganov et al., 2006(Vaganov et al., , 2011)).Yet applications generating "pseudo ring-width" series from global climate models have so far been hindered by the complexity of forward models needed to accomplish this task.Recently, Tolwinski-Ward et al. (2011a, b) introduced a simplified model version (the so-called Vaganov-Shashkin Lite (VSL) model) that skilfully reproduced climate-driven variability in North American tree-ring width chronologies.This simple and efficient model of non-linear climatic controls on tree-ring width requires only latitude, monthly temperature, and monthly precipitation as inputs and has the potential to simulate tree rings across a wide range of environments, species, and spatial scales.The VSL model may thus theoretically serve as a link between the proxy data and climate variables in any reconstruction methodology that can support use of a forward model.
The goals of this study are to (a) to assess the VSL model performance by examining the relations between simulated and observed growth at 2287 globally distributed sites, (b) identify optimal growth parameters found during the model calibration, and (c) to evaluate the potential of the VSL model as an observation operator for data-assimilation-based reconstructions of climate from tree-ring width.This study is the first report examining the VSL model's application on a global scale, and thereby across a wide range of species and site conditions, using a gridded data set based on direct observations from meteorological stations.All raw tree-ring width measurements from the International Tree Ring Data Bank (ITRDB) are processed and analysed in a consistent manner.
The paper is organised as follows: we provide a brief introduction to the climate data, the structure and parameterisations of the VSL model, and the tree-ring chronologies.The main results and discussion detailing the coherence between observed and simulated growth and associated climatic controls across space and over time are then presented.Finally, based on our findings, we discuss possible pathways of how tree-ring data could be used in data assimilation and climate reconstruction contexts.

Climate data
The climate data used to simulate tree-ring growth are monthly CRU TS3.1 mean near-surface temperature and CRU TS3.1.1 precipitation data sets from the Climatic Research Unit (CRU; Mitchell and Jones, 2005).These data cover the global land surface (except for Antarctica) at a 0.5 • × 0.5 • spatial resolution from 1901 to 2009 (updated from Mitchell and Jones, 2005; http://badc.nerc.ac.uk/data/ cru) and are based on direct observations from meteorological stations.It should be kept in mind that the temporal and spatial density of the underlying observations is variable.For instance, comparison of the available network of tree rings (see Fig. 1 for locations) with instrumental station records (see Jones et al. (1997) and Frank et al. (2007) for gridded temperature data) shows that both the instrumental and treering networks are dense for the European Alps.On the other hand, tree-ring chronologies are much more prevalent in the boreal zone, such as in the vast areas of Siberia and northern Canada.However, both instrumental and tree-ring networks are both very sparse for most tropical regions.One exception is India, where instrumental data dominate over the available tree-ring network.Note that sites of tree-ring chronologies are, on average, located 164 m higher than the elevation of the corresponding CRU grid cell.Correction for elevational differences requires knowledge of the site-specific climatological lapse rate and is therefore difficult to achieve (e.g.tests assuming a constant moist-adiabatic lapse rate on a site-by-site basis did not improve fits between tree-ring and model data) and was not further considered.

Forward model description: the VSL model
We use the Vaganov-Shashkin Lite (VSL) model and a Bayesian parameter estimation approach (Tolwinski-Ward et al., 2013) to produce synthetic tree-ring chronologies (Tolwinski-Ward et al., 2011a, b;version 2.3 accessible from http://www.ncdc.noaa.gov/paleo/softlib/softlib.html).The VSL model only requires monthly mean temperature, accumulated precipitation, and latitude to model growth.A further 12 adjustable parameters are required related to climatic limitations on growth, parameterisations for soil moisture availability, and the calendar periods when annual increment of wood is responsive to climate.The net effect of the dominating non-linear climatic controls on tree growth are implemented in terms of the principle of limiting factors and threshold growth response functions.
For a specific site and given month i, growth is calculated from the minimum of the growth responses to temperature (g T ) and moisture (g M ), modulated by the response to insolation (g E ): (1) Annual ring width is then defined as the sum of the monthly growth increments G(i) over a specified growth period I (Tolwinski-Ward et al., 2011a): In Eq. (1), g E is the ratio of mean monthly day length L relative to that in the summer solstice month at the corresponding latitude, computed using standard trigonometric approximations: The partial growth rates g T and g M are defined as ramp functions, with growth parameters representing climate thresholds below which temperature T and moisture M are not high enough for growth (T 1 , M 1 ), and thresholds above which temperature or moisture sustain non-limited growth (T 2 , M 2 ).Values of these partial growth responses are between zero and one (Tolwinski-Ward et al., 2011a): . (5) Site-specific growth parameters (T 1 , T 2 , M 1 , and M 2 ) are necessary due to diverse environmental conditions, site ecologies, and species represented in this global network.Accordingly, the VSL growth response function parameters were determined for each chronology via Bayesian parameter estimation using temperature, precipitation, site latitude, and ring-width data during the calibration interval.This scheme assumes uniform priors for the growth response parameters, and independent, normally distributed errors for the ring-width model.The median of each of the growth parameters in the resulting posterior probability distribution was finally taken as the "calibrated" growth response values for a particular site.A more detailed description of the approach can be found in Tolwinski-Ward et al. (2013).The Bayesian approach is computationally efficient, is able to incorporate expert knowledge about the likely values for parameters, and tends to guard against model overfitting (Tolwinski-Ward et al., 2013).
In addition to the above-described Bayesian approach to estimate the optimal growth parameters, we also estimated the optimal growth parameters for every site using a simple optimisation procedure as described in Tolwinski-Ward et al. (2011a).Growth parameters in this optimisation approach were sampled uniformly across intervals, and the growth parameter set producing the simulation that correlated most significantly with the corresponding observed series was then used to simulate the tree-ring series (note that throughout this paper we use Pearson's correlation coefficients).Results from both growth parameter estimation approaches are generally similar.There was some evidence that the optimisation approach sometimes overfit the model due to the lack of independence between the fitting procedure and the treering-VSL relationships.Nevertheless, the comparable results support the validity of either of these approaches with some possible advantages for the Bayesian parameter estimation (Tolwinski-Ward et al., 2011a, 2013).
Soil moisture M(i) is calculated from monthly temperature and accumulated precipitation data, T (i) and P (i), via the empirical Leaky Bucket model of hydrology from the National Oceanic and Atmospheric Administration's Climate Prediction Center (CPC), allowing for sub-monthly updates of soil moisture to account for the non-linearity of the soil moisture response (Huang et al., 1996;Tolwinski-Ward et al., 2011a; codes available from http://www.ncdc.noaa.gov/paleo/softlib/softlib.html).
In addition to the four parameters controlling the simulated growth response to temperature and soil moisture, all other parameters are taken from published studies (e.g.Huang et al., 1996;van den Dool et al., 2003;Fan and van den Dool, 2004;Evans et al., 2006;Vaganov et al., 2006; Table 1.VSL model parameters.

Parameter description Parameter Value
Temperature response parameters Threshold temp.for g T > 0 Tolwinski -Ward et al., 2011a;Zhang et al., 2011;Touchan et al., 2012), as is outlined in Table 1, which have demonstrated the applicability of VSL for a wide range of environmental conditions.
As growth period, I , we use a 16-month interval, starting in the previous September and ending in December of the modelled year (previous March to current June for the Southern Hemisphere).This integration interval accounts for some persistence in tree-ring growth (leading to autocorrelation) and showed the best overall performance in tests from a subset of our data set (which are all sites in Africa, New Zealand, Tasmania, Chile, Mongolia, southeastern Asia, Austria, Norway, Florida, and California).This interval is furthermore consistent with the outcomes from Tolwinski-Ward et al. (2011a) in their simulations of North American tree-ring series and with the importance of autumn phenology controlling interannual variability of net ecosystem productivity (Wu et al., 2013).

Tree-ring chronologies
A total of 2287 standardised tree-ring chronologies were used for our primary analyses.To obtain these series, we considered all available raw tree-ring width measurements from 2918 sites from 163 different tree species from the ITRDB (Grissino-Mayer and Fritts, 1997; http://www.ncdc.noaa.gov/paleo/treering.html).Prior to further analysis, we attempted to identify and correct data and metadata errors including repetitive measurements, feet-meter encoding, inappropriate decimal points and hyphenations, incorrect series labelling, or misplaced series positions (see Table S1 in the Supplement).
To remove the biological age trend and other lowerfrequency signals that may be driven by stand dynamics or disturbances, we processed all data sets in the program ARSTAN (Cook, 1985) using standard dendroclimatological methods (i.e.detrending and standardisation to dimensionless growth indices).We tested different standardisation methods on a subset of 97 locations from all six continents and selected our detrending approach based on this subset.Methods tested include (a) a negative exponential curve and linear regression fits, (b) a smoothing spline fit with a 50 % frequency response cut-off (which is the wavelength at which 50 % of the amplitude of a signal is retained) equal to threefourths of a series length, and (c) a smoothing spline with a 200 yr frequency response cut-off.These tests and Pearson's correlation analyses between the differently standardised series and the monthly climate parameters led us to favour a hierarchical approach: fitting negative exponential curves and linear regression curves of any slope were applied as firstorder detrending options; but if these two methods were not suitable (e.g. because values of a fitted curve approached zero), a smoothing spline was fit with a 50 % cut-off frequency at 75 % of each series length.The detrending methods applied allow retention of climate-related signals in any given chronology of the order of the median segment length of the constituent tree-ring series (Cook et al., 1995).
Combining multiple measurements into a single ringwidth chronology should then ideally represent a coherent climate signal of a particular species at a site (Cook et al., 1990).Standard chronologies (TRW ITRDB ) of all these sites were developed by a biweight robust mean estimate (Cook et al., 1990) and the variance was stabilised to minimise artefacts from changing sample size through time (Osborn et al., 1997;Frank et al., 2007).Further, we required that the sample depth of a chronology (number of samples from trees used for every point in time) is always ≥ 8 and that the 1901-1970 interval was fully covered by tree-ring data.These criteria resulted in the 2287 chronologies for further analyses.

Aggregation of observed tree-ring chronologies
In addition to comparing individual observed and VSLmodelled chronologies, we further analysed whether the signal-to-noise ratio can be enhanced by aggregating several neighbouring observed series into one "super-chronology", which subsequently is modelled with VSL.
Aggregation of chronologies is performed using a search radius, R of 200 and 600 km, centred at each land grid cell of the CRU climate data set.The former radius approximates the grid resolution of many global climate models and the latter reflects a typical decorrelation length scale in large-scale precipitation and temperature fields.We do not use greatly expanded search radii like Cook et al. (2010) or Ljungqvist et al. (2012), for instance, whose work was motivated by the construction of a gridded climate reconstruction from spatially sparsely distributed proxy series under considerations of long-range teleconnections (Jones et al., 1997).For the larger search radius, we introduced the condition that at least one chronology must be within 200 km of a grid cell such that both search radii result in the same spatial coverage.
In contrast to modelling individual chronologies, modelling aggregated chronologies requires a first-order separation of the climatic influences on growth into "moisturelimited" and "temperature-limited" chronologies to enhance the signal-to-noise ratio.Note that this aggregation approach means that changes in the limitations over time cannot be accounted for as might apply to the debate over the socalled "divergence phenomena" (D'Arrigo et al., 2008;Esper and Frank, 2009).Yet, site aggregation is particularly important in mountainous terrain with more complex climatic conditions and rapidly changing conditions within small spatial scales that locally impact tree growth (Salzer and Kipfmueller, 2005).Accordingly, we differentiated the sites based upon their principal climatic drivers and thus separately aggregated all temperature-and moisture-limited chronologies located within each grid cell node's search radius.Specifically, we classified temperature-and moisturesensitive chronologies using the long-term averaged monthly growth responses g T and g M from the VSL model for all months where g > 0. Temperature-sensitive chronologies were defined as sites that have more temperature-than moisture-limited months (i.e. a greater number of months where g T > g M ).All other series are moisture-limited (equal or more months with g T > g M ).Aggregated tree-ring width (ATRW) series for temperature-and moisture-limited series were computed for every grid cell with one or more temperature and/or moisture-limited tree-ring chronologies located within the given search radius.For the 200 km (600 km) search radii, 10 % (7 %) of all land grid cells contain only temperature-limited series, 7 % (5 %) are only moisture-limited, and 11 % (37 %) of the cells include both temperature-and moisture-limited chronologies.
After identifying the series to be aggregated, we defined the aggregated tree-ring width ATRW ITRDB as weighted means of the corresponding TRW ITRDB series.Weights were assigned according to four factors (Eqs.6-9), and multiplication of the four relative weights gives the weight value we apply to each individual tree-ring series within the search radius.The weighting factors are (a) mean Rbar of each chronology, (b) mean EPS of each chronology, (c) distances d between the site location and the grid node, and (d) the smaller p value of two correlations between CRUTS3.1 temperature or precipitation and TRW ITRDB : for 12-month averages (January-December in the Northern Hemisphere, July-June in the Southern Hemisphere, respectively) and for 5-month averages for growing season interval (May-September and November-March, respectively).These two seasonal windows consider the summer months when climate generally dominates tree growth (Briffa et al., 2002a) and also possibilities that growth is driven by a much wider seasonal window.
Weights are determined by to following functions (see Fig. S1), which vary between 0 and 1: There is no a priori reason to choose any particular weighing function, but we determined the functions based on theoretical considerations (e.g.see Ljungqvist et al., 2012, who used the Gaussian distance) and the parameter distribution as shown in Table 1.We re-calibrated the VSL simulations for all grid cells containing at least one aggregated tree-ring series.The calibration was performed using monthly temperature and precipitation series from CRU and with the same VSL model parameter set-up described in Table 1.

Experimental and analysis design
The following experimental design is used in our study.First, we modelled (following Sect.2.2) each individual observed chronology with data from the closest CRU grid cell as climatic input.These simulated chronologies are termed TRW VSL .In addition, we also modelled the aggregated chronologies using climate data from the centre of the search radius (closest grid cell), but calibrated using the aggregated observed chronologies ATRW ITRDB .These chronologies are termed ATRW VSL .All observed and modelled chronologies were normalised (i.e.given a mean of zero and standard deviation of unity) over the 1901-1970 period.
For both TRW VSL and ATRW VSL , we utilised 1901-1970 as the main calibration period for the growth parameters.However, for assessing result consistency, we also performed split-sample validation experiments in which we used the 1901-1935 period for calibration and the 1936-1970 period for validation and vice versa.

Quality indicators for tree-ring chronologies
In order to assess the chronology or site attributes that may affect the robustness of the climate signal in TRW ITRDB , we computed various characteristics of all 2287 ring-width chronologies for the 1901-1970 common period (Fig. 2; see Fig. 1 for locations).The average inter-series correlation between all series in a chronology (Rbar) indicates common variance (Cook et al., 1990), while the expressed population signal (EPS) quantifies the degree to which a chronology matches a hypothetical population chronology (Wigley  (Cook et al., 1990); EPS is the expressed population signal and the vertical dashed line denotes the commonly cited 0.85 EPS criterion (Wigley et al., 1984;Cook et al., 1990).Mean 1901Mean -1970 Rbar and EPS values were calculated using a 30 yr window that lags 15 yr; MSL is the mean segment length; and ASD is the average sample depth (48 bins were distinguished in each variable).et al., 1984; a threshold value of 0.85 is routinely, but arbitrarily, cited in dendrochronological literature).Rbar and EPS statistics were calculated over moving 30 yr windows with 15 yr overlap.A total of 93 % of all mean EPS values in our network are above 0.85 over the 1901-1970 time period.Rbar values roughly follow a normal distribution with a mean of 0.39.These statistics suggest that the vast majority of chronologies have a relatively robust signal that should allow climatic drivers to be inferred.On average tree-ring chronologies are composed of 30 series with a mean segment length of 170 yr.However, these distributions have relatively long tails reflecting sites that are exceptionally well replicated or composed of long-lived tree species (e.g.Bristlecone pine).The mean segment length and our detrending approach suggest that we should be able to test inter-annual to (at least) multi-decadal variability in our modelling efforts.The geographical distribution of the sites is dominated by low altitudes and the mid-latitudes, but there are also a considerable number of high-latitude and high-elevation sites particularly in the Northern Hemisphere.

Site-specific Bayesian parameter estimation
The Bayesian approach provides estimates of the optimal parameters T 1 , T 2 , M 1 , and M 2 for each site, with descriptive statistics of the parameter estimation for split-sample validation given in Table 2. Notably, we find great consistency among the temperature and moisture parameters that were calculated independently and separately for the 1901-1935 and 1936-1970 periods.We find for example that the mean temperature for growth onset is 4.6 and 4.8 • C for the early and late calibration periods.Similarly, temperature no longer limited growth above 16.3 and 16.5 • C again for the independent early and late calibration periods.Figure 3   ther illustrates the estimated temperature and moisture parameters T 1 , T 2 , M 1 , and M 2 for each site and different tree genera.In terms of the currently debated (see Discussion) growth-onset threshold temperature (T 1 ), we find most estimates fall within 4-6 • C with values for all sites and species between 1.80 and 8.41 • C (see also Table 2).Despite the generally narrow range of the parameters for all study sites, some evidence for species-dependent values also exists.For spruce, our estimates of T 1 and T 2 are generally higher than for pine trees, while M 2 is lower.We also find differences among the pine species.For ponderosa and pinyon pines, T 1 is about 1 to 2 • C lower than for other pine species (Fig. 3e, compare also to Fig. 3a in which these two species clearly influence the shape of the pine histogram towards cooler threshold temperatures).
Although these parameter values were computed independently, we found clear associations among the four growth parameters.For instance, T 1 and T 2 values and M 1 and M 2 Table 2. Statistics of the Bayesian estimation of the site-by-site tuned VSL growth response parameters T 1 , T 2 , M 1 , and M 2 for the 1901-1935and 1936-1970calibration intervals. Calib. 1: 1901-1935 Calib.values are significantly correlated with each other (r = 0.60 and r = 0.47, both p < 0.01).The often negative correlation between temperature and precipitation at inter-annual time scales is expressed by the inversely related temperature and moisture growth parameters, with the relationship between T 2 and M 2 (r = −0.63,p < 0.01) being the most noteworthy.Moreover, temperature-limited sites tend to have higher T 1 and T 2 values in comparison to the moisturelimited sites, whereas higher M 1 and M 2 values are characteristic of moisture-limited sites (please refer to Sect.2.4 for the definition of temperature-and moisture-limited sites).

Relationships between actual and simulated tree-ring growth
The VSL model was able to simulate tree-ring series for 2271 out of the 2287 sites, suggesting general applicability to diverse environments and species.Of the remaining 16 sites, nine were high-latitude or high-elevation sites in Nepal (3), Siberia (1), and Canada/Alaska (5) for which no growth was simulated because temperatures never reached the threshold temperature T 1 for growth initiation (hence g T = 0).The remaining seven sites were (sub-)tropical sites in Florida (1), Mexico (4), Argentina (1), and Indonesia (1), for which neither temperature-nor moisture-limited growth (hence g T = 1 and g M = 1).Similarity between the 2271 actual TRW ITRDB and the corresponding modelled TRW VSL chronologies was assessed by Pearson's correlation coefficients during the 1901-1970 period (Fig. 1).The average of Pearson's correlation coefficients between all TRW ITRDB and TRW VSL is 0.29 (standard deviation = 0.22).Approximately 41 % of the correlation coefficients between TRW VSL and TRW ITRDB are statistically significant at the 95 % confidence level (r > 0.35), taking into account the reduction of degrees of freedom due to autocorrelation (Mitchell Jr. et al., 1966).Note, however, that the 1901-1970 period includes the calibration interval.Thus, analysing the correlation coefficients between TRW ITRDB and TRW VSL in the split-sample validation, we find significant (p = 0.05; r >∼ 0.44) correlations during all four calibration/validation intervals (calibration 1901-1935 and validation 1936-1970 and vice versa) for 13.3 % of all the series.Significant pairwise relationships between the TRW ITRDB and TRW VSL series are found across the globe on all continents.In particular, geographic areas with tendencies for higher model-data agreement are found in central Europe, Scandinavia, eastern/central Siberia and in the United States with a cluster of high correlations in California where r > 0.8.Weak relationships are characteristic of the Africa sites, but are also evident in New Zealand, Tasmania, Alaska, a NW/SE band from central Canada to NE USA, the northern part of the British Isles, and parts of southern Europe (Fig. 1).
To illustrate the monthly growth response values, g T and g M , we plot both the year-to-year variability and mean values over all simulated years  for one high-altitude (Fig. 4a) and one low-altitude site (Fig. 4b month-wise minimum of either temperature or precipitation, Fig. 4a clearly shows that growth at the high-altitude site is temperature-limited (g T < g M ).The lower-elevation site (Fig. 4b), on the other hand, is temperature-limited only briefly at the beginning and at the end of the growing season, while precipitation limits growth during the summer months (g M < g T ).Importantly, tree growth at the highelevation site is negatively correlated with growth at the lower-elevation site.This negative correlation holds for both the observed and modelled chronologies with coefficients of −0.13 and −0.44, respectively.These relationships demonstrate the VSL model's ability to incorporate joint influences of temperature and precipitation on growth, and simultaneously suggest care is required when developing regional tree-ring time series (see Sect. 3.4.2below).Our simulated growth responses for two selected sites in Switzerland appears to be widely applicable across mountainous forest regions with large climatic gradients such as the southwestern USA (Salzer and Kipfmueller, 2005;Tolwinski-Ward et al., 2011a) and Mongolia (Poulter et al., 2013).Oneyear lagged autocorrelation coefficients (AR1) were assessed from TRW ITRDB and TRW VSL series at all sites.Mean coefficients of 0.46 and 0.42, respectively, over all series clearly reveal persistence; i.e.21 and 18 %, respectively, of the variance in the width of the ring in the current year can be predicted based upon simulated growth during the past year.

VSL simulations for aggregated TRW series
The aggregation generally enhances the correlation between ATRW VSL and ATRW ITRDB .For instance, the average of all coefficients for R = 200 km during the common 1901-1970 interval is 0.34 for the temperature-limited series and 0.38 for the moisture-limited series compared with a mean of 0.29 for the non-aggregated tree-ring series.The spatial patterns of correlation coefficients for temperature-and moisture-limited ATRW for both search radii are shown in Fig. 5. Notably, ATRW VSL and ATRW ITRDB show spatial coherency and capture the main climate signals and are thus suitable for data assimilation approaches (see Discussion).The 600 km search radius improves the relationships for both temperature and precipitation.Coherent regionalscale correlation patterns for temperature-limited series can be found in Scandinavia, western/central Siberia, and along the western and eastern United States (Fig. 5a and c).Regional correlation patterns for moisture-sensitive series encompass the United States and central Canada (Fig. 5b and  d).Regions such as central Europe, Turkey, central Asia, and the Andes are characterised by both temperature and moisture limitations, which is an expression of a complex terrain with small-scale climatic features and large climatic gradients.High correlation coefficients between ATRW VSL and ATRW ITRDB are present in all those regions, with a clear tendency towards temperature-limited, high-elevation sites.
Poor agreement for the temperature-limited series can be  found in Tasmania, New Zealand, western-central Canada, and in the eastern-central US, while moisture characteristics are poorly captured in Patagonia, northern Canada and in the European Alps.
Based on the split-sample validation experiment, Fig. 6 shows a scatterplot of the correlation coefficients between the ATRW VSL and ATRW ITRDB , 200 km search radius, during different calibration/validation intervals (1901-1935/1936-1970, and 1936-1970/1901-1935). Points along the 1 : 1 line are sites for which VSL produces simulations that are stable in time, whereas points below this line represent sites for which the model yields a higher correlation during the calibration period.We find model skill (r > 0.44) in Pearson's correlation coefficients in 36 % of the series during the calibration period, of which 54 % are also significant in the validation period.Correlations are stable for moisturelimited sites in large areas of the United States and central Canada, the lower-elevation sites in west-central Europe, Turkey, Mongolia, the western Himalayas, and in mid-Chile, whereas temperature-limited sites dominate Scandinavia, Siberia, the eastern Himalayas, central Europe, southern South America, and parts of the western United States.Correlations are not particularly skilful for New Zealand, Tasmania, and Patagonia.

Growth parameters
The site-by-site parameterisation for the moisture and temperature thresholds revealed interesting global-scale patterns.Notably, the parameter estimation of the VSL model provides novel insights into the currently debated temperature threshold below which growth may no longer take place.Our estimations of growth-onset temperature (ca.4-6 • C) is in agreement with the temperature range found by Tolwinski-Ward et al. (2013), who used a four-parameter beta distribution to calibrate the model parameters on a subset of 277 se-ries in the United States.More striking is the convergence of findings with assessments of growing season temperatures at tree line based upon in situ loggers or remote sensing and interpolated climatic data (Körner and Paulsen, 2004;Körner, 2012), as well as observational studies of wood formation (Rossi et al., 2007).The wide range (1-9 • C) of the priors considered in our approach allowed ample room for deviation from the previous literature, so we interpret the agreement from these various lines of evidence as strong support for growth onset thresholds near 5 • C.
At the level of individual species or individual sites, our results compare well with findings of Vaganov et al. (2006), who, based on the VS model, reported minimum and optimal temperatures in the northern Taiga of 4 and 16 • C for larch, 6 and 18 • C for spruce, and 5 and 18 • C for Scots pine, respectively.The differences in growth thresholds might reflect the better adaption of pines to lower-temperature conditions and drought.A comparison of the site values of the threshold temperature T 1 at the corresponding elevations for different tree genera is shown in Fig. S2 in the Supplement.Species-related difference in moisture parameters may be due to needle structure as well as frequency and density of stomata which influence water loss due to evapotranspiration (Vaganov et al., 2006).Our results contribute to understanding species specific growth responses as recently documented in an analysis of 36 species from across the European continent (Babst et al., 2013).
At the same time, uncertainties arise from the parameterisation of M and estimation of M 1 and from the missing representation of winter precipitation stored as snowpack in the CPC Leaky Bucket model and its impact on the timing of snowmelt (Tolwinski-Ward et al., 2011a).Temperaturelimited sites tend to have higher T 1 and T 2 than moisturelimited sites and vice versa for M 1 and M 2 .This may be an expression of different climatic conditions, but may also be related to edaphic conditions and plant physiological differences.Figure S3 of the Supplement shows scatterplots of joint posterior relationships at each site between parameters (a) T 1 and T 2 , (b) M 1 and M 2 , (c) T 2 and M 2 , and (d) T 1 and M 1 .Application of the VSL model in more controlled frameworks, such as in stands with mixed species compositions, will be useful to isolate and define species-specific growth characteristics while simultaneously excluding the possibilities that non-uniform species distributions along elevational gradients contribute to systemic differences in the obtained model parameters.

Relationships between actual and simulated growth
The comparison between observed (ATRW VSL ) and modelled (ATRW ITRDB ) tree-ring width series showed generally good agreement and robust calibration/validation statistics.We found significant model skill in independent verifications across all forested continents -tendencies that increased when aggregating the tree-ring sites.

P. Breitenmoser et al.: Forward modelling of tree-ring width
The broad-scale pattern of temperature-limited growth for high latitudes (e.g.Taymyr Peninsula, Russia), and high altitudes (e.g. the European Alps) extends previous analyses of subsets of the ITRDB network (e.g.Briffa et al., 2002a;Wettstein et al., 2011;Babst et al., 2013).Our findings furthermore largely match the distribution of potential climatic constraints to plant growth in Nemani et al. (2003) and Beer et al. (2010) as well as the climate classification based on Köppen-Geiger by Kottek et al. (2006).
However, also prominent in our analysis were locations (e.g.Tasmania and New Zealand) where VSL simulations revealed little skill.To explore possible reasons for this we conducted initial simulations where climate variability was forced from a reanalysis product (20CR; Compo et al., 2011) and from the individual climate stations themselves.We found some evidence for improved skill which suggests that (a) the quality of the instrumental data sets and (b) strong orographic effects can affect our simulations.Caution is particularly warranted where quality and quantity instrumental stations rapidly decreases back in time (e.g. of the Mediterranean, Asia, and Southern Hemisphere locations with low population density and/or strong orographic effects).
The fact that VSL reproduces first-order autocorrelation structures (AR1) similar to the actual tree-ring chronologies highlights the importance of including last year's growing season for growth simulation in the VSL model in order to retain potentially useful climate information contained in tree rings (see Table 1).Similarly, Wettstein et al. (2011) report median first-order autocorrelation coefficients of 0.47 for 762 ITRDB-derived tree-ring width series located poleward of roughly 30 to 40 • northern latitude, while Tingley et al. (2012) found first-order autocorrelation coefficients of approximately 0.6 for the tree-ring width series used in Mann et al. (2008).This reasonable agreement between the autocorrelation structure in actual versus simulated tree-ring data suggests that mechanistic modelling may help to reduce the effect of biases in the spectral characteristics of tree-ring chronologies on climate reconstructions (Meko, 1981;Cook et al., 1999;Franke et al., 2013).The global applicability of a 16-month seasonal window as applied herein will require further assessments for individual regions and species (for instance as in Wu et al., 2013).It is also conceivable to apply differential weighting to previous years' growth, but such schemes would be best implemented when understanding of lagged physiological processes and the roles of carbon reserves are further advanced (e.g.Babst et al., 2013;Zielis et al., 2013).

Use of VSL in data assimilation
The favourable calibration-validation statistics of the VSL model and its ability to skilfully simulate growth for a wide range of environments suggest that the model could possibly be used for data assimilation approaches.In the follow-ing section, several pathways are briefly outlined (see also Brönnimann et al., 2013) and connections are made to how such assimilation approaches could be implemented using the VSL model.
Data assimilation, in general, combines information from observations (here: tree-ring width, termed y) with a numerical model simulation (termed x b ) to obtain a physically consistent estimate of the climate state, x.It can be formulated as a cost function problem of the form (10) where B and R represent the error covariance matrices of the model simulation and of the tree-ring widths, respectively.H is the observation operator which simulates the tree rings from the model state.Provided that the model state encompasses all required climatic variables, this could be the VSL model.How VSL is used in the approach then depends on how Eq. ( 10) is minimised.Following Brönnimann et al. (2013), we distinguish "covariance-based approaches", "analogue approaches", and "nudging approaches".
Covariance-based approaches use observations (here: tree rings) to correct the model simulations based on the model error covariance matrix.Formally, these approaches (e.g.ensemble Kalman filter, 4D-VAR) require the matrix H, which is the Jacobian of H (requiring slight adaptations to VSL).A further constraint of this family of approaches is the state vector, which easily becomes excessively large.Bhend et al. (2012) proposed an off-line assimilation approach which does not require the entire model state to be analysed.
Analogue approaches minimise the cost function by selecting among existing simulations rather than by correcting them.Hence, the cost function becomes (11) It is minimised by evaluating it for all available x, which can be an ensemble of simulations (e.g.particle filter; Goosse et al., 2010) or different slices of a long control simulation (e.g.proxy surrogate reconstruction; Franke et al., 2011).The analogue approach does not pose any restrictions on the size of the state vector or on the observation operator.Hence, VSL or also its parent model VS (provided, again, that all relevant climate variables are in the model state) can directly be used, either in the form of individual chronologies (an example is given in Brönnimann et al., 2013) or aggregated chronologies.Moreover, in contrast to many implementations of the ensemble Kalman filter, R can be non-diagonal.The results from our study, i.e. the comparison of observed and modelled tree-ring width, can directly be used to estimate R. Furthermore, our results can be used to assess biases relative to observations, which need to be accounted for, e.g. as an additional term in the observation operator.
In the third approach, nudging, the cost function is not solved explicitly.Rather small increments are added to the tendency equations.These increments depend on a target, which might be a "classical" climate reconstruction.In this case, VSL may guide the aggregation of tree rings in order to increase the signal-to-noise ratio in classical reconstruction approaches.
Our analysis has advanced the notion that the VSL model has the potential to be used as an observation operator in palaeoclimatic data assimilation, albeit in different form depending on the approach chosen.
The VSL model within data assimilation approaches motivates future analyses on the "divergence problem" by explicitly being able to take changes of limitations over time into account.

Conclusions
We investigated relationships between climate and tree-ring data on global scale using the non-linear, process-based Vaganov-Shashkin Lite (VSL) forward model of tree-ring width formation (Tolwinski-Ward et al., 2011a).We investigated relations between actual tree-ring growth and growth simulated based on climatic influences.A network of 2287 globally distributed tree-width chronologies has shown to provide a large and high-quality sample space for these analyses.Bayesian estimation of the VSL growth response parameters yielded values that are stable with respect to the choice of calibration interval for a wide range of species, environmental conditions, and time intervals, showing the model's general applicability for worldwide studies using CRU climate input data.A benefit from this parameterisation approach was also a novel global assessment of the growth-onset threshold temperature of approximately 4-6 • C for most sites and species, thus providing new evidence for on-going debates regarding the lower temperatures at which growth may take place (Anchukaitis et al., 2012;Mann et al., 2012;Körner, 2012).Moreover, our results demonstrate the VSL model skill across a wide range of environments, species, and different time intervals.Spatial aggregation of tree-ring chronologies based upon chronology quality, climate response, and proximity reduced non-climatic noise contained in the tree-ring data and resulted in an improved relationship between actual and modelled tree growth.We further demonstrate that the potential of the VSL model can be exploited, for instance, in data assimilation approaches to reconstruct the climate of the past.

Fig. 2 .
Fig.2.Binned chronology characteristics for 1901-1970.Rbar is the average inter-series correlation between all series from different trees(Cook et al., 1990); EPS is the expressed population signal and the vertical dashed line denotes the commonly cited 0.85 EPS criterion(Wigley et al., 1984;Cook et al., 1990).Mean 1901Mean  -1970     Rbar and EPS values were calculated using a 30 yr window that lags 15 yr; MSL is the mean segment length; and ASD is the average sample depth (48 bins were distinguished in each variable). fur-
) in Switzerland.Simulations (TRW VSL ) as well as the observed chronologies (VSL ITRDB ) for both sites are shown in Fig. 4c.Because the simulated growth response is controlled by the www.clim-past.net/10etal.: Forward modelling of tree-ring width

Fig. 5 .
Fig. 5. World map displaying Pearson's correlation coefficients between ATRW VSL and ATRW ITRDB for a search radius of 200 km (a, b) and 600 km (c, d) at temperature-limited (a, c) and moisturelimited (b, d) grid cells.Note that in panels (c) and (d) at least one chronology has to be within the 200 km search radius to make the correlations comparable.