Figures
Abstract
Harboring a large number of endemic species, the Alps and the Western Carpathians are considered as major centers of biodiversity. Nonetheless, the general opinion until the turn of the millennium was that both Central European mountain regions did not provide suitable habitat during the Last Glacial Maximum, but were colonized later from southern refuges. However, recent molecular genetic studies provide new evidence for peripheral Alpine refuges. We studied the phylogeography of the calciphilous land snail O. dolium across its distribution in the Alps and the Western Carpathians to assess the amount of intraspecific differentiation and to detect potential glacial refuges. A partial sequence of the mitochondrial COI was analyzed in 373 specimens from 135 sampling sites, and for a subset of individuals, partial sequences of the mitochondrial 16S and the nuclear histone H3 and H4 were sequenced. A molecular clock analysis was combined with a reconstruction of the species’ geographic range history to estimate how its lineages spread in the course of time. In order to obtain further information on the species’ past distribution, we also screened its extensive Pleistocene fossil record. The reconstruction of geographic range history suggests that O. dolium is of Western Carpathian origin and diversified already around the Miocene-Pliocene boundary. The fossil record supports the species’ presence at more than 40 sites during the last glacial and earlier cold periods, most of them in the Western Carpathians and the Pannonian Basin. The populations of O. dolium display a high genetic diversity with maximum intraspecific p-distances of 18.4% (COI) and 14.4% (16S). The existence of various diverged clades suggests the survival in several geographically separated refuges. Moreover, the sequence patterns provide evidence of multiple migrations between the Alps and the Western Carpathians. The results indicate that the Southern Calcareous Alps were probably colonized only during the Holocene.
Citation: Harl J, Duda M, Kruckenhauser L, Sattmann H, Haring E (2014) In Search of Glacial Refuges of the Land Snail Orcula dolium (Pulmonata, Orculidae) - An Integrative Approach Using DNA Sequence and Fossil Data. PLoS ONE 9(5): e96012. https://doi.org/10.1371/journal.pone.0096012
Editor: Donald James Colgan, Australian Museum, Australia
Received: December 13, 2013; Accepted: April 1, 2014; Published: May 7, 2014
Copyright: © 2014 Harl et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Funding: This project was funded by the Austrian Science Fund (FWF, Proj.-No. 19592-B17) and the “Freunde des Naturhistorischen Museums Wien”. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
Introduction
The Pleistocene climate changes shaped the phylogeographic patterns of various organisms [1]. In particular, the severe cooling starting with the end of the Early Pleistocene (about 900 kya) was the starting point for massive glaciations in the northern hemisphere [2]. Mountainous regions such as the Central European Alps were heavily affected due to shifts in temperature and humidity, and the expansion of glaciers, resulting for many taxa in the fragmentation of populations, complete or local extinction, and the loss of variation due to genetic bottlenecks [1]. Thereby, the Last Glacial Maximum (LGM; 30–18 kya [3]) is of most relevance in respect to the current distribution of Central European species. Although the existence of glacial refuges at the periphery of the Alps was discussed more than half a century ago [4], the general opinion until end of the 20th century was that glacial refuges were located mainly in southern regions [1]. However, molecular genetic analyses and fossil data revealed the existence of northern refuges in the Western Carpathians [5], [6] and the Pannonian Basin [7]. Several peripheral Alpine refuges were proposed for silicophilous mountain plants [8] and calciphilous land snails such as Arianta arbustorum [9]–[12], Carychium minimum, Carychium tridentatum [13], Trochulus oreinos [14] and Trochulus villosus [15].
Most molecular genetic studies investigating glacial refuges were based on the assumption that populations diverged during isolation in geographically separated areas, and that populations of former refuge areas are now characterized by high genetic diversity and the presence of rare (private) alleles [16]. Nevertheless, a major obstacle in identifying signatures of Pleistocene refuges is that phylogenetic signals are blurred because of migration and intermixture of previously separated populations. Therefore, species showing low dispersal capabilities and specific habitat requirements, which is the case in most land snails, might be suited to infer past distributional patterns.
In the present study, we investigate the phylogeographic patterns of the land snail Orcula dolium (Draparnaud, 1801). The species inhabits all major limestone areas of the Alps and the neighboring Western Carpathians [17], [18]. O. dolium is usually associated with mountainous forest habitats or rocky landscapes with patches of vegetation, but it is also found on rocky slopes at high altitudes up to 2160 m above sea level (asl) [17] (and data of present study). Within the Alps, the east-west orientated Central Alps, consisting mainly of silicate rock, represent the largest distributional barrier for the species; they separate the Northern Calcareous Alps from the Southern Calcareous Alps. The Vienna Basin, closing the Pannonian Basin to the north, constitutes another distributional gap because it separates the populations of the Northern Calcareous Alps (including the Wienerwald) and the Western Carpathians. Recent investigations of loess profiles from the Pannonian Basin show that during the Late Pleistocene O. dolium also occurred in the lowlands of the region [19], although the species seems to have vanished from it during the Holocene. The margins of the Northern Calcareous Alps and the Western Carpathians, enclosing the Vienna Basin east and west, harbor the majority of the described 23 subspecies [20]. As these areas were only partially glaciated during the LGM, they also come into consideration as potential refuges for the species.
We perform a comprehensive phylogeographic study of O. dolium analyzing mt and nc markers to detect potential glacial refuges and to assess the amount of intraspecific variability. Furthermore, we test whether the described subspecies are differentiated genetically in the mitochondrial (mt) and nuclear (nc) markers. Our sample covers almost the entire range of the species. We include data of the Pleistocene fossil record to obtain insights into the species’ past distribution. In order to estimate which potential areas were inhabited by ancestral populations of O. dolium and to trace the distribution patterns of the mt lineages throughout time, a molecular clock analysis is performed and combined with a phylogeographic range reconstruction.
Methods
Study Area and Sampling
Specimens were collected from a large part of the species’ distribution, including several Alpine and Western Carpathian type localities. O. dolium is not protected by conservation laws of the countries where the collections were performed. Thus, in general, permissions were not necessary. For protected areas in Austria, permissions were provided by federal states authorities. Permit numbers: RU5-BE-64/011-2013 (Lower Austria), FA13C-53 Sch 6/6–2007 (Styria) and N10-117-2008 (Upper Austria). Most samples of O. dolium were collected in the Northern Calcareous Alps (Austria and Germany), a lesser fraction in the Western Carpathians (Slovakia), the Southern Calcareous Alps (Austria, Slovenia and Italy) and the Western Alps (Switzerland), totaling 373 specimens of 135 sites (Table 1). The habitats include wooded areas in the lowland, mountainous vegetated areas and rocky slopes in the alpine zones, with an altitudinal range from 120 m to 2160 m asl. Elevation and position were determined via GPS. At every sample locality, if available, a minimum of three living specimens was collected, prepared for DNA analyses and stored in 80% ethanol following the protocol of [21]. Orcula conica (Rossmässler, 1837) (ID: 3899; Trögener Klamm, Carinthia) was used as outgroup. Selected type specimens of O. dolium subspecies included in the present study are shown in Fig. S1. For the inference of substitution rates used in the molecular clock analysis, six specimens of Orculella bulgarica (Hesse, 1915) and four specimens of Orculella aragonica (Westerlund, 1897) were included. Samples of Orculella bulgarica were collected by Barna Páll-Gergely (Shinshu University, Matsumoto, Japan; Turkish samples) and Alexander Reischuetz (Greek samples). DNA samples of O. aragonica were obtained from Benjamín Gómez-Moliner (Universidad del País Vasco, Vitoria, Spain). Voucher specimens of the first three taxa are deposited in the Natural History Museum Vienna, the whereabouts of the O. aragonica vouchers are provided in [22]. Every individual sample was assigned an ID consisting of a unique specimen number and a locality tag. The latter encodes the Alpine geographic mountain region as classified in the SOIUSA system [23], with localities of each region numbered from west to east (Table 1). Due to their geographic vicinity to the Northern Calcareous Alps, the sites located in the Fischbach Alps, Eastern Styrian Prealps and Lavanttal Alps as well as those of the Wienerwald, which are classified to the Central Alps in the SOIUSA system, are treated as Northern Calcareous Alpine mountain regions here. A map providing an overview of the mountain regions investigated is shown in Fig. 1. The Slovakian sites were each assigned to one of the geological areas defined for the Carpathians [24]. To illustrate distribution patterns, the haplotypes in the phylogenetic trees and the histone network are marked by different colors, corresponding to the SOIUSA mountain regions as shown in Fig. 1.
The colors correspond to those used in the mt trees (Figs. 1 and 4) and in the histone H4/H3 network (Fig. 5). The outlines of the Northern Calcareous Alps (NCA), the Southern Calcareous Alps (SCA), the Western Alps (WA) and the Central Alps (CA) are framed in black. The names of the mountain regions and abbreviations are provided in the figure.
Molecular Markers and Primer Design
DNA was extracted from 373 O. dolium specimens with the QIAgen Blood and Tissue Kit, using a piece of foot tissue. A partial region of the mt cytochrome c oxidase I gene (COI) was sequenced in all 373 specimens. From a subset of 54 individuals (including representatives of the major mt clades, type localities, peripheral geographic regions and contact zones of distinct mt clades), two additional markers were amplified: a section of the mt 16S gene, as well as a nc sequence comprising almost the entire sequences of the histone 4 and 3 genes and the complete internal spacer region (H4/H3). The two histone genes are orientated in opposite direction and are separated by a non-coding spacer region, an arrangement which is probably not universal in gastropods but was verified so far particularly for species of the informal group of Orthurethra (sensu [25]). The COI and 16S sequences were also amplified in the eleven outgroup specimens. The COI forward primer COIfolmerFwd [14] is a variant of the standard primer LCO1490 [26]; H2198-Alb [10] was used as reverse primer. New 16S primers were designed for the amplification of a fragment of approximately 850 bp. The forward primers, 16SLOrc1_fw and 16SLOrc2_fw, bind about 50 bp away from the 5′-end of the 16S rRNA gene and the reverse primer 16SLOrc_rv is situated in a conserved region about 850 bp downstream. The H4/H3 primers, OrcH4_left1 and OrcH4_left2 (positioned at the 3′-end of the H4 gene) were designed based on alignments of H3 and H4 sequences from GenBank. The reverse primer OrcH3_right1 (at the 3′-end of the H3 gene) was published by [27]. Internal primers for sequencing (OrcH3S_right3 and OrcH4S_left3) were positioned close to the spacer in the coding sequences to obtain the complete 1100 bp fragment with two sequencing runs. The PCR primers (OrcH4_left1, OrcH4_left2 and OrcH3_right1) cover a wider spectrum of Orthurethra taxa, while the internal primers (OrcH3S_right3 and OrcH4S_left3) were especially adapted to the genus Orcula. All primers are listed in Table 2.
PCR Amplification and Cloning
COI and 16S fragments were amplified with the Roche Taq DNA polymerase for direct sequencing. The PCR started with a denaturation step for 3 min at 94°C, followed by 35 cycles with 30 s at 94°C, 30 s at the particular annealing temperature (Table 2), and 1 min at 72°C, followed by a final extension for 7 min at 72°C. The PCR for the nc H4/H3 fragments was performed with the standard protocol of the Finnzymes Phusion polymerase, which has proofreading activity. PCR started with a denaturation step for 30 s at 98°C, followed by 35 cycles with 10 s at 98°C, 10 s at the particular annealing temperature (see Table 2), and 30 s at 72°C, followed by a final extension for 7 min at 72°C. Purification and sequencing of the PCR products (both directions) were performed at LGC Genomics (Berlin, Germany), using the PCR primers for sequencing, except for the H4/H3 section, which was sequenced with the internal primers OrcH3S_right3 and OrcH4S_left3. The cloning of PCR products was performed for the H4/H3 primer design phase and for those individuals that proved to be heterozygous with respect to insertions and deletions (indels), resulting in varying spacer lengths and thus impeding direct sequencing. Nine of the 14 specimens with heterozygous sequences yielded fragments differing in spacer length. The fragments were purified using the QIAquick Gel Extraction Kit (QIAGEN), extended by A-endings with the DyNAzyme II DNA polymerase (Finnzymes) and cloned with the TOPO TA cloning kit (Invitrogen). Two to four clones were sequenced until the two main variants were obtained. In two cases (samples WIW2_392 and MLF1_3939) more than two length variants were obtained, which is not completely unexpected in multi-copy gene families. Sequencing was also performed at LGC Genomics (Berlin) using M13 universal primers. All sequences are deposited in GenBank under the following accession numbers: KC568830–KC569204, KJ656162–KJ656172 (COI); KC569205–KC569260, KJ656173–KJ656183 (16S); KC569261–KC569327 (H4/H3).
Sequence Analysis and Phylogenetic Tree Calculation
Sequences were edited using Bioedit 7.1.3 [28]. When directly sequenced H4/H3 fragments provided ambiguous positions, the respective sites were filled with the corresponding IUPAC codes. The complete COI data set comprises 374 sequences (including a single specimen of O. conica as outgroup). The alignment of the 655 nucleotide sites was straightforward because there were no indels. Statistical analyses were performed using all COI sequences. Haplotype and nucleotide diversity were calculated with DnaSP 5.10 [29], and uncorrected genetic p-distances between the clades of O. dolium and O. conica were calculated with MEGA 5.2 [30]. For phylogenetic tree calculation, identical haplotypes from the same geographic areas (SOIUSA codes) were collapsed resulting in a total of 220 COI sequences of O. dolium (197 haplotypes). Prior to the phylogenetic tree inference, a search for the best fitting substitution model was performed with jModelTest 0.1 [31]. A Bayesian inference (BI) was calculated with MrBayes 3.2.2 [32], [33] for 5×106 generations (samplefreq = 100; nruns = 2; nchains = 4), applying the parameters obtained from the model test (GTR+G+I; nst = 6, rates = invgamma). Tracer 1.5 [34] was used to assess whether the two runs had converged and when the stationary phase was reached. The first 25% of the trees were discarded as burnin and a 50% majority rule consensus tree was calculated from the remaining trees.
The 16S sequences (55 specimens including O. conica) were aligned with ClustalX [35] using default parameters. The original alignment contained 879 positions of which all 61 gap positions were excluded with TrimAl 1.3 [36], implemented in the Phylemon 2.0 web tools [37], using the “no gap” option. Another 84 sites were excluded by performing the “strict” option, leaving 734 positions in the final alignment for the phylogenetic tree analyses. Of those, 217 sites were variable, compared to 331 in the original and 298 sites in the “gaps excluded” alignment. The “no gap” alignment was also used for calculation of uncorrected p-distances. The BI was performed with two data partitions (COI and 16S “strict”), using the substitution models suggested by jModelTest 0.1 [31] (COI: HKY+I+G: nst = 2, rates = invgamma; 16S “strict”: HKY+I+G: nst = 2, rates = invgamma), and allowing MrBayes to evaluate the model priors independently. A Maximum Likelihood (ML) tree was calculated with MEGA 5.2 [30], applying the sequence evolution model GTR+G (5 rate categories)+I and performing 1000 bootstrap replicates with SPR (Subtree-Pruning-Redrafting) as heuristic method for tree inference. Support values of the ML analysis are provided for the combined COI and 16S tree.
A Median-Joining network was calculated from the H4/H3 data with Network 4.6.0.0 (Fluxus Technology Ltd.) using the default settings. Two networks were calculated: one without gaps and one keeping the gap positions. The networks were then post-processed to reduce unnecessary median vectors using the MP option. As the gap positions contained valuable information, we show the network including the gap sites.
Molecular Clock Analysis and Geographic Range Reconstruction
A reconstruction of the historical biogeography of O. dolium was performed with Lagrange 2.0.1 [38]. Lagrange 2.0.1 uses a dispersal-extinction-cladogenesis (DEC) modeling, which allows analyzing the ML values of rate transitions as a function of time. The calculations were based on a molecular clock dated linearized BI tree calculated with BEAST 1.7.5 [39]. The fossil record of O. dolium comprises only material of the Holocene and the Middle and Late Pleistocene. This might be due to the fact that during these time periods climate conditions promoted the accumulation of loess sediments containing high numbers of gastropod shells. In contrast, fossil records of the Early Pleistocene and the Pliocene are very scarce and lack many land snail species endemic to Central Europe, among them O. dolium. Consequently, dating the stem of the tree with the species first occurrence in the fossil record was not a reasonable option. Hence, substitution rates were inferred from a comparison of two other Orculidae, Orculella aragonica and O. bulgarica. The latter is widespread in South-Eastern Europe and Western Asia, whereas its sister species, O. aragonica, is distributed only in the Iberian Peninsula [22]. The earliest fossil record of O. aragonica/O. bulgarica dates from about 1.8 mya at the Almenara-Casablanca karst complex (Castellón, Eastern Spain), which contains Miocene to Early Pleistocene sediments [40], [41]. We assumed that the earliest record coincides with the time period when the ancestral lineages of the two species split from each other. Calculation of substitution rates and molecular clock analysis were both performed with the same COI and 16S (“strict”, 696 positions) alignments (including the sequences of O. conica and the two Orculella species). For each partition, TN93+G was applied both to the measurement of substitution rates and the inference of the molecular clock dated tree. Molecular clock tests were performed independently for the COI and 16S alignments with MEGA 5.2 under the TN93 model, using a discrete Gamma (G) distribution to model differences in evolutionary rates among sites. The null hypotheses of equal evolutionary rates throughout the trees were rejected at a 5% significance level. Substitution rates (COI: 0.02333, 16S: 0.02528; substitutions/ma) were assigned in the prior settings of BEAUti 1.7.5 (part of the BEAST package), and uncorrelated relaxed lognormal molecular clocks were implemented for both sequence partitions. “Speciation: Yule Process” was chosen as tree prior. The BEAST analysis was performed with four independent runs with 107 generations each (sample freq.: 1000). Tracer 1.5 [34] was used to assess whether the four runs had converged. The four independent runs were then combined with LogCombiner 1.7.5 (part of the BEAST package). Subsequently, 25% of the trees were discarded as burnin and a 50% majority rule consensus tree was calculated from the remaining 3×105 trees. Median node heights and 95% highest posterior density (95% HPD) intervals are provided for major nodes in the results section. The rate-calibrated linearized tree was then prepared for Lagrange configurator 20130526, together with a range matrix in which each lineage was assigned either to the Northern Calcareous Alps, the Southern Calcareous Alps, the Western Alps or the Western Carpathians. Migration was permitted between all regions, but lower probabilities (‘0.5’ instead of ‘1.0’) were assigned in the dispersal constraints for migration between areas not being immediately adjacent (Southern Calcareous Alps and Western Carpathians; Southern Calcareous Alps and Western Alps; Western Carpathians and Western Alps). The ancestors were allowed to occupy a maximum of two geographic areas. Alternative analyses were run with unlimited range sizes (allowing the taxa to inhabit all four geographic areas) and assigning the same probabilities (‘1.0’) for migration between the four areas in the dispersal constraints.
Literature Search for Fossil Records of O. dolium
We screened various papers for fossil records of O. dolium. The dating of the Alpine sites seems to be rather tentative because all but one site were dated only using reconciliation with vertebrate fossils of the same layers. Most reliable is the stratigraphic dating of recently investigated loess profiles of the Pannonian Basin, including measurements of carbonate content variations, low-field magnetic susceptibility and radiocarbon dating using macro charcoal fragments. Most of the Western Carpathian records are based on investigations of soil profiles as well. The publications with positive records of O. dolium are listed according to the countries investigated: Austria [42], [43], Croatia [44], [45], Czech Republic [46]–[48], Germany [49], Hungary [19], [50], Serbia [51]–[53] and Slovakia [47], [48]. Since the evaluation of the fossil record of Orcula dolium was based on literature data only, no permits were required for that part of the study. We assigned each record to a time period of either cold or warm climate phases of Middle and Late Pleistocene according to the Quaternary divisions of the North European climate cycles of Zagwijn [54]: the Weichselian (115–11 ka ago), Saalian (350–130 ka ago) and Elsterian (475–370 ka ago) were considered as cold climate stages (glacials), compared to the warm stages (interglacials), Holocene (11 ka ago - present), Weichselian-Saalian interglacial ( = Eemian; 130–115 ka ago) and Saalian-Elsterian interglacial ( = Holsteinian; 370–350 ka ago). However, it has to be noted that massive glaciations of the Alps and northern Europe occurred only at some times of the glacials, most recent during the LGM (30–18 kya). The maps showing the distribution of genetic clades and fossil records were prepared using ArcMap Desktop 10.0 and manually edited in Adobe Photoshop CS4 version 12.
Results
Mitochondrial Clades
Among the 373 individuals of O. dolium investigated, 197 COI haplotypes are observed. The phylogenetic trees calculated with different algorithms show similar topologies. Ten major clades, some of them divided into distinct sub-clades, are highly supported in all analyses. Six of the major clades occur exclusively in the Western Carpathians (2, 3, 5, 7, 9 and 10), three clades are distributed solely in the Alpine region (4, 6 and 8). The final clade, clade 1, is subdivided into four sub-clades themselves restricted to particular regions (Alps: 1A; Western Carpathians: 1B, 1C and 1D) (Figs. 2, 3 and 4). Sub-clade 1A is distributed throughout the Northern Calcareous Alps but occurs also in a distinct area of the Southern Calcareous Alps. Similarly, clade 6 is distributed in separate areas of the Northern Calcareous Alps and the Southern Calcareous Alps. Clade 4 also has a disjunct distribution, but in an east-west orientation, occurring in the eastern-most part of the Northern Calcareous Alps (4B) and in the Western Alps (4A). In contrast, clade 8 is geographically restricted to a small area in the Northern Calcareous Alps. The BI tree based on COI and 16S reveals the same highly supported clades, sub-clades and overall topology as the COI tree, but some nodes are better supported. Nevertheless, the topology remains partly ambiguous, e.g., the relationships between clades 1, 2 and 3 as well as between 4, 5, and 6 are unresolved trichotomies.
The central tree provides an overview. Posterior probabilities are given at the main nodes. Black dots indicate nodes with maximum support. A specimen of O. conica is used as outgroup. More details are provided in the partial trees with different colors corresponding to the geographic mountain regions depicted in Fig. 1. The scale bars indicate the expected number of substitutions per site according to the model of sequence evolution applied.
The symbols indicate the presence of the respective mt clades. If specimens of different mt clades are present at the same localities, the respective clade symbols are shown encircled with a line pointing towards the locality. A picture of the shell of a type specimen of O. dolium (Syntype NHMW 14765/1820.26.61/2; size: 6.8 mm in height) is provided in the left upper corner.
Posterior probabilities and ML bootstrap values are provided for all nodes above sub-clade level. The scale bar indicates the expected number of substitutions per site according to the model of sequence evolution applied. The dots indicate nodes with posterior probabilities of 1.0 and bootstrap values of more than 98% (black) or 95% to 97% (grey) in the ML analyses, respectively.
Nuclear Data Set
The H4/H3 data is shown as a phylogenetic network (Fig. 5). It is not possible to incorporate the division of the mt clades into the H4/H3 network because the pattern differs considerably from that obtained in the mt trees - specimens of the same mt clades often provide very different H4/H3 alleles. The general scheme of the network exhibits three frequent alleles (HT1, HT2 and HT3), each encircled by several similar haplotypes, which differ by a few substitutions or indels only. Additionally, there are several unique alleles, particularly in specimens from the Western Carpathians (e.g. in 1996_STR2, 1376_MLF5, 3952_MLF6, 3915_MLF3 and 1373_MLF7). In the Alpine populations, HT1 is the most frequent allele and, apart from a unique Western Alpine allele (differing by 2 substitutions from HT1), it is the only one found in the Western Alps and the Southern Calcareous Alps. In the Western Carpathians, HT1 is detected in a single specimen (STR1_3909) only. However, HT1 is encircled by several slightly differing haplotypes of both the Northern Calcareous Alps and the Western Carpathians. In contrast, HT2 and HT3 are each found exclusively in either the Northern Calcareous Alps or the Western Carpathians, respectively. Thus, all populations having HT2 additionally feature HT1 or haplotypes slightly deviating from HT1. Within the Alps, the populations of the eastern-most margin of the Northern Calcareous Alps exhibit the largest allele diversity.
Bars indicate substitutions within the H3, H4 and the spacer region. The loss or gain of bases (indels) in the spacer region is displayed by arrows pointing in the respective direction and the numbers indicate in how many bases the haplotypes/alleles differ from each other. The three most common alleles are named as HT1, HT2 and HT3. The size of the circles corresponds to the number of sequences providing the same allele. The colors correspond to one of the mountain regions defined in Fig. 1. To facilitate the comparison of nc and mt data, the specimen labels and the clade affiliation are indicated next to the haplotype circles.
Genetic Diversity
Intraspecific distances are extraordinarily high for the COI gene, with up to 18.4% mean p-distance between the clades (Table 3). The maximum p-distance within the Alpine populations is 16.9%, compared to 18.3% in the Western Carpathian ones. The mean p-distance between O. dolium and O. conica is only slightly higher at 18.4% (Table 3). The respective haplotype and nucleotide diversities of the COI clades are high in all populations (Table 4). Distances were also calculated for the 16S sequences (Table 5). The largest mean p-distance between clades of O. dolium is 14.4%, the distance between O. dolium and O. conica is 19.3%, indicating that the conserved parts of the 16S suffered less from saturation than COI. Regarding the nc H4/H3 sequences, the Western Carpathian populations almost consistently show larger p-distances (max.: H3 1.5%, H4 2.7%, Spacer 2.6%; clones excluded: H3 0.3%, H4 1.6%, Spacer 1.0%) than the Alpine ones (max.: H3 0.6%, H4 1.9%, Spacer 2.1%; clones excluded: H3 0.3%, H4 1.2%, Spacer 1.5%). The max. p-distances within O. dolium (H3 1.5%, H4 2.7%, Spacer 2.6%; clones excluded: H3 0.3%, H4 1.9%, Spacer 1.8%) are comparable to the mean distances between O. dolium and O. conica (H3 0.67%, H4 1.30%, Spacer 2.38%).
Molecular Clock Analysis and Reconstruction of the Phylogeographic Range Evolution
The linearized tree resulting from the BEAST 1.7.5 analysis was combined with a biogeographical range reconstruction using Lagrange v.20130526 (Fig. 6). For the nodes marking major splits of mt lineages, the node ages and the 95% HPD intervals (in mya) are provided. The alternative ancestral subdivision/inheritance scenarios with likelihoods of 15% or more are also indicated in the tree. As a main result, the analyses suggest that O. dolium originated 6.92 to 4.13 mya (95% HPD interval) around the Miocene-Pliocene boundary. However, the analyses support that the broad diversification into numerous lineages happened during the Pleistocene. According to the reconstruction of the species’ geographic range history (ancestors allowed to occupy a maximum of two geographic areas/dispersal constraint for areas not immediately adjacent: ‘0.5’), the ancestral O. dolium was distributed in the Western Carpathians (Maximum Likelihood of ancestral stage at cladogenesis event: 0.83). The alternative analysis run without range restrictions (ancestors allowed to occupy all four geographic areas/dispersal constraint for areas not immediately adjacent: ‘0.5’) predicted similar ancestral ranges but the most likely ancestral range scenarios generally obtained lower likelihood values. Moreover, alternative ancestral ranges predicted for several nodes comprised differing geographic areas (Fig. S2). Analyses run with the same range constraints as in the two previous analyses, but with differing dispersal constraint (migrations between all areas are equally likely) resulted in identical ancestral ranges for almost all nodes, only the ML values differed slightly (by a maximum of 10 percent) (data not shown).
The map shows the distribution areas of O. dolium in the four Alpine and Carpathian mountain areas sampled (encoded by different colors). Small black dots represent localities sampled in the present study. The linearized molecular clock dated BI tree shows the relationships of selected mt lineages (COI/16S data) of O. dolium. Black and grey dots indicate nodes with high posterior probabilities (see figure for values). The colored symbols at the branch tips indicate the geographic origin of the haplotypes. The ancestors were allowed to occupy a maximum of two geographic areas. At the cladogenesis events (nodes), all alternative ancestral subdivision/inheritance scenarios with likelihoods of 15% or more are indicated, together with the respective likelihoods, and separated by an “or”. When scenarios for cladogenesis events involve two ancestral areas, the symbol for the likely ancestral area/−s is/are provided left to each of the two branches. For nodes representing major splits, node ages and 95% posterior HPD intervals are indicated. A time scale in mya is given below.
Middle and Late Pleistocene Distribution Derived from Fossil Record
The information about the fossil distribution of O. dolium is displayed in Fig. 7. Separate maps are presented for cold and warm Pleistocene climate stages because several localities provide records of both glacial and interglacial periods. The literature reports O. dolium from about 100 fossil sites, most of which are located in the Western Carpathians, including 40 sites in Slovakia alone [47], [48]. Four Czech localities provide Middle Pleistocene to Holocene fossils, among them the northern-most confirmed site, which is located about 30 km north of Prague [47]. Around 20 sites are located in the Western Carpathians of Northern Hungary, with Holocene, Weichselian, Eemian and Elsterian deposits [19], [48], [50]. Some of these represent the earliest fossil records of the species, particularly the sites in the Hungarian Bükk mountains and the surroundings of Budapest. Deposits from the Slovenská Skala in South-Eastern Slovakia also date back to the Elsterian and Holsteinian about 400 ka ago [48]. Fewer data are available from Austria and Germany. Most of the 18 sites of the Northern Calcareous Alps are dated to the Holocene or the Weichselian. Only a single Alpine site (Vienna, Austria) was tentatively classified to the end of the Middle Pleistocene by [42]; it is not included in the maps because an assignment to a cold or warm climate stage is not possible. The Nuβloch site (Baden-Württemberg, Germany) is the earliest fossil occurrence in the foreland of the Western Alps in the Early and Middle Weichselian, with an increasing frequency in the Late Weichselian layers [49]. Numerous fossil records of O. dolium were reported from loess deposits of the Pannonian Basin [19], [52], [55], most of which are located along the Danube in southern Hungary and Northern Serbia, where the rivers Sava and Danube join. Most Pannonian sites provide Weichselian deposits, although the species was already present close to Novi Sad in Serbia at the end of the Saalian [51]. In Eastern Croatia, fossils of O. dolium date back from the Weichselian to the Saalian (over 200 ka ago) [44], [45].
The LGM glacier line indicates the maximum extent of glaciers during the Weichselian (35 to 19 ka). The grey-shaded areas correspond to the current distribution of O. dolium.
Discussion
Origin and Diversification of O. dolium
Zimmermann [56] hypothesized that O. dolium originated in the Northern Calcareous Alps like several congeneric species. Frank [42] took up the same position and stated that O. dolium emerged in the Northern Calcareous Alps in the Early Pleistocene. However, these assumptions regarding the species’ origin are rather tentative because no Pliocene and Early Pleistocene fossils are known from Central Europe. The earliest record of O. dolium from South-Eastern Slovakia and North-Western Hungary dates back to the Elsterian (475–370 kya) and even earlier [19], [48], whereas almost all records from Alpine sites were assigned to the Weichselian (115–11 ka ago) or the Holocene; only a single record (Austria, Vienna) was vaguely dated to the Late Middle Pleistocene [42]. In general, the gastropod record before the Middle Pleistocene is extremely scarce because sediments predominantly consist of red clay, being inappropriate for the fossilization of gastropod shells [19]. Besides, environmental dynamics in the mountainous regions, especially in the preferred limestone habitats, offer few opportunities for shell fossilization. Hence, in the present study, hypotheses regarding the origin and diversification of O. dolium are mainly based on molecular genetic data.
Contradicting the assumption of Zimmermann [56] and Frank [42], the variability of the mt and nc markers and the geographic distribution of haplotypes support an origin of O. dolium in the Western Carpathians. Seven out of ten clades occur in the Western Carpathians, including two highly diverged clades (9 and 10) which split from the basal nodes of the trees. The populations of the Alps and the Western Carpathians are not reciprocally monophyletic - lineages of both areas derive from three nodes in the mt trees each. The number of H4/H3 alleles highly diverged from the three main alleles HT1, HT2 and HT3 is larger in the Western Carpathians as well. The geographic range reconstruction supports a scenario in which the most recent common ancestor was distributed in the Western Carpathians (ML: 0.83) around 6.26 mya (95% HPD: 5.15 to 7.42); the Northern Calcareous Alps were probably settled later (Fig. 6). A scenario in which the MRCA was distributed in both the Western Carpathians and the Northern Calcareous Alps obtained only low support (ML: 0.12). The distributional patterns of mt clades/variants can best be explained by multiple (probably two or three) migrations between Alps and Western Carpathians with single specimens or populations carrying unique or similar mt variants, respectively. Since the Alpine clades are embedded within the diversity of the Western Carpathians, a predominant east-west migration route is most probable. Alternative to scenarios with multiple migrations, the Alpine diversity could have resulted from a single migration involving multiple individuals carrying strongly differing mt variants. However, the geographically distinct distributions of the Alpine mt clades suggest that the lineages evolved independently from each other. Moreover, the molecular clock analysis indicated that the Alpine mt clades separated from their closest related Western Carpathian lineages during different time periods. The results of the analysis suggested that the Alpine mt sub-clade 1A descended from the Western Carpathians rather recently during the Middle Pleistocene, whereas clade 8 and the cluster including clades 4 and 6 separated from their closest related Western Carpathians lineages probably during the Early Pleistocene already. The nc H4/H3 sequence patterns also support at least two independent migration events, as two highly diverged, geographically more or less separated clusters (HT1, HT2 and similar variants) were found.
Pleistocene Refuges and Postglacial Expansion Routes
One of the present study’s main objectives is the detection of potential glacial refuges of O. dolium. The four major limestone areas currently inhabited by O. dolium (Western Alps, Northern Calcareous Alps, Southern Calcareous Alps and Western Carpathians) are treated separately, as is the Pannonian Basin, in which the species is apparently not found nowadays.
The extensive fossil record of the Western Carpathians, with data from both glacials and interglacials, confirms that the area was permanently settled, at least during the second half of the Middle Pleistocene and the Late Pleistocene [48]. In particular the extensive record of Weichselian (115–11 ka ago) fossils provides evidence that O. dolium was widely distributed during the last glacial (Fig. 7). Moreover, despite the comparably small sample size, the Western Carpathian populations show a large genetic diversity with complex distribution patterns. There is no clear geographic structure regarding the distribution of mt clades and nc alleles, and the data do not indicate a serious loss of genetic diversity due to genetic bottlenecks. Unlike in the Alpine region, the loss of potential habitat presumably was less significant in the Western Carpathians, which were not affected by massive glaciations during Pleistocene cold stages. The scattered distribution of limestone bedrock in the Western Carpathians is another factor, which may have triggered diversification and preservation of numerous genetic lineages.
The eastern part of the Northern Calcareous Alps potentially provided the largest Alpine refuge for calciphilous taxa because it continuously offered non-glaciated limestone areas. Patterns of endemism and comparative phylogeographic analyses in Alpine plants [8] provide additional evidence for refuges in this area. Similarly, the Northern Calcareous Alps harbor a number of endemic species such as Trochulus oreinos and Cylindrus obtusus which probably originated in that region [14]. Haplotypes of all four Alpine mt clades are found here, with populations located somewhat separated from each other, and therefore suggesting several smaller refuges. Moreover, the respective populations show a high diversity in nc H4/H3 alleles. The most common mt clade (1A) is distributed from Lower Austria (Gutenstein Alps) in the east to Tirol (Lechtal Alps) in the west, spanning a distance of 400 km (Fig. 3). A distinct mt clade (4B) is present at the eastern edge of the Alps (Wienerwald), a region which is geologically somewhat isolated from the Northern Calcareous Alps due to the predominance of siliciclastics in the intermediate region. The large nc diversity with several H4/H3 alleles, each strongly diverged from the most common alleles HT1 and HT2, and the presence of Weichselian fossil deposits support the assumption that the Wienerwald served as a long-term refuge. The Gutenstein Alps at the northeastern margin of the Northern Calcareous Alps may have provided a refuge area for the population exhibiting mt clade 6. Another potential refuge was located close to the former LGM glacier line in Salzburg and Upper Austria - mt clade 8 is distributed exclusively in this region, and fossil deposits indicate the local presence of the species during the Weichselian [42]. Some of the specimens carrying haplotypes of mt clade 8 possess the nc allele HT2 or similar ones, which elsewhere occur in the Wienerwald only. This can either be explained by past gene flow between the two currently separated populations or by ancestral polymorphism, i.e., the persistence of ancestral histone variants in both areas.
The geographic range reconstruction suggests that the Western Alps were settled from the Northern Calcareous Alps during the Middle Pleistocene 1.05 to 0.66 mya (95% HPD interval). Samples from the Western Alpine sites all form a single mt sub-clade (4A) and have the nc haplotype HT1 or similar types. Although material is available from three sample sites only, the Western Alpine populations show higher distances in the mt sequences (max. uncorrected p-distance COI: 3.8%) than any other Alpine clade. The unique presence of the highly diverged mt clade 4A provides support for a refuge in the Western Alps. Since the Western Alps were almost completely covered by ice during several glacials, populations may have outlasted the glacial periods in several smaller refuges at the Western Alpine margins of Switzerland and France as was proposed for Trochulus villosus [15] or Carychium tridentatum [13]. Fossils in Early to Late Weichselian deposits of the Western Alpine foreland (Nuβloch, Baden-Württemberg, Germany) clearly support that assumption, at least for the last glacial period [49]. We have no molecular data from the German and French areas, but a common ancestry of the western populations is supported by similar conchological traits (collection material of the Natural History Museum, Vienna and Naturmuseum Senckenberg, Frankfurt am Main; Harl et al. in prep.).
The Southern Calcareous Alps were almost completely covered by glaciers during the LGM and hence provided only a small potential refuge for calciphilous taxa in the eastern-most part: a small non-glaciated area in the Karawanks [57]. The phylogenetic data argue against a permanent refuge of O. dolium in the Southern Calcareous Alps and suggest that clades 1A and 6 probably reached this region rather recently. In the Southern Calcareous Alps, sub-clade 1A is found at a single site (GAI1) only, with haplotypes similar to those in populations of the Tirolian and Bavarian Northern Calcareous Alps. Since both areas were fully glaciated, clade 1A might have crossed the Central Alps in this area after the LGM. Apart from this single locality, all other sites of the Southern Calcareous Alps possessed haplotypes of clade 6, which is distributed in the Northern Calcareous Alps as well. The similarity of clade 6 haplotypes rather indicates a very recent expansion during the Late Pleistocene or Holocene. However, the presence of a distinct variant of clade 6 at a single site in the Karawanks (KWN1) could be an indication for a Southern Calcareous Alpine refuge.
The fossil record indicates a more or less continuous presence of O. dolium in the Pannonian Basin at least from the Saalian onwards (over 200 ka ago). The habitats of the Pannonian populations were probably patchily distributed forests near rivers, which were present in the area even during the LGM [7]. The occurrence of trees at the respective sites is additionally supported by the co-occurrence of other woodland species in the same loess strata, such as Semilimax semilimax, Ena montana or Aegopinella ressmanni [45], [58]. One might ask whether riparian drift from Alpine or Western Carpathian regions could account for the presence of O. dolium in the fossil record of the southern part of the Pannonian Basin. However, the high abundance of fossil O. dolium in the Pannonian Basin indicates a local source. The contemporary absence of the species is probably the result of anthropogenically induced loss of suitable habitat. Ložek [48] stated that deforestation and dehydration are probably the reasons why these areas lack several gastropod species which were still widely distributed during the Eemian. Thus, the expansion of agricultural areas is a reasonable explanation for the decline of O. dolium populations in the Pannonian Basin during the Holocene.
Genetic Differentiation and Taxonomic Considerations
The intraspecific distances measured for the mt genes are among the highest found in pulmonate species (uncorrected p-distances: COI, 18.4%; 16S, 14.4%). By comparison, the genetically diverse helicid taxa Theba pisana and Arianta arbustorum show COI distances of 13.6% [59] and 15% [10], respectively. Regarding the non-protein coding 16S, higher intraspecific distances were found in the clausiliid species Charpentieria itala with 24.2% [60]. 16S divergence is also high in the bradybaenid Euhadra quaesita with 14.1% [61]. Regarding the nc sequences analyzed, the largest p-distance measured within the protein coding H4 and H3 sequences is 0.8%. The highest distance observed in the non-protein coding spacer region is 1.8%. [59] reported p-distances of 0.5% in the non-coding ITS1 sequences of Theba pisana, whereas the intraspecific sequence divergence within Arion subfuscus, a species extremely variable in its mtDNA, is only 0.3% for the ITS1 sequence [62].
Considering the large genetic variability found within the populations of O. dolium, the question arises whether some of the lineages might even represent distinct species. However, there are no indications of hybridization barriers that would suggest splitting the groups into different species. The data indicate gene flow between clades, as suggested by the fact that the specimens displaying the main histone gene variants each feature haplotypes of very distant mt clades (see Fig. 5). However, whether genetic groups correspond to currently accepted subspecies remains a problematic issue. More than 20 subspecies have been described for O. dolium, equally divided between the Alps and the Western Carpathians [20]. Most were characterized by minor differences in shell shape and the formation of the aperture folds. Our study includes specimens from several type localities of Slovakian subspecies, namely of O. d. titan Brancsik, 1887, O. d. brancsikii Clessin, 1887, O. d. minima Brancsik, 1887, and O. d. cebratica Westerlund, 1887 (Table 1 and Fig. S1). However, none of the clades can be definitively attributed to one of these subspecies. For instance, populations of the very slender, large-shelled O. d. brancsikii share the same mt haplotypes (clade 1C), and the nc haplotype HT3, with compact, small O. d. minima morphs. The very large-shelled O. d. titan possesses a diverged mt haplotype within clade 9 but shows the most common nc haplotype HT1, which otherwise is found in the Alps only. The Slovakian O. d. cebratica (clade 5) clusters with the Alpine mt clades 4 and 6 but displays a distinct nc haplotype (HT3). Despite the extremely high conchological and genetic variability in the Western Carpathians, we could not detect a clear correspondence between genetic haplotype groups and subspecies ranges defined by [17]. The geographically isolated and conchologically aberrant endemic of the Wienerwald, O. dolium infima, is genetically differentiated from other populations in the mt trees (sub-clade 4B). In contrast, individuals from the type localities of O. dolium gracilior Zimmermann, 1932, O. dolium pseudogularis Wagner, 1912, O. dolium edita Ehrmann, 1933, and O. dolium raxae Gittenberger, 1978, all possess haplotypes of the homogeneous mt clade 1 and exhibit the nc haplotypes HT1 or its derivatives, none of them forming a distinct sub-clade. The latter two taxa, O. d. edita and O. d. raxae, initially were of special interest for our study because they were reported to occur only at high altitudes [17], [56]. Nonetheless, specimens from the corresponding localities do not differ from the common types of the surrounding lowlands in the markers analyzed. However, for final taxonomic decisions, the most important question is whether the presumed morphological distinctness of the various subspecies can be confirmed by morphometric investigations. Such analyses are under way to quantify morphological differences and to focus on the effect of altitude on shape formation and shell size (Harl et al., in prep.).
Supporting Information
Figure S1.
Pictures of selected types of O. dolium subspecies. Specimens collected at the type localities of the respective subspecies were investigated in the present study. The types shown shall rather be considered as examples for the species’ variability than defined discrete entities of morphologically separated populations. Some of the types indeed represent extreme morphs but transitional forms are found in most populations. The pictures were already published by Harl et al. (2011) together with data on all other currently known subspecies. In the following we provide the collection data of the specimens shown: A: O. d. dolium (syntype NHMW 14765/1820.26.61/2), B: O. d. edita (syntype LML ALT/5319/1), C: O. d. raxae (syntype LML ALT/5354/1), D: O. d. pseudogularis (syntype NHMW 56158), E: O. d. gracilior (syntype LML ALT/5343), F: O. d. infima (syntype LML/ALT5353/1), G: O. d. brancsikii (syntype ? NHMW J. N. 22075), H: O. d. titan (syntype NHMW 68377 (5448)/3), I: O. d. cebratica (syntype MNHG Wstld2090), J: O. d. minima (syntype 27044/2). Abbreviations for Museums: NHMW (Naturhistorisches Museum, Wien), LML (Oberösterreichisches Landesmuseum, Linz) and MNHG (Naturhistoriska Museum, Göteborg). The scale bar indicates 5 mm.
https://doi.org/10.1371/journal.pone.0096012.s001
(TIF)
Figure S2.
Reconstruction of the geographic range evolution. The map shows the distribution areas of O. dolium in the four Alpine and Carpathian mountain areas sampled (encoded by different colors). Small black dots represent localities sampled in the present study. The linearized molecular clock dated BI tree shows the relationships of selected mt lineages (COI/16S data) of O. dolium. Black and grey dots indicate nodes with high posterior probabilities (see figure for values). The colored symbols at the branch tips indicate the geographic origin of the haplotypes. The ancestors were allowed to occupy all four geographic areas. At the cladogenesis events (nodes), all alternative ancestral subdivision/inheritance scenarios with likelihoods of 15% or more are indicated, together with the respective likelihoods, and separated by an “or”. When scenarios for cladogenesis events involve two or more ancestral areas, the symbol for the likely ancestral area/−s is/are provided left to each of the two branches. For nodes representing major splits, node ages and 95% posterior HPD intervals are indicated. A time scale in mya is given below.
https://doi.org/10.1371/journal.pone.0096012.s002
(TIF)
Acknowledgments
We thank Simone Latkolik and Barbara Däubl for technical assistance in the lab. Alexander Reischuetz and Ira Richling provided essential samples from the Western Carpathians and the Western Alps, respectively. Benjamín Gómez-Moliner provided important DNA samples of Spanish Orculella used for the molecular clock analysis. Alexander Reischuetz, Peter Reischuetz and Barna Páll-Gergely contributed valuable information and suggestions in conversations, and we are especially grateful to Wilhelm Pinsker for discussions and critical comments on the manuscript.
Author Contributions
Conceived and designed the experiments: JH MD LK HS EH. Performed the experiments: JH LK. Analyzed the data: JH LK EH. Contributed reagents/materials/analysis tools: JH MD LK HS EH. Wrote the paper: JH EH. Participated in the design and coordination of the study: MD LK HS EH. Carried out the molecular genetic studies, sequence analyses, phylogenetic analysis: JH. Created the graphics: JH. Participated in the sequence analyses: LK. Contributed in drafting the manuscript: EH. Drafted the manuscript: JH. Participated in the field work: JH MD LK HS EH. Approved the final manuscript: JH MD LK HS EH.
References
- 1. Hewitt G (2004) Genetic consequences of climatic oscillations in the Quaternary. Philos Trans R Soc London Ser B Biol Sci 359: 183–195.
- 2. Muttoni G, Carcano C, Garzanti E, Ghielmi M, Piccin A, et al. (2003) Onset of major Pleistocene glaciations in the Alps. Geology 31: 989–992.
- 3. Ivy-Ochs S, Kerschner H, Reuther A, Preusser F, Heine K, et al. (2008) Chronology of the last glacial cycle in the European Alps. J Quat Sci 23: 559–573.
- 4. Holdhaus K (1954) Die Spuren der Eiszeit in der Tierwelt Europas. Abhandlungen der Zool Gesellschaft Wien 18: 1–493.
- 5. Benke M, Brändle M, Albrecht C, Wilke T (2009) Pleistocene phylogeography and phylogenetic concordance in cold-adapted spring snails (Bythinella spp.). Mol Ecol 18: 890–903.
- 6. Sommer RS, Zachos FE (2009) Fossil evidence and phylogeography of temperate species: “glacial refugia” and post-glacial recolonization. J Biogeogr 36: 2013–2020.
- 7. Willis KJ, van Andel TH (2004) Trees or no trees? The environments of central and eastern Europe during the Last Glaciation. Quat Sci Rev 23: 2369–2387.
- 8. Tribsch A, Schönswetter P (2003) Patterns of endemism and comparative phylogeography confirm palaeoenvironmental evidence for Pleistocene refugia in the Eastern Alps. Taxon 52: 477–497.
- 9. Haase M, Bisenberger A (2003) Allozymic differentiation in the land snail Arianta arbustorum (Stylommatophora, Helicidae): historical inferences. J Zool Syst Evol Res 41: 175–185.
- 10. Gittenberger E, Piel WH, Groenenberg DSJ (2004) The Pleistocene glaciations and the evolutionary history of the polytypic snail species Arianta arbustorum (Gastropoda, Pulmonata, Helicidae). Mol Phylogenet Evol 30: 64–73.
- 11. Haase M, Misof B (2009) Dynamic gastropods: stable shell polymorphism despite gene flow in the land snail Arianta arbustorum. J Zool Syst Evol Res 47: 105–114.
- 12. Haase M, Esch S, Misof B (2013) Local adaptation, refugial isolation and secondary contact of alpine populations of the land snail Arianta arbustorum. J Molluscan Stud 79: 241–248.
- 13. Weigand AM, Pfenninger M, Jochum A, Klussmann-Kolb A (2012) Alpine crossroads or origin of genetic diversity? Comparative phylogeography of two sympatric microgastropod species. PLoS One 7: e37089.
- 14. Duda M, Sattmann H, Haring E, Bartel D, Winkler H, et al. (2011) Genetic differentiation and shell morphology of Trochulus oreinos (Wagner, 1915) and T. hispidus (Linnaeus, 1758) (Pulmonata: Hygromiidae) in the northeastern Alps. J Molluscan Stud 77: 30–40.
- 15. Dépraz A, Hausser J, Pfenninger M, Cordellier M (2008) Postglacial recolonization at a snail’s pace (Trochulus villosus): confronting competing refugia hypotheses using model selection. Mol Ecol 17: 2449–2462.
- 16. Provan J, Bennett KD (2008) Phylogeographic insights into cryptic glacial refugia. Trends Ecol Evol 23: 564–571.
- 17. Klemm W (1974) Die Verbreitung der rezenten Land-Gehäuse-Schnecken in Österreich. Denkschriften der Österrerreichischen Akad der Wissenschaften (mathematisch-naturwissenschaftliche Klasse) 117: 1–503.
- 18.
Lisický MJ (1991) Mollusca Slovenska (The Slovak molluscs). Bratislava: VEDA vydavateľstvo Slovenskej akadémie vied.
- 19. Fűköh L, Krolopp E, Sümegi P (1995) Quaternary malacostratigraphy in Hungary. Malacol Newsl 1: 1–219.
- 20. Harl J, Sattmann H, Schileyko A (2011) Types of the extant taxa of the landsnail genus Orcula Held 1837 (Gastropoda: Stylommatophora: Orculidae). Arch für Molluskenkd 140: 175–199.
- 21. Kruckenhauser L, Harl J, Sattmann H (2011) Optimized drowning procedures for pulmonate land snails allowing subsequent DNA analysis and anatomical dissections. Ann Naturhist Mus Wien B 112: 173–175.
- 22.
Arrébola J, Razkin O, Gómez-Moliner B, Páll-Gergely B (2012) Redescription of Orculella aragonica (Westerlund 1897), an Iberian species different from O. bulgarica (Hesse 1915) (Gastropoda: Pulmonata: Orculidae). North West J Zool 8.
- 23.
Marazzi S (2005) Atlante Orografico delle Alpi. SOIUSA.
- 24. Hovorka D, Méres Š, Ivan P (1993) Pre-Alpine Western Carpathians basement complexes: lithology and geodynamic setting. Mittleilungen der Österreichischen Geol Gesellschaft 86: 33–44.
- 25. Armbruster GFJ, Bohme M, Bernhard D, Schlegel M (2005) The H3/H4 histone gene cluster of land snails (Gastropoda: Stylommatophora): ts/tv ratio, GC3 drive and signals in Stylommatophoran phylogeny. J Mollus Stud 71: 339–348.
- 26. Folmer O, Black M, Heah W, Lutz R, Vrijenhoek R (1994) DNA primers for amplification of mitochondrial cytochrome C oxidase subunit I from diverse metazoan invertebrates. Mol Mar Biol Biotechnol 3: 294–299.
- 27.
Cadahía L, Harl J, Duda M, Sattmann H, Kruckenhauser L, et al.. (2014) New data on the phylogeny of Ariantinae (Pulmonata, Helicidae) and the systematic position of Cylindrus obtusus based on nuclear and mitochondrial DNA marker sequences. J Zool Syst Evol Res: in press.
- 28. Hall TA (1999) BioEdit: a user-friendly biological sequences alignment editor and analysis program for Windows 95/98/NT. Nucleic Acids Symp Ser 41: 95–98.
- 29. Librado P, Rozas J (2009) DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics 25: 1451–1452.
- 30. Tamura K, Peterson D, Peterson N, Stecher G, Nei M, et al. (2011) MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Mol Biol Evol 10: 2731–2739.
- 31. Posada D (2008) jModelTest: phylogenetic model averaging. Mol Biol Evol 25: 1253–1256.
- 32. Huelsenbeck JP, Ronquist F (2001) MRBAYES: Bayesian inference for phylogeny. Bioinformatics 17: 754–755.
- 33. Ronquist F, Huelsenbeck JP (2003) MrBayes 3: Bayesian phylogenetic inference under mixed models. Bioinformatics 19: 1572–1574.
- 34.
Rambaut A, Drummond AJ (2009) Tracer v1.5. Available: from BEAST Software website: http://tree.bio.ed.ac.uk/software/tracer/.
- 35. Larkin MA, Blackshields G (2007) ClustalW and ClustalX version 2. Bioinformatics 23: 2947–2948.
- 36. Capella-Gutiérrez S, Silla-Martínez JM, Gabaldón T (2009) trimAl: a tool for automated alignment trimming in large-scale phylogenetic analyses. Bioinformatics 25: 1972–1973.
- 37. Sánchez R, Serra F, Tárraga J, Medina I, Carbonell J, et al. (2011) Phylemon 2.0: a suite of web-tools for molecular evolution, phylogenetics, phylogenomics and hypotheses testing. Nucleic Acids Res 39: W470–W474.
- 38. Ree RH, Smith SA (2008) Maximum likelihood inference of geographic range evolution by dispersal, local extinction, and cladogenesis. Syst Biol 57: 4–14.
- 39. Drummond AJ, Rambaut A (2007) BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evol Biol 7: 214.
- 40. Robles F, Martínez-Ortí A (2009) Moluscos continentales de los alrededores de Molina de Aragón (Guadalajara, España), con notas sobre Orculella bulgarica (Hesse, 1915) (Gastropoda, Orculidae). Iberus 27: 99–105.
- 41. Agustí J, Santos-Cubedo A, Furió M, De Marfá R, Blain H-A, et al. (2011) The late Neogene-early Quaternary small vertebrate succession from the Almenara-Casablanca karst complex (Castellón, Eastern Spain): Chronologic and paleoclimatic context. Quat Int 243: 183–191.
- 42.
Frank C (2006) Plio-pleistozäne und holozäne Mollusken Österreichs. Vienna.
- 43. Frank C, Terhorst B, Damm B, Thiel C, Frechen M, et al. (2011) Pleistocene loess deposits and mollusc assemblages in the Eastern Pre-Alps. Quat Sci J 60: 126–136.
- 44. Hupuczi J, Molnár D, Galović L, Sümegi P (2010) Preliminary malacological investigation of the loess profile at Šarengrad, Croatia. Cent Eur J Geosci 2: 57–63.
- 45. Molnár D, Hupuczi J, Galović L, Sümegi P (2010) Preliminary malacological investigation on the loess profile at Zmajevac, Croatia. Cent Eur J Geosci 2: 52–56.
- 46. Ložek V (2001) Molluscan fauna from the loess series of Bohemia and Moravia. Quat Int 76–77: 141–156.
- 47. Ložek V (2006) Last Glacial paleoenvironments of the West Carpathians in the light of fossil malacofauna. J Geol Sci 26: 73–84.
- 48. Ložek V (1964) Quartärmollusken der Tschechoslowakei. Rozpr Ústředního ústavu Geol 31: 1–374.
- 49. Moine O, Rousseau D-D, Antoine P (2005) Terrestrial molluscan records of Weichselian Lower to Middle Pleniglacial climatic changes from the Nussloch loess series (Rhine Valley, Germany): the impact of local factors. Boreas 34: 363–380.
- 50. Sümegi P, Gulyás S, Persaits G, Gergely Páll D, Molnár D (2011) The loess-paleosol sequence of Basaharc (Hungary) revisited: Mollusc-based paleoecological results for the Middle and Upper Pleistocene. Quat Int 240: 181–192.
- 51. Marković SB, Oches EA, Gaudenyi T, Jovanovic M, Hambach U, et al. (2004) Palaeoclimate record in the Late Pleistocene loess-paleosol sequence at Miseluk (Vojvodina, Serbia). Quaternaire 15: 361–368.
- 52. Mitrović B (2007) Pleistocene malacofauna of the Požarevac Danube area (NE Serbia). Geološki Anal Balk poluostrva 68: 81–89.
- 53. Marković SB, McCoy WB, Oches EA, Savić S, Gaudenyi T, et al. (2005) Paleoclimate record in the Upper Pleistocene loess-paleosol sequence at Petrovaradin brickyard (Vojvodina, Serbia). Geol Carpathica 56: 545–552.
- 54. Zagwijn WH (1992) The beginning of the ice age in Europe and its major subdivisions. Quat Sci Rev 11: 583–591.
- 55. Marković SB, Bokhorst MP, Vandenberghe J, McCoy WD, Oches EA, et al. (2008) Late Pleistocene loess-palaeosol sequences in the Vojvodina region, north Serbia. J Quat Sci 23: 73–84.
- 56. Zimmermann S (1932) Über die Verbreitung und die Formen des Genus Orcula Held in den Ostalpen. Arch für Naturgeschichte 1: 1–57.
- 57. Van Husen D (1997) LGM and late-glacial fluctuations in the Eastern Alps. Quat Int 38–39: 109–118.
- 58. Sümegi P, Krolopp E (2002) Quatermalacological analyses for modeling of the Upper Weichselian palaeoenvironmental changes in the Carpathian Basin. Quat Int 91: 53–63.
- 59. Greve C, Hutterer R, Groh K, Haase M, Misof B (2010) Evolutionary diversification of the genus Theba (Gastropoda: Helicidae) in space and time: A land snail conquering islands and continents. Mol Phylogenet Evol 57: 572–584.
- 60. Scheel BM, Hausdorf B (2012) Survival and differentiation of subspecies of the land snail Charpentieria itala in mountain refuges in the Southern Alps. Mol Ecol 21: 3794–3808.
- 61. Watanabe Y, Chiba S (2001) High within-population mitochondrial DNA variation due to microvicariance and population mixing in the land snail Euhadra quaesita (Pulmonata: Bradybaenidae). Mol Ecol 10: 2635–2645.
- 62. Pinceel J, Jordaens K, Pfenninger M, Backeljau T (2005) Rangewide phylogeography of a terrestrial slug in Europe: evidence for Alpine refugia and rapid colonization after the Pleistocene glaciations. Mol Ecol 14: 1133–1150.