Close Home
iForest - Biogeosciences and Forestry
vol. 11, pp. 651-659
Copyright © 2018 by the Italian Society of Silviculture and Forest Ecology
doi: 10.3832/ifor2572-011

Collection: INCOTW - Sassari, Italy (2017)
“International Congress on Cork Oak Trees and Woodlands”
Guest Editors: Piermaria Corona, Sandro Dettori

Research Articles

Delineation of seed collection zones based on environmental and genetic characteristics for Quercus suber L. in Sardinia, Italy

Giovanbattista De Dato (1)Corresponding author, Angela Teani (1), Claudia Mattioni (2), Maurizio Marchi (1), Maria Cristina Monteverdi (1), Fulvio Ducci (1)


The assessment of seed zones or regions of provenance (RoP) to preserve local adaptation of tree species is an important element for the correct management of forest reproductive materials (FRM), as well as for common practices in forestry and restoration activities ([72]). “Local-is-best” is the rationale underlying the delineation of tree seed zones ([15]), whose populations are expected to share similar environmental and climatic pressures; thus, planting populations adapted to local conditions allows both to preserve genotypes that are adapted to local conditions and to avoid the introduction of poorly adapted genotypes, which may cause loss of adaptive potential ([35]). In this perspective, local populations can ensure higher fitness in terms of survival, reproduction, productivity, disease resistance and abiotic resilience, than non-local materials ([35], [15]), unless climatic conditions change drastically.

In European countries, the identification of RoPs is regulated by the European Union (EU) directive 1999/105/CE and is a mandatory tool for the management of FRM. According to the directive, the RoP for a species or sub-species is the area or group of areas subject to sufficiently uniform ecological conditions in which stands or seed sources show similar phenotypic or genetic characters, taking into account altitudinal boundaries where appropriate. In such a context, two main approaches can be generally adopted, both based on clustering methods. When genetic data are available, a hierarchical clustering on sample populations may be performed ([39]) to aggregate stands sharing similar genetic traits. Although such a method represents the most accurate and reliable, genetic data are generally hard to be obtained. Moreover, an adequate sampling scheme would be highly expensive and difficult to be implemented for wide areas. For this reason, climate (i.e., ecological or environmental variability) is generally considered a good proxy for genetic variability. Actually, it is well known that climate mainly represents the most important driver for spatial distribution of forest species and reliable time series are one of the most relevant data for forest monitoring ([22]). In this case, the delineation of RoPs is commonly based on overall estimates of environmental homogeneity, mainly climate and/or soil characteristics ([49], [18]), and on a spatial clustering procedure. Main variables can be represented by long-term (30 years at least) averages of temperature and precipitation, and/or derived climatic indices which combine their values (e.g., aridity, continentality, etc. - [47]). However, the lack of species-specific genetic information in establishing RoPs could lead to underestimate the area of species range showing significant genetic differentiations, which may potentially cause increasing levels of inbreeding among planted trees ([41]). Therefore, integrating RoP maps with genetic data is an effective strategy to obtain a sound tool for the legislation of FRMs.

The climatic context of the Sardinia island has often been described in the literature as an example of a thoroughly Mediterranean island. Recently, Canu et al. ([6]) have classified the different environments across the island in 43 ecological groups, according to well know climatic belts. However, such groups are too many for a viable and cost-effective management of Forest Genetic Resources under the 1999/105/CE directive. For this reason, in 2016 the Regional Government of the Autonomous Region of Sardinia (Italy) developed a program aimed to screen seed sources and delimit RoPs with a specific target on cork oak (Quercus suber L.) populations, in order to protect and valorise available gene pools. Cork oak is a widely distributed species across the Mediterranean basin, both on the European side (Italy, France, Spain, Portugal) and North Africa (Tunisia, Morocco, Algeria). The species is widely cultivated in the Iberian Peninsula (Portugal, Western Spain, Catalonia and Balearic Islands) and France (Provence, Corsica, Aquitaine, Languedoc). In Italy, it is present along the Tyrrhenian coastline (Tuscany and Latium), on the main islands (Sicily and Sardinia), on the ionic Calabrian side and in Apulia. Cork oak is generally sporadic, and best growing in areas with cold, moist winters and hot summers. It is commonly found at 300-600 m but can occur up to 1000 m a.s.l. It prefers acid soils and avoids limestone and carbonated substrates. Sardinian cork oak stands cover about 140.000 ha ([32]) and have great socio-economic importance in the region due to its exploitation for cork extraction and processing, and for silvopastoral activities ([9], [61]). However, the risk of losing cork forests is high because of forest fires and diseases, overgrazing, land abandonment and exploitation, with subsequent biodiversity loss and deterioration of ecosystem services ([11]), exacerbated by the projected climate changes ([23]).

Studies of the genetic and geographical structure of Q. suber have been previously reported on different populations around the Mediterranean range using different molecular markers ([44]). First molecular characterization of Sardinian germplasm was accomplished by RAPD markers in several economically and agriculturally strategic stands ([5]). However, to establish a detailed RoP map for the species, an updated and refined assessment of the genetic structure of Sardinian populations is required. To this purpose, microsatellite markers represent a valid method for the assessment of genetic diversity due to their high levels of polymorphism, their codominant nature and their transferability.

The objectives of this study are: (1) to determine the ecological regions of provenance for Cork oak in Sardinia; (2) to investigate the genetic diversity among Cork oak populations at the regional scale; (3) to identify possible areas of interest to be selected as seed stands.

Materials and methods 

RoP delineation based on environmental variables

An ensemble of current climate data and soil characteristics was used to analyse the environmental variability across Sardinia. The ClimateEU software (freely available at ⇒ https:/­/­sites.­ualberta.­ca/­~ahamann/­data/­climateeu.­html) is a standalone software which can be used to downscale 4-km raster surfaces at higher spatial resolutions according to local elevation. Parameter-elevation Regressions on Independent Slopes Model (PRISM) monthly data (0.5° resolution) can be extracted and downscaled to scale-free point data ([14]), and the value of seasonal and annual climate variables for specific locations can be calculated based on latitude, longitude and elevation. The original meteorological data is derived from the CRU-TS 3.22 database ([26]), which includes more than 2000 meteorological stations across the whole Europe and 13 stations in the region of interest (within the grid box centred at N37.5-E7.5 and N42.5-E7.5). Such dataset has proven to be highly suitable for local studies in the Mediterranean region even when long time series are required, thanks to the fair representativeness of local trends ([46]). Variable surfaces are generated using lapse-rate adjustments according to latitude and longitude ([70]), based on the difference between the grid elevation and the elevation of the location of interest.

The software was used to generate high-resolution (250 m) raster maps of 36 biologically-relevant climatic variables (such as growing and chilling degree days, heating and cooling degree days, Hargreaves moisture deficit and reference evaporation, etc. - see Tab. S1 in Supplementary material) for the 1981-2010 normal period in Sardinia. A routine was created by calculating the downscaled values for each pixel. In addition, soil data (pedology and soil characteristics) were derived from the Italian eco-pedologic map ([12]).

All the obtained maps were then geo-referenced in the ETRS89/UTM32N (EPSG: 25832) reference system and used to perform a spatial principal component analysis (sPCA). All the data were centred and scaled separately. A standardisation procedure was applied as the variables were expressed in different scales (e.g., degrees Celsius for temperatures, millimetres for precipitation, pure numbers for indices). The aim of sPCA was to include as much environmental variability information as possible in further analysis (clustering), avoiding at the same time the collinearity between variables. As PCA compresses raw data variability into few uncorrelated (new) variables, PC components are not affected by collinearity between environmental factors. Only a subset of PC components was used and subjected to unsupervised classification ([27]) summing up to 99% of the explained variance.

Hierarchical clustering was applied to PC data in a GIS environment (“i.cluster” and “i.maxlik” functions of GRASS-GIS - [24]) using as weight the eigenvector of each ecological variable for each component obtained from PCA. Pairwise Euclidean distances between pixels were calculated for each variable. The maximum number of iterations was set to 50, while the convergence value (the points at which cluster means become stable) was set to 98% ([47]). The number of spatial clusters was set based on a dendrogram computed using the PC scores of 1000 points randomly distributed across Sardinia, regardless Quercus suber distribution.

Plant material and DNA extraction

A total of 277 individuals were randomly sampled at 10 sites (Tab. 1) for genetic analyses. Leaves were collected from individual trees and kept under silica gel for long-term conservation. Total genomic DNA was extracted by grinding 50-60 mg of dry material in liquid nitrogen and using the DNeasy® Plant mini kit (QIAGEN, Hilden, Germany). The extracted DNA was diluted to 5 ng µl-1 of working solution and quantified by spectrophotometry (Eppendorf Biophotometer® D30, Eppendorf AG, Hamburg, Germany).

Tab. 1 - Main environmental characteristics of the sampling sites. (Pop): population label; (Site): official forest name; (#ind): number of sampled individuals; (N): latitude (°); (E) longitude (°); (Alt): altitude (m a.s.l.); (MAT): mean annual temperature (°C); (MAP): mean annual precipitation (mm); (AHM): annual heat:moisture index (MAT+10)/(MAP/1000); (RoP): region of provenance each population belongs to. Further details of climate for each site are reported in Tab. S2 (Supplementary material).

Microsatellite analysis

Samples were analysed with a set of 5 polymorphic nuclear markers (nuSSR - Tab. 2), previously developed for other oak species (QpZAG9, QpZAG15, QpZAG36, QpZAG110 in Q. petrea (Matt) Liebl.; QrZAG7 in Q. robur L.). As a rule, SSRs are species-specific markers, which must be developed de novo for each species, mainly because they usually occur in non-coding regions of the genome, which are not highly conserved. However, Mendelian inheritance and transferability of nuSSRs from other oak species to cork oak had been previously demonstrated ([30], [62]). The initial set of markers also included MSQ4 ([17]), QrZAG20 ([37]) and QpZAG46 ([63]), though they have been subsequently discarded for poor/none amplification and/or unreliable scoring in most samples.

Tab. 2 - nuSSR markers selected for molecular analysis. (Na): number of alleles; (Ne): effective numbers of alleles; (He): expected heterozigosity; (Fis): fixation index; (ƒna): null allele frequency.

The five markers selected in this study showed a pattern of allelic richness and diversity similar to that obtained using three to seven nuSSRs markers in cork oak by Löpez-Aljorna et al. ([42]), Da Costa ([13]) and Zucca ([73]). The latter author analysed three populations from northern, central and southern Sardinia using seven nuSSRs, of which five are in common with our study. Comparing the diversity at each locus with our results (see below), we found a double amount of alleles (22 versus 11) at locus ZAG110, while at the other loci (ZAG36, ZAG9, ZAG15 and ZAG7) the number of alleles was similar, ranging between 4 and 10. Thus, we selected a posteriori the five markers best discriminating among individuals and populations. Moreover, the unbiased probability of identity (uPI - [51]) was computed to detect the chance that two individuals drawn at random from a sample will have the same genotype at multiple loci. The estimated uPI for the combination of the five selected markers over all populations was 1.3 × 10-4.

PCRs were performed in triplexes in Eppendorf Mastercycler® pro (Eppendorf AG, Hamburg, Germany). Each sample was amplified in a 12.5 µl total volume reaction according to QIAGEN Type-it® Microsatellite PCR kit (QIAGEN, Hilden, Germany). Cycling protocol was as follows: 95 °C for 5 min; 28 cycles at 95 °C for 30 sec, 57 °C for 90 sec, 72 °C for 30 sec; 60 °C for 30 min. 1 µl of each amplification product was added to a mixture of 20 µl HIDI formamide and 0.3 µl LIZ and denatured at 95 °C for 5 min. Samples were run on ABI Prism 3130® Avant DNA sequencer (Applied Biosystems, Foster City, CA, USA). Results were analysed by GeneMapper® Software (Life Technologies, Carlsbad, CA, USA), and allelic profiles were scored by automatic binning and visual checking.

Statistical analysis

Number of observed (Na) and effective (Ne) alleles, observed (Ho) and expected (He) heterozygosity ([50]) were calculated using GENALEX 6.5 software ([52], [53]). The inbreeding coefficient (Fis) for each population and each locus was obtained computing a hierarchical AMOVA using the software ARLEQUIN ver. Statistical significance was determined by a non-parametric approach using 1000 permutations ([21]). The assessment of allelic richness by measuring allele frequencies must take into account the variation in population sizes (Tab. 1). For this reason, allelic richness (A) and the number of private alleles (Pa) were computed using the rarefaction method with the HP-RARE software ([36]). Null allele frequency for each locus was estimated using the package FreeNA ([7]). Principal coordinates analysis (PCoA) was performed with GENALEX 6.5 using a pairwise population matrix of Nei’s unbiased genetic distance ([50]). Correlations among sites for measures of genetic, environmental and geographic distance were calculated by Mantel test and partial Mantel tests using the software PAST ver. 3.19 ([25]), with 9999 permutations. Partial Mantel tests were conducted to determine the strength of the correlation between two distance matrices aimed at assessing the dependence between two matrices of distances while controlling the effect of a third distance matrix.

Genetic structure was evaluated by a Bayesian approach performed using STRUCTURE software ver. 2.3.4 ([55]). Twenty iterations for each value of K between 1 and 10 were run with a burning period of 10.000 and 105 Markov chain Monte Carlo cycles. The most likely number of clusters was determined by the Evanno’s method ([20]), using the STRUCTURE HARVESTER software ([19]). Twenty runs were averaged by CLUMPP ver. 1.1.2 web service ([33]) and graphically represented using the software DISTRUCT ver. 1.1 ([58]).

Selection of priority areas and overlay with regions of provenance

Grids for genetic parameters were generated in DIVA-GIS software ver. 7.5 ([28]), based on a grid cell of 6 minutes and applying a circular neighbourhood with a diameter of one degree. Four different classes of alleles can be defined: (i) rare alleles with local distribution (private); (ii) frequent alleles observed with local distribution; (iii) rare alleles observed in a large area; (iv) frequent alleles commonly distributed. Priority for conservation should be given to populations that retain frequent alleles with local distribution, indicating the presence of genotypes adapted to specific environments. However, the Reserve Selection method implemented in DIVA-GIS and developed by Rebelo & Siegfried ([57]) considered all these classes for the identification of the priority conservation areas. A point-to-grid analysis was applied, selecting the “Complementary” option and giving equal weight to each allele. This method allows the selection of conservation priority areas not only based on their allelic richness but also on differences/complementarity in allelic composition, by identifying the minimum number of sites needed to include most of the genetic diversity estimated by nuclear markers. Then, the procedure ranks the geographic units by indicating the first chosen area as that with the highest allelic richness; each successive area selected best aggregates the intra-specific diversity already represented within the previously selected priority area.

Finally, the selected zones (seed stand candidates) were overlaid onto the RoPs as characterised above in the same GIS environment.


Regions of provenance

Four groups were selected by hierarchical clustering as best matching between environmental variability and spatial fragmentation. In the sPCA, 99% of the cumulated variance was accounted for by five components, which were selected as input for the GRASS-GIS modules. As a result, four spatial clusters were obtained (Fig. 1). The northern and northern-eastern part of the island was included in a unique RoP, named Northern Region, RoP1. The western and south-eastern parts were included in RoP2 and RoP3, respectively. The last group mainly included the higher elevation areas in the central part of the island (Central region, RoP4). RoP1 and RoP4 were colder than RoP2 and RoP3 concerning mean annual temperature (MAT); RoP1 was the driest region and RoP4 the wettest. The average climatic values and extremes of each RoP are reported in Tab. S2 (Supplementary material).

Fig. 1 - The regions of provenance (RoP) obtained in this study based on climatic variables.

Concerning the climate characteristics of sampling sites, the average values of the 36 climatic variables considered are reported in Tab. S3. Two main clusters of sampling sites were detected by multivariate analysis (PCA) based on the selected 36 climatic variables (Fig. 2). One cluster grouped the southern sites SUL, SET, MAN and MAR, including the northern LIT, which also has the lowest altitude; the other cluster grouped the northern sites MOL, LIM, TER, FIO and BAR mainly according to the directions of the minimum temperatures, soils, the beginning of the frost-free period and the variables related to seasonal precipitation in autumn and winter. Accordingly, a good correlation between geographic and climatic distance matrices was detected (r = 0.45, p < 0.05).

Fig. 2 - Principal component analysis (PCA) of the sampling sites according to climatic variables.

Genetic diversity

In Tab. 2 genetic diversity parameters for each locus are reported. The number of alleles per locus (Na) ranged from 2.7 for locus QpZAG9 to 13.9 for QpZAG110. The highest number of effective alleles (Ne) was observed at locus QpZAG110. Expected heterozygosity (He) was higher for QpZAG110 (0.88) than the other analysed loci. Fixation index (Fis) showed no significant departure from the Hardy-Weinberg equilibrium. Null allele frequency ranged between 0 and 0.043 (Tab. 2). Global Fst for all loci and pairwise Fst were very similar when estimated from raw allele frequencies or allele frequencies corrected for null alleles (data not shown). Thus, null alleles had very little influence on this analysis, and all further tests were performed using uncorrected allele frequencies.

Tab. 3 reports the genetic diversity parameters estimated at each sampling site. The number of alleles (Na) ranged from 37 (SUL) to 28 (MAR). Allelic richness (A) showed values from 6.3 (SUL and TER) to 5.6 (MAR). The number of private alleles ranged from 0.0 (LIT, FIO, BAR, MAR) to 0.3 (SUL). Expected heterozygosity (He) showed homogeneity among the populations, varying around 0.56. Fixation indices (Fis) showed no significant departure from zero in all populations and varied between -0.03 (MAN) to 0.07 (LIM).

Tab. 3 - Genetic diversity parameters estimated at each sampling site. (Na): total number of alleles; (A): allelic richness; (Pa): number of private alleles; (He): expected heterozigosity; (Fis): fixation index.

AMOVA analysis (Tab. 4) showed a low differentiation among populations (Fst = 0.009, p = 0.49), not significantly different from 0. The molecular variance among individuals within populations was Fis = 0.024. Fis was not significantly different from 0 (p = 0.06). Molecular variance was therefore mainly located at the individual level (Fit = 0.033, p < 0.05).

Tab. 4 - Results of the analysis of molecular variance (AMOVA).

Genetic structure

Nei’s unbiased genetic distance (uNei) values among populations varied between 0.001 and 0.059 (Tab. 5). Significant differences were observed among MAR and all the other populations. SET and SUL showed some significant distances with northern and central populations. Also the distance between SET and SUL, the southernmost populations, was significant.

Tab. 5 - Pairwise population matrix of Nei’s unbiased genetic distance. (*): p<0.05; (**): p<0.01; (***): p<0.001.

Mantel test between genetic and the geographic distances (IBD = isolation by distance) showed a significant correlation coefficient (r = 0.23, p < 0.05). Differently, the correlation between genetic and climatic distances was not significant (r = 0.18, p = 0.12). These results were confirmed by the partial Mantel test: when geographical factors were controlled a lack of isolation by climatic distance (r = 0.16, p = 0.18) was observed; whereas when climatic factors were controlled, we detected significant correlations between genetic differentiation and geographical distance (r = 0.32, p < 0.05).

The analysis of population structure revealed that all populations belonged to the same deme (data not shown).

Selection of priority populations and overlay with the regions of provenance

According to the Reserve Selection Analysis, four priority areas could be delimited (Tab. 6). Priority 1 was assigned to TER in the RoP1 Northern. Priority 2 was located in SUL in the RoP2 Western. An area of Priority 3 was located in SET, in RoP3. The LIM population was reported as an area of Priority 4, located in the RoP1. These sites contained almost the entire set of alleles detected in the island in this survey, missing only one allele (ZAG9_234) out of 52, present in only two sites (MAR and MAN) with a frequency around 5%. At the opposite, the four private alleles (measured only in one site with a frequency averaging 3%) were comprised in the selected sites (Fig. S2 in Supplementary material).

Tab. 6 - Priority areas identified by Reserve Selection Analysis. The number in parentheses refers to the corresponding region of provenance (RoP).


Environmental features, such as climate and soil characteristics, play important roles in determining tree population dynamics. Thus, the delimitation of RoPs is mainly based on the spatial homogeneity of environmental parameters. Their number and size vary greatly among countries, due to different approaches and criteria used for their delineation, based on ecological units or vegetation zones, phenotypic or genetic similarities (or their combination), the economic importance of the species ([40]). In some countries, including Italy, RoP delineation is done not at the national but at the regional level. Concerning the whole Sardinian environment, four regions were detected, which minimised the ecological variability within and maximised differences among RoPs. Although previous authors identified a higher number of environmental groups as adequate to explain and characterize the whole ecological variability of the island ([6]), such studies were always referred to more general spatial analyses (i.e., pure climatological classifications), aiming at including Sardinia in a broad European context and matching the local climate with the European ones. Indeed, a high number of groups is often unsuitable for cost-effective and straightforward management of FGR at the regional level. In fact, RoPs delineation should take into account the trade-off between statistical accuracy and practical complexity ([40]), in view of seeds and FRM transferability ([67]).

Regarding the species-specific RoPs for cork oaks, Varela ([68]) proposed seven RoPs for Portugal, based on environmental parameters, such as climate, lithology and management. Similarly to other studies for different tree species ([34], [3], [65]), these regions usually overlap with well-delimited areas easily identified for practical operations. In this study, a different approach strictly based on spatialised ecological parameters was used, also taking into account elevation. As a consequence, our map shows numerous small enclosed areas of RoP4, laying on the more elevated portion of the island, surrounded by areas belonging to a different RoP.

As expected, the populations sampled in this study revealed a low genetic differentiation and diversity. Although the seed dispersal is too low to be a determining factor in homogenising oak populations, the high level of pollen flow could influence variability within populations and contribute to decrease the inter-population variability ([42]), similarly to other species with an extensive gene flow and long distance anemophilous pollination ([45], [54], [2]). In addition, in the case of cork oak, anthropogenic pressure cannot be excluded, considering the long (for at least 2000 years) exploitation of the species with likely intentional acorn transport ([64]).

Molecular characterisation of Sardinian cork oak populations has already been assessed through RAPD markers ([5]), whose results revealed a significant differentiation among populations, differently from our study. Although the sampled populations do not coincide in these two studies, a different pattern could be expected, as RAPD markers exhibit a certain lack of reproducibility due to mismatch annealing and competition for amplification among various fragments ([16], [38]). More recently, Zucca ([73]) analysed different Mediterranean cork oak provenances in a common garden by nuSSR, confirming the low diversity of Sardinian populations compared to other European provenances (on average lower by 30%) and a concurrent low differentiation. Likewise, a low differentiation among populations and diversity within population was shown by AFLP markers among Portuguese cork oaks ([8]).

Based on the lack of well-defined different gene pools, seed collection should be sampled regardless of the location. However, through Reserve Selection Analysis, we determined the minimum number of areas to be sampled for including the highest number of alleles and private/rare alleles. According to Hoban & Schlarbaum ([29]), the goal of capturing all alleles with high probability may be unrealistic, while a reasonable goal of sampling many but not all alleles in the various categories (local rare, local common, etc.) may be preferred. This could help identify key sites for seed collection activities, because resources allocated for these activities are frequently scarce, limiting the possibility to conduct massive samplings in several areas. According to the Reserve Selection Analysis carried out, the areas of high priority identified in this study hold the majority of allelic richness and private alleles for cork oak in Sardinia: each selected population complemented at best the intra-specific diversity already represented within the previously selected priority area. As a consequence, the second-ranking priority area was not necessarily that with the second highest diversity (given that a large portion of this diversity had already been enclosed in the first priority area), but rather areas showing the allelic composition best complementing to that of the already prioritised area. Indeed, the four selected sites for cork oak in this study were not those with the highest rank of allelic richness, nor the ones with the highest rank of private alleles. The four priority areas were recognised in the northern sites TER and LIM and southern sites SET and SUL, mainly following the distribution of common alleles at the local scale and selecting almost the whole allelic pattern, without missing any private alleles. The sites covered three out of four RoPs, while no population was selected in the high elevation RoP4. Reserve Selection Analysis has widely been applied at broad (continental) scale, for selecting in-situ conservation areas for endangered species ([66], [69]), but also for defining seed collection zones ([31]). Our findings support previous suggestions about ex-situ conservation and seed collection, confirming that material should come from geographically or ecologically distant populations to maximise the collection of genetic diversity ([71]).

The selection of spatially distributed populations may also better represent the existing quantitative variation across the species range. In fact, it is known that neutral molecular diversity of populations does not always reflect their actual variation of quantitative traits. On the other hand, it has been reported that many plant populations are adapted to their local habitats, i.e., populations often include genetic variants showing a higher fitness in the local environment. Our study revealed a moderate isolation by distance (IBD, the correlation between genetic and geographic distances), although the correlation between genetic and climatic distance among sites was not significant. However, we could observe some significant genetic distances between several southern and northern populations. This suggests a likely influence of environmental conditions on the fitness of these cork oak populations, more related to latitude than altitude. The hypothesis is supported by Ramírez-Valiente et al. ([56]), who asserted that, in spite of the low genetic diversity and the weak among-population divergence using molecular markers, cork oak was able to develop local adaptations in response to natural selection. Additionally, in cork oak seedlings exposed to high irradiance and drought under controlled environment, provenances from dry places exhibited higher tolerance to these environmental stressors, despite no relevant change in morphology ([1], [59]). Other studies analyzed the variability of candidate genes related to different environmental conditions, suggesting the existence of local adaptations ([48], [43]). Epigenetic regulation has also been reported as an important factor for adaptation of cork oak to specific environmental stresses ([10]).

Therefore, in order to improve this first assessment on the genetic resources of cork oak in Sardinia, common gardens and transplanting along the island in different RoPs will be necessary ([4]): these would provide tools to characterise provenances according to their genetic differentiation and plastic responses ([60]), especially to drought, the main threat in the Mediterranean region.


This study integrated the information on the population genetic structure of cork oak by nuSSR markers with environmental characteristics of Sardinia in order to delineate regions of provenance for the species. Four RoPs were selected as best matching between environmental variability and spatial fragmentation. Almost the entire genetic diversity was observed at the individual level. According to the Reserve Selection Analysis, four priority areas were identified which are proposed as seed collection sites for cork oak in Sardinia. These findings can be useful for future management of cork oak genetic resources and future plantations in the island. However, further efforts aimed to investigate the functional response of different provenances are necessary by common gardens and provenance trials along the island.


This study was funded by the Project LR 7/2007 TENDER “Multifunzionalità delle foreste a quercia da sughero” - Regione Autonoma della Sardegna. The authors want to thank the Agenzia Fo.Re.S.T.A.S., in particular, Gen. Manager dr. Casula, Dr. Sara Maltoni, forest officers, workers. Molecular analyses by sequencer were performed at the National Research Council - Institute of Agro-environmental and Forest Biology (CNR-IBAF) in Porano (TR, Italy), assisted by Marcello Cherubini.


Aranda I, Castro L, Pardos M, Gil L, Pardos JA (2005). Effects of the interaction between drought and shade on water relations, gas exchange and morphological traits in cork oak (Quercus suber L.) seedlings. Forest Ecology and Management 210: 117-129.
::CrossRef::Google Scholar::
Beatty GE, Brown JA, Cassidy EM, Finlay CMV, McKendrick L, Montgomery WI, Reid N, Tosh DG, Provan J (2015). Lack of genetic structure and evidence for long-distance dispersal in ash (Fraxinus excelsior) populations under threat from an emergent fungal pathogen: implications for restorative planting. Tree Genetics and Genomes 11 (3): 148.
::CrossRef::Google Scholar::
BFW (2017). Austrian regions of provenance. Web site.
::Online::Google Scholar::
Boshier D, Broadhurst L, Cornelius J, Gallo L, Koskela J, Loo J, Petrokofsky G, St Clair B (2015). Is local best? Examining the evidence for local adaptation in trees and its scale. Environmental Evidence 4: 311.
::CrossRef::Google Scholar::
Bullitta S, Dettori S, Manchinu M, Filigheddu MR, Piluzza G (2011). Characterization of Sardinian cork oak (Quercus suber L.) genetic resources for economically important traits. Genetic Resources and Crop Evolution 58: 1007-1020.
::CrossRef::Google Scholar::
Canu S, Rosati L, Fiori M, Motroni A, Filigheddu R, Farris E (2015). Bioclimate map of Sardinia (Italy). Journal of Maps 11: 711-718.
::CrossRef::Google Scholar::
Chapuis M-P, Estoup A (2007). Microsatellite null alleles and estimation of population differentiation. Molecular Biology and Evolution 24: 621-631.
::CrossRef::Google Scholar::
Coelho AC, Lima MB, Neves D, Cravador A (2006). Genetic diversity of two evergreen oaks (Quercus suber L. and Quercus ilex subsp. rotundifolia Lam.) in Portugal using AFLP markers. Silvae Genetica 55: 105-118.
::CrossRef::Google Scholar::
Corona P, Dettori S, Filigheddu MR, Maetzke F, Scotti R (2005). Site quality evaluation by classification tree: an application to cork quality in Sardinia. European Journal of Forest Research 124: 37-46.
::CrossRef::Google Scholar::
Correia B, Valledor L, Meijón M, Rodriguez JL, Dias MC, Santos C, Cañal MJ, Rodriguez R, Pinto G (2013). Is the interplay between epigenetic markers related to the acclimation of cork oak plants to high temperatures? PLoS One 8: e53543.
::CrossRef::Google Scholar::
Costa A, Pereira H, Madeira M (2010). Analysis of spatial patterns of oak decline in cork oak woodlands in Mediterranean conditions. Annals of Forest Science 67: 204.
::CrossRef::Google Scholar::
Costantini E, Barbetti R, Fantappiè M, Abate G, Lorenzetti R, Napoli R, Marchetti A, Rivieccio R (2014). The soil map of Italy: a hierarchy of geodatabases, from soil regions to sub-systems. In: “GlobalSoilMap: Basis of the Global Spatial Soil Information System” (Arrouays D, McKenzie N, Hempel J, Richer de Forges A, McBratney AB eds). CRC Press Taylor and Francis Group, London, UK, pp. 109-112.
::Online::Google Scholar::
Da Costa JSPN (2011). Differentiation and genetic variability in cork oak populations (Quercus suber L.). PhD thesis, Faculdade de Ciencias Departamento de Biologia Animal, University of Lisbon, Lisbon, Portugal, pp. 120.
::Online::Google Scholar::
Daly C, Halbleib M, Smith JI, Gibson WP, Doggett MK, Taylor GH, Curtis J, Pasteris PP (2008). Physiographically sensitive mapping of climatological temperature and precipitation across the conterminous United States. International Journal of Climatology 28: 2031-2064.
::CrossRef::Google Scholar::
De Kort H, Mergeay J, Vander Mijnsbrugge K, Decocq G, Maccherini S, Kehlet Bruun HH, Honnay O, Vandepitte K (2014). An evaluation of seed zone delineation using phenotypic and population genomic data on black alder. Journal of Applied Ecology 51: 1218-1227.
::CrossRef::Google Scholar::
Demeke T, Sasikumar B, Hucl P, Chibbar RN (1997). Random Amplified Polymorphic DNA (RAPD) in cereal improvement. Maydica 42: 133-142.
::Google Scholar::
Dow BD, Ashley MV, Howe HF (1995). Characterisation of highly variable (GA/CT)n microsatellites in the bur oak, Quercus macrocarpa. Theoretical and Applied Genetics 91: 137-141.
::CrossRef::Google Scholar::
Ducci F (2015). Genetic resources and forestry in the Mediterranean region in relation to global change. Annals of Silvicultural Research 39: 70-93.
::CrossRef::Google Scholar::
Earl DA, Von Holdt BM (2012). STRUCTURE HARVESTER: a website and program for visualizing STRUCTURE output and implementing the Evanno method. Conservation of Genetic Resources 4: 359-361.
::CrossRef::Google Scholar::
Evanno G, Regnaut S, Goudet J (2005). Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study. Molecular Ecology 14: 2611-2620.
::CrossRef::Google Scholar::
Excoffier L, Lischer HEL (2010). Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Molecular Ecology Resources 10: 564-567.
::CrossRef::Google Scholar::
Ferrara C, Marchi M, Fares S, Salvati L (2017). Sampling strategies for high quality time-series of climatic variables in forest resource assessment. iForest - Biogeosciences and Forestry 10: 739-745.
::CrossRef::Google Scholar::
Giorgi F, Lionello P (2008). Climate change projections for the Mediterranean region. Global and Planetary Change 63: 90-104.
::CrossRef::Google Scholar::
GRASS Development Team (2017). Geographic Resources Analysis Support System (GRASS) Software, version 7.2. Open Source Geospatial Foundation, Web site.
::Online::Google Scholar::
Hammer O, Harper DAT, Ryan PD (2001). PAST: Paleontological statistics software package for education and data analysis. Palaeontologia Electronica 4: 1-9.
::Google Scholar::
Harris I, Jones PD, Osborn TJ, Lister DH (2014). Updated high-resolution grids of monthly climatic observations - the CRU TS3.10 Dataset. International Journal of Climatology 34: 623-642.
::CrossRef::Google Scholar::
Hastie T, Tibshirani R, Friedman J (2001). 14 Unsupervised learning. In: “The Elements of Statistical Learning: Data Mining, Inference, and Prediction” (Hastie T, Tibshirani R, Friedman J eds). Springer Open Ltd, New York, USA, pp. 272-280.
::CrossRef::Google Scholar::
Hijmans RJ, Guarino L, Cruz M, Rojas E (2001). Computer tools for spatial analysis of plant genetic resources data: 1. DIVA-GIS. Plant Genetic Resources Newsletter 127: 15-19.
::Online::Google Scholar::
Hoban S, Schlarbaum S (2014). Optimal sampling of seeds from plant populations for ex-situ conservation of genetic biodiversity, considering realistic population structure. Biological Conservation 177: 90-99.
::CrossRef::Google Scholar::
Hornero J, Gallego FJ, Martinez I, Toribio M (2001). Testing the conservation of Quercus spp. microsatellites in the cork oak, Q. suber L. Silvae Genetica 50: 162-167.
::Online::Google Scholar::
Hufford KM, Veneklaas EJ, Lambers H, Krauss SL (2016). Genetic delineation of local provenance defines seed collection zones along a climate gradient. AoB PLANTS 8: plv149.
::CrossRef::Google Scholar::
INFC (2005). Secondo Inventario Nazionale delle Foreste e dei Serbatoi di Carbonio [2nd National Inventory of Forests and Carbon Sinks]. Web site. [in Italian]
::Online::Google Scholar::
Jakobsson M, Rosenberg NA (2007). CLUMPP: a cluster matching and permutation program for dealing with label switching and multimodality in analysis of population structure. Bioinformatics 23: 1801-1806.
::CrossRef::Google Scholar::
Jansson G, Danusevicius D, Grotehusmann H, Kowalczyk J, Skrppa T, Krajmerova D, Wolf H (2013). Norway Spruce (Picea abies (L.) H. Karst.). In: “Forest Tree Breeding in Europe. Current State-of-the-Art and Perspectives” (Paques LE ed). Springer Netherlands, pp. 123-176.
::CrossRef::Google Scholar::
Jones TA (2013). When local isn’t best. Evolutionary Applications 6: 1109-1118.
::CrossRef::Google Scholar::
Kalinowski ST (2005). HP-RARE 1.0: a computer program for performing rarefaction on measures of allelic richness. Molecular Ecology Notes 5: 187-189.
::CrossRef::Google Scholar::
Kampfer S, Lexer C, Glössl J, Steinkellner H (1998). Characterization of (GA)n microsatellite loci from Quercus robur. Hereditas 129: 183-186.
::CrossRef::Google Scholar::
Karp A, Edwards K, Bruford M, Vosman B, Morgante M, Seberg O, Kremer A, Boursot P, Arctander P, Tautz D, Hewitt G (1997). Molecular technologies for biodiversity evaluation: opportunities and challenges. Nature Biotechnology 15: 625-628.
::CrossRef::Google Scholar::
Kleinschmit JRG, Kownatzki D, Gregorius H-R (2004). Adaptational characteristics of autochthonous populations consequences for provenance delineation. Forest Ecology and Management 197: 213-224.
::CrossRef::Google Scholar::
Konnert M, Fady B, Gömöry D, O’Hara S, Wolter F, Ducci F, Koskela J, Bozzano M, Maaten T, Kowalczyk J (2015). Use and transfer of forest reproductive material in Europe in the context of climate change. European Forest Genetic Resources Programme (EUFORGEN), Bioversity International, Rome, Italy, pp. 75.
::Online::Google Scholar::
Kremer A, Ronce O, Robledo-Arnuncio JJ, Guillaume F, Bohrer G, Nathan R, Bridle JR, Gomulkiewicz R, Klein EK, Ritland K, Kuparinen A, Gerber S, Schueler S (2012). Long-distance gene flow and adaptation of forest trees to rapid climate change. Ecology Letters 15: 378-392.
::CrossRef::Google Scholar::
Löpez-Aljorna A, Bueno MA, Aguinagalde I, Martín JP (2007). Fingerprinting and genetic variability in cork oak (Quercus suber L.) elite trees using ISSR and SSR markers. Annals of Forestry Science 64: 773-779.
::CrossRef::Google Scholar::
Magalhães AP, Verde N, Reis F, Martins I, Costa D, Lino-Neto T, Castro PH, Tavares RM, Azevedo H (2015). RNA-Seq and gene network analysis uncover activation of an ABA-dependent signalosome during the cork oak root response to drought. Frontiers in Plant Science 6: 1195.
::CrossRef::Google Scholar::
Magri D, Fineschi S, Bellarosa R, Buonamici A, Sebastiani F, Schirone B, Simeone MC, Vendramin GG (2007). The distribution of Quercus suber chloroplast haplotypes matches the palaeogeographical history of the western Mediterranean. Molecular Ecology 16: 5259-5266.
::CrossRef::Google Scholar::
Maliouchenko O, Palmé AE, Buonamici A, Vendramin GG, Lascoux M (2007). Comparative phylogeography and population structure of European Betula species, with particular focus on B. pendula and B. pubescens. Journal of Biogeography 34: 1601-1610.
::CrossRef::Google Scholar::
Marchi M, Castaldi C, Merlini P, Nocentini S, Ducci F (2015). Stand structure and influence of climate on growth trends of a Marginal forest population of Pinus nigra spp. nigra. Annals of Silvicultural Research 39: 100-110.
::CrossRef::Google Scholar::
Marchi M, Chiavetta U, Castaldi C, Di Silvestro D, Contu F, Ducci F (2016). Regions of provenance for reproductive materials of the three main forest species of Abruzzi. Journal of Maps 12: 94-97.
::CrossRef::Google Scholar::
Modesto IS, Miguel C, Pina-Martins F, Glushkova M, Veloso M, Paulo OS, Batista D (2014). Identifying signatures of natural selection in cork oak (Quercus suber L.) genes through SNP analysis. Tree Genetics and Genomes 10: 1645-1660.
::CrossRef::Google Scholar::
Nanson A (2001). The new OECD scheme for the certification of forest reproductive materials. Silvae Genetica 50: 181-186.
::Online::Google Scholar::
Nei M (1973). Analysis of gene diversity in subdivided populations. Proceedings of the National Academy of Sciences USA 70: 3321-3323.
::CrossRef::Google Scholar::
Paetkau D, Waits LP, Clarkson PL, Craighead L, Vyse E, Ward R, Strobeck C (1998). Variation in genetic diversity across the range of North American Brown Bears. Conservation Biology 12: 418-429.
::CrossRef::Google Scholar::
Peakall R, Smouse PE (2006). GENALEX 6: genetic analysis in Excel. Population genetic software for teaching and research. Molecular Ecology 6: 288-295.
::CrossRef::Google Scholar::
Peakall R, Smouse PE (2012). GenAlEx 6.5: genetic analysis in Excel. Population genetic software for teaching and research-an update. Bioinformatics 28: 2537-2539.
::CrossRef::Google Scholar::
Piotti A, Leonardi S, Buiteveld J, Geburek T, Gerber S, Kramer K, Vettori C, Vendramin GG (2012). Comparison of pollen gene flow among four European beech (Fagus sylvatica L.) populations characterized by different management regimes. Heredity 108: 322-331.
::CrossRef::Google Scholar::
Pritchard JK, Stephens M, Donnelly P (2000). Inference of population structure using multilocus genotype data. Genetics 155: 945-959.
::Online::Google Scholar::
Ramírez-Valiente JA, Lorenzo Z, Soto A, Valladares F, Gil L, Aranda I (2009). Elucidating the role of genetic drift and natural selection in cork oak differentiation regarding drought tolerance. Molecular Ecology 18: 3803-3815.
::CrossRef::Google Scholar::
Rebelo AG, Siegfried WR (1992). Where should nature reserves be located in the Cape floristic region, South Africa? Models for the spatial configuration of a reserve network aimed at maximizing the protection of floral diversity. Conservation Biology 6: 243-252.
::CrossRef::Google Scholar::
Rosenberg NA (2004). DISTRUCT: a program for the graphical display of population structure. Molecular Ecology Notes 4: 137-138.
::CrossRef::Google Scholar::
Rzigui T, Cherif J, Zorrig W, Khaldi A, Nasr Z (2017). Adjustment of photosynthetic carbon assimilation to higher growth irradiance in three-year-old seedlings of two Tunisian provenances of Cork Oak (Quercus suber L.). iForest - Biogeosciences and Forestry 10: 618-624.
::CrossRef::Google Scholar::
Sáenz-Romero C, Lamy JB, Ducousso A, Musch B, Ehrenmann F, Delzon S, Cavers S, Chalupka W, Dagdas S, Hansen JK, Lee SJ, Liesebach M, Rau HM, Psomas A, Schneck V, Steiner W, Zimmermann NE, Kremer A (2017). Adaptive and plastic responses of Quercus petraea populations to climate across Europe. Global Change Biology 23: 2831-2847.
::CrossRef::Google Scholar::
Sedda L, Delogu G, Dettori S (2011). Forty-four years of land use changes in a Sardinian cork oak agro-silvopastoral system: a qualitative analysis. The Open Forest Science Journal 4: 57-66.
::CrossRef::Google Scholar::
Soto A, Lorenzo Z, Gil L (2003). Nuclear microsatellite markers for the identification of Quercus ilex L. and Q. suber L. hybrids. Silvae Genetica 52: 63-66.
::Online::Google Scholar::
Steinkellner H, Fluch S, Turetschek E, Lexer C, Streiff R, Kremer A, Burg K, Glössl J (1997). Identification and characterization of (GA/CT)n - microsatellite loci from Quercus petraea. Plant Molecular Biology 33: 1093-1096.
::CrossRef::Google Scholar::
Toumi L, Lumaret R (1998). Allozyme variation in cork oak (Quercus suber L.): the role of phylogeography and genetic introgression by other Mediterranean oak species and human activities. Theoretical and Applied Genetics 97: 647-656.
::CrossRef::Google Scholar::
UKFC (2017). Regions of provenance and native seed zones. Forestry Commission, UK, Web site.
::Online::Google Scholar::
Van Zonneveld M, Scheldeman X, Escribano P, Viruel MA, Van Damme P, Garcia W, Tapia C, Romero J, Sigueñas M, Hormaza J (2012). Mapping genetic diversity of cherimoya (Annona cherimola Mill.): application of spatial analysis for conservation and use of plant genetic resources. PLoS One 7: e29845.
::CrossRef::Google Scholar::
Vander Mijnsbrugge K, Bischoff A, Smith B (2010). A question of origin: where and how to collect seed for ecological restoration. Basic and Applied Ecology 11: 300-311.
::CrossRef::Google Scholar::
Varela MC (1997). Regions of provenance for Cork oak in Portugal. In: “Quercus suber Network. Report of the 3rd and 4th meetings” (Turok J, Varela MC, Hansen C eds). Sassari (Italy) 9-12 June 1996 and Almoraima (Spain) 20-22 Feb 1997. International Plant Genetic Resources Institute, Rome, Italy, pp. 37-42.
::Google Scholar::
Vinceti B, Loo J, Gaisberger H, Van Zonneveld MJ, Schueler S, Konrad H, Kadu CA, Geburek T (2013). Conservation priorities for Prunus africana defined with the aid of spatial analysis of genetic data and climatic variables. PLoS One 8: e59987.
::CrossRef::Google Scholar::
Wang T, Hamann A, Spittlehouse DL, Murdock TQ (2012). ClimateWNA-high-resolution spatial climate data for western North America. Journal of Applied Meteorology and Climatology 51: 16-29.
::CrossRef::Google Scholar::
Whitlock R, Hipperson H, Thompson DBA, Butlin RK, Burke T (2016). Consequences of in-situ strategies for the conservation of plant genetic diversity. Biological Conservation 203: 134-142.
::CrossRef::Google Scholar::
Williams MI, Dumroese RK (2013). Preparing for climate change: forestry and assisted migration. Journal of Forestry 111: 287-297.
::CrossRef::Google Scholar::
Zucca GM (2011). Molecular and phenoptypic characterization of Quercus suber L. and Pinus uncinata R. populations in the Mediterranean basin. PhD thesis, Scienze dei Sistemi Agrari e Forestali e delle Produzioni Alimentari, University of Sassari, Sassari, Italy, pp. 124.
::Google Scholar::


Paper Contents

Paper Sections

Paper Figures

Paper Tables

Supplementary Material



De Dato G, Teani A, Mattioni C, Marchi M, Monteverdi MC, Ducci F (2018).
Delineation of seed collection zones based on environmental and genetic characteristics for Quercus suber L. in Sardinia, Italy
iForest - Biogeosciences and Forestry 11: 651-659. - doi: 10.3832/ifor2572-011
First Previous Next Last
© iForest

Download Reference

Paper ID# ifor2572-011
Title Delineation of seed collection zones based on environmental and genetic characteristics for Quercus suber L. in Sardinia, Italy
Authors De Dato G, Teani A, Mattioni C, Marchi M, Monteverdi MC, Ducci F
Close Download