Phylogeography and demographic history of Lacerta lepida in the Iberian Peninsula: multiple refugia, range expansions and secondary contact zones
© Miraldo et al; licensee BioMed Central Ltd. 2011
Received: 4 October 2010
Accepted: 17 June 2011
Published: 17 June 2011
The Iberian Peninsula is recognized as an important refugial area for species survival and diversification during the climatic cycles of the Quaternary. Recent phylogeographic studies have revealed Iberia as a complex of multiple refugia. However, most of these studies have focused either on species with narrow distributions within the region or species groups that, although widely distributed, generally have a genetic structure that relates to pre-Quaternary cladogenetic events. In this study we undertake a detailed phylogeographic analysis of the lizard species, Lacerta lepida, whose distribution encompasses the entire Iberian Peninsula. We attempt to identify refugial areas, recolonization routes, zones of secondary contact and date demographic events within this species.
Results support the existence of 6 evolutionary lineages (phylogroups) with a strong association between genetic variation and geography, suggesting a history of allopatric divergence in different refugia. Diversification within phylogroups is concordant with the onset of the Pleistocene climatic oscillations. The southern regions of several phylogroups show a high incidence of ancestral alleles in contrast with high incidence of recently derived alleles in northern regions. All phylogroups show signs of recent demographic and spatial expansions. We have further identified several zones of secondary contact, with divergent mitochondrial haplotypes occurring in narrow zones of sympatry.
The concordant patterns of spatial and demographic expansions detected within phylogroups, together with the high incidence of ancestral haplotypes in southern regions of several phylogroups, suggests a pattern of contraction of populations into southern refugia during adverse climatic conditions from which subsequent northern expansions occurred. This study supports the emergent pattern of multiple refugia within Iberia but adds to it by identifying a pattern of refugia coincident with the southern distribution limits of individual evolutionary lineages. These areas are important in terms of long-term species persistence and therefore important areas for conservation.
Evidence from phylogeographic studies suggests that southern Europe, including the peninsulas of Iberia, Italy and the Balkans and areas near the Caucasus and the Caspian Sea, have functioned as refugial areas for species survival during periods of adverse climatic conditions . Recent studies also emphasize the important role that these refugial areas had in shaping the evolutionary history of species that have persisted within these regions for several ice ages. As suggested by early studies [2–4] the topographic complexity and geographic mosaic of habitats in southern refugial peninsulas have favoured the occurrence of multiple disjunct refugia, allowing the persistence of isolated populations within them during glacial periods. Within the Iberian Peninsula, complex species histories have been revealed for a variety of taxa, with some showing remarkable patterns of phylogeographic concordance [see  and references therein] involving deep genetic subdivisions, high haplotype richness and distinct hybrid zones. Not only has the Iberian Peninsula sourced the northern redistribution of species after ice ages, but it has also facilitated diversification through patterns of repeated population fragmentation, contraction, expansion and admixture. Phylogeographic analyses for the golden-striped salamander, Chioglossa lusitanica [6–8] and Schreiber's Lizard, Lacerta schreiberi [9, 10] provide good examples of the complexity that most likely typifies many species within this major peninsular glacial refugium.
However, even though the Iberian Peninsula is the best studied glacial refugium, most phylogeographic studies have focused on species that either have a narrow distribution within the region (e.g Chioglossa lusitanica , Lacerta schreiberi , Lissotriton boscai ) or involve species groups that, although distributed across the entire region, generally have a genetic structure that relates to older cladogenetic events (e.g Podarcis spp. [12, 13], Alytes spp , Oryctolagus cunniculus [15, 16]). In order to better understand the complex phylogeographic history of Iberian species, and the way they have responded to Pleistocene climatic oscillations, it is important to study species with distributions encompassing the entire Iberian Peninsula. For this purpose we have used the ocellated lizard, Lacerta lepida (Daudin, 1802), as a model to study the impact of Pleistocene climatic changes in generating and structuring intraspecific genetic diversity on this regional scale. Lacerta lepida is typically Mediterranean, with a distribution encompassing all the Iberian Peninsula, and with phylogeographic structure across the region . Several mitochondrial lineages that appear to have non-overlapping geographic ranges were recently described, suggesting a history of allopatric differentiation in multiple refugia during the Plio-Pleistocene . Here we assess the broader phylogeographic patterns within Lacerta lepida with the specific aims to i) clarify the distribution of mtDNA phylogroups; ii) identify refugial areas within these phylogroups during the glacial periods; iii) date the main demographic and evolutionary events within Lacerta lepida; and finally iv) identify secondary contact zones between the different phylogroups.
It is generally accepted that phylogeographic histories recovered using only mtDNA as a marker are constrained to reveal one genealogy that reflects the maternally inherited natural history of an organism. Relationships among phylogroups inferred through mtDNA might be discordant with inferences made based on nuclear genes (Harrison 1991; Avise 2000) and discordances have been illustrated in several studies [e.g. [18–24]]. Particularly within the Iberian Peninsula, several studies have emphasized the importance of using different types of markers to fully recover the complex evolutionary and demographic scenarios that most likely characterize the species that have persisted there across the Quaternary. For example, in Lacerta schreiberi evidence for gene flow and ancestral introgression between apparently allopatric mtDNA lineages was only detected by the use of nuclear markers . Here we use both mtDNA and nDNA derived genealogies, since their contrasting molecular and population properties (principally uniparental versus biparental mode of inheritance and contrasting population sizes) are mutually informative when there have been opportunities for secondary contact, gene flow and hybridization between diverging populations.
Lizards were captured under licence during the years 2005, 2006 and 2007, sampling most of the distribution of Lacerta lepida in Portugal, Spain, and France. Sampling intensity was concentrated in regions known to contain high genetic divergence within western and south-eastern part of Iberia . Lizards were captured using tomahawk traps or by hand, and tissue samples were taken by clipping 1 cm of the tail tip that was subsequently preserved in 100% ethanol. After tissue sampling, animals were immediately released back into the wild in the place of capture. Geographic coordinates of sampling sites were recorded with a GPS.
DNA extraction, amplification and sequencing
Total genomic DNA was extracted from ethanol-preserved muscle tissue using a salt extraction protocol [25, 26]. A fragment of 627 base pairs (bp) of the mitochondrial DNA (mtDNA) cytochrome b (cytb) gene was amplified using the truncated version of primer L14841  (CYTBF, 5'-CCA TCC AAC ATC TCA GCA TGA TGA AA-3') and the modified version of primer MV16  (CYTBR, 5'- AAA TAG GAA GTA TCA CTC TGG TTT-3') to increase specificity for Lacerta lepida, using information from published sequences in GenBank. Polymerase chain reactions (PCRs) were performed in a total volume of 25 μl, containing 0.5 U of Taq polymerase (BIOTAQ™), 4 mM of MgCl2, 0.4 mM of each nucleotide (Bioline), 0.4 μM of each primer, 2 μl of 10 × NH4 reaction buffer (Bioline) and approximately 50 ng of DNA. PCR amplifications were conducted as follows: DNA was initially denaturated at 94°C for 3 min followed by 35 cycles of denaturation at 94°C for 45 s, annealing at 51°C for 45 s and extension at 72°C for 45 s, plus a final extension step at 72°C for 5 min. Negative controls (no DNA) were included for all amplifications. PCR products were visualized on a 2% agarose gel and purified by filtration through QIAquick® columns (Qiagen) following the manufacturer's recommendations. Purified PCR products were sequenced in both directions using the above primers, Taq polymerase BigDye Terminator v3.1™ (Applied Biosystems) and 30-90 ng of PCR product.
Intron 7 of the β-fibrinogen gene (β-Fibint7) has been successfully used as a nuclear marker in several vertebrate phylogeographic and phylogenetic studies [e.g. [13, 29–32]]. Specifically it was recently employed for a phylogenetic study of the genus Lacerta  where it revealed sufficient variation within Lacerta lepida to be useful for phylogeographical inference. Initial β-fibint7 amplifications were performed using primers FIB-B17U (5'- GGA GAA AAC AGG ACA ATG ACA ATT CAC - 3') and FIB-B17L (5' - TCC CCA GTA GTA TCT GCC ATT AGG GTT - 3')  and the conditions described in Paulo et al. . However, due to low amplification and sequencing success a nested PCR approach as suggested by Sequeira et al.  was subsequently adopted. A fragment of 788bp was first amplified from genomic DNA using primers FIB-B17U and FIB-B17L (PCRa). The product of this reaction (1 μl) was then used as a template for a subsequent PCR of 691bp (PCRb) using primer BFXF  and BFX8 . Both amplifications were performed in a total volume of 25 μl, and included reagents in the same concentrations as those specified for cytb gene fragment. PCR cycle conditions were the same as described for cytb fragment but the primer annealing temperatures were 55°C and 56°C, for PCRa and PCRb respectively. Negative controls (no DNA) were included for all amplifications. Purified PCRb products were then sequenced with primers BFBX and BFX8 using identical sequencing conditions as for the mtDNA cytb sequencing. All PCRs and sequencing reactions were performed in a MJ Research thermocycler (PTC-240 DNA Engine Tetrad 2) and sequences were obtained using an ABI 3700 capillary sequencer.
DNA sequences were aligned by eye using BioEdit Sequence Alignment Editor 7.01 . β-fibint7 alleles of heterozygous individuals were inferred using Phase version 2.1 [34, 35]. Ten replicate Phase runs were conducted using 1000 burnin steps and 1000 iterations. The Phase probability threshold was first set to 0.90, but not all genotypes were resolved at this threshold. Garrick et al. , suggest that lowering the Phase threshold to 0.60 reduces the number of unresolved haplotypes with little or no increase in the number of false positives. In contrast, omitting unresolved haplotypes may have a substantial impact on downstream analyses, such as the estimation of Tajima's D and Fu's Fs statistics, as unresolved haplotypes usually represent rare alleles . We therefore reduced the Phase threshold to 0.60 to maximise recovery of haplotypes. Several tests implemented in the software RDP3 [Recombination Detection Program, ] were used to detect evidence for recombination in β-fibint7: RDP , GENECONV , Maximum Chi Square [40, 41], Chimaera  and Sister Scanning . Due to the small size of the fragment used in the analyses (315 bp), the window size for the recombination detection methods was set to 20bp.
Phylogenetic analysis and haplotype network construction
Phylogenetic relationships among cytb haplotypes were inferred using Bayesian inference as implemented in Mrbayes v3.1.2 [43, 44]. Data was divided into two partitions: (1st+2nd) and (3rd) codon positions. The most appropriate substitution model for each partition was selected using Modeltest v3.7  and the Akaike Information Criterion (AIC) . The model selected for both partitions was the same (GTR + I) but with a different proportion of invariable sites (pinv(1st+2nd) = 0.79 and pinv(3rd) = 0.09). The different partitions were allowed to evolve at different rates, with unlinked topology and unlinked parameters for the nucleotide substitution models. Four Markov chains were run for 10 million generations. To avoid local optima we used two independent runs, and to improve swapping of states between heated and cold chains the heating parameter was decreased to 0.05. Trees were sampled every 1000 generations and the first 100 trees were discarded as burn-in. Posterior probabilities were obtained from the 50% majority-rule consensus tree of the retained trees.
Intraspecific gene genealogies were inferred using the median-joining (MJ)  and the statistical parsimony (SP)  network construction approaches. The MJ network was computed with the program Network 4.5.0 (http://www.fluxus-engineering.com) and the SP network was inferred using the program Tcs 1.21 . For the MJ approach the parameter ε was set to 0, preventing less parsimonious pathways from being included in the analysis. The SP network was inferred with a parsimony confidence limit of 95%, allowing the inclusion of less parsimonious alternatives that fall within this confidence limit. Ambiguities within networks were resolved where possible following the criteria of Crandall & Templeton . The phylogenetic and β-Fibint7 network analyses were rooted using gene sequences of the African sister species Lacerta pater (Lataste, 1880) [17, 51] (Genebank accession numbers: AF: 378967 and EU: 365413 for cytb and β-Fibint7 genes respectively).
Estimation of divergence times
Divergence times within and between phylogroups were estimated using Beast version 1.4.8  and the cytb dataset. Beast performs Bayesian statistical inferences of parameters, such as divergence times, by using Markov Chain Monte Carlo (MCMC) as a framework. Input files were generated with Beauti version 1.4.8 . The nucleotide substitution model and its parameter values were selected according to the results of Modeltest version 3.7 , with upper and lower bounds around the values defined as 120% and 80% respectively . According to the Bayesian Information Criterion (BIC) and the hierarchical Likelihood Ratio Tests (hLRT's), the model of nucleotide substitution identified as the best fit to the data is the HKY model  with a gamma distribution (Γ) for substitution rates across sites (shape parameter, α = 0.2889) and no category of invariable sites. An uncorrelated lognormal relaxed molecular clock was used with mean mutation rate of 0.01 mutations/site/million years and standard deviation of 0.0015, assuming a normal distribution, as prior information. The mutation rate used was based on the provisional calibration of a reptile-specific molecular clock presented in Paulo et al.  who used published cytochrome b data sets from two studies of Canary island reptiles (Gallotia spp.) [56, 57] calibrated with the geological age of El Hierro. The standard deviation included in our study encompasses a faster mutation rate of 0.0115 and a slower mutation rate of 0.0085 mutations/site/million years. This takes into account potential underestimation of the mutation rate due to violation of the assumption of immediate island colonization in the clock calibration of Paulo et al. , or overestimation due to the longer generation time and larger body size of Lacerta spp. when compared to Gallotia spp. . Two runs were executed for 106 generations, sampling every 500 generations and discarding the first 10% as burn in. Results of the two runs were displayed and combined in Tracer v1.4  to check for stationarity and ensure that ESSs were above 200. For all analyses one sequence of Lacerta pater (Genebank accession number: AF378967) was included as an outgroup.
In a second approach to estimate divergence times within phylogroups the method of Saillard et al.  was employed where each extant haplotype descending from the most recent common ancestor (MRCA) represents a time interval between the present and the MRCA. Average distances from the MRCA within each phylogroup are calculated from the number of mutational steps separating each haplotype from the MRCA. The absolute timing of divergence is then calculated by multiplying the observed values by the average mutational changes per lineage per million years (Myr). Three mutation rates were used: 0.01 mutations/site/million years, representing the average mutation rate; a faster mutation rate of 0.0115, and a slower mutation rate of 0.0085.
Neutrality tests and demographic analyses
In order to detect departures from a constant population size under the neutral model, Tajima's D , Ramos-Onsins & Rozas' R 2  and Fu's Fs  tests were applied to both mtDNA and nDNA datasets. Both R 2 and Fs statistics have been shown to be the best statistical tests to detect population growth (R 2 has been suggested to behave better for small sample sizes whereas Fs is better for bigger ones) . Population expansions have also been shown to leave particular signatures in the distribution of pairwise sequence differences [64, 65]. We capitalized upon this by employing statistics based on the mismatch distribution to test for demographic expansions. The observed distribution of pairwise differences between haplotypes within phylogroups was compared with the expected results under a sudden-demographic and a spatial-demographic expansion model. Statistically significant differences between observed and simulated expected distributions were evaluated with the sum of the square deviations (SSD) and Harpending's raggedness index (hg) [66, 67]. For the β-Fibint7 gene we used two datasets: the first including only haplotypes inferred by Phase with probability higher than 0.90, and the second including all haplotypes inferred with the less stringent Phase probability threshold of 0.60. Tests were performed with Arlequin version 3.11  for Tajima's D, Fu's Fs, SDD and hg, and with DNAsp version 4.50  for R 2 and expected values for the mismatch distribution. The demographic history of each phylogroup was also inferred using a coalescent-based approach. The model used to infer past population dynamics was the Bayesian Skyline Plot (BSP, ) as implemented in Beast version 1.4.8 . For each phylogroup we carried out two independent runs of 10 million generations each. Trees and parameters were sampled every 1000 iterations, with a burn in of 10%. Results of each run were visualized using Tracer v1.4  to ensure stationarity and convergence had been reached, and that effective sample sizes (ESS) were higher than 200.
Geographical distribution of alleles and inference of refugial areas
Using predictions from coalescent theory (haplotypes at the tips of a tree are younger than interior haplotypes to which they are connected) ancestral and derived haplotypes within each phylogroup were identified, thus obtaining a temporal framework for inferring haplotype origin within phylogroups. The null hypothesis of random geographic distribution of haplotypes was tested statistically using Geodis version 2.5  with significance testing by permuting the data 106 times. When non-random associations of haplotypes with geography were detected, the geographic distribution of ancestral versus derived haplotypes (interior and tip haplotypes, respectively) was further explored to identify possible refugial areas and the directionality of previously detected demographic and spatial expansions. Coalescent theory predicts refugial areas to exhibit a high frequency of ancestral haplotypes , and therefore the average number of mutations of a haplotype from the MRCA is expected to be significantly lower in refugial areas than would be expected by chance. In contrast, recently colonized areas are predicted to show a high incidence of derived haplotypes , and therefore the average number of mutations of a haplotype from the MRCA is expected to be significantly higher than would be expected by chance. To test the hypothesis that refugial areas are located in the south of a phylogroup's distribution, with the north having been more recently colonized, we estimated the average distances of alleles from the MRCA within both the north and south of the geographic range of a phylogroup. We then tested whether these values were significantly high or low by randomizing the geographic states of alleles to generate a null distribution of mean haplotype distance from the MRCA within both the south and north of a phylogroup's distribution.
Average pairwise genetic distances between Lacerta lepid a mitochondrial phylogroups.
uncorrected p distances (%)
HKY + Γ corrected (%)
Phylogenetic and network relationships among haplotypes
The root of the network is located along the branch that connects the very divergent lineages L and N, allowing for the inference of the most recent common ancestor (MRCA) for each phylogroup (Figure 3). Although this identification is straightforward for clade N (haplotype 126), the networks reveal two probable ancestral haplotypes within clade L (haplotypes 127 and 104), which are connected to haplotype 126 through a loop. SP and MJ networks constructed with 0-fold degenerate sites only, thus reducing homoplasy within the data set (data not shown) result in the collapse of this reticulation, and haplotype 126 (lineage N) connects unambiguously to haplotype 127 (lineage L), supporting 127 as the ancestral haplotype within phylogroup L. Divergent mitochondrial haplotypes from different phylogroups occurring in sympatry were detected in several populations. The admixed populations were: populations 6, 23, 24 and 25 (phylogroups L2 and L4); populations 15, 17, 18, 19, 20, 21 22 and 96 (phylogroups L4 and L5); population 33 (phylogroups L4 and N); population 44 (phylogroups L1 and L4) and population 56, 58, 78 and 82 (phylogroups L1, L3 and L5) (Additional file 1, Table S1).
Results from mismatch distribution and neutrality tests for cytb mtDNA phylogroups and for the β-fibint 7 nuclear gene.
Spatial genetic structure
The nuclear data failed to recover the phylogroups detected by the mtDNA dataset. Nevertheless, when the nuclear dataset is grouped according to the mtDNA phylogroups previously identified, structure in the distribution of alleles can be detected. For this analysis each of the mtDNA phylogroups was considered as a geographic region and Geodis was used to test for geographical structure amongst the nuclear genetic variation. Nuclear haplotypes are shared among some of the mtDNA phylogroups (with the exception of phylogroup N): all haplotypes from L3 (B7 and B20) occur either in L1 (B7) or in L5 (B20); all haplotypes from L2 (B4, B6, B13 and B14) occur in L4 (with B13 also occurring in L5); all haplotypes from L5 (B5, B13 and B20) occur either in L4 (B5 and B13), either in L3 (B20) or in L2 (B13) and finally haplotypes from L1 (B7 and B17) occur either in L3 (B7) or in L4 (B17). Additionally all phylogroups share haplotype B1 (the most common and one of the most ancestral haplotypes) and all phylogroups apart from L3 and L5 also share haplotype B15 (the most ancestral haplotype). It is important to note that only geographically close phylogroups share nuclear haplotypes. It is also possible to identify private haplotypes within the geography of several mtDNA phylogroups: B2 was only detected in phylogroup L1; haplotypes B3, B9, B10, B14, B16 and B18 in L4 and haplotypes B8, B11, B12 and B19 in phylogroup N.
Divergence time estimates in million years (Ma) from the most recent common ancestor (mrca) of all group members from each Lacerta lepida phylogroup estimated using the method of Saillard et al. (2000) using three different mutation rates; and divergence time estimates for the main nodes recovered in the phylogenetic analysis using a mean mutation rate of 2% in Beast.
(mean ± s.d.)
0.64 ± 0.12
0.75 ± 0.14
0.55 ± 0.10
0.45 ± 0.21
0.53 ± 0.24
0.39 ± 0.16
0.59 ± 0.38
0.63 ± 0.40
0.47 ± 0.27
0.92 ± 0.30
1.08 ± 0.35
0.80 ± 0.24
0.59 ± 0.19
0.70 ± 0.22
0.51 ± 0.15
0.85 ± 0.34
0.99 ± 0.41
0.74 ± 0.28
Neutrality tests and demographic analyses
Geographical distribution of alleles and inference of refugial areas
Detailed analysis of the distribution of mitochondrial genetic variation within Lacerta lepida across the Iberian Peninsula revealed a complex phylogeographic history for the species. The cytb genealogy clearly defines 6 geographic phylogroups with generally non-overlapping geographic distributions. The strong association of mtDNA genetic variation with geography suggests a history of allopatric divergence in different refugia within the Iberian Peninsula, a pattern that has been described for several taxa within the region [see , and references therein]. Although this pattern of differentiation of distinct evolutionary units in allopatry was less evident from the analysis of the nuclear data, the distribution of nuclear haplotypes is not in conflict with the mtDNA phylogroups. Failure of the nuclear gene genealogy to reveal concordant genetic structure with the mitochondrial genealogy can be expected if one takes into account the fact that nuclear genes take on average four times longer to reach monophyly than mitochondrial ones . In fact, most of the intraspecific differentiation within Lacerta lepida is of relatively recent origin, with the majority of phylogroups (L1 to L5) estimated to have diversified within the Pleistocene, increasing the probability of mitochondrial lineages not being reciprocally monophyletic for nuclear markers. Nevertheless, phylogroup N is estimated to have diverged from the remaining phylogroups during the Miocene (5.6-13.7 Ma) representing a much older cladogenetic event within the group. This older split therefore allows for more time for lineage sorting at the nuclear level. This is evident in the composition of nuclear haplotypes of phylogroup N, which are almost all private with the exception of the ancestral haplotypes (which are shared among almost all phylogroups). Thus, although mtDNA lineages have not reached monophyly at the nuclear level, some level of nuclear differentiation between mtDNA lineages is detected by the existence of private alleles, and this is most pronounced for phylogroups that represent older cladogenetic events (L1, L4 and N). Incomplete lineage sorting seems therefore a plausible explanation for the discrepancies in mitochondrial and nuclear gene genealogies. Recent phylogeographic and phylogenetic studies focusing on lizards in the Iberian Peninsula reveal similar discrepant patterns between mtDNA and nuclear genealogies [10, 13]. In the case of Podarcis wall lizards Pinho et al.  showed that, despite a complete lack of monophyly at the nuclear level, most of the species show historical reproductive isolation and that the lack of monophyly is mainly a result of shared ancestral polymorphism.
Although it is plausible that incomplete lineage sorting is responsible for the failure of nuclear genealogy to recover the mitochondrial phylogroups in Lacerta lepida, current gene flow, mainly male-biased, at zones of secondary contact can also be invoked to explain the patterns observed. Evidence for this comes from the observation that geographically close phylogroups share more derived haplotypes. For example allele B20 occurs only in the very divergent phylogroups L3 and L5 near the zone of contact, but it was not detected in the ancestral phylogroup L4, suggesting that nuclear gene flow may be occurring between the L3 and L5 mtDNA lineages. The same is true for allele B4 which, although more widespread within phylogroup L4, is also found in individuals of L2 near the zone of contact, but was not detected in phylogroup L3, which is phylogenetically closer to phylogroup L2.
Historical biogeography of Lacerta lepida
Divergence within Lacerta lepida is estimated to have started in the Miocene, approximately 9.4 Mya (5.58-13.66), with divergence of phylogroup N from phylogroup L. Within phylogroup L estimated divergence times are much younger and inferred to have been initiated in the Plio-Pleistocene, approximately 1.96 Mya (1.13-2.91 Mya). Interestingly, although divergences between the mitochondrial phylogroups have a Plio-Pleistocene (within group L) or a late Miocene (for group N) origin, haplotype diversities within phylogroups indicate a strong influence of later Pleistocene events, with diversification within phylogroups starting between 0.92 and 0.45 Ma. The importance of the Pleistocene climatic oscillations in promoting species differentiation in the Iberian Peninsula has been emphasized by previous studies [see , for a recent review] and this clearly also seems to be the case for Lacerta lepida.
Phylogroup N is distributed across the Betic Mountain range in south-western Spain and its distribution approximately coincides with the described subspecies Lacerta lepida nevadensis. The existence of clear significant morphological differentiation between L. l. nevadensis and the remaining subspecies of Lacerta lepida [74–76], together with high levels of allozyme differentiation  and mitochondrial genetic differentiation (our study and ), suggests that the two lineages (N and L) might in fact reflect two different species. Paulo et al.  have inferred that divergence between phylogroups L and N were initiated by overseas dispersal between what was then the Iberian mainland and the Betic Massif that at that time existed as an island between Iberia and North Africa. Under this scenario subsequent contact between the phylogroups would have been initiated after the merging of the Betic Massif with Iberian mainland, due to the closing of the Betic corridor 7.8-7.6 Ma [see [77, 78] for a detailed explanation of the kinematics of the western Mediterranean basin]. The Betic Mountain range is a region of high endemism for plants and animals, and its importance as a refugium for other taxa has been highlighted before [see ]. Interestingly, some of these taxa show similar divergence times and distribution as phylogroup N (e.g. Salamandra salamandra longirostris  and Alytes dickhilleni ).
Phylogroups L2 and L3
The monophyletic group composed of phylogroups L2 and L3 is estimated to have started diverging from the remaining phylogroups in the early Pleistocene, approximately 1.5 Mya (0.82-2.27). Interestingly these two phylogroups have a disjunct distribution, with phylogroup L2 occupying the south of Portugal whereas L3 occupies the north-western parts of the Iberian Peninsula. The intervening region between phylogroups L2 and L3 is occupied by phylogroup L5. A vicariant event during the Middle Pleistocene (0.82-2.27 Mya) triggering divergence between the L2-L3 lineage and the remaining populations of Lacerta lepida seems probable. It is noteworthy that most phylogeographic studies within Iberia reveal similar phylogenetic breaks associated with the same period (e.g. Chioglossa lusitanica [6, 7], Oryctolagus cuniculus , Lacerta schreiberi [9, 82]), suggesting a common vicariant history. Such a vicariant event was most likely climatically-mediated as no apparent geographical barrier exists within western Iberian Peninsula.
The western Algarve region in southern Portugal has been indicated as the evolutionary centre for several species and also as a key refugial area [83, 84]. The region harboured relicts of temperate forests during the Last Glacial Maximum , probably providing suitable conditions for species survival through glacial periods. The high genetic distance and disjunct distribution found between phylogroups L2 and L3 is most likely the consequence of fragmentation of a once more continuous range in response to a cooling climate. Subsequent climate amelioration during an interglacial period probably resulted in this distribution gap being colonized by L5, expanding from another nearby refugium. A similar distribution pattern of genetic variation is found within Discoglossus galganoi across Portugal, with two closely related phylogroups distributed in the south and north and a less related phylogroup bisecting their distribution . The distribution of Discoglossus galganoi phylogroups and the evolutionary relationship between them is remarkably similar to that found within Lacerta lepida.
Phylogroups L1, L4 and L5
The refugia for Lacerta lepida phylogroups L1, L4 and L5 were probably located in the south-eastern side of the Guadalquivir basin. Support for this comes from the distribution of ancestral haplotypes 127, 140 and 8 within phylogroup L4 that are found only in this region (localities 10, 12, 32, 36, 37, 38, 39, 40, 41). The two very divergent haplotypes within L4 (60 and 61) are also restricted to this region. The widespread distribution of the most frequently sampled haplotype 1 suggests that the spatial and demographic expansion detected within L4 was of a "leading edge" type [87, 88], with few individuals rapidly colonizing adjacent regions, leading to a decrease in genetic diversity on the newly colonized areas. The different phylogroups are most likely the result of three different expansions from the southern refugia dominated by different ancestral haplotypes, followed by further divergence in allopatry.
Refugial areas, range expansions, and secondary contact
The lower mutation rate and slower fixation rate of nuclear genes, when compared to the mitochondrial genome, mean that nuclear genealogies may be more indicative of older demographic events . Therefore, the analysis of the distribution of ancestral nuclear alleles can be useful in the identification of refugial areas. Haplotype B15 is the root of the nuclear gene network, representing the most ancestral allele in the dataset. The distribution of B15 is primarily in the southern region of Iberia, where it occurs with higher frequency than in northern regions (70% of samples with haplotype B15 are from southern latitudes). Furthermore the southern region of Iberia (to the south of river Tagus) shows higher nuclear haplotype richness, with 90% of the identified haplotypes occurring in this region from which 67% do not occur further north. The high nuclear haplotype richness found in the southern regions of Iberia, together with the high frequency of B15 in this region suggest that southern populations are older and the source of subsequent northern expansions.
The geographical and genealogical relationships of mitochondrial haplotypes within phylogroups provide further indications as to the probable refugial areas from where range expansions subsequently occurred. The same pattern of older (interior) haplotypes being typically found within southern sampling sites is also evident for the mitochondrial dataset. For example within phylogroup L2, the ancestral haplotype 96 and the related interior haplotype 93 are most frequent in southern sites 25, 63, 64, 66, 69 and 70, coincident with the Algarve and southern Alentejo region in Portugal. Within phylogroup L3 the ancestral haplotype 51 and the related interior haplotype 46 occur either in the southern region of the phylogroup distribution, to the south of Douro river (sites 82 and 78), or in sites located at the centre of the phylogroup distribution (sites 87, 89 and 92). Within phylogroup L4, the ancestral haplotype 127 and the closely related haplotypes 140 and 8 occur again only in the south-eastern limit of this phylogroup (sites 10, 12, 32, 36, 37, 38, 39, 40 and 41), near the source of the Guadalquivir river in Jaen, Spain. Finally the ancestral haplotype of phylogroup L5 (haplotype 13) and the related interior haplotypes 14 and 17 are most frequent within southern sites of this phylogroup (sites 22, 95, 96 and 99), along the river Tagus basin. These patterns are statistically significant for phylogroups L3 and L5 (Figure 8), which also reveal an excess of derived haplotypes in the north of their geographic ranges. These repeated pattern of ancestral haplotypes occurring more frequently towards the southern range of phylogroups, and in some cases derived haplotypes being more frequent within more northern sites, suggests that southern regions within phylogroup ranges have probably functioned as refugia during past adverse climatic conditions, from which subsequent northern expansions occurred.
The evolutionary history of Lacerta lepida is consistent with a history of population fragmentation and allopatric divergence in southern refugia, followed by recent demographic and spatial range expansions. Demographic expansions are supported by the neutrality and mismatch distribution tests that reveal signatures of demographic and spatial expansions in almost all phylogroups and by the BSP analysis, which show a trend of population expansion in all phylogroups between 0.1 to 0.15 million years ago (Figure 7). Spatial expansions are supported by the latitudinal variation observed for the distributions of ancestral and derived haplotypes. This history of diversification in allopatry followed by range expansions has lead to the establishment of at least four secondary contact zones where very divergent mitochondrial haplotypes belonging to different phylogroups were found in sympatry (Figure 4). Further analyses of these zones may provide insights into the mechanisms involved in speciation and divergence within this lizard species.
Mitochondrial and nuclear gene genealogies in Lacerta lepida provide evidence for a history of isolation and divergence in allopatry resulting in the diversification of six genetically and geographically distinct lineages. Although diversification within the group is largely concordant with the onset of the major glaciations at the beginning of the Pleistocene (approximately 2 Mya), an earlier event associated with the Miocene was also identified. This event (approximately 9 Mya), which marks the divergence of lineages N and L, seems to be associated with geological events related to the evolution of the Mediterranean basin. Signatures of recent demographic and spatial expansions were apparent in all phylogroups, with several phylogroups having established zones of secondary contact. Analyses of the distribution of ancestral and derived alleles within each phylogroup, and inferences related to the biogeography of Lacerta lepida, allowed the identification of probable refugia within the Iberian Peninsula, suggesting southern refugial areas within phylogroups. Our work highlights the importance of these areas for the long-term conservation and management of diversity with Lacerta lepida across its geographic range.
This work was supported by a grant from FCT, Portugal (POCI/BIA-BDE/59288/2004). AM was supported by a PhD grant from FCT (SFRH/BD/16996/2004). We thank S. Carranza (University of Barcelona) and J. Pleguezuelos (University of Granada) for providing several lizard samples. We thank Juan Pablo de la Vega, Luis Garcia, Gabriel Marin, Eugenia Zarza-Franco, Rui Osório, Pedro Figueiredo, Nuno Valente and Matthew Osborne for helping with fieldwork. We also thank Mário Miraldo and António Firmino for assembling the lizards' traps. All samples were collected with appropriate licenses in each country (Instituto de Conservação da Natureza in Portugal and the Regional Ministry of Environment of the Government of Andalucia, Asturias, Galicia, Extremadura, Castilla y Leon, Murcia, Castilla la Mancha and Madrid in Spain).
- Hewitt GM: Genetic consequences of climatic oscillations in the Quaternary. Philosophical Transactions of the Royal Society B: Biological Sciences. 2004, 359: 183-195. 10.1098/rstb.2003.1388.View Article
- Cooper SJB, Hewitt GM: Nuclear DNA sequence divergence between parapatric subspecies of the grasshopper Chorthippus parallelus. Insect Molecular Biology. 1993, 2: 185-194. 10.1111/j.1365-2583.1993.tb00138.x.View ArticlePubMed
- Hewitt GM: Some genetic consequences of ice ages, and their role, in divergence and speciation. Biological Journal of the Linnean Society. 1996, 58: 247-276.View Article
- Hewitt GM: After the Ice: Parallelus meets Erythropus in the Pyrenees. Hybrid zones and the evolutionary process. Edited by: Harrison RG. 1993, New York: Oxford University Press
- Gomez A, Lunt DH: Refugia within refugia: patterns of phylogeographic concordance in the Iberian Peninsula. Phylogeography of Southern European Refugia. Edited by: Weiss S, Ferrand N. 2007, Dordrecht: Springer
- Alexandrino J, Froufe E, Arntzen JW, Ferrand N: Genetic subdivision, glacial refugia and postglacial recolonization in the golden-striped salamander, Chioglossa lusitanica (Amphibia: Urodela). Molecular Ecology. 2000, 9: 771-781. 10.1046/j.1365-294x.2000.00931.x.View ArticlePubMed
- Alexandrino J, Arntzen JW, Ferrand N: Nested clade analysis and the genetic evidence for population expansion in the phylogeography of the golden-striped salamander, Chioglossa lusitanica (Amphibia: Urodela). Heredity. 2002, 88: 66-74. 10.1038/sj.hdy.6800010.View ArticlePubMed
- Sequeira F, Alexandrino J, Rocha S, Arntzen JW, Ferrand N: Genetic exchange across a hybrid zone within the Iberian endemic golden-striped salamander, Chioglossa lusitanica. Molecular Ecology. 2005, 14: 245-254.View ArticlePubMed
- Paulo OS, Dias C, Bruford MW, Jordan WC, Nichols RA: The persistence of Pliocene populations through the Pleistocene climatic cycles: evidence from the phylogeography of an Iberian lizard. Proceedings of the Royal Society B: Biological Sciences. 2001, 268: 1625-1630. 10.1098/rspb.2001.1706.View ArticlePubMedPubMed Central
- Godinho R, Crespo EG, Ferrand N: The limits of mtDNA phylogeography: complex patterns of population history in a highly structured Iberian lizard are only revealed by the use of nuclear markers. Molecular Ecology. 2008, 17: 4670-4683. 10.1111/j.1365-294X.2008.03929.x.View ArticlePubMed
- Martínez-Solano I, Teixeira J, Buckley D, Garcia-Paris M: Mitochondrial DNA phylogeography of Lissotriton boscai (Caudata, Salamandridae): evidence for old, multiple refugia in an Iberian endemic. Molecular Ecology. 2006, 15: 3375-3388. 10.1111/j.1365-294X.2006.03013.x.View ArticlePubMed
- Harris DJ, Sá-Sousa P: Molecular phylogenetics of Iberian Wall lizards (Podarcis): Is Podarcis hispanica a species complex?. Molecular Phylogenetics and Evolution. 2002, 23: 75-81. 10.1006/mpev.2001.1079.View ArticlePubMed
- Pinho C, Harris DJ, Ferrand N: Non-equilibrium estimates of gene flow inferred from nuclear genealogies suggest that Iberian and North African wall lizards (Podarcis spp.) are an assemblage of incipient species. BMC Evolutionary Biology. 2008, 8: 63-10.1186/1471-2148-8-63.View ArticlePubMedPubMed Central
- Martínez-Solano I, Gonçalves HA, Arntzen JW, García-París M: Phylogenetic relationships and biogeography of midwife toads (Discoglossidae: Alytes). Journal of Biogeography. 2004, 31: 603-618. 10.1046/j.1365-2699.2003.01033.x.View Article
- Branco M, Ferrand N, Monnerot M: Phylogeography of the European rabbit (Oryctolagus cuniculus) in the Iberian Peninsula inferred from RFLP analysis of the cytochrome b gene. Heredity. 2000, 85: 307-317. 10.1046/j.1365-2540.2000.00756.x.View ArticlePubMed
- Branco M, Monnerot M, Ferrand N, Templeton AR: Postglacial dispersal of the european rabbit (Oryctolagus cuniculus) on the Iberian Peninsula reconstructed from nested clade and mismatch distribution analyses of mitochondrial DNA variation. Evolution. 2002, 56: 792-803.View ArticlePubMed
- Paulo OS, Pinheiro J, Miraldo A, Bruford MW, Jordan WC, Nichols RA: The role of vicariance vs. dispersal in shaping genetic patterns in ocellated lizard species in the western Mediterranean. Molecular Ecology. 2008, 17: 1535-1551. 10.1111/j.1365-294X.2008.03706.x.View ArticlePubMed
- Ujvari B, Dowton M, Madsen T: Population genetic structure, gene flow and sex-biased dispersal in frillneck lizards (Chlamydosaurus kingii). Molecular Ecology. 2008, 17: 3557-3564.PubMed
- Lindell J, Méndez-de la Cruz FR, Murphy RW: Deep biogeographical history and cytonuclear discordance in the black-tailed brush lizard (Urosaurus nigricaudus) of Baja California. Biological Journal of the Linnean Society. 2008, 94: 89-104. 10.1111/j.1095-8312.2008.00976.x.View Article
- Dowling DK, Friberg U, Lindell J: Evolutionary implications of non-neutral mitochondrial genetic variation. Trends in Ecology & Evolution. 2008, 23: 546-554. 10.1016/j.tree.2008.05.011.View Article
- Zink RM, Barrowclough GF: Mitochondrial DNA under siege in avian phylogeography. Molecular Ecology. 2008, 17: 2107-2121. 10.1111/j.1365-294X.2008.03737.x.View ArticlePubMed
- Leaché AD, McGuire JA: Phylogenetic relationships of horned lizards (Phrynosoma) based on nuclear and mitochondrial data: Evidence for a misleading mitochondrial gene tree. Molecular Phylogenetics and Evolution. 2006, 39: 628-644. 10.1016/j.ympev.2005.12.016.View ArticlePubMed
- McGuire JA, Linkem CW, Koo MS, Hutchison DW, Kristopher A, David L, Orange I, Lemos-Espinal J, Riddle BR, Jaeger JR: Mitochondrial introgression and incomplete lineage sorting through space and time: phylogenetics of Crotaphytid lizards. Evolution. 2007, 61: 2879-2897. 10.1111/j.1558-5646.2007.00239.x.View ArticlePubMed
- Thorpe RS, Surget-Groba Y, Johansson H: The relative importance of ecology and geographic isolation for speciation in anoles. Philosophical Transactions of the Royal Society B: Biological Sciences. 2008, 363: 3071-3081. 10.1098/rstb.2008.0077.View Article
- Aljanabi SM, Martinez I: Universal and rapid salt-extraction of high quality genomic DNA for PCR-based techniques. Nucleic Acids Research. 1997, 25: 4692-4693. 10.1093/nar/25.22.4692.View ArticlePubMedPubMed Central
- Sunnucks P, Hales DF: Numerous transposed sequences of mitochondrial cytochrome oxidase I-II in aphids of the genus Sitobion (Hemiptera: Aphididae). Molecular Biology and Evolution. 1996, 13: 510-524.View ArticlePubMed
- Kocher TD, Thomas WK, Meyer A, Edwards SV, Paabo S, Villablanca FX, Wilson AC: Dynamics of mitochondrial-DNA evolution in animals - amplification and sequencing with conserved primers. Proceedings of the National Academy of Sciences of the United States of America. 1989, 86: 6196-6200. 10.1073/pnas.86.16.6196.View ArticlePubMedPubMed Central
- Moritz C, Schneider CJ, Wake DB: Evolutionary relationships within the Ensatina-Eschscholtzii complex confirm the ring species interpretation. Systematic Biology. 1992, 41: 273-291.View Article
- Prychtko TM, Moore WM: The utility of DNA sequences of an intron from the beta-fibrinogen gene in phylogenetic analysis of woodpeckers (Aves: Picidae). Molecular Phylogenetics and Evolution. 1997, 8: 193-204. 10.1006/mpev.1997.0420.View Article
- Dolman G, Phillips B: Single copy nuclear DNA markers characterized for comparative phylogeography in Australian wet tropics rainforest skinks. Molecular Ecology Notes. 2004, 4: 185-187. 10.1111/j.1471-8286.2004.00609.x.View Article
- Godinho R, Mendonca B, Crespo EG, Ferrand N: Genealogy of the nuclear beta-fibrinogen locus in a highly structured lizard species: comparison with mtDNA and evidence for intragenic recombination in the hybrid zone. Heredity. 2006, 96: 454-463. 10.1038/sj.hdy.6800823.View ArticlePubMed
- Sequeira F, Ferrand N, Harris DJ: Assessing the phylogenetic signal of the nuclear β-Fibrinogen intron 7 in salamandrids (Amphibia: Salamandridae). Amphibia-Reptilia. 2006, 27: 409-418. 10.1163/156853806778190114.View Article
- Hall TA: BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows95/98/NT. Nucleic Acids Symposium Series. 1999, 41: 95-98.
- Stephens M, Smith NJ, Donnelly P: A new statistical method for haplotype reconstruction from population data. The American Journal of Human Genetics. 2001, 68: 978-989. 10.1086/319501.View ArticlePubMed
- Stephens M, Scheet P: Accounting for decay of linkage disequilibrium in haplotype inference and missing-data imputation. The American Journal of Human Genetics. 2005, 76: 449-462. 10.1086/428594.View ArticlePubMed
- Garrick R, Sunnucks P, Dyer R: Nuclear gene phylogeography using PHASE: dealing with unresolved genotypes, lost alleles, and systematic bias in parameter estimation. BMC Evolutionary Biology. 10: 118-
- Martin DP, Williamson C, Posada D: RDP2: recombination detection and analysis from sequence alignments. Bioinformatics. 2005, 21: 260-262. 10.1093/bioinformatics/bth490.View ArticlePubMed
- Martin DP, Rybicki E: RDP: detection of recombination amongst aligned sequences. Bioinformatics. 2000, 16: 562-563. 10.1093/bioinformatics/16.6.562.View ArticlePubMed
- Padidam M, Sawyer S, Fauquet CM: Possible emergence of new geminiviruses by frequent recombination. Virology. 1999, 265: 218-225. 10.1006/viro.1999.0056.View ArticlePubMed
- Posada D, Crandall KA: Evaluation of methods for detecting recombination from DNA sequences: computer simulations. Proc Natl Acad Sci USA. 2001, 98: 13757-13762. 10.1073/pnas.241370698.View ArticlePubMedPubMed Central
- Smith JM: Analyzing the mosaic structure of genes. Journal of Molecular Evolution. 1992, 34: 126-129.PubMed
- Gibbs MJ, Armstrong JS, Gibbs AJ: Sister-Scanning: a Monte Carlo procedure for assessing signals in recombinant sequences. Bioinformatics. 2000, 16: 573-582. 10.1093/bioinformatics/16.7.573.View ArticlePubMed
- Huelsenbeck JP, Ronquist F: MRBAYES: Bayesian inference of phylogenetic trees. Bioinformatics. 2001, 17: 754-755. 10.1093/bioinformatics/17.8.754.View ArticlePubMed
- Ronquist F, Huelsenbeck JP: MrBayes 3: Bayesian phylogenetic inference under mixed models. Bioinformatics. 2003, 19: 1572-1574. 10.1093/bioinformatics/btg180.View ArticlePubMed
- Posada D, Crandall KA: MODELTEST: testing the model of DNA substitution. Bioinformatics. 1998, 14: 817-818. 10.1093/bioinformatics/14.9.817.View ArticlePubMed
- Akaike H: Information theory and an extension of the maximum likelihood principle. Proceedings of the 2nd International Symposium on Information Theory; Akadémia Kiado, Budapest. Edited by: Petrov BN, Csaki F. 1973, 267-281.
- Bandelt HJ, Forster P, Rohl A: Median-joining networks for inferring intraspecific phylogenies. Molecular Biology and Evolution. 1999, 16: 37-48.View ArticlePubMed
- Templeton AR, Crandall KA, Sing CF: A cladistic analysis of phenotypic associations with haplotypes inferred from restriction endonuclease mapping and DNA sequence data. III. Cladogram estimation. Genetics. 1992, 132: 619-633.PubMedPubMed Central
- Clement M, Posada D, Crandall KA: TCS: a computer program to estimate gene genealogies. Molecular Ecology. 2000, 9: 1657-1659. 10.1046/j.1365-294x.2000.01020.x.View ArticlePubMed
- Crandall KA, Templeton AR: Empirical tests of some predictions from coalescent theory with applications to intraspecific phylogeny reconstruction. Genetics. 1993, 134: 959-969.PubMedPubMed Central
- Perera A, Harris DJ: Genetic variability in the Ocellated lizard Timon tangitanus in Morocco. African Zoology. 45: 321-329.
- Drummond AJ, Rambaut A: BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evolutionary Biology. 2007, 7: 214-10.1186/1471-2148-7-214.View ArticlePubMedPubMed Central
- Rambaut A, Drummond A: BEAUTi v1.4.2. Bayesian Evolutionary Analysis Utility. 2007
- Emerson BC: Alarm bells for the molecular clock? No support for Ho et al.'s model of time-dependent molecular rate estimates. Systematic Biology. 2007, 56: 337-345. 10.1080/10635150701258795.View ArticlePubMed
- Hasegawa M, Kishino H, Yano TA: Dating of the human-ape splitting by a molecular clock of mitochondrial DNA. Journal of Molecular Evolution. 1985, 22: 160-174. 10.1007/BF02101694.View ArticlePubMed
- Thorpe RS, McGregor DP, Cumming AM, Jordan WC: DNA evolution and colonization sequence of island lizards in relation to geological history: mtDNA, RFLP, Cytochrome B, Cytochrome Oxidase, 12s RRNA sequence, and Nuclear RAPD analysis. Evolution. 1994, 48: 230-240. 10.2307/2410090.View Article
- González P, Pinto F, Nogales M, Jiménez-asensio J, Hernández M, Cabrera VM: Phylogenetic relationships of the Canary Islands endemic lizard genus Gallotia (Sauria: Lacertidae), inferred from mitochondrial DNA sequences. Molecular Phylogenetics and Evolution. 1996, 6: 63-71. 10.1006/mpev.1996.0058.View ArticlePubMed
- Martin AP, Palumbi SR: Body size, metabolic rate, generation time, and the molecular clock. Proceedings of the National Academy of Sciences of the United States of America. 1993, 90: 4087-4091. 10.1073/pnas.90.9.4087.View ArticlePubMedPubMed Central
- Rambaut A, Drummond A: Tracer v1.4. 2007, [http://beast.bio.ed.ac.uk/Tracer]
- Saillard J, Forster P, Lynnerup N, Bandelt HJ, Nørby S: mtDNA variation among Greenland Eskimos: the edge of the Beringian expansion. The American Journal of Human Genetics. 2000, 67: 718-726. 10.1086/303038.View ArticlePubMed
- Tajima F: Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics. 1989, 123: 585-595.PubMedPubMed Central
- Ramos-Onsins SE, Rozas J: Statistical properties of new neutrality tests against population growth. Molecular Biology and Evolution. 2002, 19: 2092-2100.View ArticlePubMed
- Fu YX: Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics. 1997, 147: 915-925.PubMedPubMed Central
- Slatkin M, Hudson RR: Pairwise comparisons of mitochondrial DNA sequences in stable and exponentially growing populations. Genetics. 1991, 129: 555-562.PubMedPubMed Central
- Rogers AR, Harpending H: Population growth makes waves in the distribution of pairwise genetic differences. Molecular Biology and Evolution. 1992, 9: 552-569.PubMed
- Harpending HC, Sherry ST, Rogers AR, Stoneking M: The genetic structure of ancient human populations. Current Anthropology. 1993, 34: 483-496. 10.1086/204195.View Article
- Harpending HC: Signature of ancient population growth in a low-resolution mitochondrial DNA mismatch distribution. Human Biology. 1994, 66: 591-600.PubMed
- Excoffier L, Laval G, Schneider S: Arlequin ver. 3.0: An integrated software package for population genetics data analysis. Evolutionary Bioinformatics Online. 2005, 47-50.
- Rozas J, Sanchez-DelBarrio JC, Messeguer X, Rozas R: DnaSP, DNA polymorphism analyses by the coalescent and other methods. Bioinformatics. 2003, 19: 2496-2497. 10.1093/bioinformatics/btg359.View ArticlePubMed
- Drummond AJ, Rambaut A, Shapiro B, Pybus OG: Bayesian coalescent inference of past population dynamics from molecular sequences. Molecular Biology and Evolution. 2005, 22: 1185-1192. 10.1093/molbev/msi103.View ArticlePubMed
- Posada D, Crandall KA, Templeton AR: GeoDis: a program for the cladistic nested analysis of the geographical distribution of genetic haplotypes. Molecular Ecology. 2000, 9: 487-488. 10.1046/j.1365-294x.2000.00887.x.View ArticlePubMed
- Wakeley J: Coalescent theory: An introduction. 2008, Greenwood Village, Colorado: Roberts & Company Publishers
- Moore WS: Inferring phylogenies from mtDNA variation: mitochondrial-gene trees versus nuclear-gene trees. Evolution. 1995, 49: 718-726. 10.2307/2410325.View Article
- Mateo JA, Castroviejo J: Variation morphologique et revision taxonomique de l'espece Lacerta lepida Daudin, 1802 (Sauria, Lacertidae). Bulletin du Museé de Histoire Naturele de Paris. 1990, 12: 691-706.
- Mateo JA, López-Jurado LF: Variaciones en el color de los lagartos ocelados; aproximacion a la distribuicion de Lacerta lepida nevadensis Buchholz 1963. Revista Espanola de Herpetologia. 1994, 8: 29-35.
- Mateo JA, López-Jurado LF, Guillaume CP: Variabilité électrophorétique et morphologique des lézards ocellés (Lacertidae): un complexe d'espèces de part et d'autre du détroit de Gibraltar. Comptes Rendus de L'Academie des Sciences Serie iii-Sciences de la Vie-Life Sciences. 1996, 319: 737-746.
- Rosenbaum G, Lister GS, Duboz C: Relative motions of Africa, Iberia and Europe during Alpine orogeny. Tectonophysics. 2002, 359: 117-129. 10.1016/S0040-1951(02)00442-0.View Article
- Rosenbaum G, Lister GS, Duboz C: Reconstruction of the tectonic evolution of the western Mediterranean since the Oligocene. Journal of the Virtual Explorer. 2002, 8: 107-130.
- Garcia-Paris M, Alcobendas M, Alberch P: Influence of the Guadalquivir river basin on mitochondrial DNA evolution of Salamandra salamandra (Caudata: Salamandridae) from southern Spain. Copeia. 1998, 1998: 173-176. 10.2307/1447714.View Article
- Arntzen JW, Garcia-Paris M: Morphological and allozyme studies of midwife toads (genus Alytes), including the description of two new taxa from Spain. Contributions to zoology. 1995, 65: 5-34.
- Biju-Duval C, Ennafaa H, Dennebouy N, Monnerot M, Mignotte F, Soriguer R, El Gaied A, El Hili A, Monoulou JC: Mitochondrial DNA evolution in Lagomorphs: origin of systematic heteroplasmy and organization of diversity in European rabbits. Journal of Molecular Evolution. 1991, 33: 92-102. 10.1007/BF02100200.View Article
- Paulo OS, Jordan WC, Bruford MW, Nichols RA: Using nested clade analysis to assess the history of colonization and the persistence of populations of an Iberian Lizard. Molecular Ecology. 2002, 11: 809-819. 10.1046/j.1365-294X.2002.01484.x.View ArticlePubMed
- Mesquita N, Hanfling B, Carvalho GR, Coelho MM: Phylogeography of the cyprinid Squalius aradensis and implications for conservation of the endemic freshwater fauna of southern Portugal. Molecular Ecology. 2005, 14: 1939-1954. 10.1111/j.1365-294X.2005.02569.x.View ArticlePubMed
- Fritz U, Barata M, Busack SD, Fritzsch G, Castilho R: Impact of mountain chains, sea straits and peripheral populations on genetic and taxonomic structure of a freshwater turtle, Mauremys leprosa (Reptilia, Testudines, Geoemydidae). Zoologica Scripta. 2006, 35: 97-108. 10.1111/j.1463-6409.2005.00218.x.View Article
- Zagwijn WH: Migration of vegetation during the Quaternary in Europe. Courier Forschungsinstitut Senckenberg. 1992, 153: 9-20.
- Martínez-Solano I: Phylogeography of Iberian Discoglossus (Anura: Discoglossidae). Journal of Zoological Systematics & Evolutionary Research. 2004, 42: 298-305. 10.1111/j.1439-0469.2004.00257.x.View Article
- Hewitt GM: The subdivision of species by hybrid zones. Speciation and its consequences. Edited by: Otte D, Endler J. 1989, Sunderland, Massachusetts: Sinauer Associates, 85-110.
- Hewitt GM: Postglacial distribution and species substructure: lessons from pollen, insects and hybrid zones. Evolutionary patterns and processes. Edited by: Lees DR, Edwards D. 1993, London: Academic Press, 14: 97-123.
- Avise JC: Molecular Markers, Natural History and Evolution. 2004, Sunderland, Massachusetts: Sinauer Associates, 2
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.