Patterns of phylogeography and vicariance of Chamaerops humilis L. (Palmae)

1132 http://journals.tubitak.gov.tr/botany/ Turkish Journal of Botany Turk J Bot (2014) 38: 1132-1146 © TÜBİTAK doi:10.3906/bot-1404-38 Patterns of phylogeography and vicariance of Chamaerops humilis L. (Palmae) Juan Luis GARCÍA-CASTAÑO1,2,*, Anass TERRAB1,2, María Ángeles ORTIZ1, Tod Falor STUESSY2, Salvador TALAVERA1 1Department of Plant Biology and Ecology, University of Seville, Seville, Spain 2Department of Systematic and Evolutionary Botany, Biodiversity Centre, University of Vienna, Vienna, Austria * Correspondence: [email protected] 1. Introduction Dramatic changes in the past 30 million years (Ma) in the West Mediterranean region have taken place in geology and climate, resulting in significant biogeographic impacts. The closing of the Tethys Sea, the isolation of the Mediterranean from the Atlantic Ocean, and the more recent Pleistocene glaciations are only 3 of the more prominent events that have affected plant distributions in this region (Thompson, 2005). It is not surprising, therefore, that a number of investigations have attempted to understand better the biogeographic and evolutionary implications in different species within this dynamic region (Ortiz et al., 2009). One of the most precise means of revealing the biogeographic history of plants within any region is the use of molecular markers at the population level. The spatial distribution of alleles for which phylogenetic relationships can be estimated is referred to as phylogeography (Avise, 2000). Markers often employed in these studies are haplotypes from the chloroplast genome as well as microsatellites (plastid and/or nuclear) and amplified fragment length polymorphisms (AFLPs). A number of studies have investigated phylogeographic patterns and processes in the West Mediterranean region. They have focused primarily on relationships established during the Tertiary in genera such as Erophaca (Casimiro- Soriguer et al., 2010), or Quercus and Pinus (Magri et al., 2007). Several studies have also examined environmental impacts and population changes during and after the Pleistocene, such as in Abies (Terrab et al., 2007), Hypochaeris (Ortiz et al., 2009), or Juniperus (Terrab et al., 2008). Another genus worth investigating, due to its interesting Euro-African distribution and rich fossil record, is Chamaerops L. (Palmae). This genus contains only a single species, Chamaerops humilis L., native to the West Mediterranean (Figure S1, online). It is 1 of 2 palms native to Europe and it is also the palm that occurs at the most northern latitudes (del Cañizo, 2011). This palm species was also reported during classical times (Gaius Plinius Secundus, better known as Pliny the Elder, 1st century AD), suggesting possible ancient traditional human use. Chamaerops belongs to tribe Trachycarpeae of subfamily Coryphoideae (Dransfield et al., 2008), and Abstract: Chamaerops humilis L. is 1 of 2 native palms occurring in Europe and the only native palm in the West Mediterranean region. Our aims were: (1) to describe its phylogeographic structure; (2) to infer a biogeographic scenario to explain its present distribution; and (3) to assess changes in its distribution from the last interglacial period. Twenty-two populations were sampled. An amplified fragment length polymorphism analysis produced 226 fragments, which allowed recognition of 4 groups of populations: (1) E Iberian Peninsula plus Al Hoceima (NE Morocco), the Balearics, and Sardinia; (2) France, the Italian Peninsula, and Sicily; (3) SW Iberian Peninsula and NW Morocco; and (4) S Morocco (var. argentea André; the first 3 groups are currently included in var. humilis). The phylogenetic information and molecular clock-related estimates, combined with geological and fossil history from the Eocene to present, suggest that C.humilis occurred in Central Europe in the Tertiary, reaching (1) Spain and (2) Italy, with expansion from here across North Africa to (3) N Morocco and S Spain, and (4) S Morocco. Climatic changes may also help explain the fragmented current distribution of this species. The groups of populations are sufficiently genetically distinct to recommend conservation of at least some populations in each region. Key words: Amplified fragment length polymorphisms, fossils, Mediterranean, palms, phylogeography, Quaternary, Tertiary Received: 08.04.2014 Accepted: 12.08.2014 Published Online: 17.11.2014 Printed: 28.11.2014 Research Article GARCÍA-CASTAÑO et al. / Turk J Bot 1133 molecular phylogenetic studies reveal that the closest generic relatives are Guihaia J.Dransf., S.K.Lee & F.N.Wei; Trachycarpus H.Wendl.; Rhapidophyllum H.Wendl. & Drude; Maxburretia Furtado; and Rhapis L.f. ex Aiton (Baker et al., 2009). It is small, multistemmed, and rhizomatous, with palmate leaves (do Amaral Franco, 1980; Herrera, 1989; del Cañizo, 2011) and of variable life form, leaf colour, and presence of thorns (Maire, 1980). The reproductive system of this dioecious species is interesting; it has leaf-mediated pollination, whereby the leaves release scents attractive to coleopterans (Herrera, 1989; Düfay and Anstett, 2004), and it has mammal-mediated dispersal of its fleshy fruits, mainly by foxes or badgers (Fedriani and Delibes, 2011). The present distribution of C. humilis is restricted to the West Mediterranean (Médail and Quézel, 1996), both on the European (where it is the only native palm species) and African sides (Figure1). In Europe, it occurs on the coasts of S Portugal, S and E Spain, SE France (rediscovered by Médail and Quézel, 1996), and W Italy. In Africa, it occurs mainly in Morocco but also in N Algeria and N Tunisia. In Libya it apparently has become extinct (Jafri and El-Gadi, 1977). The species has its widest distribution in Morocco, where 2 varieties have been distinguished: a form with green leaves in the north (var. humilis) and a form with glaucous leaves in the south (var. argentea André), which is often found at very high altitudes of up to 2000 m (Maire, 1980; Govaerts and Dransfield, 2005). The species also occurs on most of the large islands of the West Mediterranean, i.e. the Balearics, Sardinia, Sicily, and Malta (where it apparently has become extinct). In Corsica it does not seem to be autochthonous (do Amaral Franco, 1980). The oldest fossil records of Chamaerops are from the Eocene in Europe (Palamarev, 1989) and it may be considered as one of the most ancient fleshy-fruited species in the Mediterranean basin (Palamarev, 1989; Herrera, 1995). It appears, therefore, to have reduced its distribution from paratropical Eocene forests. Nowadays, with other species such as Laurus nobilis, it constitutes the Mediterranean ‘Laurocerasus-belt’, which is surrounded by species derived from more recent geological periods (Mai, 1989). The purposes of this paper, therefore, are: (1) to describe the phylogeographic structure of the species throughout its distributional range in the West Mediterranean region; (2) to infer a biogeographic scenario to explain its present distribution; and (3) to assess the most recent changes in the distribution of this species from the last interglacial period by use of niche modelling analyses. 2. Materials and methods 2.1. Amplified fragment length polymorphism 2.1.1. Sampling The 22 sampled populations are shown in Figure 1 and listed in Table 1; they cover most of the distributional area of the Figure 1. (A) Distribution map of Chamaerops humilis, populations selected for sampling and Structure results, and (B) colour-coded population assignment for the different Bayesian-inferred groups obtained for the different Ks set (black colour for not clearly assignable populations, NAPs; optimal value of K = 4, according to the Akaike and Schwarz criteria, is highlighted); see Section 2 and Tables 1 and 2. GARCÍA-CASTAÑO et al. / Turk J Bot 1134 Table 1. Population samples of Chamaerops humilis and Trachycarpus fortunei. Origin (population, country, locality/region, and coordinates), sample size (n), Fragtot, Fragpoly (%), number of phenotypes (nph), HD (mean and SE), and DW are shown. Population Country Locality/region Coordinates (m a.s.l.) n Fragtot Fragpoly (%) nph HD (SE) DW Continental Europe Portugal 1 Portugal Vila Real de Santo António (Algarve) 37°14′20.5″N, 7°26′16.0″W (2 m) 10 156 114 (72.1%) 10 0.272 (0.015) 268.8 Spain 2Spain Hinojos (Huelva) 37°17′04.0″N, 6°23′56.0″W (55 m) 10 161 98 (60.9%) 10 0.227 (0.014) 183.9 3Spain Pruna (Seville) 36°59′47.1″N, 5°12′33.3″W (698 m) 10 162 105 (64.8%) 10 0.254 (0.014) 296.4 4Spain Tarifa (Cádiz) 36°03′23.8″N, 5°32′55.6″W (319 m) 10 162 131 (80.9%) 10 0.296 (0.015) 233.9 5Spain Cabo de Gata (Almería) 36°43′45.4″N, 2°11′25.4″W (27 m) 10 141 100 (70.9%) 10 0.249 (0.015) 177.3 6Spain Peñón de Ifach (Alicante) 38°38′20.7″N, 0°04′27.2″E (17 m) 10 156 113 (72.4%) 10 0.265 (0.015) 187.0 7Spain Garraf (Barcelona) 41°15′13.3″N, 1°48′57.6″E (40 m) 10 153 106 (69.3%) 10 0.252 (0.015) 193.9 France 8 France Cap Taillat 43°10′13.1″N, 6°38′38.3″E (36 m) 10 158 98 (62.0%) 10 0.237 (0.015) 188.5 Italy 9Italy Circeo (Lazio) 41°13′54.6″N, 13°02′28.3″E (10 m) 10 155 88 (56.8%) 10 0.209 (0.015) 363.1 Mediterranean islands 10 Spain Cabo Formentor (Mallorca, the Balearics) 39°56′39.2″N, 3°10′16.0″E (114 m) 10 153 95 (62.1%) 7 0.224 (0.015) 253.8 11 Italy Castelsardo (Sardinia) 40°53′29.8″N, 8°39′15.6″E (78 m) 10 159 112 (70.4%) 10 0.258 (0.015) 221.6 12 Italy Palermo (N Sicily) 38°12′05.6″N, 13°16′03.9″E (17 m) 10 157 95 (60.5%) 10 0.237 (0.015) 303.5 13 Italy Syracuse (S Sicily) 36°59′59.6″N, 15°15′58.0″E (9 m) 10 156 95 (60.1%) 10 0.217 (0.014) 250.0 NW Africa 14 Morocco Cape Spartel 35°47′22.6″N, 5°53′33.4″W (223 m) 10 165 97 (58.8%) 10 0.238 (0.015) 242.3 15 Morocco Oued-Lau (Tétouan) 35°31′51.2″N, 5°11′55.1″W (37 m) 10 161 103 (64.0%) 10 0.246 (0.014) 238.9 16 Morocco Sidi Yahya du Rharb (the Mamora) 34°12′14.9″N, 6°16′58.0″W (100 m) 10 159 96 (60.4%) 7 0.229 (0.014) 251.0 17 Morocco Al Hoceima 35°13′16.1″N, 3°58′54.1″W (330 m) 10 163 117 (71.8%) 9 0.269 (0.014) 206.3 18 Morocco Tazzeka 34°09′11.9″N, 4°18′36.8″W (367 m) 10 162 126 (77.8%) 10 0.297 (0.014) 262.9 19 Morocco Beni Mellal 32°18′49.3″N, 6°16′16.0″W (1536 m) 10 158 92 (58.2%) 9 0.217 (0.014) 238.5 20 Morocco Essaouira 31°30′15.4″N, 9°41′43.6″W (100 m) 10 162 102 (63.0%) 9 0.247 (0.015) 250.3 21 Morocco Tizi-n-Test 30°48′09.2″N, 8°23′33.6″W (1110 m) 10 157 105 (66.9%) 9 0.229 (0.015) 237.6 22 Morocco Tafraoute 29°45′24.6″N, 8°50′03.3″W (1610 m) 10 149 96 (64.4%) 10 0.236 (0.014) 217.8 TOTAL Chamaerops humilis L. 220 208 201 (96.7%) 210 Trachycarpus fortunei (Hook.) H.Wendl. Garden Alamillo Park nursery (Seville, Spain) 37°25′22.1″N, 5°59′40.0″W (8 m) 10 106 61 (57.5%) 8 GARCÍA-CASTAÑO et al. / Turk J Bot 1135 species, including the largest islands. In every population we collected leaf material from 10 different stems that were placed far enough apart to avoid resampling the same genetic individual; the species is rhizomatous, which required caution and care in making the field collection. Material was dried and preserved in silica gel and stored prior to DNA extraction. Vouchers of sampled populations have been deposited in the herbarium of the University of Seville (SEV). Trachycarpus fortunei (Hook.) H.Wendl., a species belonging to the sister clade (Baker et al., 2009), was chosen as the outgroup representative (10 individuals were analysed as well). 2.1.2. DNA extraction and AFLP protocol Genomic DNA was extracted according to Doyle and Doyle (1987) following the CTAB protocol: after precipitation with isopropanol and centrifugation, the DNA pellet was washed with 70% ethanol and subsequently dried in a vacuum centrifuge. Finally, DNA was resuspended in TE buffer and treated with RNAse at 37 °C for 30 min. The AFLP procedure followed established protocols (Vos et al., 1995; Blears et al., 1998). Genomic DNA (~0.5 µg) was digested with 2 restriction endonucleases, EcoRI and MseI, and the resulting fragments were ligated to double-stranded EcoRI and MseI adaptors at 37°C for 2 h; the ligated DNA fragments were diluted with TE0.1 buffer. Preselective and selective amplifications were performed in a thermal cycler (Gene Amp PCR System 9700, PEApplied Biosystems), and the polymerase chain reaction (PCR) protocols also followed Vos et al. (1995). The number of fragments was reduced ~16-fold (4 × 4) through the use of preselective primers based on the sequences of EcoRI and MseI adaptors with the addition of nucleotides (EcoRI-A and MseI-C), matching the nucleotide downstream from the restriction sites. Preselective PCR products were diluted with TE0.1 buffer and then amplified with primers having 2 additional selective bases (3 in total), of which the first base was the same as that added in the previous amplification; this resulted in a further ~256-fold (16 × 16) reduction of the number of fragments, for ~4096-fold in total (16 × 256). Moreover, the EcoRI-based primers were labelled with fluorescent dyes. For the selective PCR, the following 3 primer combinations were used (selected, based on amplification effectiveness and variability, from a trial of 32): EcoRI (FAM)-ACT/MseI-CAG, EcoRI (VIC)- AGG/MseI-CAC, and EcoRI (NED)-ACC/MseI-CAC. The fluorescence labelled selective amplification products were separated on an automated sequencer (3130xl Genetic Analyzer, Applied Biosystems|Hitachi). Raw data were aligned with internal 500 ROX size standard using ABI PRISM GeneScan 3.7 (Applied Biosystems, 1989–2000). The resulting GeneScan files were imported into Genographer ver. 1.6.0 (Benham et al., 1999) and the fragments scored using the ‘thumbnail’ option, which allowed comparison of each fragment signal over all samples; the result of the scoring process was exported as a presence/absence (1/0) matrix. Genotyping and scoring errors are known to be very common in AFLP procedures and estimation of repeatability is recommended (Bonin et al., 2007). Repeatability was tested in 32 individuals (14.5%, 1–4 individuals per population) and was calculated as the ratio of the total number of mismatches (band presence vs. band absence) in the replicated individuals divided by the number of fragments found; this gave us an average value of 0.1 for the whole data set. Assuming that genotyping and scoring errors were randomly distributed and no bias at the population level would occur, we could follow the recommendation of 0.1-limit error (Bonin et al., 2007). Other authors such as Skrede et al. (2006) have also followed this criterion as being more than adequate to assess population-level variation. 2.1.3. Genetic diversity and population structure Within-population genetic diversity was assessed as the total number of AFLP fragments (Fragtot), the number and percentage of polymorphic fragments [Fragpoly (%)], and the number of different AFLP phenotypes (nph), all of them corrected with the 0.1 scoring error obtained (i.e. considering as monomorphic those fragments that were present in at least 90% of the individuals and ignoring those fragments present only at a maximum of 10% at a population level). We also calculated the average gene diversity index per population (HD = 1 – Σxi 2, where xi is the population frequency of ‘alleles’ 1 and 0 at locus i; the average gene diversity is the resulting mean across all loci). The ‘frequency-down-weighted marker values’ index (DW) was calculated as the ratio of means. These estimates were obtained using FAMD ver. 1.108β (Schlüter and Harris, 2006), Arlequin ver. 3.0.1 (Excoffier et al., 2005), and AFLPdat (Ehrich, 2006). FAMD ver. 1.108β was also used to gather individuals into population groups. Finally, the number of fragments exclusively shared between pairs of populations was assessed and also corrected with the scoring error obtained (i.e. at least present in more than 10% of the individuals of each population). Pairwise FST values based on Euclidean distances for the between-groups comparison were obtained using FAMD ver. 1.108β. The matrix was exported to SplitsTree ver. 4.5 (Huson and Bryant, 2006) where a Neighbour-Net was constructed to establish the different relationships between sampled populations. The support for each node was calculated separately with FAMD ver. 1.108β, using 1000bootstrap replicates based on chord distances for a neighbour-joining tree, by applying the Bayesian method to calculate null allele frequencies and with the amongpopulation nonuniform prior assumption. GARCÍA-CASTAÑO et al. / Turk J Bot 1136 Bayesian analyses using Structure ver. 2.2 (Pritchard et al., 2000; Falush et al., 2003) were applied, setting a burn-in period of 1 × 105, and run separately with lengths of 5 × 105 MCMC (5 iterations) and 1 × 105 MCMC (10 iterations) for K = 1–8. As very long distances are usually not traversed by pollen (Carrión, 2002), the simpler no-admixture model was applied, starting at prior population information, under the assumption of allele frequencies in each population being correlated. Average groups were calculated with CLUMPP ver. 1.1 (Jakobsson and Rosenberg, 2007) and graphically displayed using Distruct ver. 1.1 (Rosenberg, 2004). Apopulation was assigned to a specific group if it demonstrated at least 70% association. At this population level, the optimal K was selected according to the Akaike information criterion (AIC: Akaike, 1974; AICc and AICu: McQuarrie and Tsai, 1998) and the Bayesian or Schwarz information criterion (BIC: Schwarz, 1978). Additionally, principal coordinate analyses (PCoAs) based on Euclidean distances were applied to assess relationships among the individuals and populations (pairwise ΦST values), using FAMD ver. 1.108β. To test for isolation-by-distance, pairwise FST values for populations were compared with their geographical distance using Mantel bilateral tests based on Spearman correlations (1000 permutations; XLSTAT ver. 2007.6, Addinsoft). Finally, genetic differentiation among populations and different groupings were assessed by applying analyses of molecular variance (AMOVAs), based on Euclidean distances, using the FCT values obtained with Arlequin ver. 3.0.1. For groups of populations, the number of private shared fragments was calculated and a bio-neighbour- joining tree was constructed with SplitsTree ver. 4.5 based on ρST(2) values, which were used to compensate for unequal sample sizes (El Mousadik and Petit, 1996); the bootstrap support for each node was calculated in the same way as explained above. 2.1.4. Dating analysis The use of a molecular clock based on plastid (atpB, matK, rbcL, rps16 intron, trnL–trnF) or nuclear (18S, ITS, ms –malate synthase gene) sequences was rejected, because identities between these sequences, previously analysed by Baker et al. (2009), were higher than 92% (normally 99%– 100%) between C. humilis and T. fortunei. An attempt was made, therefore, to use a molecular clock based on AFLP data (Kropf et al., 2009). As has been pointed out by some authors (e.g., Ehrich et al., 2009), this clock is not free from analytical problems. DN72 genetic dissimilarity coefficient values were calculated with POPGENE ver. 1.32 (Yeh and Boyle, 1997). To check the feasibility of using this clock, we tested a range of FIS values (0–1), rates (fast and slow rate alternatives, i.e. 0.027 Ma/DN72 and 0.517 Ma/DN72), and generation time values, which were then calibrated at the Late Eocene. This period was chosen as fossil evidence showed that Chamaerops existed in Europe during the Tertiary, extending back as far as the Eocene. Fossils (stems, roots, seeds, leaves, or cuticles) assigned to C.humilis or very similar morphological species have been found from the Eocene (56–34 Ma) in Germany; the Lower Oligocene (34–28 Ma) in France; the Upper Oligocene (28–23Ma) of Bulgaria; the Miocene (23–5.3 Ma) in Poland, Moldavia, and Romania; the Lower Miocene (23–16 Ma) in Sardinia; the Upper Miocene (11– 5.3 Ma) in Greece; and the Pliocene (5.3–2.6 Ma) in Spain (Pons, 1981; Palamarev, 1989; Velitzelos and Gregor, 1990; Givulescu and Barbu, 1999; Palamarev et al., 2005; Kovar- Eder et al., 2006; Dransfield et al., 2008; Carrión, 2012; and references cited therein). To our knowledge, there are no fossil records from N Africa. 2.2. Niche modelling Niche modelling analyses were conducted for 3 temporal climatic frames: the last interglacial period (~120,000– 140,000 years BP; Otto-Bliesner et al., 2006), the last glacial maximum (~21,000 years BP; Braconnot et al., 2007), and the present (~1950–2000; Hijmans et al., 2005). Data for the last glacial maximum were obtained from 2 atmospheric models, from the Community Climate System Model (CCSM) and from the Model for Interdisciplinary Research on Climate (MIROC). Bioclimatic data were obtained from the WorldClim Global Climate Data web page (http:// www.wordclim.org) and all 19 WorldClim variables were used. Soil data (FAO/IIASA/ISRIC/ISSCAS/JRC, 2012) were obtained from the Harmonized World Soil Database (http://webarchive.iiasa.ac.at/Research/LUC/External- World-soil-database/HWSD_Data/HWSD_RASTER.zip), with selection of the most precise mapping unit. Data sets were rescaled to 2.5 arc-minutes with Quantum GIS ver. 1.8.0-Lisboa in order to make comparable the different scenarios. Present location data were obtained from vouchers deposited in SEV, which included our extensive sampling; additionally, we added locations from Algeria (Royal Botanic Gardens Herbarium, Kew: K000208642, and K000208639, K000208640, K000208641, K000208643, and K000208644) obtained from the Global Biodiversity Information Facility (GBIF) data set (accessed 19 June 2013). Finally, predictive analyses, based on the whole set of variables, were performed using the maximum entropy method with the software MaxEnt (Phillips et al., 2006) and results were shown for 2 thresholds: maximum training sensitivity plus specificity logistic threshold and 10-percentile training presence logistic threshold. 3. Results 3.1. AFLP profile A total of 226 fragments were obtained ranging from 60 to 447 bp (208 in Chamaerops), 4 of them being GARCÍA-CASTAÑO et al. / Turk J Bot 1137 monomorphic; 7 fragments were present in 90% of the Chamaerops individuals (for further details, see Table 1). For different primers the results were as follows: (1) EcoRI (NED)-ACC/MseI-CAC, 51 fragments, 60–405 bp, 1 fragment monomorphic in Chamaerops; (2) EcoRI (VIC)-AGG/MseI-CAC, 126 fragments, 63–447 bp, 4 fragments monomorphic in Chamaerops; (3) EcoRI (FAM)-ACT/MseI-CAG, 49 fragments, 65–395 bp, 2 fragments monomorphic in Chamaerops. At the 90% level, 18 fragments were exclusive to Trachycarpus fortunei and 120 to Chamaerops. AFLPdat (after including 0.1 scoring error correction) showed the existence of 218 different phenotypes out of the 230 samples analysed, with 7–10 phenotypes per population. 3.2. Genetic diversity The results of the estimated genetic diversity patterns are shown in Table 1. The number of total fragments varied from 141 in Cabo de Gata (Spain) to 165 in Cape Spartel (Morocco); the percentage of polymorphic fragments ranged from 56.8% in Circeo (Italy) to 80.9% in Tarifa (Spain); the average gene diversity (HD) ranged from 0.209 in Circeo (Italy) to 0.296 in Tarifa (Spain) and 0.297 in Tazzeka (Morocco); and the rare fragments index (DW) varied from 177.3 in Cabo de Gata (Spain) to 363.1 in Circeo (Italy). 3.3. Population structure The Bayesian analysis conducted with Structure, in conjunction with the Akaike and Schwarz criteria used to select the optimum number of clusters, pointed to the existence of 4 main groups (see Table 2 and Figure 1). These were defined from K 2 to 4 in the following order: (1) E Spain (Cabo de Gata, Peñón de Ifach, and Garraf; pops. 5, 6, and 7), the Balearics (pop. 10), Sardinia (pop. 11), and Al Hoceima (pop. 17; NE Morocco), hereafter labelled as the E Spain-Sardinia-NE Morocco group of populations; (2) France (pop. 8), continental Italy (pop. 9), and Sicily (Palermo and Syracuse; pops. 12 and 13), hereafter labelled as the France-Italy group of populations; (3) SW Iberian Peninsula (Portugal: Vila Real de Santo António; pop. 1; and Spain: Hinojos and Pruna; pops. 2 and 3), and NW Morocco (Cape Spartel and Tétouan; pops. 14 and 15), hereafter labelled as the S Spain-N Morocco group of populations; and (4) S Morocco: Beni Mellal (pop. 19), Essaouira (pop. 20), Tizi-n-Test (pop. 21), and Tafraoute (pop. 22), hereafter labelled as the SMorocco group of populations. Three populations showed an intermediate position between these main groups: Tarifa, pop. 4 (50% of groups 1 and 3), the Mamora, pop. 16 (54% of group 3 and 46% of group 4), and Tazzeka, pop. 18 (52% of group 3, 40% of group 1, and 8% of group 4), named nonassignable populations (NAPs) hereafter. The Neighbour-Net obtained with the FST pairwise distances (Figure 2) is consistent with the grouping calculated with Bayesian analysis (Table 2 and Figure 1). Bootstrap values obtained were 91%BS (E Spain-Sardinia- NE Morocco group of populations), 98% BS (France- Italy group of populations), 60% BS (S Spain-NMorocco group of populations), and 75% BS (SMorocco group of populations). Tarifa, Tazzeka, and the Mamora populations also occupied an intermediate position among groups, as in the Bayesian analysis (Table 2; Figures 1 and 2). In Figure 3A we present PCoA results for principal component (PCo) 1 (10.2%) and PCo 2 (5.5%) for the individuals analysed. The Mamora individuals appear at an intermediate position between the S Morocco and S Spain-N Morocco groups of populations, and those Figure 2. Neighbour-Net for the populations of Chamaerops humilis sampled, including Trachycarpus fortunei as an outgroup representative, based on the Euclidean distances. Bootstrap values obtained on the chord distances for the neighbour-joining tree; see Table 1 for codes. GARCÍA-CASTAÑO et al. / Turk J Bot 1138 from Tarifa and Tazzeka are dispersed in the centre of the graph. These patterns are clearer when the average values are shown (Figure 3B). If the PCoA is based on the population, not the individual, distance matrix, the percentage of variance explained is much higher (PCo 1 = 38.6% and PCo 2 = 26.9%), but the relative position of the populations remains almost the same as in Figure 3B (results not shown). Table 2. Akaike and Schwarz criteria values for the different Bayesian-inferred groups obtained for the different Ks set (optimal value of K= 4 is highlighted). K Group SS (d.f.) Pop. SS (d.f.) Ind. SS (d.f.) Total SS (d.f.) AIC AICc AICu BIC 1 0.000 (0) 2294.209 (21) 5495.000 (198) 7789.209 (219) 1409.04 4.59 4.60 784.71 2 595.009 (2) 1699.200 (19) 5495.000 (198) 7789.209 (219) 1393.56 4.52 4.53 772.62 3 937.489 (3) 1356.720 (18) 5495.000 (198) 7789.209 (219) 1384.83 4.48 4.49 767.29 4 1196.307 (4) 1097.902 (17) 5495.000 (198) 7789.209 (219) 1378.36 4.45 4.47 764.21 5 1072.849 (4) 1221.360 (17) 5495.000 (198) 7789.209 (219) 1384.44 4.48 4.50 773.68 6 1199.046 (4) 1095.164 (17) 5495.000 (198) 7789.209 (219) 1382.27 4.47 4.50 774.90 7 1068.203 (4) 1226.006 (17) 5495.000 (198) 7789.209 (219) 1388.59 4.50 4.53 784.62 8 875.551 (4) 1418.658 (17) 5495.000 (198) 7789.209 (219) 1396.81 4.53 4.58 796.23 Figure 3. PCoA for the first 2 axes, showing the different Bayesian-inferred groups of populations as well as the not clearly assignable populations (NAPs) for (A)individuals and (B) their mean values per population; see Table 1 for codes. GARCÍA-CASTAÑO et al. / Turk J Bot 1139 The Mantel test applied between the geographic and the genetic distances was positively significant (Spearman’s ρ = 0.629, P = 0.001). Moreover, the number of exclusive shared fragments between pairs of populations, which reveal very closely related populations, was null or very low, with only 5 matches of 1 fragment between: pop. 13 and pop. 15, 13 and 22, 18 and 21, 10 and 16, and 19 and 20. The different AMOVAs, considering no grouping, i.e. globally, and different groupings (excluding the 3 NAPs), showed that the significantly best grouping was the same as that defined by Bayesian inference, i.e. E Spain- Sardinia-NE Morocco, France-Italy, S Spain-NMorocco, and S Morocco (Table 3). Of the ρST(2) values calculated for the groups of population assortment, S Spain-N Morocco and S Morocco had the closest relationship (ρST(2) = 0.083), the latter showing the largest separation from France-Italy (ρST(2) = 0.128) and E Spain-Sardinia-NE Morocco (ρST(2)=0.125). In the bio-neighbour-joining tree (Figure 4), E Spain- Sardinia-NE Morocco appeared as the basally branching group in the Chamaerops core, followed by France-Italy and, finally, by S Morocco and S Spain-N Morocco. Table 3. Results of the AMOVAs, considering different groupings and excluding not clearly assignable populations (NAPs: pops. 4, 16, and 18) according to Bayesian-inferred groups for K = 4. The best (lowest) Akaike and Schwarz criteria result is highlighted; see Section 2 for details. FCT (95% CI) Group % (d.f.) Pop. % (d.f.) Ind. % (d.f.) AIC AICc AICu BIC No grouping -0.00 (0) 24.52 (18) 75.48 (171) 1377.31 4.44 4.45 752.98 var. humilis (pops. 1–18)/ var. argentea (pops. 19–22) 0.103 (0.079–0.130) 10.35 (1) 19.24 (17) 70.41 (171) 1367.42 4.40 4.41 746.48 N Mediterranean, islands included (pops. 1–13)/ S Mediterranean (pops. 14–22) 0.044 (0.031–0.057) 4.36 (1) 21.83 (17) 73.81 (171) 1370.99 4.41 4.43 750.05 Iberia (pops. 1–7)/France-Italy (pops. 8–9)/ Islands (pops. 10–13)/Morocco (pops. 14–22) 0.046 (0.033–0.060) 4.63 (3) 20.76 (15) 74.61 (171) 1365.05 4.39 4.41 750.90 Iberia (pops. 1–7)/France-Italy (pops. 8–9)/ Islands (pops. 10–13)/N Morocco (pops. 14–17)/ S Morocco (pops. 19–22) 0.072 (0.054–0.091) 7.19 (4) 18.31 (14) 74.50 (171) 1358.84 4.36 4.39 748.08 S Spain-N Morocco (pops. 1–3 and 14–15)/ E Spain-Sardinia-NE Morocco (pops. 5–7,10–11 and 17)/ France-Italy (pops. 8–9 and 12–13)/ S Morocco (pops. 19–22) 0.173 (0.150–0.196) 17.28 (3) 10.05 (15) 72.66 (171) 1343.77 4.29 4.31 729.62 Figure 4. Bio-neighbour-joining tree for the Bayesian-inferred groups of populations (excluding NAPs) of Chamaerops humilis sampled, including Trachycarpus fortunei as an outgroup representative, based on the ρST(2) distances. Bootstrap values obtained on the chord distances for the neighbour-joining tree. Divergence times obtained from an AFLP-based molecular clock (Kropf et al., 2009). GARCÍA-CASTAÑO et al. / Turk J Bot 1140 3.4. Dating analysis Plausible generation times arose when FIS and rate values were the highest among those proposed, i.e. FIS = 1 and 0.571 Ma/DN72; smaller values resulted in higher generation times, of more than 1000 years at slow rates. When calibrating the first divergence between T. fortunei and C.humilis in the Late Eocene (at least 34 Ma), which resulted in a generation time of 309.45 years, the remaining nodes were then set to 8.32 Ma, 7.58 Ma, and 5.83 Ma for the E Spain-Sardinia-NE Morocco, France-Italy, and SMorocco and S Spain-N Morocco groups of populations respectively (Figure 4), i.e. in the Miocene. 3.5. Niche modelling Results from niche modelling revealed a relatively narrow distribution of C. humilis during the last interglacial (Figure 5A), having been even narrower in the northern shore of the Mediterranean in the last glacial maximum (Figures 5B and 5C), although the degree of this narrowing was highly dependent on the atmospheric model chosen. The 4 groups of populations currently recognised also showed congruent distributional areas during this period. For the S Spain-N Morocco and the S Morocco groups, range alteration was not important (increasing even the potential range of this species in N Africa), but for the continental areas of the France-Italy and E Spain- Sardinia-NE Morocco groups, reduction might have been significant. For the islands sampled, there were extensions of suitable areas during the 3 time frames, and, remarkably, Corsica showed a very low likelihood of occurrence. There was probability of occurrence in Crete, as well as in some areas of Greece, Turkey, or the Levant, with the last ranges of suitable areas suffering minor alteration during the last glacial maximum (Figures 5A–5D). 4. Discussion To infer the biogeographic processes that might explain the present distribution and patterns of genetic diversity among populations of Chamaerops humilis, several factors need to be taken into account. The phylogenetic context is important, because relationships to related genera might offer insights on location of origin. Likewise, fossil evidence is significant, because of the unusual morphology of this palm and relative abundance of dated fossil material. The historical climatic and geological changes in the West Mediterranean region, one of the most dynamic areas of the world, must also be taken into account. Finally, levels of genetic diversity within and among groups of populations as revealed in the phylogeographic analysis might suggest regions of ancestral refugia and/or recent migrations. These different sources of data, therefore, have been considered in the development of the biogeographic scenarios presented herein. 4.1. Outgroup relationships and fossil evidence To infer the ancient distribution of Chamaerops humilis, it is helpful to understand the phylogenetic relationships between it and close generic relatives. Chamaerops belongs A C B D Figure 5. Suitable areas for Chamaerops humilis estimated in 3 periods: last interglacial (LIG, ~120,000–140,000 years BP; A); last glacial maximum (LGM, ~21,000 years BP), B for CCSM and C for MIROC atmospheric models; and present period (D, ~1950–2000). Results shown were obtained from niche modelling analyses considering both bioclim and soil data, for 2 thresholds: maximum training sensitivity plus specificity logistic threshold (in dark grey) and 10-percentile training presence logistic threshold (in black). LIG and LGM are represented with current coastlines. A C B D A C B D A C B D GARCÍA-CASTAÑO et al. / Turk J Bot 1 Figure S1. Chamaerops humilis individual in a typical Mediterranean Quercus L. forest from pop. 3 (Pruna, Spain). Supplementary material