Most Southern Scots Pine Populations Are Locally Adapted to Drought for Tree Height Growth

Most populations of Scots pine in Spain are locally adapted to drought, with only a few populations at the southernmost part of the distribution range showing maladaptations to the current climate. Increasing tree heights are predicted for most of the studied populations by the year 2070, under the RCP 8.5 scenario. These results are probably linked to the capacity of this species to acclimatize to new climates. The impact of climate change on tree growth depends on many processes, including the capacity of individuals to respond to changes in the environment. Pines are often locally adapted to their environments, leading to differences among populations. Generally, populations at the margins of the species’ ranges show lower performances in fitness-related traits than core populations. Therefore, under expected changes in climate, populations at the southern part of the species’ ranges could be at a higher risk of maladaptation. Here, we hypothesize that southern Scots pine populations are locally adapted to current climate, and that expected changes in climate may lead to a decrease in tree performance. We used Scots pine tree height growth data from 15-year-old individuals, measured in six common gardens in Spain, where plants from 16 Spanish provenances had been planted. We analyzed tree height growth, accounting for the climate of the planting sites, and the climate of the original population to assess local adaptation, using linear mixed-effect models. We found that: (1) drought drove differences among populations in tree height growth; (2) most populations were locally adapted to drought; (3) tree height was predicted to increase for most of the studied populations by the year 2070 (a concentration of RCP 8.5). Most populations of Scots pine in Spain were locally adapted to drought. This result suggests that marginal populations, despite inhabiting limiting environments, can be adapted to the local current conditions. In addition, the local adaptation and acclimation capacity of populations can help margin populations to keep pace with climate change. Our results highlight the importance of analyzing, case-by-case, populations’ capacities to cope with climate change.


Introduction
The impact of climate change on tree growth is complex, because it depends on many biotic and abiotic processes, and hence, it remains difficult to predict. Northern and western European forests may benefit from warmer temperatures, with positive effects on tree growth [1][2][3][4], but at the same time, tree mortality has increased during the last years [5]. Contrarily, southern and eastern European forests can experience reductions in their growth rates as a consequence of more frequent and intense episodes of drought [1,6]. In natural populations, tree growth strongly depends on climate and stand structure, and on processes related to competition and facilitation among individuals. This complexity

Climate Data
We used two sources of gridded climate data for the present and future conditions: (1) We used the Gonzalo-Jiménez [42] Spanish climatic model-the average climate was calculated for the period between 1961 and 1999-covering Spain at 30 arc sec resolution (~1 × 1 km), to calibrate our models. We used a total of 51 temperature and precipitation-related climate variables to accurately assess the main climatic drivers promoting among-population differentiation and plasticity (Table 1). We named the climatic variables of the planting site and the population origin _s or _p, respectively (e.g., AHM_s and AHM_p for the annual heat-moisture index of the planting site or at the population origin, respectively). Table 1. Temperature-and precipitation-related climate variables considered in the present study. Summer heat-moisture index SHM = MAXWMT/(PREC.sum/1000), based upon [10] jan, feb, mar, apr, may, jun, jul, ago, sep, oct, nov, dec are the abbreviations for the months of the year, from January to December. annl means annual. win, spr, sum and aut are the abbreviations for winter, spring, summer and autumn, respectively. (2) We used the WorldClim climate database [43,44]-the average climate was calculated for the period between 1970 and 2000-at 30 arc sec resolution (~1 × 1 km), to make spatial predictions of tree height growth for the current and future climates. The climate values for the planting sites and population origins for the current climate were represented by the average climate, as calculated for the period 1970-2000. For future climates, we averaged five global circulation models (GISS-E2-R, HadGEM2-AO, IPSL-CM5A-LR, MIROC-ESM, and NorESM1-M) available in WorldClim [44] for the Representation Concentration Pathway (RCP) 8.5 trajectory for the year 2070. RCPs are greenhouse gas concentration (not emissions) trajectories adopted by the IPCC in the last Assessment Report (AR5) [45]. The RCP 8.5 trajectory is the most aggressive trajectory, and it predicts an increase in the mean temperature of between 2.0 • C and 3.7 • C during the 21st century. Hence, using this scenario for our predictions will show the worst expected outcome. The climate values for the planting sites were represented by the average climate predicted for the year 2070, while the climate values for the populations origins were represented by the average climate of the period 1970-2000.

Selection of Climate Variables
To identify linear and non-linear co-variations between tree height growth and the climate variables, both between the planting sites and at the populations' origins, we used Pearson and Spearman correlation analyses. Based on the correlation results, we set a cutoff of ρ ≥ |0.5| to select the climate variables for the planting sites, and of ρ ≥ |0.08|, to select the climate variables for the populations' origins. For further analyses, all climate variables were standardized: i.e., the mean was subtracted from each value and divided by the standard deviation. We selected 29 climate variables for the planting sites, and five variables for the populations. Climate variables for the planting sites included temperature-and precipitation-related variables, while climate variables for the populations were mostly precipitation-related (Table A1 in Appendix A).

Linear Mixed-Effect Models of Tree Height Growth
We fitted linear mixed-effect models of tree height growth that accounted for both the climate of the planting site and the climate of the population origin [10,12,13]. To build the best supported linear mixed-effect model, we followed two steps: (1) We fitted a battery of linear mixed-effect models of pairwise combinations between one climate variable of the planting site, and another climate variable of the population origin from the set of previous selected climate variables to find the best combination. These models included the structure of the experimental design in the random part (two random effects, blocks nested into planting sites, and populations), and the linear and quadratic terms; and the linear interaction term of these climate variables in the fixed part [12] (Equation (1)). We selected the model with the lowest Akaike Information Criterion value (AIC) [46].
where H ijk is tree height growth of the ith individual of the jth population in the kth planting site. α s is the set of n parameters associated with the fixed effects of the model, clim_p ij is the climate at the population of origin of the ith individual of the jth population; clim_s ik is the climate at the planting site of the ith individual in the kth planting site. βs are the random effects. ε ijk is the residual distribution of the ith individual of the jth population in the kth planting site following a Gaussian distribution.
(2) We built the best supported linear mixed-effect model with the best combination of climate variables selected in the previous step. First, the random part remained invariable, and included the experimental design structure described above (Equation (1)). Second, we defined a full model that included the linear and quadratic terms of both the climate at the planting site and the climate at the population origin, and the linear interaction term of both terms. The model selection of the predictor variables in the fixed part was conducted by using a hierarchical backward selection procedure. We used the AIC criteria, following the rule that net increments of lower than two units of AIC associated with the elimination of any parameter in the full model determined the exclusion of the parameter from the final model [47,48]. We started with the selection of the two-variable interaction (round 1) and then tested the quadratic effects of both climate variables (rounds 2 & 3) and so on downwards the main effects of each predictor (round 4). Fixed effects were tested, using the maximum likelihood (ML), and random effects were tested using the restricted maximum-likelihood method (REML). We computed the percentage of variance explained by the fixed effects of the best supported model, MR 2 , and the percentage of variance explained by the fixed and random effects together, CR 2 , [49,50]. The goodness-of-fit of the best model was assessed by examining the predicted vs. observed values. We used normal error distribution with an identity link. We used R (version 3.2.3, 10 December 2015) run in the Linux-GNU operating system to perform all of the analyses, and the "lme4" and "lmerTest" packages [50,51].

Assessment of the Local Adaptation of Tree Height Growth for Spanish Scots Pine Populations under the Current Climate
Using the best supported linear mixed-effect model, we calculated for each population the tree height growth at its local environment (H L ), and the optimum tree height growth (H OPT ) attainable at the same environment by a non-local population. To calculate H L , we replaced the climatic values of the planting site by those of the population climate of origin (i.e., clim_s = clim_p). The estimation of H OPT was done in two steps. First, we identified the value of the climate of the population origin, providing the optimum tree height growth (CL OPT ), using the best supported mixed-effect model. To do that, we calculated the first-order partial derivate of the best supported linear mixed-effect model with respect to the climate variable of the population origin and settled it to zero [14]. Second, we used CL OPT in our best supported mixed-effect model to obtain H OPT .
The amount of local adaptation of tree height growth for each population was then calculated as follows: A positive value of LA H indicates that a non-local population would outperform the local one, suggesting that the local population is not locally adapted. Values of LA H that are equal to or close to zero indicate local adaptation. The higher the value of LA H is, the higher the degree of maladaptation is.
We computed the climate lags (CL H ) for each population associated with LA H values, to evaluate the climatic causes of maladaptation; i.e., if the local population does not reach the optimum tree height growth, it could be because it is currently living in a drier or warmer climate (negative values of CL H ) or in a wetter or cooler climate (positive values of CL H ).

Spatial Predictions of the Local Adaptation of Tree Height Growth for Southern Scots Pine Populations under the Impacts of the Current and Future Climate
We predicted the optimum tree height growth (H OPT ) and local tree height growth (H L ) at age 15, based on our models. Our predictions were performed across the species' range in Spain [52]. Then, we estimated LA H as the difference among them. Moreover, we predicted tree height growth by the year 2070 under the RCP 8.5. scenario (H L -FUT) [53]. Differences between H L -FUT and H L , (D FUT-PRES ), can inform us about the climate change impact on Scots pine populations. Positive values of D FUT-PRES mean that Scots pine populations will increase the tree height growth, while negative values would indicate a decrease. We used the 'raster' package in R for all of the spatial computations.

Linear Mixed-Effect Models of Tree Height Growth
First, from the battery of models fitted, the model with the lowest AIC included the spring precipitation of the planting site, PREC.spr_s, and the summer heat moisture index of the population origin, SHM_p, (Table A2). Second, the best supported linear mixed-effect model contained the linear term of the spring precipitation of the planting site, PREC.spr_s, and the linear and quadratic terms of the summer heat moisture index of the population origin, SHM_p, and the linear interaction between PREC.spr_s and SHM_p (Tables 2 and 3). The fixed effects of the model explained 62.36% (MR 2 ) of the variance, and 75.14% (CR 2 ) of the variance was explained by the fixed and random effects together. An examination of the residuals indicated that the main assumptions of linear mixed-effect models were met ( Figure A2). The spring precipitation of the planting sites (PREC.spr_s) had a positive effect on tree height growth ( Figure 1). The summer heat moisture index (SHM_p) drove the among-population differentiation in the studied populations. Specifically, populations that originated in sites with intermediate values of the summer heat moisture index, SHM, displayed the highest tree height growth, while populations that originated in sites with either high or low values of SHM presented lower tree height growth. This pattern varied across the PREC.spr_s gradient, due to the interaction term between PREC.spr_s and SHM_p (Table 3 and Figure 1). Asterisks indicate statistical significance at a p-value of lower than: 0.001 (***); 0.01 (**); 0.05 (*); 0.1 (.) or 1 (n.s., not significant). PREC.spr_s is the spring precipitation at the planting site. SHM_p is the summer heat moisture index at the population origin.

Assessment of the Local Adaptation of Tree Height Growth for Spanish Scots Pine Populations under the Current Climate
Most of the populations, except two, underperformed at their local environment, as shown by LAH values that were different from zero (Table 4). However, the LAH values were small, and they were specifically below 5 cm in 11 out of the 16 populations tested ( Table 4). Most of the populations were currently growing under drier conditions, as we found negative values in climatic lags CLH. This may prevent the optimum tree height growth from being reached. Two out of the 16 populations were growing under wetter conditions, as we found positive values in CLH (Table 4).

Assessment of the Local Adaptation of Tree Height Growth for Spanish Scots Pine Populations under the Current Climate
Most of the populations, except two, underperformed at their local environment, as shown by LA H values that were different from zero (Table 4). However, the LA H values were small, and they were specifically below 5 cm in 11 out of the 16 populations tested ( Table 4). Most of the populations were currently growing under drier conditions, as we found negative values in climatic lags CL H . This may prevent the optimum tree height growth from being reached. Two out of the 16 populations were growing under wetter conditions, as we found positive values in CL H (Table 4).

Spatial Predictions of the Local Adaptation of Tree Height Growth for Southern Scots Pine Populations under the Impacts of the Current and Future Climate
Most southern Scots pine populations are locally adapted (LA H values close to zero) or slightly maladapted (LA H values of around 5; gray to light-yellow colors in Figure 2) to current climates. Just a few populations of Scots pine at the southernmost part of the distribution presented large values of LA H (yellow to red colors in Figure 2).

Spatial Predictions of the Local Adaptation of Tree Height Growth for Southern Scots Pine Populations under the Impacts of the Current and Future Climate
Most southern Scots pine populations are locally adapted (LAH values close to zero) or slightly maladapted (LAH values of around 5; gray to light-yellow colors in Figure 2) to current climates. Just a few populations of Scots pine at the southernmost part of the distribution presented large values of LAH (yellow to red colors in Figure 2). Our model predicted an increase in tree height growth in year 2070, particularly in the north (Pyrenees), and in high-elevation areas (Figure 3). The model predicted a decrease in tree height growth in the center and in the north-west. Populations at the southernmost part of the distribution were predicted to decrease moderately (Figure 3). Our model predicted an increase in tree height growth in year 2070, particularly in the north (Pyrenees), and in high-elevation areas (Figure 3). The model predicted a decrease in tree height growth in the center and in the north-west. Populations at the southernmost part of the distribution were predicted to decrease moderately (Figure 3).

The Main Climatic Drivers Shaping Among-Population Differentiation and Phenotypic Plasticity Responses of Tree Height Growth
The main climatic driver shaping genetic differentiation among the populations, in terms of tree height growth, was the summer moisture index (SHM_p), a proxy of drought. This result suggests that drought has probably been a selective factor, which is in agreement with previous studies [39,54]. For instance, our model predicted that populations that originated under more stressful conditions, i.e., with higher values of the summer moisture index calculated for the period of years between 1961 and 1999 (Figure 1), would present a slower rate of growth, which is a suggested strategy for increasing drought adaptation by reducing aboveground biomass [39,55]. Other adaptive responses to reduce cavitation risks, induced either by drought or frost stresses, have been described for the same populations. For example, populations from dry sites show large average tracheid lumen diameters, to assure hydraulic conductivity, whilst populations from sites with frequent freeze-thaw events present small tracheid lumen diameters [53]. Also, it has been suggested that southern populations adapt their leaf conductance and stomata control for less transpiration through the needles [54].
Phenotypic plastic responses were driven by spring precipitation (PREC.spr_s). We found a common pattern where tree height growth decreases under drier conditions (Table 3 & Figure 1). This result agrees with the findings of a previous study, where annual precipitation drove phenotypic plastic responses in tree height-diameter allometry for the same species [40]. Now, this study has allowed us to identify that tree height growth is largely driven by the rainfall that falls during the spring season. This could be related to the fact that maximum cambial activity for this species has been described to take place between mid-March and end of August [55], suggesting that Scots pine populations make better use of the rainfall in the early growing season than in the late season (autumn).
Both results highlight that water availability drives the phenotype variations in tree height growth in Scots pine populations, which is consistent with the limiting role that water plays in Mediterranean ecosystems [56].

Most Southern Scots Pine Populations are Locally Adapted to the Current Climate
The finding of local adaptations for most of the southern Scots pine populations was in agreement with our expectations (Table 4 & Figure 2). One possible explanation could be related to the species' demographic history. During the last glacial maximum, Scots pine had many refuge areas

The Main Climatic Drivers Shaping Among-Population Differentiation and Phenotypic Plasticity Responses of Tree Height Growth
The main climatic driver shaping genetic differentiation among the populations, in terms of tree height growth, was the summer moisture index (SHM_p), a proxy of drought. This result suggests that drought has probably been a selective factor, which is in agreement with previous studies [39,54]. For instance, our model predicted that populations that originated under more stressful conditions, i.e., with higher values of the summer moisture index calculated for the period of years between 1961 and 1999 (Figure 1), would present a slower rate of growth, which is a suggested strategy for increasing drought adaptation by reducing aboveground biomass [39,55]. Other adaptive responses to reduce cavitation risks, induced either by drought or frost stresses, have been described for the same populations. For example, populations from dry sites show large average tracheid lumen diameters, to assure hydraulic conductivity, whilst populations from sites with frequent freeze-thaw events present small tracheid lumen diameters [53]. Also, it has been suggested that southern populations adapt their leaf conductance and stomata control for less transpiration through the needles [54].
Phenotypic plastic responses were driven by spring precipitation (PREC.spr_s). We found a common pattern where tree height growth decreases under drier conditions (Table 3 & Figure 1). This result agrees with the findings of a previous study, where annual precipitation drove phenotypic plastic responses in tree height-diameter allometry for the same species [40]. Now, this study has allowed us to identify that tree height growth is largely driven by the rainfall that falls during the spring season. This could be related to the fact that maximum cambial activity for this species has been described to take place between mid-March and end of August [55], suggesting that Scots pine populations make better use of the rainfall in the early growing season than in the late season (autumn).
Both results highlight that water availability drives the phenotype variations in tree height growth in Scots pine populations, which is consistent with the limiting role that water plays in Mediterranean ecosystems [56].

Most Southern Scots Pine Populations are Locally Adapted to the Current Climate
The finding of local adaptations for most of the southern Scots pine populations was in agreement with our expectations (Table 4 & Figure 2). One possible explanation could be related to the species' demographic history. During the last glacial maximum, Scots pine had many refuge areas scattered throughout southern Europe [57,58]. These populations may have adapted to the particular conditions of this part of the range during the postglacial migrations [24,33]. Furthermore, the fragmented distribution of the species in the southern part of the range, and the low levels of gene flow among populations might have favored local adaptation [31]. However, we found maladaptation in the southernmost populations (Figure 2). This result could be explained accordingly by the gene flow asymmetry theory across species' ranges [14]: maladapted populations might have received alleles from core and northern populations that are pre-adapted to cooler or wetter conditions than those found at the southernmost part of the range, limiting local adaptation [10,14].
Nevertheless, our findings of local adaptation may change if other fitness-related traits were considered. For instance, fitness-related traits can present trade-offs across the species' ranges [21,23]. Typically, demographic compensation has been found in several marginal populations [21][22][23]. We could expect the same for Scots pine, given that the southernmost Spanish populations show low early recruitment as a consequence of high seedling mortality, low seed production, and high predation rates [56]. Hence, a wider perspective of Scots pine adaptation patterns using other traits than those related with growth is desirable.

The Importance of Considering Genetic and Plastic Effects for Evaluating Tree Height Growth for Future Climates
Our predictions should be interpreted as the combination of adaptation and plasticity effects in tree height growth. The contribution of genetic effects to explain tree growth variation in our data was lower than that of phenotypic plasticity, although it was of similar magnitude (∆AIC = 28.69 and ∆AIC = 33.94, respectively, Table 1). This result highlights that both adaptations to local conditions and the capacity to adjust to the environment through plastic responses are crucial for understanding the performance of Scots pine populations for future climates.
Future predictions for an average tree, showing an increase in tree height, suggest that plasticity can help with tree acclimation to the new climate, and it can then compensate for the imprint of local adaptation, which is a common characteristic in trees [57]. Our results are partially in agreement with [27], where moderate increases in tree growth for Scots pine populations are predicted in Spain under the RCP 8.5 scenario for 2050, although tree growth does decrease, especially in southern low-elevation sites, in 2070. In the northern and high-elevation areas (Pyrenees), warmer conditions would benefit tree growth by lengthening the growing season, thus increasing net photosynthesis, stomatal conductance, and specific hydraulic conductivity [58]. However, opposite results of tree growth decrease have been reported for this species in Spain [35,59], as lower precipitation regimes combined with higher temperatures would presumably lead to a general pattern of more frequent drought events, higher evapotranspiration, and reduced soil water availability in the future. However, none of these studies have accounted for the plasticity capacity of the species, and this could explain the differences in the forecasts. Moreover, few populations in the north-west and in the central part of the species distribution in Spain have been predicted to decrease in tree height growth. This prediction could be reflected by the sharp decrease in spring precipitation that is expected for year 2070 (RCP 8.5.) in these areas ( Figure A3). Our forecast could be related to site characteristics (such as soil depth, nutrient availability, geomorphological features, etc.) that can make water availability difficult.
Finally, caution is needed when interpreting our predictions beyond the year 2070, as we do not know whether plasticity will help with acclimation to new environments in the far future, and evolutionary adaptation may then become necessary. Contrary to our expectations, in terms of the species' growth, the forecast for southern Scots pine populations is generally favorable. Nonetheless, our results are based on common gardens. In natural populations, other factors, such as above-and below-ground competition, herbivory, management practices, etc., can strongly modify our expectations of Scots pine growth and adaptation to climate change.

Conclusions
Most of the populations of Scots pine in Spain were locally adapted to drought. This result suggests that marginal populations, despite inhabiting limiting environments, can become adapted to their current local conditions. In addition, the local adaptation and acclimation capacity of populations can help marginal populations to keep pace with climate change. Our results highlight the importance of analyzing populations' capacity to cope with climate change, on a case-by-case basis. Chambel-who measured, gathered, and stored the data used in this research article.

Conflicts of Interest:
The authors declare no conflict of interest. T.#, TMC.# and TMF.# mean minimum, maximum and mean average monthly (#) temperatures ( • C); PREC.ann means annual precipitation (mm); PREC.win, PREC.spr, PREC.sum, PREC.aut, for winter, spring, summer and autumn, respectively; TM means mean annual temperature ( • C); MINCMT means mean of the minimum temperatures from the coldest month ( • C); FP means frost period (month); DD5 means degree days-period over 5 • C; TD means continentality; AHM means annual heat moisture index; and SHM means summer heat moisture index.  Figure A1. Network of provenance common gardens in Spain (GENFORED): Planting sites (green crosses), populations (pink circles). The species distribution range is shown in gray (two isolated populations are represented by "+" in gray).