Molecular genetic analysis of two native desert palm genera, Washingtonia and Brahea, from the Baja California Peninsula and Guadalupe Island

Abstract The complex geological and ecological processes that have generated high levels of biodiversity and endemism in the Baja California Peninsula have been the subject of intensive study. However, relatively little is known about phylogeography of the iconic endemic palm species of this region. We therefore analyzed a total of 2,294 bp of chloroplast and 738 bp of nuclear sequence data in 169 samples of five native palm species from Baja California, Sonora and Guadalupe Island. We found that Washingtonia and Brahea palms had low levels of genetic diversity and were highly structured, with the majority of species and major geographic regions being characterized by distinct haplotypes. We also found strong support for currently recognized species in Washingtonia, but our results were less clear cut for Brahea due to haplotype sharing. Furthermore, patterns of population structure were broadly consistent with historical vicariant events such as the inundation of the Isthmus of La Paz, the formation of the Sea of Cortez, and the more recent colonization and isolation of Guadalupe Island's palms. Our findings contribute toward a growing appreciation of the complexity of plant responses to past geological changes and also provide valuable baseline genetic data on relict American palm species.


| INTRODUCTION
Spatial patterns of species abundance and distribution are of a fundamental importance as they contribute toward virtually every facet of ecological knowledge from conservation to the understanding of ecosystem function and structure. Current species distributions arise from the interplay of various evolutionary, ecological and geological forces as well as from stochastic effects acting on species and populations. Multiple factors have been implicated in shaping species distributions including contemporary climate conditions (Currie et al., 2004;Svenning & Skov, 2007), Quaternary climatic oscillations (Hewitt, 1996(Hewitt, , 2000(Hewitt, , 2004, topography and habitat heterogeneity (Kerr & Packer, 1997), dispersal capacity (Svenning, Normand, & Skov, 2008), and biotic interactions (Araújo & Luoto, 2007;Kissling, Rahbek, & Bohning-Gaese, 2007).
The biological diversity of the Baja California Peninsula has long been recognized as unique and significant (Grismer, 2000;Lindell, Ngo, & Murphy, 2006;Riemann & Ezcurra, 2005;Wiggins, 1980). The Baja California is also considered a "natural laboratory" in which multiple geological and ecological events are thought to have triggered cycles of vicariant episodes followed by dispersion along climatic gradients (Mantooth, Hafner, Bryson, & Riddle, 2013;Riddle, Hafner, Alexander, & Jaeger, 2000). The major geological processes that have been evoked in explaining patterns of biotic diversification in the region have been classified on the basis of their timing and strength (Dolby, Bennett, Lira-Noriega, Wilder, & Munguia-Vega, 2015). The continental rifting of the Peninsula away from mainland Mexico and the volcanic upwelling of Guadalupe Island are the most ancient (>5 million years ago, MYA) and resulted in ecological speciation, island endemism and disjunctions in mainland-peninsular and peninsular-island sister species distributions (Aleixandre, Hernandez Montoya, & Mila, 2013;Grismer, 2000;León De La Luz, Rebman, & Oberbauer, 2003;Riddle et al., 2000;Rosas-Escobar, Gernandt, Pinero, & Garcillan, 2011). More recently (1-3 MYA) and more locally, complex physical interactions between the land and the sea, such as the inundation of the Isthmus of La Paz and the formation of the temporal mid-peninsular seaway, resulted in a north-south genetic discontinuity within the peninsula and the generation of regions of high local endemism (Dolby et al., 2015; León-De la Luz & Breceda, 2006;McCauley, Cortés-Palomec, & Oyama, 2010;Riddle et al., 2000). This complex interplay was followed by progressive aridification of the Baja Peninsula and adjacent areas after the Last Glacial Maximum (LGM) (Hafner & Riddle, 1997;Lindell et al., 2006;Riddle et al., 2000). Following drastic changes in precipitation, a number of endemic species are thought to have become either locally extinct or restricted to isolated permanent water bodies such as oases and canyons (Grismer & McGuire, 1993;Wehncke, López-Medellín, & Ezcurra, 2010).
Palms (Arecaceae) are one of the most distinctive of all plant families and are emblematic of the tropics (Tregear, Rival, & Pintaud, 2011). Although the majority of palm species are restricted to the wet tropics, several groups occur in dry open savannas and desert oases (Tomlinson, 2006). Palm populations on the Baja California Peninsula are considered to be relicts of historically more widespread and continuous populations that are now largely confined to sites where permanent water exists either above or below the ground (Axelrod, 1950;Cornett, 1985;Cornett, Glen, & Stewart, 1986;Felger & Joyal, 1999;Grismer & McGuire, 1993). This is evident in the fossil record, which shows that during the late Cretaceous, palms were common across North America and extended much further than their current geographic distribution (Couvreur, Forest, & Baker, 2011;Harley, 2006). Two native palm genera (Washingtonia spp. and Brahea spp.) are currently present on the Baja Peninsula (Minnich, Franco-Vizcaino, & Salazar-Ceseña, 2011). The genus Washingtonia is represented by two species, W. robusta and W. filifera, but the taxonomical relationships within this genus are poorly resolved (Felger & Broyles, 2007;Felger & Joyal, 1999;Henderson, Galeano, & Bernal, 1995;McClintock, 1978).

W. filifera is distributed from southeastern California and western
Arizona to northern Baja California (McClintock, 1978). By contrast, W. robusta is locally distributed near Cataviña and Sierra Asamblea and is more common in the southern half of the Baja peninsula (Minnich et al., 2011). It also occurs on the Mexican mainland, where it is restricted to a few riparian canyons at the southern edge of the Sonoran Desert and appears to be relictual, being geographically and ecologically the most narrowly distributed palm species in the region (Felger & Joyal, 1999). There is still unresolved debate over the relictual status of Washingtonia palms, with some authors having suggested that their current distribution is of recent invasive origin (Cornett, 1987b(Cornett, , 2008 or might be linked to the movement of native human populations (Hicks, 1989;Minnich et al., 2011). One the other hand, many characteristics of Washingtonia palms such as their patchy geographical distributions, low numbers, and affinity to tropical climates would support fossil evidence in suggesting that their current distribution may be relictual (Axelrod, 1950;Felger & Joyal, 1999).
The Brahea complex comprises nine species, two of which are endemic and restricted to Baja California (B. brandegeei and B. armata) and one (B. edulis) to Guadalupe Island. The remaining Brahea species are absent from the Baja California Peninsula but can be found on the Mexican mainland, with some species (B. nitida and B. dulcis) extending their range into Central America. Brahea is the least studied genus of American palms and the relationships between and within its species are not clearly described (Henderson et al., 1995). The geographical range of B. brandegeei extends from the northern Baja California Sur (Sierra San Francisco) to Sierra La Laguna at the southernmost tip of the peninsula (Felger, Johnson, & Wilson, 2001;Minnich et al., 2011).
B. armata extends from just south of the United States-Mexico border southward to just north of the state line of Baja California Sur in the central peninsula (Franco-Vizcaino, Lopez-Beltran, & Salazar-Cesena, 2007;Wiggins, 1980). Nonetheless, the exact limits of the ranges of B. armata and B. brandegeei are unclear because of taxonomic uncertainties (Felger & Joyal, 1999;Felger et al., 2001;Henderson et al., 1995). The Guadalupe fan palm or B. edulis is distributed across the northwestern side of Guadalupe Island (Garcillán, Vega, & Martorell, 2012;León De La Luz et al., 2003) with small patches found along the bottoms of arroyos throughout the island (Oberbauer, 2005).
In spite of the ecological, ethno-biological and economical importance of Washingtonia and Brahea palms, many aspects of their biology including species delimitations, genetic relationships between species and levels of intraspecific genetic diversity remain to be elucidated.
Previous studies have examined certain aspects of the fundamental ecology of these species such as growth rates, mechanisms of seed dispersal and herbivore-plant interactions (Bullock & Heath, 2006;Franco-Vizcaino et al., 2007;Wehncke & López-Medellín, 2015;Wehncke, López-Medellín, Wall, & Ezcurra, 2013;Wehncke et al., 2010). The geographic distributions of certain species based on morphological characteristics such as stem size, leaf shape and color, fruit shape and size, have also been described (Cornett 2008;Felger & Joyal, 1999;Oberbauer, 2005;Minnich et al., 2011). However, because morphological differences can sometimes be ambiguous and the majority of Washingtonia and Brahea populations on the peninsula are extremely hard to access (Felger & Joyal, 1999) no single study has effectively elucidated the distribution of all of the palms found in the Baja California Peninsula and Guadalupe Island.
In recent decades, molecular genetic studies have proven instrumental in resolving ongoing debates over the biogeographical origins of diverse taxa, including the Aracaceae (Baker & Couvreur, 2013;Couvreur et al., 2011;Kissling, Eiserhardt, et al., 2012;Kondo et al., 2012). However, the majority of studies of palms have focused on commercially important species (Barrett, Bacon, Antonelli, Cano, & Hofmann, 2016) and we are not aware of any published genetic studies of Washingtonia or Brahea species, with the exception of a single allozyme study of W. filifera in California (McClenaghan & Beauchamp, 1986).
Here we used Washingtonia and Brahea palms as a case study to explore how tropical species may have responded to the complex geological and ecological history of the mostly arid Baja California Peninsula and adjacent Guadalupe Island. These genera share similar life histories and habitat requirements yet both occupy the northern limit of the Arecaceae distributional range, are endemic and locally rare. To provide the first baseline genetic data on Washingtonia and Brahea palms from this region, we used a classical Sanger sequencing approach. By combining chloroplast and nuclear sequence data, we first set out to determine species delimitations on the Baja Peninsula.
We then characterized patterns of haplotype diversity and distribution to address a number of questions: (i) Do sympatrically occurring

| Sampling
A total of 176 samples were collected from a total of five palm species (Washingtonia and Brahea) identified on the basis of their morphology in the field. Specimens were collected from virtually all accessible oases throughout the Baja California Peninsula as well as from Guadalupe Island and two sites in Sonora, Mexico ( Figure 1). We avoided sampling immediately adjacent trees in order not to include close relatives in our dataset. Although our total sample size is modest, it is important to note that the species in question are rare and endemic, and normally grow in inaccessible canyons. Many sites could only be reached by foot or by riding a mule for several hours and in some cases days. Additionally, some stands were very small and were represented by only 5-10 individuals.  Henderson et al., 1995) have argued that B. brandegeei may be synonymous to B. elegans from the Mexican mainland, we also included one herbarium specimen of this species.
More details on sampled species, sierras and sites, and the number of amplified chloroplast and nuclear sequences can be found in Figure 1, Tables 1 and S1. The genetic data were generated from silica-gel dried leaf samples collected in the field during 2015 and 2016.
Five loci (one nuclear and four chloroplast) were deemed suitable for further study and these were sequenced in all of the samples using the PCR conditions described above. The presence of amplified target DNA fragments was verified on 1.5% agarose gels stained with ethidium bromide. All positive PCR products were sent to Arizona GeneCore for sequencing in both directions. The resulting sequence chromatograms were assembled and manually edited using BioEdit 7.0.0 (Hall, 1999). The chloroplast sequences were then concatenated into a single sequence, as were the nuclear sequences, allowing subsequent analyses to be carried out separately for these two types of loci.

| Molecular data analysis
Our sample size was somewhat restricted due to the difficulty of collecting samples from these locally rare and geographically restricted endemic species. We therefore tested if our final sample size was sufficient to recover most of the genetic diversity within the sampled populations. For this analysis, we constructed haplotype accumulation curves separately for the nuclear and chloroplast sequence data using the function haploAccum in the R package SPIDER (R Development Core Team 2008) with a total of 10,000 permutations. Molecular diversity statistics including the number of variable sites, the number of haplotypes, nucleotide and haplotype diversities, and the number of shared haplotypes were then computed using DnaSP 5.1 (Librado & Rozas, 2009) and ARLEQUIN 3.5.1.3 (Excoffier & Lischer, 2010).
To assess levels of genetic divergence between putative taxonomic species and a priori defined geographical regions ( Figure 1, Table 1 and see below), we quantified pairwise differentiation using F st (Weir & Cockerham, 1984) and Φ st (Kimura, 1980). The former looks only at haplotype frequency differences (Excoffier, Smouse, & Quattro, 1992) while the latter also incorporates haplotype sequence similarity. We also used hierarchical analysis of molecular variance (AMOVA) to partition the genetic variation within Washingtonia and Brahea using ARLEQUIN. We used AMOVA to estimate relative support for three alternative scenarios ( Figure 2). Each scenario specified a hierarchical partitioning of populations within regions/groups, and the nesting scheme that maximized the among-group variance component was considered the best divergence scenario. The first scenario was based on the current taxonomical delimitation of the species (i.e. two groups were specified for Washingtonia, corresponding to W. filifera and W. robusta and three groups were specified for Brahea corresponding to we tested the hypothesis that relict populations persisted within each of the sierras, which should be reflected in each sierra harboring a genetically distinctive population. Eight groups were specified for Washingtonia and ten for Brahea, each corresponding to a different sierra (For more details see Figure 2). AMOVA and genetic differentiation analyses were not carried out for the Brahea chloroplast data, as only three haplotypes were present and almost all of the individuals shared the same haplotype. The single B. elegans sample was also excluded from these analyses. For each AMOVA, statistical significance was determined using 10,000 permutations of the dataset.
Finally, we used a Bayesian approach to estimate the number of genetic clusters present within the data as implemented in BAPS 6 (Corander, Marttinen, Siren, & Tang, 2008;Corander, Siren, & Arjas, 2008;Corander, Waldmann, & Sillanpaa, 2003). For this analysis, individuals were first clustered without any previous knowledge of their sampling locations. This analysis was run using the mixture model to determine the most probable number of clusters (K) within the data.
K was set from 1 to 10 and five replicates were performed for each value of K. An admixture analysis was then carried out with 50 reference individuals and 1,000 iterations. Finally, we spatially clustered the individuals using Voronoi tessellation and geo-referenced sample location data. All of the spatial analyses were conducted separately for Washingtonia and Brahea.

| Ecological niche modeling
To better understand the current and past population dynamics of Washingtonia and Brahea palms along the Baja Peninsula we performed Ecological Niche Modeling (ENM) using the maximum entropy method implemented in Maxent 3.4.0 (Phillips, 2017;Phillips, Anderson, & Schapire, 2006). Maxent attempts to estimate a probability distribution of species occurrence that is closest to uniform while still subject to environmental constraints (Elith et al., 2011).
Information on the geographic distribution of Washingtonia and Brahea was based on a total of 49 and 44 unique presence records respectively. We combined our field data, information of natural occurrence of Washingtonia in the USA described in Cornett et al. (1986) and georeferenced occurrence data from www.gbif.org and www.bajaflora. org (Accessed March 2017). Each model included seven bioclimatic T A B L E 1 Numbers of Washingtonia and Brahea samples collected and successfully sequenced at nuclear and plastid loci. Also shown are the number of population sampled for each species and sampling sites are included for each population and species. The sampled populations (sierras) correspond to those shown in Figure 1 Species Sierra (Population)  variables that we considered might have the strongest effect on the distribution of these species (mean temperature of wettest quarter, mean temperature of driest quarter, mean temperature of warmest quarter, mean temperature of coldest quarter, precipitation of wettest quarter, precipitation of driest quarter, and precipitation of coldest quarter). By restricting this analysis to a subset of the most promising variables, we avoided potential overfitting of the Maxent model (Peterson & Nakazawa, 2008). To gain a better perspective on temporal changes in the distributions of Washingotnia and Brahea palms, we predicted suitable climate envelopes for each genus using current ecological conditions ) and the following three historical time periods: mid Holocene (~6,000 years BP), the last glacial maximum (~22,000 years BP) and the last interglacial (~130,000 years BP).

Number of collected samples per population
All of the bioclimatic data were downloaded from WorldClim (www. worldclim.org) as a set of raster layers with a spatial resolution of 30 s  (Swets, 1988). Jackknife analyses were performed to determine the relative contribution of each of environmental variable to the model of projected current distribution (Phillips, 2017).
In order to test the null hypothesis that the two Baja palm's species occupy identical climatic environments ("niches"), we performed a niche identity test using ENMtools 1.4.4 (Warren, Glor, & Turelli, 2010) with 1,000 pseudoreplicates. We then compared a distribution of niche similarities obtained from pairs of pseudoniches based on randomly sampled occurrence points with the actual niche overlap between Washingtonia and Brahea (Warren et al., 2010). The latter was quantified using Schoener's D (Schoener, 1968) and the standardized Hellinger distance (I) (Warren, Glor, & Turelli, 2008). Schoener's D assumes that the suitability scores are proportional to species abundance, whereas Hellinger's-based I quantifies the probability distributions of two ecological niche models. Both similarity metrics range from 0 (no niche overlap) to 1 (identical niches) (Warren et al., 2010).

| Molecular markers
We evaluated a total of ten chloroplast and nuclear loci that have been proposed for molecular studies of the Arecaceae (Bacon et al., 2008;CBOL Plant Working Group 2009;Jeanson et al., 2011). Of these, we selected five informative loci: matK, CISP 4, trnT-trnD, trnG-trnS, and trnF-trnL for further analysis. The other loci either lacked sequence variation in a representative sample of individuals (rbcL and trnC-rpoB), failed to PCR amplify (trnfM-trnS) or amplified more than one PCR product (ITS2 and CISP 8).  Figures 3 and   4). Levels of sequence diversity were generally low within the morphologically defined species, with individuals carrying either a single common haplotype or a handful of closely related haplotypes at low frequency ( Table 2). As shown by the haplotype accumulation curves, our sample size was sufficiently large to encompass the majority of the genetic diversity present within the studied species ( Figure S1).

| Relationships among and within species
Parsimony networks were used to visualize relationships among the samples based on the concatenated chloroplast and nuclear sequences (Figures 3a and 4a  whereas another haplotype was only found in samples from the SL population (Figure 3c). Correspondingly greater haplotype diversity was found at the nuclear sequence data, with a total of 11 haplotypes detected, although most of these were present at low frequencies.
There was also considerably more haplotype sharing among the Brahea species, with only B. edulis being distinguished by a unique nuclear haplotype (Figure 4a). Spatial distributions of the nuclear haplotypes for Brahea did not reveal any clear geographical patterns along the Peninsula, but confirmed the uniqueness of Guadalupe Island palms (Figure 4c).
Pairwise differentiation statistics provided additional quantitative support for the patterns described above (Tables 3 and 4, Table S2).
For example, within Washingtonia and for both the chloroplast and nuclear data, pairwise F st and Φ st values were large and highly significant for the comparisons between W. filifera (SJ) and the peninsular and mainland W. robusta populations. In contrast, none of the comparisons among the sierras of the Baja Peninsula were significant (Table 3a,b).
Within Brahea, we focused on the nuclear data, as only three chloroplast haplotypes were present, one of which was found in 78% of the samples. For the nuclear sequence data, Brahea samples were grouped into ten a priori defined geographical regions (See Section 2,  (Table S2).
Next, we used AMOVA to assess how sequence variation was partitioned among species and a priori defined geographical regions ( Figure 2 and  plastid (94%) data (Scenario 2 in Figure 2a). For Brahea palms, the proportion of variance attributable to differences among the ten sierras was greatest, at around 66% (Scenario 3 in Figure 2b), followed by among taxonomically designated species (62%) and B. brandegeei being restricted to SL (50%) ( Table S4).

| Bayesian cluster analysis
The above results were strongly supported by Bayesian cluster analyses of the combined chloroplast and nuclear sequence datasets  Figure 1 and Table 1 for more details).   Figure 1 and Table 1 for more details).

T A B L E 4 Estimates of genetic differentiation among the ten populations of
Specifically, all of the Washingtonia samples clustered together, as did all of the Brahea samples with the exception of B. edulis from Guadalupe Island, which formed a third cluster.
Finally, we used information on the sampling sites to conduct spatial clustering analyses within BAPS separately for Washingtonia and Brahea based on the chloroplast and nuclear data respectively ( Figure 5).
For Washingonia, Voronoi tessellation graphs of both the chloroplast ( Figure 5a) and nuclear (Figure 5c) sequence data revealed two geographically separate genetic clusters, corresponding to W. filifera in the north (SJ) and W. robusta throughout the rest of the Baja Peninsula. Based on the nuclear data, an additional cluster was recovered corresponding to samples from the mainland Washingtonia populations (Figure 5c). For Brahea, the chloroplast data revealed three clusters (Figure 5b) that exactly mirror the spatial haplotype distributions shown in Figure 2c. By contrast, the nuclear data indicated the presence of only two genetic clusters (Figure 5d), one corresponding to the peninsular and mainland populations combined, and the second to Guadalupe Island.  (Table S5).

| DISCUSSION
We have generated baseline molecular data on species boundaries and the population structure of native palms of Baja California and Guadalupe Island, which provide new insights into the histories of these unique species. Our initial prediction regarding the levels of

| Species delimitation
One of the main goals of our work was to establish species boundaries for Washingtonia and Brahea palms in the study area. The distribution of Washingtonia inferred from our sequence data is broadly in line with previously published morphological information (Felger & Joyal, 1999;Minnich et al., 2011). Despite alleged hybrids having been reported on the basis of morphological data at Cataviña (Cornett, 1987b) and hybridization between W. robusta and W. filifera being common under cultivation (Hodel, 2014), we found no shared haplotypes between W. robusta and W. filifera and the distributional limits of these two species inferred from the sequence data also did not overlap, the southern limit of W. filifera being Sierra Juarez and the northern limit of W. robusta being Cataviña (Figures 3 and 4). Overall, given the lack of shared haplotypes and strong sequence divergence between W. filifera and W. robusta, our data support the classification of these two taxa as separate, geographically delimited, species.
In contrast, our sequence data revealed mixed support for previous There are several possible explanations for the observed differences between the chloroplast and nuclear sequence data within Brahea. First, introgression of chloroplast DNA cannot be discounted, but this seems unlikely given the degree of geographic isolation of Guadalupe Island and its geological origin (Moran, 1996). Second, incomplete linage sorting of shared ancestral polymorphisms could be important, especially given the relatively low levels of observed sequence polymorphism, which could reflect slow chloroplast mutation rates in palms (Wilson, Gaut, & Clegg, 1990). Third, the chloroplast genome reflects seed-mediated gene flow (Hamrick & Nason, 2000) suggesting that a contrasting pattern at the two types of marker could arise if partial geographical isolation interrupted seed-mediated gene flow but not pollen-mediated genetic exchange (Hamilton & Miller, 2002;Petit et al., 2005). Regardless of the exact explanation, incongruent phylogenetic signals between plastid and nuclear data are not uncommon in plants (Pillon et al., 2013), suggesting that an interesting avenue for further research in palms would be to develop a larger and more informative marker panels, for instance using new approaches such as genotyping by sequencing (Narum, Buerkle, Davey, Miller, & Hohenlohe, 2013).

| Patterns of haplotype diversity and structure
Although we expected that sympatrically occurring Washingtonia and Brahea palms would exhibit similar patterns of haplotype structuring along the Baja Peninsula we found the opposite (Figures 3 and 4). For example, as inferred with AMOVA and BAPS, Washingtonia palms on the Baja Peninsula were clustered into two groups comprising palms from SJ and the rest of the peninsula. Contrastingly, AMOVA analysis of the nuclear sequences in Brahea favored a scenario in which each sierra was genetically distinct (Figure 2b, Scenario 3). These differences may be explained by the slightly distinct habitat requirements and differences in the drought and freezing tolerances of Washingtonia and Brahea (Cornett, 1987a,b;Franco-Vizcaino et al., 2007;Minnich et al., 2011). In accordance, ENM indicated that although the strongest ecological predictor of the contemporary distribution of Washingtonia was the mean temperature of the warmest quarter. Predicted temporal changes in the distribution of these two genera over the past 130,000 years also indicated that Washingotnia and Brahea probably responded differently to historical changes in climate (Figures 6 and   7). For example, the contemporary distribution of Washingtonia suggests that its range may have contracted over recent historical time whereas Brahea has maintained a similar distribution for at least 6,000 years.
Dispersal is an essential element that shapes patterns of genetic structure in the Arecaceae (Eiserhardt, Svenning, Kissling, & Balslev, 2011). Although many vertebrates like bats, coyotes and foxes have been proposed as possible seed dispersal agents for Washingtonia and Brahea, the effect of these agents on palm genetic diversity and structure has not been properly tested (Cornett, 2008;Wehncke et al., 2010). In the case of Brahea, it also appears that water pulses may be more significant in terms of seed dispersal than vertebrates (Wehncke et al. 2010;Wehncke, López-Medellín, & Ezcurra, 2009). A comparative study that would quantify the differences in dispersal capacity between Washingtonia and Brahea would definitively help to further our understanding of this topic.
As expected, levels of haplotype diversity were low in Washingtonia for both the chloroplast and nuclear sequence data, whereas nuclear diversity in peninsular Brahea was somewhat higher (Figures 3 and 4). The levels of diversity that we found were broadly León de la Luz, & Cota-Sánchez, 2016;Nason, Hamrick, & Fleming, 2002). This scenario was also supported by ENM analysis (Figure 6c) (Minnich et al., 2011); ENM results also supported this explanation ( Figure 6). However, W. robusta is naturally rare in the northern peninsula, has higher water requirements and seems to be less tolerant to low temperatures than Brahea (Minnich et al., 2011).
We therefore cannot discount the possibility that these observed patterns could have been partly affected by selection favouring the retention of specific haplotypes through linkage disequilibrium to functional loci. One way to test for this would be to genotype functional genetic markers putatively involved in adaptation (Luikart, England, Tallmon, Jordan, & Taberlet, 2003).
The pattern of sequence diversity found in Brahea may be explained by the presence of refugia not in the south, but in the middle part of the peninsula (Garrick, Nason, Meadows, & Dyer, 2009).
However, ENM analysis contradicted this explanation and showed that the predicted distribution of Brahea during the LGM was concentrated across relatively large areas in SL and SSPM (Figure 7). Another possibility that may explain the high genetic diversity of Brahea palms in the middle of the peninsula would be that this region represents the contact zone of recolonization between palms of the southern (SL) and northern (SSPM) refugia during the LGM. This hypothesis would also help to explain why it has been hard to established precise species boundaries between B. brandegeei from the southern and B. armata from the northern peninsula using morphological characteristics (Felger & Joyal, 1999;Minnich et al., 2011).
Moving to the Guadalupe Island palms, we found comparably much lower levels of chloroplast and nuclear diversity. Only one chloroplast and one nuclear haplotype was recovered for 17 and 14 B. edulis individuals respectively. This pattern is in contrast to results found for other tree species of the Island (e.g., Pinus radiata, Karhu, Vogl, Moran, Bell, and Savolainen (2005) and Cupressus guadalupensis, Rosas-Escobar et al. (2011). Although founder effects are often invoked to explain low levels of genetic variation, many other factors may influence the genetic diversity of endemic island plants, including extrinsic ones like island size, distance to the mainland and the amount of available habitat. Another group of factors comprise intrinsic properties such as mutation rate, dispersal capacity, diversity of the population of provenance, and the capacity to respond to human disturbance (Stuessy, Takayama, Lopez-Sepulveda, & Crawford, 2014). When all of these factors are taken into account, the low genetic diversity detected for B. edulis is not surprising. For example, although the total area of Guadalupe Island is about 240.00 km 2 , B. edulis occupies only 0.8 km 2 at the foggy northern extreme of the Island (Garcillán et al., 2012;Oberbauer, 2005). This restricted available habitat reduces the potential establishment of new populations and limits effective population size (Frankham, 1996). Another important factor that may explain the low diversity of B. edulis is the low diversity of Brahea populations from the Baja Peninsula. For example, only two chloroplast haplotypes were detected in over 60 peninsular Brahea individuals. Additionally, the long generation time of Brahea (Bullock & Heath, 2006) and low mutation rates in the Arecaceae (Wilson et al., 1990) may both limit the rate at which new genetic variation is acquired by B. edulis (Stuessy et al., 2014).

| Congruence with historical scenarios
During the last six million years, the Baja California Peninsula has undergone many geological changes that have profoundly affected the distribution and genetic structure of plant and animals species (Dolby et al., 2015;Riddle et al., 2000). Most notably, the formation of the Gulf of California and separation of the Baja California Peninsula from the Mexican mainland, which began approximately 6 million years ago, resulted in the complete or partial isolation of many plant and animal populations on either side of the break (Dolby et al., 2015;Garrick et al., 2009;Gutiérrez Flores et al., 2016;Riddle et al., 2000).
Currently, the Gulf of California represents a major barrier to gene flow even for certain organisms with moderate dispersal capabilities (Grismer, 2000;Riddle et al., 2000), so in some respects it is unsurprising that it also appears to impeded gene flow in sessile organisms like plants (Cavender-Bares et al., 2015). Consistent with this, we found a strong pattern of haplotype differentiation between mainland and peninsular W. robusta at both nuclear and chloroplast sequences. Our results from the chloroplast data for Brahea were similar, but as only one individual from the Mexican mainland was analyzed, we could not draw any firm conclusions. Nonetheless, although we lacked samples from the Mexican mainland for Brahea, we were able to include an endemic Brahea palm species from Guadalupe Island. This volcanic island has never been in contact with the Baja Peninsula and as a consequence almost 16% of its native plant species are endemics (Batiza, 1977;León De La Luz et al., 2003). Deep nuclear isolation of Guadalupe Island's Brahea from its peninsular sister species is consistent with a strong effect of dispersal barriers like sea channels on the population structure of Brahea palms (Karhu et al., 2005;Rosas-Escobar et al., 2011).
Historically, the Baja has also undergone a series of submersions and uplifts (Lindell et al., 2006;Riddle et al., 2000). During these partial submersions, the Peninsula was fragmented into one or more islands and the Pacific Ocean was connected to the Gulf of California.
More specifically, inundation of the Isthmus of La Paz isolated the Cabo region (SL) from the rest of the Baja Peninsula, whereas the mid-peninsular seaway temporally isolated the southern and northern parts of the Baja Peninsula (Lindell et al., 2006;Riddle et al., 2000).
Nonetheless, the exact locations of the proposed seaways as well as the timing of their formation are still under debate (Dolby et al., 2015).
Regardless of the precise time and location estimates, these breaks are thought to be responsible for the observed patterns of genetic structuring in a variety of taxa (Garrick et al., 2009 Washingtonia's nuclear and plastid sequence data did not provide any support for either of these breaks, with clustering results indicating the presence of a single genetically homogeneous population along the peninsula. By contrast, for Brahea, the chloroplast data strongly supported the presence of two genetically differentiated clades along the peninsula (Figures 3a,c, and 5b), providing support for historical marine intrusion into lowland areas of the Isthmus of La Paz.