A DNA barcode library for the butterflies of North America

Introduction

DNA barcoding is an effective tool for addressing the widely recognized need for an improved understanding of biodiversity. By employing sequence diversity in short, standardized gene regions, such as the 648 base pair segment of the 5′ region of mitochondrial cytochrome c oxidase 1 (CO1) employed for the animal kingdom (Hebert et al., 2003), DNA barcoding allows both the identification of specimens and the discovery of new species. Since its introduction, this marker has been adopted in fields ranging from population genetics (Hajibabaei et al., 2007), phylogenetics (Hajibabaei et al., 2007), and phylogeography (Dapporto et al., 2019) to ecology (Valentini, Pompanon & Taberlet, 2009) and conservation (Dincă et al., 2018). It has also gained application in contexts ranging from the detection of marketplace fraud (Galimberti et al., 2013) to the suppression of illegal trade in endangered species (Rehman et al., 2015).

DNA barcoding enables the identification of specimens without morphological analysis by querying their CO1 sequences against a reference library of DNA barcodes obtained from carefully identified vouchers. These reference sequences are curated in the Barcode of Life Datasystem (BOLD) (boldsystems.org) (Ratnasingham & Hebert, 2007), an informatics platform that also hosts collateral data such as specimen images and collection data. Aside from identifying specimens, DNA barcoding can help to delineate species boundaries (Cariou et al., 2020), an important task since species play a central role in biodiversity assessments and conservation actions. Rather than simply presuming that the current taxonomic system is valid, DNA barcoding provides a basis for testing this assertion. Prior studies have shown that closely allied congeneric species of Lepidoptera typically show more than 2% divergence (Hebert et al., 2003; Huemer et al., 2014; Dincă et al., 2015). Although some sister species do show lower divergence (Burns et al., 2007; Cong et al., 2016), cases where species share barcode sequences can reflect over-splitting (Vila et al., 2010) or introgressive hybridization (Zakharov et al., 2009; Cong et al., 2017). Conversely, when members of a putative species show high sequence divergence, this often signals the presence of overlooked species (Hebert et al., 2004; Burns et al., 2007; Dincă et al., 2013). DNA barcode data has also enabled the development of algorithms that employ sequence information for species delimitation. The latter methods cluster specimens into Molecular Operational Taxonomic Units (MOTUs) that have been shown to correspond closely with recognized species in groups with well-established taxonomy (Ratnasingham & Hebert, 2013).

Most studies have tested the capacity of DNA barcodes to discriminate species when viewed from a local or regional context, and only a few have examined resolution at a continental scale (e.g.,  Kerr et al., 2007; Lukhtanov et al., 2009; Bergsten et al., 2012; Huemer et al., 2014; Zahiri et al., 2017; Dincă et al., 2021). The latter studies are important because they can reveal cases of low interspecific divergences, potentially reducing the effectiveness of DNA barcoding for species delimitation. Prior studies have shown the general effectiveness of DNA barcoding for butterflies (Lukhtanov et al., 2009; Dincă et al., 2011; Dincă et al., 2015; Lavinia et al., 2017; Dincă et al., 2021), but have also exposed discordances with current taxonomy including probable cases of synonymy (e.g.,  Vila et al., 2010) and frequent instances of overlooked species (e.g.,  Hebert et al., 2004; Burns et al., 2008).

The butterfly fauna of North America has seen more intensive morphological study than any other comparably diverse insect lineage on this continent (Warren et al., 2012). Despite this attention, there remains uncertainty in the status of many taxa, often reflecting the subjectivity inherent in decisions on species boundaries based on morphology alone. The present study assembles a comprehensive DNA barcode library for the butterfly fauna of North America, delivering an identification system for most of these species while testing the current taxonomy. This work also provides an overview of patterns of genetic diversity and offers insights on mechanisms responsible for shaping the genetic diversity of the butterflies of this continent.

Methods

Sampling

This study sought to recover DNA barcodes for the butterfly fauna of Canada and USA. BOLD hosts a checklist for 846 species (CL-NABUT) derived from Pelham’s list (Warren et al., 2012) with a few changes based on recent publications. Table S1 provides a condensed version of this checklist for the 648 species with persistent populations in North America, excluding those introduced by humans.

The sampling program aimed to capture geographic and phylogenetic diversity for each species in continental North America (i.e., islands beyond the continental shelf were excluded). Specimens from Canada (6,935) and USA (7,037) were sequenced when possible, but this left some gaps which were filled by analyzing specimens from Central America (602), South America (27), Europe (7), Asia (4) and unvouchered (11) records from GenBank (Fig. 1). Overall 14,626 specimens were analyzed, and associated metadata are available on BOLD (v4.boldsystems.org) in the public dataset “DS-USCANLEP” (dx.doi.org/10.5883/DS-USCANLEP). From this total, 10,425 vouchers are held in public natural history collections, 3,864 in private collections, 259 derive from GenBank, and 78 were unvouchered. Permission from all the institutional and private collections was obtained to access and study the records. The Centre for Biodiversity Genomics made the largest contribution (4,474 records), followed by the Canadian National Collection (1,771) and the Smithsonian’s National Museum of Natural History (1,010).

iNEXT (INterpolation and EXTrapolation) (Hsieh, Ma & Chao, 2016), was employed to estimate sampling completeness using R (R Studio Team, 2016). This approach implements the Chao1 diversity estimator (Chao, 1984) to generate accumulation curves that can be used to estimate the total haplotype diversity in a species. This asymptotic value was compared with the observed haplotype diversity to quantify sampling completeness, an approach used to estimate coverage for European butterflies (Dincă et al., 2021). Because the present study targeted species resident in North America, levels of genetic diversity for introduced species and tropical strays (for which sampling was limited) were likely to be underestimated. As a result, estimates of sampling completeness excluded species whose distribution primarily falls outside North America. Moreover, species represented by fewer than six specimens were also excluded, reducing consideration to 402 of the 648 species (Table S1). For each species in the barcode library, we recorded the number of specimens (N), the number of observed haplotypes (H), the fraction of haplotype diversity retrieved (R), and the number of additional haplotypes which remain to be sampled (L).

Results

Genetic diversity

Maximum intraspecific distances averaged 0.97% (range 0–8.4%) while the nearest neighbor distance averaged 3.7% (range 0–14.3%), almost 4-fold higher than the maximum intraspecific distance (Fig. 2, Table S2). A barcode gap was present in 573 of 755 species (75.9%) represented by at least two individuals (Fig. 2, Table S2).

The Bayesian trees revealed that 604 of the 755 species (80%) represented by two or more individuals formed monophyletic groups, while the other 151 (20%) were either paraphyletic or polyphyletic (Figs. S7–S12, S13 and Table S2). Species with just a single barcode sequence were necessarily excluded from this analysis, but 55 of the 59 possessed barcodes distinct from their nearest neighbor. While discrimination between paraphyly and polyphyly is not essential for specimen identification, it is critical to distinguish those species characterized by overlapping phylogenetic branches from those species sharing barcodes with their nearest neighbor(s). Twenty-six species (3.4%) fell in the first category and 125 (16.6%) species in the latter (Fig. S13).

Use of a fixed distance threshold (2.5%) exposed 79 cases (9.7%) where intraspecific distance was above the set threshold and 324 (39.8%) cases where interspecific divergence was below it. MOTU delineation revealed a variable number of entities depending on the method employed. BIN analysis was performed on all but one species (the sole sequence for Calephelis rawsoni did not qualify for analysis) and revealed 772 BINs, comprised of 540 Matches (66.4%), 55 Splits (6.8%), 181 Merges (22.3%), and 37 Mixtures (4.6%). By comparison, GMYC analysis generated 862 taxonomic units, partitioned in 527 Matches (64.7%), 63 Splits (7.8%), 150 Merges (18.4%), and 74 Mixtures (9.1%) (Figs. S7–S12). Overall, the three analyses provided concordant support for 369 species recognized by the current taxonomy, but they also revealed 24 species split into two or more units, 124 grouped with one or more nearest neighbor(s), and 34 in both previous categories. The performance of the three methods is compared in Table S2.

Discussion

This study generated a DNA barcode reference library for 96% of the North American butterfly fauna (814 of 846 species) with an average of 18 records per species. This level of sampling captured 67% of the estimated haplotype diversity, but at least 4,705 haplotypes await detection. The present library was effective in assigning newly encountered specimens to either a species or to a small number of closely allied species. Because butterflies are key bioindicators (Syaripuddin, Sing & Wilson, 2015) and umbrella species (New, 1997), this library facilitates their use in tracking the impacts of habitat loss, fragmentation, and climate change. Along with the reference barcodes available for other taxonomic groups (Ratnasingham & Hebert, 2007), the present work also provides a basis for exploring species interactions.

The sequence records analyzed in this study were mainly obtained during two ways. The first approach involved the collection of fresh material and lasted about two decades (2000–2019), while the latter spanned three years (2015–2017) and consisted of the targeted analysis of specimens held in two major natural history collections—the Smithsonian’s National Museum of Natural History and the Canadian National Collection. Previous studies on Lepidoptera, largely regional in scale, indicated that most species are separated from their nearest neighbor by a barcode gap (e.g., Janzen et al., 2005; Hajibabaei et al., 2006; Lavinia et al., 2017). Increased geographical coverage should lead to higher intraspecific genetic variation as a result of the isolation-by-distance effect (Wright, 1943), and reduced interspecific divergences as more species are analyzed (Avise, 2000). Both factors should reduce the difference between intra- and inter-specific divergence. It needs emphasis that exceptions to this model were observed in butterflies and that narrow suture zones with high diversity alternate with extended regions with little variation (e.g., Gompert et al., 2010; Dapporto et al., 2019; Platania et al., 2020). This continental-scale study revealed that average maximum CO1 divergences within species was nearly four times less than the average minimum distance to the nearest neighbor. Reflecting this fact, 76% of the species displayed a barcode gap ensuring their unambiguous identification. Identification of specimens using the criterion of monophyly raised identification success as 80% of North American butterfly species formed a monophyletic cluster. These results approximate those obtained in a study of European butterflies where monophyly was met for 94% of species at a regional level (Iberia) (χ2 = 23.1; p-value <0.001) versus 85% for the continent (Dincă et al., 2015) (χ2 = 2.71; p-value = 0.099). As the European study only considered approximately 60% of the fauna, the proportion of monophyletic species is likely to decline with increased coverage. In a study examining how increased geographic scale affected the barcode gap in Asian butterflies, Lukhtanov et al. (2009) found that the barcode gap decreased with distance, but that species resolution was nearly unaltered. Interestingly, they showed that monophyly, unlike the barcode gap, was less affected by geographic coverage. The evidence for reduced identification success in the present study is likely explained, at least in part, by the 5-fold higher sampling effort (18 specimens/species) than in Lukhtanov et al. (2009) (3 specimens/species). Although the presence of a barcode gap or monophyly are sufficient conditions to ensure the correct assignment of a sequence to its correct species, these are not essential criteria. For example, a species whose component sequences are paraphyletic or polyphyletic can meet neither criteria, but can be perfectly diagnosable (Ross, Murugan & Li, 2008; Bergsten et al., 2012).

Because the capacity of DNA barcoding (Ratnasingham & Hebert, 2007) to deliver a correct identification ultimately depends upon the presence of diagnostic (or a diagnostic combination of) characters, sequence sharing by species indicates that their discrimination will be compromised. This study revealed that 15% of North American butterfly species share their DNA barcodes with another species. Work on European butterflies revealed 3% barcode sharing at a regional level (Iberia), and 7% at a continental scale (Dincă et al., 2015). Although the latter value is about half that observed for North America (χ2 = 12.6; p-value <0.001), it was inferred based on 60% of the European fauna and more comprehensive taxonomic and geographic coverage will almost certainly increase the incidence of barcode sharing for European butterflies. Because the level of barcode sharing does not only depends on geographic coverage, it is important to ascertain the factors underlying this pattern. First, both the incomplete sorting of ancestral polymorphisms and introgression can lead to barcode sharing, particularly between recently diverged species. The former factor should be less important for mitochondrial than nuclear genes because their lower effective population size facilitates the loss of ancestral polymorphisms. The impact of introgression is more controversial. Although Haldane’s rule predicts that the heterogametic sex (females in Lepidoptera) of hybrid individuals is not likely to pass on mitochondrial DNA (Haldane, 1922), introgression has often been reported in butterflies (e.g., Sperling, 1993; Gompert et al., 2006; Wahlberg et al., 2009). This could be explained by a more general hypothesis, broadly applicable to all organisms, suggesting that low purifying selection on introgressed mitochondrial genes favours the transfer of these elements (over nuclear ones) across species boundaries ( Harrison, 1993). Another feature that facilitates contact between closely related lineages, increasing the likelihood for introgression, is the high dispersal capability of some butterflies (Stevens, Turlure & Baguette, 2010). Second, although butterflies possess a relatively well-established taxonomy, recent studies have exposed taxonomic uncertainty including cases of over-split species (e.g., Vila et al., 2010). Such unrecognized cases of synonymy can produce barcode sharing. Third, although the specimens examined in this study were identified by specialists, DNA barcode results revealed a number of misidentifications. In other cases, the discrimination of morphologically similar species is so difficult that diagnostic characters might have been misinterpreted inflating the incidence of barcode sharing. A full investigation of the roles played by these three factors is beyond the scope of this study, however, an exhaustive study of 42,000 specimens representing nearly 5,000 species of European Lepidoptera showed that just 40% of non-monophyletic species were generated by biological factors (i.e., introgression, incomplete lineage sorting) while 60% reflected methodological problems such as misidentifications (Mutanen et al., 2016).

The incidence of barcode sharing in North American butterflies varies nearly 5-fold with latitude, being far higher among species found in the North/alpine (42%) than at Mid-latitude (25%) or South (9%) locales. A similar pattern was evident for 1541 North American noctuoid moth species as the Canadian fauna showed 10% barcode sharing (Zahiri et al., 2014), versus 7% for the continent (Zahiri et al., 2017).

Although the main scope of this work was to build a comprehensive barcode library for the identification of North American butterflies, effort was made to capture diversity from regions where contact zones, hybrid zones, and phylogeographic breaks congregate (Swenson & Howard, 2005). This strategy aimed to maximize the recovery of haplotype diversity, yet ascertainment bias is inevitable and likely impacted both geography (Yang, Ma & Kreft, 2013) and taxonomy (Troudet et al., 2017). While under-sampling can exaggerate the sequence divergence between species, comprehensive knowledge of intraspecific diversity will decrease interspecific distances and expose barcode sharing (e.g., Wiemers & Fiedler, 2007; Dasmahapatra et al., 2010). Our sample sizes were similar (40 specimens/species) for two regions (North/alpine, Mid-latitude), but lower in the South (16 specimens/species). However, iNEXT indicated that similar proportions of CO1 diversity were recovered from each bioregion (66%—North/alpine, 64%—Mid-latitude, 71%—South). Although iNEXT brings statistical rigor to such estimates, its accuracy closely depends on sampling quality. For instance, under-sampling biodiversity hotspots would give the illusion of low diversity and inflate estimates of sampling completeness. This result suggests that the observed differences in barcode sharing are not an artefact of varied sampling coverage. This pattern could well reflect the different exposure of species in each bioregion to the impacts of Quaternary glaciation (Hewitt, 1996). Species in northern regions could have experienced recurrent cycles of isolation in glacial refugia, leading to subsequent opportunities for secondary contact and mitochondrial exchange (Hewitt, 2000). Moreover, the small size of the populations at the leading edge of the species distribution might aid the fixation of introgressed mitochondria (Kingman, 1982). Another consequence of rapid expansion into deglaciated habitats would be low density at the leading edge. This situation can create difficulties in finding a conspecific mate and can weaken isolating mechanisms (e.g., Shelly & Bailey, 1992; Alatalo et al., 1998; Willis, Ryan & Rosenthal, 2011), leading to heterospecific mating and introgression (e.g., Wirtz, 1999; Randler, 2002). Not only scarcity of conspecifics could favor hybrid formation, but it could also enhance their persistence because of decreased competition with parental populations (Arnold, 1997).

Aside from this latitudinal pattern, barcode sharing varied nearly 4-fold among butterfly families (χ2 = 25.13; p-value <0.001), from 10% in Hesperiidae to 38% in Pieridae. Interestingly, a similar pattern was also observed at lower taxonomic rank, among genera, where the incidence of barcode sharing showed a weak and non-significant correlation with the number of species in a genus suggesting taxonomic localization. The genus Colias (Pieridae) was particularly impacted as all 22 species shared at least one of their barcode sequences with another species (Table S2), perhaps reflecting their recent radiation (Chew & Watt, 2006; Wheat & Watt, 2008). Not only this pattern explains the introgression of haplotypes between hybridizing species such as C. eurytheme/C. phildice in the eastern USA (Gerould, 1946; Jahner, Shapiro & Forister, 2011), and C. eurytheme/C. eriphyle in the west (Taylor Jr, 1972), but it also lays the foundation for operational issues such as misidentifications reflecting the unsettled taxonomy of the genus (Wheat & Watt, 2008).

DNA barcoding combined with species delimitation methods enables rapid, cost-effective surveys of biodiversity. While this is particularly beneficial for poorly-studied taxonomic groups, it can also disclose overlooked diversity in well-studied taxa. BIN analysis indicated that ∼11% of North American butterfly species (6.8% Splits, 4.6% Mixture) were split into two or more entities. The incidence of such cases mirrors values (9–12%) reported in prior studies on Lepidoptera (Huemer et al., 2014; Zahiri et al., 2014). GMYC analysis showed an even higher discordance with current taxonomy, showing that about 17% of species (7.8% Splits, 9.1% Mixture) potentially involve overlooked diversity. Interestingly, when the same approach was applied to 60% of European butterfly species, there was even higher discordance (28%) (Dincă et al., 2015). It is probable that the varied habitats in the Mediterranean basin (Blondel et al., 2010), coupled with the presence of southern refugia during the Pleistocene (Schmitt, 2007), created more genetic structure and/or speciation in Europe (Vodă et al., 2016; Dapporto et al., 2019; Scalercio et al., 2020). Employing a fixed divergence threshold (2.5%), about 10% of North American butterfly species exceeded this criterion. Similar values (8–12% with a 2% threshold) have been reported in other studies on Lepidoptera (Huemer et al., 2014; Zahiri et al., 2014). Based on these results, it is likely that a considerable number of cryptic species await description or that Evolutionary Significant Units (ESUs) within species deserve protection (Avise, 1989). Detailed studies (e.g., nuclear genetic, morphological, ecological) (DeSalle, Egan & Siddall, 2005) should be undertaken on these lineages (Polasky & Solow, 1999).

Conclusion

This study has generated one of the first continental-scale DNA barcode libraries for an entire taxonomic group. Beyond providing an identification system for most (>96%) North American butterflies, it creates the foundation needed to test the current classification. This library also delivers an overview of large-scale patterns of genetic diversity revealing cases of evolutionary interest such as potential hybridization and of importance to biodiversity conservation such as cryptic diversity and evolutionary significant units. As such, this study provides a basis for improving understanding of the mechanisms that have shaped genetic diversity in the North American butterfly fauna.

Supplemental Information