Geographical and within-population variation of constitutive chemical defences in a Mediterranean oak ( Quercus ilex )

Aim of study: to assess whether constitutive levels of total phenols and condensed tannins acting as chemical defences in Quercus ilex vary across regions, populations and families. Area of study: acorns from six Spanish regions with natural presence of Q. ilex were collected for later sowing in the greenhouse of the University of Extremadura at Plasencia, western Spain. Material and methods: 1176 acorns from 12 populations (2 per region) were sown according to a nested design (6 regions x 2 popu-lations/region x 7 families/population). After five months, 588 emerged seedlings were harvested for phenol extraction from their leaves. Quantification of total phenols through the Folin-Ciocalteu method and condensed tannins by Porter's reagent method was carried out. Main results: total phenolics and condensed tannins correlated positively. Significant differences were found for the concentration of total phenolics and condensed tannins among regions and among families, but not among populations within regions. The lowest levels of constitutive defences were found in the northern Iberian Peninsula. Heritability values (total phenols 0.37±0.08 and condensed tannins 0.48±0.36) were lower than those obtained of total tannins in a previous study. Research highlights: constitutive levels of heritable chemical defences in holm oak significantly depend on their geographical origin and family.


Introduction
Chemical defences in plants have a high diversity due to their coevolution with natural enemies in different environments (Freeman & Beattie, 2008;Agrawal & Heil, 2012). This chemical diversity has been extensively studied in antagonistic interactions with herbivores. These organisms can vary the defensive response of plants according to their origin and impact (Thompson, 2005;Be-cerra et al., 2009). In addition, the pattern of constitutive chemical defences seems to be imprinted by phylogeographical signals (Eckert et al., 2008). In species with wide geographical distribution this diversity is even greater and its variability in response to different sources of stress affects plant survival (Solla et al., 2015;Ivetić et al., 2016;Ivetić & Devetaković, 2017;Martín et al., 2019). However, and despite its importance, little is known to what extent intraspecific variation in plant defenses is related to geographical factors (but see Solla et al., 2016;Moreira et al., 2018;López-Goldar et al., 2019). Geographical factors such as latitude or altitude have been shown to influence this defensive diversity across many plant taxa (Abdala-Roberts et al., 2016;Stevens et al., 2016;Moreira et al., 2017;Abdala-Roberts et al., 2018;Bogdziewicz et al., 2019;Moreira et al., 2020). Regarding latitude, a classic paradigm in ecology argues that at lower latitudes the herbivore pressure and also the production of plant defences increase (Rasmann & Agrawal, 2011;Pearse & Hipp, 2012;Moreira et al., 2014;Anstett et al., 2016). But this paradigm in ecology focused mainly on explaining clinal patterns among species and less on intraspecific variation (Hahn & Maron, 2016). In fact, recent studies of intraspecific variation of secondary metabolites (Stevens et al., 2016;Moreira et al., 2017;Moreira et al., 2018) have observed a positive relationship between latitude and defensive chemical levels of plants, suggesting that these patterns within species may not follow the same trends. A wide geographical distribution can also encompass greater differences in the edaphic variables (pH, moisture, temperature, porosity...) and climatic conditions. The potential of these environmental factors to modify basal phenolic levels has also been reported (Kraus et al., 2004;Salminen & Karonen, 2011;Moreira et al., 2014;Abdala-Roberts et al., 2016;Abdala-Roberts et al., 2018). The constitutive resistance levels correspond to defensive traits of plants that are relatively stable or whose variation depends on internal factors (Schultz, 1988). They are genetically controlled and show extensive variation between genotypes and provenances (Dicke & Hilker, 2003;Solla et al., 2016;Gallardo et al., 2019).
Recent studies have approached the interspecific diversity in the constitutive defences of several species of the genus Quercus , or the geographical component in the defensive chemical variability of species such as Quercus robur . But there are few known studies on the intraspecific variation of these secondary metabolites in holm oak (Quercus ilex; but see Solla et al., 2016) despite the declining status this species experiences in certain valuable habitats such as dehesas. The diversity in the response of plants to stress-causing agents allows the selection of those individuals, populations and/or provenances of greater tolerance. To achieve this, it is essential to study the chemical pattern of constitutive defences in the species and identify later what compounds could grant greater tolerance.
The presence and variations of chemical defences in holm oak according to the region of origin are still quite unknown (Solla et al., 2016). Advances in this field are needed to develop selection programmes of forest material for reproduction, because of their involvement in the adaptation and response to the environmental stress factors and, therefore, against future scenarios of changing environmental conditions (Isabel et al., 2020).
The holm oak is the most abundant forest species in the Iberian Peninsula and the dominant one in "dehesas", which constitute a High Natural and Cultural Value System (HNCV) thus considered by the European Union. Dehesas are socio-ecological systems exploited by humans throughout history and they constitute the economic base of the Western Iberian rural world (Escribano & Pulido, 1998;Campos, 2004;Pereira et al., 2004). Nevertheless, they are experiencing symptoms of decay in recent decades due to various factors (Brasier, 1992(Brasier, , 1996Sánchez et al., 2006;Corcobado et al., 2013) whose consequences are becoming increasingly intense. Although the importance of tree decay varies among regions, over time, the persistence of forests and dehesas could be at risk, with direct negative implications on resilience and associated biodiversity (Pulido et al., 2001;Plieninger et al., 2015;Duque-Lazo et al., 2018).
This study addresses the geographical variability in chemical defences of a widely distributed species in the Mediterranean basin, the holm oak. The main chemical defences in this species, as in other oaks, are based on phenolic compounds, especially high-molecular-weight tannins, both hydrolyzables (Ht, derivatives of gallic and ellagic acids) and condensed (Ct, also known as "proanthocyanidins"). These compounds reduce the digestibility of plant tissues in herbivores (Haslam, 2007), cause metabolic alterations in insects (Barbehenn et al., 2009) and inhibit the activity of root pathogens (Oliva et al., 1999;Kraus et al., 2003). Other studies have shown a potent antimicrobial activity (Scalbert, 1991) and suggest that pre-infection levels in plant tissues might condition the interaction between certain pathogens such as Phytophthora and their hosts (Conrad et al., 2017), although this resistance mechanism is still poorly known in the Mediterranean context. In this study, the constitutive chemical defences of holm oak in six regions of Spain are analysed to answer the following questions: i. Does the secondary metabolite concentration in holm oak show significant geographical and local variation?
ii. Does variation in defences depend on latitude or soil pH?
iii. To what extent this chemical variation is heritable?
iv. Which provenances are likely candidates for selection?

Plant material
Acorns from six Spanish regions with natural presence of holm oaks (located in National Parks and their surroundings) were collected between November and December 2015. Two natural populations were selected at random in each region (Fig. 1). The studied morphologies were always from Q. ilex subsp. rotundifolia in Cabañeros (C), Guadarrama (G), Monfragüe (M), Ordesa (O) and Sierra Nevada (S). In Picos de Europa (P), acorns were collected from a population subsp. rotundifolia and another one of Q. ilex subsp. ilex. Within each region, mature acorns without signs of infestation were collected. In order to evaluate their viability, acorns were subject to flotation before sowing. To determine whether differences are genetically-based, acorns from these regions were grown in a common garden at the University of Extremadura greenhouses in Plasencia (UTM Zone 29N X: 748862; Y: 4435709; 395 meters above sea level).

Experimental design
A nested type design with a total of 588 seedlings (6 regions x 2 populations/region x 7 families/population x 7 seedlings/family) was implemented at the greenhouse. Before sowing, acorns were placed in trays with sterilized sand for germination. Once the radicles emerged, the collected acorns were planted at the end of December in 28-cell plastic trays, 450 ml in volume, containing sand and peat (1:1, pH: 5.5). In order to maximize seedling availability, 2 acorns were planted per alveolus, selecting later the first emerged. The emergence took place in February and in general seedlings grew successfully.

Chemical analysis
Five months later, when the seedlings were at least 20 cm in height, three genotypes were harvested per each mother tree studied, removing the substrate carefully to avoid plant damage and subsequently processing the sample for freezing at -80ºC. Frozen seedlings were weighed and 5-6 leaves per plant were extracted for determination of phenolic content as described by Gallardo et al. (2019). These leaves were lyophilized using a Telstar LyoQuest lyophilizer (temperature -55ºC and 0.001-0.002 mbar pressure) and ground to a fine particle size. Once the milling was completed, tissue sampling was stored at -80ºC until analysis. The phenolic content of the plants was extracted from lyophilized material with 70% (v/v) aqueous methanol for 60 minutes in an ultrasonic bath at room temperature. The crude extracts were centrifuged at 10000 rpm for 5 minutes at 4ºC and the supernatant was collected and stored at -80ºC. The total phenolics content (Tp) of the extract was determined by the Folin-Ciocalteu method (Makkar, 2003). Crude extracts were mixed with 50 volumes of 10% Folin-Ciocalteu reagent (Merck KGaA, Darmstadt, Germany) and 40 volumes of 7.5% (w/v) sodium carbonate. In the control tube, the extract volume was replaced by methanol. The mixture was stirred gently and maintained in the dark and at room temperature for 45 minutes. After incubation, the absorbance was measured at 725 nm. Gallic acid (Sigma-Aldrich, Merck KGaA, Darmstadt, Germany) was used as standard and results were expressed as milligrams of gallic acid equivalents (G.A.E.) per grams of lyophilized sample.
The butanol-HCl assay (Porter et al., 1986) was used to quantify condensed tannins (Ct) using procyanidin B2 (Sigma) as a reference compound. Briefly, crude extracts were mixed with 100 volumes of n-butanol/acetone 1:1 (46% each) plus HCL (1.85%) and ferric ammonium sulphate (0.04%). In the control tube, the extract volume was replaced by methanol. Samples were heated at 70º. After 45 minutes of incubation, the samples were cooled and the absorbance at 550 nm was measured, with final results expressed as milligrams of procyanidin B equivalents (PB.E.) per grams of lyophilized sample.

Statistical analyses
The effects of region, population (nested in region) and family (nested in population) were analysed through a general linear mixed model (GLMM). Region, population and family were used as random factors and the Tp content (expressed in mg GAE/g sample) and Ct content (expressed in mg PBE/g sample) were used as dependent variables. Seedling weight was added to this model as a covariate in the Tp and Ct analysis to control for phenological effects on plant size due to differences in germination time. Data were analysed to check normality (by Kolmogorov-Smirnov test) and homoscedasticity (through Levene´s test). The relationships between both dependent variables as well as the influence of the germination time and seedling size (measured from their weight during the harvest of samples) on the constitutive defences were analysed with Pearson´s linear correlation. For the analysis of the effect of soil pH on the variables studied, the pH data were obtained from the Geochemistry Database of the Geological and Mining Institute of Spain (IGME, http://info.igme.es/Geoquimica/) and from the study of "Heavy Metals, Organic Matter and other Parameters of Agricultural Soils and Pastures in Spain" (Rodríguez-Martín et al., 2009). Thus, the regions were classified as acidic (less than 7 pH) and alkaline (greater than 7 pH).
Heritability of chemical traits (ĥ 2 ) was calculated following the procedure described by Solla et al. (2016) for the same species, as the additive genetic variance (VA) divided by the phenotypic variance (VP). The ĥ 2 estimate was corrected for a generalized selfing rate of 1% to 3% in native holm oak stands (Ortego et al., 2014), using an r coefficient of 0.27 because half-sibs have a quarter of their alleles in common (¼VA). The following equation [Eq. 1] was used: where  2 f(gen) is the variance component among families (genotypes),  2 r the among-regions variance and  2 e the error variance. Standard errors were obtained following Jayaraman (1999). The quantitative genetic differentiation (Q ST ) was estimated using the additive genetic variance in the denominator, as described by Gilbert & Whitlock (2015), but with the r coefficient from Ortego et al. (2014) for more precision [Eq. 2]: Qst = ( 2 r) (2 * 3.7 *  2 f(gen)) + ( 2 r)

Total phenolics and condensed tannins depend on region and family
Highly significant differences were observed in the concentration of Tp and Ct according to the region (p<0.001 for both compounds) and families within populations (p<0.05 for both compounds, Table 1). By contrast, there were no significant effects of the population. There was a significant relationship (r=0.589, p<0.001, r 2 =0.347) between the Tp obtained by the Folin-Ciocalteu method and the accumulation of Ct measured by the Porter method. The percentages of variance explained by region, population, and family for Tp were 21.37%, 1.61% and 30.72%, respectively. For Ct values were 17.21%, 0.004% and 32.00%, respectively (% variance calculated from results of GLMM).

Variation of defences with latitude and soil
Concentrations of Tp and Ct showed negative correlations with latitude, though they did not reach statistical significance due to low sample size (N=6 regions). Correlation values were r=-0.681, p=0.136 and adjusted r 2 =0.464 in Tp, and r=-0.236, p=0.653 and adjusted r 2 =0.055 in Ct. These chemical levels against the increasing soil pH with the same sample size (N=6 regions), also showed a non-significant negative correlation (r=-0.584, p=0.223 and adjusted r 2 =0.341 in Tp, and r=-0.467, p=0.351 and adjusted r 2 =0.212 in Ct).
Using populations within regions (N=12) instead of using the region as the unit for association with geographic variables, the negative trend is reinforced. The correlation with populations´ latitude was significant in Tp (r=-0.644, p=0.024 and adjusted r 2 =0.414), but not in Ct (r=-0.181, p=0.573 and adjusted r 2 =0.033). The correlation with soil pH of populations was not significant in Tp (r=-0.564, p=0.056 and adjusted r 2 =0.318) or Ct (r=-0.421, p=0.174 and adjusted r 2 =0.177).
To analyze the north-south geographical differences in defence concentration, we grouped all the popula-tions from the northern regions (P, O and G regions) and all the populations from the southern regions (M, C and S regions). Significantly higher averages were obtained in the defensive levels of the southern area versus the northern one for Tp (F 1,10 =5.72, p<0.05, Fig.2), but not for Ct (F 1,10 =1.38, p=0.26, Fig.3 To analyze differences in the defences according to soil type, two major groups were distinguished: populations from the acidic soil regions of G, M and C, and populations from the alkaline soil regions of P, O and S. Higher averages were obtained in the defensive levels on acidic soils than on alkaline ones, but they were weakly significant (F 1,10 =4.65, p=0.05) in Tp (Fig.2) and not significant (F 1,10 =2.15, p=0.17) in Ct (Fig.3

Heritability of total phenols and condensed tannins
Defences showed moderate narrow-sense heritability across populations (ĥ 2 ) both for Tp (0.37±0.08) and Ct (0.48±0.36) but the standard error for Ct was much larger than for Tp. The quantitative genetic differentiation (Q ST ) was 0.10 and 0.08, respectively (Table 2).

Discussion
This study shows significant patterns of variation among regions and families in constitutive and heritable secondary metabolites along a large portion of holm oak's range. Generally, higher levels of constitutive defences were observed at lower latitudes. This supports the hypothesis that at lower latitude, plants invest more in defences, also at the intraspecific level for the case of holm oak in the studied regions. Furthermore, as expected, Tp and Ct showed a positive correlation. Their narrow-sense heritability was moderate, lower than those obtained for the total tannin content in the same species (0.83; Solla et al., 2016) and similar to the broad-sense heritability of Ct in wild stands of the genus Populus (0.59; McKown et al., 2014). Despite the relatively low sample size, the results of this study offer a better understanding of the chemical defences of the holm oak for possible selection programmes.

Defences in holm oak show significant geographical and local variation
Significant differences were found for the concentration of total phenolics and condensed tannins among regions and among families, but not among populations within regions. The non-significant differences between the populations within each region indicate that our esti-mates characterize defence levels of the region of origin. Previous studies have analysed total phenols and condensed tannins in acorns of Q. ilex and Quercus suber (Cadahía et al., 1993) and in leaves of several oak species (Faeth, 1986;Dawra et al., 1988) without considering geographical variation (but see Solla et al., 2016). This work analyses for the first time the geographical variation of total phenols and condensed tannins in Q. ilex leaves. Nonetheless, a previous study explored the genetic variation in the concentration of total tannins in the same species but only in southern peninsular regions: Solla et al. (2016) showed that total tannin content did not differ significantly among regions, but it differed significantly among families. Our results confirm that, as compared with Solla et al. (2016), the geographical variation of chemical defences is specific to the phenolic group studied. In addition, the differences among regions increase when the sampling latitudinal range is extended. Moreira et al. (2018) analysed the variation of condensed tannins in another oak species, Q. robur, in 38 populations along an 18º latitudinal gradient, also obtaining significant differences. This disparity in the estimates of the different phenolic defence groups shows the need to continue studying chemical patterns in Quercus. Probably, a greater number of sampled regions of the species, including populations in the geographical limits of the natural distribution, as well as relict stands, could give a better characterization of the differences in the patterns of chemical defences. Indeed, the high genetic diversity of the holm oak, the postglacial development of several European Quercus refuges (Urli, 2013) and the capacity of plant secondary metabolites to evolve rapidly (Moore et al., 2014) are factors that make the intraspecific variation in holm oak an interesting and beneficial resource in extreme conditions.
Regarding the evaluated phenolic groups, we found a positive correlation between Tp and Ct, as expected. Therefore, both traits could be used as regional indicators of holm oak leaf defences. However, the proportion of Ct within Tp varied among regions and families. Considering the not only defensive role of these phenolic compounds that also perform other ecophysiological functions (Moore et al., 2014), it is expected that the Tp and Ct relationship may be under selection by biotic and abiotic factors and is fit as necessary to adapt to new environments. Finally, for the Picos de Europa region, our results for Tp and Ct were consistent with those found previously in a similar study performed using the same methodology (Gallardo et al., 2019).

Effects of latitude and soil pH on defences
Though our sample sizes were limited, we found that the lowest levels of constitutive defences correspond to regions in the northern Iberian Peninsula. Theory predicts that chemical defences should decrease as latitude increases (Rasmann & Agrawal, 2011;Abdala-Roberts et al., 2016). The results of this work are consistent with such prediction. However, recent studies have found a positive relationship between latitude and phenolic compounds in species such as Q. robur  or Pinus pinaster (López-Goldar et al., 2019). Our average values showed an increase in chemical defences with decreasing latitude. This trend could be linked to underlying factors associated with latitudinal variation such as climate and soil . But as previously seen, the results are specific to the phenolic group evaluated and the species studied. In addition, new studies would be necessary to know if these results remain exposed to biotic and abiotic stress agents. Following our results, Cabañeros and Monfragüe regions stood out especially above the rest for their high defensive levels and could be of interest for genotype selection.
On the other hand, increasing pH also showed a negative trend on phenolic defences. Higher averages were obtained on acidic soils (G, M and C) than on alkaline ones (P, O and S), but they were not significant due to low sample size. Soil pH influences the availability of nutrients for plants and may therefore be the cause of deficiency or toxicity (Benton, 2003). Tannins may induce allelopathic responses at the ecosystem level, including toxicity against pathogens (Kraus et al., 2003). One of the most widely spread and aggressive pathogens in recent decades is Phytopthora cinnamomi, the soil-borne oomycete that causes root rot disease and is identified as the primary causal invasive agent of Quercus decline in the Iberian Peninsula (Brasier, 1992(Brasier, , 1996Sánchez et al., 2006;Corcobado et al., 2014). This pathogen develops better in acidic than basic soils (Blaker & MacDonald, 1983;Dixon, 1984;Duvenhage & Kotzé, 1991). In fact, calcium mineral nutrition raises the soil pH and increases the tolerance of Q. ilex to Phytophthora root disease (Serrano et al., 2013). Assuming the important function that tannins exert on soil pathogens (Scalbert, 1991;Stong et al., 2013), plants in basic soil would not need to produce high levels of these compounds because the pathogen does not develop well at high pH. However, in acidic soils, this greater production of phenolic compounds could be the defensive response to a higher pressure by the pathogen.

Heritability of constitutive defences
The heritability of total phenols and condensed tannins in Q. ilex leaves was moderate. In the context of tree improvement determining the narrow-sense heritability is a useful tool to know how much of the observed variation of a certain quantitative character is due to the additive genetic component. In our study the narrow-sense heritability was moderate and likely to allow selection (Eriksson et al., 2013;Alcaide et al., 2019aAlcaide et al., , 2019b. Nevertheless, these results must be taken with caution due to the low number of families evaluated (Jensen & Barr, 1971;Visscher & Goddard, 2015). Heritability values were lower than those obtained in a previous work in the same species, but for the total tannin content (Tt) with fewer families sampled (0.37 in Tp and 0.48 in Ct versus 0.83 in Tt; Solla et al., 2016). However, this could enable shortterm evolutionary change in response to selection and makes the selection of material from regions with higher constitutive levels meaningful in future species improvement programmes. Genetic differentiation by quantitative traits (Q ST ) was 0.08 for Ct and 0.10 in Tp. Q ST -F ST estimates are one of the best empirical tests to detect signals of greater genetic divergence than that expected by neutral processes alone. Thus, by relating those chemical traits with significant Q ST -F ST to environmental factors, inferences could be made on what causes the patterns of genetic divergence. In general, there is a lack of knowledge regarding this geographical variation in chemical defences (but see, for Pinus pinaster, López-Goldar et al., 2019). Solla et al. (2016) studied the Q ST in total tannins, observing a higher genetic differentiation among Q. ilex regions (0.12). Despite our Q ST estimates lack F ST to compare with, the results suggest that genetic divergence of defensive traits may be due to adaptive variation to environmental differences among regions. Therefore, this study offers a first insight into the characterization of variations that Q. ilex secondary metabolites undergo in the Mediterranean area. Nevertheless, these results should be taken with caution owing to the relatively small size of sampled families. In addition, it would be very interesting to study the fixation index (F ST ), a measure of gene flow between populations, to find out if there is divergent selection in the species and thus, confirm whether the different families are favoured under different environmental conditions (Merilä & Crnokrak, 2001).

Applicability for selection programmes
The lowest phenol concentrations were found in the northern region, while the highest values were found in Cabañeros (highest Tp) and Monfragüe (highest Ct), in the southern region and on acidic soils. Likewise, the heritability of the evaluated characters was moderate but enough to show genetic control and opportunities for artificial selection (Eriksson et al., 2013;Alcaide et al., 2019a). The important variability observed among chemical defences patterns from regions may be a source of useful ecotypes in areas under high levels of biotic stress and, especially, in forests affected by the infection with pathogens of the genus Phytophthora, which are known to suffer from tannin inhibition (Stong et al., 2013). This suppressive ability of tannins is a key factor, since it is expected that plants with elevated tannin levels will develop better under unfavourable environmental conditions.
The basic unit to distinguish the forest reproductive materials (FRM) is the region of provenance, and in Spain these have been established for Q. ilex according to their distribution and patterns of variation (Alía et al., 2009). Southern areas with higher levels of constitutive defences seem to be candidates for possible selection of FRM against widespread pathogens such as P. cinnamo-mi. In particular, Cabañeros and Monfragüe consistently showed the highest levels of defensive compounds. However, in a widely distributed species such as the holm oak, other factors should be taken into account because individuals that could be good candidates for their high defensive chemical levels, may not adapt suitably to other types of soil, climatic conditions or phenological conditions, among others. Confirmatory studies are needed in order to establish a sound basis for future selection. Selection and breeding programmes that use this FRM often focus on productive traits, much less defensive traits, which often go unnoticed. However, in recent decades and due to the accelerated change in climatic conditions, the proliferation of pathogenic organisms has increased (Delgado-Baquerizo et al., 2020). It is urgent to deepen the search and selection of individuals with outstanding defensive traits, as is already done in the European chestnut tree (Castanea sativa) against the ink disease, also produced by Phytophthora cinnamomi. Thus, this work contributes to the improvement of breeding programmes, considering the defensive traits among the selection factors and deepening the knowledge of intraspecific variation of secondary metabolites in Quercus ilex.