We apply a molecular approach to quantify the level of hymenopteran parasitoids infestation in the larvae of the marsh fritillary (Euphydryas aurinia), a declining butterfly species, in western Bohemia, Czech Republic, in two subsequent years. We used the novel primer HymR157 in combination with known universal 28SD1F to establish a PCR detection system which amplifies hymenopteran parasitoids, but not the lepidopteran host. In the 14 sampled E. aurinia colonies, the infestation rates per individuum were 33.3% and 40.2%; whereas per sampled larval colony, these were on average 38.5% (range 0–100) and 40.1% (0–78). The per-colony infestation rates correlated with the numbers of larval webs censused per colony the year prior to sampling the parasitoids, pointing to a time lag in parasitoid infestation rates. The levels of the hymenopteran parasitoid prevalence are thus relatively high, supporting the importance of parasitoids for the population dynamics of the threatened host. The detection primers we developed can detect a range of hymenopteran parasitoids on other butterfly hosts.

butterfly ecology, Braconidae, Lepidoptera, Marsh Fritillary, molecular detection, Nymphalidae, population dynamics

Hymenopteran parasitoids are one of the most diverse groups of animals in terrestrial ecosystems and play a key role in the natural regulation of their host populations (La Salle and Gauld 1991; Forbes et al. 2018). The impact of parasitoids on their hosts can vary depending on their ecological specialisation, but in general they are known to cause significant levels of mortality in their hosts (Hawkins 1994). High levels of parasitism may also pose a potential threat to many threatened butterfly species, especially to various specialists in fragmented landscapes (cf. Anton et al. 2007). Species of the genus Euphydryas Scudder, 1872 (Nymphalidae: Melitaeini) represent a suitable system for studying host-parasitoid interactions, as egg clutches and webs with gregarious caterpillars can easily be detected in the field (Stamp 1981; Hula et al. 2004; Johansson et al. 2019). Previous studies using rearing have shown that the level of parasitism varied annually and depended mainly on weather conditions, which play a key role in the synchronisation between larval development and the emergence of parasitoids (Porter 1979, 1983).

The Marsh Fritillary, Euphydryas aurinia (Rottemburg, 1775) is an EU-protected butterfly, declining in many European countries (van Swaay et al. 2010), including the Czech Republic (Hejda et al. 2017). It belongs to a genetically polymorphic group of closely related taxa, the “E. aurinia complex” (cf. Korb et al. 2016), with a wide Palearctic distribution and regional habitat and host plant specificity (e.g., Munguira et al. 1997; Singer et al. 2002; Junker et al. 2010; Korb et al. 2016). In Central and Western Europe, its main habitats are oligotrophic grasslands, and the most frequently used host plant is Succisa pratensis Moench (Dipsacaceae) (Warren 1994; Anthes et al. 2003; Konvička et al. 2003; Meister et al. 2015). The butterfly is monovoltine, with flight period in late spring/early summer when the mated females oviposit on leaf rosettes of the host plant. The larvae feed gregariously in silken webs on the plants until overwintering. They enter hibernation with the host plants’ senescence in mid-September and resume feeding solitarily in April. In the Czech Republic, the distribution is restricted to the western part of the country (Fig. 1), where it forms three distinct metapopulation clusters inhabiting ≈90 separate oligotrophic meadow patches interconnected by dispersal (Zimmermann et al. 2011; Junker et al. 2021; Tájek et al. 2023). This system is monitored annually by counting larval nests (cf. Ojanen et al. 2013) and displays remarkable within-site and inter-annual dynamics with booms and bursts (John et al. in rev.).

Figure 1. 

The distribution of Euphydryas aurinia in the Czech Republic, historic records included, based on Beneš et al. (2002), with actualisations. The inset in upper right corner shows the position of the country in Europe.

Like many other insects, E. aurinia hosts numerous hymenopteran parasitoids (Wahlberg et al. 2001; Eliasson and Shaw 2003; Stefanescu et al. 2009). The braconids of the genus Cotesia (Cameron, 1891), gregarious endoparasitoids of Lepidoptera, can be considered the most important and numerous. Their adult females oviposit into haemolymph of lepidopteran caterpillars; their larvae feed internally, break through the cuticle in the prepupal larval instar, and form silky external cocoons, in which they pupate, and from which the adult wasps hatch (Kester and Barbosa 1991; Pakarinen 2011). It was long believed that the main hymenopteran parasitoid of European Melitaeini butterflies was Cotesia melitaearum (Wilkinson, 1937), remarkable for its plurivoltine development, in which successive wasp broods oviposit on successive caterpillar instars (Shaw et al. 2009; Pakarinen 2011). A molecular approach complicated the matter by revealing that C. melitaearum is a complex of several cryptic species (Kankare et al. 2005c; Stefanescu et al. 2009). Regardless, the following hymenopteran parasitoids of E. aurinia have been so far reported from the Czech Republic: Braconidae – Cotesia melitaearum (Wilkinson, 1937), C. tibialis (Curtis, 1830); Pteromalidae – Pteromalus puparum (Linnaeus, 1758); Ichneumonidae – Ichneumon emancipatus (Wesmael, 1845), Ichneumon gracilicornis (Gravenhorst, 1849) (Shaw et al. 2009; Yu et al. 2012).

Diverse methods were used so far to study the parasitoids of E. aurinia, and related butterflies, ranging from field counts of hymenopteran cocoons (Ford and Ford 1930; Porter 1983), captive rearing (Eliasson and Shaw 2003), field experiments with captive-reared material (Stamp 1981) to population genetic studies targeting parasitoid adults (Lei and Hanski 1997; Van Nouhuys and Lei 2004). However, the question pivotal to the butterfly population dynamics and conservation, that of infestation rates relative to population cycle and state of the butterfly colonies, seems to be little explored. This is probably due to the work requirements for rearing both butterflies and parasitoids (cf. Klapwijk and Lewis 2014), combined with destructivity of such methods for field populations. To quantify the parasitism rates in the Czech Republic populations of E. aurinia, we developed a molecular method, allowing rapid and low-cost detection of Hymenoptera parasitoids’ incidence.

DNA-based methods are increasingly used in studies of parasitoid-hosts interactions (Zhu et al. 2019; Jeffs et al. 2021). The protocols so far developed for Lepidoptera/Hymenoptera systems mainly focused on COI locus, commonly referred as barcode (Folmer et al. 1994; Hebert et al. 2004), with the hope that the solid barcoding databases will assist species` identification (Toro-Delgado et al. 2022). Particularly good results were obtained via a reversal approach when the adult parasitoids were screened for host DNA shortly after their emergence (Rougerie et al. 2011) or in a species-poor natural system (high Arctic: Wirta et al. 2014). However, use of COI-based primers may be unreliable without subsequent sequencing, because deeply phylogenetically conserved bases are few and too far between in COI to place a group specific primer. Therefore, it was necessary to find a novel primer or primer pair which would amplify Hymenoptera but not Lepidoptera in the mixed samples containing the DNA of known lepidopteran host and unidentified hymenopteran parasitoids. We found such a potential primer in the nuclear region encoding the 28S ribosomal DNA. The novel primer, together with a primer published by Larsen (1992), targets part of the 28S gene and aims to amplify only Hymenoptera.

In this paper, we quantify Hymenoptera parasitoids infestation rates in a selection of the Czech Republic populations of Euphydryas aurinia and relate the infestation rates to the stage of the butterfly population cycle. Additionally, we document utility of our primers’ combination for rapid Hymenoptera infestation assessment in butterflies.

We sampled E. aurinia caterpillars in western Bohemia (Fig. 1) in late August and early September. We sampled two caterpillars per larval web (105 webs in total) from 13 sites in 2019 and four per web (90 webs in total) from 9 sites in 2020.

While sampling the caterpillars, we recorded the following: Julian date, to account for infestation changes during larval period; longest and shortest dimension of the larval web (cm); sward height, i.e., visually estimated height of surrounding vegetation in 2.5 metre radius circles around each larval web sampled; host plant density, expressed as the number of Succisa flowerheads in the circle; and webs density, expressed as the number of larval webs in the circle.

The material was stored in 96% ethanol, the DNA was extracted using the Tissue Genomic DNA Mini Kit (GenAid Biotech, Taiwan) following the manufacturer’s protocol.

We targeted a part of the 28S D1 region. We used primer 28SD1F (GGGGAGGAAAAGAAACTAAC; Larsen et al. 1992) in combination with a new primer HymR157 (TGGCCCCATTCAAGATGG) with a resulting product of 164–167 bp. For the primer design, we assembled a library of target sequences (Hymenoptera parasitoids) and non-target sequences (Lepidoptera) from sequences available in GenBank (primarily PopSet 300390962, Heraty et al. 2011). We aligned the sequences in GENEIOUS PRIME 2020.2.4 (https://www.geneious.com) software and used AMPLICON software (Jarman 2003) to identify sections with concentrated nucleotides that consistently differ between the target and non-target groups. We then manually screened these promising sections and used general rules of thumb to design candidate primers. We aimed for at least three differences between target and non-target groups in the first five positions at the 3’ end of the primer for reliable specificity. We then tested the candidate primer in an in-silico PCR in GENEIOUS and optimized melting temperature in the primer pairs by extending or shortening the primers.

Each in vitro PCR reaction contained 6.5 µl of Combi PPP Mastermix (Top-Bio, Czech Republic), 4.5 µl of H₂O, 0.5 µl of both reverse and forward primer, and 1 µl of DNA template. The cycling conditions of PCR were as follows: 94 °C of initial denaturation (5 mins), 30 cycles at 94 °C denaturation (40 secs), 50 °C annealing (30 secs), and 72 °C elongation (1 min), with the final elongation at 72 °C (5 mins). The presence/absence of PCR products was checked using agarose electrophoresis (1.5% gel, 150V, 30 mins). The samples with a band on the gel were assumed as positive; i.e., individuals infested by parasitoids (Fig. 2).

Figure 2. 

Electrophoresis gels used to assess whether the primers used can discriminate lepidopteran hosts and hymenopteran parasitoids a various adult hymenopteran parasitoids (PCRs are positive) and positive (P) and negative (N) samples of E. aurinia b four species of adult butterflies; PCRs are negative. The adult specimens of Hymenoptera and Lepidoptera were identified by M. Rindoš, M. Konvička, and Z. Faltýnek Fric.

To test the utility of the primer used, we also carried out control reactions which contained DNA extracted from various adult hymenopteran parasitoids and butterflies (Fig. 2) to assure that we amplified only the potential parasitoids and not the butterfly. Some of the obtained PCR products (n=4) of parasitoids from positive E. aurinia samples were sequenced in SEQme (Czech Republic) to confirm hymenopteran origin. We checked the identity of the obtained sequences by using BLAST (nBLAST algorithm) (https://blast.ncbi.nlm.nih.gov, Altschul et al. 1990).

Positive results were recalculated to infestation levels per larval web (1/0 factor, infested or not) and per site (% of larvae sampled). To relate the per-web infestation level to larval web properties, we carried out logistic regressions (binomial error distribution; in R 4.2.3, R Core Team 2018) with infestation (1/0) as the dependent variable; and the longest and shortest web ground projections, surrounding vegetation height, Succisa density, and number of larval webs as predictors. We used the information theory approach, comparing the fitted regression Akaike information criteria (AIC) with AIC of the null model, y~+1, and considered models with ΔAIC > ≈2.0 as fitting the data.

To relate the per site percentual infestation to larval counts at the sites, we used data from annual monitoring of the sites, ongoing since 2001 (Hula et al. 2004). Over this time, larval counts were obtained for roughly ¾ of site x year combinations (John et al. in rev., Suppl. material 1).

Out of 210 (year 2019) and 358 (year 2020) E. aurinia caterpillars assayed for hymenopteran DNA, we obtained 70 and 144 positive results, respectively; i.e., the total infestation rates were 33.3% and 40.2% per individuum. On a per-site basis, this translates to mean±SD / median / range 38.5±29.89 / 40 / 0–100 per cent in 2019, and 40.1±26.51 / 50 / 0–77.5 per cent in 2020.

The sequences of the positive samples were 124–126 bp long. The most similar sequences in GenBank according to nblast algorithm are those of Cotesia glomerata (Linnaeus, 1758), with the query identity 95.2–96.03%. The next similar sequences did not even reach a match of 93%.

According to the logistic regressions, none of the recorded properties of larval webs were related to infestation of the web (Table 1).

At the level of individual colonies, the infestation rates were highly variable (Fig. 3a). The per-site infestation levels did not correlate with larval web counts from the same year (Spearman’s r_s = 0.11, t_(17df) = 0.46, p = 0.65) or with web counts in the subsequent year (r_s = 0.31, t_(16df) = 1.32, p = 0.21), but did correlate positively with the larval webs’ counts from the previous year (i.e., 2018 web counts for 2019 sampling, and 2019 web counts for 2020 sampling: r_s = 0.46, t_(16df) = 2.07, p = 0.054) (Fig. 3b).

Figure 3. 

Per-site hymenopteran parasitoids infestation rates in colonies of the butterfly Euphydryas aurinia in two consecutive years (2019–20), with information of caterpillar web counts in the respective colonies in 2018–2021 (above), and illustration of the relationship between hymenopteran parasitoids infestation rates and E. aurinia caterpillar web counts in the previous year (below).

Because the absolute values of web counts highly varied among the sites (Fig. 3a), depending, e.g., on the site areas, we also recalculated web counts into percentage fractions of a ten-year (2011–2021) maximum for the given site, and recalculated the correlations, with no significant results (identical year percentage web counts: r_s = –0.23, t_(17df) = –0.99, p = 0.37; previous year percentage web counts: r_s = 0.31, t_(16df) = 1.31, p = 0.21, subsequent year percentage web counts: r_s = 0.10, t_(16df) = 0.40, p = 0.69).

The novel combination of primers 28SD1F (Larsen 1992) and HymR157 allowed us labour-efficient detection of high infestation rate in the declining E. aurinia butterfly by Hymenoptera parasitoids. Sequencing a selection of the obtained PCR products suggested that some of the parasitoids belong to the genus Cotesia Cameron, 1891. This is supported by the fact that other hymenopteran parasitoids known from our region, Ichneumon spp. and Pteromalus spp., attack pre-pupation larvae and pupae, respectively, whereas we worked with pre-diapause larvae. The gene 28S D1 is currently little represented in nucleotide databases for genus Cotesia, which, together with the unsettled taxonomy of Cotesia wasps infecting Melitaeinae butterflies (cf. Kankare et al. 2005a, b), precludes specific identification at this moment. This highlights the need to continue building comprehensive reference libraries for species identification. Such libraries are currently well developed for the COI gene, but supplementing it with a marker in another locus could improve discrimination power in some complex cases. Still, our approach allowed quantifying hymenopteran infestation rates in the declining butterfly, revealing that per-colony infestation rate is affected by larval webs density in the previous season.

Although the role of parasitoids on population fluctuations of E. aurinia, and related species, had been proposed almost a century ago (Ford and Ford 1930), relatively few authors quantified the natural infestation rates, with widely varying results. Infestation rates of <5% were reported, e.g., for the American congeneric species E. editha (Boisduval, 1852) (Singer and Erlich 1979) and E. chalcedona (Doubleday, 1847) (Lincoln et al. 1982). Similar results were obtained for E. aurinia from Sweden, with rates 2.6% (Eliasson and Shaw 2003), and Spain, where the rates were 2.4–5.1% (Stefanescu et al. 2009). The latter study in fact covered a newly recognised species, E. beckeri (Lederer, 1853), feeding on Lonicera spp. The rates detected by us, 33.3% and 40.2% in two consecutive years, are more comparable to the situations reported for American E. phaeton (Drury, 1773) (up to ≈10% in larvae prior to diapause) (Stamp 1981), E. maturna (Linnaeus, 1758) in the Czech Republic (69%) (Dolek et al. 2006) and Sweden (32%) (Eliasson and Shaw 2003). For E. aurinia, values higher than ours (≈90%) were reported by Ford and Ford (1930) from Britain during peaks of cyclic fluctuation of the butterfly population, whereas Klapwijk and Lewis (2014) reported a high range of infestation rates within individual webs (4–83%) from Britain. This all points to a high variation among sites, seasons, parasitoids’ generations, and Euphydryas species in hymenopteran infestation rates. In detailed studies of the closely related Melitaeinae model species, Melitaea cinxia (Linnaeus, 1758), this variation was attributed to spatial positions of butterfly colonies, competition among parasitoids and hyperparasitism, and annual variation in weather (e.g., Lei and Hanski 1997; Van Nouhuys and Lei 2004). Arguably, some of the reported variation may also be due to the diversity of methods applied by various authors, ranging from field counts of infested caterpillars (Lei and Hanski 1997), through captive rearing (Stamp 1981; Eliasson and Shaw 2003; Klapwijk and Lewis 2014), to molecular methods as used here. Possibly, the molecular detection reveals higher infestation rates than rearing, because some of the infested larvae may die prior to the parasitoid emergence. Causes of this mortality are then interpreted as “unknown” (e.g., in Eliasson and Shaw 2003).

Klapwijk and Lewis (2014) observed that the probability of caterpillar web infestation increased in webs isolated from other E. aurinia larval webs. We did not detect any relationship to the webs or surrounding vegetation parameters, but these results may be biased, because – for conservation concerns – we sampled the caterpillars solely from colonies containing a high number of larval webs during the sampling. The mean±SD webs’ number of sampled colonies were 116±97.3 and 60±66.6, whereas the numbers across all colonies were 53±57.8 (n = 44) and 19±38.0 (n = 70) in 2019 and 2020, respectively (John et al. in rev.). Conservation concerns also prevented quantifying the proportion of infestation per larval web, which would require killing all the caterpillars (cf. Klapwijk and Lewis 2014).

It also should be noted that our approach did not distinguish parasitoids from hyperparasitoids; i.e., the insects developing within parasitoids and thus killing them (Nair et al. 2016). A hyperparasitoid, however, can infest only a larva already infested by a parasitoid, and hence hyperparasitoids presence does not affect our findings on hymenopteran infestation rates.

With all the limitations, we found that the infestation rate per site positively correlated with per-site caterpillar webs’ numbers of the previous year. This is fully expected if the parasitoids need a rich resource supply (i.e., high host density) to multiply in a butterfly colony, depleting the hosts’ numbers in the process (Ford and Ford 1930; Frazer 1954; Porter 1981, 1983). A time delay in parasitoids infestation, and higher likelihood of infestation of larger and more connected host colonies, were found by Lei and Hanski (1997) in the metapopulation system of M. cinxia and its parasitoids. Although the inter-annual abundance changes of Melitaeini butterflies’ colonies are likely influenced by numerous other factors, including variation in weather (Brunbjerg et al. 2017) or site vegetation management (Johansson et al. 2019; Tájek et al. 2023), natural enemies’ pressure certainly plays a significant role.

The DNA-based method detected high hymenopteran parasitoids’ infestation rates in colonies of the declining butterfly, Euphydryas aurinia. These rates, however, widely varied among the butterfly colonies and between two study years, likely interfering with, and possibly driving, the inter-annual butterfly abundance changes within colonies, as well as the metapopulation dynamics of the butterfly, described in detail by John et al. (in review). Our primer combination seems to be promising for wider use in detecting infestation of butterflies by Hymenoptera parasitoids. In our control tests, it amplified various Hymenoptera but not several butterfly species tested (Fig. 2). However, potential users should first test that the detection works also in their system to avoid false positive and false negative results.

Part of the development of the primers took place during the stay of JH in the lab of Graham Stone, University of Edinburgh. The work was funded by the Technology Agency of the Czech Republic (SS01010526). JH was supported by Czech Science Foundation grant no. 20-30690S.