Geographic variation in growth, survival, and susceptibility to the processionary moth (Thaumetopoea pityocampa Dennis & Schiff.) of Pinus halepensis Mill. and P. brutia Ten.: results from common gardens in Morocco

Abstract

Key message

Adaptation of Pinus brutia and Pinus halepensis to harsh Moroccan environments varied considerably among populations, with variation following marked geographic clines.

Context

Mediterranean pines are drought tolerant species that play relevant ecological and economic roles in North Africa. Mediterranean pines harbor huge intraspecific variation in adaptive traits.

Aims

Exploring the relative performance of different seed sources and analyzing the genotype by environment (G × E) interaction becomes, thus, essential for maximizing the efficiency of any reforestation program.

Methods

We present the results of a comprehensive quantitative analysis of growth, survival, and resistance to the processionary moth of 57 populations of P. halepensis Mill. and P. brutia Ten. tested under severe drought conditions in two contrasting Moroccan test sites differing in water regime, soil, and physiographic conditions.

Results

The results indicated the existence of considerable variation between species and among populations within species in growth, survival, and susceptibility to the processionary moth (Thaumetopoea pityocampa Dennis & Schiff.), but low plasticity across sites, despite the apparently large environmental differences between them. On average, P. brutia performed better than P. halepensis. The G × E interaction was significant for all traits and ages, but an in-depth analysis revealed that the main source of the interaction was due to large differences in residual variances across sites, while the relative population ranks remained highly consistent across sites.

Conclusion

Growth and survival of Aleppo pine populations followed well-defined geographical clines, with populations from the Eastern range of the distribution, particularly those of Greece, showing outstanding performance. The two local Aleppo pine populations were, however, also among the best-growing and surviving populations. Superiority of Eastern Aleppo pine populations was not considered high enough to recommend their use in reforestation programs. Aiming to preserve local genetic resources, local seed sources are preferred for reforestation in the dry areas of Morocco. Planting middle to high elevation Turkish populations of the non-native P. brutia is also an alternative, especially in more mesic sites. The reduced G × E suggests that the above recommendations may be valid, even in the context of future climate change.

1 Introduction

Forest tree species, and particularly Mediterranean pines, usually harbor extremely large intraspecific genetic diversity (Fady 2012; Soto et al. 2010). Current intraspecific variability is essentially the result of demographic processes that have occurred during species redistribution after historical bottlenecks, modulated by adaptive responses to the new local environmental conditions, which are particularly heterogeneous in the Mediterranean basin (Alrababah et al. 2010; Grivet et al. 2009). Population isolation together with high frequency of recurrent extreme events (e.g., fire, drought, biotic outbreaks) has also contributed to the evolutionary history of Mediterranean pines, increasing stochastic differentiation among populations by means of genetic drift (Budde et al. 2017; Tapias et al. 2004). On the other hand, long-distance gene flow via pollen dispersal in these outcrossing wind-pollinated species is likely to have counteracted population differentiation and helped to maintain high levels of within-population diversity (Savolainen et al. 2007). Human exploitation of genetic resources of Mediterranean pines since ancient time has largely eroded within-population variation, whereas large-scale planting over the last century has reshuffled the current genetic structure (Richardson et al. 2007). As a result, Mediterranean pines show complex geographic patterns of intraspecific genetic variation which vary among the different species depending on the relative relevance of all these evolutionary forces (Fady 2012; Grivet et al. 2013; Soto et al. 2010).

Pinus halepensis Mill. (Aleppo pine) and Pinus brutia Ten. (Turkish red pine or Brutia pine) are closely related taxa of the Mediterranean Basin that were once considered as subspecies of the complex halepensis/brutia group but now are recognized as independent species (Chambel et al. 2013; Phillips and Gladfelter 1991). The two species are photophilous fast-growing colonizers especially adapted to low moisture conditions. P. halepensis is widely distributed across the Mediterranean region from Morocco to Syria, being the most important forest tree species of North Africa. P. brutia spreads naturally across the Eastern Mediterranean basin, being the most abundant conifer species in Turkey. Both species encounter each other in some areas of Greece and Turkey, where they naturally hybridize, with unidirectional introgression from P. halepensis to P. brutia (Bucci et al. 1998). Between them, P. halepensis and P. brutia cover more than 7 million ha in the Mediterranean basin (Chambel et al. 2013). Despite their phylogenetic proximity, the two species showed marked divergence for several key adaptive traits (Grivet et al. 2013), resulting in contrasting performance when planted together in different environmental conditions (Bariteau and Pommery 1992; Fisheri et al. 1986; Kaundun et al. 1997; Weinstein 1989). Aleppo pine is generally more tolerant of water stress, while Brutia pine is comparatively more tolerant of frost damage (Chambel et al. 2013).

Both Aleppo and Brutia pines show marked population differentiation but relatively lower variability within populations (Grivet et al. 2009), with both genetic identity and diversity commonly following geographical clines (Grivet et al. 2009; Kurt et al. 2011b). Genetic diversity of Aleppo pine tends to diminish from East to West following the suggested pathway of postglacial recolonization (Gomez et al. 2000; Grivet et al. 2009). In Brutia pine, altitude has been shown to be an important driver of population differentiation (Kurt et al. 2011b). Common garden tests have also recurrently shown large population differentiation in different life history traits in both Aleppo (Calamassi et al. 2001; Climent et al. 2008; Santos del Blanco et al. 2013; Voltas et al. 2015) and Brutia pines (Bariteau and Pommery 1992; Isik et al. 2000; Kurt et al. 2011b). In accordance with the assumed relevance of local adaptation processes, variation among populations frequently shows marked geographical, altitudinal, or climatic clinal patterns (Dangasuk and Panetsos 2004; Kurt et al. 2011b). For example, longitudinal and latitudinal gradients have been reported for several adaptive traits in Aleppo pine, with northern and eastern populations prioritizing growth at the expense of reproduction (Climent et al. 2008) and water use efficiency (Voltas et al. 2008; Voltas et al. 2015). In Brutia pine, quantitative differentiation for height growth among populations has been shown to follow marked altitudinal clines (Kurt et al. 2011b), with populations from middle to high elevations tending to grow faster than populations from the lower peripheral areas (Isik et al. 2000). Differentiation in quantitative traits among P. brutia populations, commonly greater than that for neutral markers, has also been attributed to local adaptation processes (Kurt et al. 2011a). Based on these results and in order to conserve local genetic resources, the abovementioned authors suggested that the selection of material for reforestation should be based on the geographic origin along altitudinal gradients (Kurt et al. 2011a). It should be noted, however, that most of the previous work studying among-population variation in P. brutia has focused on Turkish populations, but the distribution of the species is considerably wider than that (Kurt et al. 2011b). Differentiation between populations was also observed in other areas, again following marked geographical clines (Dangasuk and Panetsos 2004).

Frequently, the observed differences among populations in adaptive traits are environmentally dependent, and the relative performance of the different seed sources rank differently depending on the testing site, resulting in a strong genotype by environment interaction (G × E) (Atzmon et al. 2004; Isik et al. 2000; Santos del Blanco et al. 2013; Taïbi et al. 2015). Formal quantifications of the practical relevance of the observed G × E patterns are, however, scarce (Isik et al. 2000; Taïbi et al. 2015). Studying in-depth the G × E patterns is essential in order to understand which are the main environmental drivers responsible for the differences in population performance, something vital not only for delimiting seed transfer zones (e.g., De la Mata and Zas 2010) but also for mitigating the impact of climate change (e.g., Taïbi et al. 2014).

One interesting point that is receiving increasing attention is that the local origin is not always the best adapted one to a given environment (Leimu and Fischer 2008). Many factors can explain this lack of local adaptation, including a strong contribution to population differentiation of evolutionary forces other than natural selection such as genetic drift associated to founder effects or other bottlenecks (Lenormand et al. 2009), limited variability within local populations (Schiller and Atzmon 2009), and changes in environmental limiting factors across time and space (Jaramillo-Correa et al. 2014). In particular, climate change may cause important shifts in environmental conditions that may affect the relative performance of the local and non-local populations (e.g., St Clair and Howe 2007). Thus, the optimal suitability of local seed sources for reforestation plans under climate change is being increasingly challenged, while assisted migration is slowly gaining some credit as a potential tool for mitigating its effects (Aitken and Whitlock 2013).

Exploring the performance of different seed sources from the whole distribution area under current conditions and quantifying and analyzing the patterns of G × E across the environmental variation of a particular region are, thus, essential for maximizing the efficiency of future reforestation programs. This is particularly important for traits related to tolerance to abiotic (e.g., extreme drought or frost events) and biotic (e.g., pest outbreaks) stress factors. Failing to account for the intraspecific genetic variation in tolerance to abiotic and biotic stresses can result in dramatic losses associated with seed transfer. Well-known examples include the extremely high mortality of Pinus pinaster Aiton in Aquitaine (France) in the 80s due to the use of a Portuguese provenance intolerant to frost (Benito-Garzón et al. 2013) or the massive failure of plantations of P. pinaster in Provence (France) due to the use of a seed source (Corsican provenance) highly susceptible to the bast scale Matsucoccus feytaudi Duc. (Di Matteo and Voltas 2016). The processionary moth (Thaumetopoea pityocampa Dennis & Schiff.) is the main insect defoliator of Mediterranean pines and can cause important growth losses and overall decline (Jacquet et al. 2012). Despite its major relevance as a factor causing decline in Mediterranean pines—especially under current global warming (Robinet et al. 2007)—little effort has been spent to date in exploring the extent to which intraspecific variation of Mediterranean pines also expands to variation in their susceptibility to this insect (Bariteau and Pommery 1992).

Specifically, in Morocco, where Mediterranean pines (mainly P. pinaster and P. halepensis) cover extensive areas and provide a wide array of ecosystem services, little information about the intraspecific adaptive variation and its patterns is available (Destremau 1974). Additionally, Moroccan populations of Aleppo pine are known to have low diversity levels as a consequence of ancient bottlenecks (Salim et al. 2010). This reduced variation may compromise the suitability of local origins for use across (spatial or temporal) environmental heterogeneity within the country (Boulli et al. 2001). In recent decades, extensive reforestation efforts with pine trees (mainly Aleppo pine) have been carried out in the country to restore denuded land (Wahid et al. 2010). To identify Mediterranean pine genotypes adapted to the severe and stressful conditions of these areas, several provenance tests of drought resistant pine species (P. halepensis and P. brutia) were established in the early 1990s (Chambel et al. 2013). This paper presents a comprehensive quantitative analysis of growth, survival, and susceptibility to the processionary moth in two provenance trials of the above two species across 21 years after planting. The study of age-age correlations and variance trends over time together with the analysis of the sources of G × E interaction and their practical implications will help to better understand the patterns of adaptive variation to harsh Moroccan environments. Given that the main goal of reforestation efforts in Morocco is to provide ecosystem services such as water regulation and soil protection rather than to produce timber, special attention should be given to adaptive traits such as survival and resistance to biotic stress, as well as to the conservation of local genetic resources.

2 Material and methods

2.1 Plant material

Genetic material collected from 57 populations (8 of Brutia pine and 49 of Aleppo pine) and from seven different countries (Italy, Tunisia, Spain, Turkey, France, Greece, and Morocco), covering most of the natural range of the species under study, was planted in two field trials located in different environments (Table S1). The trial at Chatba included 56 populations, while the one at Izarene comprised 52, with 51 common to both trials (Table S2 and Fig. 1).

Fig. 1

Natural distribution of P. halepensis (blue area) and P. brutia (green area) and geographical location of the sampled populations of the two species (circles: P. halepensis; squares: P. brutia) and of the two provenance trials in Morocco (red stars). The distributions maps were obtained from EUFORGEN 2009, www.euforgen.org

For each pine population, open-pollinated seed was collected at the site of origin from 20 to 30 mother trees that were located at least 100 m apart to minimize relatedness among them. Collected seeds from each origin were mixed and sown in 1.256-cm³ (plastic bag) containers in the dar Essalam (Rabat) nursery, using standard containerized nursery practices. After growing in the nursery for a year, the seedlings were planted following conventional planting procedures, in February and March of 1992.

2.2 Test sites

Although the two test sites (named Chatba and Izarene) fall within the same climatic region (Peel et al. 2007), they differ widely in environmental conditions, Chatba being much drier than Izarene (Table S1). Chatba has a semi-arid climate (around 400 mm of annual precipitation) with mild winters and very dry summers (average summer precipitation < 20 mm). Humidity from the Atlantic Ocean (50 km from the sea), however, alleviates the water deficit. The site is completely flat while the soil is deep with few stones on the surface. The soil is dominated by clay and presents a small ferruginous nodule surmounting a friable limestone tuff at an average depth of 1.3 m. Annual grasses dominate the understory layer and are likely to pose very little competition to the experimental trees (Fig. S1A). The Izarene site has a subhumid Mediterranean climate with temperate winters. Annual precipitation in Izarene (900 mm per year) is more than double that in Chatba, but summers remain extremely dry (around 40 mm summer precipitation), and the Atlantic influence is lower (Table S1). The soil is derived from schists covered with colluvial clay sand, being stonier and shallower than that in Chatba. The Izarene site presents a pronounced slope, and the secondary shrub vegetation caused major competition to the experimental trees (Fig. S1B). Altogether, climatic water stress is supposedly greater in Chatba, but edaphic, geomorphic, and other environmental characteristics (e.g., competing vegetation, soil preparation, etc.) in Izarene are less favorable for pine growth.

Before planting, the Chatba site was tilled with a disc ripper. In Izarene, because of the sloping terrain, the preparation of the soil was restricted to the opening of holes of 50 × 50 cm following the contour lines. Both trials followed a complete randomized block design with 9 and 10 replications in Chatba and Izarene, respectively. Experimental plots in Chatba included nine trees (squares of 3 by 3) (Table S1). Three out of the 10 replicates in Izarene had experimental plots of nine trees as in Chatba, while the remaining seven included five aligned trees. In total, 7903 trees were planted in the two test sites. Spacing was 3 × 3 m, except in a small number of situations in the Izarene trial where spacing was reduced up to 2 × 2 m because of the irregularities of the terrain.

2.3 Measurements

All trees were assessed at both sites for survival, diameter at breast height (DBH), and total height ca. every 3–5 years during the period 1992–2013. Height and survival were recorded in the two sites at 2, 7, 10, 13, 15, 18, and 22 years of age, with an extra assessment in Chatba at the age of 4. Diameter at breast height was measured in both sites at ages 13, 15, 18, and 22.

In 2004, the processionary moth, probably the most destructive insect herbivore of Mediterranean pines, caused severe defoliation in the two test sites. In order to evaluate the degree of population variation in susceptibility to this important forest pest, all trees in the two test sites were assessed for the defoliation caused by the insect using a visual subjective score (0, 0–20, 20–40, 40–60, 60–80, 80–100% of the crown defoliated).

2.4 Climate data

For each pine population, climate variables of the region of origin were obtained from the Global Climatic Model “WordClim” Database for the period 1960–1990 (www.wordclim.org, spatial resolution of 1 km²). Six bioclimatic variables were used in the present study to characterize the climate at the region of origin: annual mean temperature (T), maximum temperature of the warmest month (Tmax), minimum temperature of the coldest month (Tmin), temperature range (TAR = Tmax–Tmin), annual precipitation (P), and annual summer precipitation (Ps). These variables have been previously shown to be relevant drivers of population differentiation in Mediterranean pines (Climent et al. 2008; Voltas et al. 2015).

2.5 Statistical analyses

2.5.1 Joint analysis across sites

Joint analysis of growth and defoliation data from the two sites was carried out for each assessment age by fitting the following general linear mixed model in SAS:

$$ {Y}_{ijklm}=\mu +{S}_i+{SP}_j+{POP}_k\left({SP}_j\right)+{S}_ix\kern0.5em {SP}_j+{S}_ix\kern0.5em {POP}_k\left({SP}_j\right)+{b}_l\left({S}_i\right)+{b}_l\left({S}_i\right)x\kern0.5em {POP}_k\left({SP}_j\right)+{\varepsilon}_{ijklm} $$

(1)

were Y_ijklm is the height or diameter of each single tree, μ is the overall mean, S_i, SP_j, and POP_k(SP_j) are the fixed effects (in uppercase) of the sites i, the species j, and the population k within the species j, respectively, S_i x SP_j and S_i x POP_k(SP_j) are the fixed interactions between the corresponding effects, and b_l(S_i) and b_l(S_i) x POP_k(SP_j) are the random effects (in lowercase) of blocks within sites and plots within blocks, respectively. This model analyzed species and provenance variation with the appropriate degrees of freedom considering the multi-tree plot design, in which trees within plots are not completely independent measures. The model was fitted with the PROC MIXED procedure of the SAS System (Littell et al. 2006) for each trait and assessment age (2, 7, 10, 13, 15, 18, and 22) separately using individual spatially unadjusted raw tree data. Heterogeneous residual variance models were used when needed. Only the 51 populations common to both sites were included in these analyses. Best linear unbiased estimators (BLUEs) of fixed effects were derived from these models for characterizing the performance of each population and species in each site.

Survival (0/1) at each assessment age was analyzed across sites with a generalized mixed model (PROC GLIMMIX in SAS; (Littell et al. 2006)) similar to the mixed model described above, assuming a binomial distribution and a logit link function.

2.5.2 Analysis of the G × E

To explore the causes of the G × E interaction, we used a likelihood approach similar to that used in De la Mata and Zas (2010). Firstly, full mixed models were fitted considering site as a fixed factor and populations and blocks within sites as random factors (species were not considered in these models; the 51 common populations of the two species were treated as different levels of a single factor):

$$ {\overline{Y}}_{ijk}=\mu +{S}_i+{pop}_j+{S}_i\kern0.5em x\kern0.5em {pop}_j+{b}_k\left({S}_i\right)+{\varepsilon}_{ijk} $$

(2)

In these full models, which were fitted independently for each assessment age, plot means $\left({\overline{Y}}_{ijk}\right)$ were used to avoid computing power issues, and all possible causes of G × E interaction were considered and analyzed, that is, heterogeneity of residual variance, which reflects within-site microenvironmental variation; heterogeneity of genetic variance related to scalar effects across sites; and deviations from perfect correlation among sites, which can be associated to genotype rank changes across environments, i.e., cross-over interaction (Yang 2002). Then, reduced models restricting the R and G variance-covariance to constraint the different sources of the interaction were fitted (see details in Table S4). Hypothesis testing regarding the constraints imposed on the population/residual covariance structure was done by comparing the restricted log-likelihoods (RLL) of the constrained model and the unconstrained model (see Table S4 for details). Under the null hypothesis that the full covariance model is not different from the reduced covariance model, the log-likelihood ratio LLR = − 2(RLL_{reduced model} − RLL_{full model}) is distributed approximately as χ² with degrees of freedom equal to the difference among the number of covariance parameters specifying the full model and the reduced model (Fry 2004).

The analysis of the G × E patterns in the case of survival was conducted with the same likelihood procedure having previously transformed the plot mean values (arcsin transformation) to achieve normality.

2.5.3 Age trend analyses

Age trends of population variation were explored by plotting variance components and population correlations across sites (as estimated from the previous full models for each assessment age, model [2]). Plots of population ranks against age for each assessed trait were also constructed.

2.5.4 Performance of the local origins and geographical and climatic clines

The superiority of the local Aleppo pine population (Chatba and Izarene, see Table S2) in relation to the average of all the remainder populations of the species was tested by means of specific contrasts in the PROC MIXED (growth traits) and PROC GLIMMIX (survival) (Littell et al. 2006), using single site analysis at age 22 (same nomenclature as before):

$$ {Y}_{ijk}=\mu +{POP}_i+{b}_j+{POP}_i\;x\kern1em {b}_j+{\varepsilon}_{ijk} $$

(3)

Simple correlations were calculated to explore the relationships among the performance of individual populations in each site with the geographic and bioclimatic variables of the place of the origin of each population. These analyses were carried out independently for the populations of P. halepensis and P. brutia using the data from age 22 for growth and survival and from age 13 for processionary moth defoliation.

Data availability

The phenotypic data have been deposited in DIGITAL.CSIC repository: http://hdl.handle.net/10261/163608; DOI https://doi.org/10.20350/digitalCSIC/8535 (Sbay and Zas 2018).

3 Results

3.1 Intra- and interspecific variation

Pine establishment was excellent in both sites with 86.1 and 79.0% of the trees surviving 21 years after planting (age 22) in Chatba and Izarene, respectively. At age 22, surviving trees averaged 10.6 ± 0.25 and 10.4 ± 0.4 m in height and 15.0 ± 0.06 and 15.6 ± 0.08 cm in diameter (mean ± SE), in Chatba and Izarene, respectively, with a standing basal area of 15.7 and 14.6 m² ha⁻¹. At this point, average growth and survival across species was not statistically different between sites (Table 1).

Table 1 Summary of the mixed models for the analysis of height, diameter and survival at age 22 after planting, and defoliation by the processionary moth at age 13, showing the degrees of freedom (DF), F ratios (F), and associated probability levels (P > F) for the fixed factors site, species, populations within species, and their interactions

Survival and growth at age 22 were highly variable among species and among populations within species (Table 1). P. brutia was consistently the species with the higher growth and survival rates on both sites (Fig. S2). Survival of the tested populations varied from 69.1 to 97.5% in Chatba and from 52.1 to 96.9% in Izarene (Fig. 2d). Although P. brutia populations showed the highest average survival rates on both sites, some P. halepensis populations reached even higher survival rates, especially in the Izarene test site (Fig. 2d and Table S3).

Fig. 2

Environment-Environment bi-plots representing the height (a) and diameter (b) growth, the defoliation by the processionary moth (c), and the survival (d) of the 51 common populations in the two test sites. Growth and survival was assessed at age 22 and processionary defoliation at age 13. Symbols denote the species (squares: P. brutia; circles: P. halepensis) and colors represent the country of origin. The dotted diagonal line denotes the 1:1 relationship (i.e., lack of plasticity). Codes for outstanding populations are also included (see Table S2)

Height growth varied 1.6- and 1.5-fold among populations, while diameter growth varied 1.9- and 1.8-fold in Chatba and Izarene, respectively (Fig. 2a, b and Table S3). Again, despite the mean growth superiority of P. brutia, some Aleppo pine populations from Greece (e.g., Halkida, Kassandra, Elea, Poligiros, Istaea) stood out as the best growers of all (Fig. 2). Rank correlations between population growth and survival at age 22 were positive and significant in both test sites (r = 0.36, p = 0.006 and r = 0.45, p < 0.001 for height and diameter in Chatba (N = 56) and r = 0.53, p < 0.001 and r = 0.44, p = 0.001 for height and diameter in Izarene (N = 52), respectively).

Defoliation due to the processionary moth in 2014 was high in both sites, with needle loss averaging 24.6 and 37.4% in Chatba and Izarene, respectively. No significant differences between sites but strong differences between species and among populations within species were observed in the damage caused by the processionary moth (Table 1). Variation among populations was especially evident in Izarene, where the range of variation was much larger than in Chatba (Fig. 2c and Table S3). All the P. brutia populations were consistently the least affected by the insect, while Turkish and Tunisian P. halepensis populations were the most affected (Fig. 2c).

3.2 Age trends and age to age correlations

Although average height growth and survival at age 22 was not statistically different between sites (Table 1), significant differences were observed at earlier ages (Table S5), with slightly greater growth and survival in Chatba than in Izarene. Judging from the F ratios of the mixed models, variation among species and among populations within species became larger with age (Table S5) with population and residual variances largely increasing as trees aged in the two test sites (Fig. 3). However, while population variances were quite similar and followed similar trends across sites, residual variances were drastically higher in Izarene for the three traits, and in the case of height growth, differences in residual variances between sites drastically increased with age (Fig. 3). Population correlations between sites for growth traits (height and diameter) were consistently high across ages (r around 0.90; Fig. 3). Regarding survival, correlation across sites was very low at early ages and tended to increase up to moderate values (r = 0.40) as trees aged (Fig. 3). Estimates of site to site correlation for survival were not very accurate because of large standard errors.

Fig. 3

Trends of variance components for height (a), diameter (b) and survival (c) across ages. Note the similar population variances but the large differences in the residual variances between the two test sites

A high correspondence among ages in population relative performance can also be observed when comparing the population ranks across ages (Fig. 4). With the exception of the very early age assessment of height and survival, population rank changes were not very important and tended to diminish with age for the three assessed traits (Fig. 4).

Fig. 4

Population rank changes as the trees aged in the two provenance trials (Chatba (left panels) and Izarene (right panels) for height (up), diameter (middle), and survival (bottom). Measurement ages are indicted in the x-axis. The color of the lines identifies the species. The local (Aleppo pine) population within each site is also identified

3.3 Genotype by environment interaction

Variation between species and among populations within species differed across sites as revealed by the corresponding significant interaction for all traits and ages (Tables 1 and S5). Likelihood ratio tests identified heterogeneity of residual variances as the main cause of the interaction in all cases (Tables 2 and S6). Population variance, on the contrary, was consistently homogeneous across sites for all traits (Table 2) and ages (Table S6). Models assuming a perfect correlation between sites (r = 1) fitted significantly worse than the full models (Tables 2 and S6), suggesting that deviations from perfect correlations are also contributing to the G × E interaction. However, judging from the likelihood of the corresponding reduced models, the relative contribution of this last source of interaction was comparatively much lower than that of the heterogeneity of residual variances. The high correspondence between population means across sites also reflects small deviations from perfect correlations across sites (see Fig. 2).

Table 2 Likelihood ratio tests for sources of the genotype by environment interaction (G × E) for height, diameter and survival at age 22, and processionary moth incidence at age 13

3.4 Performance of the local populations

The local populations (Chatba and Izarene, see Table S2) were among the best growing P. halepensis populations, with the Chatba population ranking eighth for diameter and height in Chatba and the Izarene population ranking ninth for both traits in Izarene at age 22. The superiority of the local P. halepensis population compared to the average of the remaining P. halepensis populations was statistically confirmed by specific contrasts in the Chatba site (F_{1, 376} = 61.8, p < 0.001 for height and F_{1, 372} = 71.6, p < 0.001 for diameter) but not in the Izarene test site (F_{1, 372} = 1.5, p = 0.219 for height and F_{1, 372} = 1.1, p = 0.303 for diameter). In the case of survival, although the local populations ranked above average in both sites, their survival was not statistically different than the mean survival of the other P. halepensis populations (F_{1, 376} = 1.4, p = 0.228 in Chatba and F_{1, 372} = 0.0, p = 0.956 in Izarene). Regarding the incidence of the processionary moth, no significant differences were observed between the local and non-local populations.

3.5 Correlation with biogeographical characteristics of origin

At age 22, height and diameter of P. halepensis populations at both sites were positively correlated with the longitude and negatively correlated with the altitude of their seed source origins (Fig. 5). Populations from the eastern Mediterranean Basin and those from lower altitudes tended to grow faster at both sites. Additionally, in the Izarene site, a significant latitudinal cline was also observed, with northern populations showing the highest growth (Fig. 5). In both test sites, growth at age 22 and annual temperature range were also significantly correlated. Populations from the coldest climates with high temperature oscillations at monthly and annual scales (i.e., continentality) grew slower than those from milder and warmer climates (Table S7). Neither geographic nor climatic variables were correlated with survival of P. halepensis populations (Table S7).

Fig. 5

Relationship between population height at the two test sites (superior panels a, b, and c for Chatba and inferior panels d, e, and f for Izarene) and the geographical characteristics at origin (altitude (a, d), latitude (b, e), and longitude (c, f)). Regression lines and Pearson correlation coefficients are shown for those significant (p < 0.05) within species relationships (continuous line for P. halepensis and dashed line for P. brutia populations)

Variation in processionary moth defoliation among Aleppo pine populations was positively correlated with the longitude and the annual mean temperature of the region of origin at both sites (Table S7). Populations from the eastern part of the species range and those originating from warmer climates tended to be more severely affected by the insect. In Izarene, where variation of population susceptibility to the insect was more pronounced, incidence of the processionary moth was also negatively correlated with altitude and with annual temperature range (Table S7).

Regarding P. brutia populations, highly significant correlation was observed in Izarene at age 22 between growth and altitude of the seed collection zone, with populations originating from higher altitudes growing faster (Fig. 5). Moreover, similarly to what was observed for P. halepensis, a positive relationship between growth at age 22 and longitude was also recorded, although this correlation was only significant for DBH in Chatba (Table S8). Finally, neither survival nor susceptibility to the processionary moth was correlated to geographic or climatic variables (Table S8).

4 Discussion

4.1 Among and within species variation

Growth and survival varied between the two pine species, with differences increasing as trees aged. The studied Brutia pine populations were, on average, the ones that better adapted to the Moroccan test site conditions. Brutia pine and Aleppo pine are known to vary in their ability to tolerate abiotic stress (e.g., cold and water stress) (Climent et al. 2009; Schiller 2000), and these differences commonly result in different survival and growth rates when planted together under stressful environmental conditions (Bariteau and Pommery 1992; Fisheri et al. 1986; Parlmberg 1975; Weinstein 1989). Relative performance of these two species, however, varies depending on the most critical environmental factor. Patterns observed in the present study are consistent with the species responses associated with moderate to severe summer water stress in the two test sites. P. halepensis is typically assumed to be more drought tolerant than P. brutia (Chambel et al. 2013; Schiller 2000), so the better average growth and survival of Brutia pine populations may therefore suggest that the experimental conditions were not, in fact, so harsh.

In agreement with the well-known strong population differentiation for neutral markers (Grivet et al. 2009; Kurt et al. 2011b), and for a number of adaptive traits (reviewed in Chambel et al. 2013), within species variation in growth and survival was also found to be extensive for both the Aleppo and the Brutia pine. Previous results suggested that the relative performance of the populations depended on the environmental conditions of the testing sites (Atzmon et al. 2004; Schiller and Atzmon 2009; Taïbi et al. 2015). In the case of Aleppo pine, results presented here are consistent with the variation observed in Israel for moderate to severe drought stress conditions (Schiller and Atzmon 2009; Weinstein 1989). As observed in those studies, Aleppo pine populations from Greece were the ones that performed best, showing outstanding growth and survival. This pattern was reflected in a geographical cline in which growth and survival were positively correlated with the longitude of population seed zones. These results, showing the superior performance of Eastern populations of Aleppo pine, are quite consistent across different studies (Climent et al. 2008; Voltas et al. 2015; Weinstein 1989) and have been explained in terms of the interaction between demographic processes associated with genetic bottlenecks resulting from severe reduction in population sizes during the glacial ages (Morgante et al. 1998) and the subsequent east to west colonization and adaptation to new environmental conditions (Gómez et al. 2005; Grivet et al. 2009). The result of this interaction is the specialization into divergent adaptive syndromes in which populations from the East show higher diversity levels and greater growth but are less water conservative, while populations from the West are less variable and more water conservative, resulting in lower growth potential (Climent et al. 2008; Voltas et al. 2008; Voltas et al. 2015). These results do not suggest, however, that Eastern populations of Aleppo pine are poorly adapted to drought stress. The highly heterozygous Greek Aleppo pine populations have been also shown to be well adapted to the extremely dry conditions of the Yatir desert in Israel (Schiller and Atzmon 2009). The physiological ability of those populations to tolerate severe water deficit (Calamassi et al. 2001; David-Schwartz et al. 2016) might be related to their higher diversity levels. This has led some authors to suggest heterozygosity as a selection criterion for selecting seed sources for reforestation in highly stressful environments (Schiller and Atzmon 2009).

It is worth mentioning that, although the success of the Aleppo pine populations in the stressful Moroccan conditions tended to diminish from East to West, some Moroccan populations (i.e., western range), including the local population of each test site, were among the best growers and survivors. This result, which was especially pronounced in the drier (although not necessary the more stressful) site, suggests that adaptive processes to the local environmental conditions have contributed to the differentiation among populations (Grivet et al. 2013; Grivet et al. 2011). It should be noted, however, that local adaptation has been tested here with just two populations in two test sites. Formal robust test of local adaptation would require a more complete design of reciprocal transplants (Kawecki and Ebert 2004), so results presented here should be managed with care. Nevertheless, given the relatively good performance of the local origins and in order to preserve local gene pools, using local seed sources is recommended as the most conservative strategy for reforestation, especially when timber production is not the main motivation of reforestation plans, as is the case in Morocco. Furthermore, the best growing Aleppo pine populations from the East of the Mediterranean basin were among the most affected by the processionary moth. Importing seeds from these origins might therefore lead to undesired outcomes related to their interaction with local herbivores, as has been found previously for other species (Di Matteo and Voltas 2016).

Local environmental conditions can, however, drastically shift due to climate change. Increase of temperature and aridity in the near future may alter the relative performance of pine populations (Petit et al. 2005; Taïbi et al. 2015). Nevertheless, the scarce practical relevance of the genotype by environment interaction (see discussion below), as well as the general patterns of variation in growth and survival observed across a wide range of environmental conditions (Climent et al. 2008; Schiller and Atzmon 2009; Voltas et al. 2008; Voltas et al. 2015; Weinstein 1989) suggest that the expected shifts in relative performance due to climate change will not be high. Given that premise and the slight superiority of the non-native seed origins, there is not yet enough evidence for implementing programs of assisted population migration, as it has been previously suggested (Taïbi et al. 2015; Taïbi et al. 2014).

Regarding Brutia pine, significant variation in growth and survival was also observed, but the level of population differentiation was lower, probably because of the reduced number of populations tested (eight populations). In agreement with previous findings (Isik et al. 2000; Weinstein 1989), the results showed that the Eastern Turkish populations performed well at both sites, while low-elevation populations from the edge of the species range in Greece performed worse. According to previous studies, altitude is an important driver of population differentiation (Dangasuk and Panetsos 2004; Kurt et al. 2011b; Weinstein 1989): low-elevation populations grew faster and survived better at dry conditions in Israel (Weinstein 1989), whereas middle to high elevation populations performed in average better across a large array of test sites in Turkey (Isik et al. 2000). Our results are in agreement with those reported by Weinstein (1989) and suggest that the more favorable environmental conditions prevailing in the middle to high elevations in Turkey are favoring fast-growing genotypes, whereas in peripheral areas characterized by harsher conditions, more conservative and slow-growing forms prevailed. However, the reduced number of P. brutia populations tested in the present study argues for caution in interpreting these results. P. brutia populations from middle to high elevations in Turkey not only showed good survival and growth rate but were also the most tolerant to the processionary moth, so they seem to be good alternatives for reforestation purposes in Morocco. Because hybridization between P. brutia and P. halepensis is unidirectional (from P. halepensis to P. brutia but not reciprocally) (Bucci et al. 1998), planting the exotic P. brutia will likely have no impact on the conservation of native P. halepensis genetic resources.

4.2 Population plasticity

In spite of the apparent environmental contrast between the two sites, our results consistently indicated low levels of plasticity, as has been previously reported for both Aleppo pine (e.g., Taïbi et al. 2015) and Brutia pine (e.g., Isik et al. 2000). Furthermore, at young ages, growth and survival was slightly higher in Chatba, the drier site, in which average annual precipitation during the course of the experiment was less than half that in Izarene. Summer drought was, however, extremely intense in both sites, with prolonged periods (up to 3–4 months) of almost no rain. Indeed, the two sites are classified within the same climatic region (Peel et al. 2007), and drought is assumed to be the main limiting factor in both sites. In Izarene, however, water deficit during the summer might have been exacerbated by the shallow and stony soils, the steeper slope, and the stronger competition with the accompanying vegetation. In other words, the stronger water deficit in Chatba seems to be compensated by its more favorable conditions for pine growth and survival in terms of physiographic, edaphic, and competition factors. This may explain the better growth and survival at this site at early ages. The good establishment and growth of Aleppo pine and Brutia pine in the two test sites in spite of the extreme summer droughts suggest, in any case, that these two species, or at least some of their populations, are well adapted to water stress.

Strong age to age correlation from age 7–10 onwards was also observed in both test sites, with population ranks becoming relatively stable after this age. Important population rank changes were recorded however when comparing population performance at younger ages. Many factors, including the post-planting stress that usually occurs in artificial reforestations, the incidence of environmental maternal effects that are likely stronger at earlier ages (Bischoff and Mueller-Schaerer 2010), or the ontogenetic shifts in the trade-offs between different life history traits as trees mature (Climent et al. 2008), may explain these patterns.

Results also indicated that among-population variation remained similar across sites, with little changes in population ranks across sites. Although the statistical analyses indicated that the G × E interaction was significant, an in-depth analysis of the sources of this interaction revealed a high correspondence between the relative population performances across the two test sites, with very high site to site correlation estimates. This result apparently contradicts previous findings in which a highly significant G × E interaction was observed for growth and survival and assumed to be of high practical relevance (Atzmon et al. 2004; Isik et al. 2000; Klein et al. 2013; Santos del Blanco et al. 2013; Taïbi et al. 2015; Taïbi et al. 2014). It should be noted, however, that despite the high site to site correlation estimates, some differential patterns can be observed when comparing the plasticity of populations coming from different climates at origin. Populations from colder origins with less summer drought (e.g., Brutia pine from Turkey and Aleppo pine from France, Italy, and Spain) performed better on the more mesic site (Izarene), while populations from warmer and more drought-prone origins (e.g., Aleppo pine from Greece, southern Tunisia and certain origins of Morocco) performed better in the drier one (Chatba). Although it is difficult to infer any trend based on results from only two sites, these findings suggest, again, that adaptive processes are likely to have contributed to the within species differentiation processes.

Most of the previous studies analyzing the G × E interaction in these pine species did not formally explore their causes and practical relevance (but see Isik et al. 2000; Taïbi et al. 2015). In our case, a likelihood analysis of the G × E interaction (Yang 2002) revealed that the heterogeneity of residual variances—reflecting differences in within-site environmental heterogeneity—was, by far, the main source of interaction, with residual variance being much higher in Izarene. This only affects the statistical efficiency in the two test sites, as the greater the residual variance the lower the ability for differentiating the performance among populations. However, heterogeneity of residual variance has no relevant implications regarding the relative performance of the tested provenances in the two environments (rank changes) and suggests low practical relevance of the G × E interaction. This result is likely to have arisen due to the harsher conditions in Izarene, where the poorer soils and the stronger inter- and intraspecific competition have probably increased the within-site environmental heterogeneity, resulting in larger small-scale variability that the statistical model was unable to account for. Accounting for microenvironmental heterogeneity by spatial adjustments will, thus, probably improve the efficiency of the statistical analyses (Zas 2006). Unfortunately, coordinates of each individual tree in the test sites were not available. The reduced practical relevance of the G × E found here should be thus confirmed by further studies exploring a wider range of environmental conditions with stronger effects on pine performance.

5 Conclusions

Analysis of growth and survival in the test sites revealed interesting patterns of adaptive variation. Overall, results indicate (i) large variation among and within species, with (ii) a high consistency in the relative performance across sites despite large differences in their statistical discrimination power resulting from contrasting residual variances, (iii) a strong geographical variation within species for growth, survival, and defoliation due to processionary moth in the analyzed populations, showing marked geographical clines, and (iv) low levels of plasticity across sites despite their environmental differences. The local populations were well adapted to their respective location but were outperformed by other Aleppo and Brutia populations from the Eastern Mediterranean region. However, the superiority of Eastern populations of Aleppo pine is not sufficiently high to justify their recommendation for afforestation in Morocco. In order to conserve local genetic resources, local seed sources and middle to high elevation P. brutia populations from Turkey are preferred. The low relevance of the GxE interaction suggests that this recommendation would hold in spite of future climate change.