Pre-Quaternary divergence and subsequent radiation explain longitudinal patterns of genetic and morphological variation in the striped skink, Heremites vittatus
© The Author(s). 2017
Received: 3 February 2017
Accepted: 12 May 2017
Published: 9 June 2017
Many animal and plant species in the Middle East and northern Africa have a predominantly longitudinal distribution, extending from Iran and Turkey along the eastern Mediterranean coast into northern Africa. These species are potentially characterized by longitudinal patterns of biological diversity, but little is known about the underlying biogeographic mechanisms and evolutionary timescales. We examined these questions in the striped skink, Heremites vittatus, one such species with a roughly longitudinal distribution across the Middle East and northern Africa, by analyzing range-wide patterns of mitochondrial DNA (mtDNA) sequence and multi-trait morphological variation.
The striped skink exhibits a basic longitudinal organization of mtDNA diversity, with three major mitochondrial lineages inhabiting northern Africa, the eastern Mediterranean coast, and Turkey/Iran. Remarkably, these lineages are of pre-Quaternary origin, and are characterized by p-distances of 9–10%. In addition, within each of these lineages a more recent Quaternary genetic diversification was observed, as evidenced by deep subclades and high haplotype diversity especially in the Turkish/Iranian and eastern Mediterranean lineages. Consistent with the genetic variation, our morphological analysis revealed that the majority of morphological traits show significant mean differences between specimens from northern Africa, the eastern Mediterranean coast, and Turkey/Iran, suggesting lineage-specific trait evolution. In addition, a subset of traits exhibits clinal variation along the eastern Mediterranean coast, potentially indicating selection gradients at the geographic transition from northern Africa to Anatolia. The existence of allopatric, morphologically and genetically divergent lineages suggests that Heremites vittatus might represent a complex with several taxa.
Our work demonstrates that early divergence events in the Pliocene, likely driven by both climatic and geological factors, established the longitudinal patterns and distribution of Heremites vittatus. Subsequent radiation during the Pleistocene generated the genetic and morphological diversity observed today. Our study provides further evidence that longitudinal diversity patterns and species distributions in the Middle East and northern Africa were shaped by complex evolutionary processes, involving the region’s intricate geological history, climatic oscillations, and the presence of the Sahara.
KeywordsPhylogeography Western Palearctic Cline Latitude Longitude Climate oscillation Heremites vittatus Intraspecific variation
Intraspecific patterns of biological diversity are often the result of geological and climatic processes in the past. In the western Palearctic, the best understood example is the widespread presence of latitudinal gradients in genetic diversity in European taxa [1–4]. During interglacial periods in the Pleistocene, these species rapidly and repeatedly extended from Mediterranean refugia into central and northern Europe. Genetic bottlenecks during the cyclical expansions led to allele loss and decreased heterozygosity, and ultimately a reduction in genetic diversity in northern areas. While latitudinal gradients can be found in many species in Europe, the geographical shape of southern Europe made longitudinal movements more infrequent. Longitudinal patterns of genetic diversity are thus often confined to the different Mediterranean peninsulas, and can be difficult to disentangle from latitudinal gradients [5, 6].
Compared to continental Europe, relatively little is known about the structure and evolution of biological diversity in other regions of the Palearctic, notably in the Middle East and northern Africa [7–12]. Several species in northern Africa display post-glacial latitudinal range expansions, often into Europe despite the potential barrier of the Mediterranean Sea . Remarkably, many species in this region share a roughly longitudinal distribution extending from Iran and Turkey along the eastern Mediterranean coast into northern Africa . Unlike in Europe, strong sea barriers are absent in this region, which potentially facilitated longitudinal movements. In addition, the presence of the Sahara as a southern barrier for many species may have further supported longitudinal migration. European species with more pronounced longitudinal patterns of genetic diversity often originated in Asia Minor and expanded westwards along the north Mediterranean coast. As such, Asia Minor has become known as a center of genetic diversity in the Middle East [14–16]. However, whether species that expanded from Asia Minor into northern Africa instead of Europe are also characterized by longitudinal patterns remains poorly understood. More generally, few studies have addressed the biogeographic mechanisms and underlying evolutionary timescales in species with longitudinal distributions in this region [5, 6, 17, 18]. Here, we examined these questions in the striped skink, Heremites vittatus, a scincid lizard with a predominantly longitudinal distribution across the Middle East and northern Africa.
Lizards of the genus Mabuya sensu lato (s.l.) are some of the most widely distributed skinks with a circumtropical distribution. For a long time, these skinks were collectively allocated to the genus Mabuya s.l., until the Afro-Malagasy Mabuya taxa were assigned to the genus Euprepis [19, 20], and later re-assigned to Trachylepis . The Middle Eastern Mabuya s.l. species were preliminarily included into Trachylepis, although they form a separate radiation . To account for this, the genus Heremites was recently revalidated for these species, including H. vittatus, H. auratus, and H. septemtaeniatus . Of these, H. vittatus has the largest distribution range, occurring in Algeria, Tunisia, Libya, Egypt, Israel, Jordan, Lebanon, Cyprus, Syria, Turkey, Iraq and Iran [24–29]. Despite its wide distribution, the striped skink has been regarded as a monotypic species . In several parts of this large distribution range, e.g., in Turkey, Iran, and Lebanon, morphometric studies uncovered considerable morphological variation within local populations and differences between populations on a regional scale [31–36]. Part of this variation is likely due to adaptation to local environmental conditions, rather than phylogenetic divergence. For example, in a population in southeastern Turkey, uniform skinks are more abundant in habitats with low grass cover, while striped specimens are more abundant in habitats with high grass cover, indicative of disruptive selection by visual predators . A distribution model for H. vittatus demonstrated that the species is predominantly found in areas with high winter precipitation (>300 mm), and rainy winters may be the driving factor behind its distribution . Populations in northern Africa inhabit wetlands and oases , while those in Iran and Turkey live in mountainous areas . Due to its poor dispersal skills and its dependence on humid habitats in the arid Middle Eastern environment, the striped skink is a sensitive indicator of geographic processes that have driven the distribution and intraspecific evolution of animals in this region .
In this work, we combined investigations of mtDNA (cytochrome b) sequence and multivariate morphological variation across the species’ distribution range, and provide the first comprehensive analysis of the biogeography and evolutionary history of the striped skink. We chose to study both genetic and morphological variation because evolutionary processes can affect genetics and morphology in different ways and at different times during the period of divergence. Unlike putatively neutral molecular markers, morphological traits are under varying degrees of selection; thus, it may be expected that morphological traits can show both clinal variation as a result of local adaptation, as well as trait conservation as a result of deep divergence events . Together, insights from genetic and morphological variation can thus provide for a richer and more integrated understanding of the underlying evolutionary processes [39, 40].
P-distances between major haplogroups. Below diagonal, mean p-distances between MHGs. Above diagonal, standard error of the mean (SEM) for mean p-distances between MHGs. On diagonal, mean p-distance ± SEM within MHGs. Standard errors were determined from 1000 Bootstrap replicates
MHG1 (Cyprus, Israel, Jordan)
MHG3 (Jordan, Lebanon, Syria)
0.062 ± 0.008
MHG1 (Cyprus, Israel, Jordan)
0.012 ± 0.003
0.010 ± 0.004
MHG3 (Jordan, Lebanon, Syria)
0.010 ± 0.004
0.010 ± 0.003
0.003 ± 0.002
0.002 ± 0.002
Phylogeography and divergence times
Within-group mean and range of traits with significant mean differences between groups. Trait differences between groups were tested for significance with Dunnett’s modified Tukey-Kramer (DTK) test
W-E; W-C; C-E
# middorsals (MDN)
**; NS; ****
# midventrals (MVN)
****; ****; NS
# longitudinal scale rows (LSN)
***; NS; ****
# subdigital lamellae under 4th toe (SDLN)
****; ****; NS
Normalized width of tail at tail base (HTc)
****; **; NS
Normalized height of tail at tail base (VTc)
***; **; NS
Normalized head length (HLc)
**; NS; *
Normalized head width (HBc)
*; NS; NS
Normalized length of hind leg (LHc)
**; NS; *
Normalized length of front leg (LFc)
NS; NS; *
Normalized length of tibia (LTc)
****; ****; **
Normalized length of forearm (LAc)
***; **; *
Normalized inter-nare distance (INAc)
*; NS; *
Alternatively, the effect of latitude/longitude on the position in the morphospace could be explained by clinal variation within geographical groups. While the western (Algeria, Tunisia, Libya) and eastern (Turkey, Iran) groups have a roughly longitudinal distribution, the central group (Egypt, Israel, Jordan, Lebanon, Cyprus, Syria) is predominantly distributed along the latitudinal dimension (Fig. 7a). To account for these different geographical extensions, we regressed traits against longitude for the western and eastern group, and against latitude for the central group. 10 of 19 traits (53%) showed evidence of clinal variation in one geographical group. No trait showed evidence for clinal variation in more than one group. Of the 10 traits, seven showed clinal variation in the central group. Two traits showed clinal variation in the western group, and one trait in the eastern group (Fig. 7f).
We next investigated how mean differences between groups overlap with clinal variation within groups. Of 19 traits, three (16%) showed evidence of neither mean nor clinal variation (Fig. 7b). Six traits (32%) showed evidence of only mean differences (Fig. 7c), while three traits (16%) showed evidence of only clinal variation. Seven traits (37%) showed evidence of both mean and clinal variation. Of the latter seven traits, five showed clinal variation in the central group (Fig. 7d-e). Taken together, our results suggest that the majority of morphological traits in H. vittatus exhibit mean differences between geographical groups with either no evidence for clinal variation or clinal variation along the eastern Mediterranean coast.
Our mitochondrial phylogeography suggests that Heremites vittatus diverged from H. auratus in the late Miocene (13.3 mya; range 7.9–19.7 mya), consistent with a previous estimate that dated this split to 10.54 mya . Heremites vittatus likely emerged as a distinct species in southwest Asia, where the Heremites radiation originated . In the mid Pliocene (3.6 mya; range 2.3–5.2 mya), the ancestors of the eastern Mediterranean and northern African lineages split off from the ancestral Irano-Anatolian lineage and became allopatric. This divergence event coincided with prominent geological events in Anatolia during the late Neogene. The Anatolian mountain ranges, e.g., the Amanus and Taurus mountains, created geographical barriers between the Anatolian and eastern Mediterranean regions in the South. The vast lake system in central Anatolia delimited the distribution into western Anatolia. Other emerging dispersal barriers subsequently supported this divergence, notably the Amik plain in southern Turkey. As a result, vicariant cladogenesis, frequent among other animals in this geologically complex region (e.g. green lizards), may explain the divergence of the Irano-Anatolian lineage and the ancestors of the eastern Mediterranean and northern African lineages [14, 43, 45, 46].
During much of the Pliocene, the climate was relatively warm and humid in the Palearctic , thus potentially facilitating the colonization of northern Africa from the eastern Mediterranean. In the late Pliocene, however, the climate suddenly aridified and became colder [41, 48, 49]. Around this time, the northern African lineages split off (2.5 mya; range 1.6–3.6 mya). Consistent with the preference of the modern species for more humid habitats in northern Africa [24, 37], this climatic turn may have restricted the species to the remaining humid regions, and isolated them from populations inhabiting the eastern Mediterranean coast. Similar divergence events coinciding with the transition from hotter/wetter to colder/drier climate in the late Pliocene/early Pleistocene have previously been observed in other animal species of the Mediterranean (e.g., tree frogs, Typhlops vermicularis) [8, 11, 50]. Notably, H. vittatus is the only Heremites species in northern Africa today, suggesting that other Heremites species never migrated into northern Africa, or were unable to withstand the drastic climatic changes at the Pliocene/Pleistocene transition.
During the Pleistocene, the Sahara increasingly aridified, which restricted more humid habitats to a few coastal areas (e.g., Cyrenaica) and oases [41, 51]. This aridification may have led to the further fragmentation of the northern African lineage, as evidenced by the deep split between the Tunisian and Libyan haplotypes in the Cyrenaica (1.9 mya; range 1.0–2.9 mya). Although these climatic cycles must have influenced distributions during the Pleistocene, the mitochondrial lineages remained separate and apparently did not mix. The extant, allopatric populations in Libya, Tunisia, and Algeria are thus potential remnants of a once more extended distribution in northern Africa. Similar fragmentation processes in the Pleistocene likely affected other animal species in northern Africa currently limited to the Cyrenaica and coastal areas of Tunisia and Algeria, such as lizards of the genus Ophisops .
In the late Pleistocene, the Irano-Anatolian and eastern Mediterranean lineages further radiated, establishing the mitochondrial diversity observed today. The Egyptian, Cypriot-Israeli-Jordanian, Jordanian-Lebanese-Syrian, Lebanese, and Syrian major haplogroups (MHGs) all evolved into monophyletic lineages during this time. Based on our dataset, most of these lineages occupy non-overlapping regions, suggesting that these haplotypes are not admixed. One interesting exception is the co-occurrence of the Jordanian-Lebanese-Syrian MHG with the Cypriot-Israeli-Jordanian MHG at one locality in Jordan (Al Himma), raising the possibility of admixture between these lineages in this region. In our dataset, the monophyletic lineage on the island of Cyprus is sister to the lineages in Israel and Jordan and diverged approximately 0.5 mya (range 0.3–0.8 mya). During the Pleistocene, Cyprus was likely never connected to the mainland through a land bridge [25, 52], so that the species must have reached Cyprus via overseas dispersal, which has also been suggested for other species in Cyprus, e.g., Hyla savignyi [8, 25, 53]. However, these hypotheses may not be correct if the true mainland source of the Cypriot population was not included in our dataset.
Although our study shows again the usefulness of fast-evolving mtDNA genes to dissect more recent intraspecific divergence in Palearctic reptiles (e.g., [54–58]), we cannot rule out the possibility of discordant evolutionary histories of the mitochondrial and nuclear genome [59–61]. The mtDNA lineage evolution described in this study should thus be corroborated with nuclear sequence variation before additional inferences (e.g. for taxonomic purposes) are made on the evolutionary history of this complex.
Our multidimensional scaling analysis suggests that morphological variation in Heremites vittatus is linked to latitude and longitude. At the level of individual traits, most show mean differences between three morphological groups (western, central, eastern) coinciding with the putative ranges of the three old mitochondrial lineages. We note that while these results are suggestive, proof of concordant genetic and morphological variation would require data collection from identical specimens. Interestingly, we further find that these traits show either (1) no evidence for within-group clines, or (2) clinal latitudinal variation in the central (eastern Mediterranean) group.
The first category, mean differences with no evidence for within-group clines, argues for a close association of trait variation with the three old mtDNA lineages. For example, the western (northern African) group has an average of 18.2 subdigital lamellae, while the central (eastern Mediterranean) and eastern (Turkish/Iranian) group have on average 16.4 and 16.1, respectively, with no evidence for a clinal increase of lamellae in the central group (Fig. 7c). Thus, although they are likely not closely related the eastern and central group share a similar number of subdigital lamellae, which may reflect an ancestral state that was independently maintained in both lineages. By contrast, the higher number of subdigital lamellae in the western group may be an adaptation to the xeric environmental conditions across northern Africa because sand-dwelling lizards often possess pedal specializations that evolved as (convergent) adaptations to sandy habitats . Mean differences between geographical groups can thus be explained by lineage-specific ancestral or derived trait states.
The second category, mean differences between groups with latitudinal clines within the central (eastern Mediterranean) group, suggests that some traits are affected by selection gradients in addition to lineage-specific trait states. One such trait is relative head length. Lizards in the western group have on average longer heads relative to their body length than the eastern group, with no evidence of clinal variation in either group. By contrast, relative head length in the geographically intermediate central group decreases gradually from the western (northern African) to eastern (Turkish/Iranian) trait states: the farther north specimens were collected along the eastern Mediterranean coast, the shorter their heads are relative to body length. The large extension (roughly 1000 km) of this cline suggests that the central group is likely not a hybridization zone of the western and eastern group, but rather a distinct population affected by gradual environmental change. Remarkably, we found little evidence that similar clines exist in the western or eastern group, pointing to unique environmental changes at the transition from northern Africa to Anatolia [45, 46].
Taken together, our morphological analysis reveals that morphological variation in Heremites vittatus has a clear geographical organization, and is likely both affected by lineage-specific trait conservation and local trait adaptation.
Many species in the Middle East and northern Africa exhibit predominantly longitudinal distributions, and are potentially characterized by longitudinal patterns of biological diversity, but the underlying biogeographic mechanisms and evolutionary timescales remain largely unclear. In this paper, we provide the first comprehensive assessment of the mitochondrial phylogeography and multivariate morphological variation of the striped skink, Heremites vittatus, a species with a predominantly longitudinal distribution across the Middle East and northern Africa. We suggest that H. vittatus evolved as a separate species in southwest Asia in the late Miocene. We show that the species diverged into three main mitochondrial lineages roughly 2.5–3.6 mya in the eastern, central, and western regions of the distribution range, introducing a longitudinal organization of genetic diversity. This divergence was likely driven by a combination of geological changes in Anatolia and climatic changes in northern Africa. The morphological variation reflects this pattern of mitochondrial divergence: the majority of morphological traits show significant differences between these regions, supporting the idea of independently evolving lineages. We also find that subsequent local radiation in the Pleistocene, likely driven by climatic oscillations and locally emerging geographical barriers, led to the mitochondrial diversity observed today. Morphologically, clinal variation along the eastern Mediterranean coast in a subset of traits suggests that ongoing selection pressures and local adaptation may play an important role in shaping populations at the transition from northern Africa to Anatolia. Our findings of allopatric, morphologically and genetically divergent lineages raise the possibility that Heremites vittatus represents a complex with several undescribed taxa. However, taxonomic decisions should be based on nuclear, mitochondrial, and morphological data collected in the context of future studies on the inter- and intraspecific systematics and taxonomy of Heremites. In conclusion, our study suggests that longitudinal patterns and species distributions in the Middle East and northern Africa may be the result of complex evolutionary processes, driven by the region’s geological and climatic history, geographical setup, and the presence of the Sahara.
DNA sequence analysis
The DNA sampling consists of 46 H. vittatus specimens (Additional file 3), including one sequence from GenBank, covering most of the distribution range (Fig. 1a). We also sequenced five H. auratus specimens. We included Trachylepis quinquetaeniata from the African lineage of Mabuya s.l. as an outgroup, because this lineage is one of the closest relatives of Middle Eastern Heremites [23, 63–65]. We did not include in the analyses sequences of H. vittatus from GenBank that only partially overlap the alignment generated in this study (but see Additional file 2).
Small pieces of tissue were sampled either from museum specimens (toe or tongue) of various ages and storage conditions, or were provided by colleagues based on field collections (part of tail or liver). Samples were stored in 70% ethanol or EDTA buffer, and genomic DNA was extracted with standard protocols . A portion of the mitochondrial cytochrome b gene was amplified by polymerase chain reaction (PCR) with the primers mt-c-emys (5′-CCG GAT CAA ACA AYC CAA CAG G-3′)  and mt-fs-h (5′-CCA GTA GAA CAC CCA TTC ATC ATC ATT GGC CAA CTA-3′) . This resulted in an amplicon with a length of 394 bp, corresponding to positions 14,762–15,155 in the mitochondrial genome of Eumeces egregius (positions 633–1026 in the cytochrome b gene of E. egregius) .
DNA sequences were aligned with the ClustalW algorithm implemented in BioEdit , and SNPs were verified by checking the sequence chromatograms. No insertions or deletions were found in the alignment. We deposited new sequences in GenBank under accession numbers MF101923-MF101972 (see Additional file 3). We used MEGA5.2.2  to determine nucleotide diversity and to calculate uncorrected p-distances based on pairwise deletion of ambiguous sites, and DnaSP  to identify haplotypes and determine haplotype diversity. We then divided the dataset by codon position, and ran PartitionFinder 1.1.1  to find the best partition scheme and associated substitution models. The best partition scheme included the three codon positions each as a separate partition (lnL = −1930.14357, BIC = 4553.24905498), with the following models: K80 + G (position 1), HKY + I (position 2), TrN + G (position 3).
We used BEAST v1.8.1 and associated software tools  to obtain time-calibrated estimates of phylogenetic divergence. We used the partitioned dataset, and set the substitution models available in BEAUTi equivalent to the PartitionFinder model selection for the partitions separately (partition 1: HKY with gamma site heterogeneity; partition 2: HKY with invariant site heterogeneity; partition 3: TN93 with gamma site heterogeneity). In accordance with previous recommendations , we used an uncorrelated lognormal relaxed clock model with a constant coalescent tree prior and a random starting tree for all three partitions. Since dated fossils are not available for any Heremites species or related genera (Trachylepis, Chioninia, Eutropis, Mabuya sensu stricto), we instead generated time-calibrated divergence times through estimates of the substitution rate. Substitution rates for cytochrome b in other skinks (Eumeces, Chalcides, Scincus), geckos (Hemidactylus), and lacertid lizards (Lacertini) range from 1.15–1.35% per million years [64, 76–78]. We used a normally distributed prior of the substitution rate (ucld) with an initial mean of 1.25% and a standard deviation of 0.5% for all partitions in BEAUTi, and then optimized these values by checking for convergence and high ESS scores in Tracer v1.6. We then ran five separate rounds of each 107 MCMC iterations and logged parameters every 1000 steps, which generated 50,000 trees. We combined these runs in LogCombiner v1.8.1 with a burn-in of 10% of each run, generated a consensus tree with TreeAnnotator v1.8.1, and produced the final tree with FigTree v1.4.2.
In addition, we calculated a maximum likelihood (ML) tree with the program RAxML 7.0.4  using the rapid hill climbing algorithm . The dataset was partitioned into the three codon positions corresponding to the PartitionFinder partitioning strategy, and run under the suggested  GTR + G substitution model in RAxML.
Sequence evolution within a species does not necessarily follow dichotomous patterns, because co-existence of ancestral and derived haplotypes may introduce reticulate patterns of evolutionary relationships. To account for this phenomenon, we computed a haplotype network with the software NETWORK v4.611  (http://www.fluxus-engineering.com), using the median-joining (MJ) algorithm with default settings (epsilon = 0). We assigned haplotypes to major haplogroups (MHGs) with the Bayesian implementation of the Poisson tree processes (PTP) model , a method to delimit putative species on a phylogenetic tree, using the BEAST tree (Fig. 4) and the PTP default settings (100,000 MCMC generations, thinning 100, burn-in 0.1, seed 123).
We examined 146 H. vittatus specimens from across the distribution range of the species (Additional file 4, Fig. 1b). Specimens were only included in the dataset if the collection locality could be reliably geo-referenced. We did not distinguish among sexes of the examined specimens. For bilateral traits, the left character state was recorded. We measured 60 traits and excluded traits that showed no or only sporadic variation, leaving two binary, eleven meristic and 14 mensural traits for further analyses (Additional file 4). All mensural characters, except for snout-vent length, were first normalized for size by dividing through the snout-vent length (TL, HT, VT, ALL, HL, TH, LH, LF, LT, LA) or the head length (SL, HB, INA), respectively. Data were then ln(x + 1)-transformed. We first conducted a non-metric multidimensional scaling analysis (NMS), which is a statistical approach to reduce the dimensionality of the data set by estimating the similarity of every pairwise comparison of samples based on a similarity coefficient. We used the Gower index (= range-normalized Manhattan index) for meristic and mensural traits, and the Jaccard index for binary characters as defaulted in the program PAST . We retrieved the first three NMS axes and performed a multivariate multiple regression to test for the effect of geographical location on the 3D position in morphospace. We summarized the results with a type 2 MANOVA using Pillai’s test to obtain separate η2 point estimates for the effect of latitude and longitude, and conducted a Bootstrap analysis with 105 repeats to obtain 95% confidence intervals for the η2 point estimates. We restricted subsequent analyses to five meristic and 14 mensural traits that showed substantial variation, and assigned specimens to three geographical groups (western, central, eastern). We conducted pairwise multiple comparison tests for significant mean differences between these groups using Dunnett’s modified Tukey-Kramer test for uneven sample sizes and heterogeneous variance as implemented in the “DTK” R package. We also regressed individual traits separately for each group against the longitude (western and eastern groups) and latitude (central group). All statistical analyses were done in R software v3.2.4 .
The following colleagues generously provided tissue samples or access to museum specimens: Markus Auer, Wolfgang Böhme, Salvador Carranza, Michael Franzen, Philippe Geniez, Gunther Köhler, Nicolà Lutzmann, Petros Lymberakis, Eskandar Rastegar-Pouyani, Boaz Shacham, Franz Tiedemann. Data science specialists Steven Worthington and Simo Goshev at the Institute for Quantitative Social Science, Harvard University, provided statistical support. We thank two anonymous Referees and the Editor for their invaluable suggestions, and Severin Uebbing for comments on an earlier version of the manuscript. David J. Sparrow (Paphos) kindly provided the photographs of Heremites vittatus in Fig. 1.
Availability of data and materials
The DNA sequences generated during the current study are deposited in GenBank under accession numbers MF101923-MF101972.
FB designed the study, generated and analyzed the data, made the figures and wrote the paper. AS helped design the study, analyze the sequence data and write the paper. HSG helped generate the sequence data. MW helped design the study, analyze the sequence data, and write the paper. All authors read and approved the final manuscript.
The authors declare that they have no competing interests.
Consent for publication
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
- Taberlet P, Fumagalli L, Wust-Saucy A-G, Cosson J-F. Comparative phylogeography and postglacial colonization routes in Europe. Mol Ecol. 1998;7:453–64.View ArticlePubMedGoogle Scholar
- Hewitt GM. Post-glacial re-colonization of European biota. Biol J Linn Soc. 1999;68:87–112.View ArticleGoogle Scholar
- Petit RJ, Aguinagalde I, de Beaulieu J-L, Bittkau C, Brewer S, Cheddadi R, Ennos R, Fineschi S, Grivet D, Lascoux M, et al. Glacial refugia: Hotspots but not melting pots of genetic diversity. Science. 2003;300:1563–5.Google Scholar
- Hewitt G. The genetic legacy of the quaternary ice ages. Nature. 2000;405:907–13.View ArticlePubMedGoogle Scholar
- Stewart JR, Lister AM, Barnes I, Dalén L. Refugia revisited: individualistic responses of species in space and time. Proc R Soc B Biol Sci. 2010;277:661–71.View ArticleGoogle Scholar
- Conord C, Gurevitch J, Fady B. Large-scale longitudinal gradients of genetic diversity: a meta-analysis across six phyla in the Mediterranean basin. Ecol Evol. 2012;2:2600–14.View ArticlePubMedPubMed CentralGoogle Scholar
- Kapli P, Lymberakis P, Crochet PA, Geniez P, Brito JC, Almutairi M, Ahmadzadeh F, Schmitz A, Wilms T, Pouyani NR, et al. Historical biogeography of the lacertid lizard Mesalina in North Africa and the Middle East. J Biogeogr. 2015;42:267–79.Google Scholar
- Gvoždik V, Moravec J, Klütsch C, Kotlik P. Phylogeography of the Middle Eastern tree frogs (Hyla, Hylidae, Amphibia) as inferred from nuclear and mitochondrial DNA variation, with a description of a new species. Mol Phylogenet Evol. 2010;55:1146–66.Google Scholar
- Migliore J, Baumel A, Juin M, Médail F. From Mediterranean shores to central Saharan mountains: key phylogeographical insights from the genus Myrtus. J Biogeogr. 2012;39:942–56.View ArticleGoogle Scholar
- Carranza S, Arnold EN, Pleguezuelos JM. Phylogeny, biogeography, and evolution of two Mediterranean snakes, Malpolon monspessulanus and Hemorrhois hippocrepis (Squamata, Colubridae), using mtDNA sequences. Mol Phylogenet Evol. 2006;40:532–46.View ArticlePubMedGoogle Scholar
- Kornilios P, Ilgaz Ç, Kumlutaş Y, Lymberakis P, Moravec J, Sindaco R, Rastegar-Pouyani N, Afroosheh M, Giokas S, Fraguedakis-Tsolis S, et al. Neogene climatic oscillations shape the biogeography and evolutionary history of the Eurasian blindsnake. Mol Phylogenet Evol. 2012;62:856–73.Google Scholar
- Kyriazi P, Poulakakis N, Parmakelis A, Crochet PA, Moravec J, Rastegar-Pouyani N, Tsigenopoulos CS, Magoulas A, Mylonas M, Lymberakis P. Mitochondrial DNA reveals the genealogical history of the snake-eyed lizards (Ophisops elegans and O. occidentalis) (Sauria: Lacertidae). Mol Phylogenet Evol. 2008;49:795–805.Google Scholar
- Husemann M, Schmitt T, Zachos FE, Ulrich W, Habel JC. Palaearctic biogeography revisited: evidence for the existence of a North African refugium for Western Palaearctic biota. J Biogeogr. 2014;41:81–94.View ArticleGoogle Scholar
- Bilgin R. Back to the suture: The distribution of intraspecific genetic diversity in and around Anatolia. Int J Mol Sci. 2011;12:4080–103.View ArticlePubMedPubMed CentralGoogle Scholar
- Rokas A, Atkinson RJ, Webster LMI, Csokas G, Stone GN. Out of Anatolia: longitudinal gradients in genetic diversity support an eastern origin for a circum-Mediterranean oak gallwasp Andriscus quercustozae. Mol Ecol. 2003;12:2153–74.View ArticlePubMedGoogle Scholar
- Kornilios P, Ilgaz Ç, Kumlutaş Y, Giokas S, Fraguedakis-Tsolis S, Chondropoulos B. The role of Anatolian refugia in herpetofaunal diversity: an mtDNA analysis of Typhlops vermicularis Merrem, 1820 (Squamata, Typhlopidae). Amphibia-Reptilia. 2011;32:351–63.Google Scholar
- Feliner GN. Patterns and processes in plant phylogeography in the Mediterranean Basin. A review. Perspect Plant Ecol Evol Syst. 2014;16:265–78.View ArticleGoogle Scholar
- Fady B, Conord C. Macroecological patterns of species and genetic diversity in vascular plants of the Mediterranean basin. Divers Distrib. 2010;16:53–64.View ArticleGoogle Scholar
- Mausfeld P, Vences M, Schmitz A, Veith M. First data on the molecular phylogeography of scincid lizards of the genus Mabuya. Mol Phylogenet Evol. 2000;17:11–4.View ArticlePubMedGoogle Scholar
- Mausfeld P, Schmitz A, Böhme W, Misof B, Vrcibradic D, Rocha CFD. Phylogenetic affinities of Mabuya atlantica Schmidt, 1945, endemic to the Atlantic Ocean archipelago of Fernando de Noronha (Brazil): Necessity of partitioning the genus Mabuya Fitzinger, 1826 (Scincidae: Lygosominae). Zool Anz. 2002;241:281–93.View ArticleGoogle Scholar
- Bauer AM. On the identity of Lacerta punctata Linnaeus 1758, the type specimen of the genus Euprepis Wagler 1830, and the generic assignment of Afro-Malagasy skinks. Afr J Herpetol. 2003;52:1–7.Google Scholar
- Mausfeld P, Schmitz A. Molecular phylogeography, intraspecific variation and speciation of the Asian scincid lizard genus Eutropis Fitzinger, 1843 (Squamata: Reptilia: Scincidae): taxonomic and biogeographic implications. Org Divers Evol. 2003;3:161–71.Google Scholar
- Karin BR, Metallinou M, Weinell JL, Jackman TR, Bauer AM. Resolving the higher-order phylogenetic relationships of the circumtropical Mabuya group (Squamata: Scincidae): An out-of-Asia diversification. Mol Phylogenet Evol. 2016;102:220–32.View ArticlePubMedGoogle Scholar
- Schleich HH, Kästle W, Kabisch K. Amphibians and reptiles of North Africa. Koenigstein: Koeltz Scientific Publishers; 1996.Google Scholar
- Baier F, Sparrow DJ, Wiedl H-J. The amphibians and reptiles of Cyprus. 2nd ed. Frankfurt/M: Edition Chimaira; 2013.Google Scholar
- Anderson SC. The lizards of Iran. Contrib Herpetol. 1999;15:1–442.Google Scholar
- Disi AM, Modrý D, Necas P, Rifai L. Amphibians and reptiles of the Hashemite Kingdom of Jordan. Frankfurt am Main: Edition Chimaira; 2001.Google Scholar
- Leviton AE, Anderson SC, Adler K, Minton SA. Handbook to the Middle East Amphibians and Reptiles, vol. 8. Athens: Society for the Study of Amphibians and Reptiles, Contributions to Herpetology; 1992.Google Scholar
- Bader T, Riegler C, Grillitsch H. The herpetofauna of the Island of Rhodes (Dodecanese, Greece). Herpetozoa. 2009;21:147–69.Google Scholar
- Heremites vittatus, The Reptile Database. http://reptile-database.reptarium.cz/species?genus=Heremites&species=vittatus. Accessed 20 Jan 2017.
- Rastegar-Pouyani N. Intra- and interspecific geographic variation in the Iranian Scincid lizards of the genus Trachylepis Fitzinger 1843 (Sauria: Scincidae). Iran J Anim Biosyst. 2006;2:1–11.Google Scholar
- Budak A. Türkiye’de Mabuya vittata (Scincidae, Lacertilia)‘nın bireysel ve coğrafik variasyonu üzerinde çalışmalar [Studies on individual and geographic variation of Mabuya vittata (Scincidae, Lacertilia) in Turkey]. Ege Üniv Fen Fakültesi İlmi Raporlar Serisi. 1973;162:1–25.Google Scholar
- Tok CV, Göçmen B, Mermer A. Kuzey Kibris Türk Cumhuriyeti Mabuya vittata (Seritli Kertenkele) (Sauria: Scincidae) örnekleri hakkinda [On the specimens of Mabuya vittata (Banded Lizard) (Sauria: Scincidae) from the Turkish Republic of Northern Cyprus]. Turk J Zool. 1999;23:583–9.Google Scholar
- Nassar F, Challita M, Sadek R, Hraoui-Bloquet S. Sexual dimorphism and female reproductive cycle in the scincid lizard Trachylepis vittata (Olivier, 1804) in Lebanon (Reptilia: Scincidae). Zool Middle East. 2013;59:297–301.View ArticleGoogle Scholar
- Rastegar-Pouyani N, Fattahi R. Sexual dimorphism in Trachylepis vittata (Olivier, 1804) (Sauria: Scincidae) in the Zagros Mountains, western Iran. Turk J Zool. 2015;39:59–65.View ArticleGoogle Scholar
- Van der Winden J, Strijbosch H, Bogaerts S. Habitat related disruptive pattern distribution in the polymorphic lizard Mabuya vittata. Acta Oecol. 1995;16:423–30.Google Scholar
- Fattahi R, Ficetola GF, Rastegar-Pouyani N, Avcı A, Kumlutaş Y, Ilgaz Ç, Hosseinian Yousefkhani SS. Modelling the potential distribution of the Bridled Skink, Trachylepis vittata (Olivier, 1804), in the Middle East. Zool Middle East. 2014;60:208–16.Google Scholar
- Camargo A, Sinervo B, Sites JW Jr. Lizards as model organisms for linking phylogeographic and speciation studies. Mol Ecol. 2010;19:3250–70.View ArticlePubMedGoogle Scholar
- De Queiroz K. Species concepts and species delimitation. Syst Biol. 2007;56:879–86.View ArticlePubMedGoogle Scholar
- Leaché AD, Koo MS, Spencer CL, Papenfuss TJ, Fisher RN, McGuire JA. Quantifying ecological, morphological, and genetic variation to delimit species in the coast horned lizard species complex (Phrynosoma). Proc Natl Acad Sci. 2009;106:12418–23.Google Scholar
- Le Houérou HN. Climate, flora and fauna changes in the Sahara over the past 500 million years. J Arid Environ. 1997;37:619–47.View ArticleGoogle Scholar
- Kornilios P, Kyriazi P, Poulakakis N, Kumlutaş Y, Ilgaz Ç, Mylonas M, Lymberakis P. Phylogeography of the ocellated skink Chalcides ocellatus (Squamata, Scincidae), with the use of mtDNA sequences: A hitch-hiker’s guide to the Mediterranean. Mol Phylogenet Evol. 2010;54:445–56.Google Scholar
- Ahmadzadeh F, Flecks M, Rödder D, Böhme W, Ilgaz Ç, Harris DJ, Engler JO, Üzüm N, Carretero MA. Multiple dispersal out of Anatolia: biogeography and evolution of oriental green lizards. Biol J Linn Soc. 2013;110:398–408.Google Scholar
- Güçlü Ö, Candan K, Kankiliç T, Kumlutaş Y, Durmuş SH, Poulakakis N, Ilgaz Ç. Phylogeny of Trachylepis sp. (Reptilia) from Turkey inferred from mtDNA sequences. Mitochondrial DNA. 2014;25:456–63.Google Scholar
- Kosswig C. Zoogeography of the Near East. Syst Zool. 1955;4:49–73.View ArticleGoogle Scholar
- Blondel J, Aronson J. Biology and wildlife of the Mediterranean Region. Oxford: Oxford University Press; 1999.Google Scholar
- Willis KJ, Kleczkowski A, Crowhurst SJ. 124,000-year periodicity in terrestrial vegetation change during the late Pliocene epoch. Nature. 1999;397:685–8.View ArticleGoogle Scholar
- Webb T, Bartlein PJ. Global Changes During the Last 3 Million Years: Climatic Controls and Biotic Responses. Annu Rev Ecol Syst. 1992;23:141–73.View ArticleGoogle Scholar
- Zachos J, Pagani M, Sloan L, Thomas E, Billups K. Trends, rhythms, and aberrations in global climate 65 Ma to present. Science. 2001;292:686–93.View ArticlePubMedGoogle Scholar
- Bennett KD. Milankovitch cycles and their effects on species in ecological and evolutionary time. Paleobiology. 1990;16:11–21.View ArticleGoogle Scholar
- Caujapé-Castells J, Jansen RK, Membrives N, Pedrola-Monfort J, Montserrat JM, Ardanuy A. Historical biogeography of Androcymbium Willd. (Colchicaceae) in Africa: evidence from cpDNA RFLPs. Bot J Linn Soc. 2001;136:379–92.Google Scholar
- Hadjisterkotis E, Masala B, Reese DS. The origin and extinction of the large endemic Pleistocene mammals of Cyprus. Biogeographia. 2000;XXI:593–606.Google Scholar
- Poulakakis N, Kapli P, Kardamaki A, Skourtanioti E, Göcmen B, Ilgaz C, Kumlutas Y, Avci A, Lymberakis P. Comparative phylogeography of six herpetofauna species in Cyprus: late Miocene to Pleistocene colonization routes. Biol J Linn Soc. 2013;108:619–35.Google Scholar
- Avise JC, Arnold J, Ball RM, Bermingham E, Lamb T, Neigel JE, Reeb CA, Saunders NC. Intraspecific phylogeography: The mitochondrial DNA bridge between population genetics and systematics. Annu Rev Ecol Syst. 1987;18:489–522.Google Scholar
- Pouyani ER, Noureini SK, Pouyani NR, Joger U, Wink M. Molecular phylogeny and intraspecific differentiation of the Eremias velox complex of the Iranian Plateau and Central Asia (Sauria, Lacertidae). J Zool Syst Evol Res. 2012;50:220–9.Google Scholar
- Joger U, Fritz U, Guicking D, Kalyabina-Hauf S, Nagy ZT, Wink M. Phylogeography of western Palaearctic reptiles – Spatial and temporal speciation patterns. Zool Anz. 2007;246:293–313.View ArticleGoogle Scholar
- Guicking D, Lawson R, Joger U, Wink M. Evolution and phylogeny of the genus Natrix (Serpentes: Colubridae). Biol J Linn Soc. 2006;87:127–43.View ArticleGoogle Scholar
- Rastegar-Pouyani E, Rastegar-Pouyani N, Noureini SK, Joger U, Wink M. Molecular phylogeny of the Eremias persica complex of the Iranian plateau (Reptilia: Lacertidae), based on mtDNA sequences. Zool J Linnean Soc. 2010;158:641–60.Google Scholar
- Ferchaud AL, Eudeline R, Arnal V, Cheylan M, Pottier G, Leblois R, Crochet PA. Congruent signals of population history but radically different patterns of genetic diversity between mitochondrial and nuclear markers in a mountain lizard. Mol Ecol. 2015;24:192–207.Google Scholar
- Rato C, Carranza S, Perera A, Carretero MA, Harris DJ. Conflicting patterns of nucleotide diversity between mtDNA and nDNA in the Moorish gecko, Tarentola mauritanica. Mol Phylogenet Evol. 2010;56:962–71.Google Scholar
- Rato C, Carranza S, Harris DJ. When selection deceives phylogeographic interpretation: The case of the Mediterranean house gecko, Hemidactylus turcicus (Linnaeus, 1758). Mol Phylogenet Evol. 2011;58:365–73.Google Scholar
- Luke C. Convergent evolution of lizard toe fringes. Biol J Linn Soc. 1986;27:1–16.View ArticleGoogle Scholar
- Carranza S, Arnold EN. Investigating the origin of transoceanic distributions: mtDNA shows Mabuya lizards (Reptilia, Scincidae) crossed the Atlantic twice. Syst Biodivers. 2003;1:275–82.View ArticleGoogle Scholar
- Miralles A, Carranza S. Systematics and biogeography of the Neotropical genus Mabuya, with special emphasis on the Amazonian skink Mabuya nigropunctata (Reptilia, Scincidae). Mol Phylogenet Evol. 2010;54:857–69.Google Scholar
- Miralles A, Fuenmayor GR, Bonillo C, Schargel WE, Barros T, García-Perez JE, Barrio-Amorós CL. Molecular systematics of Caribbean skinks of the genus Mabuya (Reptilia, Scincidae), with descriptions of two new species from Venezuela. Zool J Linnean Soc. 2009;156:598–616.Google Scholar
- Sambrook J, Russell DW. Molecular Cloning. A laboratory manual. New York: Cold Spring Harbor Laboratory Press; 2001.Google Scholar
- Fritz U, Siroký P, Kami H, Wink M. Environmentally caused dwarfism or a valid species--Is Testudo weissingeri Bour, 1996 a distinct evolutionary lineage? New evidence from mitochondrial and nuclear genomic markers. Mol Phylogenet Evol. 2005;37:389–401.Google Scholar
- Rastegar-Pouyani E. Molecular phylogeography and evolution of the lacertid genus Eremias of the Iranian Plateau and Central Asia. (Reptilia, Lacertidae). Dissertation. Heidelberg: University of Heidelberg; 2007.Google Scholar
- Kumazawa Y, Nishida M. Complete mitochondrial DNA sequences of the Green Turtle and Blue-Tailed Mole Skink: Statistical evidence for archosaurian affinity of turtles. Mol Biol Evol. 1999;16:784–92.Google Scholar
- Hall TA. BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT. Nucleic Acids Symp Ser. 1999;41:95–8.Google Scholar
- Tamura K, Peterson D, Peterson N, Stecher G, Nei M, Kumar S. MEGA5: Molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Mol Biol Evol. 2011;28:2731–9.View ArticlePubMedPubMed CentralGoogle Scholar
- Librado P, Rozas J. DnaSP v5: A software for comprehensive analysis of DNA polymorphism data. Bioinformatics. 2009;25:1451–2.View ArticlePubMedGoogle Scholar
- Lanfear R, Calcott B, Ho SY, Guindon S. Partitionfinder: combined selection of partitioning schemes and substitution models for phylogenetic analyses. Mol Biol Evol. 2012;29:1695–701.View ArticlePubMedGoogle Scholar
- Drummond AJ, Suchard MA, Xie D, Rambaut A. Bayesian phylogenetics with BEAUti and the BEAST 1.7. Mol Biol Evol. 2012;29:1969–73.Google Scholar
- Drummond AJ, Ho SYW, Phillips MJ, Rambaut A. Relaxed Phylogenetics and Dating with Confidence. PLoS Biol. 2006;4:e88.View ArticlePubMedPubMed CentralGoogle Scholar
- Carranza S, Arnold EN, Geniez P, Roca J, Mateo JA. Radiation, multiple dispersal and parallelism in the skinks, Chalcides and Sphenops (Squamata: Scincidae), with comments on Scincus and Scincopus and the age of the Sahara Desert. Mol Phylogenet Evol. 2008;46:1071–94.View ArticlePubMedGoogle Scholar
- Carranza S, Arnold EN. Systematics, biogeography, and evolution of Hemidactylus geckos (Reptilia: Gekkonidae) elucidated using mitochondrial DNA sequences. Mol Phylogenet Evol. 2006;38:531–45.View ArticlePubMedGoogle Scholar
- Carranza S, Arnold EN, Amat F. DNA phylogeny of Lacerta (Iberolacerta) and other lacertine lizards (Reptilia: Lacertidae): did competition cause longterm mountain restriction? Syst Biodivers. 2004;2:57–77.Google Scholar
- Stamatakis A. RAxML-VI-HPC: maximum likelihood-based phylogenetic analyses with thousands of taxa and mixed models. Bioinformatics. 2006;22:2688–90.View ArticlePubMedGoogle Scholar
- Stamatakis A, Blagojevic F, Nikolopoulos DS, Antonopoulos CD. Exploring New Search Algorithms and Hardware for Phylogenetics: RAxML Meets the IBM Cell. J VLSI Signal Proc Syst Signal Image Video Technol. 2007;48:271–86.View ArticleGoogle Scholar
- Bandelt H-J, Forster P, Röhl A. Median-joining networks for inferring intraspecific phylogenies. Mol Biol Evol. 1999;16:37–48.View ArticlePubMedGoogle Scholar
- Zhang J, Kapli P, Pavlidis P, Stamatakis A. A general species delimitation method with applications to phylogenetic placements. Bioinformatics. 2013;29:2869–76.View ArticlePubMedPubMed CentralGoogle Scholar
- Hammer Ø, Harper DAT, Ryan PD. PAST: Paleontological statistics software package for education and data analysis. Palaeontol Electron. 2001;4:1–9.Google Scholar
- R Core Team. R: A language and environment for statistical computing. Vienna: R Foundation for Statistical Computing; 2016.Google Scholar