Browse Subject Areas

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Multiple Glacial Refugia of the Low-Dispersal Ground Beetle Carabus irregularis: Molecular Data Support Predictions of Species Distribution Models

  • Katharina Homburg ,

    Affiliation Leuphana University Lüneburg, Institute of Ecology, Lüneburg, Germany

  • Claudia Drees,

    Affiliations Leuphana University Lüneburg, Institute of Ecology, Lüneburg, Germany, University of Hamburg, Biocentre Grindel and Zoological Museum, Behavioural Biology, Hamburg, Germany

  • Martin M. Gossner,

    Affiliation Technische Universität München, Department of Ecology and Ecosystem Management, Freising-Weihenstephan, Germany

  • László Rakosy,

    Affiliation Babeş-Bolyai University Cluj-Napoca, Faculty of Biology and Geology, Cluj-Napoca, Romania

  • Al Vrezec,

    Affiliation National Institute of Biology, Department of Freshwater and Terrestrial Ecosystems Research, Ljubljana, Slovenia

  • Thorsten Assmann

    Affiliation Leuphana University Lüneburg, Institute of Ecology, Lüneburg, Germany


  • Katharina Homburg, 
  • Claudia Drees, 
  • Martin M. Gossner, 
  • László Rakosy, 
  • Al Vrezec, 
  • Thorsten Assmann


Classical glacial refugia such as the southern European peninsulas were important for species survival during glacial periods and acted as sources of post-glacial colonisation processes. Only recently, some studies have provided evidence for glacial refugia north of the southern European peninsulas. In the present study, we combined species distribution models (SDMs) with phylogeographic analyses (using mitochondrial DNA = mtDNA) to investigate if the cold-adapted, stenotopic and flightless ground beetle species, Carabus irregularis, survived the Last Glacial Maximum (LGM) in classical and/or other refugia. SDMs (for both a western European and for a Carpathian subgroup) were calculated with MAXENT on the basis of 645 species records to predict current and past distribution patterns. Two mtDNA loci (CO1 and ND5, concatenated sequence length: 1785 bp) were analyzed from 91 C. irregularis specimens to reconstruct the phylogeography of Central and eastern European populations and to estimate divergence times of the given lineages. Strong intra-specific genetic differentiation (inter-clade ΦST values ranged from 0.92 to 0.99) implied long-term isolation of major clades and subsclades. The high divergence between the nominate subspecies and the Carpathian subspecies C. i. montandoni points to two independent species rather than subspecies (K-2P distance 0.042 ± 0.004; supposed divergence of the maternal lineages dated back 1.6 to 2.5 million years BP) differing not only morphologically but also genetically and ecologically from each other. The SDMs also inferred classical as well as other refugia for C. irregularis, especially north of the Alps, in southeastern Europe and in the Carpathians. The coincidences between the results of both methods confirm the assumption of multiple glacial refugia for the studied species and the usefulness of combining methodological approaches for the understanding of the history of low-dispersal insect species.


The glaciers and permafrost grounds of the Quaternary ice ages shaped distribution patterns in biodiversity, as they caused species extinctions as well as the retraction of species' distribution ranges. Glacial refugia were important for species survival in glacial and interglacial periods and hosted many species which had recently been widespread across Europe [1]. Temperate and warm-adapted species in particular are assumed to have colonised central and northern parts of Europe after the Last Glacial Maximum (LGM), originating from classical refugia such as the southern European peninsulas (Iberian, Apennine and Balkan) [1], [2]. Cold-adapted species that had recently inhabited alpine and arctic habitats are assumed to have survived in the margins of southern European mountain chains [3]. In recent decades, palaeontological, palynological and phylogeographic studies have provided evidence for glacial refugia far north of the southern European peninsulas and the Alps (so-called cryptic refugia) [4], [5]. In contrast to the species which were widespread across Europe [6], [7], species with restricted and disjunct ranges in the European low and high mountain chains are poorly studied with regard to spatial patterns in genetic variation. However, some mountain plant species inhabiting the Alps and geographically relatively close mountain systems (e.g. the Central European mountain ranges, the Carpathians) show high genetic differentiation [8], [9]. Zoological studies have covered arthropod species with high dispersal abilities (e.g. due to flight activity or ballooning): Butterflies [10] and spiders [11] show low level differentiations and indicate young and even post-glacial splits between populations from the Alps and adjoining high and low mountain ranges [12]. In contrast, low-dispersal species are likely to show stronger genetic differentiation patterns due to the absence of gene flow (e.g. flightless carabids; [13], [14]); however, the population histories of these species have been studied only insufficiently to date.

New study methods such as species distribution models (SDMs) have become useful tools to predict potential habitat distribution of species under current climate conditions and – when analysed in conjunction with palaeoclimatic data – to project species' past distributions [15]. Many studies detected classical as well as cryptic refugia using past predictive models (e.g. [16], [17]). Unaffected by incomplete sampling, the extinction of genetic variants and large-scale range shifts of species, past predictive models can contribute to the localization and sizing of species' glacial refugia [18]. Consequently, past predictive models and phylogeographic analyses can be used complementarily to develop reasonable scenarios of species' range retractions during glaciations and of range expansions during interglacial stages.

In this study we applied a combined approach using modelling and of phylogeographic analysis of two mitochondrial DNA loci (CO1 and ND5) to investigate the history of the flightless ground beetle species, which occurs in a disjunct range across high- and low-altitude mountains of Central and eastern Europe: Carabus irregularis. We aimed to test (i) whether the cold-adapted species survived the LGM in classical and/or cryptic refugia and (ii) whether the low-dispersal species shows old phylogenetic groups and high intra-specific variation. Since the CO1 locus has been successfully used in barcoding [19], we also checked (iii) whether the phylogeny of Carabus irregularis is consistent with the intra-specific taxonomy of the species.

We addressed the following specific study questions: (1) Where did C. irregularis find suitable climate conditions (refugia) to survive the LGM? (2) Does the species' genetic differentiation pattern provide molecular evidence for low dispersal and for the detected potential glacial refugia? (3) Is the phylogeny of C. irregularis congruent with its taxonomic subspecies delineation?

Materials and Methods

Study species

Carabus (Platycarabus) irregularis Fabricius, 1792 occurs in cold and mesic habitats in montane to subalpine altitudes (e.g. beech forests and high mountain meadows). Suitable habitats are predominantly on limestone and are rich in snails, the preferred food of C. irregularis. The disjunct distribution of C. irregularis ranges from the Teutoburg Forest in the north to the Dinaric Mountains in the south and from the Carpathians in the east to the Jura Mountains in the west (Fig. 1) [20]. In the taxonomic literature, three subspecies are well-accepted even though different names are sometimes used [21]. In this study we follow Casale and Kryzhanovskij [22] and use the following binomials: C. i. irregularis s.str. in Central Europe including the northern Alps, C. i. bucephalus Kraatz, 1879 in the southeastern Alps and the northern Balkan Peninsula, and C. i. montandoni Buysson, 1882 in the Carpathians (Fig. 1). While the distribution ranges of C. i. irregularis and C. i. bucephalus pass into each other in the southeastern Alps, C. montandoni is geographically separated from the other two subspecies (see subspecies affiliation in Table S1). The species is declining in some regions [23], [24] and already extinct in some parts of its former distribution range (e.g. in Belgium [25]).

Figure 1. Species records of Carabus irregularis used for species distribution modelling and populations sampled for genetic analyses.

SDMs included species records displayed as small crosses. Genetic analyses included populations of the three subspecies – C. i. irregularis (circles), C. i. bucephalus (squares) and C. i. montandoni (triangles). Genetic analyses resulted in two major clades: a Carpathian clade (Clade A = blue colours) and a European clade (Clade B = green colours). The clades split into five subclades: an eastern Carpathian (Ae = dark blue), a southern Carpathian (As = light blue), a northern European (Bn = dark green), a southern European (Bs = olive green) and a purely Slovenian clade (Bsl = light green).

Species distribution modelling

In order to predict the potential habitat distribution of the study species under current and past climate conditions, we compiled 1,005 data points of observed occurrences of C. irregularis from the literature and databases (e.g. [26][33]), from museum and private entomological collections as well as our own field surveys (see complete list of records in Table S1). Since spatial clumping of species records due to different sampling intensities may bias SDMs, redundant localities were excluded by the model software, and the final locality input included 645 records. In addition to the geographic separation, there are also differences in the altitudinal distribution pattern of the two western European subspecies (C. i. irregularis and C. i. bucephalus) and the Carpathian subspecies (C. i. montandoni). According to the compiled species records, C. i. irregularis and C. i. bucephalus occur in high as well as low altitudes in the western part of its distribution range and on average in lower altitudes (altitudes as low as 250 m a.s.l.; mean altitude of western records = 774 m a.s.l.; standard error (SE) = 13.71; confidence interval (CI) = 26.87; see Fig. S1). From the eastern part of the range (the Carpathian mountain system), however, the species was scarcely recorded from altitudes below 600 m a.s.l. (mean altitude of Carpathian records = 912 m a.s.l.; SE = 45.50; CI = 89.19; see Fig. S1). Due to geographically separated distribution ranges and non-overlapping confidence intervals of the altitudinal distribution of the subspecies, we could not exclude different habitat preferences of the Carpathian subspecies. Thus, species records were split into two sets before modelling: one included occurrences of the Carpathian subspecies, the other comprised occurrences of the two other subspecies from the western part of the distribution range (Central Europe and Balkan). Climate data (monthly temperature and rainfall values) was generated according to Hijmans et al. [34] to 19 bioclimatic variables (BIO1-BIO19; Table S2) which representing factors of climate which are particularly significant in determining species' distributions [18]. For carabids, temperature and humidity variables have been demonstrated to be the most important environmental factors influencing habitat selection [35].

Current climate data (i.e. period 1950–2000) was downloaded with a spatial resolution of 30 arc s (∼1×1 km) from the WorldClim database ( Past (Last Glacial Maximum = LGM; 21,000 BP) climate data downscaled from the MIROC general circulation models to the required spatial resolution (2.5 arc min) was provided by the Paleoclimate Modelling Intercomparison Project Phase II ( Since our SDMs were to include only predictors with high explanatory power and ecological relevance [36], [37], we included only the variables which are most likely to influence the occurrence of C. irregularis. Therefore, we extracted all 19 bioclimatic variables at coordinates of species records and performed a principal component analysis (PCA). To reduce multicollinearity, we compiled a correlation matrix and excluded all strongly correlated variables (Pearson's r2>0.75; [38]). Final models included five uncorrelated predictor variables: ‘Isothermality’ (BIO3), ‘temperature seasonality’ (BIO4), ‘minimum temperature of coldest month’ (BIO6), ‘precipitation of wettest month’ (BIO13) and ‘precipitation of driest month’ (BIO14) – climatic measures reflecting habitat conditions which are likely to affect the occurrence of C. irregularis [30], [31]. The bioclimatic layers were cut to cover the recently documented distribution range of the species sufficiently but not too generous, since large layer extents could have a negative effect on SDM results [39]. The layers used ranged from 30 to 60 decimal degrees North and from 0 to 30 decimal degrees East.

We used the maximum entropy approach [40], which repeatedly performed better than other methods using species presence data in comparative studies [41], [42]. Three sets of SDMs – one with the western record set (hereafter western distribution model), one with only Carpathian records (hereafter Carpathian distribution model) and one model including all records together (hereafter entire distribution model) were generated with MAXENT 3.3.3 (, with automatically sampled random background points and species records split into 75% model training and 25% model evaluation data. Model quality was evaluated by a bootstrap validation and the Area Under Curve method (AUC>0.7; [43]). The logistic output of 100 model replicates was transformed to species maps showing the potential distribution of C. irregularis under current climate conditions and under LGM conditions. Species maps were processed in DIVA-GIS v7.1.7 ( We compared the Maxent outputs of the three model sets using the software ENMtools [44]. Niche overlap was quantified by two metrics of niche similarity introduced by Warren et al. (2008) [45]: D (Schoener's statistic for niche overlap; [46]) and I (derived from Hellinger distance). Both measures range from 0 (no niche overlap) to 1 (identical niches). In addition, a niche identity test was run to evaluate niche similarity between the western and the Carpathian distribution model [44].

Sample collection and molecular methods

In total 91 specimens of the study species were collected from 23 localities across its entire distribution range (Table 1, Fig. 1). In spite of conscientious efforts and communications with numerous coleopterologists, no fresh beetles were available from certain parts of the distribution range (possibly due to the recent decline mentioned above). Adult beetles were collected from hibernation in autumn or caught by pitfall trapping during the activity periods. The specimens were preserved at –80°C or in pure ethanol at –20°C until DNA extraction.

Table 1. Populations of Carabus irregularis sampled for genetic analyses and occurrence data used for species distribution models.

Mitochondrial genomes within the genus Carabus can be influenced by horizontal gene flow (secondary contact after speciation) [47][49]. In natural habitats, C. irregularis seems to produce hybrids with related species (members of the same subgenus Platycarabus: C. creutzeri and C. depressus) [50]. Therefore, we used six individuals of the latter species from three different sites as outgroups for rooting phylogenetic trees and for checking the reliability of mitochondrial DNA (mtDNA) as a phylogenetic marker. In comparative studies, the two different mitochondrial genes studied, CO1 and ND5, have outperformed nuclear markers in terms of giving a reliable estimates in timing of splits and phylogenetic reconstruction because of their maternal heritability [51]. Moreover, CO1 is a valuable barcoding marker in animals [19] which makes it possible to study the performance of such a marker system on an infraspecific level.

Genomic DNA was extracted from one femur of each specimen by CTAB (Cetyldimethylethyl-Ammoniumbromid) lysis [52], where we reached a DNA concentration of min. 200 ng/ µl. Two regions of mtDNA, Cytochrome c Oxidase I (CO1) and NADH Dehydrogenase 5 (ND5), were amplified and sequenced using the primers Jerry (C1-J-2183) 5′-CAA CAT TTA TTT TGA TTT TTT G-3′ and Pat (L2-N-3014) 5′-TCC AAT GCA CTA ATC TGC CAT ATT A-3′ [53] for CO1 and ND5-His (V1.06-1) 5′-CCT GTT TCT GCT TTA GTT CA-3′ and ND5-Phe (V1.04-4) 5′-GTC ATA CTC TAA ATA TAA GCT A-3′ [54] for ND5. Polymerase chain reactions (PCRs) were carried out on a Thermal Cycler TGradient (Biometra GmbH, Göttingen, Germany). The PCR mix of 20 µl contained 1 µl of extracted DNA, 2.5 µl of 10× DreamTaq Green Buffer and 0.4 µl dNTP mixture (10 mM each; both Fermentas GmbH, St. Leon-Rot, Germany), 0.1 µl of each primer (50 pmol/ µl; Biomers GmbH, Ulm, Germany), 0.2 µl DreamTaq Polymerase (5 U/ µl; Fermentas GmbH, St. Leon-Rot, Germany) and 16.2 µl of DNA-free water. Cycling conditions for CO1 started with an initial denaturation for 5 min at 94 °C, followed by 35 cycles of denaturation for 45 at 94°C s, annealing for 45 s at 50°C and extension for 1 min at 72°C. Final elongation was performed for 8 min at 72°C. For ND5 the following PCR program was used: initial denaturation at 94 °C for 5 min, 35 cycles (1 min at 94°, 1 min. at 50°C, 2 min at 70°C), then final elongation for 7 min at 70°C.

PCR products were visualised on 2% agarose gels, dyed with Roti-Safe Gelstain (Carl Roth GmbH, Karlsruhe, Germany) and then purified with the GeneJET PCR Purification Kit (Fermentas GmbH, St. Leon-Rot, Germany). The manufacturer's protocol was slightly modified: only 50 µl elution buffer were used for final elution to ensure a minimum DNA concentration of 100 ng/ µl. Purified PCR products were again checked on agarose gels. Sequencing (both forward and reverse strands) of PCR products was carried out at the laboratory of LGC Genomics (Berlin, Germany) using an Automatic Sequencer 3730 xl (Applied Biosystems, Foster City, USA).

Sequence alignment and phylogenetic analysis

Forward and reverse sequences were assembled and manually corrected using the software Geneious v4.8.5 [55]. Sequences were aligned and checked for reading frame errors in MEGA5 [56]. Prior to phylogenetic analysis, MrModeltest v2.3 [57] was used to identify the best nucleotide substitution model for each mtDNA fragment (CO1 and ND5) ranked by Akaike Information Criterion (AIC). For each gene, we performed an independent run of 60×106 generations, sampling one tree every 1,000 generations. Based on the selected evolutionary models, we investigated phylogenetic relationships within the species C. irregularis (including two outgroups) using Bayesian inference (BI), maximum likelihood (ML) and maximum parsimony (MP). The BI tree was calculated by MrBayes v3.2.1 [58] running 2×106 generations, while reaching an average standard deviation of split frequencies<0.01. The ML method was conducted in RAxML v7.3.0 [59] performing 10,000 thorough bootstrap replicates. MP analysis was carried out in MEGA5 [56] using 5,000 bootstrap replicates. Divergence times were estimated by using strict molecular clocks in BEAST v1.7.1 [60] considering gene-specific mtDNA substitution rates for the two studied loci (1.45% My−1 for CO1, and 1.59% My−1 for ND5) based on rates determined for other carabid beetle species [51], [61]. Population size through time was applied by a Yule speciation model, and Markov-Chain-Monte-Carlo (MCMC) chain length was set to 200×106, whereas other parameters remained in default settings.

Intra-specific divergence and genetic diversity

Measures of DNA sequence variation within different hierarchical levels (populations, subclades and clades) and in all sequences together (total) were computed for concatenated sequences and for CO1 and ND5 separately using DnaSP v5.10.01 [62]. Variation measures included the number of polymorphic/segregating sites (S), number of haplotypes (H), haplotype diversity (h) and nucleotide diversity (π). In addition, we performed an analysis of molecular variance (AMOVA) using ARLEQUIN v3.5 [63] to estimate the level of genetic differentiation (ΦST) among populations, subclades and clades. We used the median-joining approach [64] implemented in the NETWORK Software v4.6.1.0 ( to illustrate phylogenetic and geographic patterns in haplotype diversity. Intra-specific genetic divergence in CO1 – both within and between major clades – and inter-specific genetic divergence in the studied outgroup species of the subgenus Platycarabus were calculated by pairwise comparisons (4,005 intra-specific pairs and 554 inter-specific pairs) using Kimura's two parameter (K-2P) model in MEGA.


Current and past predicted distribution

High average training AUC values over all replicate runs for all models – the entire (AUC = 0.92±0.001), the western (AUC = 0.92±0.001) and the Carpathian distribution model (AUC = 0.98±0.002), reflected the high accuracy of our models. Niche overlap was high between the entire distribution model and the western distribution model (D = 0.85; I = 0.98) and low between the entire and the Carpathian distribution model was low (D = 0.39; I = 0.71). Lowest niche overlap was found between the western and the Carpathian distribution model (D = 0.29; I = 0.58). The potential (a) current and (b) past (LGM) distribution of C. irregularis computed with the entire distribution model is documented in Fig. S2. Fig. 2 shows maps of the potential (a and b) current and (c and d) past (LGM) distribution of C. irregularis calculated with the western and the Carpathian record set, respectively. The SDMs for the western European and the Carpathian populations differed markedly from each other.

Figure 2. Western European distribution model for (a) current and (c) past climate conditions and Carpathian distribution model for (b) current and (d) past climate conditions.

Past distribution is shown for the Last Glacial Maximum (LGM, 21,000 BP). Warmer colors show areas with better predicted conditions (green: p>0.3, yellow: p>0.5, red: p>0.7).

The western distribution model for current climate showed the most suitable conditions (p>0.7) in the Jura and Vosges Mountains in eastern France, parts of the Alps (in Switzerland, Germany, Austria and Slovenia), the Southern Carpathians (Transylvanian Alps) and the Apuseni Mountains in Romania as well as in the western Balkan (Dinaric Alps). Moderate suitability (0.5<p<0.7) was found in the lower mountain regions of southern Germany (e.g. Black Forest, Bavarian Forest) and in the low mountain range of northwestern and Central Germany (i.e. Teutoburg Forest, Weser Mountains and Harz) and also in the northern Carpathians (Tatra Mountains). Low suitability (0.3<p<0.5) was attributed to the eastern Carpathians (Romania, Ukraine), the Ardennes (France, Belgium, Luxembourg), large parts of Germany and also up to Denmark. Areas with high suitability rates were always surrounded by areas with lower rates, resulting in a disjunct distribution of suitable habitats. Thus, the model was widely consistent with the documented actual distribution of C. irregularis. However, the study species is not known from the regions with low suitability rates, with the exception of the Romanian and Ukrainian Carpathians.

The Carpathian distribution model for current climate conditions (Fig. 2b) displayed the most suitable conditions in parts of the Alps and several high-altitude regions in the Carpathian mountain system.

The past predictive model of C. irregularis using the western European data set showed the most suitable conditions (p>0.7) during the LGM along the edges of the large Alpine glacier: in the southwestern Alps, the Vosges and Jura and in the Ligurian part of the Apennines. Well suited climatic conditions were also shown in a large area on the Balkan Peninsula (large parts of the Dinaric Alps) and in a small part of the northwestern German low mountain range. Moderately suitable regions were found in the central Apennines, in the eastern parts of the Alps and in the western Balkan. In the Carpathians, models for LGM conditions showed only small regions with low suitability rates (see arrows in Fig. 2c).

The past distribution model for the Carpathian data set (Fig. 2d) illustrated large areas with highly suitable climate conditions in the southern parts of the Eastern Alps and also high suitability rates for the Tatra Mountains and for some parts of the Southern and the Serbian Carpathians.

Sequence characteristics

For 97 individuals (including outgroups), the concatenated mtDNA sequence matrix included 1785 bp: partitioned into 786 bp at the CO1 and 999 bp at the ND5 locus. Sequences were submitted to GenBank (; see GenBank Accession numbers in Table S3). From these mtDNA sequences, we detected 31 haplotypes in C. irregularis and 5 in the outgroups (2 C. depressus, 3 C. creutzeri). In total, 144 sites were variable (67 in CO1, 77 in ND5) and 127 sites were parsimony-informative (CO1: 56, ND5: 71), with no significant difference between the two studied loci (Table S4). Both mitochondrial loci were heavily biased toward high A+T-contents of 71.8% (CO1) and 80.6% (ND5).

Phylogenetic analyses and divergence time estimation

For both loci, the best evolutionary model identified by the AIC was the general time reversible model with Gamma distributed rate variation across sites and a proportion of invariable sites (GTR+G+I) [65]. Since there has been concern about estimation problems in this model type due to interactions between the proportion of invariable sites and the Gamma distribution [59], we also checked the model without +I and chose this simpler model for the ML analysis in RAxML.

BI, MP and ML yielded highly congruent phylogenetic trees with the same major nodes and branching order. The majority-rule consensus tree (Fig. 3) displays the Bayesian posterior probability (BPP), parsimony bootstrap percentages (PB) and maximum likelihood bootstrap percentages (MLB), where BBP ≥ 0.5 and bootstrap values (PB and MLB)≥50% represent well-supported nodes [66]. In addition, nodal ages and 95% confidence intervals are represented for each node.

Figure 3. Majority-rule consensus tree for Carabus irregularis.

The tree shows results for the Bayesian posterior probability (BPP)/the parsimony bootstrap percentages (PB)/the maximum likelihood bootstrap percentages (MLB) for well-supported nodes. Divergence times are displayed with median values and 95% confidence intervals for nodal ages (bars) in million years BP (MY; scale at the bottom). See Table 1 for population abbreviations and Figure 1 for clade abbreviations and colour codes.

In the phylogenetic tree, all basal nodes had very high support (BPP≥0.98; PB≥93; MLB≥65) and branches represent two major clades with the following groups of populations: (A) an eastern clade including two Carpathian subclades and (B) a western clade. The eastern clade covers (Ae) an eastern Carpathian subclade (populations Car I and Car II) and (As) a Southern Carpathian subclade (Car III, Car IV and Car V). The western clade comprises (Bn) a northern/central German subclade (Teu, Sue, Dei, Hz I, Hz II, Hai I, Hai II and Hai III), (Bs) a subclade including populations from Southern Germany (Bla, Sw I, Sw II, Sw III), the French Jura mountains (Jur I, Jur II) and Slovenia (Slo I) and (Bsl) a purely Slovenian subclade (Slo II, Slo III; Fig. 3). Assumed divergence times spanning the last 2 million years (Fig. 3) indicated that major clades and subclades within C. irregularis were formed long before the last glacial period. The separation of C. irregularis from other species of the subgenus Platycarabus was found to have taken place earlier than the separation of the species used as outgroups (C. creutzeri and C. depressus). Additionally, the split into the two major lineages within C. irregularis appears to be older than the split of the two species C. creutzeri and C. depressus. Only one of the three subspecies was monophyletic (C. i. montandoni), while C. i. irregularis and C. i. bucephalus are paraphyletic.

Intra-specific diversity and geographic distribution of genetic diversity

The overall haplotype diversity (h = 0.882 ± 0.0006) and nucleotide diversity (π = 0.017±0.002) were high, with clade A showing higher haplotype and nucleotide diversity (h = 0.901; π = 0.016) than clade B (h = 0.838; π = 0.009; Table 2). In general, haplotype and nucleotide diversity tended to be higher in populations from the southern part of the distribution range (Fig. S3). Within clades, subclade As showed higher diversity than Ae and the subclades Bsl and Bs are more diverse than subclade Bn (see Table 3 for a summary of computed diversity measures of different hierarchical levels). The mean intra-specific variation in CO1 (K-2P; d = 0.027 ± 0.004) was also very high. The inter-clade K-2P divergence was 0.042 ± 0.004 and the frequency distribution of intra-specific genetic divergence in CO1 within and between major clades showed a clear pattern (Fig. 4). On the population level, K-2P distance values ranged from 0.0 to 0.016 between populations within clade B to 0.022 between populations within clade A. Between populations of clade A and B, genetic distance ranged from 0.040 to 0.044 (Table S5).

Figure 4. Frequency distribution of intra-specific genetic divergence in CO1 – within (light grey) and between major clades (dark grey).

Pairwise distances (4005 intra-specific comparisons within C. irregularis) were calculated using Kimura's two parameter (K-2P) model.

Table 2. Diversity statistics based on concatinated CO1+ND5 sequences (1785 bp) of C. irregularis populations (sampled with more than one specimen) and in subclades and clades (as determined by the phylogenetic tree).

Table 3. Results of Analysis of Molecular Variance (AMOVA) for different hierarchical levels: subclades, clades and subspecies.

In addition, we detected a significant differentiation structure among all populations (ΦST = 0.935, p<0.001) and no shared haplotypes between the clades, subclades and subspecies but between populations within clades and subclades. Additionally, AMOVAs for different hierarchical groups (subclades, clades, subspecies) revealed significant genetic variance among these hierarchical groups, among populations and also within populations of the hierarchical groups. Most of the diversity was observed among subclades (88.1%), among clades (72.3%) and among subspecies (63.7%), while low percentage of variance was detected within populations (3.1 to 4.8%; Table 3). Differentiation between the single populations of clade B ranged from lower values in geographically closer populations (ΦST = 0.10) to higher values in geographically distant populations (ΦST = 0.99). Populations of clade A also showed very high ΦST values (>0.9; Table S5), even though the sampled populations were geographically closer than those from clade B. Inter-clade ΦST values comparing populations of clade A and clade B ranged from 0.92 to 0.99. The haplotype network (Fig. 5) supports these findings and illustrates a similar pattern: high genetic distance between the two geographically distant (∼1,500 km) major clades, but also large genetic distance and a high number of substitutions between the geographically closer (∼300 km) Carphathian subclades (Ae and As). In contrast, populations within subclade Bs showed very low genetic distance, although geographic distance is quite large (between Bla/Sw I-III/Jur I-II and SloI ∼900 km).

Figure 5. Median-joining network of mtDNA haplotypes based on concatenated CO1+ND5 sequences.

Circles represent haplotypes; circle size and numbers in brackets indicate the haplotype frequency within our samples. Small black circles/nodes indicate intermediate haplotypes between observed haplotypes. Haplotype circles are filled corresponding to the colour code for the subclades also used in Figure 1 and 3.


Potential glacial refugia

Our past predictive models showed several independent mountainous regions across Central and southern Europe with highly suitable climate conditions for C. irregularis which could have acted as glacial refugia during the LGM. Some of these potential refugia (e.g. in the Dinaric Alps) conformed to classical refugia on the southern European peninsulas [1], [67], [68] and other refugia at the margins of the Alps and in neighbouring mountain ranges (e.g. the Jura, the Black Forest), some of which have quite recently been recognised as refugia [69][71]. In addition, our model pointed to refugial areas far north of the southern European peninsulas and the Alps. So-called cryptic refugia have been assumed for a growing group of other species (plants as well as animals) [5], [72].

While the southern European peninsulas are typical glacial refugia of current lowland plant and animal species [2], [73], mountain species potentially persisted in refugia along the border of the Alps when higher altitudes were covered with ice sheets [69], [74]. According to our results, C. irregularis seems to have survived in refugia of both low- and high-altitude species. Our models projected relatively large and consistent areas of suitable and moderately suitable conditions during the Last Glacial Maximum (comprising the maximum extension of the continental glaciers), suggesting that C. irregularis was relatively widespread during glaciations. In general, numerous mountain species (even though mostly defined as species with their main distribution above the timber line) are assumed to have survived in several spatially restricted and isolated refugia rather than in single extensive and consistent refuge areas ([75] and references therein; [76]). Consequently, C. irregularis, a cold-adapted montane to subalpine species, can be expected to have survived the last ice age in multiple spatially restricted refugia within the projected potential distribution area in the Alps and some neighbouring mountain ranges. It appears less important whether there is a significant difference between the historical and the recent distribution range. The fact that the study species tended to have cryptic refugia north of the Alps seems to be of greater importance. While cold-tolerant species, e.g. some mammalian herbivores, have long been believed to have survived in Central European steppe tundra surrounding the glacial ice-sheets, cryptic refugia were surprisingly also found for more temperate animal and tree species [72], [77]. Topographically sheltered sites are assumed to have provided suitable climate conditions and stable microclimates for species with differing habitat requirements [3], [72].

Species distribution modelling is a relatively new but already frequently used method for inferring species potential distribution ranges from occurrence data. The method has been criticised for including abiotic factors only and not taking biotic parameters such as competition and dispersal into account, which might explain discrepancies between the documented range and the projected species range [15], [16]. Nonetheless, the high model quality of our simulations (AUC = 0.92 and 0.98) represented high consistence between actual and predicted occurrence. Another uncertainty, especially of past predictive modelling, is the accuracy of past climate data. We used data processed by the downscaling method, which is currently the best large-scale data available but still requires further validation [18]. Additionally, the modelling approach is based on niche conservatism and does not consider changes in habitat requirements. Nevertheless, past predictive models are useful tools to visualise species distribution patterns and potential glacial refugia. More traditional approaches to identifying refugia are pollen and (sub-)fossil analyses (e.g. [78]). Since no subfossils are known for Platycarabus species (cf. [79]), we chose a molecular method (mtDNA analysis) to supplement our findings from modelling and to review assumptions concerning potential glacial refugia.

Genetic evidence for glacial refugia

Phylogenetic analyses (using BI, MP and ML approaches) revealed two well-supported and spatially congruent major clades: one Central European (B) with three subclades and one Carpathian clade (A) with two subclades. The ages of all major phylogenetic splits predating the last ice age, the absence of shared haplotypes and high genetic intra-specific differentiation between clades as well as subclades implied that the ancestors of phylogenetic (sub)clades in C. irregularis survived many glacial-interglacial cycles isolated from each other. The Central European clade consisted of three independent and genetically distinct Central European subclades: one subclade with a more southern distribution (Bs: with populations from southern Germany, France and Slovenia), one Slovenian subclade (Bsl) and another one including haplotypes from the northern part of the species' distribution range (Bn: central to northern Germany).

Thus, the assumption of multiple refugia of C. irregularis in Central Europe inferred from past predictive models was supported by our phylogenetic results, with at least two refugia located close to the Alps or on the Balkan Peninsula and one northern refugium, potentially in central Germany. However, phylogenetic analyses indicated the existence of independent Carpathian refugia, which were hardly evident from the western European distribution model for LGM conditions (Fig. 2c, arrows). The past predictive model using the eastern records showed several potential refugia in the Carpathian mountain system (Fig 2d). These potential refugia coincide with the finding of two subclades of clade A one from the eastern (Ae) and one from the southern Carpathians (As), with divergence times dating back more than one million years BP.

As well as refugia in northern Central Europe, Carpathian refugia are not accepted as classical refugial areas and phylogeographic studies including the Carpathian region are still rare. However, since the Carpathians have faced a different glacial history with more local glaciers and larger areas of suitable alpine habitats than the well-studied Alps [9], [80], investigation of the southeastern European mountain range appears likely to be very interesting. Some recent studies (on mountainous plant and animal species) detected low-level differentiation between the Carpathians and the Alps and assumed the Carpathians to have been colonised by adjacent refugia e.g. in the western Alps [8] or the Balkans [71]. In contrast, an increasing number of phylogeographic studies have found very old and divergent phylogenetic lineages, suggesting a major refugium for plants [9], [80] as well as for vertebrates [81][84] in the Carpathian region. In congruence with our phylogeographic results for C. irregularis, alpine plants, namely Hypochaeris uniflor and Campanula alpina, showed distinct and strongly differentiated genetic groups within the Carpathians, also pointing to long-term isolation and restricted gene flow between several areas within the mountain system [9], [80], [85]. While Ronikier et al. [80] discovered the strongest differentiation between western and eastern Carpathian samples, results of Mráz et al. (2007) [85] resembled our findings in that they also showed high differentiation between southern and eastern Carpathian populations. Examples from animals also revealed genetically distinct and long-term isolated phylogenetic groups in the Carpathians [81][84].

DNA barcoding and systematic delineation

The phylogenetic results of our study agreed only in part with subspecies systematics: Whereas the Carpathian subspecies C. i. montandoni was monophyletic (clade A), the other two subspecies (C. i. irregularis and C. i. bucephalus) were paraphyletic (clade B) and we suggest to rank the populations belonging to clade B as one subspecies. The high proportion of polymorphic sites within ND5 (7.70%) and also within the CO1 loci (8.52%), which is commonly used as a barcoding marker [19] – also for ground beetles [86], questions the conventional taxonomy. DNA barcoding using the CO1 locus for rapid species delineation [87], implies that inter-clade K-2P distances higher than 2.8 to 3.4% indicate species rather than subspecies or populations. Following these thresholds, the divergence of 0.042 ± 0.004 between clade A (taxonomic subspecies C. i. montandoni) and B (irregularis s. str.) points to independent species rather than subspecies of the same species. Although the commonly used Kimura's two-parameter substitution model (K-2P) has been repeatedly criticised [88], [89], estimates using other substitution models showed similar results for species delineation [89]. For C. irregularis, the calculation of other (simpler as well as more complex) genetic distance metrics shows a consistent pattern with K-2P (p distance between clade A and B = 0.042; Tamura-Nei distance between clade A and B = 0.044). However, due to the fact that substitution rates can vary even within the same species [90] it is necessary to exercise caution when defining a species from nucleotide rates alone.

Using morphological and molecular data Casale et al. [91] reconstructed the phylogeny of the subgenus Platycarabus. Following their phylogeny C. irregularis shows the strongest differentiation from all other species within the subgenus (basal split) and C. depressus and C. creutzeri are the second major diverged species pair. Surprisingly, the split between the two main lineages of C. irregularis we found is older than the divergence between the latter species pair. This underlines both (i) the possible species status of the Carpathian populations and the (ii) need to incorporate a population-based approach to understand the phylogeny of these highly differentiated ground beetles with their ancient lineages, even below or at the species level.


Past predictive modelling and phylogenetic analyses acted as supplements and imply that C. irregularis survived the last glacial periods in long-term isolated classical refugia on the edges of the Alps as well as in other refugia (in Central Europe and some parts of the Carpathian Mountains). Altogether our results indicate that the Carpathian subspecies of C. irregularis differs not only in terms of its geographical distribution but also genetically and (due to differing habitat preferences as revealed from the poorly overlapping SDMs) ecologically from rest of the species. Thus we conclude that C. irregularis comprises at least two evolutionarily significant units (ESUs, sensu Moritz, [92]) indicating that the species has a very interesting history, which should be investigated in more detail by further phylogeographic analyses – also in the framework of the Platycarabus group.

Supporting Information

Figure S1.

Means and 95% confidence intervals of altitudinal distribution of the western and the Carpathian subspecies of Carabus irregularis .


Figure S2.

Entire distribution model for (a) current and (c) past climate conditions. Past distribution is shown for the Last Glacial Maximum (LGM, 21,000 BP). Warmer colors show areas with better predicted conditions (green: p>0.3, yellow: p>0.5, red: p>0.7).


Figure S3.

Nucleotide diversity of 19 studied Carabus irregularis populations. Black circles symbolise values higher than the global mean nucleotide diversity; white circles stand for lower than mean diversity. Circle size is proportional to the difference from the mean diversity. Three populations (Dei, Jur I, Car IV) are omitted due to insufficient sample sizes.


Table S1.

Species records of Carabus irregularis and subspecies affiliation (according to Turin et al. (2003)) used for species distribution modelling.


Table S2.

Climate data (monthly temperature and rainfall values) generated according to Hijmans et al. (2005) [34] to 19 bioclimatic variables.


Table S3.

Population and specimen codes and GenBank Accession numbers of the analysed CO1 and ND5 sequences.


Table S4.

Summary of statistical analyses on measures of DNA sequence variation within clades, subclades and in all sequences together (total) for CO1 and ND5 separately: N = sample size, S = number of polymorphic/segregating sites, Eta = total number of mutations, k = average pairwise nucleotide difference per sequence, H = number of haplotypes, h = haplotype diversity, π = nucleotide diversity, θG = mutation parameter per sequence).


Table S5.

Genetic differentiation (ΦST; in the upper part of the matrix) and intra-specific genetic divergence (K-2P distances; in the lower part of the matrix) of Carabus irregularis. K-2P distances greater than the threshold value assumed for species delineation by Wiemers and Fiedler [87] are printed in bold.



We are very grateful to Oleg Aleksandrowicz, Wouter Dekoninck, Božidar Drovenik, Michael Andreas Fritze, Matthias Hartmann, Günter Hofmann (Hessische Koleopterologen), Pavel Hrusa, Andrej Kapla, Damir Kovac (Forschungsinstitut Senckenberg), Roman Lohaj, Dirk Maes, Werner Marggi, Wolfgang Paill, Alexandr Putchkov, Lucija Šerić-Jelaska, Mieczysław Stachowiak, Christiane Schilling (Nds. Landesmuseum Hannover) and Friedrich Weber who provided us with useful information on and records of C. irregularis. In addition, we thank Nora Drescher, Nils Homburg, Andrej Kapla, Bernard Lassalle, Jacqueline Loos, Tamar Marcus, Maja Rabe, Florian Schäfer, Reinhard Schäfer and Friedrich Weber for their help collecting samples, Marianne Peters for technical assistance in the lab and Henrik von Wehrden and Patric Brandt for their help in species distribution modelling.

Some of the beetles were collected during the field campaign of the PopGeneLand project (AS 309/3-1) funded by the DFG Priority Program 1374 "Infrastructure-Biodiversity-Exploratories". We also acknowledge the field work permits issued by the responsible environmental authorities.

Author Contributions

Commented on previous versions of the manuscript: CD MMG LR AV TA. Conceived and designed the experiments: KH CD TA. Performed the experiments: KH CD TA. Analyzed the data: KH CD. Contributed reagents/materials/analysis tools: KH CD MMG LR AV TA. Wrote the paper: KH.


  1. 1. Hewitt GM (2000) The genetic legacy of the Quaternary ice ages. Nature 405: 907–913.
  2. 2. Taberlet P, Fumagalli L, Wust-Saucy A-G, Cosson J-F (1998) Comparative phylogeography and postglacial colonization routes in Europe. Molecular Ecology 7: 453–464.
  3. 3. Stewart JR, Lister AM, Barnes I, Dalen L (2010) Refugia revisited: individualistic responses of species in space and time. Proc Biol Sci 277: 661–671.
  4. 4. Bhagwat SA, Willis KJ (2008) Species persistence in northerly glacial refugia of Europe: a matter of chance or biogeographical traits? Journal of Biogeography 35: 464–482.
  5. 5. Provan J, Bennett KD (2008) Phylogeographic insights into cryptic glacial refugia. Trends in Ecology & Evolution 23: 564–571.
  6. 6. Valtuena FJ, Preston CD, Kadereit JW (2012) Phylogeography of a Tertiary relict plant, Meconopsis cambrica (Papaveraceae), implies the existence of northern refugia for a temperate herb. Molecular Ecology 21: 1423–1437.
  7. 7. Knopp T, Merilä J (2009) The postglacial recolonization of Northern Europe by Rana arvalis as revealed by microsatellite and mitochondrial DNA analyses. Heredity 102: 174–181.
  8. 8. Kropf M, Comes HP, Kadereit JW (2012) Past, present and future of mountain species of the French Massif Central - the case of Soldanella alpina L. subsp. alpina (Primulaceae) and a review of other plant and animal studies. Journal of Biogeography 39: 799–812.
  9. 9. Ronikier M (2011) Biogeography of high-mountain plants in the Carpathians: An emerging phylogeographical perspective. Taxon 60: 373–389.
  10. 10. Schmitt T, Haubrich K (2008) The genetic structure of the mountain forest butterfly Erebia euryale unravels the late Pleistocene and postglacial history of the mountain coniferous forest biome in Europe. Molecular Ecology 17: 2194–2207.
  11. 11. Muster C, Berendonk TU (2006) Divergence and diversity: lessons from an arctic-alpine distribution (Pardosa saltuaria group, Lycosidae). Molecular Ecology 15: 2921–2933.
  12. 12. Schmitt T (2009) Biogeographical and evolutionary importance of the European high mountain systems. Frontiers in Zoology 6: 9.
  13. 13. Matern A, Drees C, Vogler AP, Assmann T (2010) Linking genetics and ecology: reconstructing the history of relict populations of an endangered semi-aquatic beetle. In: Habel JC, Assmann T, editors. Relict Species, Phylogeography and Conservation Biology.Heidelberg:Springer. pp. 253-266.
  14. 14. Schmidt J, Opgenoorth L, Höll S, Bastrop R (2012) Into the Himalayan Exile: The Phylogeography of the Ground Beetle Ethira clade Supports the Tibetan Origin of Forest-Dwelling Himalayan Species Groups. PLoS ONE 7: e45482.
  15. 15. Elith J, Leathwick JR (2009) Species distribution models: Ecological explanation and prediction across space and time. Annual Reviews in Ecology, Evolution and Systematics 40: 677–697.
  16. 16. Rebelo H, Froufe E, Brito JC, Russo D, Cistrone L, et al. (2012) Postglacial colonization of Europe by the barbastelle bat: agreement between molecular data and past predictive modelling. Molecular Ecology 21: 2761–2774.
  17. 17. Vega R, Fløjgaard C, Lira-Noriega A, Nakazawa Y, Svenning J-C, et al.. (2010) Northern glacial refugia for the pygmy shrew Sorex minutus in Europe revealed by phylogeographic analyses and species distribution modelling. Ecography: 260-271.
  19. 19. Hebert PDN, Ratnasingham S, deWaard JR (2003) Barcoding animal life: cytochrome c oxidase subunit 1 divergences among closely related species. Proceedings of the Royal Society B 270: 96–99.
  20. 20. Turin H, Penev L, Casale A (2003) The genus Carabus in Europe: A synthesis. Fauna Europaea Evertebrata. Sofia: Pensoft.
  21. 21. Bousquet Y, Brezina B, Davies A, Farkac J, Smetana A (2003) Tribe Carabini Latreille, 1802. In: Löbl I, Smetana A, editors.Catalogue of Palaearctic Coleoptera Vol 1: Archostemata, Myxophaga, Adephaga. Stenstrup: Apollo Books.pp. 118-201.
  22. 22. Casale A, Kryzhanovskij OL (2003) Key to the adults. In: Turin H, Penev L, Casale A, editors. The genus Carabus in Europe - A synthesis. Pensoft: Sofia.pp. 73-124.
  23. 23. Lorenz WMT (2003) Rote Liste gefährdeter Lauf- und Sandlaufkäfer (Coleoptera Carabidae s. l.) Bayerns. Schriftenreihe Bayerisches Landesamt für Umweltschutz 166: 102–111.
  24. 24. Assmann T, Dormann W, Främbs S, Gürlich S, Handke K, et al. (2003) Rote Liste der in Niedersachsen und Bremen gefährdeten Sandlaufkäfer und Laufkäfer (Coleoptera: Cicindelidae et Carabidae) mit Gesamtartenverzeichnis, 1. Fassung vom 1.6.2002. Informationsdienst Naturschutz Niedersachsen 23: 70–95.
  25. 25. Desender K, Dekoninck W, Maes D (2008) Een nieuwe verspreidingsatlas van de loopkevers en zandloopkevers (Carabidae) in België. Brussel: Instituut voor Natuur- en Bosonderzoek.
  26. 26. Franzen M (1995) Nachweise von Carabus irregularis Fabricius, 1792 aus Rheinland-Pfalz (Coleoptera: Carabidae). Fauna Flora Rheinland-Pfalz 8: 5–15.
  27. 27. Gebert J (2006) Die Sandlaufkäfer und Laufkäfer von Sachsen, Beiträge zur Insektenfauna Sachsens, Teil 1, Band 4 (Cicindelini-Loricerini). Dresden: Entomologische Nachrichten und Berichte. Beiheft 10: 1–180.
  28. 28. Hejda R (2011) Map of distribution of Carabus irregularis in the Czech Republic. In: Zicham O, editor. Biological Library - BioLib Retrieved on 2011-10-24 Available: http://wwwbiolibcz/en/taxonmap/id348/.
  29. 29. Sokolár F (1909) Carabus irregularis F. . Entomologische Rundschau 1909: 87–89.
  30. 30. Weber F (1966) Zur Verbreitung von Carabus irregularis Fabr. im Teutoburger Wald und Eggegebirge (Westfalen). Entomologische Blätter für Biologie und Systematik der Käfer (Krefeld) 62: 1–5.
  31. 31. Hartmann M (1998) Die Verbreitung von Carabus irregularis F., C. linnei PANZER und C. sylvestris PANZER in Thürningen (Coleoptera, Carabidae). Thüringer Faunistische Abhandlungen V: 147–152.
  32. 32. Köhler F, Flechtner G (2007) Coleoptera (Käfer). In: Dorow WHO, Kopelke J-P, editors. Naturwaldreservate in Hessen 7/22 Hohestein Zoologische Untersuchungen 1994-1996, Teil 2.Wiesbaden.pp. 101-192.
  33. 33. Flechtner G (2000) Coleoptera (Käfer). In: Flechtner G, Dorow WHO, Kopelke J-P, editors. Naturwaldreservate in Hessen 5/22 - Niddahänge östlich Rudingshain Zoologische Untersuchungen 1990-1992 Mitt Hess Landesforstverwaltung.Wiesbaden.pp. 1-349.
  35. 35. Thiele HU (1977) Carabid beetles in their environments. Berlin: Springer.
  36. 36. Araújo MB, Guisan A (2006) Five (or so) challenges for species distribution modelling. Journal of Biogeography 33: 1677–1688.
  37. 37. Rödder D, Dambach J (2010) Bioclimatic Models as predictive GIS tools for the identification of potential refugia and possible migration pathways. In: Habel J, Assmann T, editors.Relict Species, Phylogeography and Conservation Biology ..Heidelberg : Springer. pp. 373-384.
  39. 39. Anderson RP, Raza A (2010) The effect of the extent of the study region on GIS models of species geographic distributions and estimates of niche evolution: preliminary tests with montane rodents (genus Nephelomys) in Venezuela. Journal of Biogeography 37: 1378–1393.
  40. 40. Phillips SJ, Anderson RP, Schapire RE (2006) Maximum entropy modeling of species geographic distributions Ecological Modelling. 190: 231–259.
  43. 43. Fielding AH, Bell JF (1997) A review of methods for the assessment of prediction errors in conservation presence/absence models. Environmental Conservation 24: 38–49.
  45. 45. Warren DL, Glor RE, Turelli M (2008) Environmental niche equivalency versus conservatism: quantitative approaches to niche evolution. Evolution 62: 2868–2883.
  46. 46. Schoener TW (1968) Anolis lizards of Bimini: resource partitioning in a complex fauna. Ecology 49: 704–726.
  47. 47. Sota T, Vogler AP (2001) Incongruence of mitochondrial and nuclear gene trees in the carabid beetles Ohomopterus. Systematic Biology 50: 39–59.
  48. 48. Deuve T, Cruaud A, Genson G, Rasplus JY (2012) Molecular systematics and evolutionary history of the genus Carabus (Col. Carabidae). Molecular Phylogenetics and Evolution 65: 259–275.
  49. 49. Brouat C, Meusnier S, Veyrier R, Streiff R (2006) Haldane's rule in Carabus: interspecific mating between C. punctatoauratus and C. splendens using experimental tests and molecular markers. Entomologia Experimentalis et Applicata 120: 189–194.
  50. 50. Breuning S (1937) Monographie der Gattung Carabus L. (VII. Teil). Bestimmungs-Tabellen der europäischen Coleopteren. Troppau, Schlesien: Verlag Edmund Reitter's Nachf. Emmerich Reitter.pp. 1361-1610.
  51. 51. Andujar C, Serrano J, Gomez-Zurita J (2012) Winding up the molecular clock in the genus Carabus (Coleoptera: Carabidae): assessment of methodological decisions on rate and node age estimation. BMC Evolutionary Biology 12: 40.
  52. 52. Doyle JJ, Doyle JL (1987) A rapid DNA isolation procedure for small quantities of fresh leaf tissue. Phytoch Bull 19: 11–15.
  53. 53. Simon C, Frati F, Beckenbach A, Crespi B, Liu H, et al. (1994) Evolution, Weighting, and phylogenetic utility of mitochondrial gene-sequences and a compilation of conserved polymerase chain-reaction primers. Annals of the Entomological Society of America 87: 651–701.
  54. 54. Su Z-H, Tominaga O, Okamoto M, Osawa S (1998) Origin and diversification of hindwingless Damaster ground beetles within the Japanese islands as deduced from mitochondrial ND5 gene sequences (Coleoptera, Carabidae). Molecular Biology and Evolution 15: 1026–1039.
  55. 55. Drummond A, Ashton B, Cheung M, Heled J, Kearse M, et al.. (2009) Geneious v4.7, Available:
  57. 57. Nylander JAA (2004) MrModeltest v2. Program distributed by the author Evolutionary Biology Centre, Uppsala University.
  58. 58. Huelsenbeck J, Ronquist F (2001) MRBAYES: Bayesian inference of phylogenetic trees. Bioinformatics 17: 754–755.
  61. 61. Osawa S, Su Z-H, Imura Y (2004) Molecular phylogeny and evolution of carabid ground beetles: Springer.
  63. 63. Excoffier L, Smouse PE, Quattro JM (1992) Analysis of molecular variance inferred from metric distances among DNA haplotypes: application to human mitochondrial-DNA restriction data. Genetics 131: 479–491.
  64. 64. Bandelt H-J, Forster P, Röhl A (1999) Median-joining networks for inferring intraspecific phylogenies. Molecular Biology and Evolution 16: 37–48.
  65. 65. Rodríguez F, Oliver JL, Marín A, Medina JR (1990) The general stochastic model of nucleotide substitution. Journal of Theoretical Biology 142: 485–501.
  66. 66. Hillis DM, Bull JJ (1993) An empirical test of bootstrapping as a method for assessing confidence in phylogenetic analysis. Systematic Biology 42: 182–192.
  67. 67. Hewitt GM (1999) Post-glacial re-colonization of European biota. Biological Journal of the Linnean Society 68: 87–112.
  68. 68. Hewitt GM (1996) Some genetic consequences of ice ages, and their role in divergence and speciation. Biological Journal of the Linnean Society 58: 247–276.
  69. 69. Schönswetter P, Stehlik I, Holderegger R, Tribsch A (2005) Molecular evidence for glacial refugia of mountain plants in the European Alps. Molecular Ecology 14: 3547–3555.
  70. 70. Schmitt T, Hewitt GM, Müller P (2006) Disjunct distributions during glacial and interglacial periods in mountain butterflies: Erebia epiphron as an example. Journal of Evolutionary Biology 19: 108–113.
  71. 71. Pauls SU, Lumbsch HT, Haase P (2006) Phylogeography of the montane caddisfly Drusus discolor: evidence for multiple refugia and periglacial survival. Molecular Ecology 15: 2153–2169.
  72. 72. Stewart JR, Lister AM (2001) Cryptic refugia and the origins of the modern biota. Trends in Ecology & Evolution 16: 608–613.
  73. 73. Petit RJ, Aguinagalde I, Beaulieu J-Ld, Bittkau C, Brewer S, et al. (2003) Glacial refugia: Hotspots but not melting pots of genetic diversity. Science 300: 1563–1565.
  74. 74. Schmitt T, Augenstein B, Finger A (2008) The influence of changes in viticulture management on the butterfly (Lepidoptera) diversity in a wine growing region of southwestern Germany. European Journal of Entomology 105: 249–255.
  75. 75. Holderegger R, Thiel-Egenter C (2009) A discussion of different types of glacial refugia used in mountain biogeography and phylogeography. Journal of Biogeography 36: 476–480.
  76. 76. Schoville SD, Roderick GK, Kavanaugh DH (2012) Testing the ‘Pleistocene species pump’ in alpine habitats: lineage diversification of flightless ground beetles (Coleoptera: Carabidae: Nebria) in relation to altitudinal zonation. Biological Journal of the Linnean Society 107: 95–111.
  77. 77. Willis KJ, Rudner E, Sümegi P (2000) The Full-Glacial Forests of Central and Southeastern Europe. Quaternary Research 53: 203–213.
  78. 78. Elias SA (1994) Quaternary insects and their environments. Washington,D. C.:Smithsonian Institute.
  79. 79. Breuning S (1932-1937) Monographie der Gattung Carabus L. - Bestimmungstabellen der europäischen Coleopteren. Troppau: E. Reitter.
  80. 80. Ronikier M, Cieslak E, Korbecka G (2008) High genetic differentiation in the alpine plant Campanula alpina Jacq. (Campanulaceae): evidence for glacial survival in several Carpathian regions and long-term isolation between the Carpathians and the Alps. Molecular Ecology 17: 1763–1775.
  81. 81. Kotlik P, Deffontaine V, Mascheretti S, Zima J, Michaux JR, et al. (2006) A northern glacial refugium for bank voles (Clethrionomys glareolus). Proceedings of the National Academy of Sciences of the United States of America 103: 14860–14864.
  82. 82. Deffontaine V, Libois R, Kotlik P, Sommer R, Nieberding C, et al. (2005) Beyond the Mediterranean peninsulas: evidence of central European glacial refugia for a temperate forest mammal species, the bank vole (Clethrionomys glareolus). Molecular Ecology 14: 1727–1739.
  83. 83. Babik W, Branicki W, Crnobrnja-Isailovic J, Cogalniceanu D, Sas I, et al. (2005) Phylogeography of two European newt species - discordance between mtDNA and morphology. Molecular Ecology 14: 2475–2491.
  84. 84. Babik W, Branicki W, Sandera M, Litvinchuk S, Borkin LJ, et al. (2004) Mitochondrial phylogeography of the moor frog, Rana arvalis. Molecular Ecology 13: 1469–1480.
  85. 85. Mráz P, Gaudeul M, Rioux D, Gielly L, Choler P, et al. (2007) Genetic structure of Hypochaeris uniflora (Asteraceae) suggests vicariance in the Carpathians and rapid post-glacial colonization of the Alps from an eastern Alpine refugium. Journal of Biogeography 34: 2100–2114.
  86. 86. Maddison DR (2012) Phylogeny of Bembidion and related ground beetles (Coleoptera: Carabidae: Trechinae: Bembidiini: Bembidiina). Molecular Phylogenetics and Evolution 63: 533–576.
  87. 87. Wiemers M, Fiedler K (2007) Does the DNA barcoding gap exist? – a case study in blue butterflies (Lepidoptera: Lycaenidae). Frontiers in Zoology 4.
  88. 88. Srivathsana A, Meier R (2012) On the inappropriate use of Kimura-2-parameter (K2P) divergences in the DNA-barcoding literature. Cladistics 28: 190–194.
  89. 89. Collins RA, Boykin LM, Cruickshank RH, Armstrong KF (2012) Barcoding's next top model: an evaluation of nucleotide substitution models for specimen identification. Methods in Ecology and Evolution 3: 457–465.
  90. 90. Caccone A, Valerio S (2001) Molecular biogeography of cave life: a study using mitochondrial DNA from bathysciine beetles. Evolution 55: 122–130.
  91. 91. Casale A, Prüser F, Arndt E, Mossakowski D (1998) Phylogenetic relationships in the subgenus Platycarabus Morawitz, 1886 (Coleoptera: Carabidae: Carabini). In: Ball GE, Casale A, Vigna Taglianti A, editors.Phylogeny and classification of Caraboidea (Coleoptera: Adephaga).Torino: Museo Regionale di Scienze Naturali.pp. 429-447.
  92. 92. Moritz C (1994) Defining 'Evolutionarily Significant Units' for conservation. Trends in Ecology & Evolution 9: 373–375.