Climate warming differently affects Larix decidua ring formation at each end of a French Alps elevational gradient

Abstract

Key message

The ongoing global warming tends to reduce Larix decidua radial growth and ring wood density at the bottom of an elevational gradient in the French Alps, while it has a less marked effect at the top.

Context

The ongoing climate warming affects most tree species across their biogeographic distribution range. The bottom and the top of an elevational gradient are appropriate environments to observe the effect of climate warming in contrasted temperature conditions.

Aims

To retrospectively study the effect of the ongoing climate warming on Larix decidua trees located in warm and cold conditions, i.e., at the bottom (1200 m) and the top (2300 m) of an elevational gradient, respectively.

Methods

Dendroecological analysis of two groups of Larix decidua trees located at very low (1200 m) and very high (2300 m) elevations. Construction of climatic response curves for annual ring variables. Decomposition of the variance of the response models into their high and low frequency components.

Results

At Briançon, in the French Alps, the increase of the daily maximum temperature during the last 50 years is already almost 3 °C for the March–October period. The results reveal contrasted behaviors at both ends of the investigated elevational gradient. Basal area increment, ring width, and their earlywood and latewood components increase or level off at the top of the gradient, while they all strongly decrease at the bottom. At the bottom, the low frequency warming effect explains the recent ring width decrease, with no influence of the high frequency temperature variation. At the top, both the low frequency warming effect and the high frequency temperature effect explain each about half of the ring width variation encompassed by the model. Latewood density displays the same trend as ring width, while earlywood density decreases at both ends of the gradient.

Conclusion

Such opposed trends of the radial growth time trends between the top (2300 m) and the bottom (1200 m) of the gradient reflects the respectively favorable and unfavorable effects of climate warming at such ends of this Larix decidua elevational gradient. We propose that the strong ring width decrease observed at the warmer bottom announces a coming dieback. The corresponding wood density decrease will affect wood hydraulic properties in a way that is unclear.

1 Introduction

Broad scale forest die-off triggered by more frequent and intense drought associated to the global warming has been widely documented during the recent years (Allen et al. 2015). In general, mountain forests are recognized to be more fragile and more sensitive to disturbances than lowland forests (Tse-ring et al. 2010). The maximum warming associated to climate change is expected to occur in artic regions and in mountain regions at medium to high elevations (Kotlarski et al. 2015; Wang et al. 2014) and thus may jeopardize the future of mountain forests. In some areas of the Alps, temperature already augmented twice as fast as the average increase in the northern hemisphere during the previous century (Rebetez and Reinhard 2008). This current warming rate is greater at upper elevation (Mountain Research Initiative EDW Working Group et al. 2015).

Several recent studies reported climate change–related declines in several European mountain tree species, i.e., Scots pine (Pinus sylvestris L.) in South-Western and North-Eastern Germany (Taeger et al. 2013); mountain pine (Pinus mugo ssp. Uncinata (DC.) Domin) in Western Tatra Mountains, Slovakia (Solár 2013); Norway spruce (Picea abies (L.) H. Karst.) in Austria (Vospernik and Nothdurft 2018); and silver fir (Abies alba Mill.) in Spain, Italy, Romania, and France (Gazol et al. 2015; Cailleret et al. 2014). These declines are generally observed at the bottom and could correspond to a contraction of the warm edge of these species elevational distributions.

Unlike these species, and to our knowledge, Larix decidua (Larix decidua Mill.) showed no particular climate change–related decline so far. Larch is a light-demanding pioneer species with a distribution that benefited from agricultural activities and especially grazing (Heinz H. Ellenberg 1988; Schulze et al. 2007) and growing well on soils with a high water supply (Fourchy 1952; Keller et al. 1997). But one modeling study shows that its distribution is threatened by climate change (Dyderski et al. 2018). Indeed, Larix decidua is regulating water loss during atmospheric drought more efficiently than Picea abies (Peters et al. 2019). It is able to preserve a comparatively high transpiration rate even during periods of moderate water deficit, thanks to its high capacity of water uptake (Anfodillo et al. 1998; Badalotti et al. 2000). However, this species is also known to be nearly as sensitive to soil water deficit as Norway spruce (Vanoni et al. 2016; Schuster and Oberhuber 2013; Lévesque et al. 2013; Charra-Vaskou et al. 2012). Decreased water availability and increased atmospheric demand constrain Larix decidua from low to higher than expected elevations in dry Alpine valleys (Obojes et al. 2018; Schuster and Oberhuber 2013; Eilmann and Rigling 2012): low elevation Larix decidua stands are sensitive to water stress (M. Saulnier et al. 2019). The bottom and the top of the elevational distribution of a mountain tree species are appropriate environments to observe the effect of climate warming in much contrasted temperature conditions. The climate warming tends to impose more intense drought stress to the trees located at the lower end of an elevational gradient, while it tends to improve the growing conditions at the top. Global warming could increase tree radial growth at high elevation, while it would decrease it at low elevation (Primicia et al. 2015; Sidor et al. 2015; Marqués et al. 2016), as already observed for Larix decidua in an Inner Alpine dry valley (Obojes et al. 2018).

In the southwestern part of its Alpine natural area, Larix decidua elevational distribution ranges from as low as 1200 m to as high as 2500 m, making it one of the mountain species with the largest temperature variation along its elevational gradients (Fourchy 1952). Larix decidua is an important forest tree species throughout the Alps and in some regions of central Europe. As one of the very few deciduous conifer species, it confers a very distinctive appearance to the Alpine landscapes (Garbarino et al. 2013). It produces a durable, high quality wood and plays an important ecological, economic, and cultural role all over its natural area (Riou-Nivert 2001). In this study, we retrospectively investigated the annual ring response of adult Larix decidua trees to climatic variation during the last 50 years in the French Alps. We selected three as low and two as high as possible sampling plots in the region of Briançon, France, where we found the longest Larix decidua elevational gradients.

We hypothesize that during the last 50 years, radial growth increased at the top, while it decreased at the bottom of the Larix decidua elevational gradient. We completed the information brought by basal area increment and ring width with the early- and latewood components, and the corresponding density variables of the annual radial growths. Ring density is a proxy of the proportion of cell wall and lumen area in a wood volume, and provides information on wood hydraulic properties (Dalla-Salda et al. 2011). The large anatomical variation between early- and latewood (Zobel 1989) suggests that both parts of the ring play very different hydraulic roles (Zimmermann 1983), as shown in Douglas fir for vulnerability to cavitation (Guillermina Dalla-Salda et al. 2014). Therefore, it seems relevant and interesting to investigate Larix decidua’s response to climate not only at the ring level, but also at the early- and latewood levels. We studied the relationships between these annual ring variables and the corresponding climatic time series from the closest Météo-France weather station. We postulate that the global warming explains the radial growth increase at the top and the radial growth decrease at the bottom of the gradient. Finally, we discuss to what extent the observed variations are possible predictors of future changes.

2 Material and methods

2.1 Sample collection

We collected the wood samples in November 2016 at two of the three low elevation plots and in June 2017 for the two high elevation plots, coupled with the sampling of the third low elevation plot (Rozenberg et al. 2020). All plots are located in the forêt communale of St-Martin-de-Queyrières, near Briançon, Hautes-Alpes, France. The forest is owned by the village of St-Martin-de-Queyrières and managed by the Office National des Forêts (French National Forest Service) since the middle of the nineteenth century. It is a North, North-West facing slope ranging from 1200 to 2500 m and covered mostly by pure uneven-aged larch forests with some patches of Abies alba, Pinus sylvestris, Pinus cembra, and Pinus uncinata. The forest is mostly naturally regenerated with occasional plantations of local origins in case of regeneration failure. The soil types found along the slope are more or less complex succession and associations of Calisol, Eutric Brunisol, Colluviosol, and Regosol. The forest is subjected to low intensity periodic thinning done according to the management plan (Table 1).

Table 1 Location of sampling sites

We collected 144 5-mm breast height increment cores with a power-driven increment borer, following a perpendicular to the slope orientation to avoid compression wood. The trees were all adult, dominant trees with no visible major defect. The freshly collected increment cores were stored in polycarbonate honeycomb boxes, brought back to the laboratory, where they were oven-dried and sawn to a 2-mm uniform thickness. Then, they were resin-extracted with a 48-h pentane bath. They were dried again and then X-rayed using the indirect procedure invented by Polge in 1966 and continuously improved since then (Polge 1978; Mothe et al. 1998). Finally the X-ray films were scanned at 4000 dpi to acquire the microdensity profiles using the WINDENDRO software (Windendro 2008e Regent instruments Canada nc.; http://www.regent.qc.ca/products/dendro/DENDRO.html).

The microdensity profiles were visually checked and compared with the sawn samples then cross-dated using pointer years. Some profiles had to be discarded because ring boundaries were blurred by a very strong grain angle. When the increment core missed the pith, we estimated the number of missing rings using the curvature of the innermost rings of the increment core. Then, tree age was estimated by arbitrarily adding 10 years to the number of rings counted from pith. We assumed that 10 was the average number of years taken by the trees to reach the increment core collection height, breast height.

2.2 Annual radial growth variables

The raw microdensity profiles were imported into R (R Core Team 2018), and the following conventional annual radial growth variables were calculated using specially designed R functions, i.e., ring width, basal area increment (ring area, BAI), earlywood width, latewood width, earlywood density, latewood density, and standard deviation of ring microdensity (ring standard deviation). Earlywood-latewood boundary was calculated using the conventional extreme average method based on the mean of the extreme ring density values (Rozenberg et al. 2002). Basal area increment was calculated as${\mathrm{BAI}}_i=\pi \times \left({r}_i^2-{r}_{i-1}^2\right)$, where BAI_i is the basal area increment of the ith ring from the pith, and ${r}_i^2$ and ${r}_{i-1}^2$ are the radius of the ith and (i − 1)th rings from the pith. The radii are calculated as cumulated ring width. BAI is a dendroecological approach based on the assumption that it expresses a constant growth distributed all around the growing stem (Biondi and Qeadan 2008).

2.3 Detrending

Our objective is to ensure that the main trends in our results and conclusions are not biased but robust. Besides the conventional ring variables based on the earlywood-latewood model (Ruiz Diaz Britez et al. 2014), we introduced BAI. Endogenous and exogenous factors may affect radial variation of any annual radial growth variable. To study relationships with climate, we need to keep as far as possible the climatic signal in our time series, but eliminate other possibly confounding sources of variation. We especially need to eliminate the low frequency cambial age effect, which is generally assumed to be the major confounding effect, but keep the low frequency climatic effects like the possible uniform effect of climate change and global warming (Cook et al. 1990). We tested separately each non-climatic effect available (cambial age, Larix decidua bud moth defoliation, and heartwood-sapwood), and we adjusted from the effects that proved to affect significantly the radial variation of the ring variables. Larvae of larch bud moth (Zeiraphera griseana) episodically defoliate alpine stands of Larix decidua during considerable outbreaks striking most of the time at 8–10-year intervals, with a maximum at 1800 m altitude (Baltensweiler et al. 2008; Dormont et al. 2006). The defoliation reduces photosynthesis, decreases carbohydrates production, and lowers annual ring growth. Such reduction of annual ring growth leaves a typical mark in annual ring time series (Battipaglia et al. 2014; Mélanie Saulnier et al. 2017; Cerrato et al. 2019). To separate bud moth outbreak-induced growth decreases from those produced by climatic factors or by other biotic factors, we used a 15-year moving window approach (Büntgen et al. 2009). We tested the heartwood-sapwood effect with a one-way ANOVA. We found that only the cambial age significantly affected radial variation. A one-way ANOVA showed that the cambial age effect was significantly different between the elevational classes (Fig. 4). Then, we adjusted linearly all the ring variables from the cambial age effect estimated at both elevational levels by regional curve standardization (RCS) procedure (Esper et al. 2003). We compared three methods, i.e., raw (no detrending), residual RCS, and ratio RCS (Esper et al. 2003). The terms residuals and ratio refer to the terminologies used in Esper et al. (2003). We obtained age-related curve for both elevations (Fig. 4). For each individual time series (or individual tree), we obtained an adjusted value of tree-ring growth by subtracting (residuals) or dividing (ratio) the observed value by age-depended estimation for each year. We obtained the residual and ratio RCS chronologies by time-depended averaging of individual series. Figure 5 shows raw, residual, and ratio detrended time series for ring width at low and high elevations. We compared the results obtained with the three types of times series, the two adjusted and the one non-adjusted. For ring width, earlywood width, and latewood width, we found very similar results that lead to mostly identical conclusions with the three methods (Fig. 6). For BAI, the results obtained with the raw and residual RCS time series were similar, while the ratio RCS method appeared to be obviously not adapted. The possible bias associated to the ratio RCS applied to specific ring variables has already been reported (Cook and Peters 1997). For ring, earlywood, and latewood density, the results obtained with the raw and residual RCS time series were similar, while the ratio RCS time series were apart and probably biased. The corresponding results are not shown in Section 3, since the objective of the article is not to work on method development but on larch response to climate. A small selection of these methodological results is provided in the Appendix section.

The cambial age effect shown in Fig. 4 is high near the pith and gradually decreases with the ring number counted from the pith (Zobel and Sprague 1995). In this study, we used the last 50 rings counted from the bark (from 1967 to 2016), in trees that are generally much older than 70 years (estimated age is between 89 and 106). It means that in our time series, the cambial age effect is certainly lower than the climatic effect. In Section 3, we present the results based on the residual RCS detrended time series for both the conventional ring variables and BAI, knowing that in this study, both raw and residual RCS methods lead to the same conclusions.

2.4 Climatic data

Climatic data are daily values of minimum and maximum temperatures and precipitations coming from the closest Météo-France weather station (French national weather service). This weather station was first located in Briançon (N 44.905856, W 6.633046, 1306 m) and started to record data in 1966. In 2003, a new station was created and started to record climatic data in Villard-St-Pancrace (N 44.864861, W 6.61541, 1276 m), 4.8 km south of the Briançon weather station. From 2004 to 2010, both weather stations have been operating together. In 2010, the Briançon weather station stopped operating. Only the Villard-St-Pancrace weather station continued to record climatic data since then. The climatic time series during the Briançon-Villard-S-Pancrace common recording periods show that there is nearly no difference between both sites for maximum and minimum temperatures and very little for the precipitations (Fig. 9).

We grouped the data from the two weather stations to create a unique climatic database going from 1966 to 2017. We used this database to study the annual variation of the three climatic variables available from 1966 to 2017. For each climate variable, we tested a large number of intra-annual periods, and we selected the annual periods for which the climate-time relationships were maximum (Fig. 8). The corresponding models describe the effect of climate change on the three climatic variables available, i.e., minimum and maximum temperature, and precipitation in the Briançon region.

2.5 Variation of ring variables against time

The mean adjusted annual radial growth variables and BAI were plotted against time (years) at each elevational level for the period 1945–2016. We estimated the Spearman correlation coefficient (rho) among the ring variables for the same 1945–2016 period. We chose the Spearman non-parametric rank correlation instead of the Pearson, because the first method is more flexible for data that are not necessarily normally distributed, or for relationships which are not strictly linear, such as with time series data.

We tested the relationships between climate and annual radial growth variables at each elevational level using also the Spearman correlation coefficient, using for the ring variable time series from 1966 to 2016, corresponding to the climatic variable time series. While we found no significant relationships with precipitation, we found significant relationships with both temperature variables and stronger relationships with maximum than with minimum temperature. Thus, we decided to do the dendroecological analysis with maximum temperature only.

We created two maximum temperature time series going from 1966 to 2017 for the two elevational classes (low and high elevation) by linearly adjusting the temperature time series of the elevational effect as in Latreille et al. (2017). The parameters of the elevational temperature model are based on temperature data recorded from 2008 to 2016 in four INRA field weather stations installed in four elevational plots located a few kilometers away in the forest of Villard-St-Pancrace, at 1350, 1700, 2000, and 2300 m elevation.

2.6 Data analysis

We tested the variation with time (annual periods) of minimum and maximum temperature from 1966 to 2017 with linear regression analysis for varying annual periods. The duration and dates of the tested annual periods varied from 2 to 12 months, starting from January 1 to December 31, respectively, with a 1-week step variation from one tested period to the other. For each climatic variable, we selected the period maximizing the coefficient of correlation of the time-trend relationship (“maximizing period”). We estimated the confidence interval at 95% of the beginning and ending dates of the maximizing period.

The relationships between annual radial growth variables and climate were tested at each elevational level using simple linear regression. The effect of low and high frequency climate variation was tested using the analysis of variance (ANOVA):

$$ {x}_{ijk}=\mu +{y}_i+{z}_j+{\varepsilon}_{ijk} $$

where x_ijk is the annual radial growth variable, μ is the mean, y_j is the year variable, z_j is the annual climate variable, and ε_ijk is the residual. The part of the variation explained by the linear year (y_j) effect is the low frequency warming effect (Fig. 1), while the part of the variation explained by the annual climatic variable (maximum temperature here) is the remaining high frequency effect.

Fig. 1

Chart of the strongest significant linear relationships between minimum and maximum temperature, and precipitation (rainfall and snow) across time (year) for the maximizing annual periods based on the Météo-France weather data of Briançon (1306 m) and Villard-St-Pancrace (1276 m). The equations and adjusted R-squares of the corresponding linear regressions are printed below each chart

We used this ANOVA model to decompose the model variance into its two components, the global warming (GW) effect and the high frequency temperature (HFT) effect. We used the R gvlma package to test and validate the linear model assumptions.

3 Results

3.1 Climate change

There are significant linear time trends for the three climatic variables. These relationships are quite strong for minimum temperature (Adj. R-square = 0.48), strong for maximum temperature (Adj. R-square = 0.62), and low (but significant) for precipitation (Adj. R-square = 0.09). The periods maximizing the relationships are very similar for minimum and maximum temperature (from the beginning of March to the end or the middle of October, respectively), while it is very different for precipitation (from the beginning of January to the beginning of April). They correspond to a 0.0281-°C and 0.058-°C annual increase for minimum and maximum temperature, respectively, and to a significant but very low 0.0154 mm annual decrease for precipitation.

3.2 Annual radial growth variable time series

From 1945 to around 1985, adjusted BAI increases in a parallel way at both low and high elevations, with BAI growth being faster at the low elevation. Then after 1985, low elevation BAI starts to level off and finally decreases after 2000. During the same period, high elevation BAI continues to increase. Adjusted ring width is mostly stable at the high elevation except for a marked decrease during the early 1970s, while at the low elevation, it steadily decreases from the beginning (1945) to the end (2016) of the study period. Earlywood width and latewood width show mostly similar trends (Fig. 6). The high and low elevation lines cross during the 1990s: low elevation radial growth becomes lesser than high elevation after this period. Both curves merge and cross between 1995 and 2005. Earlywood density is mostly stable from 1945 to the beginning of the 1990s at both elevations, with the low elevation earlywood being generally higher than high elevation. Then, both high and low elevation earlywood density strongly decreases from the beginning of the 1990s to 2005, where they seem to mostly stabilize again until 2016. Latewood density show a different trend, with high elevation latewood density being mostly stable or slightly increasing, and low elevation latewood density showing a trend very similar to low elevation earlywood density. Ring density variation is intermediate between that of earlywood and latewood densities (Fig. 6). From 1945 to approximately 1985 (earlywood density, Fig. 2), 1990 (ring density, Fig. 6), and 2000 (latewood density, Fig. 2), density is lower at the high elevation than at the low elevation. This tendency strongly changes afterwards, temporarily for earlywood density, more clearly for ring density, and considerably for latewood density (Figs. 2 and 6).

Fig. 2

Average variation across time of four adjusted annual radial growth variables at the low (red line) and high (blue line) elevation. Light colors represent the 95% confidence interval around the mean. The red and blue numbers are the number of trees used in the calculation of the mean at each year and elevational level (same color code)

3.3 Relationships among the ring variables

At the low and high elevations, ring width, and BAI are quite strongly positively correlated (0.82 and 0.74, respectively). Both variables are positively correlated with latewood density, with this correlation being stronger at the low elevation, especially for ring width. The sign of the moderate correlation between earlywood density and ring width changes with elevation: it is positive at the low and negative at the high elevation. The same trend is observed for BAI but with a much weaker relationship at the low elevation (Table 2).

Table 2 Spearman’s correlation coefficients among the four selected ring variables at both elevation levels. All correlation coefficients are significant at the 0.001 (***), 0.01 (**), or 0.05 (*) level, except when indicated (NS)

3.3.1 Relationships between annual radial growth variables and annual temperature

Table 3 and Fig. 3 show the strongest significant relationships that we found between the annual radial growth variables and mean maximum daily temperature of the period maximizing the relationship (the annual period maximizing the coefficient of correlation of the time-trend relationship). We found similar but weaker relationships with minimum temperature, and no significant relationship with precipitations. Table 3 shows the results for all the ring variables, while Fig. 3 shows the same results for a selection of four variables (basal area increment, ring width, earlywood density, and latewood density).

Table 3 Correlations (Spearman’s correlation coefficient) between the ring variables and the mean maximum temperature of periods selected for maximizing the ring-climate relationship.

Fig. 3

A selection of relationships between four ring variables (basal area increment, ring width, earlywood, and latewood densities) and mean maximum temperature of the selected maximizing annual periods at both low and high elevation. Rho is the Spearman correlation coefficient and p is the associated probability. p = 0 corresponds to significance levels below 10⁻¹². The gradient of gray of the year points corresponds to between-year time variation from 1967 (dark gray) to 2016 (light gray)

Table 3 shows that it is possible to find very strongly significant relationships with maximum temperature for all the ring variables at both elevations. There are significant relationships for a large number of periods and sometimes for all the tested period. This is the case, for example, for latewood density at the low elevation (Fig. 7). The maximizing periods are different between variables and elevations. At the low elevation, all the relationships are negative and generally quite strong to very strong (Spearman’s coefficient of correlation ranging between − 0.72 and − 0.86 except for basal area increment, − 0.50). At the high elevation, all the relationships are quite strong to strong (from 0.60 to 0.90), stronger than at the low elevation, with six positive relationships and two negative (for mean ring density and earlywood density).

Figure 3 shows the corresponding relationships for the four selected variables. Patterns of earlywood width and latewood width are very similar to that of ring width, while the pattern of ring density is very similar to that of earlywood density (not shown). The pattern of ring density standard deviation is very similar to that of latewood density (not shown). All the relationships are linear.

3.4 Variance decomposition based on the ANOVA model for four selected variables

The results in Table 4 show the variance decomposition based on the ANOVA model for the low frequency “global warming” effect (the part of the variation explained by the linear relationship between annual temperature and year) and the remaining high frequency “annual temperature” effect. There are strong differences between elevations and ring variables. The global warming effect explains most part of the variation at the low elevation for ring width and early- and latewood density, and at the high elevation for BAI and earlywood density. For ring width and latewood density at the high elevation, and for BAI at the low elevation, both global warming and high frequency effects explain analogous parts of the total variance.

Table 4 Variance decomposition based on the ANOVA model. The values are the proportions of the total variance explained by each effect (low frequency “Global warming” and “high frequency temperature”) and by the complete model (“Total,” equal to model R-square)

The total R-square of the models in Table 4 are of the same magnitude order than those of the models presented in Fig. 3, with small differences due to the fact that not exactly the same annual values or the same models are used in both analyses. For example, for ring width at the high elevation, R-square is 0.40 in Fig. 3, while it is 0.42 in Table 4. Nevertheless, both analyses show the same trends with similar between-traits and between-elevation level differences.

4 Discussion

The climate dataset of the Briançon, Villard-St-Pancrace Météo-France weather stations shows that the average temperature of a large part of the year including the Larix decidua growing season (March–October) linearly increased between 1967 and 2017 (Fig. 1). Two warm years (1966 and 1969) preceded the cold period from 1970 to 1980. The temperature of these two warm years is nevertheless at the same level than the temperature of the coldest years after 1998. The succession of cold and relatively cold years from 1970 to 1980 reduced the radial growth of the corresponding annual rings in the limiting conditions of the high elevation plots, but not in the warmer conditions of the bottom of the gradient (Fig. 2). For the whole study period, the maximum temperature increase exceeded 3 °C (0.58 °C by 10 years during the last 50 years). According to the French national weather service Météo-France, the average increase of maximum temperature is 0.32 °C/10 years for France (Impacts du changement climatique: Atmosphère, Températures et Précipitations 2020). This result is consistent with the faster warming observed in the Western Alps than in lower elevation regions for maximum temperature (Rebetez and Reinhard 2008). For minimum temperature, the values are approximately the same: 0.28 °C/10 years at Briançon, Villard-St-Pancrace, and 0.29 °C/10 years for France. A similar trend was recently observed in Switzerland (Vitasse et al. 2018).

The variation of the annual ring time series during the same period (Fig. 2) shows very different temporal trends at the low (under 1300 m) and high (about 2300 m) elevations. At the low elevation, ring width decreases from more than 2.5 mm to less than 1 mm, while at the high elevation, it remains stable around 1.5 mm. BAI better represents the overall tree growth potential than ring width and its components (Biondi and Qeadan 2008), which nevertheless remain pertinent at least for comparison with previously published results.

Earlywood and latewood densities display different trends reflecting distinct determinisms (Fig. 2). The pattern of mean ring density (Fig. 6) is very similar to that of earlywood density. Before 1990, mean earlywood density was 0.47 kg/dm³ (standard deviation 0.016) at the high elevation and 0.50 (0.016) at the low elevation. After the 2000s, it decreased to respectively 0.39 (0.016) and 0.42 (0.016) kg/dm³. The relatively uniform latewood density trend suggests no global warming effect on this trait at the high elevation. In contrast, the strong latewood density decreases at the low elevation from the 1980s to the beginning of the 2000s (Fig. 2, from 0.82 (0.016) to 0.68 (0.016) kg/dm³) indicates a major global warming effect in the warmer conditions of the bottom of the Larix decidua distribution. Another obvious consequence is the alteration of the elevational trend of wood density. During the period preceding the observed warming, wood density tends to be lower at the higher than at the lower elevation. This well-known trend (Zobel 1989; Rossi et al. 2014; Apiolaza 2012) is reversed afterwards, especially for latewood density (Fig. 2; see also Fig. 6).

The strong and continuous radial growth decrease at the low elevation for all variables (ring, early- and latewood width, and from 1990 for BAI) is consistent with predictions (Primicia et al. 2015; Sidor et al. 2015; M. Saulnier et al. 2019) and with observations of different species including Larix decidua (Obojes et al. 2018; Jiao et al. 2019; Zhu et al. 2018; Marqués et al. 2016). However, we find contradictory trends for the high elevation radial growth between BAI and ring width, even if partly explained by the geometrical relationship between both traits. While ring width remains rather stable all across the study period, high elevation BAI increases and thus follows the predicted (Primicia et al. 2015; Sidor et al. 2015; M. Saulnier et al. 2019) and sometimes observed (Bošeľa et al. 2014; Zhu et al. 2018; Yin et al. 2018) global warming–related growth trend.

A radial growth decrease often precedes forest tree mortality (Cailleret et al. 2017, 2019). Therefore, the growth decrease observed at the low elevation could be interpreted as a warning signal for coming diebacks related with increased drought stress during the growing season. In our case, a wood density decrease also accompanies the observed low elevation growth decline. In general, a radial growth drop-off is related to a density increase (Zobel 1989). Nevertheless, the summer 2003 heat wave and drought in France induced in Douglas fir a mean density decrease of the 2003 ring explained by a smaller proportion of lower density latewood (Martinez Meier et al. 2008). To our knowledge, few results are available concerning variation with time of ring density variables in forest mountain species. Wood density of four dominant European species decreased by 8 to 12% during the last 100 years (Pretzsch et al. 2018). Several studies suggest that the global warming tends to decrease minimum ring and/or earlywood density (Camarero and Gutiérrez 2017; Fonti and Babushkina 2016). In our case, earlywood (at both elevations) and latewood (at the low elevation only) densities strongly decreased during the study period, with most part of the decrease observed during a 20 to 25-year period between the mid-1980s and the mid-2000s. Then, the three ring density variables seem to level off after 2005. Only latewood density at the high elevation remains more or less stable all along the study period (Fig. 2). Water stress is known to affect ring density: for example, a major heat wave induced the formation of a narrower, less dense ring with incomplete latewood in Douglas fir (Martinez Meier et al. 2008). Wood density is mainly related with the proportion of cell wall and lumen area in the xylem. It is generally accepted that the variation of wood density has an effect on conduction properties and on resistance to drought (Dalla-Salda et al. 2011). During the growing season, stomata of anisohydric species like Larix decidua remain open even when temperature increases up to critical levels. Because of its high water uptake capacity, Larix decidua can preserve a rather high transpiration rate even during periods of limited water deficit (Anfodillo et al. 1998). Transpiration peak, water availability fall, and carbohydrates rerouting to vital functions rather than to tracheid cell wall thickening could explain the observed density decrease.

The distribution of the tested periods for the relationship between latewood density and temperature at the low elevation is shown in Fig. 7. The maximum correlation coefficient in this distribution corresponds to the selected response curve. The heat map in Fig. 7 and all the heat maps for the other ring variables (not shown) reveal that there is a wide range of periods generating significant correlations nearly as high as the maximizing one. A variation of more or less 50 to 100 days of the beginning and ending dates of the maximizing periods does not decrease much the value of the correlation coefficient. In other words, the average temperature of a wide range of periods explains more or less the same, large part of the inter-annual variation of ring characteristics. With Fig. 3, we show the selected temperature response curves for BAI, ring width, earlywood density, and latewood density. The earlywood width trend is very similar to that of ring width, and the ring density trend is very similar to the earlywood density one (Fig. 6). Other studies of climate-growth relationships based on monthly mean temperature seem to have always found lower values of correlation coefficients than in this study (M. Saulnier et al. 2019; Danek et al. 2018; Carrer and Urbinati 2004; Rolland et al. 1998). All the relationships are obviously linear, except maybe the latewood density low elevation relationship, which shows a relatively clear negative sigmoidal trend, not tested in this study.

In Table 4, we split the explained variation for the ring variables into two components, i.e., the low frequency year effect, which is equivalent to the linear global warming effect highlighted in Fig. 1, and the high frequency residual effect, equivalent to the part of the between-year temperature variation not explained by this linear global warming effect. We show that the weight of these two components is very different between both elevation levels. The global warming effect alone explains most part of the annual ring linear variation for all the traits except BAI at the low elevation, while it explains much less at the high elevation. More generally, it means that the global warming alone explains all or nearly all of the strong ring width and density linear decrease at the warm bottom of the Larix decidua elevational gradient, as well as the earlywood density linear decrease and the BAI linear increase at the top. It means that the global warming is the main driver of larch radial growth decrease during the last 50 years at the warm bottom of the study site. For ring width, latewood width, and latewood density at the high elevation, about half or more of the variation is explained by the high frequency year effect. In this cold-limiting environment, we hypothesize that the high frequency inter-annual effect of climate partly masks the potentially favorable effect of the global warming. Studies at larger geographical scales are necessary to estimate the geographical extension and the variation of such trends.

At the species level, wood density is generally positively related with resistance to drought (Hacke et al. 2001). Nevertheless, the same relationship becomes more complex at the intra-species level and when wood density is decomposed into its annual ring components. Early-, latewood, and even the transition between early- and latewood may play different roles in resistance to drought (Dalla-Salda et al. 2011; Dalla-Salda et al. 2014). Wood microdensity depicts the proportion of cell wall in a small (between micrometers and millimeters) wood volume. A wood density decrease can be explained by a cell wall width decrease, a lumen area (and a tracheid diameter) increase, or a combination of both. For example, according to a multi-site multi-species study in the northern hemisphere, tracheid size has the main effect on earlywood density, while cell wall width determines latewood density (Bjorklund et al. 2017). The consequences on wood hydraulic properties are different. Thus, it is difficult at this stage to infer precisely the effect of a wood density decrease on Larix decidua hydraulic properties and on its resistance to drought. Should we expect coming to a hydraulic failure at the bottom of the Larix decidua elevational distribution? Direct studies using resistance to drought traits like resistance to cavitation could bring light to the consequences of the observed wood density decrease.

The main environmental variable consistently decreasing with elevation is temperature. However, many other environmental variables vary with elevation. Among them, soil properties sometimes play a very important role on water availability and could explain part of the observed behaviors.

Can we do something to lower the risk of sub-vitality or even mortality? At short term, reducing competition is a way to mitigate the effect of drought (Gleason et al. 2017). At longer term, a possible solution is to stimulate substitution by assisted migration of more drought-adapted origins or species. Up to which elevation will the lower end of the Larix decidua distribution climb? Larger geographical studies along continuous elevational gradients are necessary to identify the critical limit, remembering that this elevational threshold will ascent with global warming.

5 Conclusion

The warming trend observed in the Briançon region is faster than in lowland areas and faster for maximum than for minimum temperature. It explains most part of the growth and density decrease observed at the warm bottom of the studied Larix decidua elevational gradient. Such recent, fast, and marked decrease can be interpreted as a warning signal for coming diebacks related with increased drought stress during the growing season. A parallel growth increase is observed for BAI at the cold end of the Larix decidua elevational gradient, where global warming is supposed to improve growing conditions, but not for ring width and its earlywood and latewood components. BAI and ring width seem to bring relevant complementary information.

Appendix

Fig. 4

Cambial age time series for ring width

Fig. 5

Comparison of the raw and detrended time series for ring width and BAI

Fig. 6

Same as Fig. 2 of the article for the other four ring variables (earlywood and latewood width, mean ring density, and ring density standard deviation). Average variation across time of four adjusted annual ring variables at the low (red line) and high (blue line) elevation. Light colors represent the 95% confidence interval around the mean. The red and blue numbers are the number of trees used in the calculation of the mean at each year and altitudinal level (same color code)

Fig. 7

Heat map and level curves of the correlation coefficients of the response curves obtained for the tested periods for latewood density (Df) at the low elevation level

Fig. 8

Same as Fig. 3 in the article for the other four ring variables (earlywood and latewood width, mean ring density, and ring density standard deviation) and mean maximum temperature of the selected maximizing annual periods at both low and high elevation (periods’ table and the relationships for the other ring variables are in the Appendix section). Rho is the Spearman correlation coefficient and p is the associated probability. p = 0 corresponds to significance levels below 10–12. The gradient of gray of the year points corresponds to between-year time variation from 1967 (dark gray) to 2016 (light gray)

Fig. 9

Météo-France Briançon and Villard-Saint-Pancrace climatic time series during the 2004-2010 period