Boreal tree-rings are influenced by temperature up to two years prior to their formation: a trade-off between growth and reproduction?

Large spatial and between-tree variability has recently been observed in the response of boreal forests to ongoing climate change, spanning from growth stimulation by increasing temperatures to drought limitation. To predict future responses of boreal forests, it is necessary to disentangle the drivers modulating the temperature-growth interaction. To address this issue, we established two inventory plots (at a treeline and closed-canopy forest) and assembled site chronologies in Picea glauca stands at the transition between boreal forest and tundra in Northern Quebec, Canada. In addition to site chronologies, we established a set of chronologies containing, for each year, exclusive subsets of tree-rings with specific cambial age (young/old), tree dimensions (small/large) and tree social status (dominant/suppressed). All chronologies were correlated with climatic data to identify the course of climatic conditions driving variability in tree-ring widths. Our results show that the growth of P. glauca correlates significantly with summer temperature in tree-ring formation years and during up to two prior summers. Tree-ring width is positively influenced by summer temperatures in tree-ring formation year and two years prior to tree-ring formation. In addition, climate-growth correlations indicate a negative effect of summer temperature one year before tree-ring formation at the closed-canopy forest site. The pattern of climate-growth correlations is tightly synchronized with previously published patterns of climate-reproduction correlations of P. glauca, suggesting a growth-reproduction trade-off as a possible factor modulating the response of boreal forests to summer temperatures. Climatic signal does not differ between pairs of chronologies based on subsets of cambial ages, stem dimensions or tree competition status at the treeline site. However, the response to summer temperatures one year before tree-ring formation is significant only in mature (old, large and dominant) individuals at the closed-canopy site. The inverse pattern of temperature-growth correlations during a sequence of three years challenges predictions of how boreal forests respond to climate change.


Introduction
Boreal forests represent an important carbon sink storing approximately 32% of the global forest stock (Pan et al 2011). Their responses to recent climate change (ACIA 2004) are of crucial environmental, cultural and socio-economic importance, as they often positively feed back on to the climatic system (Dunn et al 2007, Babst et al 2019. However, the trends and intensity of these plastic changes are often nonlinear often non-linear, vary geographically and depend on forest stand structure (Lloyd and Bunn 2007). Indeed, although the arctic-boreal ecotone is considered the ecological limit of trees because of their strong temperature limitation (Körner 2012), previous studies often failed to identify a uniform pattern(s) of climate-growth interaction. Responses to climate change tend to vary along geographic gradients (Wilmking and Juday 2005), temporally (Wilmking et al 2004, Lloyd and Bunn 2007, D'Arrigo et al 2008, Ohse et al 2012 or as a result of forest stand conditions and local microsites (Nicklen et al 2019). Our understanding of the specific factors and mechanisms modulating the climate-growth interaction of boreal forests thus remains limited.
One possible reason behind spatial and temporal variations in climatic signal are tree micro-site conditions, which alter the climate-growth interactions of entire forest stands (Carrer 2011, Galván et al 2014, Buras et al 2016. Specifically in boreal ecosystems, the climatic signal of site chronologies seems to be driven by the proportions of trees with contrasting ecological preferences (Wilmking et al 2004, Wilmking and Juday 2005, Pisaric et al 2007 given by their dimensions and social status (Trouillier et al 2019). Sensitivity to climatic conditions (Carrer and Urbinati 2004, Mérian and Lebourgeois 2011, Konter et al 2016 and biotic disturbances (Caccianiga et al 2008) is often size-or age-specific, complicating the extraction of environmental information from site chronologies without dendrometric data for individual trees. It has been shown (e.g. Buras et al 2016) that site chronologies from uneven-aged boreal forests might not be a reliable proxy archive for palaeoclimatic reconstructions or for constructing prediction models (Soja et al 2007, Miller andSmith 2012). This issue can be overcome by analysing individual tree size or age cohorts separately (Trouillier et al 2019), which also provides more detailed insight into ecology of entire forest stands.
In addition, the disruption of the expected positive temperature signal in tree-ring chronologies of boreal species might be a consequence of tight trade-offs between individual plant sinks of carbohydrates. Net primary productivity is low in boreal forests (Gillman et al 2015) and, therefore, trees must balance the allocation of limited amounts of carbohydrates among individual survival strategies. For instance, the tradeoff between growth and reproduction results in synchronization of narrow tree-rings (pointer years) with years of extensive seed production (Juday et al 2003, Drobyshev et al 2010, Hacket-Pain et al 2018, Nicklen et al 2019. Alternatively or in addition, it may result in significant correlations between tree-ring widths and climatic conditions during key phases of cone development (Wilmking et al 2004, Ohse et al 2012. A good understanding of the factors modulating climate-growth responses of boreal ecosystems is a substantial prerequisite for predictions of future environmental changes in Arctic and sub-Arctic regions (Soja et al 2007). We aimed to address this issue by establishing two inventory plots in forest patches dominated by white spruce (Picea glauca [Moench] Voss.) in a remote region of Northern Quebec, representing a transition between boreal forest and tundra vegetation. We constructed site tree-ring chronologies for the plots and performed climategrowth correlation analysis, unusually considering months for up to three years before tree-ring formation to account for a possible delayed effect of carbon allocation to cone production (Roland et al 2014). We performed analyses separately for tree-rings of young/ old cambial ages, small/large stems, and suppressed/ dominant individuals (sensu Trouillier et al 2019). We hypothesized that (i) the pattern of climate-growth correlations would be inverse to the pattern(s) of climate-reproduction correlations (Juday et al 2003, Krebs et al 2012, Roland et al 2014, and (ii) as a result of a lower level of non-climatic noise (e.g. competition), chronologies based on old tree-rings and large or dominant trees would have a stronger climatic signal compared to ones based on young tree-rings and small or suppressed individuals.

Study area
The study was performed in the vicinity of Whapmagoostui-Kuujjuarapik village (further abbreviated as WK; N55°15′, W77°45′; formerly known as Poste-dela-Baleine) in Northern Quebec, Canada ( figure 1(a)). The area of the municipality is bordered by Hudson Bay of the Arctic Ocean to the west and the Great Whale River to the south. The climate is sub-Arctic with a mean annual temperature of −4.1°C and mean total annual precipitation 639 mm over the period of 1935-2017 (figure 1(c); compiled data from Environment Canada and CRU TS4.01, see below). There is a recently increasing trend in both mean annual temperature and mean summer temperature, most apparent since 1990s (figure 1(d)). The climate is continental in winter (as a result of freezing of Hudson Bay), with low temperatures and a lack of precipitation, and maritime in summer, with mild temperatures and sufficient precipitation.
The vegetation belongs to the forest sub-zone of the forest-tundra, characterized by high spatial heterogeneity of lichen-dominated granitic outcrops, grassdominated sandy areas (nearby the coast) and scattered spruce or shrub stands in valleys or on raised beach deposits ( figure 1(b (Payette 1983(Payette , 1993). An increase in tree cover density (Payette and Filion 1985), a switch of the dominant growth form from stunted to erect (Caccianiga and Payette 2006) and an expansion of treelines northwards and towards the coast Filion 1985, Laliberté and, all resulting from climate change, have been observed over the last decades along the eastern coast of Hudson Bay.

Data collection
We established two inventory plots in pure P. glauca stands near the village of WK during the summer of 2018 ( figure 1(b)). The WK1 site represents a forest with strongly scattered small trees and a low canopy cover situated at the upper margin of an isolated forest patch (~treeline sensu Körner, 2012), surrounded by a boulder field and lichen-dominated rock outcrops. Shrubs (namely Alnus viridis var. crispa (Aiton) House, Salix planifolia Pursh and Juniperus communis L.) cover the forest floor at the WK1 site. By contrast, the site WK2 was established in the interior of a large old growth forest with larger dimensions of individual trees and dense canopy closure (closed-canopy forest). The WK2 site mostly lacks shrubs with only a sparse cover of A. viridis var. crispa along the borders of the plot.
At both inventory plots, we marked all living and dead trees taller than 0.5 m with metal tags and trigonometrically measured their positions. To determine the social status (dominance) of each tree, we measured its basic dendrometric parameters (stem diameter at 1.0 m above the ground, tree height, crown diameter and crown base height) and noted information about its health (vitality, crown breaks, defoliation). In total, we measured 97 and 44 trees at WK1 and WK2, respectively. In addition, we recorded the position of each shrub with a crown diameter greater than 2 m (A. viridis var. crispa and S. planifolia) or 1 m (J. communis) together with the height of the highest branch. Finally, we extracted one or two increment cores from all the trees (both living and dead) with a stem diameter greater than 5 cm. We took the cores from opposite sides of the trunk using a Pressler corer at the height of 1 m measured from the stem base. In total, we cored 82 trees at WK1 and 34 trees at WK2 (224 cores in total).

Processing of increment cores and inventory data
Tree-ring cores were processed using traditional dendrochronological approaches (Stokes and Smiley 1996). All cores were glued to wooden supports and the visibility of individual tree-rings was enhanced by sanding the wood surface with progressively finer sandpaper. Afterwards, tree-ring widths of individual cores were measured to the precision of 0.01 mm using TimeTable and PAST4 software (Knibbe 2004). The accuracy of cross-dating was verified by visual comparison (i) of the measured time series with each other (first the cores from single trees, later between trees), (ii) with the pattern of narrow rings described by Caccianiga et al (2008) for our study region, and also (iii) with a chronology of P. glauca from the nearby Castle peninsula (N56°10′, W76°33′), available from the ITRDB database (CANA159, Jacoby, D'Arrigo & Buckley). We applied the same approach also to crossdating samples from dead trees. After the cross-dating, the two cores from each tree were averaged together to obtain growth curves on a tree level. For off-pith cores, the distance between the innermost ring and pith position (pith offset) was estimated using the graphical method (Applequist 1958). The pith offset was divided by mean tree-ring width of the ten innermost tree-rings to produce the number of tree-rings missed near the pith. The estimated number of missing tree-rings was used to assign a cambial age (i.e. the number of years since the formation of the first tree-ring) to each tree-ring of each sample. In addition, the time series of tree-ring widths were converted into cumulative series of increasing stem diameter by cumulating the widths of tree-rings formed before each year (+pith offset). The reconstruction of past tree dimensions was combined with information about their spatial distribution to calculate time series of a modified Hegyi (1974) competition index. This competition index is defined as: where n is the number of neighbour trees within a 11 m radius around the tree under consideration, d iY is the diameter of ith neighbour tree in year Y (cm), d Y is the diameter of the tree under consideration in year Y (cm) and dist i is the distance between the tree under consideration and the ith neighbour tree (m). We modified the original Hegyi index by including the parameter 'Outside', which represents the proportion of a circle area with a radius of 11 m around each tree that lies outside the plot borders. By using this modification, we reduced the risk of underestimating the level of competition for trees situated near the plot borders.
For each site we built a set of seven site chronologies using two different detrending approaches and three types of autocorrelation treatment. We considered splines with 50% frequency cut-off at 30 years and 80 years for detrending (Cook and Peters 1981), and for each of those options, we assembled standard (all autocorrelation retained), residual (all autocorrelation removed) and 'ARSTAN' chronologies (autocorrelation synchronous among series retained, random autocorrelation removed) (Cook 1985). The last chronology was based on raw measurements without detrending and autocorrelation treatment. To select the chronology which was the most relevant for further analyses, we compared the structure of autocorrelation and lagged responses to environmental forcing in raw chronologies with other types of chronologies using ARMA modelling (Box et al 2008). Generally, raw chronologies were ARMA(2, 2) and ARMA (1, 1) processes for the WK1 and WK2 sites, respectively. Standard chronologies based on a 80 year spline had the most similar structure of temporal autocorrelation and lagged responses (table S1 is available online at stacks.iop.org/ERL/14/124024/mmedia) and were therefore used in further analyses.
In addition to site chronologies based on all samples and tree-rings, we also built chronologies including only specific subsets of tree-rings. These subsets were defined based on threshold values of cambial age (below or above 60 years), stem diameter (below or above 15 cm) and the tree competition index (above or below 8 m −1 ). By this, we assembled six additional chronologies for each site (12 in total). This approach was motivated by possibly different growth trends and climatic signal of juvenile/mature, small/ large and suppressed/dominant individuals (Brienen et  Climate-growth analysis Mean monthly temperatures and monthly sums of precipitation for climate-growth analysis were obtained for the station located at WK airport from Environment Canada. This dataset extends back to 1926 for temperature and to 1935 for precipitation; however, it contains numerous gaps (e.g. 1945-1946, 2000s). We filled them using CRU TS4.01 data (Harris et al 2014), which tend to tightly correlate with WK station data for both temperature (mean r = 0.97) and precipitation (0.84). There is mostly non-significant persistence (temporal autocorrelation) in the climatic data (figure S1).
Bootstrapped correlation coefficients were calculated between site chronologies and monthly resolved climatic data over the period of 1935-2017. We included a period of 39 months into the analysis, starting from June three years before tree-ring formation (pppJune) to August of tree-ring formation. The same design of climate-growth analysis was applied also for all chronologies based on subsets of tree-rings. By considering the climate of more than a single year before tree-ring formation in climate-growth correlations, we aimed to quantify the trade-off between secondary growth and seed production. In P. glauca, cone development represents a multi-year process of carbon reserves accumulation (three and two years before seed dispersal), cone initiation (one year before seed dispersal), and seed growth and maturation (the year of seed dispersal; Roland et al 2014). Our effort to quantify the effect of cone production on tree growth was limited by the absence of cone production monitoring in the remote WK region. However, to indirectly assess this interaction, we compared patterns of climate-growth correlations with climate-reproduction models developed for P. glauca in other parts of the North American boreal region, namely Alaska (Juday et al 2003, Roland et al 2014 and the Yukon territory (Krebs et al 2012). Further in the text, those models are referred to as ROL, JUD and KRE, respectively.

Forest structure
The spatial distribution of individual trees and shrubs at both sampling sites is shown in figure 2. Mean tree age equalled 73 and 136 years at the sites WK1 and WK2, respectively. The oldest measured tree-rings date back to 1844 and 1821 at the sites WK1 and WK2, respectively. Individual trees exhibited very slow growth in both plots with mean tree-ring width 0.85 mm and 0.63 mm, respectively. In addition, trees had a larger stem diameter (18.8 cm) and were influenced by stronger competition (14.2 m −1 ) at the WK2 site compared to the WK1 site (14.0 cm and 5.1 m −1 ). For detailed statistics of additional dendrometric parameters, see table S2. Generally, we identified three periods in which temperature had a significant effect on the site chronologies (figure 4). Positive correlations occur during summer of the tree-ring formation year (Jun-Aug) and spring-summer two years before tree-ring formation (ppMar-ppMay). In addition, the negative correlations appear between chronology WK2 and summer temperature in the year preceding tree-ring formation (pJun-pAug). In contrast to temperature, correlations with precipitation are weaker and mostly non-significant, with the exception of negative effect of precipitation in June of tree-ring formation year at WK2 site.
The pattern of temperature signal revealed by climate-growth correlations is synchronized but inverse to that of cone production ( figure 4). Specifically, all the climate-reproduction models under assessment (KRE, JUD, ROL) consistently show (i) a positive effect of temperature during summer one year before cone production (pApr-pAug) and (ii) a negative effect of summer temperature two years before cone production (ppApr-ppAug). In case of tree-ring formation  summer, temperature-growth correlations (positive) are synchronized with temperature-reproduction correlations (negative) derived by the ROL model. However, the opposite response of growth and reproduction to temperature in tree-ring formation year is not supported by the other two models, indicating a positive (JUD) or no (KRE) effect of temperature on cone production.
The pattern of temperature-growth correlations is similar regardless of the detrending option and autocorrelation treatment. The significance of correlations during the three periods described above is not sensitive to the type of chronology at the WK2 site, but the climatic signal of pJun-pAug temperature becomes stronger and significant in WK1 chronologies retaining only high-frequency variability (i.e. all chronologies based on a 30 year spline and residual chronology based on 80 year spline; figure S2).
Chronologies of tree-ring subsets with specific cambial ages, stem dimensions and tree competition status and their climate-growth correlations Individual pairs of standard chronologies based on tree-ring subsets correlate significantly (p<0.01) with each other and their extremes are synchronized (figure S3). Coherency statistics (EPS, Rbar and SNR; sensu Wigley et al 1984) of chronologies based on old, large and dominant trees are higher compared to those for young, small and suppressed trees, and in the cases of Rbar and SNR even higher than respective values for the site chronologies (table S3). Climate-growth correlations exhibit a similar pattern as those for site chronologies in cases of a positive temperature effect during growing seasons of tree-ring formation and two years before tree-ring formation ( figure 5, S4-S5). The negative effect of previous year's summer temperature observed in the WK2 site chronology is not significant in cases of tree-ring subsets including young and small individuals (figure 5).

Discussion
Climatic drivers of Picea glauca growth at the boreal forest-tundra transition Climate-growth correlations of both site chronologies indicate that the growth of P. glauca is driven by a pattern of temperature over three consecutive summers. Specifically, significant correlations appear during summer seasons in tree-ring formation year and up to two years before tree-ring formation. Positive correlations between tree-ring width and spring/ summer temperatures in tree-ring formation year Figure 4. Climate-growth correlations for site standard chronologies WK1 (red) and WK2 (blue) considering the period from June three years before tree-ring formation (pppJun) to August of the tree-ring formation year (Aug). Full bars indicate significant correlations (p<0.05). Rectangles below correlation charts summarize three statistical models of climatic drivers of Picea glauca cone production, with green, orange and white indicating positive, negative and no effects, respectively. The full span of rectangles indicates the periods considered in specific models. The lower/higher intensity of colours for precipitation/temperature-reproduction correlations derived by Roland et al (2014) model reflects their smaller/greater effect size (β-estimates below/above 0.6).
have repeatedly been reported for both boreal and alpine treelines (D'Arrigo et al 2014, 2009, Körner 2012, Babst et al 2019, which represent the cold limits of the existence of trees. Higher spring/summer temperatures promote the formation of wider treerings through the earlier onset, longer duration and higher intensity of cambial activity (Rossi et al 2008).
The significance of correlations between tree-ring width chronologies and summer temperature one year before tree-ring formation (pJun-pAug) is sensitive to the detrending approach at the WK1 treeline site but universal at the WK2 closed-canopy site. Negative correlations in this period seem to be a phenomenon that is highly typical of boreal forests dominated by P. glauca (Wilmking et al 2004, Tardif et al 2008, Andreu-Hayles et al 2011, Porter and Pisaric 2011, Ohse et al 2012, Trouillier et al 2019. In dry continental regions it (i) appears synchronously with positive correlations with precipitation, and (ii) becomes stronger during dry periods, so it was often interpreted as a response to drought stress (Ohse et al 2012). However, when not frozen over, Hudson Bay gives the region a maritime summer climate with large amounts of precipitation ( figure 1(c)). This is the reason why boreal forests around its shoreline are reported not to be limited by drought (Tardif et al 2008, Bhiry et al 2011and our figure 4 for WK region).
We suggest that processes other than drought-limitation must be indirectly involved in the negative effect of previous-year temperature on secondary tree growth. Indirect drivers are supported also by the fact that negative correlations appear only in the case of tree-ring width chronologies and not in other dendrochronological proxies (e.g. maximum latewood density; Barber et al 2000;Andreu-Hayles et al 2011). The tightly synchronized but inverse pattern of climate-growth and climate-reproduction correlations at the closed-canopy site suggests the existence of a trade-off between secondary growth and reproduction. Cone production in P. glauca consists of two crucial physiological processes (sensu Roland et al 2014) in periods for which we observed significant temperature-growth correlations, namely storage of carbohydrates for reproduction (ppMar-ppMay) and cone initialization (pJun-pAug). Although the lack of nonstructural carbohydrates is not a dominant limiting factor of cambial activity (Körner 2015), previous studies have indicated that cone production in P. glauca reduces secondary growth, resulting in the formation of unusually narrow tree-rings (Nicklen et al 2019) or even extreme pointer years (Juday et al 2003). Given the high frequency of seed years (once per 2-4 years) reported from the nearby region of Central Quebec (Rossi et al 2012a), cone production has the capacity to indirectly drive the climatic signal of summers one and two years prior to tree-ring formation. The lower values of correlations of standard chronologies for the WK1 (treeline) site compared to site WK2 (closedcanopy forest) support this interpretation, as trees at both northern and alpine treelines have been documented to have a limited frequency and intensity of sexual reproduction (  Effect of age, size and social status on climatic signal The pattern of climate-growth correlations does not differ between individual age, size and competition classes at the WK1 treeline site. By contrast, only larger, older and dominant individuals show significant climatic signal of temperature during summer one year before tree-ring formation at the WK2 closed-canopy site. In previous studies, the greater climatic sensitivity of larger trees was associated with stronger exposure of crowns to irradiation (Mérian and Lebourgeois 2011), stronger hydraulic constraints on water transport (Carrer and Urbinati 2004), higher water use efficiency (Brienen et al 2017a) or a weaker effect of noise in dominant trees caused by competition (Fritts 1976). In addition, the greater sensitivity of large trees might be associated with the intensity of the growth-reproduction trade-off. Although P. glauca is fertile already from the juvenile stage (Sutton 1969), individuals with trunk diameters smaller than 10 cm (Krebs et al 2012) and those growing in dense stands (Rossi et al 2012b) produce substantially fewer cones and seeds. This trade-off is therefore not as tight in young and small trees, resulting in non-significant climatic correlations in year of cone initiation (one year before tree-ring formation), when carbohydrates stored for cone formation become depleted. This further supports the above described indirect influence of reproduction on correlations between site chronologies and temperature one year before treering formation.
Contrary to climate-growth correlations one year before tree-ring formation, both sites and all size/ dimension/competition classes show similar positive correlations with temperature two years before treering formation. Summers two and three years before cone production are important for the accumulation of carbohydrates later utilized in cone mass growth (Roland et al 2014). Although the amounts of cones produced at the treeline (Johnson et al 2017, Trant et al 2018 and by small and suppressed individuals are low (Krebs et al 2012, Rossi et al 2012b, juvenile individuals of P. glauca are able to produce single cones from the age of four and substantial quantities of cones from the age of 30-35 years (Sutton 1969). Thus, trees producing cones in their juvenile period were probably partly included also in our chronologies based on subsets of tree-rings of young, small and suppressed trees. In such trees, the amount of carbohydrates is limited by their small leaf area and low levels of irradiation caused by their under-canopy position. This increases the importance of growth-reproduction trade-off in periods when carbohydrates for reproduction are accumulated, that is, summer two years before treering formation.
Effect of detrending on the pattern of climatic signal The most prominent variability among different types of chronologies appeared in the significance of correlations with summer temperatures one year before tree-ring formation (figure S2). Retaining only high-frequency variability in the series (by using a 30 year spline and/or pre-whitening series) resulted in strong negative correlations at both sites. Conversely, the effect of temperature one year before tree-ring formation decreased (even to below the significance level at the WK1 site) when chronologies contained also mid-frequency variability (80 year spline). Moreover, at both sites, all variants of chronologies had stronger and statistically significant correlations with previous year's summer temperature if the climatic data were also detrended before the analysis (not shown). This suggests that the process responsible for temperature-growth correlations during the summer before tree-ring formation probably varies from year to year without any long-term trend. This further supports the assumption that the effect of previous year's summer temperature (containing a strong trend; see figure 1(d)) on tree growth is probably indirectly mediated through some tree-internal physiological process.

Conclusions
The course of temperature during three consecutive summers drives intra-annual variability in tree-ring widths of P. glauca in maritime parts of the boreal forest ecosystem. As expected for a cold maritime environment, tree growth is directly stimulated by spring and summer temperatures in the year of treering formation. However, the opposite correlations between secondary tree growth and summer temperature one and two years prior to tree-ring formation observed the most prominently at closed-canopy site highly resemble inverse patterns of climate-reproduction correlations. This suggests that climatically driven reproduction effort indirectly affects tree growth performance. A simultaneously positive and negative effect of summer temperature on tree growth poses a considerable challenge for both predictions of forest responses to climate change and for climatic reconstructions based on tree-ring chronologies of boreal species.
Trees of different size, age and social status growing in specific microsites share a similar pattern of climate-growth correlations but maximize the climatic signal in different periods. More specifically, exclusively old, large and dominant trees at the closed-canopy site respond significantly to summer temperature in the year before tree-ring formation. This illustrates that age/size-related variability in climatic sensitivity must be considered in palaeoclimatological studies in order to maximize the climatic signal and reduce non-climatic noise.