- Research article
- Open Access
Founder events and pre-glacial divergences shape the genetic structure of European Collembola species
BMC Evolutionary Biologyvolume 16, Article number: 148 (2016)
Climate oscillations in the Cenozoic reduced species richness and genetic diversity of terrestrial and aquatic animals and plants in central and northern Europe. The most abundant arthropods in temperate soils are Collembola that live in almost any soil-related habitat. Extant species show little morphological variation to Eocene fossils, suggesting persistence of species in stable habitats for millions of years. Collembola are able to evade adverse climatic conditions by moving into deeper soil layers and are tolerant to frost and draught. If these adaptations sufficed for surviving glacial periods remains open and needs to be investigated in a phylogeographic context, i.e. investigating spatial structure on molecular level. We investigated the molecular variation of three common species of Collembola at a pan-European scale to identify glacial refuges and post-glacial colonization patterns with three genetic markers.
All genes revealed remarkable genetic structure between but not within populations, suggesting density dependent processes for establishment of populations (founder-takes-all principle), which is common for European animals and plants. In contrast to the post-glacial recolonization patterns of many aboveground organisms, divergence times of most geographic lineages indicate preservation of genetic structure since the Miocene.
Collembola survived severe climatic changes including those during Quatenary glaciation and kept high genetic variance across Europe. Likely the buffering of temperature oscilliations in soil and the ability to evade adverse climatic conditions due to cold-tolerance and horizontal migration enabled Collembola to evade strong selective pressure of abiotic forces.
Collembola are small wingless hexapods that have been among the first arthropods on land, oldest fossils date to the Early Devonian ~400 mya [1, 2]. Fossils from Baltic amber (55–35 mya) have been assigned to extant species [3–6], indicating persistence with little morphological modification over very long periods of time. They are worldwide distributed and are ubiquitous in leaf litter with more than 90 % of the individuals inhabiting the upper 10 cm of the soil . They reach high local density in temperate forests (> 105 ind./m2; ) and significantly contribute to decomposition processes, soil respiration and nutrient cycling [9, 10].
The above-and belowground food web is intimately linked  but also differs in many respects. Mobility of soil arthropods is limited due to the porous structure of soils. Further, abiotic constraints, such as temperature fluctuations and drought, are less severe in soil, and soil-living animals can evade adverse climatic conditions by moving deeper into the soil [12, 13]. Indeed, freezing avoidance but also frost and draught tolerance is widespread in soil invertebrates, including Collembola [14–17]. Consequently, Collembola are little affected by low temperature conditions and in temperate ecosystems they typically remain active during winter [18–20]. This suggests that Collembola in soil may have suffered less from Quaternary climate changes than animals above the ground.
Collembola form part of the decomposer system and predominantly feed on dead organic matter and associated microorganisms [21, 22], resources which likely were at least temporarily available during Quaternary glaciations of central Europe. Therefore, present day populations of central European Collembola may well descend from relict populations that survived glacial periods, rather than from south European populations that recolonized empty habitats when glaciers retreated northwards. Effects of the recurrent and severe abiotic disturbances during Quaternary ice-ages on the genetic structure of soil-living animals has been little investigated [23, 24], disregarding one of the most species-rich terrestrial animal communities.
Few studies investigated the genetic structure of Collembola in Europe and the existing studies focused on large surface-dwelling rather than soil-dwelling species [25–28]. Generally, there is very limited information on phylogeographic patterns of soil-living invertebrates across Europe. The only study available analyzing the genetic structure of soil-living arthropods on a wide geographic area investigated the oribatid mite species Steganacarus magnus . In agreement with the above considerations the results suggest that S. magnus survived Quaternary glaciations in central Europe in cryptic refugia and indicates that phylogeographic patterns of aboveground animals cannot easily be transferred to those of the belowground system. Collembola resemble oribatid mites in various respects but typically are faster reproducing, more mobile and less frost tolerant [19, 29–33].
The ability to survive Quaternary glaciations and the recolonization potential from glacial refuges likely differs between Collembola and oribatid mites. We investigated the genetic structure of three species of Collembola that are widely distributed and common in the northern hemisphere across thirteen European countries, including potential refuge areas south of the Alps. Ceratophysella denticulata (Bagnall, 1941) (Hypogasturidae) and Folsomia quadrioculata (Tullberg, 1871) (Isotomidae) are large species (up to 2.5 mm) that live in the uppermost soil layer (hemiedaphic) and reproduce sexually. In contrast, the third species, Isotomiella minor (Schaeffer, 1896) (Isotomidae), is small (up to 1.3 mm), lives deeper in soil (euedaphic) and reproduces via parthenogenesis. The investigated species are generally described as palearctic or holarctic in their distribution ranges, only I. minor is considered cosmopolitan, but its taxonomic status needs further attention [34–36]. In fact, all three species are more widely distributed but the part of introduced and indigenous populations in other regions of the world remains uncertain . Ceratophysella denticulata feeds on bacteria and fungi but presumably also on other soil animals such as nematodes [21, 38]. While C. denticulata functions as secondary decomposer, F. quadrioculata is a typical primary decomposer feeding predominantly on litter . In contrast, I. minor presumably mainly feeds on bacteria [39, 40]. Due to its parthenogenetic mode of reproduction, I. minor likely recovers more quickly from disturbances and colonizes new habitats faster than sexual species . Species of Hypogasturidae and Isotomidae have been found in Baltic amber (55–35 mya) [3–5, 41, 42] suggesting that they have been present in Europe for millions of years.
To explore divergences of lineages at a wide timespan, we analyzed one mitochondrial and two nuclear markers with different substitution rates. The mitochondrial COI gene is widely used for barcoding but previous studies demonstrated high intraspecific variance in Collembola [43, 44]. In order to resolve divergences deeper in time we added two nuclear markers that evolve at a lower rate than COI. The D3-D5 loop of 28S rDNA is a useful species marker for soil-living animals  but also exhibits intra-individual variance in Collembola [46, 47]. Further, we included the nuclear protein-coding gene Histone 3 with an expected mutation rate intermediate between COI and the 28S rDNA fragment (D3-D5 loop). For estimating divergence times of genetic lineages in absence of fossil records of the investigated species, we applied the general mutation rate for arthropods of COI, ranging from 1.5–2.3 % sequence divergence per million years [48, 49], the mutation rates of the two nuclear genes are unknown and were estimated by BEAST  in a combined alignment.
Assuming that soil buffered Quaternary temperature extremes, Collembola presumably survived in local patches in central Europe and expanded from these patches. Accordingly, endemic haplotypes are present in central and northern Europe. Further, we expected different genetic lineages to coexist locally due to expansions from isolated central European refuge populations. Due to more extensive survival during Quaternary glaciations of Collembola in southern Europe, the local molecular variance and genetic distances within populations from southern Europe should be higher than those within central and northern European populations.
Sampling of animals and DNA extraction
Leaf litter including humus layers from about two square meters of deciduous and coniferous forests were collected in 19 locations in 13 countries in Europe, including northwest and middle-west Russia, and the Pleistocene refuge areas of the Balkan (Montenegro, Serbia, Croatia, Greece), Italy and Spain, covering a geographic sampling range of 4,500 km on the east–west and ~4,000 km on the north–south axis, extending other studies on European Collembola by including countries from the Mediterranean to Scandinavia and western Eurasia. Animals were extracted by heat , collected in 96 % EtOH and stored at−20 °C until further analyses. For species identification specimens were sorted under a dissecting microscope and determined by light microscopy following Hopkin . If possible, five individuals per species and sampling location were sequenced (Table 1). In preliminary analyses we sequenced 12–20 individuals per species from 8 sampling locations and genetic variance within locations was very low, with uniformly a single genetic lineage dominating sampling locations. After these preliminary analyses we decided to restrict the sample size per sampling location to five individuals because for this study the dominant haplotypes are important to analyse genetic structure. Analysing more individuals per sampling location would include rare haplotypes and result in uneven sample sizes due to differences in natural abundance, which varied strongly from very few to hundreds of individuals per sampling location for all species.
Genomic DNA was extracted from single individuals of C. denticulata (n = 54), F. quadrioculata (n = 56) and I. minor (n = 55) using the DNeasy® Blood and Tissue Kit (Qiagen, Hilden, Germany) following the manufacturer’s protocol for animal tissue. Purified DNA was eluted in 30 μl buffer AE and stored at−20 °C until further preparation. PCRs of the nuclear Histone 3 and 28S rDNA (D3-D5 region) and the mitochondrial COI gene were performed in 25 μl volumes containing 12.5 μl SuperHot Taq Mastermix (Genaxxon Bioscience GmbH, Ulm, Germany) with 1.5 μl of each primer (10 pM), 4.5 μl H2O, 2 μl MgCl2 (25 mM) and 3 μl template DNA. Primers and PCR programs are given in Additional file 2: Table S1. Positive PCR products were purified with the QIAquick PCR Purification Kit (Qiagen, Hilden, Germany) following the manufacturer’s protocol and sent for sequencing to the Göttingen Genome Laboratory (Institute for Microbiology and Genetics, Georg August University of Göttingen). All sequences are available at NCBI GenBank (KF684371-KF684865, Additional file 2: Table S2). DNA was extracted from entire specimens but secondary vouchers (same morphological species from the same population) were deposited at our collections at J.F. Blumenbach Institute of Zoology and Anthropology, Georg August University Göttingen, Germany. For molecular clock analyses, COI and 28S rDNA sequences of outgroup taxa from the same family, congeneric species, and the same species from additional, mainly non-European, locations were downloaded from NCBI and BOLD databanks (Additional file 2: Table S3).
Phylogenetic analyses, divergence time estimation and population structure
Sequences were edited and corrected with Sequencher 4.10 (Gene Codes Corporation, USA), coding nucleotide sequences (H3 and COI) were translated into protein sequences using the standard and invertebrate mitochondrial codes implemented in Sequencher. For each species, nucleotide and protein sequences, were aligned separately and combined (as concatenated nucleotide matrix of three genes) with Clustal W  implemented in BioEdit 7.0.1 .
The best fit model of sequence evolution for each alignment (COI, 28S, H3, combined matrix) was inferred with to the hLRT in TOPALi v2.5  using the PHYML algorithm. Phylogenetic trees were calculated with RAxML v7.0.3  and MrBayes v3.1.2 . For Maximum likelihood analyses the model of sequence evolution was GTR + I + G (all four alignments) and 10,000 bootstrap replicates were calculated. For Bayesian Inference lset parameters were nst = 6, rates = invgamma (all four alignments), the MCMC (Markov Chain Monte Carlo) chains were run for ten million generations that were sampled every 1,000th generation. For the 10,000 sampled generations a burnin of 2,500 was used, eliminating the first 25 % of the remaining generations. In the absence of fossil or biogeographic calibration points for the investigated species in Europe, a strict molecular clock for the COI nucleotide alignment was applied in BEAST v1.7.4 . For constructing trees we used the Yule Process  as preliminary analyses indicated quicker convergence and higher probabilities and likelihoods than coalescent tree priors. However, topologies with different tree priors did not vary. The Yule Process also is more appropriate for the genetically highly diverged lineages as the substitution rate among branches is more variable than with coalescent priors. We used the widely adopted mutation rate for COI in arthropods of 2.3 % pairwise sequence divergence per million years, corresponding to a rate of 0.0115 [48, 49].
Convergence of the MCMC chain after 600 million generations (sampled every 60,000th generation) with a burnin of 25 % was confirmed using Tracer v1.4 . Divergence estimates were calculated with three datasets: (1) all COI and 28S sequences of this study and NCBI combined, with strict clock settings for COI (fixed rate of 0.0115) and estimated rates for 28S, (2) all COI sequences obtained in this study with a strict clock (fixed rate of 0.0115), and (3) all COI sequences of this study extended with sequences from non-European countries obtained from NCBI and BOLD databanks with a strict clock (fixed rate of 0.0115). The extended datasets (2) and (3) included additional sequences from Antarctica, Australia, Canada, Chile, New Zealand, South Africa and northern France and a more detailed outgroup sampling for better estimation of the substitution rate. Outgroups covered additional species within the respective genera and families, i.e. Hypogastrura (five species), Xenylla grisea and Gomphiocephalus hodgsoni for Hypogastruridae and Folsomia (three species), Parisotoma notabilis, Anurophorus septentrionalis and Isotoma (three species) for Isotomidae (Additional file 2: Table S3). For F. quadrioculata and I. minor, non-European 28S sequences of the D3-D5 region were not available and the combined datasets were extended only with the above mentioned outgroup taxa. The outgroup settings of the two isotomid species were identical, except for five additional sequences of the parthenogenetic species Parisotoma notabilis, to account for variance in the substitution rate due to the reproductive mode of I. minor.
Median-joining haplotype networks for the nucleotide datasets of COI, H3, 28S and the concatenated dataset were generated for all three species using the program Network 4.6 (Fluxus Technology, Suffolk, Great Britain). Molecular variance (AMOVA) within and between populations and isolation by distance (Mantel test) of all three genes (uncorrected p-distances) were analyzed separately in ARLEQUIN  with 20,000 permutations. To infer the molecular divergence times of C. denticulata, F. quadrioculata and I. minor, we generated a phylogenetic tree of 10 families with 23 genera and eight outgroup taxa including the taxa studied by D’Haese , using 28 COI sequences available at NCBI and a strict molecular clock in BEAST as described above.
Intraspecific variance of the three genes was similar in the three species. For C. denticulata the nucleotide alignments of COI, H3 and 28S had 240, 75 and 45 variable positions, respectively. Respective values for F. quadrioculata were 234, 65 and 2 and for I. minor 253, 50 and 15.
Isolation by distance was rejected as not being significant for all datasets and all species; however, genetic variance among populations was high, explaining 96–97 % (C. denticulata), 87–88 % (F. quadrioculata) and 92–99 % (I. minor) of the genetic variance of COI, H3 and the concatenated dataset (COI and H3; Table 2). Genetic distances between populations were very high (Table 3), ranging for COI between 14 and 20 % in C. denticulata and I. minor and 11–17 % in F. quadrioculata. Within population distances were zero or very low (<2 %). The only exceptions were one population of C. denticulata (Croatia: H3, 5.5 %), two populations of F. quadrioculata (Croatia: COI, 8 % and Greece: COI, 9.5 %; H3, 4.3 %) and two populations of I. minor (Russia: COI, 5.8 % and Germany: COI, 4.9 %) in which distinct COI and H3 lineages co-occurred. Genetic distances between and within regions (COI and H3) were highest in southern Europe and decreased towards the north, but were lowest within central Europe (Table 4). Between several populations in central and northern Europe genetic distances were lower than average (Table 5). Distances were low for H3 but high for COI between populations from Norway, Germany and France (F. quadrioculata; 14 %), from Montenegro and Italy (C. denticulata, 16 %), and from Norway and Spain (I. minor; 14 %).
The nuclear 28S rDNA (D3-D5 region) had the lowest genetic variance with little (6.2 % and 1.8 % in C. denticulata and I. minor, respectively) or almost no variation between populations (0.2 % for F. quadrioculata). Accordingly, the network analysis separated only haplotypes of C. denticulata into distinct geographic clusters (Additional file 1: Figure S1).
Phylogenetic trees based on the combined matrix gave the best resolution and statistical support for internal and terminal nodes (Figs. 1, 2 and 3). Topologies of ML and BI trees were similar with nearly all sampling locations clustering in separate clades with high support (pp: 0.97-1; bp: 95–100). Populations from southern European countries comprised a number of isolated COI and H3 lineages which derived early in each of the trees. Central and northern European populations were more homogeneous with shorter and more derived branches. Common patterns in each of the three species were the presence of one or two closely related phylogeographic clades with a wide distribution range across central and north-east Europe (blue and green lineages). Sampling locations in the Mediterranean were not covered by phylogeographic clades but rather by genetically distinct and phylogenetically isolated lineages (orange lineages). Interestingly, both isotomid species had phylogenetically related lineages with disjunct distributions. In F. quadrioculata, individuals from Spain and Croatia (black lineage) were monophyletic in the single gene phylogenetic trees and carried identical 28S sequences. In I. minor, individuals from Norway and Spain had distinct mitochondrial haplotypes, but nearly identical H3 and 28S alleles. Coexistence of genetic lineages did not occur in each of the three species except for one population from Croatia of F. quadrioculata. Different to the other species, C. denticulata had a central-eastern phylogeographic clade (black lineage) but the geographic range and spatial coherence of this lineage remains open due to the limited number of sampling locations.
Haplotype diversity of the mitochondrial gene was generally high in the three investigated species (Table 6) but most haplotypes differed in only a few positions. Haplotype networks of the COI and H3 gene were complex in each of the species as the majority of populations were separated by >10 mutation steps (Additional file 1: Figures S1-S3). For simplifying the dataset, we grouped the COI haplotypes and H3 and 28S alleles into lineages according to clades with very high [posterior probabilities (pp): 1, bootstrap (bs): 100] or high node support (pp: 0.97-99, bs: 89–100) in the phylogenetic trees based on single gene alignments (Additional file 1: Figures S4–S5).
Estimation of divergence times
Trees for molecular divergence estimates calculated with COI in BEAST differed only slightly from those generated with the combined matrix (Additional file 1: Figure S6). All populations and clades that comprised larger geographic areas were recovered, only positions of populations from Italy (C. denticulata), Greece (F. quadrioculata), Croatia (I.minor) and Montenegro (C. denticulata and I. minor) differed from the trees calculated with the combined matrices. Divergence estimates of major lineages fell into the Miocene (23.03-5.3 mya) in each of the three species, considerably predating Quaternary glaciations (1.8 my to present). Results were consistent among the three datasets, demonstrating robustness of the molecular clock analyses (Fig. 4).
Three major radiation events can be extrapolated from the datasets of the three species (Fig. 4, Additional file 1: Figure S6). First, large phylogeographic clades of lineages now present in central and north-east Europe separated from south European lineages during the Early and Middle Miocene (23.03-11.6 mya); second, clades in central and northern Europe separated during the Late Miocene (11.6-5.4 mya), and third, populations, i.e. sampling locations, radiated during the Pleistocene (0.011-1.8 mya). Divergence estimates of the COI datasets were consistent, divergences in dataset (2) without additional outgroups were even more conserved, i. e. lineages in general were younger. In the two isotomid species (F. quadrioculata and I. minor) divergence estimates of larger phylogeographic clades were younger in dataset (1) (combined COI and 28S) compared to datasets (2)–(3) (COI only), with predominantly Middle Miocene origin (16–11.5 mya). However, this probably is an artefact of the 28S gene which is strongly conserved among genera of Folsomia (H. von Saltzwedel, unpublished data) and probably other Isotomidae. In general, divergences of clades and populations were more recent in F. quadrioculata than in the other two species but still concentrated in the Miocene, independent of the dataset analyzed.
The phylogeographic dataset can only estimate the molecular age of European lineages of the respective species, which ranged between 22.5 and 25.8 my for C. denticulata, between 12.8 and 20.6 my for F. quadrioculata and between 17.2 and 24.9 my for I. minor, depending on the dataset analyzed. Estimated divergence times of the three Collembola species based on a phylogeny of 28 COI sequences can infer the ages of the respective genera, which were about twice as old with 48 ± 13 my for C. denticulata, 56 ± 14 my for F. quadrioculata and 49 ± 15 my for I. minor (Additional file 1: Figure S7). The age of the species must be between these two age estimates. As both analyses do not conflict with the fossil record, we take this as support for our molecular clock analyses and the Miocene origin of large phylogeographic clades in European Collembola.
This study investigated the genetic structure of Collembola at a pan-European scale. The results suggest three largely consistent phylogeographic patterns of the three investigated species. First, populations of Collembola in Europe are highly structured with high genetic variance between, but low variance within populations. Second, molecular divergence estimates and genetic distances suggest that common ancestors of present day populations in central and southern Europe probably persisted in isolated populations for millions of years. This is in agreement with our hypothesis that Quaternary climate changes were less severe for soil animals than for aboveground animals. Third, Miocene divergences dominate among geographic lineages indicating that diversification events correlate with wet and warm climate, and a biome change in central and eastern Europe from forest to grassland. The presence of the three species of Collembola in Europe during the Miocene is supported by their estimated Eocene origin that also correlates with the fossil record from Baltic amber [3–6].
In each of the three species, the majority of sampling locations were dominated by a single haplotype, suggesting that only few early colonizing individuals founded the populations which grew and expanded rapidly thereby preventing invasion of other lineages. This is conform to the founder-takes-it-all process  resulting in low genetic variance within but high variance between populations. Molecular divergence estimates of populations substantially predated Quaternary glaciations and originated in the Miocene. Long-term persistence of distinct genetic lineages without mixing of populations presumably is related to low mobility and high local abundances of springtails.
Strong molecular differentiation between populations indicates the existence of cryptic species. Indeed, recent molecular studies suggest that due to the widespread occurrence of cryptic species the number of Collembola species may be magnitudes higher than current estimates based on morphology [64, 65]. However, without evidence of differences in physiological, feeding, or mating preferences, we refrain from discussing potential cryptic species but rather consider them as species with high intraspecific genetic variance in COI. High genetic variation in COI have been reported from other arthropods in particular those living in soil, including Collembola [28, 43, 44, 65], oribatid mites [24, 66, 67] and Opilionida . In the sexual species C. denticulata genetic variance in nuclear genes was high but genetic distances within the slowly evolving 28S gene still were below the intra-lineage threshold previously estimated for this species .
Genetic distances between and within regions were high in the south compared to distances between and within the central and northern European regions. This suggests that Collembola follow the “southern richness and northern purity” scenario of recolonization of central and northern Europe . Lower distances within the central region suggest more recent expansion of populations or recent colonization of genetic lineages from east Europe or central Asia. However, molecular divergence time estimates of the three Collembola species indicate pronounced radiations into geographic lineages during the Miocene (20–5 mya). These divergences are considerably older than for most aboveground animals, suggesting that extinctions of populations during the Quaternary and recolonization of central Europe thereafter  had only minor effects on the population structure of present day Collembola populations.
The Miocene was warm and humid resembling subtropical climate with central Europe separating later into marine and continental climatic zones along west–east gradients with seasonality increasing to the east . For little sclerotized arthropods, susceptible to desiccation such as Collembola, the warm climate and high precipitation during the Miocene potentially facilitated large scale expansion. Accordingly, our divergence estimates indicate that southern and central European lineages separated during the Mid Miocene Climate Optimum (MMCO, 17–14 mya) when precipitation and temperature were higher than in the Early Miocene . Increasing seasonality and decreasing precipitation during the Miocene also changed the vegetation across Europe favoring the spread of deciduous forests, open woodlands and grassland vegetation in eastern and southern Europe, which expanded considerably during the Late Miocene (Messinian, 7.2-5.3 mya) [71, 73]. Interestingly, clades of C. denticulata and F. quadrioculata that cover the geographic regions from south-east, central-east and central Europe (Croatia, Russia, Serbia Austria, Germany, Poland) were dated to be of Late Miocene origin. The lower than average genetic distance in these clades suggests that a single or a few lineages expanded into these areas, possibly due to fragmentation of forest habitats into open woodlands, causing local extinctions or adaptations. In contrast, lineages of the parthenogenetic species I. minor were less separated into geographic regions, but rather split continuously from a ubiquitous and widespread ancestral lineage. This is also reflected by the equal distances of haplotypes in the network of the nuclear marker (28S), indicating that no region is more ancestral to others and suggesting separation of these lineages from a common ancestor. Accordingly, the two closely related (28S and H3) but geographically separated populations of Norway and eastern Spain likely represent relict populations that separated from an ancestral lineage about 12 mya and survived in isolated patches. Related lineages either went extinct or were not sampled. Further, divergences into geographic lineages mainly occurred during the Early and Middle Miocene and Pleistocene in I. minor, with only one divergence within the central European clade during the Late Miocene. In contrast, divergences of central European lineages of C. denticulata and F. quadrioculata predominantly took place several million years later in the Late Miocene and Early Pliocene.
Collembola are fast colonizers of newly formed and young habitats such as glacier forelands , but the geographic origin of these founder populations are unknown. The distinctness of populations, including the most closely related populations from the east of Europe suggests a regional recruitment of colonizers from local source populations that survived in glacial refugia and consequently diverged from more distant populations for millions of years. A more detailed sampling likely would reveal even more divergent haplotypes providing a more coherent picture of distribution ranges of genetic lineages. Further, a more detailed sampling in north-east and central-east Europe may uncover if north European Collembola descended from geographically distant source populations in Eurasia or from genetically distinct source populations from within Europe.
Southern European populations were dominated by isolated lineages, supporting the existence of south European refuge areas but also suggesting that the Alps function as major dispersal barrier. Accordingly, regional extinctions in central Europe and rapid recolonization of empty habitats from the surrounding, rather than from southern refugia, likely resulted in the strong phylogeographic structure of populations north of the Alps. Strong genetic differentiation, little or no mixing between populations and ancient diversifications resemble the pattern reported for epigeic springtail species of the genus Lepidocyrtus in the Mediterranean basin  and the soil-living oribatid mite species S. magnus . This indicates that evolutionary changes in soil are much slower and follow different trajectories as in aboveground invertebrates. However, coexistence of genetically distinct lineages in populations of the oribatid mite species S. magnus was more common than in the studied Collembola. This suggests that different factors regulate the establishment of populations among belowground invertebrates and that competition and early habitat colonization is more important for springtails than for oribatid mites [32, 74, 75].
Similar to distribution patterns of European springtail species , we show that genetic structure within Collembola species also deviate from common aboveground patterns. Overall, the results suggest that populations of Collembola across Europe are highly structured and that diversifications into regional lineages predominantly occurred in the Miocene with these lineages still persisting. Divergence dates and geographic structure suggest that soil microarthropods are less affected by global and long-term climatic changes than aboveground animals and plants. Buffering of adverse climatic conditions, the presence of organic matter colonized by bacteria and fungi in soil, large population sizes and the ability to tolerate freezing conditions likely contributed to the lower responsiveness of soil microarthropods. This suggests that in soil evolutionary processes are slowed down compared to the aboveground system, rendering soil-living detritivore microarthropods “living fossils” with communities in central Europe reflecting Miocene rather than Quaternary diversification events.
Hirst S, Maulik S. On some arthropod remains from the Rhynie Chert (Old Red Sandstone). Geol Mag. 1926;63:69–71.
Whalley P, Jarzembowski E. A new assessment of Rhyniella, the earliest known insect, from the Devonian of Rhyne. Scotland Nat. 1981;291:317.
Rapoport E. The geographical distribution of Neotropical and Antarctic Collembola. Pacific Insects Monogr 1971:99–118.
Hädicke C, Haug C, Haug J. Adding to the few: a tomocerid collembolan from Baltic amber. Paleodiversity. 2013;6:149–56.
Zawischa D. Einschlüsse im Bernstein. Arbeitskr Paläontologie Hann. 1993;21:1–32.
Weitschat W, Wichard W. Diplurans and Springtails - Insecta: Diplura and Collembola. In: Atlas of Plants and Animals in Baltic Amber. München: Dr. Friedrich Pfeil; 2002. p. 86–7.
Fountain MT, Hopkin SP. Biodiversity of Collembola in urban soils and the use of Folsomia candida to assess soil “quality”. Ecotoxicology. 2004;13:555–72.
Petersen H, Luxton M. A comparative analysis of soil fauna populations and their role in decomposition processes. Oikos. 1982;39:288.
Seastedt TR. The role of microarthropods in decomposition and mineralization processes. Annu Rev Entomol. 1984;29:25–46.
Wolters V. Soil invertebrates - effects on nutrient turnover and soil structure - a review. Zeitschrift für Pflanzenernährung und Bodenkd. 1991;154:389–402.
Wardle DA, Bardgett RD, Klironomos JN, Setälä H, van der Putten WH, Wall DH. Ecological linkages between aboveground and belowground biota. Science. 2004;304:1629–33.
Healey I. An ecological study of temperatures in a Welsh moorland soil. J Anim Ecol. 1967;36:425–34.
Gass F, Gillet S, Ponge J. The use of directional traps for the assessment of short-term phenanthrene effects upon soil springtail communities. Environ Pollut. 2006;2:364–70.
Ohlsson L, Verhoef HA. Effects of diet composition on cold adaptation in temperate Collembola. Comp Biochem Physiol. 1988;91:475–9.
Worland MR. Factors that influence freezing in the sub-Antarctic springtail Tullbergia antarctica. J Insect Physiol. 2005;51:881–94.
Bahrndorff S, Holmstrup M, Petersen H, Loeschcke V. Geographic variation for climatic stress resistance traits in the springtail Orchesella cincta. J Insect Physiol. 2006;52:951–9.
Konestabo HS, Michelsen A, Holmstrup M. Responses of springtail and mite populations to prolonged periods of soil freeze-thaw cycles in a sub-arctic ecosystem. Appl Soil Ecol. 2007;36:136–46.
Hopkin SP. Biology of the Springtails. Oxford: Oxford University Press; 1997.
Block W. Cold hardiness in invertebrate poikilotherms. Comp Biochem Physiol Physiol. 1982;73:581–93.
Zettel J, Zettel U, Egger B. Jumping technique and climbing behaviour of the collembolan Ceratophysella sigillata (Collembola: Hypogasturidae). Eur J Entomol. 2000;97:41–5.
Chahartaghi M, Langel R, Scheu S, Ruess L. Feeding guilds in Collembola based on nitrogen stable isotope ratios. Soil Biol Biochem. 2005;37:1718–25.
Maraun M, Martens H, Migge S, Theenhaus A, Scheu S. Adding to “the enigma of soil animal diversity”: fungal feeders and saprophagous soil invertebrates prefer similar food substrates. Eur J Soil Biol. 2003;39:85–95.
Garrick RC, Sands CJ, Rowell DM, Hillis DM, Sunnucks P. Catchments catch all: long-term population history of a giant springtail from the southeast Australian highlands--a multigene approach. Mol Ecol. 2007;16:1865–82.
Rosenberger M, Maraun M, Scheu S, Schaefer I. Pre- and post-glacial diversifications shape genetic complexity of soil-living microarthropod species. Pedobiologia. 2013;56:79–87.
Frati F, Dell’Ampio E, Casasanta S, Carapelli A, Fanciulli P. Large amounts of genetic divergence among Italian species of the genus Orchesella (Insecta, collembola) and the relationships of two new species. Mol Phylogenet Evol. 2000;17:456–61.
van der Wurff AWG, Gols R, Ernsting G, van Straalen NM. Population genetic structure of Orchesella cincta (Collembola; Hexapoda) in NW Europe, as revealed by microsatellite markers. Pedobiologia. 2005;49:167–74.
van der Wurff AWG, Isaaks JA, Ernsting G, Van Straalen NM. Population substructures in the soil invertebrate Orchesella cincta, as revealed by microsatellite and TE-AFLP markers. Mol Ecol. 2003;12:1349–59.
Cicconardi F, Nardi F, Emerson BC, Frati F, Fanciulli P. Deep phylogeographic divisions and long-term persistence of forest invertebrates (Hexapoda: Collembola) in the north-western Mediterranean basin. Mol Ecol. 2010;19:386–400.
Webb NR, Block W. Aspects of cold hardiness in Steganacarus magnus (Acari: Cryptostigmata). Exp Appl Acarol. 1993;17:741–8.
Salomone N, Emerson BC, Hewitt GM, Bernini F. Phylogenetic relationships among the Canary Island Steganacaridae (Acari, Oribatida) inferred from mitochondrial DNA sequence data. Mol Ecol. 2002;11:79–89.
Siepel H. Life-history tactics of soil microarthropods. Biol Fert Soils. 1994;263–278.
Lindberg N, Bengtsson J. Population responses of oribatid mites and collembolans after drought. Appl Soil Ecol. 2005;28:163–74.
Petersen H. General aspects of collembolan ecology at the turn of the millennium. Pedobiologia. 2002;46:246–60.
Gao Y, Potapov M. Isotomiella (Isotomidae: Collembola) of China. Ann la Soc Entomol Fr. 2011;47:350–6.
Bedos A, Deharveng L. The Isotomiella of Thailand (Collembola: Isotomidae), with description of five new species. Entomol Scand. 1994;25:451–60.
Deharveng L, Olieveira E. Isotomiella (Collembola: Isotomidae) d´Amazonie: les espèces du groupe delamarei. Ann la Société Entomol Fr. 1990;26:185–201.
Greenslade P, Convey P. Exotic Collembola on subantarctic islands: Pathways, origins and biology. Biol Invasions. 2012;14:405–17.
Heidemann K, Hennies A, Schakowske J, Blumenberg L, Ruess L, Scheu S, Maraun M. Free-living nematodes as prey for higher trophic levels of forest soil food webs. Oikos. 2014;123:1199–211.
Ponge JF. Food resources and diets of soil animals in a small area of Scots pine litter. Geoderma. 1991;49:33–62.
Ferlian O, Klarner B, Langeneckert AE, Scheu S. Trophic niche differentiation and utilisation of food resources in collembolans based on complementary analyses of fatty acids and stable isotopes. Soil Biol Biochem. 2015;82:28–35.
Lawrence P. Ten species of Collembola from Baltic amber. Pr Muz Ziemi. 1985;37:101–4.
Christiansen K, Pike E. Cretaceous Collembola (Arthropoda, Hexapoda) from the Upper Cretaceous of Canada. Cretac Res. 2002;23:165–88.
Hogg ID, Hebert PDN. Biological identification of springtails (Hexapoda: Collembola) from the Canadian Arctic, using mitochondrial DNA barcodes. Can J Zool. 2004;82:749–54.
Torricelli G, Carapelli A, Convey P, Nardi F. High divergence across the whole mitochondrial genome in the “pan-Antarctic” springtail Friesea grisea: evidence for cryptic species? Gene. 2010;449:30–40.
Maraun M, Heethoff M, Schneider K, Scheu S, Weigmann G, Cianciolo J, Thomas RH, Norton RA. Molecular phylogeny of oribatid mites (Oribatida, Acari): evidence for multiple radiations of parthenogenetic lineages. Exp Appl Acar. 2004;33:183–201.
Porco D, Bedos A, Greenslade P. Challenging species delimitation in Collembola: cryptic diversity among common springtails unveiled by DNA barcoding. Invertebr Syst. 2012;26:470–7.
Porco D, Potapov M, Bedos A, Busmachiu G, Weiner WM, Hamra-Kroua S, Deharveng L. Cryptic diversity in the ubiquist species Parisotoma notabilis (Collembola, Isotomidae): a long-used chimeric species? PLoS One. 2012;2:e46056.
Brower AV. Rapid morphological radiation and convergence among races of the butterfly Heliconius erato inferred from patterns of mitochondrial DNA evolution. Proc Natl Acad Sci U S A. 1994;91:6491–5.
Avise JC. Molecular Markers, Natural History and Evolution. New York: Chapman & Hall; 1994.
Drummond AJ, Rambaut A. BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evol Biol. 2007;7:214.
Kempson D, Lloyd M, Ghellardi R. A new extractor for woodland litter. Pedobiologia. 1963;3:1–21.
Hopkin SP. A Key to the Collembola (Springtails) of Britain and Ireland, Field Studies Council. 2007.
Thompson JD, Higgins DG, Gibson TJ. CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice. Nucleic Acids Res. 1994;22:4673–80.
Hall TA. BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT. Nucleic Acids Symposium Series. 1999:95–98.
Milne I, Lindner D, Bayer M, Husmeier D, McGuire G, Marshall DF, Wright F. TOPALi v2: a rich graphical interface for evolutionary analyses of multiple alignments on HPC clusters and multi-core desktops. Bioinformatics. 2009;25:126–7.
Stamatakis A, Ott M, Ludwig T. RAxML-OMP: an efficient program for phylogenetic inference on SMPs. In: Malyshkin V, editor. Parallel Computing Technologies (Lecture Notes in Computer Science, vol. 3606). Berlin, Heidelberg: Springer; 2005. p. 288–302
Ronquist F, Huelsenbeck JP. MrBayes 3: Bayesian phylogenetic inference under mixed models. Bioinformatics. 2003;19:1572–4.
Drummond AJ, Suchard MA, Xie D, Rambaut A. Bayesian phylogenetics with BEAUti and the BEAST 1.7. Mol Biol Evol. 2012;29:3–6.
Gernhard T, Hartmann K, Steel M. Stochastic properties of generalised Yule models, with biodiversity applications. J Math Biol. 2008;57:713–35.
Rambaut A, Drummond AJ. Tracer v1.4. 2007; http://tree.bio.ed.ac.uk/software/tracer/.
Excoffier L, Lischer HEL. Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Mol Ecol Resour. 2010;10:564–7.
D’Haese CA. Were the first springtails semi-aquatic? A phylogenetic approach by means of 28S rDNA and optimization alignment. Proc R Soc B Biol Sci. 2002;269:1143–51.
Waters JM, Fraser CI, Hewitt GM. Founder takes all: density-dependent processes structure biodiversity. Trends Ecol Evol. 2013;28:78–85.
Emerson BC, Cicconardi F, Fanciulli P, Shaw PJA. Phylogeny, phylogeography, phylobetadiversity and the molecular analysis of biological communities. Philos Trans R Soc London B. 2011;366:2391–402.
Cicconardi F, Fanciulli P, Emerson BC. Collembola, the biological species concept and the underestimation of global species richness. Mol Ecol. 2013;22:5382–96.
Heethoff M, Domes K, Laumann M, Maraun M, Norton RA, Scheu S. High genetic divergences indicate ancient separation of parthenogenetic lineages of the oribatid mite Platynothrus peltifer (Acari, Oribatida). J Evol Biol. 2007;20:392–402.
Schäffer S, Koblmüller S, Pfingstl T, Sturmbauer C, Krisper G. Contrasting mitochondrial DNA diversity estimates in Austrian Scutovertex minutus and S. sculptus (Acari, Oribatida, Brachypylina, Scutoverticidae). Pedobiologia. 2010;53:203–11.
Boyer SL, Baker JM, Giribet G. Deep genetic divergences in Aoraki denticulata (Arachnida, Opiliones, Cyphophthalmi): a widespread “mite harvestman” defies DNA taxonomy. Mol Ecol. 2007;16:4999–5016.
Hewitt GM, Ibrahim KM. Inferring glacial refugia and historical migrations with molecular phylogenies. In Integrating Ecology and Evolution in a Spatial Context. Edited by Silvertown J, Antonovics J. Oxford: Blackwell Science Ltd.;2001:271–294.
Hewitt G. Post-glacial re-colonization of European biota. Biol J Linn Soc. 1999;68:87–112.
Bruch AA, Uhl D, Mosbrugger V. Miocene climate in Europe — patterns and evolution. Palaeogeogr Palaeoclimatol Palaeoecol. 2007;253:1–7.
Zachos JC, Dickens GR, Zeebe RE. An early Cenozoic perspective on greenhouse warming and carbon-cycle dynamics. Nature. 2008;451:279–83.
Mosbrugger V, Utescher T, Dilcher DL. Cenozoic continental climatic evolution of Central Europe. Proc Natl Acad Sci U S A. 2005;102:14964–9.
Hågvar S. Primary succession of springtails (Collembola) in a norwegian glacier foreland. Arctic, Antarct Alp Res. 2010;42:422–9.
Ingimarsdóttir M, Caruso T, Ripa J, Magnúsdóttir OB, Migliorini M, Hedlund K. Primary assembly of soil communities: disentangling the effect of dispersal and local environment. Oecologia. 2012;170:745–54.
Ulrich W, Fiera C. Environmental correlates of species richness of European springtails (Hexapoda: Collembola). Acta Oecologia. 2009;35:45–52.
We thank Patrick Pachl for help with sampling of the Balkan region. We are grateful to Jens Bast, Julian Bock, Elena Corral Hernández, Guido Humpert, Sven Marhan, Jordi Moya-Laraño and Lucija Šerić Jelaska for sending litter and soil material.
This study was funded by the DFG (Deutsche Forschungsgemeinschaft).
Availability of data and materials
DNA sequences: Genbank accessions KF684371-KF684865.
HvS, SS and IS conceived of the project and took the lead in design, discussion, synthesis, coordination and writing. All authors brought distinctive expertise to the collaboration and contributed importantly to the ideas represented, as well as through the drafting and revising of the article. All authors read and approved the final manuscript.
The authors declare that they have no competing interests.
Ethics approval and consent to participate
Sampling sites were outside Nature Reserve Areas and no permissions for soil samples were required. The field study did not involve any endangered or protected species.
Map of all sampling points and approcimate sketches of extent of permanent ice shields (white filled areas) and borders of polar desert (red, dotted line) and permafrost (black, dashed line) climate during the Last Glacial Maximum (~20,000 years ago) that strongly shaped genetic and species richness above the ground in central and northern Europe (after Hewitt ). Figure S2. Haplotype Networks of Ceratophysella denticulata. Figure S3. Haplotype Networks of Folsomia quadrioculata. Figure S4. Haplotype Networks of Isotomiella minor. Figure S5. Bayesian phylogeny based on COI. Figure S6. Bayesian phylogeny based on H3. Figure S7. Molecular divergence estimates of dataset (2). Figure S8. Molecular phylogeny and divergence estimates of 28 species of Collembola based on COI calculated with BEAST. Figure S9. Molecular divergence estimates in C. denticulata of a) dataset (1) and b) dataset (3). Outgroups: C. = Ceratophysella, H. = Hypogastrura, X = Xenylla, G = Gomphiocephalus. Figure S10. Molecular divergence estimates in F. quadrioculata of a) dataset (1) and b) dataset (3). Outgroups: F. = Folsomia, I. = Isotoma, C. = Cryptopygus, A. = Anurophorus. Figure S11. Molecular divergence estimates in I. minor of a) dataset (1) and b) dataset (3). Outgroups: F. = Folsomia, I. = Isotoma, C. = Cryptopygus, A. = Anurophorus. (PDF 920 kb)
a) Primer sequence of and b) PCR conditions for the ribosomal D3-D5 region of 28S rDNA, and the protein coding Histone 3 (H3) and COI genes. Table S2. DNA sequences available at GenBank. Table S3. Summary of additional locations and outgroup taxa included in the extended datasets (1) and (3) of the molecular clock analyses. Accession numbers from NCBI and BOLD databanks are listed. For some isotomiid species, the complete D3 fragment of 28S was not available, missing positions were filled with “?” in the alignments. (XLSX 17 kb)