DNA barcoding reference libraries of Italian Plecoptera: a gap analysis

Abstract

Plecoptera is a small order of aquatic insects that is considered one of the most endangered groups of insects due to increasingly altered freshwater ecosystems. Plecoptera nymphs can be challenging to identify at the species level because diagnostic characters for most species are lacking, while adults can be difficult to collect routinely given their short life span. Species identification by DNA barcoding is increasingly used as an alternative to morphological identification, but gaps and inaccuracies in reference databases needed for taxonomic assignment can undermine the utility of barcoding in real case studies. Here we aim to: i) quantify the number of Italian species of Plecoptera with barcodes from specimens collected worldwide and from Italy; ii) perform a regional assessment of DNA barcoding coverage; iii) calculate the intraspecific distance among available sequences to evaluate the potential presence of errors and cryptic species. As reference databases, we used both a non-curated database (BOLD) and a curated database (MIDORI2) to test the effect of sequence selection on the availability of reference sequences. We found that 67.6% and 51.8% of the Italian Plecoptera species were represented in BOLD and MIDORI2. Most of the available sequences were obtained from specimens collected outside Italy, with only 21.8% and 13.5% of the Italian species having sequences from specimens collected in Italy. Endemisms were poorly represented, and intraspecific distances for some species were high, which suggest cryptic diversity or erroneous assignments. Our results support the growing need to increase international cooperation among barcode initiatives and to promote the integration between molecular biologists and zoologists to exploit the full potential of DNA barcoding to protect biodiversity.

Keywords:

• freshwater insects
• metabarcoding
• BOLD
• MIDORI2
• endemisms

Introduction

Quantifying spatial and temporal trends of species distribution is essential to preserve biodiversity and to plan effective conservation strategies. Of metazoans, insects show astonishing diversity, with nearly one million species described worldwide and a few million yet to be discovered (Stork 2018). Assessing the diversity of insects can be challenging because of different impediments, including the considerable expertise needed for morphological identification, the lack of identification keys and cryptic taxa (Miura 2005; Laini et al. 2023). This is especially true for many freshwater insects that spend most of their life in aquatic habitats as larvae with a relatively short lifespan as terrestrial adults (e.g., Plecoptera, Ephemeroptera, Diptera, and Trichoptera). The larvae of several freshwater insect groups are often difficult to identify at the species level either because of diagnostic characters for both the first instars and mature larvae are lacking or because labour-intensive protocols are required (Zhou et al. 2010; Mynott et al. 2011). Moreover, the short lifespan of adults can complicate their collection which, thus, challenges the rapid assessment of biodiversity in aquatic environments (Cordero et al. 2017). In this context, genetic tools for the identification of preimaginal stages could boost our capacity to assess biodiversity.

DNA barcoding has emerged as a valuable tool for identifying species because it overcomes most of the challenges related to morphological identification (Hebert et al. 2004). This technique allows species-level assignments of taxonomically unknown specimens by comparing a short DNA sequence to those of known species stored in a reference database. The quality of the reference database, including both the quality of the DNA sequences and the accuracy of taxonomic assignments, is crucial for barcoding to be effective. However, several gaps and weaknesses have been identified in existing reference databases. Gaps, i.e., lack of reference sequences for species within a specific taxon, are particularly relevant for some taxa that are difficult to identify morphologically (e.g., Diptera) (Weigand et al. 2019). Moreover, within a specific taxon, easy-to-identify species are represented by numerous reference sequences, while others are poorly represented or even lack reference sequences (Csabai et al. 2023). Beyond the quantity of barcoded species, sequence quality is important for molecular identification. Cheng et al. (2023) recently showed that specimen misidentification, sample confusion, and contamination are not rare in two of the most comprehensive reference barcode databases. Problems associated with the quality of sequences can be mitigated using curated reference databases. Curated databases include reference sequences selected according to sequence quality and taxonomic assignment accuracy (Machida et al. 2017; Leray et al. 2022; Magoga et al. 2022). Furthermore, using reference sequences of specimens collected from the same area of the specimens where barcoding is needed may enhance correct taxonomic assignments (Li et al. 2022).

Plecoptera, or stoneflies, is a small order of hemimetabolous insects with 17 families and more than 3900 species (Fenoglio et al. 2021; DeWalt et al. 2023). In the last few decades, as nearly 43 new Plecoptera species per year have been described globally, it is likely that the total number of described species will significantly increase in the future (DeWalt & Ower 2019). Stonefly nymphs are generally stenoecious, while adults’ flight ability is poor. These combined factors explain the large number of endemisms within the order (Fochetti & Tierno de Figueroa 2008a). Stoneflies spend most of their life as immature stages: nymphs live from one to several years and development occurs through a number of moults that vary from 12 to 30. Alternatively, adults have a relatively short lifespan (ranging from days to a few weeks) (Fochetti & Tierno de Figueroa 2008a; Feeley et al. 2016). Species identification occurs mostly based on adult characteristics, such as wing venation, genitalia, or other sexual elements. Species identification for nymphs can be problematic because a number of them have not yet been described (see Fochetti & Tierno de Figueroa 2008a). Stoneflies are considered one of the most endangered groups of insects due to the increasingly altered lotic systems, which leads to many stonefly species being reduced to small and isolated populations, or even becoming extinct (Fenoglio et al. 2010; Lee et al. 2022). Tools for the identification of nymphs can enhance our capacity to protect plecopteran diversity via the discovery of new species and an increased capacity to quantify spatial and temporal patterns of diversity.

Here, we aimed to quantify the DNA barcode coverage of the Italian species of Plecoptera in a non-curated (BOLD, Ratnasingham & Hebert 2007) and a curated (MIDORI2, Leray et al. 2022) reference database, which are both widely used for molecular species identification. Our objectives were to: i) quantify the number of Italian species of Plecoptera with barcodes from specimens collected worldwide and from Italy; ii) perform a regional assessment of DNA barcoding coverage; iii) calculate the intraspecific distance among the available sequences to evaluate the potential presence of errors and cryptic species.

Material and methods

Italy hosts 170 Plecoptera species belonging to 23 genera and seven families. Nearly 30% of the species are endemic, and their geographical distribution is poorly known (Fochetti 2020). The list of Italian Plecoptera was obtained from Fochetti (2020) and has been updated with the recently described species Zwicknia bifrons (Newman, 1938), Protonemura bispina Vinçon, Ravizza & Reding, 2021, and Protonemura pennina Vinçon, Ravizza & Reding, 2021 by Reding (2021) and Vinçon et al. (2021). Synonyms were retrieved from Fochetti and Tierno de Figueroa (2008b) and checked against the Plecoptera species file (DeWalt et al. 2023). Some changes in the taxonomy of BOLD and MIDORI2 were made to homogenise the taxonomic assignments in both datasets. In particular, Dictyogenus alpinum (Pictet, 1841), Isogenus nubecula Newman and Siphonoperla italica (Aubert, 1953) were used instead of Dictyogenus alpinus (Pictet, 1841), Isogenus nubeculum) and Siphonoperla torrentium italica (Aubert, 1953). Moreover, both species Leuctra vinconi Ravizza & Ravizza Dematteis, 1993 and subspecies Leuctra vinconi vinconi Ravizza & Ravizza Dematteis, 1993 were employed to retrieve sequences for L. vinconi. Similarly, both species Capnia vidua Klapálek, 1904 and subspecies Capnia vidua vidua Klapálek, 1904 were utilised to retrieve sequences for C. vidua. The sequences of valid species and synonyms were retrieved from BOLD (Ratnasingham & Hebert 2007) with the function bold_specimens() implemented in the bold (Dubois & Chamberlain 2023) package for the R programming language for statistical computing (R Core Team 2017). The sequences retrieved from BOLD (accessed on 2023-08-09) included also those mined from GenBank. MIDORI2 version MIDORI2_UNIQ_NUC_GB255_CO1_RAW of the GenBank255 database was retrieved from website www.reference-midori.info. The MIDORI2 sequences were filtered to keep only plecopteran species.

Gaps were removed from the sequences of both databases with the function RemoveGaps() of the Biostrings package (Pagès et al. 2023). Number of sequences, as well as mean, minimum and maximum sequence lengths, were calculated for each species. The mean, minimum and maximum ratios between the number of non-ambiguous nucleotides and the total number of nucleotides were calculated for each species as a measure of sequences quality.

The number of sequences of specimens collected in Italy was also quantified. For BOLD, this was done by setting the geo argument to Italy in the bold_specimens() function. For MIDORI2, country information was retrieved from sequence accession numbers using the function entrez_fetch() of the rentrez package (Winter 2017).

For each species, sequences were aligned with the function AlignSeqs() of the R package DECIPHER (Wright 2016) with default settings. Intraspecific distances were calculated from aligned sequences using the Kimura’s 2-parameters distance implemented in the function dist.dna() of the R package ape (Paradis et al. 2019). Intraspecific distances served as an indication of the presence of identification errors or cryptic species. The code used for the data analysis is available in the Supplementary Material.

Results

Species-level analysis

More species with a barcode were found in BOLD (n = 115) than in MIDORI2 (n = 88). Species with at least a barcode in BOLD belonged to valid species excepted S. italica whose sequences were attributed to S. torrentium italica. Of the 115 species found in BOLD, 107 had public records, while eight had private records (Leuctra apenninicola Ravizza 1988, Leuctra fochettii Vinçon & Graf 2017, Leuctra grafi Vinçon & Vitecek 2017, Leuctra insubrica Aubert 1949, Leuctra ravizzai Ravizza Dematteis & Vinçon 1994, Nemoura fulviceps Klapálek 1902, Protonemura austriaca Theischinger 1976, Protonemura bipartita Consiglio 1962). On average, 20.6 sequences per species (range 1–117) were found in BOLD, with a mean length of 645.9 (285–1531) and an average ratio between non-ambiguous and total number of bases of 0.998 (0.669–1.000). The number of sequences from specimens collected in Italy was 325. These sequences were representative of 37 species (21.8% of the Italian Plecoptera species). Fourteen (10 public and 4 private) of the 55 endemic species was represented in BOLD. Further calculations were done by referring to public data only.

In MIDORI 2, 88 species were found, with an average of 12.4 sequences (range 1–66) per species with a mean length of 666.1 (range 300–1545) and an average ratio between non-ambiguous and total number of bases of 0.997 (0.907–1.000). The number of sequences from specimens collected in Italy was 147. These sequences were representative of 23 species (13.5% of the Italian species of Plecoptera). Finally, three out of the 55 endemic species were represented in MIDORI2.

Regional-level assessment

Regions of Northern Italy generally had higher species richness than Southern regions (Figure 1(a)). DNA barcoding coverage did not mirror species richness for both BOLD (Figure 1(b)) and MIDORI2 (Figure 1(c)), with a coverage ranging from 0.42 to 0.97 (mean = 0.74) in BOLD and from 0.42 to 0.89 (mean = 0.67) in MIDORI2. Endemic species (Figure 2(a)) showed lesser DNA barcoding coverage in both BOLD (Figure 2(b)) and MIDORI2 (Figure 2(c)). For endemic species, DNA barcoding coverage ranged from 0.1 to 0.67 in BOLD (mean = 0.28) and from 0 to 0.5 in MIDORI2 (mean = 0.15).

Figure 1. Number of plecopteran species in Italian regions and barcode coverage in BOLD and MIDORI2.

Figure 2. Number of plecopteran endemisms in Italian regions and barcode coverage in BOLD and MIDORI2. Regions coloured in grey have 0 endemic species with a barcode.

Intraspecific distances

The intraspecific distance of species with more than one unique sequence varied among species and reference databases. In BOLD, Protonemura nitida (Pictet 1836) (0.10) showed the highest mean intraspecific distance, followed by Protonemura nimborella (Mosely, 1930) (0.10) and Leuctra fusca (Linnaeus, 1758) (0.08) (Figure 3). In MIDORI 2, P. nimborella and Perlodes intricatus (Pictet, 1841) had the highest mean intraspecific distance (0.13), followed by P. nitida (0.12) (Figure 3).

Figure 3. Intraspecific distance of the sequences of Italian stoneflies stored in BOLD and MIDORI2. The dashed red line represents the likely maximal value for intraspecific divergence of 0.035 according to (Gattolliat et al. 2016).

Genus-level analysis

All the genera had at least one species with a barcode, although the number of species with a barcode varied within genus and between reference databases (Figure 4). Species-poor genera had generally a complete coverage (Amphinemura Ris, 1902, Besdolus Ricker, 1952, Brachyptera Newport, 1848 in BOLD, Capnia Pictet, 1841, Capnopsis Morton, 1896, Chloroperla Newman, 1836, Dictyogenus Klapálek, 1904, Dinocras Klapálek, 1907 in BOLD, Isogenus Newman, Nemurella Kempny, 1898, Perla Geoffroy, 1762, Perlodes Banks, 1903, Rhabdiopteryx Klapálek, 1902 in BOLD, Siphonoperla Zwick, 1967, Tyrrhenolecutra Consiglio, 1957, Xanthoperla Zwick, 1967). The genera with the largest number of species without a public barcode were Leuctra Stephens, 1836 (24 in BOLD, 33 in MIDORI2), Protonemura Kempny, 1898 (20 in BOLD, 22 in MIDORI2), Isoperla Banks, 1906 (7 in BOLD, 9 in MIDORI2) and Nemoura Latreille, 1796 (7 in BOLD, 9 in MIDORI2).

Figure 4. Number of Italian Plecoptera species with and without barcodes in a) BOLD and b) MIDORI2, divided by genera.

Discussion

The quality of reference databases is crucial for an accurate taxonomic assignments of sequences in DNA barcoding. Using an updated list of Italian stoneflies to query curated and non-curated reference databases, we found gaps and inaccuracies that can bias species identification by molecular methods. A regional assessment showed that gaps were not spatially structured, although regions of northern Italy hosted a higher species richness than other regions due to the abundance of mountain reliefs that host cold and well oxygenated running waters. Gaps were particularly severe for endemic species, which suggest potential problems in using molecular species identification for conservation actions and reinforce the need to address and coordinate efforts of national and supra-national barcoding initiatives (Csabai et al. 2023).

DNA barcoding coverage of Italian Plecoptera in BOLD and MIDORI2

Our DNA barcoding coverage assessment for Italian Plecoptera highlighted both gaps and inaccuracies in species labelling. Coverage depended on the employed reference database, with 67.6% of species covered in BOLD and 51.8% species covered in MIDORI2. Most of the missing sequences were included in a few genera, like Leuctra and Protonemura although less speciose genera such as Brachyptera and Taeniopteryx still had nearly 50% of species without barcodes. The range of sequence quality measured as the number of non-ambiguous bases over the sequence length was smaller in MIDORI2 compare to BOLD, indicating higher sequence quality in MIDORI2. Our results thus suggest that non-curated databases like BOLD can have more coverage than curated datasets like MIDORI2, but taxonomic assignments in non-curated datasets for some species could be hindered by a lower sequence quality (Leray et al. 2022). This is especially important if species are represented by only low quality-sequences. Our sequence quality assessment did not take into account other strategies that could have increased our capacity to discriminate between good and poor quality (e.g., translation of sequences in amino acids).

Inaccuracies were few and related mostly to recent changes in taxonomy. S. italica is currently ranked as a valid species and was previously reported as a subspecies of Siphonoperla torrentium (Pictet, 1841) (Weiss et al., 2012). However, S. italica is labelled as S. torrentium italica in both BOLD and MIDORI2, which denote a temporal mismatch between the species descriptions by taxonomists and taxonomic assignments of old sequences in reference databases. Other temporal mismatches are more subtle and can be difficult to spot. For example, although Leuctra biellensis Festa, 1942 has been recently reinstated as a valid species (Vinçon et al. 2018), our query did not retrieve any sequence assigned to this species. The sequence of one of the specimens reviewed by Vinçon et al. (2018) and identified as L. biellensis is reported in BOLD as Leuctra nigra (Olivier, 1811). Similarly, P. bispina and P. pennina, the new species of Protonemura described by Vinçon et al. (2021), are reported in BOLD as Protonemura auberti Illies, 1954. Problems due to synonyms and lack of updated information being associated with reference sequences are widespread and affect different orders of aquatic insects including Diptera (Gadawski et al. 2022), Heteroptera, and Coleoptera (Csabai et al. 2023).

Intraspecific distances

We found a wide range of intraspecific distances for both databases, with some species exceeding values of 0.2 (e.g., P. nitida, P. nimborella). Previous works on intraspecific distances on Plecoptera have shown limited intraspecific divergence, although values exceeding 0.035 had been found for Leuctra nigra (Gattolliat et al., 2016), Isoperla adunca Jewett, 1962 (Gill et al. 2015) and Riekoperla alpina McLellan, 1971 (Mynott et al. 2011) among others. High intraspecific distance values can be the result of errors in different barcoding process steps, such as wrong identification, mislabelling of samples and cross-contamination (Pentinsaari et al. 2020). This was probably the case of one record of P. nimborella (accession number GBMH8738–13) that clustered with Leuctra inermis (Kempny, 1899) when directly searched in BOLD. Problems are not limited to errors within Plecoptera; a sequence assigned to Isoperla rivulorum (Pictet, 1841) (accession number MK584518) is probably a sequence of the mollusc Ancylus fluviatilis O.F.Müller, 1774. Besides obvious errors, marked genetic distances can indicate the presence of cryptic species (Morinière et al. 2017; Suttinun et al. 2022).

Endemisms

Our analysis showed that endemic species are poorly represented in both BOLD and MIDORI2. We suggest that this situation is due to the combination of two factors. On the one hand, most endemic species have been described in the past century when DNA barcoding was not yet established to support species description (e.g., Leuctra annae Consiglio 1975, Leuctra archimedis Consiglio, 1968, Protonemura lagrecai (Aubert 1954)). On the other hand, barcoding initiatives are scarce in Italy as suggested by our results on the country of specimen collection. Our estimates of barcodes from individuals collected in Italy could be underestimated because of lack of metadata associated with records. For example, Siphonoperla torrentium italica has been collected in Italy (Weiss et al. 2012), but the collection country was not reported either in BOLD or in GenBank. However, the scarcity of barcoding initiatives of Mediterranean countries, Italy included, compared to Northern European countries has also been highlighted by Csabai et al. (2023) for Heteroptera and Coleoptera. Poor DNA barcoding coverage is even more important given that southern European countries host rich and diversified fauna, including Plecoptera (Tierno de Figueroa et al. 2013). To apply DNA barcoding as a tool to protect endemic or endangered species, more initiatives are needed in Italy to use genetic approaches to infer spatial and temporal patterns of species distribution.

Implication of the results for stoneflies conservation

Identifying gaps in current reference databases can contribute to address future efforts in selecting species for DNA barcoding to boost biodiversity conservation (Hobern & Adamowicz 2021). This is even more important if we consider the growing importance of genetics in species identification. For example, DNA metabarcoding provides high-resolution taxa lists from bulk samples (Taberlet et al. 2012) or environmental DNA (eDNA) (Blackman et al. 2019; Laini et al. 2020). DNA metabarcoding is a cost-effective technique that can be used for large-scale monitoring programmes where morphological identification is logistically or financially impractical (Liu et al. 2020). Well-curated reference databases are crucial for metabarcoding to be effective in real case studies (Leese et al. 2016; Keck et al. 2017).

The importance of Plecoptera goes beyond the entomological and zoological contexts because these organisms are widely used in biogeographic studies and especially in biomonitoring. Stonefly nymphs are very sensitive to environmental conditions and their use as bioindicators is extremely diffuse. In particular, their presence is considered pivotal in many indices throughout the world (Bonada et al. 2006; Guareschi et al. 2017), and for this reason their correct taxonomic identification is critical. Moreover, current Italian biomonitoring protocols require genus- or even family-level identification. We argue that species-level identification can boost our capacity to infer causes of habitat degradation as subtle changes in community composition by environmental and anthropogenic causes. For example, the species within the genus Taeniopteryx have different responses to desiccation and, thus, to flow intermittence (López-Rodríguez & de Figueroa 2006). Hence, species-specific life history information can be used along with trait-based approaches to detect the effect of flow intermittence on invertebrate communities (Piano et al. 2020; Soria et al. 2020). Furthermore, several genera show interspecific variability in thermal requirements (e.g., Leuctra) (Tierno de Figueroa et al. 2010). Species-level identification within these genera can help to track the global warming effects on freshwater insects. Our results support a growing awareness that strengthening the international cooperation among barcoding initiatives and promoting the integration between molecular biologists and zoologists are required to exploit the full potential of DNA barcoding to protect biodiversity.

Acknowledgments

The authors thank Manuel Jesús López-Rodríguez, José Manuel Tierno de Figueroa, Gilles Vinçon, Jean-Paul G. Reding, and two anonymous reviewers for their useful suggestions.

Disclosure statement

No potential conflict of interest was reported by the authors.

Data availability statement

The data used for this work are available from the corresponding author on reasonable request.

Supplementary material

Supplemental data for this article can be accessed online at https://doi.org/10.1080/24750263.2023.2298977