Contrasting patterns of tree species mixture effects on wood δ13C along an environmental gradient

Establishing mixed-species stands is frequently proposed as a strategy to adapt forests to the increasing risk of water scarcity, yet contrasted results have been reported regarding mixing effects on tree drought exposure. To investigate the drivers behind the spatial and temporal variation in water-related mixing effects, we analysed the δ13C variation in 22-year tree ring chronologies for beech and pine trees sampled from 17 pure and mixed pine–beech stands across a large gradient of environmental conditions throughout Europe. In the pure stands, average δ13C values were lower for beech (−27.9‰ to −22.2‰) than for pine (−26.0‰ to −21.1‰), irrespective of site conditions. Decreasing SPEI values (calculated over June to September) were associated with an increase in δ13C for both species, but their effect was influenced by stand basal area for pine and site water availability for beech. Mixing did not change the temporal constancy of δ13C nor the tree reaction to a drought event, for any of the species. While the mixing effect (Δ δ13C = δ13C pure stands − δ13C mixed stands) was on average positive for beech and non-significant for pine across the whole gradient, this effect strongly differed between sites. For both species, mixing was not significant at extremely dry sites and positive at dry sites; on moderately wet sites, mixing was positive for beech and negative for pine; at sites with permanent water supply, no general patterns emerge for any of the species. The pattern of mixing effect along the gradient of water availability was not linear but showed threshold points, highlighting the need to investigate such relation for other combinations of tree species.


Introduction
Both the frequency and intensity of drought events are expected to increase in the Northern Hemisphere in the upcoming decades as a result of climate change (IPCC 2013). The impacts on tree functioning due to water shortages are an important concern to foresters. Such impact includes hydraulic failure leading to tree mortality (McDowell 2011), carbon starvation due to stomatal closure (McDowell 2011), increased fine root mortality or changes in fine root dynamics (increase in fine root productivity, Meier and Leuschner 2008). In addition, drought-induced stress increases the risk of wildfires and insect attacks that may lead to tree mortality (Schlesinger et al. 2016).
Favouring mixed-species stands has been proposed as a strategy to adapt forests to the increasing risk of water scarcity. However, quite contrasting results have been reported regarding the effects of mixing tree species on soil water availability and drought exposure of the trees. For instance, Lebourgeois et al. (2013), Pretzsch et al. (2013) and Anderegg et al. (2018) found an improved drought response in more diverse plots while Bosela et al. (2018) and Vanhellemont et al. (2019) found the opposite effect. Complementarity, which is a major determinant of mixing effects (Ammer 2019), is expected to depend both on the pool of species in presence and on the set of environmental conditions (Forrester and Bauhus 2016). A prerequisite Communicated by Rüdiger Grote.

Electronic supplementary material
The online version of this article (https ://doi.org/10.1007/s1034 2-019-01224 -z) contains supplementary material, which is available to authorized users. for synergistic interactions to occur is therefore to admix species with distinct traits in terms of physiology, phenology or morphology. When water is the resource of interest, difference in rooting pattern (Schume et al. 2004), canopy interception (Augusto et al. 2015), transpiration and water use efficiency (Gebauer et al. 2012) and phenology (Forrester 2015) between associated species are of special relevance. For a given set of species, the environmental conditions are the other major determinant of tree species mixture effects. They first influence the expression of the species traits. For example, soil constraints such as hypoxia (Kozlowski 1986), compacity (Greacen andSands 1980) or bedrock depth (Balneaves and De La Mare 1989) are well known to strongly limit the rooting depth of most tree species, leveling off any potential root stratification. A second way environmental conditions may alter complementarity is through their impact on the level of the target resource over space and time. Complementarity effects are thus expected to increase along a gradient of decreasing resource availability if the species mixture improves the availability, uptake or use efficiency of the limiting resource (Forrester 2014). Hence, while mechanisms influencing water availability are of potential relevance for mixing effects in drier sites, we would not expect to find strong species mixture effects related to water in wet sites. Variations in species mixture effects can also occur due to temporal changes in environmental conditions (del Río et al. 2014;Forrester 2014), for example in the case of an annual drought (Grossiord et al. 2014;Lebourgeois et al. 2013;Pretzsch et al. 2013). It is worth noting that temporal and spatial variation in species mixture effect might not be independent of each other. For instance on sites that are moist on average, complementarity with respect to water, while absent in normal years, might operate when water becomes limiting due to a drought event. It is also important to note that, although complementarity in water-related processes would probably be absent on such sites in average conditions, complementarity interactions could be at play for other resources (e.g. light). This means there is room to switch between different types of complementary interactions, depending on the temporal fluctuations in environmental conditions. However, the possible tradeoffs between those mechanisms as well as their respective temporal dynamics remain largely unknown.
In this context, this paper aims to investigate the drivers behind the spatial and temporal variation in the waterrelated mixing effects. For that, we use a well-documented network of pure and mixed stands of pine and beech trees distributed over a large gradient of climate conditions in Europe (Pretzsch et al. 2015;del Río et al. 2017;Dirnberger et al. 2017;Heym et al. 2017). Fagus sylvatica L. and Pinus sylvestris L. are two species which present contrasting traits relevant to water-related processes: potential root system shape (heart shaped for beech and tap root for pine), shade tolerance (high vs. low, respectively), stomatal density (around 200/mm −2 vs. 84/mm −2 ) and strategies in regard to drought resistance (anisohydric tendency vs. isohydric) (Cochard 1992;Martínez-Vilalta et al. 2004;Pflug et al. 2018;Schäfer et al. 2017;Yang et al. 2016). Based on the same network of plots, the mixture with pine and beech at the community level has been shown to have both a higher productivity (Pretzsch et al. 2015) and a higher temporal stability of productivity ) than expected from the corresponding pure stands. Complementarity with respect to water resources has been proposed as a possible mechanism to explain those species mixture effects.
We used the δ 13 C signal as an indicator of water use efficiency (Farquhar and Richards 1984) and analysed its variation in 22-year tree ring chronologies from trees of each species sampled in pure and mixed pine-beech stands.
We hypothesized that the spatial and temporal patterns in wood δ 13 C values are influenced by water availability for both species, with higher δ 13 C under harsher conditions. We further expected to observe lower wood δ 13 C values in the mixtures compared to the pure stands, with the mixing effect increasing under drought stress or as the sites are dryer. Finally, we anticipated an overall increased stability of the δ 13 C signal in the mixtures.

Study area and site/stand characteristics
The tree data used in this study come from 17 pine-beech triplets established under the COST Action FP1206 EuMIX-FOR (European Network on Mixed Forests). This network covers a large gradient of environmental conditions within the overlapping natural ranges of pine and beech (Fig. 1). For this study, 17 triplets were used. Each triplet consists of pure stands of pine and beech and a mixed stand of both species located in similar conditions. The stands are mostly even-aged and mono-layered, and are typical of the (near-) mature development stage. No silvicultural activities had been conducted in the stands during the preceding decade. A standard protocol for tree data collection (diameters, heights of trees and crown bases) and tree coring was applied. The full measurement protocol was described in detail by Heym et al. (2017). Site characteristics for each stand are presented in Table S1.
The triplets (Fig. 1) cover a wide range of environmental conditions and stand productivity classes (Table S1 and Fig. 2). Elevations range from 20-1339 m; mean annual precipitation (P) from 556 to 1175 mm; mean annual temperature (T) from 6 to 10.5 °C; and the de Martonne index (M = annual precipitation (mm)/(mean annual temperature (°C) + 10); de Martonne 1926) from 29 to 67. The Site Index (height of the 100 largest diameter trees of that species per hectare in pure stands at age 50 years; Forrester et al. 2017a, b) ranges from 11.7 to 27.6 m for beech and from 9.5 to 26.9 m for pine.
Stand age and basal area are provided in Table S2. In the mixtures, the percentage of basal area represented by pine ranged from 33% to 74%; total basal area ranged from 30 to 79 m 2 ha −1 ; the total number of trees per hectare from 248 to 2,421; and stand age from 40 to 130 years.

Isotope data
In this study, we used the wood carbon isotopic composition, which reflects the balance between the CO 2 diffusion and carboxylation rates, as it is proportional to the ratio of internal CO 2 to atmospheric CO 2 concentrations (Farquhar et al. 1982). Various environmental factors such as light, temperature and water availability do have an influence on δ 13 C, and their relative importance depends on Fig. 1 Distribution of the triplets across Europe and natural distribution range of European beech and Scots pine according to EUFORGEN (www.eufor gen.org) environmental conditions. Under limiting soil water conditions, higher δ 13 C values primarily result from the stronger decrease in stomatal conductance compared to the one in CO 2 assimilation (Farquhar et al. 1989). Where water is non-limiting, other factors than drought, including light availability, may have an influence on the wood isotopic signature. For instance, change in light competition arising from stand development could have an impact on the wood carbon isotopic composition over the 20-year study period, yet we expect this effect to be limited for such predominantly mature stands.
For each site, five cores from five different trees per species and per stand type were randomly selected from 10 to 20 dominant trees that had been cored per species in each stand. The cores were sampled to the pith at a height of 1.3 m. For each core and for each annual ring within the period 1993-2014, we used a scalpel and a stereomicroscope to sample the last third of the ring. We only took the last third of each ring in order to avoid carry-over effects on wood isotopic composition (Michelot et al. 2012). While this does not match specifically the latewood, it ensures that the carbon isotope signature we obtained for each year characterized the functioning of the selected trees during the second part of the growing season. The five samples were then pooled per species and stand resulting in four samples per triplet and year (total number of samples: 4 × 17 [number of triplets] × 22 [years] = 1.496).
At the INRA Silvatech platform (Nancy, France), the pooled samples were ground to a fine powder in a ball mill (MM400, Retsch). The 13 C/ 12 C ratio was measured with a mass spectrometer (Isoprime 100 (Isoprime Ltd., Cheadle Hulme, UK) coupled with an elemental analyser (Elementar vario, ISOTOPE cube, Elementar Analysen Systeme GmbH, Hanau, Germany)). The standard deviation for the analysis of standard saccharose was 0.12‰.
The isotopic composition (δ 13 C) relative to the standard Vienna Pee Dee Belemnite scale was calculated as follows (Eq. 1): with R standard being the isotopic ratio ( 13 C/ 12 C) of a belemnite fossil from the Pee Dee Formation, corresponding to the international standard (IAEA 1995).
Isotopic composition of pooled samples was corrected to take into account the change in the isotopic composition of atmospheric CO 2 due to industrialization, Eq. 2 (McCarroll and Loader 2004): with δ 13 C plant being the isotopic ratio of the plant, δ 13 C atm being the isotopic ratio of the atmosphere and − 6.4‰ corresponding to a preindustrial standard value.
(1) δ 13 C plant = (R sample ∕R standard ) − 1 × 1000 (2) δ 13 C = δ 13 C plant − δ 13 C atm + 6.4 Fig. 2 Components of the water balance (precipitation + maximum soil water available (SWA) − potential evapotranspiration (PET)) for each site. Climatic components of the water balance are calculated over the vegetation period (March-September) and averaged over the 1950-2014 period. Sites are ordered by increasing value of mean water resources ((P + SWA) − PET) Isotopic composition was also corrected to reflect rising atmospheric CO 2 concentrations since 1850 (reference period) following the method suggested by McCarrol et al. (2009). This nonlinear method aims to extract low frequency variations in δ 13 C series based on a theoretical plant's reaction to rising atmospheric CO 2 (Lévesque et al. 2013). Corrected δ 13 C series are hereafter referred to as δ 13 C cor .

Climate data
We used the 0.25°-gridded E-OBS dataset from EU-FP6 ENSEMBLES project (Haylock et al. 2008;van den Besselaar et al. 2011). From this dataset, we obtained series of daily minimum, maximum and mean temperatures, cumulative daily precipitation and daily average sea level pressure for the period 1950-2014. Monthly potential evapotranspiration (PET) was derived from these data following the modified Hargreaves equation (Choisnel et al. 1992, Droogers andAllen 2002), which provides estimations that are close to those obtained from the Penman-Monteith equation (Beguería et al. 2014).
We used both long-term (averaged over the 1950-2014 period) and short-term (inter-annual) water balance indices. As a long-term index of water availability, we used a simplified water balance calculated over the vegetation period (March-September), WB VP , defined as: total precipitation over the vegetation period (P) − total potential evapotranspiration over the vegetation period (PET) + maximum soil water available (SWA). Maximum soil water available (data from Forrester et al. 2017a, b) was calculated from soil depth, soil water holding capacity estimated from soil texture and the amount of stones in the soil. To take into account short-term variations in water availability during the period 1993-2014, we used the Standardized Precipitation Evapotranspiration Index (SPEI). SPEI (Vicente-Serrano et al. 2010) is a (monthly) multi-scalar index that can be calculated over different time scales, and which can be used to determine the onset, duration and magnitude of drought conditions with respect to normal conditions. The average SPEI value over 1993-2014 was zero for each site. Positive values indicate above-normal wet conditions, whereas negative values identify dry periods. SPEI values between − 0.67 and 0.67 are considered normal, values between − 0.67 and − 1.28 indicate moderate drought, and values < − 1.28 indicate severe drought (Isbell et al. 2015). SPEI was calculated over the second half of the vegetation period (June-September) which corresponds to the part of the tree rings that were sampled. Calculation was made using the SPEI package in R software (R Core Team 2014; Beguería et al. 2014).

Spatial and temporal patterns of δ 13 C in pure stands
Welch's t-tests were used on the annual δ 13 C cor values to assess the tree species effect (pine vs beech) within each site. Linear models were then used to test for the effects of site and stand variables (Table S4) on the average differences between beech and pine δ 13 C chronologies per triplet; an lntransformation of the response variable was used to reduce the heteroscedasticity of the residuals.
Next, we fitted a linear mixed effect model for each species separately on δ 13 C cor time series, considering site as a random factor in order to identify the major determinants of spatial and temporal variation in δ 13 C in pure plots for each species: where β is the vector of the fixed effects parameters, E is the matrix of the predictors of the fixed effects, S is the random factor characterized by the inter-site variance 2 site and ε is the error term. A series of climatic variables, site and stand attributes (see full list in Table S4) and their interactions were used as candidate variables for fixed effects. To avoid missing any key explanatory variable, we first used various selection procedures (Lasso, Elastic Net and stepwise forward selection with Akaike Information Criterion (AIC) and Bayesian Information Criterion (BIC)) and retained all variables selected by at least one method. Variance Inflation Factors (VIF) were calculated in order to measure the degree of multi-collinearity of the variables. The low VIF values (< 4) indicate that multi-collinearity was not a problem (O'Brien 2007). Starting from the model with the full set, the variables with the lower predictive power were then progressively removed based on the likelihood ratio test (González de Andrés et al. 2018;Zuur et al. 2009).
To investigate whether the two species had a similar temporal response to environmental fluctuations (synchronism), we calculated correlation coefficients between the beech and pine time series of non-corrected δ 13 C per site, following del . Values can range from − 1 (complete asynchrony of species response to environmental fluctuations) to + 1 (complete synchrony).

Temporal constancy
Two types of indices were used to analyse the temporal constancy of carbon isotope series: (1) the Temporal Stability Index (TS-Eq. 4) and (2) sensitivity (Eq. 5). TS is an indicator of the dispersion of corrected δ 13 C values with regard to the mean isotopic composition, while sensitivity is an indicator of the year-to-year variability of the time series. Increased TS indicates higher stability due to either a higher mean value or a lower standard deviation. High sensitivity values indicate important year-to-year variation in δ 13 C values due to high sensitivity to external (e.g. climatic) parameters.
where µ is the mean of the corrected δ 13 C series; σ is its standard deviation; n is the number of year; and S i+1 = (δ i+1 − δ i ) and δ i are the corrected isotope values (Saurer et al. 1997). TS and sensitivity are calculated per site, species and stand type (pure vs. mixed). We used mixed effects models with site as a random intercept to test the species mixture effect on these two indices, as follows: where TC j is the temporal constancy index (temporal stability or sensitivity) of one species at site j; a 0 and a 1 are the fixed parameters of the model; a oj is the random parameter associated with site; stand composition is a dummy variable with two levels (pure/mixed); and ε is the error term. Site-and stand-level variables (e.g. BA, altitude, WB VP ) that could potentially influence TC were included as additional predictors in the model (Table S4).

Spatial and temporal patterns of mixing
For each species, the mixing effect was quantified by Δ δ 13 C, the difference in δ 13 C between pure and mixed stands per year and site. For this specific analysis, non-corrected δ 13 C values were used because computing the difference in δ 13 C values from a same site cleaned the signal of long-term trends.
For each species, we first tested whether the Δ δ 13 C time series significantly deviated from the expected null value under a no mixing effect hypothesis, using one-sample t-tests for each site separately as well as across all sites.
The possible effects of average site conditions, stand structure and inter-annual climate variability on the spatial and temporal patterns of the mixing effect were then analysed by fitting linear mixed models on the Δ δ 13 C series of each species, following the same procedure as previously explained for the δ 13 C cor time series from the pure stands.
We also investigated the effect of species mixture on species asynchrony by analysing the relationship between correlation coefficients of beech and pine time series in pure and mixed stands.

Resilience and resistance to drought in pure and mixed stands
Resistance and resilience were calculated to analyse tree reaction to a drought event (Lloret et al. 2011). Resistance was calculated as the difference in δ 13 C values between a reference year, i.e. without water limitation preceding the drought, and a drought year with water shortage (δ 13 C (ref, before) − δ 13 C (drought) ); resilience was computed as the difference in δ 13 C values between two reference years preceding and following the drought year, respectively (δ 13 C (ref, before) − δ 13 C (ref, after) ). To select the drought years, we first identified the year with the lowest SPEI (June-September) values during the 1994-2013 period. In order to avoid carry-over effects, we then checked that the SPEI value for the previous year had been normal or moderately wet/dry (i.e. within the interval [− 1.28 to 1.28]). If this specification was not met, we shifted to the year with the second lowest SPEI value and started over. The reference years were selected as the wettest years preceding and following the corresponding drought year, and the associated SPEI values were further checked to ensure they were > − 0.67. Selected years and their associated SPEI and P − PET values are presented in Table S5.
The resistance index should be negative since δ 13 C cor values are expected to rise during drought events; the more negative the resistance index (low resistance), the higher the drought effect (stress). A resilience value not significantly different from zero indicates that trees have a high capacity to return to pre-drought δ 13 C cor levels after being subject to a drought event (high resilience); a negative value indicates low tree resilience.
For each index, we tested for species and stand composition effects separately according to the following mixed model: where R j is either resistance or resilience at site j; a 0 and a 1 are the fixed parameters of the model; a oj is the random parameter associated with site j; and ε is the error term. Effect is a dummy variable with two levels (pine vs beech or pure vs mixed for the species identity and stand composition effects, respectively). Site-and stand-level variables (e.g. BA, altitude, WB VP ) that could potentially influence R were also included as additional predictors in the model (Table S4).
All statistical analyses were conducted with the R software, version 3.4.1 (R Core Team 2014). Mixed models were fitted with the package "nlme" (Pinheiro et al. 2017).

Spatial and temporal patterns of δ 13 C in pure stands
In general, pine reached higher (less negative) mean isotope values than beech (Fig. 3). Average beech and pine δ 13 C cor values across all sites were − 25.3‰ and − 23.7‰, respectively. The range of site average values was [− 27.9 to − 22.2‰] for beech and [− 26.0 to − 21.1‰] for pine.
Average δ 13 C cor values were always lower for beech than for pine, for all sites. Differences were always statistically significant, with the exception of Pol1 (Table S6). The magnitude of the difference between species was site dependent with values ranging from 0.25‰ (Pol1) to 3.23‰ (Lit2). None of the variables tested in the linear model (Table S4) had a significant effect on this difference (data not shown).
There was a global coherence between the beech and pine time series ( Figure S1) that was confirmed by the high correlation coefficients values observed in most sites (Table S7). Sites with a lower level of correlation between the two species' time series were well distributed along the gradient of average water availability (sites with WB VP ranging from −249 to 632 mm, while water availability across the whole gradient ranged from −366 to 632 mm).
Linear mixed effect models adjusted on each species δ 13 C cor time series highlighted a significant, negative effect of SPEI June-September on δ 13 C cor for both beech and pine (Table 1). In pine stands, this SPEI effect depended on basal area only, while in beech stands, it depended on basal area, WB VP and slope (Table 1). In pine stands, higher basal area was associated with a more negative slope of the δ 13 C cor /SPEI relationship. In beech stands, higher basal area, lower WB VP and less steeper stand slopes were associated with a more negative slope of the δ 13 C cor /SPEI relationship.
Only beech showed a significant influence of WB VP on δ 13 C cor . This influence followed the expected pattern, i.e. higher values of average water resources were related to lower δ 13 C cor values.

Temporal constancy
Temporal stability and sensitivity values are shown in Table S3. There was no significant difference in temporal stability nor sensitivity between pure and mixed stands for any of the species (Table 2).
Climatic variables (average WB VP and average temperature over the period 1993-2014) did not explain the variability of temporal stability between sites; for pine, temporal stability increased as stands became older. No climatic or site characteristics were found to impact sensitivity (Table 2).

Spatial and temporal patterns of species mixture effect
For beech, seven sites showed no significant species mixture effect on δ 13 C, while two showed a negative effect and eight showed a significant positive effect (lower δ 13 C values in mixed stands compared to pure ones) (Fig. 4). The mean Δ δ 13 C across all sites was significantly positive (value = 0.3183, p value = < .001).
For pine, eight sites showed a significant negative effect of species mixture, five showed a significant positive effect, and four showed no significant effect (Fig. 4). The mean Δ δ 13 C across all sites was significantly different from zero (value = − 0.1080, p value = .01).
Four sites showed positive species mixture effects for both species and two showed negative effects for both species. When species mixture effect was positive for pine, it was either also positive for beech (four sites) or neutral (no significant species mixture effect; one site), but never negative. On the other hand, when species mixture effect was positive for beech, it was positive (four sites), negative (three sites) or neutral (one site) for pine. Three sites also displayed a negative species mixture effect for pine but no significant effect for beech, and three sites displayed no significant species mixture effect for either species ( Figure S2).
Looking at the drivers behind the temporal and spatial patterns of mixture effects (Table 3), we found that beech Δ δ 13 C was not influenced by any of the site or stand characterization variables, but that SPEI had a significant positive effect. Between-site variability of pine Δ δ 13 C was linked to WB VP and to mean age of the stand (the Δ δ 13 C/WB VP relationship becoming less negative as Table 1 Parameter estimates, p values and R-squared for the linear mixed models adjusted on the 1993-2014 δ 13 C cor series in the pure stands For both models, we used site as a random intercept. Marginal R-squared (R 2 m) represents the variance explained by fixed factors. Conditional R-squared (R 2 c) represents the variance explained by both fixed and random factors (whole model). WB VP is the average water balance over the vegetation period (precipitation + potential available soil water − potential evapotranspiration) calculated over the period 1950-2014.  (Table 3). There was no significant difference in correlation coefficients between the δ 13 C time series for beech and pine in pure and mixed stands ( Figure S3), thus indicating that species mixture did not change the synchrony of the two species' reactions to environmental fluctuations.

Resilience and resistance to a drought event in pure and mixed stands
In pure stands, pine displayed a negative index for resistance to drought (p value = .03). This effect was not significant for beech (p value = .07). The resilience index was not significantly different from 0 for pine but was significantly positive for beech (Table 4 and Fig. 5). No effects of site or of stand characterization variables were found except for the beech resilience index, which was significantly influenced by stand age (estimate = − 0.022; p value = .00) and site WB VP (estimate = 0.002; p value = .01) (data not shown).
There was no significant mixture effect on either resistance or resilience, for either of the species (Table 4, Fig. 5).

Spatial and temporal variation in δ 13 C in pure stands
Carbon isotope composition in tree rings was systematically higher in pine than in beech (Fig. 3), pointing to greater intrinsic water use efficiency for pine compared to beech; this is consistent with previous studies (Daux et al. 2018;Hemming et al. 1998;Szczepaneck et al. 2006). Several Fig. 4 Mean difference between δ 13 C in pure and mixed stands (Δ δ 13 C = δ 13 C pure − δ 13 C mixed ), for beech (black bullet) and pine (grey bullet). Vertical bars represent 95% confidence intervals around the mean. The dashed line indicates zero. The absence of intersection between this line and the confidence interval bars gives strong indication of a mean which is significantly different from zero. Sites are shown in increasing order of average WB VP calculated aver the 1950-2014 period Table 3 Parameter estimates, p values and R-squared for the linear mixed models adjusted on the 1993-2014 series of δ 13 C pure stands − δ 13 C mixed stands For both models, we used site as the random intercept. Marginal R-squared (R 2 m) represents the variance explained by fixed factors. Conditional R-squared (R 2 c) represents the variance explained by both fixed and random factors (whole model). WB VP is the average water balance over the vegetation period (precipitation + potential available soil water − potential evapotranspiration) calculated over the period 1950-2014  explanations for these differences are possible. Firstly, because of differences in physiological or morphological characteristics (such as higher light availability associated with lower light interception in pine stands due to a less dense canopy), the carbon uptake (A) is higher in pine resulting in lower leaf internal CO 2 concentrations, thus leading to increased δ 13 C all other things being equal. Daux et al. (2018) recently has discarded this explanation as a cause for the observed difference in isotopic composition between beech and pine because conifers usually have lower A than broadleaved trees. However, Medlyn et al. (1999) report a higher potential electron rate and maximum rate of Rubisco activity for pine compared to beech, suggesting that this general rule of lower A for conifers than for broadleaves might not hold true for pine and beech. The difference in δ 13 C levels between the two species could also be explained by lower stomatal conductance (g s ) in pine. Lower g s could originate either as a direct effect of morphological characteristics (e.g. lower stomatal density, smaller stomata), or as an indirect effect of ecological functioning. Indirect effects include (1) lower leaf area index in closed stands for pine compared to beech leading to higher evapotranspiration from the soil and the understory (Daux et al. 2018) and (2) lower access to belowground water reserves due to shallower rooting (Daux et al. 2018). These indirect effects tend to reduce water availability in pine stands. However, if indirect effects do indeed prevail, we would expect the difference in δ 13 C between species to be lower or, even, to disappear when water availability is higher. We did not observe any such pattern leading us to think that the difference in δ 13 C between species is probably due to an effect either of light interception or a difference in stomatal characteristics of beech and pine. As hypothesized, pine and beech δ 13 C cor values were significantly influenced by water availability during the last part of the vegetation period (Table 1). These results are consistent with previous studies and with physiological models of the response of isotopic discrimination and water use efficiency during carbon assimilation under soil drought conditions (Farquhar et al. 1989). Saurer et al. (2008) found that δ 13 C chronologies of pine and beech were negatively correlated with precipitation on non-water-limited sites, while González de Andrés et al. (2018) found a negative influence of water balance (P − PET over the summer) on δ 13 C cor for climatically contrasting sites. A few sites (Sp1, Sp2, Swe1, Lit1 and Lit2) did, however, not respond to SPEI, even when additional starting time and aggregation periods were considered ( Figure S4); in addition, none of the site, stand or climate characterization variables were able to explain the variability in the δ 13 C cor /SPEI relationship for this limited subset of sites.
Regarding the spatial variations in δ 13 C, we expected to find a negative relationship between carbon isotopic Table 4 Parameter estimate, standard error and p values for the models testing for species identity and species mixture effects on resilience component indices For species identity effect, beech is the value estimated for pure beech, pine − beech is the difference between pure pine and pure beech stands, and pine is the value estimated for pure pine stands. For species mixture effect, pure is the value estimated for pure stands and mixed for the species in mixed stand, and mixed − pure is the difference between mixed and pure stands. Significant effects are indicated in bold composition in tree rings and average water availability for both species (Saurer et al. 1995). While δ 13 C cor was negatively related to WP VP for beech, confirming the hypothesis of higher δ 13 C cor levels in drier sites for that species, this was not the case for pine (Table 1). We hypothesize that this is due to the nonlinear pattern of variation in average δ 13 C cor along the water availability gradient (Fig. 3), where dry sites clearly displayed higher δ 13 C cor values. It is important to note that spatial and temporal variations in carbon isotope composition are not independent of each other. Indeed, we found that for beech, due to the SPEI/ WB VP interaction, the WB VP effect on δ 13 C cor disappeared during extremely wet years (SPEI > 3), but held in other situations. Saurer et al. (1995) found a similar increase in δ 13 C cor values for pine and beech on drier sites in Switzerland. The absence of any significant influence of the interaction term "SPEI June-September × WB VP " on pine δ 13 C cor suggests the existence of local adaptation mechanisms as well as long-term genetic divergence within species; this means ecotypes vary in functional traits, as previously proposed by Weigt et al. (2015) and Härdtle et al. (2013). We also found that the inter-annual variation in δ 13 C cor was influenced by stand variables (basal area for pine and basal area and slope for beech). Basal area in both pine and beech stands (and lower WB VP in beech stands) could have been "aggravating factors" as they induce a higher sensitivity of δ 13 C cor to annual water balance (more negative δ 13 C cor /SPEI slope). The aggravating effect of basal area can probably be linked to increased competition among trees for soil water, and low WB VP is likely to be correlated with higher sensitivity to annual variations in the water balance. Surprisingly, stand slope did not have such an aggravating effect on the δ 13 C cor / SPEI relationship in the pure beech stands, probably because this effect was confounded with that of slope aspect.
The δ 13 C cor time series for both species were quite coherent within a site ( Figure S1, Table S7), which is a further indication that beech and pine trees tend to respond in a comparable way when facing similar temporal fluctuations of environmental conditions (Table 1). There were, however, a few exceptions. In most cases, those exceptions were linked to the long-term trend in the δ 13 C cor time series of one of the two species, thus decoupling pine and beech  Table 4 δ 13 C cor values. For instance, such an effect can be seen in the decreasing trend in beech δ 13 C cor series at Bel2 or in the drop in beech δ 13 C cor values around 1995-2003 at Bel1 ( Figure S1). The considerable length of those trends suggests that they are not of climatic origin but are probably rather due to changes in stand characteristics (e.g. changes in access to light, management effect), long-term changes in global circulation patterns (González de Andrés et al. 2019;Rozas et al. 2015;Sardans et al. 2017) or to tree weakening (dieback), which could have influenced the physiological functioning of the trees.
Looking now at the tree reaction to a drought event ( Table 4, Fig. 5), a significant increase in δ 13 C cor was observed for pine; surprisingly, this effect was non-significant for beech, but by a very slight margin. Inspection of Fig. 5 suggests that beech response was related to site water availability, the drought effect being generally dampened on wetter sites. Such a difference in drought reaction between dry and wet sites was not observed for pine. We attribute this to the fact that pine is a "drought-avoiding" species. Such species close their stomata quickly during water shortages to avoid damage to the conductive system (Cochard 1992;Martínez-Vilalta et al. 2004). This drought avoidance strategy is common in conifers (especially Pinus species), which tend to have lower embolism resistance than angiosperms (Martínez-Vilalta et al. 2004;Choat et al. 2012). In addition, LAI reduction through leaf shedding under drought stress has been shown to be a strategy which buffers water loss in pine trees, resulting in the formation of tree rings with increased intrinsic Water Use Efficiency (iWUE) (González de Andrés et al. 2019). Beech, on the other hand, is more anisohydric (Pflug et al. 2018;Schäfer et al. 2017). However, there was no overall significant difference between pine and beech resistances. This could indicate that, while pine reacts quicker than beech to drought, both species end up being affected in a similar way during extreme events. Betweensite variability in resistance to drought event could also be partially explained by other variables not considered in this study, such as nutrient availability. Potassium and phosphorus are known to be particularly important, as they are involved in various water-related mechanisms (Newton et al. 1986). For instance, potassium influences stomatal function, osmosis and hydraulic conductance (González de Andrés et al. 2019, Khosravifar et al. 2008, Wang et al. 2013). We did not find any lasting effect of drought on tree functioning as trees of both species were able to return to pre-drought levels of carbon isotope composition (pine) or even show lower δ 13 C cor values (beech) after the extreme event.

Temporal constancy
A species mixture effect on temporal stability or sensitivity could be caused either by a differential response (sensitivity to different parameters or differential temporality of the response) of each species to environmental changes (Loreau and de Mazancourt 2008;Hector et al. 2010), or by reduced competition in mixed stands compared to pure ones (Tilman 1999). The lack of any mixing effect in our case (Table 2) is therefore in agreement with the tight correlation we observed in pure stands between the δ 13 C time series of beech and pine trees growing in a same site (Table S7); in addition, as there was no significant difference in the correlation between the two species in mixed and pure stands, species mixture did not induce a decoupling of the species reaction to environmental fluctuations ( Figure S3).
Contrasting results have been reported with respect to diversity-stability relationship, from the total absence of a stabilizing effect in single-trophic communities (Jiang and Pu 2009) to higher stability in mixed forests . It is currently becoming more and more accepted that diversity improves stability at the community level but decreases stability, or does not affect it, at the species level. del  recently have highlighted the stabilizing/destabilizing pattern for productivity in mixed pine/ beech stands across Europe. Although water is often considered a main factor of resource-driven effects, we found no clear stabilizing/destabilizing species mixture effect on water-related processes at the species level.

Species mixture effects on δ 13 C
While the mixing effect (Δ δ 13 C) was on average positive for beech and non-significant for pine across the whole gradient, the major result or our study is that this effect strongly differed between sites, depending on the average water availability (Fig. 4).
On wet sites (i.e. sites with permanent available belowground water resources), the species mixture effect on δ 13 C should be close to zero on average, that is, if the species mixture effect is indeed mainly due to waterrelated mechanisms. If the species mixture effect differs from zero, then other mechanisms should be considered.
A key mechanism influencing δ 13 C values, and therefore Δ δ 13 C, is the access to light (Ehleringer et al. 1986;Farquhar et al. 1989). Our results showed a high variability of Δ δ 13 C in wet sites (Bel2, Swe2, Ger6 and Lit2), which could be due to a species mixture effect on light availability. Indeed,  and Barbeito et al. (2017) have shown that higher stand diversity is associated with increased crown size, which could lead to increased light absorption. In addition, Barbeito et al. (2017) showed that this effect is dependent of site conditions, the difference between pure and mixed stands being higher in more productive sites. We used delta height (the difference between the height of the cored target-species trees in mixed stands and the mean height of the mixed stand) to investigate the potential effect of light interception on Δ δ 13 C. However, light interception did not fully explain Δ δ 13 C deviation from zero on very humid sites. Indeed, some sites for which access to light did not differ between pure and mixed stands also displayed significant Δ δ 13 C deviation from zero and vice versa (we found differences in access to light but no significant Δ δ 13 C deviation from zero). On moderately wet sites (average WB VP close to zero), the same consideration should hold (species mixture effect close to zero), on the condition that the addition of a second species does not influence water availability (increased belowground competition in mixed stands compared to pure ones). If this is the case, species mixture effect could be negative for one or both species, depending on their ecophysiological characteristics. On this type of site, we found a high variability in beech Δ δ 13 C, suggesting the influence of non-water-related mechanisms as stated earlier; for pine, Δ δ 13 C was consistently negative. We could thus conclude that in moderately wet sites (theoretically non-stress sites), adding beech would induce a stress on pine, which is consistent with previous findings (González de Andrés et al. 2018). This would also be consistent with the difference in the compromise strategies of the two species between carbon uptake and water loss highlighted above for pure stands. On dry sites, the species mixture effect is expected to be positive (see Forrester and Bauhus 2016) (1) if species mixture has a positive influence on water availability, (2) if this influence is large enough to affect carbon isotope composition and (3) if potential negative species mixture effects (competition) are lower than the positive effects. We found that species mixture effect tends to be positive for both species on dry sites (SP2, CZE1, FRA1) but, as we move towards extremely dry sites (SP1, BUL1), this positive effect seems to disappear. This is probably indicative of the fact that the positive species mixture effect on water availability is not strong enough to compensate for the increasing environmental constraint. While this result questions the linear relation between complementarity and resource availability along a gradient of environmental conditions hypothesized by Forrester and Bauhus (2016), a possible decrease in the positive effects of species mixture in extremely harsh situations is not a new idea (Maestre and Cortina 2004;Tielbörger and Kadmon 2000). The species mixture effect presented in this study may therefore be indicative of a more complex structure with threshold points (Fig. 4). These threshold points could correspond to the level of average site water availability where beech starts to regulate its water consumption (Fig. 3), thus reducing competition and inducing a switch from negative to positive species mixture effect in pine. Additional sites would, however, be needed to fully test this hypothesis. The complex relationship between species mixing effects and average site water availability could also arise from interactions with other drivers of water availability, in particular those related to nutrient availability (Sardans and Peñuelas 2007).
Our models for pine confirmed the role average site water balance plays on the species mixture effect, in contrast to the lack of any effect related to inter-annual variations in SPEI (Table 3). The average site water balance effect further depended on the mean age of pine stands since, in older stands, the slope of the Δ δ 13 C/WB VP relationship was less negative. We should eliminate two unlikely causes of the age effect in mature stands such as the ones used in this study: (1) vegetation growing close to the forest floor using air with increased 12 C/ 13 C ratio due to respiration (McCarroll and Loader 2004) and (2) variation in bark refixation of respired CO 2 as bark is usually too thick in mature stands for bark refixation to play a major role (McCarroll and Loader 2004). Ontogenic changes in rooting patterns (Claus andGeorge 2005, Finér et al. 2007, Yuan and Chen 2010) could possibly be involved, but testing this hypothesis would need thorough field root investigation we did not performed. It is possible, however, that this age effect was confounded with the effect of height, considering that using mean stand height instead of age as a variable only slightly decreases the performance of the model. In contrast to pine, difference in WB VP was unable to explain the variability in mixing effects for beech. The beech model rather highlighted an annual water balance (SPEI) effect on beech Δ δ 13 C, with increased complementarity effects during wet years. Altogether, those findings are in agreement with the beech's more intense competitive strategy (González de Andrés et al. 2018).
While mixing affected the δ 13 C signal for both species in a complex way related to site water availability (Fig. 4, Table 3), it did not influence the tree response to a drought event (Table 4, Fig. 5). Considering the whole gradient, pure plots displayed higher δ 13 C levels under more intense waterlimiting conditions, and mixed stands had a similar behaviour. Such a result is consistent with the growing body of literature on the subject, which reports that species mixture does not always improve reaction to drought (Bonal et al. 2017;Grossiord et al. 2014), although this may indeed be the case in certain situations (Grossiord et al. 2015;Lebourgeois et al. 2013;Pretzsch et al. 2013).

Conclusion
We conclude from the present study that pine and beech present different levels of average δ 13 C values indicative of the compromise between CO 2 assimilation and H 2 O loss, but that the spatial and temporal variations in their δ 13 C values are similar.
We did not find any species mixture effect on tree's reaction to drought events, in accordance with the growing body of literature on this topic. However, analyses of spatial and long-term temporal variations in species mixture effects showed that mixing species leads to contrasting effects on beech and pine carbon isotope composition (a slightly positive effect for beech and no significant effect for pine) when the whole gradient of water availability is taken into account. The global pattern of species mixture effect along this gradient is consistent with some theories (Forrester and Bauhus 2016): an increasingly positive species mixture effect on drier sites until the drought constraint becomes too strong for the species mixture effect to compensate. Our study demonstrates the importance of considering the nonlinear relationships of the species mixture effect on wood isotopic composition and that the species mixture effect appears at certain threshold points. Intrinsic species characteristics concerning water-related processes play a critical role in species mixture effect, especially at moderately wet sites. A combination of the difference in the two species' CO 2 / H 2 O compromise and average environmental conditions in terms of water availability therefore determines the balance between competition and complementarity in mixed stands.