Phylogeography of SW Mediterranean firs: Different European origins for the North African Abies species

Molecular Phylogenetics and Evolution 79 (2014) 42–53

Contents lists available at ScienceDirect

Molecular Phylogenetics and Evolution journal homepage: www.elsevier.com/locate/ympev

Phylogeography of SW Mediterranean firs: Different European origins for the North African Abies species Jose M. Sánchez-Robles a,⇑, Francisco Balao a,b, Anass Terrab a, Juan L. García-Castaño a, María A. Ortiz a, Errol Vela c, Salvador Talavera a a b c

Departamento de Biología Vegetal y Ecología, Facultad de Biología, Universidad de Sevilla, Apdo. 1095, 41080 Sevilla, Spain Department of Systematic and Evolutionary Botany, Faculty Centre of Biodiversity, University of Vienna, Rennweg 14, A-1030 Vienna, Austria UMR AMAP (botAnique et bio-inforMatique de l’Architecture des Plantes), Université Montpellier-II, CIRAD TA/A51, bat. PS2, 34398 Montpellier cedex 5, France

a r t i c l e

i n f o

Article history: Received 28 November 2013 Revised 15 April 2014 Accepted 4 June 2014 Available online 24 June 2014 Keywords: Abies sect. Piceaster AFLP cpSSR Gymnosperms Mediterranean Basin Vicariance

a b s t r a c t The current distribution of Western Mediterranean Abies species is a result of complex geodynamic processes and climatic oscillations that occurred in the past. Abies sect. Piceaster offers a good study model to explore how geo-climatic oscillations might have influenced its expansion and diversification on both sides of the W Mediterranean basin. We investigated the genetic variation within and among nine populations from five Abies species by molecular markers with high and low mutation rates and contrasting inheritance (AFLP and cpSSR). Analyses revealed the opening of the Strait of Gibraltar as an effective barrier against gene flow between the Southern Iberian (A. pinsapo) and North African (A. marocana and A. tazaotana) firs. The A. pinsapo populations in Spain and likewise those of the A. marocana – A. tazaotana population complex were not differentiated, and no evidence was found to distinguish A. tazaotana at the species level. Diversification of Abies across North Africa could occur by way of at least two vicariant events from Europe, in the west, giving rise to the A. marocana – A. tazaotana complex, and in the east, giving A. numidica. Secondary contacts among species from Abies sect. Piceaster (A. pinsapo and A. numidica), and with A. alba (Abies sect. Abies) are also indicated. However, there is a closer relationship between the Algerian fir (A. numidica) and the North Mediterranean widespread A. alba, than with the Moroccan firs (A. marocana and A. tazaotana) or the Southern Iberian (A. pinsapo). We also discuss the distribution range of these taxa in its paleogeological and paleoclimatic context, and propose that part of the modern geography of the South-Western Mediterranean firs might be traced back to the Tertiary. Ó 2014 Elsevier Inc. All rights reserved.

1. Introduction The Mediterranean region has been recognized as one of the world’s hotspots where exceptional levels of biodiversity occur (Myers et al., 2000). Representing only 1.6% of Earth’s dry land, this region harbors more than 25,000 known vascular plant species. This striking diversity has been largely shaped by geo-climatic events that have affected the current distributions of W Mediterranean biota. The tectonic migration of an assemblage of continental micro-plates originally present at the margins of northeastern Spain and southern France in the Oligocene (Cavazza and Wezel, 2003; Rosenbaum et al., 2002; Rosenbaum and Lister, 2004), and which match with the current positions of Calabria (Italy), Sicily, Corsica, Sardinia, the Baetic-Rifean belt (composed by actual Baetic mountains in S of Iberian Peninsula and Rif mountains in the N of ⇑ Corresponding author. Fax: +34 954557049. E-mail address: [email protected] (J.M. Sánchez-Robles). http://dx.doi.org/10.1016/j.ympev.2014.06.005 1055-7903/Ó 2014 Elsevier Inc. All rights reserved.

Morocco), Balearic Islands and the Kabylie (Algeria) was a major determinant of dispersal and isolation for some Mediterranean species (Magri et al., 2007; Pfenninger et al., 2010). In contrast, the cycles of desiccation and transgression of the Mediterranean Sea and the Quaternary climatic oscillations enabled direct biotic interchange between the marine basins (Bocquet et al., 1978) but also created effective barriers between the African and European margins for many species (Krijgsman et al., 1999; Hewitt, 1999, 2000). These oscillations have deeply modelled the genetic structure and spatial distribution of biota creating complex phylogeographical patterns. Abies (Pinaceae) is a widely distributed monophyletic genus of Northern Hemisphere conifers (firs) comprising around 50 species (Xiang et al., 2009; Semerikova and Semerikov, 2014), and two periods were key for its diversification: the Eocene and the Miocene (Aguirre-Planter et al., 2012). Firs are distributed in three disjunct areas in the world: North America, East Asia and Southern Europe-Mediterranean. Migration routes to the southern

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

Europe-Mediterranean area are uncertain, but it is assumed that Abies populations migrated southwards from more northerly latitudes due to global climatic cooling-drying in the Eocene (Xiang et al., 2007), and a Miocene divergence between Southern Europe-Mediterranean and Asian–North American species has been estimated (Semerikova and Semerikov, 2014). The taxonomy of Abies species in Mediterranean Basin is somewhat controversial. The genus has two sections in this area: section Abies and section Piceaster (Farjon and Rushforth, 1989), but phylogenetic studies of these sections (i.e. Suyama et al., 2000; Xiang et al., 2004; Xiang et al., 2009; Semerikova and Semerikov, 2014) could not discern clearly both sections. The sect. Abies includes the widespread A. alba Miller (with a fragmented distribution ranging from the Pyrenees to the Carpathians), and a series of endemic species with restricted distributions: A. cephalonica Loudon (Greece), A. nordmanniana Spach (NE Turkey and Black Sea area), A. nebrodensis (Lojac.) Mattei (Sicily), and A. cilicica (Antoine & Kotschy) Carrière (SE Turkey). The sect. Piceaster includes only four species in the SW Mediterranean Basin: A. pinsapo Boiss. (restricted to three populations in Southern Spain), A. numidica Carr. (endemic to Babor mountains, in Algeria), and the Moroccan species A. marocana Trab. (comprising a few small populations in the Rif Mountain Range) and A. tazaotana Villar (a narrow endemic restricted to a single population on Jbel Tazaot). In sect. Abies, many studies have been reported on the genetic population structure of the widespread A. alba (Liepelt et al., 2010; Parducci et al., 1996; Vendramin et al., 1999; Ziegenhagen et al., 1998), and some surveys have also been conducted on the

43

Central and East Mediterranean fir species (Fady and Conkle, 1993; Parducci et al., 2001a,b; Scaltsoyiannes et al., 1999). In sect. Piceaster, some studies have suggested differentiation between A. pinsapo and the Moroccan firs (isozymes: Pascual et al., 1993; biometrical comparison of the needles: Se˛kiewicz et al., 2013; chloroplast microsatellites: Terrab et al., 2007). However, there are no detailed studies on all species of sect. Piceaster. Due to the restricted distribution area and location of Abies sect. Piceaster (Fig. 1), this group of firs offers a good model to explore how geo-climatic oscillations might have influenced expansion and diversification on both sides of the W Mediterranean basin. Schematically, two hypotheses can be postulated to explain the present distribution of species from sect. Piceaster in the W Mediterranean area. Abies pinsapo and North African firs (A. marocana, A. tazaotana and A. numidica) may have evolved from a single ancestor that originated in the Baetic-Rifean belt (see Fig. 1) and expanded eastwards, as proposed by Esteban et al. (2009); alternatively, North Africa may have been colonized independently on several occasions (i.e. western and eastern origins) by European firs. We investigated the genetic structure and phylogeography of Abies sect. Piceaster in the SW Mediterranean, i.e., A. numidica, A. marocana, A. tazaotana and A. pinsapo, and we also sampled one population of Abies alba from its only location area in the Iberian Peninsula, the Pyrenees. According to Liu (1971) and Parducci (2000), the migration of Abies to S Iberian Peninsula and N Africa would have occurred through this range of mountains. Our aims were to: (i) explore the distribution of genetic variation within

Fig. 1. Distribution map of the nine Abies populations sampled corresponding to the five studied species. (a) Global map; (b) Study area with sampled populations indicated by color circles and distribution area of species indicated as follows: green, A. pinsapo; brown, A. marocana; red, A. tazaotana; yellow, A. numidica; blue, A. alba; (c) Detail of species location around the Strait of Gibraltar. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

44

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

and among the disjunct areas that constitute the present ranges of these species; (ii) determine the phylogeographic relationships among sect. Piceaster species. We used Amplified Fragment Length Polymorphisms (AFLP), a very powerful fingerprinting technique for phylogenetic reconstruction in closely related organisms and recent species radiations (Meudt and Clarke, 2007) because of the low phylogenetic resolution of Abies sec. Piceaster. Additionally, chloroplast microsatellites (cpSSR) were used to complement nuclear markers. Events of longdistance gene flow by pollen could be detected because chloroplasts in gymnosperms show paternal inheritance (Birky, 1995). Furthermore, mutation rate for cpSSR is lower than for the nuclear genome (Provan et al., 1999; Wolfe et al., 1987), allowing us to have an estimate of early evolutionary events. 2. Materials and methods 2.1. Plant materials, sampling site and DNA extraction We obtained plant material from nine populations of five Abies taxa: three populations from SW Spain (A. pinsapo), four from N Morocco (A. marocana and A. tazaotana), one from NE Algeria (A. numidica) and one from NE Spain (A. alba) (Fig. 1). All sampled populations were geo-referenced by means of GPS (Garmin Geko 201, Garmin International Inc., Olathe, KS, USA) and vouchers were deposited in the Herbarium of the University of Seville, Spain (SEV). Table 1 provides details for the sampled populations. Leaves were collected from a total of 201 individuals, from 10 to 28 individuals per population, and the samples were stored in silica gel until DNA extraction. Genomic DNA was extracted using a QIAGEN DNeasy Plant Mini Kit (QIAGEN, Valencia, California, USA). 2.2. Amplified fragment length polymorphisms (AFLP) analysis The AFLP procedure followed established protocols (Vos et al., 1995; PE Applied Biosystems, 1996). Genomic DNA (50–100 ng) was digested with two restriction endonucleases (EcoRI and MseI) and ligated to double-stranded adaptors (EcoRI and MseI) at 37 °C for 2 h. Ligated DNA fragments were diluted 20-fold with TE0.1 buffer. Pre-selective primers, based on the EcoRI and MseI adaptors, were used to amplify a subset of fragments with matching nucleotides (downstream from the restriction sites). The pre-selective PCR products were diluted 20-fold with TE0.1 buffer and fragments were amplified with fluorescent-labelled selective primers. An initial screening of selective primers was performed on eight individuals from eight populations. Out of 102 primer combinations tested (with three and four nucleotides), six primer combinations were selected because of their marked polymorphism: EcoRI (6FAM)-ACA/MseI-CATA, EcoRI (6-FAM)-ACC/MseI-CATA, EcoRI (VIC)-AGG/MseI CAGG, EcoRI (VIC)-ACG/MseI CACC, EcoRI (NED)ACC/MseI-CACC and EcoRI (NED)-ACC/MseI-CAGG. To calculate the error rate of the method, replicates of the AFLP protocol (from new DNA extractions) were conducted on three individual plants for each population (27 samples, 13% of the total). 2.3. Chloroplast microsatellite (cpSSR) analysis Ten individuals per population (90 individuals in total) were genotyped using five cpSSR loci. Primer pairs were designed on the basis of the published chloroplast genome sequences of Pinus thunbergii (Wakasugi et al., 1994). The primers pt15169, pt63718 and pt71936 (GenBank accession no. U82923) were developed by Vendramin and Ziegenhagen (1997) whereas primers pt30141 (GenBank accession no. AF367975) and pt30249 (GenBank accession no. AF367976) were derived from the composed locus

pt30204 (GenBank accession no. U82922) to avoid homoplasy (Liepelt et al., 2001). Polymerase chain reaction (PCR) was conducted with a Veriti Thermal Cycler (Applied Biosystems, Foster City, CA, US) in a final volume of 10 lL using 5 lL of master mix (Bioline, London, United Kingdom), 50–100 ng of DNA template and 4 lM each forward and reverse primer (forward primers were labelled with fluorescent dye: 6-FAM for Pt71936 and Pt30141 and VIC for Pt15169, Pt63718 and Pt30249). PCR cycle profiles consisted of an initial denaturation step at 95 °C for 3 min, 35 cycles of 94 °C for 1 min, 55 °C for 1 min and 72 °C for 1 min, and a final step of 8 min at 72 °C. To calculate cpSSR error rate we tested PCR repeatability in one randomly selected sample per population (nine samples; 10% of the total), starting from the same extraction. Amplification and labelling of amplified products of both molecular markers (AFLP and cpSSR) were performed at the Centro de Investigación, Tecnología e Innovación de la Universidad de Sevilla (Seville, Spain). For each individual, 0.5 ll of 6-FAMlabelled, 0.5 ll of 6-NED-labelled and 0.5 ll of VIC-labelled selective PCR products were combined with 0.5 ll of GeneScan 500 LIZ (Applied Biosystems) and 13.5 ll of formamide. The mix was run on a capillary sequencer (ABI 3730; Applied Biosystems) at the Unidad de Genómica (Parque Científico de Madrid-Universidad Complutense, Madrid, Spain), and the GENEMARKER 1.8 software (SoftGenetics, State College, Pennsylvania, USA) was used to score amplified AFLP fragments 100–500 bp in length and cpSSR fragments. 2.4. AFLP molecular analysis AFLP scores were treated as binary (presence/absence) data and imported into R 3.0.1 software (R Development Core Team, 2013). Genetic diversity was assessed for each population using the total number of AFLP fragments scored (Fragtot), the percentage of polymorphic fragments (Fragpoly), and the number of private fragments (Fragpriv). Nei´s average gene diversity (HD) and Rarity Index, (i.e. DW or frequency-down-weighted marker value; Schönswetter and Tribsch, 2005) was calculated with the R-script AFLPdat package with modifications (Ehrich, 2006; last modified 20.10.2010). To avoid any sampling size bias, the average values and confidence intervals (bias-corrected and accelerated –BCa–; Efron and Tibshirani, 1986) of the population DW estimates were obtained by nonparametric bootstrapping at the plant individual level for each population from 1000 repetitions (Balao et al., 2010). Furthermore, allelic richness (Ar) and private-allele richness (PAr) per species were also calculated and standardized to the smallest sample size by means of the rarefaction method with the Adze software (Szpiech et al., 2008). To represent overall genetic relationships among all the 201 individuals analyzed, we constructed a dendrogram applying the neighbor-joining method (NJ) based on Nei and Li (1979) genetic distances (ape 3.0–8 package in R software; Paradis et al., 2004). Support for each node was tested by 1000 bootstrap replicates. 2.5. cpSSR molecular analysis Allele combinations across the five cpSSR loci were treated as a haplotype, considering different haplotypes as a unique combination of size variants of each cpSSR locus. A minimum spanning network (MSN) was computed in Arlequin ver. 3.5 (Excoffier et al., 2005) from a distance matrix based on the sum of squared differences of each allele to show the relationship between the different haplotypes. Haplotype variation within each population was calculated by estimating the number of haplotypes (Nh), the number of private haplotypes (Nhpriv) and the genetic diversity (He) using Arlequin and Spagedi ver 1.4 (Hardy and Vekemans, 2002). Allelic richness (Ar) and private-allele richness (PAr) per species was also

Table 1 Features of the Abies sp. populations sampled. N represents the number of individuals sampled. Fragtot: number of fragments per population, Fragpoly: percentage of polymorphic fragments, Fragpriv/PAr: number of private fragments and private allelic richness rarefacted per species, in brackets if they are fixed, HD/Ar: Nei´s average gene diversity (95% -IC) and private allelic richness rarefacted per species, DW: index of rare fragments (95% CI); Nh/Ar: number of haplotypes; Nhpriv/PAr: number of private haplotypes and private allelic richness per species; He: gene diversity and allelic richness rarefacted per species.

A. pinsapo 1 2 3 Total Rarefacted A. marocana 4 5 6 Total Rarefacted A. tazaotana 7 Rarefacted A. numidica 8

Location

Region

Coordinates

AFLP

cpSSR

N

Fragtot

Fragpoly

Fragpriv/PAr

HD/Ar

DW

Nh/Ar

Nhpriv/PAr

He/Ar

Sierra de Grazalema Sierra de las Nieves Sierra Bermeja

SW Spain SW Spain SW Spain

36°460 N/5°250 W 36°400 N/5°030 W 36°290 N/5°110 W

23 23 25 71

371 380 359 405

77.63 78.16 77.16 88.40

0 3 0 31 0.07

0.14 0.16 0.14 0.16 1.41

(0.13–0.16) (0.14–0.17) (0.13–0.16) (0.14–0.17)

3.04 3.04 2.64 2.90

(2.77–3.32) (2.78–3.24) (2.44–2.88) (2.76–3.05)

10 10 10 30

6 4 5 9

1 1 2 8 0.63

0.52 0.28 0.52 0.43 2.84

(0.34) (0.22) (0.31) (0.33)

Jbel Azilane Jbel Gharbouch N. P. Talassemtane

NW Morocco NW Morocco NW Morocco

35°110 N/5°130 W 35°070 N/5°040 W 35°070 N/5°080 W

24 28 22 74

357 349 325 381

82.35 79.37 72.31 87.40

0 2 0 5 0.02

0.15 0.15 0.12 0.15 1.38

(0.13–0.16) (0.13–0.16) (0.11–0.14) (0.13–0.16)

2.08 2.21 2.01 2.11

(1.94–2.28) (2.09–2.34) (1.93–2.10) (2.04–2.19)

10 10 10 30

8 7 6 12

2 1 1 6 0

0.38 0.36 0.34 0.37 2.12

(0.36) (0.37) (0.31) (0.35)

Jbel Tazaot

NW Morocco

35°160 N/5°060 W

24

419

89.98

15 0.04

0.20 (0.18–0.22) 1.53

3.43 (3.00–4.20)

10

5

1 0.01

0.35 (0.35) 2.00

Seraïdi (arboretum from Jbel babor)

NE Algeria

36°540 N/7°590 E

22

400

93.75

31

0.22 (0.20–0.23)

4.63 (4.21–5.08)

10

8

8

0.34 (0.31)

0.07

1.55

0.62

2.40

13(2) 0.05

0.15 (0.14–0.17) 1.41

8 0.42

0.51 (0.24) 3.20

Rarefacted A. alba 9 Rarefacted

N

Pyrenees

NE Spain

42°630 N/1°080 E

10

287

85.02

4.08 (3.50–4.72)

10

8

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

Taxa/population

45

46

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

assessed in the same way than AFLP, following a rarefaction method (see above). To test for phylogeographic signals, we compared two coefficients of genetic differentiation between populations (GST and NST; Pons and Petit, 1996). The test compared observed NST to random permuted NST (pNST = GST). Estimates of GST are based on allele frequencies and follows an infinite allele mutational model (IAM) for microsatellites, whereas NST takes into account the distance between different haplotypes, following a step-wise mutational model (SMM). Significantly higher NST than GST values indicate the presence of phylogeographic structure, i.e. different alleles are more related within population than among populations (Lynch and Crease, 1990; Pons and Petit, 1996). The test was performed using Spagedi through 20,000 permutations. NST was computed using a matrix of genetic distances based in the sum of the squared differences (RST; Hardy et al., 2003). In addition, the test was also performed considering pairwise comparison at different grouping in order to identify at which scale the phylogeographic signal occurs. 2.6. Genetic structure and genetic barriers Different approaches were used for both molecular markers (AFLP and cpSSR) to investigate the genetic structure and identify genetic barriers existing among all populations. We used analyses of molecular variance (AMOVAs; Arlequin ver 3.5; Excoffier et al., 2005) to estimate the genetic variation attributed to differences within and between predefined hierarchical groups of populations. These hierarchical groups of populations were made to test different phylogeographical scenarios. The first grouping distinguished A. pinsapo populations from the rest of the studied species. The second grouping distinguished the A. marocana – A. tazaotana populations from the rest of studied species. The third grouping distinguished three groups: (a) A. numidica and A. alba populations, (b) A. pinsapo populations and (c) A. marocana – A. tazaotana populations complex. The fourth grouping distinguished: (a) A. numidica population, (b) A. alba population, (c) A.

pinsapo populations and (d) A. marocana – A. tazaotana populations complex. The fifth grouping also distinguished four groups: (a) A. numidica and A. alba populations, (b) A. pinsapo populations, (c) A. marocana populations and (d) A. tazaotana population. Finally, a grouping considering the different populations studied by species level was also studied to confirm the current taxa delimitation. Variance components were tested for significance using an exact non-parametric test with 10,000 permutations. A Bayesian Markov Chain Monte Carlo (MCMC) clustering model was performed in Structure ver. 2.3.3 (Pritchard et al., 2000) to assess the genetic structure of the different populations. An admixture model with correlated allele frequencies and 10 independent repetitions for each K (ranging from 1 to 9) using the same patterns were employed. The number of iterations was 106, with a burn-in period of 150,000. Analyses were run on the Bioportal, University of Oslo (www.bioportal.uio.no). To determine the optimal number of clusters (K), we checked the Structure results on Structure Harvester software (Earl, 2012), which implements the Evanno et al. (2005) method using the ad hoc statistic DK. Moreover, to visualize populations spreading on a three dimensional space, a Principal Coordinate Analysis (PCoA) based on AFLP Nei–Li’s distances (Nei and Li, 1979) was performed at the individual level. Relationships among chloroplast haplotypes were also inferred using a PCoA (also based on Nei–Li’s distances). Lastly, genetic barriers (i.e. geographic areas with pronounced genetic discontinuity between the populations) were explored with the software Barrier ver. 2.2 (Manni et al., 2004). This approach combined geometry and a Monmonier’s maximum-difference algorithm to identify edges associated with high different rates based on a genetic distance matrix using Delaunay triangulation and Voronoi polygons. Briefly, Delaunay triangulation is the fastest triangulation method to connect a set of populations on a map by a set of triangles. Voronoi tessellation is obtained from the intersection of the medians of the triangles. Two populations are neighbors if the corresponding Voronoi polygons have a common edge. Once a network connecting all the populations is

A. pinsapo

A. numidica

Pops. 1-3

Pop. 8

100

95.7

A. alba

87.9

Pop. 9

100 97

A. marocana A. tazaotana Pops. 4-7

Fig. 2. Neighbor-joining dendrogram based on Nei-Li distances for AFLP markers in nine Abies populations. Numbers alongside branches indicate the bootstrap support (>60 bootstrap).

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

obtained and a corresponding distance matrix is available, the Monmonier’s algorithm computes barriers that lie on the Voronoi tessellation. Starting from the edge for which the distance value is maximum, and proceeding across adjacent edges, we crossed the edge of the triangle whose distance value is higher. To evaluate barriers robustness, population pairwise genetic distance matrix based on Jost’s D index (Jost, 2008) was bootstrapped 1000 times with the R package boot (Canty, 2002) and the 1000 data replicates were run in Barrier software. 3. Results 3.1. AFLP analyses The use of six AFLP primer combinations generated a total of 590 highly reproducible fragments (<2% error rate). Most of these fragments (587) were polymorphic and only three were monomorphic. Detailed estimates of genetic diversity analyses are shown in Table 1. Total number of fragments varied between species from 287 in A. alba to 419 in A. tazaotana, whilst the percentage of polymorphic fragments varied between 85.0% in A. alba and 93.7% in A. numidica. Private fragments varied from 5 (A. marocana) to 31 (A. pinsapo and A. numidica). Within-population Nei´s average gene diversity ranged between 0.15 (A. marocana) and 0.22 (A. numidica). This latter species also exhibited the highest DW (4.63) and A. marocana showed the lowest one (2.11). Values of allelic richness and private-allelic richness after rarefaction suggested that A. marocana is the least diverse of the species (1.38 and 0.02, respectively) and A. numidica the most diverse (1.55 and 0.07, respectively). Genetic relationships between populations are indicated in the neighbor-joining phylogram (Fig. 2). Populations were structured into three distinct and highly supported clusters. One cluster corresponded to A. pinsapo individuals (100% bootstrap). The A. maroccana – A. tazaotana individuals appeared intermingled in another cluster (97% bootstrap). The last cluster corresponded to A. alba – A. numidica individuals (100% bootstrap). However, A. alba individ-

47

uals were highly differentiated from A. numidica individuals in a different cluster (87.9% bootstrap). Intra-species relationships were generally poorly resolved (<60% bootstrap), with the exception of two A. numidica individuals, that appeared to have a highly differentiated lineage from the rest (95.7% bootstrap). 3.2. cpSSR analyses The five cpSSR loci showed high reproducibility (global error rate 3.3%). All loci were polymorphic, ranged from 10 alleles (Pt71936) to three alleles (Pt30249 and Pt15169). Alleles for loci Pt63718 and Pt30141 were three and eight, respectively. A total of 37 different haplotypes were found through the five loci combinations. Detailed estimates of genetic diversity analyses among populations are shown in Table 1. Number of haplotypes for each species varied from five in A. tazaotana to 12 in A. marocana. All species exhibited private haplotypes, ranging from one (A. tazaotana) to eight (A. pinsapo, A. numidica and A. alba). Gene diversity per population showed clinal pattern in a north–south direction varying between 0.51 in A. alba and 0.34 in A. numidica. Allelic richness per species after rarefaction showed the same pattern. A. pinsapo showed the highest private allele richness after rarefaction (0.63) and A. marocana did not show any uniqueness. The coefficient of global genetic differentiation between populations (NST = 0.65) was significantly higher (P = 0.03) than the permuted value (GST = 0.57), indicating a global phylogeographical signal (i.e. microsatellite evolution followed a Stepwise-Mutational Model, SSM). However, not all pairwise comparisons were significant for the phylogeographical test. There was no phylogeographical signal among the three A. pinsapo populations (P = 0.55), and among populations within A. marocana – A. tazaotana complex (P = 0.54). Furthermore, and congruently with its multiple connections in the MSN, A. alba population did not show a significant phylogeographical signal when was compared with A. pinsapo populations (P = 0.38) and A. numidica populations (P = 0.37), indicating a low genetic differentiation among these species

Fig. 3. Minimum Spanning Network (MSN) among the 37 haplotypes found in nine Abies populations from SW Mediterranean. Genetic distances are based on the sum of the squared differences of distances of each locus. Small black circles represent one distance unit between haplotypes and thick black lines represent higher distance values than five units.

48

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

populations. Interestingly, the pairwise comparisons between A. marocana – A. tazaotana populations complex and the different populations of the studied species were all significant (P = 0.01 for comparison with A. pinsapo populations; P = 0.01 for comparison with A. numidica population; P = 0.03 for comparison with A. alba population), indicating strong genetic differentiation. Similarly, the pairwise test for A. pinsapo populations and A. numidica population was significant (P < 0.01). The network of the 37 cpSSR haplotypes, based on a microsatellite Stepwise-Mutational Model (SSM) according to the phylogeographical test, showed a complex pattern (Fig. 3), with three intermixed groups as follows. Abies marocana – A. tazaotana complex formed a large cluster occupying an intermediate position and sharing four haplotypes (MT9, TM14, MT15 and MT18). Abies pinsapo haplotypes showed two patterns: on the one hand, a defined group (P29–P33), distant from A. marocana – A. tazaotana complex by a large genetic distance but relatively close to most of haplotypes of A. alba. On the other hand, three haplotypes (P34–P36) appeared close to several haplotypes of A. marocana – A. tazaotana complex, whilst one haplotype (MP19) was shared between A. pinsapo and A. marocana. The third group was formed by A. numidica haplotypes (N21–N27), which was linked to A. marocana – A. tazaotana complex through one A. alba haplotype (A8). However, haplotype N22 was separated by a long genetic distance from the rest of A. numidica haplotypes. This pattern also appeared in the NJ tree, where some individuals of A. numidica were differentiated from the rest with a high bootstrap support (see above). Lastly, haplotypes belonging to A. alba were linked to all species in different position in the network.

3.3. Genetic structure and genetic barrier analyses The most explanatory grouping in AMOVA for both marker types, was the fourth grouping: (a) A. numidica population, (b) A. alba population, (c) A. pinsapo populations and (d) A. marocana – A. tazaotana populations complex (43.48% for AFLP and 68.09% for cpSSR; see Table 2). This grouping explained slightly more genetic variability than the grouping considering the different populations studied at the species level (43.19% for AFLP and 64.23% for cpSSR). Bayesian clustering analyses were partially incongruent for both markers (AFLP and cpSSR). Although Structure analyses (Fig. 4) showed that the most informative representation of overall genetic structure was achieved with K = 2 for both markers (See Supplementary Fig. 1), the clusters obtained with AFLP and cpSSR markers were different. For AFLP, the first cluster was very homogeneous and consisted exclusively of the A. pinsapo individuals and the second cluster comprised individuals of the remaining species. For cpSSR, the first cluster consisted exclusively of the A. marocana – A. tazaotana individuals and the second cluster comprised individuals of the remaining species. Remarkably, nine individuals of A. pinsapo (most of them from the geographically closest Sierra Bermeja population) had a high probability of being ‘misplaced’ to the A. marocana – A. tazaotana cluster. However, the suboptimal clustering K = 4 was congruent for both markers (Fig. 4) and also with the most explanatory grouping in the AMOVA. A very homogeneous cluster was formed by A. numidica – A. alba individuals, whereas A. marocana – A. tazaotana individuals were assigned to two different groups. However these latter clusters were very

Table 2 Analysis of molecular variance (AMOVA) for AFLP and cpSSR markers employed in the studied Abies populations. For cpSSR, Sum of squares (SS), variance components and percentage of variation are given for the R-statistics. Groupings

N AFLP

A. pinsapo vs. rest of species 2

A. marocana – A. tazaotana vs. rest of species

A. numidica – A. alba vs. A. pinsapo vs. A. marocana – A. tazaotana

A. numidica vs. A. alba vs. A. pinsapo vs. A. marocana – A. tazaotana A. numidica – A. alba vs. A. pinsapo vs. A. marocana

2

3

4

4

vs. A. tazaotana Species

*

P < 0.05.

5

cpSSR

Source of variation

d.f.

SS

Among groups Among populations Among individuals Among groups Among populations Among individuals Among groups

1 7

Variance component

% of variance

R-values

3565.97 34.54 2602.64 14.75

35.88 15.32

RcT = 0.36* 1 7

192 9022.44 46.99

48.81

1 7

3192.26 27.32 2976.36 17.16

29.87 18.76

192 9022.44 46.99

51.37

2

5079.72 38.40

42.02

Among populations Among individuals Among groups Among populations Among individuals Among groups

6

1088.90 6.01

6.57

Among populations Among individuals Among groups Among populations Among individuals

192 9022.44 46.99 3 5

d.f. SS

% of variance

R-values

3144.04 65.07 3787.32 49.62

40.78 31.10

RcT = 0.31

3634

44.86

28.12

3650.63 71.59 3280.73 42.38

45.07 26.68

81

3634

44.86

28.25

2

4309.49 59.46

41.42

6

2621.88 39.21

27.32

81

3634

31.26

3 5

6507.64 104.22 423.72 3.99

68.09 2.61

81

3634.00 44.86

29.31

RcT = 0.41* 3

4313.27 42.55

31.45

81 RcT = 0.30* 1 7

RcT = 0.42

*

51.41 RcT = 0.43

*

Variance component

44.86

5376.35 39.70 792.27 4.62

43.48 5.06

192 9022.44 46.99

51.46

3

5452.37 35.65

40.97

5

716.24

5.02

5

2618.08 47.88

35.39

192 9022.44 46.99

54.01

81

3634.00 44.86

33.16

4 4

5749.01 37.55 419.61 2.40

43.19 2.76

6511.42 91.37 419.93 6.01

64.23 4.23

192 9022.44 46.99

54.05

3634

31.54

4.37

RcT = 0.43* 4 4 81

44.86

RcT = 0.45*

RcT = 0.41*

RcT = 0.68*

RcT = 0.31

RcT = 0.64*

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

49

Fig. 4. Admixture analyses of nine Abies populations from SW Mediterranean basin, performed using Structure software with K = 2 and 4, for AFLP and cpSSR markers. Each horizontal bar for each marker represents a different individual. Bars are divided in different segments with their length proportional to the estimated membership in the clusters.

hereterogeneous and did not distinguish between Moroccan fir species. Lastly, A. pinsapo individuals were also grouped into a distinct cluster with several misclassified in cpSSR. Nevertheless, Abies pinsapo individuals were not misclassified within the A. marocana – A. tazaotana cluster relative to their origin (S. Bermeja or S. Grazalema). The PCoA based on the AFLP data (Fig. 5A), was congruent with the NJ tree, showing a clear separation among A. pinsapo, A. alba – A. numidica and A. marocana – A. tazaotana, with the exception of six individuals of A. numidica that appear in an intermediate position between the A. alba – A. numidica complex and A. marocana – A. tazaotana complex. Moreover, the PCoA based on the cpSSR data (Fig. 5B), showed a slight relationship between A. pinsapo and A. marocana – A. tazaotana, based in haplotypes located in central positions of the MSN (Fig. 3). The maximal genetic differentiation indicated by Barrier was completely congruent with Structure results (K = 2). For AFLP, A. pinsapo populations were the most genetically differentiated and isolated ones from the rest of the studied taxa, while for cpSSR, the A. marocana – A. tazaotana populations complex were the most separated ones (Fig. 6). In spite of these differences, for both markers, Monmonier’s maximum-difference algorithm highlighted the important biogeographic barrier imposed by the Mediterranean Sea through the Strait of Gibraltar between A. pinsapo and the Moroccan taxa.

4. Discussion 4.1. Genetic diversity in Western-Mediterranean firs This is the first, up-to-date, study using genome-wide AFLP markers for Abies sect. Piceaster. Mediterranean firs showed similar

levels of gene diversity for these markers as that found in the Asian firs (Semerikova et al., 2011; Tang et al., 2008). Additionally, the gene diversities for each species and the strong genetic differentiation observed among species with cpSSR markers were congruent with previous studies in Mediterranean firs using the same molecular markers (Jaramillo-Correa et al., 2010; Liepelt et al., 2010; Parducci et al., 2001a; Terrab et al., 2007; Vendramin et al., 1999). However, the observed gene diversity for A. pinsapo, A. marocana and A. tazaotana in this study contrasts with the higher values observed by Terrab et al. (2007), in a study focused on these three species in which the same molecular marker was used; It is likely that the different primer pairs used by Terrab et al. (only three polymorphic loci versus five polymorphic loci used in the present study), and the different sampling size (33 individuals per populations versus 10 individuals per populations sampled in the present study) may explain this difference. The large genetic distance found among A. pinsapo haplotypes in the MSN (see Fig. 3) was also found by Terrab et al. (2007). A 21 bp indel in the locus pt15169 was responsible of this large genetic distance between haplotypes (See Supplementary Table 1). The A. numidica haplotype (N22), also separated by large genetic distance from the other haplotypes of this species, also had a large indel. And the indel presence also explained the large genetic distance among A. alba haplotypes and haplotypes of Spanish and Moroccan firs. The low gene diversity for the silver-fir (A. alba) occurring in the Pyrenees was in concordance with previous studies of this species in this region (Liepelt et al., 2010; Vendramin et al., 1999), which has been considered a glacial refugium (Hewitt, 2000; Konnert and Bergmann, 1995; Médail and Diadema, 2009). However, the low genetic diversity observed (cpSSR) across southern areas is consistent with natural range expansion towards the south (Taberlet et al., 1998).

50

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

1: Abies pinsapo 2: Abies pinsapo 3: Abies pinsapo 4: Abies marocana

A

5: Abies marocana 6: Abies marocana 7: Abies tazaotana 8: Abies numidica 9: Abies alba

0.10

0.05

0.04

2 (1

0.00

ate

0.02 0.00

) 4%

-0.02

-0.10 -0.05

rdin C oo

0.06

Coordinate 3 (7%)

0.10

0.08

-0.04 -0.06

-0.04 -0.03 -0.02 -0.01 0.00 0.01 0.02 0.03 0.04

Coordinate 1 (75%)

B

0.03

0.01

0.02

-0.03

-0.01

-0.04 -0.05 -0.010

-0.02 -0.005

0.000

0.005

0.010

-0.03

ate

0.00

2(

0.01

-0.02

12% )

0.03

0.00 -0.01

Co ord in

Coordinate 3 (4%)

0.02

Coordinate 1 (79%) Fig. 5. Principal coordinates analysis (3D-PCoA) based on Nei-Li distances for (A) AFLP at the individual level and (B) cpSSR at the haplotype level.

4.2. The origin and genetic structure of the Spanish and the Moroccan firs Although abundant Tertiary macrofossils and microfossils of the genus Abies have been found across Europe (from NE Iberian Peninsula to N Italy) throughout the Eocene (Alcalá, 1997; Cavagnetto and Anadón, 1996; Gaussen, 1964; Jimenez-Moreno, 2005; Jiménez-Moreno et al., 2010; Palamarev, 1989; Rivas-Carvallo et al., 1994; Roiron et al., 1999), the fossil pollen record only confirms the presence of Abies in the Southern Iberian Peninsula since the Miocene (Jiménez-Moreno, 2005; Jiménez-Moreno and Suc, 2007), a period during which Abies fossils were absent on the eastern coast of the Iberian Peninsula, with the exception of the Balearic Islands (See Supplementary Fig. 2A). These pollen fossil records from S Spain are located in allochthonous terranes according to Rosenbaum and Lister (2004). The most likely hypothesis is that the vicariance events that promoted speciation of A. pinsapo and A. marocana – A. tazaotana in the Baetic-Rifean mountains followed the Tertiary palaeogeographical evolution of the W Mediterranean (Linares, 2011; Cavazza and Wezel, 2003;

Rosenbaum et al., 2002). Archaic Abies (the ancestor of A. pinsapo, A. marocana and A. tazaotana) present in the NE Iberian Peninsula and S France margins in the Eocene (where currently A. alba is located) may have rafted with the Baetic-Rifean microplate during the Oligocene–Miocene. This microplate was previously fragmented (15–10 Ma), originating the Baetic range in the S Iberian Peninsula and the Rif mountains in the NW Africa (Rosenbaum and Lister, 2004). Abies pinsapo has become markedly genetically differentiated from the North African firs as is shown in the AFLP phylogram (Fig. 2) and Bayesian assignment (Fig. 4). Vicariant speciation between A. pinsapo and A. marocana – A. tazaotana could have occurred, most likely, following the final opening of the Strait of Gibraltar (around 5.33 Ma), which imposed an efficient barrier against gene flow between the Southern Iberian and Moroccan populations (Esteban et al., 2009; Pascual et al., 1993; Se˛kiewicz et al., 2013; Terrab et al., 2007). This separation has been observed in other plant species and animals (Arroyo et al., 2008; Casimiro-Soriguer et al., 2010; Escudero et al., 2008; JaramilloCorrea et al., 2010; Lavergne et al., 2013; Lumaret et al., 2002;

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

Genetic break robustness

100%

51

Genetic break robustness

63-82%

Fig. 6. Genetic discontinuity found by Monmonier’s maximum-difference algorithm for the nine Abies populations studied by AFLP and cpSSR marker. Bootstraps supports are indicated for each genetic barrier.

Ortiz et al., 2007). In addition, this lack of gene flow is supported by the restricted pollen dispersal of Abies genus in general (Mazier et al., 2008) and A. pinsapo in particular (usually less than 3 km; Arista and Talavera, 1994; Sánchez-Robles et al., unpublished). Nevertheless, cpSSR analyses indicated a certain admixture between A. pinsapo and A. marocana – A. tazaotana clusters. Four cpSSR haplotypes of A. pinsapo were very close to the A. marocanna – A. tazaotana haplotypes (one of them shared) in the MSN (Fig. 3) and several individuals of A. pinsapo were misplaced to A. marocanna – A. tazaotana cluster in the Structure results (Fig. 4). As chloroplast markers have lower mutation rates than nuclear ones (Provan et al., 1999; Wolfe et al., 1987), these genetic admixture could have resulted from historical gene flow and re-colonization before the opening of Strait of Gibraltar, rather than a result of contemporary gene flow. Lastly, Bayesian analyses indicated no strong genetic structure within A. pinsapo or within A. marocana – A. tazaotana. This lack of genetic differentiation among populations of the same species, suggests that their current geographic areas are a result of a recent fragmentation, probably due to unfavorable weather conditions, fires, or human activity (Esteban et al., 2010; Linares et al., 2009), rather than current gene flow, given the restricted pollen dispersal of the Abies genus (Mazier et al., 2008) and the geographic isolation of studied populations. 4.3. The different origins of North African firs: integrating molecular and fossil data Bayesian analysis for AFLP and cpSSR (K = 4), the NJ phylogram and the AFLP-based PCoA defined three main clusters (A. pinsapo, A. marocana – A. tazaotana and A. numidica – A. alba). It thus seems likely that the diversification of Abies in North Africa occurred through at least two vicariance events, which is congruent with the presence in different clusters of the same haplotypes and haplotype lineages between European and African species (i.e. A. pinsapo – A. marocanna – A. tazaotana and A. alba – A. numidica). As discussed below, A. pinsapo, Abies marocana and A. tazaotana may have evolved from an Abies ancestor located in the NE Iberian Peninsula and S France. This Abies ancestor could have migrated with the Baetic-Rifean microplate (Rosenbaum and Lister, 2004), originating the extant species (i.e. A. marocanna – A. tazaotana and A. pinsapo) by posterior vicariance due to the opening of the Strait of Gibraltar. This hypothesis is supported by the close and shared haplotypes between A. marocanna – A. tazaotana and A. pinsapo. The same palaeogeographical vicariance events could be applied

to A. numidica as the current geographical extent (Algerian Kabylies) exactly matches the Kabylie microplate, which was also originally fragmented from NE Iberian Peninsula during the Miocene and rafted to current position (Rosenbaum et al., 2002). An alternative or complementary route of colonization for A. numidica may have been through the Italian Peninsula during the late Pliocene–Pleistocene. A relationship between the Algerian fir (A. numidica) and the southern A. alba (from Calabria, Italy) has been demonstrated by plastidial and isozyme markers (Parducci et al., 2001b; Scaltsoyiannes et al., 1999). Sicily and Tunisia may have shared a geological connection until at least the Late Pliocene–Pleistocene (Rosenbaum et al., 2002). So A. numidica may have remained in contact with the Sicilian fir, A. nebrodensis (Parducci et al., 2001b; Parducci and Szmith, 1999) and A. alba (discontinuously distributed along the Italian Peninsula, France and N Spain). Quaternary Abies pollen records in the NW Tunisia also support the hypothesis of a past contact between these species in the late Pliocene–Pleistocene (Tiba and Reille, 1982). Similar phylogeographical pattern of Algerian populations closely linked to those of Sardinia, Corsica or Southern France have been found in other gymnosperms (Burban and Petit, 2003; Terrab et al., 2008), angiosperms (FernándezMazuecos and Vargas, 2010; Magri et al., 2007) and also in animals (Zangari et al., 2006; Pfenninger et al., 2010). Finally, A. alba appeared less genetically related to A. numidica in the MSN (Fig. 3) in spite of the Structure results (Fig. 4), whereas both taxa are genetically very close. The microsatellite mutational model used to build the network may be behind this incongruence. Although the global phylogeographical test was significant, indicating that microsatellite alleles followed a Step-wise Mutational Model (SMM), the pairwise comparisons between A. alba and A. numidica (and also A. pinsapo) were not significant. Thus, the Infinite Allele Mutational model (IAM) for microsatellites could be more appropriate than SMM, and this could have affected the relationship found between the studied species.

4.4. Secondary contacts among species Genetic associations among A. pinsapo, A. alba and A. numidica were also elucidated from cpSSR Bayesian analyses (K = 2), but absent from AFLP. As mentioned above, nuclear markers have higher mutation rates than chloroplast ones (Provan et al., 1999; Wolfe et al., 1987). For this reason, if gene flow occurrs at present, genetic associations among species would be shown also for AFLP. As this is not the case, secondary contacts may have occurred in the past.

52

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53

The genetic separation of A. pinsapo from A. alba may have occurred before the separation of the latter and A. numidica, as AFLP analyses indicated (Figs. 2 and 4). A semi-arid wedge in SE Iberian Peninsula may have acted as a barrier, preventing the gene flow between A. pinsapo and A. alba. This semi-arid wedge seems to have been present since the Miocene (Geiger, 1970; Carrión et al., 2010) and only at the end of the Pliocene and in some Pleistocene periods was this region climatically propitious for Abies. This hypothesis is supported by pollen fossil records from the Quaternary in the area (Alba-Sánchez et al., 2010; Carrión et al., 2010; see Supplementary Fig. 2B), which were not found in the Tertiary, suggesting that SE Spain harbored fir forests in this age period. Although evolutionary convergence in haplotypes cannot be discarded, the multiple occurrences of A. alba haplotypes in such distant locations within the network (MSN; Fig. 3) support the hypothesis of secondary contacts between this species and species from sect. Piceaster. The distribution range may have undergone expansion and contraction due to climatic oscillation during the Pliocene–Pleistocene (Comes and Kadereit, 1998; Kropf et al., 2003), which may have led to contacts between species (by migration following favorable weather conditions) and the reduction to actual enclaves (by extinction in adverse weather conditions). 5. Conclusion This study demonstrated different origins for firs that compose Abies sect. Piceaster. The genetic relationship among species of this section seems to be complex, indicating vicariant events and different secondary contacts with A. alba, a fir from Abies sect. Abies. Although a close genetic relationship between A. alba and A. numidica was indicated, a larger sampling of Abies alba along Italy and Southern France is necessary to clarify the relationship between these two species. Furthermore, the endangered Sicilian fir, A. nebrodensis, could provide valuable information about the influence of Sicily channel in the African-Europe Abies connection. Acknowledgments This work was funded by FEDER funds allocated to projects from the Junta de Andalucía and the Spanish Ministerio de Ciencia y Tecnología to ST (P08-RNM-03703, and CGL2009-08178 and CGL2012-32914, respectively). JMS-R was supported by a PhD fellowship from the Junta de Andalucía. We acknowledge P. Gibbs for helpful comments on the manuscript and language review. Thanks to both anonymous referees for their invaluable comments and L. Navarro-Sampedro and M. Lorenzo for helpful advice in the laboratory. We also thank to Gérard de Bélair (Univ. of Annaba) who helped to collect A. numidica samples. Finally, we thank the Servicios Generales de Biología (CITIUS) from the University of Seville for allowing the use of their facilities. Appendix A. Supplementary material Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.ympev.2014.06. 005. References Aguirre-Planter, É., Jaramillo-Correa, J.P., Gómez-Acevedo, S., Khasa, D.P., Bousquet, J., Eguiarte, L.E., 2012. Phylogeny, diversification rates and species boundaries of Mesoamerican firs (Abies, Pinaceae) in a genus-wide context. Mol. Phylogenet. Evol. 62, 263–274. Alba-Sánchez, F., López-Sáez, J.A., Benito, B.M., Linares, J.C., Nieto-Lugilde, D., LópezMerino, L., 2010. Past and present potential distribution of the Iberian Abies species: a phytogeographic approach using fossil pollen data and species distribution models. Divers. Distrib. 16, 214–228. Alcalá, B., 1997. Prospección palinológica en el Neógeno de Teruel. Teruel 85, 7–20.

Arista, M., Talavera, S., 1994. Pollen dispersal capacity and pollen viability of Abies pinsapo Boiss. Silvae Genetica 43, 2–3. Arroyo, J., Aparicio, A., Albaladejo, R.G., Muñoz, J., Braza, R., 2008. Genetic structure and population differentiation of the Mediterranean pioneer spiny broom Calicotome villosa across the Strait of Gibraltar. Biol. J. Linn. Soc. 93, 39–51. Balao, F., Valente, L.M., Vargas, P., Herrera, J., Talavera, S., 2010. Radiative evolution of polyploid races of the Iberian carnation Dianthus broteri (Caryophyllaceae). New Phytol. 187, 542–551. Birky, C.W., 1995. Uniparental inheritance of mitochondrial and chloroplast genes: mechanisms and evolution. Proc. Natl. Acad. Sci. 92, 11331–11338. Bocquet, G., Widler, B., Kiefer, H., 1978. The Messinian model. A new outlook for the floristics and systematics of the Mediterranean area. Le modele messinien. Un nouveau point de vue pour la floristique et la systematique de la region mediterraneenne.). Candollea 33, 269–287. Burban, C., Petit, R.J., 2003. Phylogeography of maritime pine inferred with organelle markers having contrasted inheritance. Mol. Ecol. 12, 1487–1495. Canty, A.J., 2002. Resampling methods in R: the boot package. R News 2, 2–7. Carrión, J.S., Fernández, S., Jiménez-Moreno, G., Fauquette, S., Gil-Romera, G., González-Sampériz, P., Finlayson, C., 2010. The historical origins of aridity and vegetation degradation in southeastern Spain. J. Arid Environ. 74, 731– 736. Casimiro-Soriguer, R., Talavera, M., Balao, F., Terrab, A., Herrera, J., Talavera, S., 2010. Phylogeny and genetic structure of Erophaca (Leguminosae), a East-West Mediterranean disjunct genus from the Tertiary. Mol. Phylogenet. Evol. 56, 441– 450. Cavagnetto, C., Anadón, P., 1996. Preliminary palynological data on floristic and climatic changes during the Middle Eocene-Early Oligocene of the eastern Ebro Basin, northeast Spain. Rev. Palaeobot. Palynol. 92, 281–305. Cavazza, W., Wezel, F.C., 2003. The Mediterranean region – a geological primer. Episodes 26, 160–168. Comes, H.P., Kadereit, J.W., 1998. The effect of quaternary climatic changes on plant distribution and evolution. Trends Plant Sci. 3, 432–438. Earl, D.A., 2012. Structure Harvester: a website and program for visualizing Structure output and implementing the Evanno method. Conserv. Genetics Resour. 4, 359–361. Efron, B., Tibshirani, R., 1986. Bootstrap methods for standard errors, confidence intervals, and other measures of statistical accuracy. Stat. Sci., 54–75. Ehrich, D., 2006. AFLPdat: a collection of R functions for convenient handling of AFLP data. Mol. Ecol. Notes 6, 603–604. Escudero, M., Vargas, P., Valcárcel, V., Luceño, M., 2008. Strait of Gibraltar: an effective gene-flow barrier for wind-pollinated Carex helodes (Cyperaceae) as revealed by DNA sequences, AFLP, and cytogenetic variation. Am. J. Bot. 95, 745–755. Esteban, L.G., de Palacios, P., García Fernández, F., Martín, J.A., 2009. Wood anatomical relationships within Abies spp. from the Mediterranean area: a phyletic approach. Investigación Agraria: Sistemas y Recursos Forestales 18, 213–225. Esteban, L.G., de Palacios, P., Aguado, L.R.L., 2010. Abies pinsapo forests in Spain and Morocco: threats and conservation. Oryx 44, 276. Evanno, G., Regnaut, S., Goudet, J., 2005. Detecting the number of clusters of individuals using the software Structure: a simulation study. Mol. Ecol. 14, 2611–2620. Excoffier, L., Laval, G., Schneider, S., 2005. Arlequin (version 3.0): an integrated software package for population genetics data analysis. Evol. Bioinform. Online 1, 47–50. Fady, B., Conkle, M.T., 1993. Allozyme variation and possible phylogenetic implications in Abies cephalonica Loudon and some related eastern Mediterranean firs. Silvae Genetica 42, 351–351. Farjon, A., Rushforth, K.D., 1989. A classification of Abies Miller (Pinaceae). Notes R. Botanic Garden, Edinburgh 46, 59–79. Fernández-Mazuecos, M., Vargas, P., 2010. Ecological rather than geographical isolation dominates Quaternary formation of Mediterranean Cistus species. Mol. Ecol. 19, 1381–1395. Gaussen, H., 1964. Les gymnospermes actuelles et fossiles. Fasc. VI, Toulouse. Geiger, F., 1970. Die Aridität in Südostspanien. Stuttgarter Geographische Studien, Stuttgart, 173 pp. Hardy, O.J., Vekemans, X., 2002. SPAGeDi: a versatile computer program to analyse spatial genetic structure at the individual or population levels. Mol. Ecol. Notes 2, 618–620. Hardy, O.J., Charbonnel, N., Fréville, H., Heuertz, M., 2003. Microsatellite allele sizes: a simple test to assess their significance on genetic differentiation. Genetics 163, 1467–1482. Hewitt, G.M., 1999. Post-glacial re-colonization of European biota. Biol. J. Linn. Soc. 68, 87–112. Hewitt, G.M., 2000. The genetic legacy of the ice ages. Nature 405, 907–913. Jaramillo-Correa, J.P., Grivet, D., Terrab, A., Kurt, Y., De-Lucas, A.I., Wahid, N., Vendramin, G.G., González-Martínez, S.C., 2010. The Strait of Gibraltar as a major biogeographic barrier in Mediterranean conifers: a comparative phylogeographic survey. Mol. Ecol. 19, 5452–5468. Jiménez-Moreno, G., 2005. Utilización del análisis polínico para la reconstrucción de la vegetación, clima y estimación de paleoaltitudes a lo largo del arco alpino europeo durante el Mioceno (21–8 Ma). Ph.D. Dissertation. University of Granada, Spain, 313p. Jiménez-Moreno, G., Suc, J.P., 2007. Middle Miocene latitudinal climatic gradient in western Europe: evidence from pollen records. Palaeogeogr. Palaeoclimatol. Palaeoecol. 253, 208–225.

J.M. Sánchez-Robles et al. / Molecular Phylogenetics and Evolution 79 (2014) 42–53 Jiménez-Moreno, G., Fauquette, S., Suc, J.P., 2010. Miocene to Pliocene vegetation reconstruction and climate estimates in the Iberian Peninsula from pollen data. Rev. Palaeobot. Palynol. 162, 403–415. Jost, L., 2008. GST and its relatives do not measure differentiation. Mol. Ecol. 17, 4015–4026. Konnert, M., Bergmann, F., 1995. The geographical distribution of genetic variation of silver fir (Abies alba, Pinaceae) in relation to its migration history. Plant Syst. Evol. 196, 19–30. Krijgsman, W., Hilgen, F.J., Raffi, I., Sierro, F.J., Wilson, D.S., 1999. Chronology, causes and progression of the Messinian salinity crisis. Nature 400, 652–655. Kropf, M., Kadereit, J.W., Comes, H.P., 2003. Differential cycles of range contraction and expansion in European high mountain plants during the Late Quaternary: insights from Pritzelago alpina (L.) O. Kuntze (Brassicaceae). Mol. Ecol. 12, 931– 949. Lavergne, S., Hampe, A., Arroyo, J., 2013. In and out of Africa: how did the Strait of Gibraltar affect plant species migration and local diversification? J. Biogeogr. 40, 24–36. Liepelt, S., Kuhlenkamp, V., Anzidei, M., Vendramin, G.G., Ziegenhagen, B., 2001. Pitfalls in determining size homoplasy of microsatellite loci. Mol. Ecol. Notes 1, 332–335. Liepelt, S., Mayland-Quellhorst, E., Lahme, M., Ziegenhagen, B., 2010. Contrasting geographical patterns of ancient and modern genetic lineages in Mediterranean Abies species. Plant Syst. Evol. 284, 141–151. Linares, J.C., 2011. Biogeography and evolution of Abies (Pinaceae) in the Mediterranean Basin: the roles of long-term climatic change and glacial refugia. J. Biogeogr. 38, 619–630. Linares, J.C., Delgado-Huertas, A., Camarero, J.J., Merino, J., Carreira, J.A., 2009. Competition and drought limit the response of water-use efficiency to rising atmospheric carbon dioxide in the Mediterranean fir Abies pinsapo. Oecologia 161, 611–624. Liu, T., 1971. A monograph of the Genus Abies. Taiwan University. Department of Forestry, Taipei, Taiwan, 608 pp. Lumaret, R., Mir, C., Michaud, H., Raynal, V., 2002. Phylogeographical variation of chloroplast DNA in holm oak (Quercus ilex L.). Mol. Ecol. 11, 2327–2336. Lynch, M., Crease, T.J., 1990. The analysis of population survey data on DNA sequence variation. Mol. Biol. Evol. 7, 377–394. Magri, D., Fineschi, S., Bellarosa, R., Buonamici, A., Sebastiani, F., Schirone, B., Simeone, M.C., Vendramin, G.G., 2007. The distribution of Quercus suber chloroplast haplotypes matches the palaeogeographical history of the western Mediterranean. Mol. Ecol. 16, 5259–5266. Manni, F., Guerard, E., Heyer, E., 2004. Geographic patterns of (genetic, morphologic, linguistic) variation: how barriers can be detected by using Monmonier’s algorithm. Hum. Biol. 76, 173–190. Mazier, F., Broström, A., Gaillard, M.J., Sugita, S., Vittoz, P., Buttler, A., 2008. Pollen productivity estimates and relevant source area of pollen for selected plant taxa in a pasture woodland landscape of the Jura Mountains (Switzerland). Veget. Hist. Archaeobot. 17, 479–495. Médail, F., Diadema, K., 2009. Glacial refugia influence plant diversity patterns in the Mediterranean Basin. J. Biogeogr. 36, 1333–1345. Meudt, H.M., Clarke, A.C., 2007. Almost forgotten or latest practice? AFLP applications, analyses and advances. Trends Plant Sci. 12, 106–117. Myers, N., Mittermeier, R.A., Mittermeier, C.G., Da Fonseca, G.A., Kent, J., 2000. Biodiversity hotspots for conservation priorities. Nature 403, 853–858. Nei, M., Li, W.H., 1979. Mathematical model for studying genetic variation in terms of restriction endonucleases. Proc. Natl. Acad. Sci. 76, 5269–5273. Ortiz, M.A., Tremetsberger, K., Talavera, S., Stuessy, T., García-Castaño, J.L., 2007. Population structure of Hypochaeris salzmanniana DC. (Asteraceae), an endemic species to the Atlantic coast on both sides of the Strait of Gibraltar, in relation to Quaternary sea level changes. Mol. Ecol. 16, 541–552. Palamarev, E., 1989. Paleobotanical evidences of the Tertiary history and origin of the Mediterranean sclerophyll dendroflora. Plant Syst. Evol. 162, 93–107. Paradis, E., Claude, J., Strimmer, K., 2004. Ape: analyses of phylogenetics and evolution in R language. Bioinformatics 20, 289–290. Parducci, L., 2000. Genetics and evolution of the Mediterranean Abies species. PhD thesis, Swedish University of Agricultural Sciences, Umea. Parducci, L., Szmidt, A.E., 1999. PCR-RFLP analysis of cpDNA in the genus Abies. Theor. Appl. Genet. 98, 802–808. Parducci, L., Szmidt, A.E., Villani, F., Wang, X.R., Cherubini, M., 1996. Genetic variation of Abies alba in Italy. Hereditas 125, 11–18. Parducci, L., Szmidt, A.E., Madaghiele, A., Anzidei, M., Vendramin, G.G., 2001a. Genetic variation at chloroplast microsatellites (cpSSRs) in Abies nebrodensis (Lojac.) Mattei and three neighboring Abies species. Theor. Appl. Genet. 102, 733–740. Parducci, L., Szmidt, A.E., Ribeiro, M.M., Drouzas, A.D., 2001b. Taxonomic position and origin of the endemic Sicilian fir Abies nebrodensis (Lojac.) Mattei based on allozyme analysis. Forest Genet. 8, 119–128. Pascual, L., Garcia, F.J., Perfectti, F., 1993. Inheritance of Isozyme. Variations in seed tissues of Abies pinsapo Boiss. Silvae genetica 42, 335–335. Pfenninger, M., Véla, E., Jesse, R., Elejalde, M.A., Liberto, F., Magnin, F., Martínez-Ortí, A., 2010. Temporal speciation pattern in the western Mediterranean genus Tudorella P. Fischer, 1885 (Gastropoda, Pomatiidae) supports the Tyrrhenian vicariance hypothesis. Mol. Phylogenet. Evol. 54, 427–436. Pons, O., Petit, R.J., 1996. Measwring and testing genetic differentiation with ordered versus unordered alleles. Genetics 144, 1237–1245.

53

Pritchard, J.K., Stephens, M., Donnelly, P., 2000. Inference of population structure using multilocus genotype data. Genetics 155, 945–959. Provan, J., Soranzo, N., Wilson, N.J., Goldstein, D.B., Powell, W., 1999. A low mutation rate for chloroplast microsatellites. Genetics 153, 943–947. R Core Team, 2013. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria, URL http://www.Rproject.org/. Rivas-Carballo, M., Alonso-Gavilán, G., Valle, M., Civis, J., 1994. Miocene palynology of the central sector of the Duero basin (spain) in relation to palaeogeography and palaeoenvironment. Rev. Palaeobot. Palynol. 82, 251–264. Roiron, P., Ferrer, J., Linan, E., Rubio, C., Diez, J.B., Popescu, S., Suc, J.P., 1999. The flora of the Rubielos de Mora lacustrine basin: climatic conditions during the Lower Miocene in the Teruel region (Spain). Comptes Rendus de l’Academie des Sciences Series IIA Earth and Planetary Science 329, 897–904. Rosenbaum, G., Lister, G.S., 2004. Formation of arcuate orogenic belts in the western Mediterranean region. Orogenic curvature: integrating paleomagnetic and structural analyses, pp. 41–56. Rosenbaum, G., Lister, G.S., Duboz, C., 2002. Reconstruction of the tectonic evolution of the Western Mediterranean since the Oligocene. J. Virt. Explorer 8, 107–130. Scaltsoyiannes, A., Tsaktsira, M., Drouzas, A.D., 1999. Allozyme differentiation in the Mediterranean firs (Abies, Pinaceae). A first comparative study with phylogenetic implications. Plant Syst. Evol. 216, 289–307. Schönswetter, P., Tribsch, A., 2005. Vicariance and dispersal in the alpine perennial Bupleurum stellatum L. (Apiaceae). Taxon, 725–732. Se˛kiewicz, K., Se˛kiewicz, M., Jasin´ska, A.K., Boratyn´ska, K., Iszkuło, G., Romo, A., Boratyn´ski, A., 2013. Morphological diversity and structure of West Mediterranean Abies species. Plant Biosyst. Int. J. Dealing Aspects Plant Biol. 147, 125–134. Semerikova, S.A., Semerikov, V.L., 2014. Molecular phylogenetic analysis of the genus Abies (Pinaceae) based on the nucleotide sequence of chloroplast DNA. Russ. J. Genet. 50, 7–19. Semerikova, S.A., Semerikov, V.L., Lascoux, M., 2011. Post-glacial history and introgression in Abies (Pinaceae) species of the Russian Far East inferred from both nuclear and cytoplasmic markers. J. Biogeogr. 38, 326–340. Suyama, Y., Yoshimaru, H., Tsumura, Y., 2000. Molecular phylogenetic position of Japanese Abies (Pinaceae) based on chloroplast DNA Sequences. Mol. Phylogenet. Evol. 16, 271–277. Szpiech, Z.A., Jakobsson, M., Rosenberg, N.A., 2008. Adze: a rarefaction approach for counting alleles private to combinations of populations. Bioinformatics 24, 2498–2504. Taberlet, P., Fumagalli, L., Wust-Saucy, A.G., Cosson, J.F., 1998. Comparative phylogeography and postglacial colonization routes in Europe. Mol. Ecol. 7, 453–464. Tang, S., Dai, W., Li, M., Zhang, Y., Geng, Y., Wang, L., Zhong, Y., 2008. Genetic diversity of relictual and endangered plant Abies ziyuanensis (Pinaceae) revealed by AFLP and SSR markers. Genetica 133, 21–30. Terrab, A., Talavera, S., Arista, M., Paun, O., Stuessy, T.F., Tremetsberger, K., 2007. Genetic diversity and geographic structure at chloroplast microsatellites (cpSSRs) in endangered West Mediterranean firs (Abies spp., Pinaceae). Taxon 56, 409–416. Terrab, A., Schönswetter, P., Talavera, S., Vela, E., Stuessy, T.F., 2008. Range-wide phylogeography of Juniperus thurifera L., a presumptive keystone species of western Mediterranean vegetation during cold stages of the Pleistocene. Mol. Phylogenet. Evol. 48, 94–102. Tiba, B.B., Reille, M., 1982. Recherches pollenanalytiques dans les montagnes de Kroumirie (Tunisie septentrionale): premiers résultats. Ecol. Mediterranea 7, 75– 86. Vendramin, G.G., Degen, B., Petit, R.J., Anzidei, M., Madaghiele, A., Ziegenhagen, B., 1999. High level of variation at Abies alba chloroplast microsatellite loci in Europe. Mol. Ecol. 8, 1117–1126. Vos, P., Hogers, R., Bleeker, M., Reijans, M., van De Lee, T., Hornes, M., Friters, A., Pot, J., Kuiper, M., Zabeau, M., 1995. AFLP: a new technique for DNA fingerprinting. Nucl. Acids Res. 23, 4407–4414. Vendramin, G.G., Ziegenhagen, B., 1997. Characterisation and inheritance of polymorphic plastid microsatellites in Abies. Genome 40, 857–864. Wakasugi, T., Tsudzuki, J., Ito, S., Shibata, M., Sugiura, M., 1994. A physical map and clone bank of the black pine (Pinus thunbergii) chloroplast genome. Plant Mol. Biol. Report. 12, 227–241. Wolfe, K.H., Li, W.H., Sharp, P.M., 1987. Rates of nucleotide substitution vary greatly among plant mitochondrial, chloroplast, and nuclear DNAs. Proc. Natl. Acad. Sci. 84, 9054–9058. Xiang, Q.P., Xiang, Q.Y., Liston, A., Zhang, X.C., 2004. Phylogenetic relationships in Abies (Pinaceae): evidence from PCR-RFLP of the nuclear ribosomal DNA internal transcribed spacer region. Bot. J. Linn. Soc. 145, 425–435. Xiang, X., Cao, M., Zhou, Z., 2007. Fossil history and modern distribution of the genus Abies (Pinaceae). Front. Forest. China 2, 355–365. Xiang, Q.P., Xiang, Q.Y., Guo, Y.Y., Zhang, X.C., 2009. Phylogeny of Abies (Pinaceae) inferred from nrITS sequence data. Taxon 58, 141–152. Zangari, F., Cimmaruta, R., Nascetti, G., 2006. Genetic relationships of the western Mediterranean painted frogs based on allozymes and mitochondrial markers: evolutionary and taxonomic inferences (Amphibia, Anura, Discoglossidae). Biol. J. Linn. Soc. 87, 515–536. Ziegenhagen, B., Scholz, F., Madaghiele, A., Vendramin, G.G., 1998. Chloroplast microsatellites as markers for paternity analysis in Abies alba. Can. J. For. Res. 28, 317–321.