Both Mature Patches and Expanding Areas of Juniperus thurifera Forests Are Vulnerable to Climate Change But for Different Reasons

Abstract

Research Highlights: Water use efficiency (WUE) varied along a gradient of Juniperus thurifera (L.) forest expansion, being higher in recently colonised areas. Background and Objectives: WUE is a classic physiological process of plants that reflects the compromise between carbon assimilation and water loss and has a profound influence on their performance in water-limited environments. Forest expansion in Mediterranean regions associated with land abandonment can influence the WUE of plants due to the existence of two opposing gradients: one of favourable–unfavourable environmental conditions and another one of increased–decreased intraspecific competition, the former increasing and the latter decreasing towards the expanding front. The main objective of this study was to elucidate how the WUE of Juniperus thurifera varied along the stages of forest expansion and to provide insight on how this variation is influenced by intraspecific competition and abiotic factors. Materials and Methods: Seventeen plots at different distances from the mature forest core were selected at three sites located in the centre of the Iberian Peninsula. For 30 individuals within each plot, we measured biometric characteristics, age, tree vigour, and C/N ratio in leaves, and the leaf carbon isotope signature (δ¹³C (‰)) as a proxy for WUE. Around each individual, we scored the percentage cover of bare soil, stoniness, conspecifics, and other woody species. Results: WUE of J. thurifera individuals varied along the forest expansion gradient, being greater for the individuals at the expanding front than for those at the mature forest. WUE was influenced by the cover of conspecifics, tree age, and C/N ratio in leaves. This pattern reveals that less favourable environmental conditions (i.e., rocky soils and higher radiation due to lower vegetation cover) and younger trees at the expanding front are associated with increased WUE. The increased cover of conspecifics decreases irradiance at the mature forest, involving milder stress conditions than at the expanding front. Conclusions: Lower WUE in mature forests due to more favourable conditions and higher WUE due to abiotic stress at expanding fronts revealed high constraints on water economy of this tree species in these two contrasting situations. Climate change scenarios bringing increased aridity are a serious threat to Juniperus thurifera forests, affecting both mature and juvenile populations although in different ways, which deserve further research to fully unveil.

1. Introduction

Plant survival in environments where water is the main limiting factor is significantly influenced by the water strategy of the species and their overall drought tolerance [1]. Plants respond to a decrease in water availability through various mechanisms such as stomatal closure to minimise water loss and reduce the risk of hydraulic dysfunction [2]. Water use efficiency (WUE) of plants partly incorporates the response to a decrease in water availability, quantifies the trade-off between carbon gain by photosynthesis, and water loss by transpiration [3], and it is positively associated with tolerance to water stress [4]. Thus, the study of this variable is key to unveiling the current strategy and the future vulnerability of plants in increasingly drier environments.

Typically, WUE has been studied considering the variation of a single abiotic factor [5], mainly soil water availability [4]. However, there are many environmental factors that affect WUE variation such as atmospheric CO₂ [6], soil nutrient availability [7], irradiance [8], number of rocks and stones in the soil [9] as well as the interactions among them [8,10]. In general, individuals subjected to more extreme conditions such as higher irradiance, lower water availability, and/or poorer soils present higher WUE [11], but the complex set of interacting factors along environmental gradients might render several alternative patterns.

Although WUE has been widely studied considering abiotic factors as the main source of variation [12], it can also be influenced by a range of biotic factors that include the ubiquitous intra- and interspecific competition [13]. Competition for the available water resources can promote an increase in WUE, especially in the case of woody species. This has been observed for Fagus sylvatica L. and Pinus sylvestris L. in mixed forests of northern Spain (interspecific competition) [14] and for monospecific stands of Quercus suber L. in the centre of Portugal (intraspecific competition) [15]. Furthermore, the interactions between both biotic and abiotic factors have also been shown to affect the WUE of plants [16].

The current increase in environmental adversity due to climate change, mainly involving an increase in mean temperatures and in the frequency and intensity of drought periods, is co-occurring with land-use changes in many biomes [17]. Rural abandonment, due to the cessation of agricultural and livestock activities, is nowadays one of the most common land-use changes in temperate zones and particularly in the Mediterranean Basin [18]. Rural abandonment is consequently promoting colonisation by nearby forests, which become denser and expand towards adjacent abandoned farmland and formerly grazed areas [19], phenomenon known as ‘forest transition’ [20]. Despite the widespread nature of this phenomenon within the Northern Hemisphere and its potential impact on many ecosystem processes and functions [21], the ecological implications of spontaneous forest expansion are still poorly understood. In particular, how forest expansion affects the water balance of a region and the hydrological processes occurring at different scales ranging from the individual to the ecosystem is both a timing and important scientific question under the global change scenarios. There is evidence suggesting that warmer and drier conditions in Mediterranean ecosystems due to climate change, coupled with the increase in plant cover associated with forest transition, could reduce water availability [22], forcing in turn an increase in the WUE of trees.

In water-limited Mediterranean environments, forest dynamics and potential expansion into abandoned land is constrained by water availability [23]. The water use efficiency of trees can thus play a key role for the colonisation of land formerly used for cattle or agriculture. Juniperus thurifera is a good example of an endemic Mediterranean forest undergoing expansion towards adjacent abandoned lands [24]. These juniper forests are low-density and open, and grow under exceptionally harsh continental conditions in terms of climate and soil conditions [25]. They constitute a unique system to study WUE variation associated with forest transition and serve to model forest dynamics under climate change scenarios.

Along the expansion of a juniper forest, intraspecific competition decreases from the mature forest to the expanding front, while unfavourable environmental conditions increase from the forest to the expanding front [26,27]. Thus, in mature forests, which are denser and act as a key source of propagules for the expansion process, intraspecific competition is greater, but abiotic conditions are more favourable (i.e., soils with more nutrients [28] and lower water stress due to reduced evaporation and heat [29]). In contrast, at the expanding front, with scattered trees and plenty of open spaces and more rocky soils [27], intraspecific competition is low, but the environmental conditions are less favourable (i.e., lower retention capacity due to rocky soils [9] and greater radiation due to lower vegetation cover [29]). It remains unknown to which extent and which factors, either biotic (i.e., intraspecific competition) or abiotic (i.e., microhabitat and microclimate), influence WUE strategies of plants along forest expansion gradients.

The main objective of this study was to elucidate how the WUE of Juniperus thurifera varied along these gradients of forest expansion by assessing the influence of intraspecific competition and abiotic factors. There are three alternative hypotheses for the pattern of WUE along the expansion gradient of these Spanish juniper forests: (i) the WUE of Spanish juniper is greater in the mature forest due to an increase in intraspecific competition as discussed in previous studies [15]; (ii) WUE is greater at the expanding front due to more adverse conditions, which have been shown to increase towards the recently colonised areas [11]; and (iii) WUE exhibits no clear pattern due either to a compensation of these two opposing gradients, to the influence of other environmental factors, or to a negligible effect of these conditions on the realised WUE.

2. Materials and Methods

2.1. Study Species and Area

Juniperus thurifera L. (Cupressaceae) is a dioecious tree species with a restricted distribution in the Western Mediterranean Basin; in particular, there are populations in Spain, Morocco, France, Algeria, and the Italian Alps. Approximately 90% of this species distribution area is located in Spain [24], where it occupies 600,000 ha of the country, from which approximately 117,000 ha are monospecific forests [30]. Generally, J. thurifera is a dominant species forming low-density forests on poor, shallow, rocky soils that can establish at high altitudes and cover broad climate ranges [31]. This species tolerates broad temperature ranges (high temperatures in summer, and low ones and frost in winter), thereby, it is adapted to water stress typical of Mediterranean-climate regions, characterised by a summer drought period and moderate or low precipitation concentrated in spring [32].

The current study was carried out in central Spain, specifically in the Guadalajara Province, at the Alto Tajo Natural Park and surrounding areas. The climate in this area is continental Mediterranean, characterised by hot and dry summers and cold and snowy winters (mean annual total rainfall ± Standard Error (SE): 477.1 ± 15.6 mm, mean annual temperature: 10.4 ± 0.2 °C, Molina de Aragón 40°50′40′′ N, 1°53′07′′ W, 1063 m a.s.l., 1951–2017 period; Agencia Estatal de Meteorología (AEMET)). Through historical and current satellite images and aerial photographs, we located areas where Juniperus thurifera forests colonised abandoned farmlands and formerly grazed lands areas (for further details see [27]). We selected three sites that presented well preserved J. thurifera forests including Maranchón (MA), Huertahernando (HU), and Ribarredonda (RI). We visually observed a gradient of forest expansion by satellite images at each site (Figure 1) and located areas of well conserved cores of mature forests and areas recently established on adjacent abandoned lands (for further details see [33]). At these three sites, we established seventeen plots at different distances form the core, mature forest (seven plots in Maranchón, five plots in Huertahernando, and five plots in Ribarredonda; Table S1). We assigned each plot as either mature or expanding front, with a transition zone between these two stages (Figure 1). Mean plot size varied among stages due to differences in the density of J. thurifera adults (i.e., mean ± SE of plot area: mature forest (0.51 ha ± 0.14), transition zone (0.81 ha ± 0.46), and expanding from (1.35 ha ± 0.47)). Mean distance among plots varied among sites due to differences in site accessibility being 1.5 km in Maranchón, 0.4 km in Huertahernando, and 0.3 km in Ribarredonda.

2.2. Sample Collection

We conducted sample collection between September and October 2017 to ensure that all trees, and therefore their leaves, had gone through a stressful water deficit during the summer drought period. Within each plot, we georeferenced all juniper trees, obtaining a minimum sample size of 50 adult trees (i.e., those whose quadratic diameter (QD) was greater or equal to 3 cm, because it was the minimum size for which we found reproductive individuals) per plot. QD (also known as equivalent diameter) was calculated as the square root of the sum of square diameter at breast height of each stem of a tree. QD is proportional to the total cross-sectional area and is the methodology used to estimate the trunk diameter for multi-stem trees [34]. Within each plot, we randomly selected 30 adult trees. For each of these selected adult trees, we measured different biometric characteristics including QD (measured at 130 cm), average crown diameter (calculated as the mean of the diameter of two perpendicular axes that passed through the axis of the trunk; for that we used a DME-distance measurer Haglöf, Långsele, Sweden) and height (using a hypsometer-Haglof Vertex IV).

Moreover, we scored the gender of all the selected individuals and visually estimated their vigour using a semiquantitative scale (0–4): (0) tree with low-density foliage and three quarters of the tree was withered; (1) half of the tree presented dry o semi-dry branches; (2) tree with several dry branches; (3) the individual presented a dry branch; and (4) tree that had no apparent damage and presented a high density foliage. To estimate tree age, we bored each individual at 50 cm from the ground using a Pressler increment borer (Haglöf, Långsele, Sweden). Then, for each tree, we obtained two wood cores from the bark to the pith, each with a length equal to the trunk radio. Additionally, for each individual, we collected 10 foliar units formed during the previous year (i.e., 2016). To set the size of the foliar units for each tree, we considered the fractal architectural pattern that was repeated along each branch. In this sense, the architectural pattern of J. thurifera leaves varied among individuals; therefore, the morphology and size of the foliar units considered could also be affected depending on individual intrinsic characteristics such as its gender and/or its age. To minimise the possible microenvironmental variability, all the selected foliar units came from branches orientated south and located in the upper part of the canopy, that is, in the areas of greatest solar exposure.

We delimitated the area of influence for each individual as a circle around it where radius was equal to the tree average crown diameter. Within the area of influence of each selected tree, we estimated the percentage of cover of bare soil, stoniness, conspecifics, and other woody species. We considered the cover of conspecific individuals as a proxy of intraspecific competition while the cover of other woody species (mainly including Genista sp., Thymus sp., and Lavandula sp.) reflected a measure of interspecific competition.

2.3. Sample Processing

To estimate the age of each individual, the extracted wood cores were air-dried, glued onto wooden mounts, and polished using sandpaper of progressively finer grain until tree rings were visible. The cores were dated with a stereomicroscope; they were scanned at 1600 d.p.i. and ring widths measured to an accuracy of 0.001 mm using the software CooRecorder v9.3 (Cybis Elektronik, Saltsjöbaden, Sweden, 2018). Cross-dating of individual series was checked using CooRecorder and COFECHA programs (Holmes, Tucson, Arizona, USA, 1983).

For each individual, we pooled the ten foliar units collected, oven-dried them at 60 °C for 72 h and finely ground them (MM300, Retsch, Haan, Germany). We weighed a subsample between three and four mg from each sample of the powdered pooled leaves (using Mettler Toledo XP6U Micro-Balance, Barcelona, Spain) and placed it into tin capsules (Sn 98 capsules, Lüdiswiss, Flawil, Switzerland) in order to analyse carbon isotope signature (δ¹³C, ‰) and C and N content. Analyses were carried out at the UC Davis Stable Isotope Facility, using a PDZ Europa ANCA-GSL elemental analyser interfaced to a PDZ Europa 20-20 isotope ratio mass spectrometer (Sercon Ltd., Cheshire, UK). Carbon isotope signature (δ¹³C, ‰) values reported were defined as:

δ¹³C = (R sample/R standard) × 1000

where R sample/R standard is the ¹³C/¹²C ratio expression for leaf samples and international standard (Vienna Pee Dee belemnite; VPDB), respectively. Subsequently, the C/N ratio was calculated for each tree by dividing total carbon content (µg mg⁻¹ sample) by total nitrogen content (µg mg⁻¹ sample).

2.4. Statistical Analyses

To study the relationship between QD and tree age among expansion stages, we fitted a full linear mixed model effect with Gaussian error distribution including age (second degree polynomial) and stage of the forest expansion gradient as fixed effects. Plot was included as random intercept. To test the effect of all the fixed components, we compared the full model (the one including all the fixed effects) with all the possible models varying in their fixed effects by including subsets of the predictors (MuMIn, package, R; [35]). We selected the model that provided the best fit to the data, the one with the lowest value of the Akaike Information Criterion (AIC; [36]). We conducted the same full linear mixed model effect with QD in logarithmic scale. Finally, we calculated marginal and conditional r² for the best model (r.squaredGLMM function, MuMIn package [35]).

We conducted Kruskal–Wallis tests to check for differences between expansion stages and sites for all the variables (QD, height, average crown diameter, tree vigour, age, C/N ratio, cover of (i) bare soil, (ii) stoniness, (iii) J. thurifera, and (iv) other woody species). When we detected significant differences, we used post-hoc Mann–Whitney U tests to determine differences between pairs of stages or sites. To select the variables to be included in further models, we first explored the correlations existing among all them through Spearman’s correlation tests.

To study the variation in δ¹³C (‰) as a proxy of WUE, we fitted a full linear mixed effect model with Gaussian error distribution including stage of the forest expansion gradient, site, J. thurifera cover, age, C/N ratio, stoniness, and gender as fixed effects. Plot was included as a random intercept.

To test the effect of all the fixed components, we compared the full model (the one including all the fixed effects) with all the possible models varying in their fixed effects by including subsets of the predictors. We selected the model that provided the best fit to the data, the one with the lowest value of the AIC (Table S2). Finally, we calculated marginal and conditional r² for the best model (r.squaredGLMM function, MuMIn package [35]).

The analyses were carried out using statistical software R version 3.6.2 [37], specifically lme4 [38], stats [37], lmerTest [39], and corrplot [40] packages.

3. Results

We did not find any differences in the relationship between QD and tree age among expansion stages (Figure 2), showing a saturated dynamic in which QD increased with age to a saturation point (between 60–70 years). At that point, QD was constant (Figure 2; logarithmic scale). Biometric characteristics (QD, height, and average crown diameter) varied along the gradient of forest expansion, being significantly greater at the mature forest and lower at the expanding front. Regarding QDD, we did not find differences among sites; however, tree height was significantly greater in Maranchón and Ribarredonda than in Huertahernando, while we found the opposite trend regarding average crown diameter (

Table 1. Mean ± SE of each variable studied for each stage of the expansion gradient and for each study site (values for each site are the mean of each variable studied for all plots located in each site). QD (quadratic diameter); different letters indicate significant differences (p < 0.05).

	Gradient Stage			Site
	Mature Forest	Transition Zone	Expanding Front	Maranchón	Huertahernando	Ribarredonda
QD (cm)	22.80 ± 1.04 a	15.80 ± 0.69 b	13.15 ± 0.56 c	16.84 ± 0.82 a	17.86 ± 0.76 a	16.35 ± 0.85 a
Height (m)	5.76 ± 0.13 a	4.42 ± 0.10 b	4.09 ± 0.09 c	4.50 ± 0.12 b	4.71 ± 0.12 a	5.04 ± 0.13 a
Average crown diameter (m)	5.35 ± 0.16 a	4.49 ± 0.12 b	4.02 ± 0.11 c	4.39 ± 0.12 b	5.14 ± 0.14 a	4.29 ± 0.15 b
Tree vigour (0–4)	2.43 ± 0.06 c	2.73 ± 0.06 b	2.93 ± 0.07 a	2.75 ± 0.06 a	2.75 ± 0.07 a	2.60 ± 0.07 a
Age (years)	44.98 ± 1.40 a	30.77 ± 1.01 b	28.81 ± 1.13 b	38.30 ± 1.26 a	28.00 ± 0.97 b	35.27 ± 1.40 a
C/N ratio	38.86 ± 0.46 a	36.85 ± 0.34 b	38.91 ± 0.32 a	37.75 ± 0.31 b	37.04 ± 0.37 b	40.02 ± 0.47 a
Bare soil (%)	10.00 ± 0.84 a	8.23 ± 0.52 a	9.82 ± 0.73 a	7.65 ± 0.61 c	9.44 ± 0.69 b	11.79 ± 0.76 a
Stoniness (%)	32.05 ± 1.03 b	42.09 ± 1.29 a	43.51 ± 1.24 a	40.10 ± 1.21 b	47.93 ± 1.07 a	29.30 ± 0.91 c
Juniperus thurifera L. (%)	18.39 ± 1.17 a	9.32 ± 0.84 b	6.37 ± 0.88 c	10.81 ± 0.88 a	10.56 ± 1.09 a	12.01 ± 1.21 a
Woody species (%)	36.55 ± 1.44 a	40.36 ± 1.28 a	40.30 ± 1.25 a	41.44 ± 1.18 b	32.08 ± 0.98 c	46.91 ± 1.45 a

). Regarding tree age, individuals at the mature forest were on average around 20 years older than those located at the expanding front and around 10 years older than trees at the transition zone. Tree age was around 40 years old on average for the three sites but there were significant differences among sites, with the individuals from Huertahernando being the youngest (ca. 10 year younger than trees from the other two sites, Table 1). Vigour of J. thurifera individuals was significantly greater at the expanding front and lower at the mature forest, while non-significant differences were found among sites (Table 1).

The cover of stones was significantly lower in the mature forest than in the transition zone and at the expanding front; among sites, stoniness was significantly greater in Huertahernando and lower in Ribarredonda (Table 1). Density of J. thurifera was greater at the mature forest and lower at the expanding front and did not vary among sites. In contrast, the cover of other woody species did not vary among stages but differed among sites, being significantly greater in Ribarredonda and lower in Huertahernando (Table 1). Finally, the cover of bare soil did not vary among stages but differed among sites, being significantly greater in Ribarredonda and lower in Maranchón (Table 1).

C/N ratio also varied among gradient stages, being greater in leaves coming from individuals located at the mature forest and at the expanding front than for those individuals located in the transition zone; the C/N ratio also varied among sites, being greater in Ribarredonda than in Maranchón and Huertahernando (Table 1).

Among all the correlations studied (Figure S1), age was positively correlated with C/N ratio (Rho = 0.26, p < 0.001) and C/N ratio was negatively correlated with total nitrogen content (Rho = −0.97, p < 0.001), but there was no correlation between C/N ratio and total carbon content (Rho = −0.01, p > 0.05). We also found that tree vigour was negatively correlated with age (Rho = −0.14, p < 0.01) and positively correlated with δ¹³C (Rho = 0.14, p < 0.01).

The best linear mixed effects model (according to the AIC) that explained the variation in δ¹³C was the one that included the stage of the forest expansion gradient, site, cover of J. thurifera, tree age, and C/N ratio as fixed effects (marginal r² = 0.42; conditional r² = 0.51). We observed that δ¹³C values varied along the forest expansion gradient, being greater at the expanding front than at the mature forest and transition zone (Figure 3a). We found that an increase in cover of J. thurifera was negatively correlated with δ¹³C values (Figure 3b). Likewise, tree age and C/N ratio were negatively and significantly correlated with δ¹³C (Figure 3c,d, respectively). Among sites, individuals from Ribarredonda showed the lowest values of δ¹³C (Table S3).

4. Discussion

Water use efficiency (WUE) varied significantly along a gradient of forest expansion being greater at the expanding front than at the mature forests and transition zone. An increase in unfavourable environmental conditions (rocky soils and lower vegetation cover with the corresponding increased irradiance) together with younger trees in recently colonised areas was associated with an increased WUE.

Differences found in WUE among different stages of the forest expansion gradient were related to differences in environmental conditions. Juniperus thurifera occupies, in general, sites with unfavourable soils, high stoniness, and extreme basicity, avoiding competition with other fast-growing tree species, which tend to occupy better sites [35]. Additionally, stones covered a larger fraction of the ground at the expanding front, which has been related to a decrease in soil water retention capacity [25] and to a lower survival of juvenile trees of J. thurifera [27]. A decrease in soil water retention capacity in stony soils could translate into less water at the tree level and affect the WUE of trees [9,41]. Furthermore, these stony soils at the expanding front suffer from greater evapotranspiration due to a lower vegetation cover, which further promotes an increase in WUE [29]. The comparison among sites revealed that WUE was highest at those sites with the greatest stoniness (Table 1 and Table S3).

An increase in the density of J. thurifera trees was related to a decreased WUE. It is known that intraspecific competition decreases the WUE of individuals for other species (i.e., Abies pinsapo; [42]). Previous studies have reported that juveniles of J. thurifera that established under the canopy of adult individuals had lower WUE and greater survival than isolated individuals likely due to better environmental conditions under canopy cover [27,43]. At the mature forest, the decreased WUE of the individuals could be related to the observed improvement of the physical soil conditions in the vicinity of trees [44]. The increased shade associated with the increased cover of J. thurifera at the mature stages could also influence the observed pattern in WUE through milder conditions of irradiation. It has been shown that shade decreased the WUE of trees, even when trees were water stressed in comparison with trees exposed to full sunlight [8]. The decreased WUE in mature forests of J. thurifera could be induced by lower light availability and therefore less irradiance exposure, in addition to other environmental factors. This lower WUE is likely to reduce the ability of trees located in mature forests to tolerate water stress [4]. It has been observed in other systems that the positive influence of improved environmental conditions under the tree canopy is counteracted by a reduced soil moisture, particularly during dry years [45]. Thus, the increased shade at mature stages involves a potential risk during dry years. Furthermore, we observed reduced tree vigour in the mature forest compared to those at the expanding front. The changing nature of J. thurifera forests and the cessation of rural activities has been linked to defoliation not related to infection [46], with water stress being the most likely causal factor for the reduced vitality of these trees [47]. However, whether this reduced vigour is the cause or effect of the observed reduced WUE remains uncertain and deserves specific research due to its relevance to model forest dynamics in general, and under climate change scenarios in particular.

The negative relationship observed between δ¹³C (a proxy of WUE) and C/N ratio in leaves (a proxy of nitrogen use efficiency; NUE) reflects the trade-off existing between WUE and NUE, which is especially important in low-resource environments [48]. It has been shown that this trade-off is influenced by water and nitrogen availabilities in the soil [49,50] and could be explained by physiological constraints at the leaf level. High NUE (i.e., high rates of photosynthesis per unit of leaf N) can occur when stomata are open and, therefore, water loss by transpiration is high (low WUE) [47]. Here, we also found a negative correlation between C/N ratio and total nitrogen content in leaves as has been shown for several species [51], but we did not find a correlation between C/N ratio and total carbon content in leaves (Figure S1). Therefore, the trade-off found here could be due to differences in soil nitrogen availability. Although higher nutrient availability in soils could be expected at the expanding front as a legacy from past land-use (i.e., farmland or agricultural activities) [52], previous studies conducted in the area did not find large differences in nutrient availability among the stages of the forest expansion gradient [27]. It is important to explore in more detail the previous land-use changes and to increase the understanding of soil nutrient dynamics and spatial heterogeneity to interpret the patterns of WUE and NUE. Farmland activities alter soil properties (i.e., soil structure, organic matter content, and soil nutrient) [50] and can have different and complex effects on the performance of trees [7].

WUE could be modulated by endogenous factors of trees such as age [53] or gender in the case of dioecious species [54], in addition to microenvironmental characteristics around each tree as discussed previously. In our study, WUE was influenced by tree age but not by gender. Contrary to previous studies showing an increased WUE with age [53], we found younger individuals of J. thurifera had greater WUE. It has been reported that adults and juveniles of this species reduce net photosynthetic rate and stomatal conductance under unfavourable climatic conditions, and also that juvenile trees have greater plasticity of cambial activity in response to environmental constrains than adult trees [55,56]. We also found a positive relationship between age and C/N ratio (Figure S1). These relationships together suggest that older trees maintain high stomatal conductance and transpiration rates (lower WUE) while younger trees reduce both rates (stomatal conductance and transpiration rates) to increase resource use efficiency under adverse environmental conditions and increase survival. A decrease in WUE at the expanding front could be expected over time as tree age, but our experimental design cannot separate the effects of age and location along the colonisation gradient; this point deserves specific attention in future studies. Moreover, we found that the relationship between QD (a surrogate for growth; [57]) and age did not vary among stages (Figure 1). Thus, the maintenance of growth and the greater vigour of younger trees growing at the expanding front might be possible by an increase in their WUE. A detailed dendrochronological analysis of trees established in different stages could reveal whether an increase in WUE at the expanding front allows for a growth increase in these recently colonised areas [58].

5. Conclusions

Our study showed an important vulnerability of Juniperus thurifera forests to future scenarios of climate change. Both expanding fronts and mature patches of this system are vulnerable, but for different reasons. At the expanding front, rocky soils, together with greater irradiance, generate stressful conditions that will be exacerbated by climate change, compromising recruitment and further expansion of these forests. Trees in mature patches showed a lower WUE, suggesting a potentially reduced capacity to accommodate increased droughts.

The limitations of our experimental design, which cannot differentiate the effect of age and the stage of the forest expansion gradient on WUE because tree age differed between the three stages, points to the need for a mechanistic study with a different experimental design to clearly separate the correlated and potentially confounding factors that occur along gradients of forest expansion. Such studies are needed to improve our understanding of the sensitivity of these highly dynamic forest ecosystems to climate change scenarios.