Wolbachia in the spittlebug Prosapia ignipectus: Variable infection frequencies, but no apparent effect on host reproductive isolation

Abstract Animals serve as hosts for complex communities of microorganisms, including endosymbionts that live inside their cells. Wolbachia bacteria are perhaps the most common endosymbionts, manipulating host reproduction to propagate. Many Wolbachia cause cytoplasmic incompatibility (CI), which results in reduced egg hatch when uninfected females mate with infected males. Wolbachia that cause intense CI spread to high and relatively stable frequencies, while strains that cause weak or no CI tend to persist at intermediate, often variable, frequencies. Wolbachia could also contribute to host reproductive isolation (RI), although current support for such contributions is limited to a few systems. To test for Wolbachia frequency variation and effects on host RI, we sampled several local Prosapia ignipectus (Fitch) (Hemiptera: Cercopidae) spittlebug populations in the northeastern United States over two years, including closely juxtaposed Maine populations with different monomorphic color forms, “black” and “lined.” We discovered a group‐B Wolbachia (wPig) infecting P. ignipectus that diverged from group‐A Wolbachia—like model wMel and wRi strains in Drosophila—6 to 46 MYA. Populations of the sister species Prosapia bicincta (Say) from Hawaii and Florida are uninfected, suggesting that P. ignipectus acquired wPig after their initial divergence. wPig frequencies were generally high and variable among sites and between years. While phenotyping wPig effects on host reproduction is not currently feasible, the wPig genome contains three divergent sets of CI loci, consistent with high wPig frequencies. Finally, Maine monomorphic black and monomorphic lined populations of P. ignipectus share both wPig and mtDNA haplotypes, implying no apparent effect of wPig on the maintenance of this morphological contact zone. We hypothesize P. ignipectus acquired wPig horizontally as observed for many Drosophila species, and that significant CI and variable transmission produce high but variable wPig frequencies.

Many Wolbachia manipulate host reproduction to propagate in host populations. For example, many strains cause cytoplasmic incompatibility (CI) that reduces the egg hatch of uninfected embryos fertilized by Wolbachia-infected sperm (Hoffmann & Turelli, 1997).
Wolbachia contribute to assortative mating and postzygotic isolation between co-occurring D. paulistorum semispecies (Miller et al., 2010), and to reinforcement of isolation between uninfected D. subquinaria and Wolbachia-infected D. recens (Jaenike et al., 2006;Shoemaker et al., 1999). In contrast, Wolbachia do not contribute to RI in the D. yakuba clade, which includes wYak-infected D. yakuba, wSan-infected D. santomea, and wTei-infected D. teissieri (Cooper et al., 2017). Thus, while some results from Drosophila strongly support contributions of Wolbachia to RI, and interest in the possibility of such effects remains high, it is unknown whether Wolbachia effects on RI are common in nature.
Prosapia ignipectus (Fitch) (Hemiptera: Cercopidae) is one of about 14 species of Prosapia and one of two commonly found in the United States, the other being its sister species P. bicincta (Say) (Hamilton, 1977). P. ignipectus occurs in southern Ontario, Canada, and the northeastern United States from Minnesota to Maine (Carvalho & Webb, 2005;Hamilton, 1977Hamilton, , 1982Peck, 1999;Thompson & Carvalho, 2016). These species vary in male genital morphology and in associations with host plants, with P. ignipectus monophagous on the late season C4 perennial grass Schizachyrium scoparium (Little bluestem) (Hamilton, 1982;Thompson, 2004) and P. bicincta polyphagous on a variety of C4 grasses, but not including Little bluestem (Fagan & Kuitert, 1969;Thompson, 2004). Both species have conspicuous dorsal coloration, standing out against their respective host plants. All P. bicincta individuals have a single narrow transverse orange line across the widest part of the pronotum and a pair of narrow orange lines across the elytra. Most P. ignipectus individuals have a solid black dorsal surface, but in Maine some P. ignipectus have P. bicincta-like coloration ( Figure 1). Notably, only 10 km separate monomorphic black and monomorphic lined P. ignipectus populations in western Maine, with little evidence of a hybrid zone and no obvious physical barriers to mixing across the boundary (Thompson & Carvalho, 2016). This morphological contact zone has persisted for at least 90 years. About 45 km southwest of this abrupt transition between aposematic color forms, three other P. ignipectus populations were found to be polymorphic with both black and lined forms-these populations are surrounded by monomorphic black populations. It has been hypothesized that Wolbachia may contribute to host RI and to preservation of the sharp Maine morphological contact zone (Thompson & Carvalho, 2016).
Here, we use collections of P. ignipectus from several sites in the northeastern United States across two years, in combination with collections of P. bicincta from Hawaii and Florida, United States, to assess modes of Wolbachia acquisition and to test for Wolbachia frequency variation through space and time. By sampling monomorphic black and lined populations and typing both Wolbachia and mtDNA haplotypes, we also test for contributions of Wolbachia to the P. ignipectus morphological contact zone. Finally, we generate whole genome Wolbachia data for phylogenetic analysis and to search for loci associated with inducing and rescuing CI LePage et al., 2017;Shropshire et al., 2018). While we cannot currently test P. ignipectus for CI in the laboratory, CI-causing Wolbachia are predicted to occur at high infection frequencies and to have specific loci associated with CI in their genomes.

| Sampling
We netted specimens from Little bluestem; sorted them by species, sex, and color form; and preserved them in 95% ethanol. The 2019 specimens (N = 4 sites) were collected on August 23. The 2020 specimens (N = 9 sites) were collected on August 9 (Silver Lake, NH), August 17 (Wonalancet, NH), and August 20 (all Maine localities) (Table S1). Collection sites were on the verges of public rights of way or privately owned land. In two cases (New Vineyard and New Portland), they correspond to sites reported in Thompson and Carvalho (2016). Specimens were collected near the height of abundance for P. ignipectus, which starts to emerge in adult form in late July and early August. We also sampled three additional spittlebug species at these sites: Lepyronia quadrangularis (Say) (N = 25), Philaenus spumarius (L.) (N = 5), and Philaenarcys killa (Hamilton) (N = 24), all of the family Aphrophoridae. Like, P. ignipectus, P. killa a is monophage on Little bluestem. L. quadrangularis is a polyphage but often abundant on Little bluestem. P. spumarius is an extreme polyphage, with a preference for forbs (herbaceous perennial dicots) but is occasionally collected from Little bluestem in the company of P. ignipectus. By screening them for Wolbachia, we tested for the possibility of horizontal Wolbachia transfer through plant interactions (Chrostek et al., 2017). Lastly, because identification of infections in sister hosts enables formal analysis of modes of Wolbachia acquisition (Conner et al., 2017;Cooper et al., 2019;Raychoudhury et al., 2009;Turelli et al., 2018), we also obtained samples of the sister species P. bicincta from Hawaii (N = 60) and Florida (N = 40) to screen for infections. P. bicincta is native to the southeastern United States (Fagan & Kuitert, 1969;Thompson & Carvalho, 2016), but has recently been introduced into the Kona Region of Big Island, Hawaii (Thorne et al., 2018).

| Wolbachia typing
We generated whole-genome Wolbachia data to type the Wolbachia infecting P. ignipectus and to search for loci associated with CI. We extracted 800ng of high molecular weight DNA (Qiagen Genomic-tip 20/G; Qiagen, Germany) from one black New Vineyard female (see below), and then input and sequenced it (Ligation Sequencing Kit, SQK-LSK109; FLO-MIN106 flow cell) for 48 hr (Oxford Nanopore Technologies). We mapped raw nanopore reads (5.8Gb of data) to all known Wolbachia sequences (NCBI taxid 953) with BLASTn and extracted reads where at least 60% of their length mapped (qcovs ≥ 60).

strains using
Prokka v.1.11 (Seemann, 2014). We used only genes present in single copy and with identical lengths in all genomes. To assess the quality of our assembly, we excluded wPig and repeated this with only wMel, wRi, wPip, and wMau.
Preliminary analysis of a few loci placed the P. ignipectus Wolbachia in group-B (see below), but we performed Bayesian analyses using the GTR + Γ + I model for sequence evolution using whole-genome data to confirm this (Höhna et al., 2016). Genes were concatenated and partitioned by codon position, with a rate multiplier, σ, assigned to each partition to accommodate variable substitution rates. We used flat, symmetrical Dirichlet priors on the stationary base frequencies, π, and the relative rate parameters, η, of F I G U R E 1 Sister species Prosapia ignipectus and Prosapia bicincta have conspicuous dorsal coloration. All P. bicincta individuals have a single narrow transverse orange line across the widest part of the pronotum and a pair of narrow orange lines across the elytra. Most P. ignipectus individuals have a solid black dorsal surface, but in Maine some P. ignipectus have P. bicincta-like coloration. P. ignipectus is monophagous on the late season C4 perennial grass Schizachyrium scoparium (Little bluestem). Little bluestem photo by Krzysztof Ziarnek, Kenraiz (CC BY-SA 4.0, https://creat iveco mmons.org/licen ses/by-sa/4.0) the GTR model (i.e., Dirichlet(1,1,1…)). As in Turelli et al. (2018), we used a Γ(2,1) hyperprior on the shape parameter, α, of the discrete-Γ model (adopting the conventional assumption that the β rate parameter equals α, so that the mean rate is 1 (Yang, 1994). The Γ model for rate variation assigns significant probability near zero when the α < 1 (accommodating invariant sites). The Γ(2,1) hyperprior on α assigns 95% probability to the interval (0.36, 4.74), allowing for small and large values. Four independent runs for each gene set produced concordant topologies. We diagnosed MCMC performance using Tracer 1.7 (Rambaut et al., 2014).

| Wolbachia and mtDNA haplotyping of black and lined color morphs
To confirm that the same Wolbachia strain infects different P. ignipectus populations and color morphs, we amplified and Sanger sequenced five protein-coding Wolbachia genes (coxA, hcpA, fbpA, ftsZ, and wsp) in both directions (Eurofins Genomics LLC, Louisville, Kentucky)(see below, Table S2). We also amplified and Sanger sequenced gatb, but sequence quality was consistently too low to include in our analyses. Samples included one infected female of each color form (black or lined), from each of the four populations (Carthage, New Portland, New Vineyard, and Strong) sampled in both years (Table S1).
To specifically assess whether Wolbachia might contribute to the morphological contact zone between New Vineyard (monomorphic black) and New Portland (monomorphic lined) P. ignipectus, we amplified and Sanger sequenced the cytochrome C oxidase I (CoI) mitochondrial locus from one male and one female from these populations, with the exception of one (New Vineyard black male) that did not produce a usable sequence. We also produced CoI sequences for one black and one lined female from the polymorphic Strong population. Covariance of Wolbachia and mtDNA haplotypes with P. ignipectus color forms would support a potential role for Wolbachia in maintaining the morphological contact zone.
We visually inspected each sequence for quality and ambiguities, and consensus sequences were used as queries for a BLASTn search and the NCBI "nr" database to confirm that orthologous genes were amplified (Altschul et al., 1990). We then used the "multiple locus query" function of the multi locus sequence typing (MLST) database to type Wolbachia (Baldo et al., 2006). Together, these data enable us to test for differentiation in Wolbachia and mtDNA between populations and color forms, including between populations monomorphic for different color forms separated by only 10 km in Maine.  Shropshire et al., 2018). Two genes (cifA/B) transgenically expressed in male D. melanogaster induce CI, while one gene (cifA) expressed in females rescues it. To identify cif loci, we used BLASTn to search for cif homologs in our whole-genome raw reads, querying the Type 1 cif pair in wMel, the Type 2 pair in wRi, the Type 3 pair in wNo, the Type 4 pair in wPip, and the Type 5 pair in wStri (Bing et al., 2020;Lindsey et al., 2018;Martinez et al., 2020). We later broadened our search for Type 1 pairs by querying wPip and wNPa pairs (Gerth & Bleidorn, 2017;Klasson et al., 2008). For each Type, we extracted raw reads that covered at least 40% of the genes. We then corrected and assembled the reads with canu 2.1.1 (Koren et al., 2017(Koren et al., , 2018Nurk et al., 2020), producing sequences with about a 1% error rate.

| Analysis of CI loci
We limit our analyses to the discovery of cif types, since we did not generate additional sequence data to further correct the long reads. The assembled genes were compared to those in Martinez et al. (2020).

| Analysis of Wolbachia frequency variation
To test for Wolbachia frequency variation, we extracted DNA from many individuals from each collection using a standard squish buffer protocol and identified Wolbachia infections using polymerase chain reaction (PCR) (Simpliamp ThermoCycler; Applied Biosystems, Singapore) (Meany et al., 2019). We amplified the Wolbachia surface protein (wsp) (Braig et al., 1998) and arthropodspecific 28S rDNA, which served as a positive control (Baldo et al., 2006) (Table S2). PCR products were visualized using 1% agarose gels. Assuming a binomial distribution, we estimated exact 95% confidence intervals for Wolbachia frequencies for each collection. We used Fisher's exact test (FET) to determine differences in frequencies among sites, between years, between sexes, and between color forms. P. bicincta populations, in combination with generally high wPig frequencies in P. ignipectus, indicate that P. ignipectus likely acquired wPig after its initial divergence from P. bicincta, although it is also possible that P. bicincta lost its resident Wolbachia following cladogenic acquisition. Because testing predictions about modes of Wolbachia acquisition requires formal analysis of Wolbachia, host nuclear, and host mtDNA phylograms and chronograms, we are unable to distinguish between introgressive and horizontal wPig transfer (Conner et al., 2017;Cooper et al., 2019;Gerth & Bleidorn, 2017;Raychoudhury et al., 2009;Turelli et al., 2018). We discuss this further below. Because Wolbachia that infect P. spumarius and wPig in P. ignipectus are both at high frequency, we also typed the Wolbachia infecting P. spumarius to determine whether a wPig-like variant infects this host species. The multiple sequence query in the MLST database supports that a different group-B strain, most closely related to the thrip species Aptinothrips rufus (ID 1945, ST 509) infects P. spumarius. Generating more sequence data will be required to resolve the phylogenetic relationships of these and other group-B strains, including Wolbachia in P. spumarius (Lis et al., 2015).

| No apparent effect of wPig on the maintenance of the morphological P. ignipectus contact zone
The Strong, Carthage, and Dixfield P. ignipectus populations Our mtDNA haplotypes are also very similar to ten P. ignipectus samples included in the Barcode of Life Database (BOLD) (Foottit et al., 2014). A single base-pair insertion present in all of our samples is absent from all ten BOLD samples. Four other sites in CoI that are polymorphic among the BOLD samples are fixed in our samples for one of the BOLD alleles. mtDNA haplotypes of P. ignipectus and P. bicincta also differ by <2% (Foottit et al., 2014).

| The wPig genome contains three divergent types of CI loci
We identified Type 1, 3, and 4 cifs in the wPig genome (Martinez et al., 2020). This specific complement of cifs is not found in any other published Wolbachia genomes, but close relatives to each wPig cif Type are. For instance, the wPig Type 1 genes are 99% identical to those in the genome of the Wolbachia infecting the gall-inducing wasp Diplolepis spinosa (Cynipidae), but less than 90% similar to any others (Martinez et al., 2020 None of the wPig cifs are truncated relative to copies with 99% identity. Additional sequencing is required to make more detailed cif comparisons.

| Pervasive wPig frequency variation
wPig varied in frequency in several ways. First, frequency varied spatially among all samples (FET, P = 0.001) (Table 1) . We interpret these results as pervasive spatial, and rare temporal and sex-specific, variation in wPig frequency. F I G U R E 3 wPig frequency varies through space and time. Circle size denotes sample size, with outline and fill color denoting sampling year and infection status, respectively. Sample means and 95% binomial confidence intervals are reported for each sample. The dashed back line denotes the geographical separation of monomorphic black and monomorphic lined Prosapia ignipectus populations

| D ISCUSS I ON
Our results suggest that wPig is a group-B Wolbachia acquired after the initial divergence of P. ignipectus from P. bicincta. Analysis of Wolbachia and mtDNA haplotypes indicates that wPig has no apparent effect on the P. ignipectus morphological contact zone in Maine. Across all samples, wPig occurs at very high frequencies, consistent with our discovery of three divergent sets of CI loci in the wPig genome. Finally, we document pervasive spatial, and rare temporal, wPig frequency variation. We discuss this in more detail below.

| Wolbachia acquisition in spittlebugs
In contrast to very high wPig frequencies in P. ignipectus, we found no evidence of Wolbachia in our sample of 100 P. bicincta. A prior report of one infected P. bicincta sample indicates that Wolbachia could infect this species (Anderson et al., 2019). If so, it must be at very low frequencies, given our credible interval here (p = 0.0 [0.0, 0.04]; N = 100). Mathematical models predict that intense CI drives Wolbachia to high frequencies, balanced by imperfect maternal transmission (Hoffmann et al., 1990;Turelli & Hoffmann, 1995); conversely, Wolbachia that do not cause strong CI tend to occur at much lower frequencies (Cooper et al., 2017;Hague, Mavengere, et al., 2020;Hamm et al., 2014;Kriesner et al., 2016). While crossing to test for CI in the laboratory is not currently feasible in this system, the presence of three sets of CI loci in the wPig genome, combined with its very high frequencies, suggests that wPig causes intense CI.
How did P. ignipectus acquire wPig? There are three possibilities: cladogenic transmission from its most recent common ancestor with its sister species, presumably P. bicincta or a close relative; by introgression from P. bicincta or another close relative; or by horizontal transmission (O'Neill et al., 1992). Given that we find no evidence for a high frequency Wolbachia in P. bicincta, cladogenic acquisition seems implausible, although we cannot fully rule it out. Without more extensive analysis of close relatives, we also cannot rule out introgression. However, opportunities for introgression with species other than P. bicincta have likely been limited. Other species of the genus Prosapia or family Cercopidae occur no further north than the US-Mexico border region, about 1,400 km from the nearest P. ignipectus populations and 3,000 km from the populations studied here.
Overall, the limited data are consistent with relatively recent noncladogenic transmission, a process that seems to be common among Drosophila species (Turelli et al., 2018). It may also be common among spittlebugs. This would be in stark contrast to obligate transovarial endosymbionts associated with amino acid nutrition in spittlebugs and other hemipterans (Koga et al., 2013).
In addition to the thrip-related Wolbachia found in P. spumarius in this study, Nakabachi et al. (2020) report that two spittlebug species, Aphrophora quadrinotata Say and Philaenus maghresignus Drosopoulos & Remane (both Aphrophoridae), harbor Wolbachia with 16S rRNA sequence that is identical to Wolbachia in two psyllid species, two whiteflies, an aphid, a planthopper, two leafhoppers, two grasshoppers, a mosquito, and a weevil. Likewise, Lis et al. (2015) report that Wolbachia they studied in P. spumarius is closely related to strains in vespids, drosophilids, whiteflies, chrysomelid beetles, and weevils based on five MLST loci. Kapantaidaki et al. (2021) also report Wolbachia infections at low levels in P. spumarius, as well as higher frequencies in Neophilaenus campestris (Fallén) (Aphrophoridae). Based on five MLST loci, their N. campestris strain is closely related to Wolbachia found in a leafhopper (Hemiptera) and cluster with Wolbachia from a planthopper, a scale insect and a psyllid (all Hemiptera), as well as two chrysomelid beetles, two butterflies, a parasitic wasp, and a mosquito. Koga et al. (2013, Table S2) report the presence of Wolbachia in the spittlebug Cosmoscarta heros (F.) (Cercopidae), in addition to A. quadrinotata and P. maghresignus.
In contrast, five specimens of Poophilus costalis (Walker) (Aphrophoridae) (Wiwatanaratanabutr, 2015), six specimens of Philaenus tesselatus Melichar (Lis et al., 2015), 37 specimens of Philaenus signatus Melichar (Kapantaidaki et al., 2021;Lis et al., 2015),  (Lis et al., 2015) were not infected. Based on limited sequence data, the emerging pattern suggests that Wolbachia infection is widespread, but far from ubiquitous among spittlebugs, and that when it does occur, it often involves Wolbachia strains similar to those infecting distantly related insects. Whole Wolbachia and host genomic data are sorely needed to test our hypothesis that horizontal Wolbachia acquisition might be common in spittlebugs.

| Little contribution of wPig to the P. ignipectus morphological contact zone
We find no evidence for differentiation in wPig or mtDNA haplotypes among P. ignipectus color forms. This includes the monomorphic black (New Vineyard) and lined (New Portland) populations that are separated by only 10 km in Maine, with no obvious barriers to dispersal or reproduction (Thompson & Carvalho, 2016). We also found no variation in wPig or mtDNA haplotypes between black and lined individuals in the polymorphic Strong population. wPig frequency also did not vary between color forms. These data indicate that wPig is unlikely to significantly contribute to the maintenance of the P. ignipectus morphological contact zone.
How common are Wolbachia effects on host RI? Obligate Wolbachia infections in co-occurring D. paulistorum semispecies contribute to assortative mating and generate hybrid inviability and male sterility (Miller et al., 2010). Wolbachia also contribute to reinforcement between Wolbachia-infected D. recens and uninfected D. subquinaria (Jaenike et al., 2006;Shoemaker et al., 1999). In contrast, Wolbachia do not contribute to premating, gametic, or postzygotic RI among the three D. yakuba-clade host species (Cooper et al., 2017). While the crossing schemes used in these Drosophila studies to dissect Wolbachia contributions to RI are not feasible in P. ignipectus and many other systems, our genetic data here lend support to our prior conjecture that Wolbachia contributions to RI observed in some Drosophila may be the exception rather than the rule (Cooper et al., 2017;Turelli et al., 2014).
With the exception of model systems like wRi in D. simulans, few estimates of the key parameters required to approximate population frequency dynamics and equilibria of Wolbachia exist (Carrington et al., 2011;Turelli & Hoffmann, 1995). wMel-like Wolbachia frequencies in the D. yakuba clade vary through space and time in west Africa (Cooper et al., 2017), due in part to effects of cold temperatures on wYak titer (Hague, Mavengere, et al., 2020). CI strength also varies in the D. yakuba clade, which may influence infection frequencies (Cooper et al., 2017;Hague, Caldwell, et al., 2020). wMel frequencies vary with latitude in D. melanogaster populations, potentially due to wMel fitness costs in the cold (Kriesner et al., 2016).
Interestingly, hot temperatures reduce wMel CI strength and transmission in transinfected Aedes aegypti used for biocontrol of human disease (Ross et al., 2017(Ross et al., , 2020, suggesting that temperature may generally influence key parameters underlying Wolbachia infection frequencies.
What underlies variable wPig frequencies in nature? High wPig frequencies and the presence of three divergent sets of cifs suggest, but do not confirm, that wPig causes strong CI. It seems plausible that some or all of these loci were horizontally acquired , but additional sequence data are required to test this.
We hypothesize that variable wPig transmission rates contribute to the frequency variation we observe, potentially due to environmental effects on titer, as observed for wYak (Hague, Mavengere, et al., 2020). Temporal variation in transmission was also observed for wRi between two samples of D. simulans collected from Ivanhoe, California, in April and November of 1993 (Carrington et al., 2011;Turelli & Hoffmann, 1995), although the relative stability of wRi frequencies in global D. simulans populations suggests that its transmission persists across a range of environmental conditions. Additional analyses of Wolbachia titer and transmission in the field, and across environmental contexts, are needed to better understand the causes of Wolbachia frequency variation. Comparing the titer and transmission of Wolbachia that occur at different frequencies in nature-for example, those that do and do not cause intense CI-would be particularly useful.

ACK N OWLED G M ENTS
We thank M. Thorne

CO N FLI C T O F I NTE R E S T
We declare no conflicts of interests.