Surprisingly high genetic divergence of the mitochondrial DNA barcode fragment (COI) within Central European woodlice species (Crustacea, Isopoda, Oniscidea).

DNA barcoding has become the most popular approach for species identification in recent years. As part of the German Barcode of Life project, the first DNA barcode library for terrestrial and freshwater isopods from Germany is presented. The analyzed barcode library included 38 terrestrial (78% of the documented species of Germany) and five freshwater (63%) species. A total of 513 new barcodes was generated and 518 DNA barcodes were analyzed. This analysis revealed surprisingly high intraspecific genetic distances for numerous species, with a maximum of 29.4% for Platyarthrushoffmannseggii Brandt, 1833. The number of BINs per species ranged from one (32 species, 68%) to a maximum of six for Trachelipusrathkii (Brandt, 1833). In spite of such high intraspecific variability, interspecific distances with values between 12.6% and 29.8% allowed a valid species assignment of all analyzed isopods. The observed high intraspecific distances presumably result from phylogeographic events, Wolbachia infections, atypical mitochondrial DNAs, heteroplasmy, or various combinations of these factors. Our study represents the first step in generating an extensive reference library of DNA barcodes for terrestrial and freshwater isopods for future molecular biodiversity assessment studies. Michael J. Raupach, Björn Rulik, Jörg Spelda.

Introduction

Isopods are a highly diverse group of invertebrates, with more than 10,300 species described to date (Boyko et al. 2008; Poore 2012). Most of these peracarid crustaceans are free-living, but a number of marine species represent bizarre ectoparasites that infest crustacean and fish species (e.g., Raupach and Thatje 2006; Williams and Boyko 2012; Hadfield et al. 2014; Smit et al. 2014). Isopods range in body length from 0.5 mm (members of the family ) up to 500 mm (species of the famous giant deep-sea isopod genus Milne-Edwards, 1879) (McClain et al. 2015). With more than 4,500 known marine species to date, isopods can be found from all shorelines of the world down to the abyssal depths of the oceans where asellote isopods have undergone a massive radiation and represent the dominant taxon (e.g., Wilson and Hessler 1987; Wägele 1989; Raupach et al. 2004, 2009). Approximately 900 isopod species colonized freshwater habitats including lakes, rivers, underground waters, and even thermal springs (e.g., Verovnik et al. 2005; Wilson 2008).

Isopods are, however, not restricted to the aquatic realms only. One group, the or woodlice, are the most successful group of crustaceans that invaded the land by far. Without doubt, these animals represent the most familiar and well-known group of isopods to humans. In contrast to other amphibious crustaceans, e.g., land crabs of the family or terrestrial hermit crabs of the genus Latreille, 1829, no developmental stage (egg, juvenile, etc.) of the requires free water and all biological activities are conducted on land (e.g., Broly et al. 2013). The have evolved a number of unique adaptations, such as the water conducting system, various forms of pleopodal lungs and the cotyledons in the marsupium (e.g., Sfenthourakis and Taiti 2015). Based on the dorsal surface of their exoskeleton, various other morphological traits as well as ecological strategies and behavior, woodlice can be roughly categorized in three main groups (Schmalfuss 1984; Hornung 2011): i) the runners, characterized with an elongate, slightly convex body and long pereopods (e.g., Latreille, 1804), ii) the clingers, with a flat broad body and short but strong pereopods (e.g., Brandt, 1833), and iii) the rollers, with a highly convex body able to roll up into a ball (pill bugs) (e.g., Brandt, 1833) (Fig. 1). Whereas their dispersion ability is rather limited, woodlice are found in almost all biomes of the world except the poles and high mountain ranges (Hornung 2011; Sfenthourakis and Taiti 2015). A hot spot of woodlice diversity is located in the Mediterranean region (Sfendourakis and Taiti 2015), and some species have been introduced to other parts of the world by humans in the past, e.g., to North America (Jass and Klausmeier 2000; Singer et al. 2012; Hornung et al. 2015) and other regions (e.g., Gruner 1966; Slabber and Chown 2002; Karasawa and Nakata 2018). Furthermore, oniscid isopods are amongst the most common and species-rich components of cave-dwelling animal groups with high numbers of troglobitic species (Sfenthourakis and Taiti 2015). In some ecosystems, e.g., European forests, woodlice perform a central role in the decomposition, being largely phytosaprophagous and often occur in very high population densities (e.g., Dias and Hassall 2005; Gongalsky et al. 2005; Hättenschwiler et al. 2005; David 2014; Špaldoňová and Frouz 2014), but also act as important prey for a broad range of predatory arthropods (Raupach 2015). Until now, more than 3,700 species of oniscid isopods have been described worldwide (Schmidt 2008; Sfenthourakis and Taiti 2015). For Germany, 49 species of terrestrial and eight species of freshwater isopods are reported so far (Grünwald 2016).

Figure 1.

Various woodlouse species of Germany A Linnaeus, 1758 B Budde-Lund, 1885 C (Brandt, 1833) D (Heller, 1858) E (Scopoli, 1763) F Strouhal, 1953 G (C. Koch, 1841) H Brandt, 1833. Scale bar: 1 mm. Photograph credits: A–G Jörg Spelda H Armin Rose.

Since its beginning almost 15 years ago, the concept of DNA barcoding for species identification has revolutionized biodiversity research (Valentini et al. 2009; Cristescu 2014). For many groups of animals, an approximate 650 base pair (bp) frag­ment of the mitochondrial cytochrome c oxidase subunit I (COI) gene was selected as marker of choice (Hebert et al. 2003a). The efficiency of DNA barcoding is based on a simple assumption: each species will most likely have similar DNA barcode sequences representing their intraspecific variability whereas the genetic variation between species exceeds the variation within species (Hebert et al. 2003a, 2003b). In this context, the German Barcode of Life initiative (GBoL; www.bolgermany.de) aims at capturing the genetic diversity of animals, fungi, and plants of Germany. Various comprehensive barcode libraries of arthropods, e.g., marine crustaceans (Raupach et al. 2015), spiders (Astrin et al. 2016), and myriapods (Spelda et al. 2011), have been generated in the past. In terms of isopods, most DNA barcoding studies focused on marine species so far (e.g., Khalaji-Pirbalouty and Raupach 2014, 2016; Raupach et al. 2015; ; Brix et al. 2018; Chew et al. 2018; Kakui et al. 2019), whereas for terrestrial and freshwater taxa almost no studies do exist (Asmyhr and Cooper 2012; Zimmermann et al. 2015, 2018a, 2018b). However, no comprehensive DNA barcode reference library has been published for these taxa until now.

In this study we present the first DNA barcode library of terrestrial and freshwater isopods with a focus on Central European species, with a specific emphasis on the . The analyzed barcode library includes 38 terrestrial (78% of the known species of Germany) and five freshwater (63%) species. In summary, 513 new barcodes were generated and a total number of 518 DNA barcodes was analyzed.

Materials and methods

Sampling of specimens

Samples used for this study were collected between 2000 and 2018 by pitfall traps, sieves, or by hand. Specimens were stored in ethanol (96%) and identified by two of the authors (MJR, JS) using a combination of keys provided in Schmölzer (1964), Gruner (1966), Hopkin (1991), and Berg and Wijnhoven (1997). In total, 513 new DNA barcodes of 46 species were generated. For our analysis we also included five DNA barcodes of the sea slater (Linnaeus, 1767) as part of a previous study (Raupach et al. 2010), generating a total data set of 518 DNA barcodes from 46 species. Five of the analyzed species, Dollfus, 1887 (n = 1, Spain), Brandt, 1833 (n = 2, France), Fabricius, 1798 (n = 2, Italy), (Budde-Lund, 1885) (n = 4, Mallorca, Spain), and Grebnitzky, 1874 (n = 1, Spain) are not recorded from Germany so far but were included for comparison. The number of analyzed specimens per species ranged from one (5 species) to a maximum of 57 for Latreille, 1804. Most isopods were collected in Germany (n = 458, 88.3%), whereas some specimens from other countries were included: Austria (3, 0.6%), Denmark (4, 0.8%), France (3, 0.6%), Italy (3, 0.6%), Luxembourg (38, 7.3%), Spain (6, 1.2%), and Switzerland (3, 0.6%).

DNA barcode amplification, sequencing, and data depository

Laboratory operations were carried out either at the Canadian Center for DNA Barcoding (CCDB), University of Guelph, following standardized protocols for COI amplification and sequencing (Ivanova et al. 2006; deWaard et al. 2008), the molecular lab rooms of the German Centre of Marine Biodiversity Research (DZMB), Senckenberg am Meer, in Wilhelmshaven, the working group Systematics and Evolutionary Biology at the Carl von Ossietzky University Oldenburg, or the Zoologisches Forschungsmuseum Alexander Koenig (ZFMK), Bonn, all located in Germany. Photographs were taken for each studied isopod before molecular work was performed. One or two legs of one body side were removed for the subsequent DNA extraction. For some very small isopods with a body length < 3 mm, e.g., specimens of Schöbl, 1860 or Leach, 1814, partial or complete specimens were used for DNA extraction. In the case of own molecular studies, DNA was extracted using the QIAmp Tissue Kit (Qiagen GmbH, Hilden, Germany) or NucleoSpin Tissue Kit (Macherey-Nagel, Düren, Germany), following the extraction protocol. Detailed information of used primers, PCR amplification and sequencing protocols are given in previous publications (see Raupach et al. 2015; Astrin et al. 2016). All purified PCR products were cycle-sequenced and sequenced in both directions at a contract sequencing facility (GATC, Konstanz, Germany), using the same primers as used in PCR. Double stranded sequences were assembled and checked for putative mitochondrial pseudogenes (numts) by analyzing the presence of stop codons, frameshifts as well as double peaks in chromatograms with the Geneious version 8.1.9 software package (Biomatters, Auckland, New Zealand) (Kearse et al. 2012). For sequence verification, BLAST searches (nBLAST, search set: others, program selection: megablast) were performed to confirm the identity of all new sequences as isopod sequences based on already published sequences (high identity values, very low E-values) (Zhang et al. 2000; Morgulis et al. 2008).

Comprehensive voucher information, taxonomic classifications, photos, DNA barcode sequences, used primer pairs and trace files including their quality are publicly accessible through the public data set “DS-BISCE” (Dataset ID: dx.doi.org/10.5883/DS-BISCE) on the Barcode of Life Data Systems (BOLD; www.boldsystems.org) (Ratnasingham and Hebert 2007). Parallel to this, all new barcode data were deposited in GenBank (accession numbers MN810569–MN810873, MT521085–MT521292).

DNA barcode analysis

Following a standardized approach of DNA barcode analysis, the BOLD workbench was used to calculate the nucleotide composition of the sequences and distributions of Kimura-2-parameter distances (K2P; Kimura 1980) within and between species (align sequences: BOLD aligner; ambiguous base/gap handling: pairwise deletion). All barcodes became subject of the Barcode Index Number (BIN) system as it is implemented in BOLD (2020–06–05). In doing so, DNA barcodes are clustered in order to produce operational taxonomic units that closely correspond to species (Ratnasingham and Hebert 2013). Using the given default settings, a recommended threshold of 2.2% was applied for a rough differentiation of intraspecific and interspecific K2P distances (Ratnasingham and Hebert 2013). It should be noted, however, that the BIN assignments on BOLD may change due to the addition of new sequences. Therefore, individual BINs can be split or merged in the light of new data (Ratnasingham and Hebert 2013).

A neighbor-joining cluster analysis (NJ; Saitou and Nei 1987) was performed for all studied species for a graphical representation of the genetic differences between sequences and clusters of sequences using MEGA 10.0.5 (Kumar et al. 2018). Again, the K2P model was chosen as the model for sequence evolution for comparison purposes with previous studies. For validation, non-parametric bootstrap support values were obtained by resampling and analyzing 1,000 replicates (Felsenstein 1985). All analyses were based on an alignment of all studied barcode sequences that was generated using MUSCLE (Edgar 2004) implemented in MEGA 10.0.5. It should be explicitly noted that this analysis is not intended to be phyloge­netic. Instead of this, the shown topology represents a graphical visualization of DNA barcode divergences/distances and putative species cluster.

Results

We analyzed 518 DNA barcode sequences of 46 isopod species. A list of species is presented in the supporting information (Suppl. material 1). Fragment lengths of the analyzed DNA barcodes ranged from 407 to 658 bp. As previously shown for arthropods, the DNA barcode region was characterized by a high AT-content: average sequence compositions were A = 24.6%, C = 18.1%, G = 21.5%, and T = 35.8%. Fourteen (30.4%) species had intraspecific distances > 2.2%, with a maximum of 29.4% for Brandt, 1833. Interspecific distances within the analyzed taxa had values between 12.6% ( Brandt, 1833; Stein, 1859) and 29.8% ( Valkanov, 1936; Budde-Lund, 1885). In total, 76 BINs were found. The number of BINs per species ranged from one (32 species, 68%) to a maximum of six ( (Brandt, 1833)). No BIN sharing between species was observed. The NJ analyses revealed non-overlapping clusters with bootstrap support values > 95% for 39 species (95%) with more than one studied specimen (Fig. 2). A more detailed topology of all analyzed specimens is presented in the supporting information (Suppl. material 2).

Figure 2.

Neighbor-joining (NJ) topology of the analyzed isopod species based on Kimura 2-parameter distances. Triangles show the relative number of individual’s sampled (height) and sequence divergence (width). Red triangles highlight terrestrial species with intraspecific maximum pairwise distances > 2.2%, whereas dark blue triangles indicate freshwater species with such distances. Numbers next to nodes represent non-parametric bootstrap values > 90% (1,000 replicates). Asterisks indicate species not recorded in Germany.

Continued.

Discussion

Our study revealed very high intraspecific distances for numerous woodlice species (Tab. 1), including abundant and well-known species as Latreille, 1804 (maximum intraspecific distances (ISD): 12.16) or (Brandt, 1833) (max ISD: 13.47). Intraspecific distance values higher than 2.2% were also shown for three of the five analyzed freshwater species (Tab. 1). The observed high variability can be caused by a number of different factors and will be discussed in the following.

Table 1.

Molecular distances based on the Kimura 2-parameter model of the analyzed specimens of the analyzed isopod species with intraspecific distances > 2.2% using the BOLD work bench. ISD = intraspecific distance. BINs are based on the barcode analysis from 05–06–2020. See methods for explanation of basis.

Species	n	BINs	Mean ISD	Max ISD
Armadillidiumvulgare (Latreille, 1804)	28	AAE6611, AAH4108, AAH4111, AAU1529	3.76	6.44
Asellusaquaticus (Linnaeus, 1758)	41	ACF1266, AEC4774, AAA1970	4.25	13.37
Oniscusasellus (Linnaeus, 1758)	33	ADM8743, ADM8116, ADK9123	2.12	5.63
Philosciaaffinis Verhoeff, 1908	3	ADM8125, AAY1058	3.63	5.44
Philosciamuscorum (Scopoli, 1763)	38	AAH4103, AAH4104	0.3	2.98
Platyarthrushoffmannseggii Brandt, 1833	33	AAV8050, AAV8051, ADK9658	9.4	29.35
Porcelliomontanus Budde-Lund, 1885	6	ADR0694, ADM7742	1.26	3.81
Porcellioscaber Latreille, 1804	57	AAC3755, AAZ0248, ABA5892, ADK8850, ADM8147	2.58	12.16
Porcelliospinicornis Say, 1818	6	ADF7011, ADI3596	3.01	5.13
Proaselluscavaticus (Leydig, 1871)	8	ADX3790, ADW6988, ADX4659	1.61	2.95
Proaselluscoxalis (Dollfus, 1892)	13	ACI1746, ACH7545	2.81	5.78
Trachelipusrathkii (Brandt, 1833)	16	AAH4102, ADK8699, ADK8533, ADM8087, ADM8088, ADF6188	6.89	16.59
Trichoniscoideshelveticus (Carl, 1908)	23	ADM7247, ADM7248, ADM7249	1.07	5.46
Trichoniscuspusillus Brandt, 1833	22	AAN7523, AAZ1993	6.8	13.47

First, phylogeographic events may generate different haplotypes and distinct mitochondrial lineages. In the case of European woodlice species, numerous studies showed complex phylogeographic patterns correlated with high variability of the studied mitochondrial markers including COI, e.g., for species of the genus Racovitza, 1908 (Bedek et al. 2019), the common sea slater (Linneaus, 1767) (Raupach et al. 2014), spp. (Klossa-Kilia et al. 2005), the common woodlouse Linnaeus, 1758 (Bilton et al. 1999), spp. (Poulakakis and Sfenthourakis 2008), Ebner, 1868 (Gentile et al. 2010), or two species of the genus Budde-Lund, 1908 (Parmakelis et al. 2008). Similar results have been also demonstrated for freshwater isopods of the genus Geoffroy, 1762 (Verovnik et al. 2004; Verovnik et al. 2005; Sworobowicz et al. 2015; Pérez-Moreno et al. 2017) and Dudich, 1925 (Ketmaier 2002; Eme et al. 2013; Kilikowska et al. 2013). Our data set revealed extremely high intraspecific distance values for the myrmecophilous isopod Brandt, 1833 (n = 33), with a maximum value of 29.4% (Tab. 1). It is a small, white, and blind oniscid isopod that is widely-distributed in Europe and strictly associated with various ant species (e.g., Mathes and Strouhal 1954; Gruner 1966; Parmentier et al. 2017). A few other species are found in the Mediterranean region, e.g., Budde-Lund, 1879 (Garci and Cruz 1986), which can be easily differentiated from . The NJ topology revealed three distinct lineages associated with three BINs (Fig. 3), but no clear correlation of the analyzed specimens to specific sampling regions. Furthermore, we found no ant-host-specific correlation of the observed lineages. For some other species distinct lineages were also detected, but no conspicuous substructures were revealed (see Suppl. material 3).

Figure 3.

Subtree of the Neighbor-joining topology based on Kimura 2-parameter distances of all analyzed specimens of Brandt, 1833 and nearest neighbor. Branches with specimen ID-number from BOLD and sample localities. Numbers next to internal nodes are non-parametric bootstrap values (in %) with values higher than 80. BIN values are based on the barcode analysis from 05-06-2020. The isopod drawing by Christian Schmidt was obtained from Raupach (2005).

Second, the presence of the inherited alpha-proteobacteria Hertig, 1936 can affect the mitochondrial variability in arthropods (e.g., Hurst and Jiggins 2005; Werren et al. 2008; Correa and Ballard 2016). These endosymbionts are transmitted vertically through host eggs and alter the biology of their host in various ways, including the induction of reproductive manipulations, such as feminization, parthenogenesis, male killing and sperm-egg incompatibility (Werren et al. 2008). If a population is infected by , patterns of mitochondrial polymorphisms will be altered by natural selection that acts on these symbionts, either increasing or decreasing the frequency distribution of haplotypes within a population (Hurst and Jiggins 2005). Previous studies documented high infection rates of within many terrestrial as well as freshwater isopod species (e.g., Bouchon et al. 1998; Rigaud et al. 1999; Cordaux et al. 2012), including numerous species that have been analyzed in this study, e.g., , , and . However, it is very difficult to distinguish demographic variation from symbiont-induced effects of mitochondrial variability (see Hurst and Jiggins 2005).

Third, the amplification and sequencing of nuclear mitochondrial pseudogenes (numts) can obscure the true mitochondrial variability within a species (Bensasson et al. 2001). Numts are nonfunctional copies of mitochondrial DNA in the nuclear genome. As consequence of reduced selection pressure, nucleotide substitutions and insertions as well as deletions may introduce stop codons and shifts in the reading frame of these inactive copies (Buhay 2009; Schizas 2012). Various studies documented such numts for a number of different crustacean taxa (e.g., Buhay 2009; Baeza and Fuentes 2013; Kim et al. 2013). For isopods, however, numts have not been reported so far, and a careful inspection of our COI sequences revealed no double peaks and translation without stop codons.

Fourth, many oniscid species, e.g., (Latreille, 1804), (De Geer, 1778), or (Scopoli, 1763), are characterized by atypical mitochondrial DNA structures that are composed of linear monomers and circular dimers, generating different mitochondrial lineages (Doublet et al. 2012, 2013). There is also a possible link between such atypical mitochondrial DNAs and heteroplasmy (i.e., the mixture of mtDNA genotypes within an organism) which has been documented for various woodlice in the past (Doublet et al. 2008, 2012). However, only a few studies are available until now, and most details still remain unclear.

Finally, distinct mitochondrial lineages that correlate with high genetic distances can give evidence for the existence of currently overseen cryptic species. Considering the previous discussed aspects, however, additional morphological and/or nuclear DNA sequence data are essential for a verification of truly distinct lineages. For freshwater and terrestrial isopods, a few studies demonstrated such integrative taxonomic approaches (McGaughran et al. 2005; Santamaria et al. 2017; Santamaria 2019). In terms of the analyzed taxa, no previous studies discussed the existence of cryptic species, and all specimens were checked and determined carefully before molecular works started.

Based on the given data we are currently unable to clarify the reasons of the observed high intraspecific variability within some of the analyzed species in detail. We suggest, however, that the detected high distances result from i) phylogeographic effects, ii) infections, iii) atypical mitochondrial DNAs and/or heteroplasmy, or, most likely, iv) various combinations of these phenomena in many cases. More specimens from different geographic regions as well as additional nuclear markers should be analyzed to verify this in detail. Despite these high intraspecific distances and multiple BINs for some species, however, high interspecific distances in combination with monophyletic lineages allow a correct determination of all studied taxa.

Conclusions

The development of new sequencing technologies changed biological science significantly. As a consequence, DNA-based approaches have become more and more popu­lar for the assessment of biodiversity and identification of specimens. Parallel analysis of thousands of specimens, bulk samples (metabarcoding) or environmental DNA (eDNA) will become routinely used techniques in modern species diversity assessment studies (e.g., Shokralla et al. 2012; Moriniere et al. 2016; Brauckmann et al. 2019; Hardulak et al. 2020). Whereas hypervariable regions of nuclear rRNA genes or other mitochondrial gene fragments may represent useful markers for such studies (e.g., Mohrbeck et al. 2015; Gillet et al. 2018; Lopez-Escardo et al. 2018), COI has become the most popular and efficient marker of choice (e.g., Andujar et al. 2018; Curry et al. 2018; Brauckmann et al. 2019; Hausmann et al. 2020). All these approaches, however, rely highly on comprehensive on-line reference libraries of correctly identified specimens (e.g., Brandon-Mong et al. 2015; Creer et al. 2016; Staats et al. 2016). Ideally, such libraries include sequence data of a species` complete distribution range that can provide additional information of phylogeographic substructures that are well-known for many species (e.g., Gentile et al. 2010; Raupach et al. 2014; Paill et al. 2021).

The necessity of DNA barcode reference libraries is also important for the modern molecular-based analysis of soil biodiversity (Taberlet et al. 2012; Orgiazzi et al. 2014; Rota et al. 2020). Reference libraries have been already published for a variety of typical soil-inhabiting taxa, e.g., earthworms (Porco et al. 2013; Pansu et al. 2015; Sun et al. 2018), mites (Young et al. 2012; Young et al. 2019), springtails (Hogg and Hebert 2004; Porco et al. 2013), spiders (Astrin et al. 2016), myriapods (Spelda et al. 2011) and ground beetles (Raupach et al. 2016; Raupach et al. 2018). In our present study we lay the foundations for a comprehensive DNA barcode data set for terrestrial and freshwater isopods of Central Europe.