A decade of stability for wMel Wolbachia in natural Aedes aegypti populations

Mosquitoes carrying Wolbachia endosymbionts are being released in many countries for arbovirus control. The wMel strain of Wolbachia blocks Aedes-borne virus transmission and can spread throughout mosquito populations by inducing cytoplasmic incompatibility. Aedes aegypti mosquitoes carrying wMel were first released into the field in Cairns, Australia, over a decade ago, and with wider releases have resulted in the near elimination of local dengue transmission. The long-term stability of Wolbachia effects is critical for ongoing disease suppression, requiring tracking of phenotypic and genomic changes in Wolbachia infections following releases. We used a combination of field surveys, phenotypic assessments, and Wolbachia genome sequencing to show that wMel has remained stable in its effects for up to a decade in Australian Ae. aegypti populations. Phenotypic comparisons of wMel-infected and uninfected mosquitoes from near-field and long-term laboratory populations suggest limited changes in the effects of wMel on mosquito fitness. Treating mosquitoes with antibiotics used to cure the wMel infection had limited effects on fitness in the next generation, supporting the use of tetracycline for generating uninfected mosquitoes without off-target effects. wMel has a temporally stable within-host density and continues to induce complete cytoplasmic incompatibility. A comparison of wMel genomes from pre-release (2010) and nine years post-release (2020) populations show few genomic differences and little divergence between release locations, consistent with the lack of phenotypic changes. These results indicate that releases of Wolbachia-infected mosquitoes for population replacement are likely to be effective for many years, but ongoing monitoring remains important to track potential evolutionary changes.


Introduction
Open field releases of Wolbachia-infected mosquitoes are becoming one of the best ways to control arbovirus transmission. Wolbachia "population replacement" programs involve the release of mosquitoes carrying a Wolbachia infection that spreads through mosquito populations and reduces their vector competence [1][2][3]. Several different Wolbachia strains from other insects have been introduced into Ae. aegypti mosquitoes through microinjection, with many of them reducing the ability of mosquitoes to transmit arboviruses including dengue, Zika and chikungunya [4][5][6]. The establishment of wMel and wAlbB Wolbachia strains in natural populations has suppressed arbovirus transmission by Ae. aegypti in release locations time after being transferred to a novel Drosophila host [30,31]. Selection experiments show that shifts in the phenotypic effects of Wolbachia are often due to host genetic changes [32]. So far there have been limited changes observed in Wolbachia genomes across several years after transinfection to novel hosts such as wCer2 in Drosophila [33] and wMelPop-CLA in Ae. aegypti [34].
In mosquitoes, there is a clear distinction between natural Wolbachia infections and novel transinfections, where the latter tends to induce deleterious effects [1]. This suggests that the effects of Wolbachia transinfections may weaken across evolutionary timescales. Since population replacement programs began only a decade ago, there is still limited information on evolutionary changes following deliberate Wolbachia invasions. In laboratory populations of Ae. aegypti, the wMelPop-PGYP strain has continued to induce complete cytoplasmic incompatibility for the last decade, but some fitness costs appear to have weakened [35]. In field populations, the wMel strain has shown stable phenotypic effects, with no evidence for changes in dengue virus blocking [21,36], cytoplasmic incompatibility [15,37] or effects on fertility [15] after a period of 1 year or more under field conditions. Whole Ae. aegypti [38] and Wolbachia [39,40] genome sequencing studies show limited genomic changes after wMel has been established in Cairns for at least 7 years.
Given the importance of tracking the long-term stability of Wolbachia infections in mosquito populations, we have now collected additional data on the phenotypic and genomic stability of wMel from the first ever population replacement releases in Cairns, Australia in 2011. Our data from up to a decade post-release show that the phenotypic effects of wMel largely remain stable in both laboratory and near-field Ae. aegypti populations. We also extend the findings of Huang et al. [39] to show limited genomic changes in Wolbachia over the span of a decade, with no divergence in wMel genomes between different release locations. Our data point to the likely long-term effectiveness of Wolbachia population replacement programs globally.

wMel remains at a high frequency in Cairns following releases
We performed ovitrapping across suburbs in Cairns in 2016 and 2018 to monitor Wolbachia infection frequencies in the Ae. aegypti population. In 2016, all suburbs where wMel releases had taken place had a wMel infection frequency above 0.96 (S1A Fig). Some pre-release suburbs had low infection frequencies (0.05-0.21), indicating spread of wMel to adjacent suburbs (e.g. to Mt Sheridan and Holloways Beach). By 2018, releases had occurred in additional locations (e.g. Trinity Beach and Redlynch) and all suburbs had wMel infection frequencies greater than 0.93 (S1B Fig), except for Redlynch where wMel was not released until 2019. Our data are broadly consistent with Wolbachia infection frequencies from an independent study which also shows that wMel has maintained high frequencies in most release locations [16].

wMel density remains stable across suburbs and laboratory generations
We measured wMel density in 4 th instar larvae reared in the laboratory from ovitraps collected across suburbs in 2018. We found no clear effect of release year (GLM: F 4,269 = 2.274, P = 0.062) or suburb (nested within release year) (F 6,269 = 2.017, P = 0.064) on Wolbachia density ( Fig 1A). These results suggest that the wMel infection has remained stable after being established in the field for different periods of time (from 1-7 years), with no clear effects of local environmental conditions on wMel in the next generation of mosquitoes.
In June 2019, we measured wMel density in adults from laboratory populations that were established from field collections at different times. We found no significant differences PLOS PATHOGENS wMel stability in Aedes aegypti between populations for females (GLM: F 4,65 = 1.767, P = 0.146, Fig 1B) or males (F 4,70 = 1.942, P = 0.113, Fig 1C), suggesting that whole-adult wMel density has not changed across different durations of laboratory rearing.

Population origin and wMel infection influence mosquito fitness
In September 2018 we performed phenotypic comparisons between long-term laboratory and near-field populations that were wMel-infected or cleared of wMel through antibiotic treatment (Fig 2). We found significant effects following Bonferroni correction (adjusted α: 0.008) of population origin on larval development time (females: F 1,44 = 19.656, P < 0.001, males:

No clear effect of wMel origin in a common host background
In the first experiment, we found that effects on fitness were driven mainly by population origin. Due to potential interactions between Wolbachia infection and background in the first experiment, we performed a second set of experiments in August 2020 to evaluate the effects of lab and field-derived wMel infections in a common host background. We found no significant differences following Bonferroni correction (adjusted α: 0.006) between populations for most traits tested, including development time (GLM: females: F 3,43 = 2.841, P = 0.049, males: F 3,43 = 3.558, P = 0.022), female wing length (F 3,74 = 1.826, P = 0.150), fecundity (F 3,433 = 1.950, P = 0.121) and adult longevity (Log-rank: females: χ 2 = 5.700, P = 0.127, males: χ 2 =

Limited effects of antibiotic treatment on fitness
To test whether tetracycline treatment could influence mosquito fitness, we performed phenotypic comparisons of wMel-infected and uninfected mosquitoes following parental tetracycline treatment (Fig 4). We found no significant effect of tetracycline treatment on any trait (all P > 0.007, adjusted α: 0.006) except for (logit transformed) egg hatch proportion where we found an interaction between treatment and Wolbachia infection type (F 1,70 = 9.261, P = 0.003). Decreased egg hatch due to tetracycline treatment was apparent only for the wMel-infected population, likely due to incomplete curing resulting in partial self-incompatibility. In this experiment, we also found no significant effect of Wolbachia infection type for any trait (all P > 0.017, adjusted α: 0.006) except for male wing length, where wMel-infected males had larger wings (F 1,69 = 12.387, P = 0.001). These results provide confidence that observed effects in comparisons between Wolbachia-infected and cured mosquitoes are due to removal of the Wolbachia infection and not off-target tetracycline-related effects on the gut microbiome or mitochondria.

wMel induces complete cytoplasmic incompatibility after eight years in the field
In June 2019, we tested the ability of wMel-infected males from different origins (Yorkeys Knob, Gordonvale and laboratory) to induce cytoplasmic incompatibility with uninfected females. Males from all wMel-infected populations induced complete cytoplasmic incompatibility, with no eggs hatching in crosses with uninfected females (Table 1). All other crosses yielded egg hatch proportions above 90% (Table 1), consistent with data from the previous experiments showing limited effects of wMel infection on egg hatch. These results show that wMel has retained complete cytoplasmic incompatibility after > 8 years in the field and one generation of laboratory rearing.

Limited changes in wMel and mitochondrial genomes across a decade in a novel Aedes aegypti host
Previous wMel genome comparisons suggest limited changes in the Wolbachia genome have occurred since field releases [39,40]. We extend this timeframe by including pre-release material from 2010 and field collections from Gordonvale and Yorkeys Knob in 2020. We also include a dataset from a previous whole genome sequencing study of individuals collected from Gordonvale in 2013 and 2018 and Yorkey's Knob in 2018 [38]. Little evidence of change was observed among the analysed wMel genomes, for which we achieved � 99.99% coverage. Within populations, each of which was represented by a sample of pooled individuals, seven loci were found to have one or more alternate alleles present at a frequency greater than our threshold value of 25% in at least one sample, all of which were SNPs or small indels. At all other positions within the genome, the wC45 F 10 (pre-release) consensus sequence was identical to that of the Gordonvale and Yorkeys Knob field collections from 2011, 2018 and 2020, and the wMel lab population (see also Huang et al. [39] and Dainty et al. [40]). Relative to these genomes, one single nucleotide deletion was observed in the wC45 F 9 pre-release genome sequence and seven single nucleotide indels were detected in the 2013 Gordonvale field collection genome sequence (S1 Table). It is probable that these eight indels represent sequencing artifacts, as they involve non-parsimonious changes in populations sampled over multiple time points, were often associated with a marked localized decrease in read depth in the variant samples and mostly represent frameshift mutations.
Of the seven loci that displayed within-sample polymorphism, two were within non-coding regions; three were within three different copies of the group II intron-encoded protein gene LtrA, each at the same position within the gene and corresponding to an Asp > Asn substitution; one was within a gene encoding a hypothetical protein and corresponds to a silent change; and one was within a tRNA-Arg gene (S2 Table). It is possible that the apparent polymorphism within the LtrA genes is due to sequence variation between gene copies, rather than the presence of alternate alleles. Multiple alleles were observed in all samples for these three loci and for the polymorphic locus within the tRNA-Arg gene. At the other three loci, allele frequencies were more variable, with only some samples displaying polymorphism, most of which were from the field populations at Gordonvale and Yorkeys Knob. Two of these loci, at positions 587,862 and 1,174,712 relative to the wMel reference genome, were previously noted by Huang et al. [39] and Dainty et al. [40]. No clear pattern was observed across different populations or time points. While it is possible that this variation is attributable to underlying fluctuation in allele frequencies within populations, or methodological differences between studies, it is more likely to be due to stochastic sampling effects, given the relatively small number of individuals pooled in some samples.
Similar to the wMel genomes, we detected very few changes between the mitochondrial genomes sequenced in this study, for which we obtained � 99.31% coverage. The genomes

PLOS PATHOGENS
wMel stability in Aedes aegypti shared 100% sequence identity across all gene-encoding regions. Nucleotide differences were observed among the consensus sequences at 24 positions within an approximately 1.8Kb region of the genome that does not contain any predicted gene sequences, is very AT rich and displayed greatly reduced read mapping depths (S1 Table). It is therefore likely that many or all of these differences represent sequencing errors. A further 16 positions within this region were found to show within-sample polymorphism, with alternate alleles present at a frequency greater than the threshold value in at least one sample (S2 Table). The lack of variation observed among the wMel genomes and among the mitochondrial genomes is consistent with a high fidelity of maternal co-transmission of both Wolbachia and mitochondria.

Discussion
In this study, we provide a comprehensive update on the genomic and phenotypic stability of the wMel Wolbachia infection, up to ten years after being released in Ae. aegypti populations. We observed few changes in the wMel genome since before releases began, with little divergence between different locations despite nearly a decade of separation. Furthermore, wMel retains complete cytoplasmic incompatibility, a stable within-host density, and limited host fitness costs. Our results suggest that any shifts in the host phenotypic effects of wMel are likely to reflect evolutionary changes within mosquito populations rather than the Wolbachia genome. This hypothesis is consistent with crossing and selection experiments showing that Wolbachia phenotypic effects depend on nuclear background [13,41]. It is also consistent with genomic data from this and previous studies showing the long-term stability of Wolbachia genomes following transinfection [33,34,39,40]. Our fitness experiments do not provide evidence for evolutionary changes in response to Wolbachia infection, with similar effects observed in mosquito populations that had been infected with wMel for many years (experiment 1) and in populations where wMel was recently introgressed (experiment 2). Comparisons of mosquito genomes prior to wMel release in Gordonvale and seven years post-release point to limited changes [38], but further work with replicated populations is required to confirm any changes due to Wolbachia infection.
Accurate phenotyping of Wolbachia infections requires careful control of the nuclear and mitochondrial background between Wolbachia-infected and uninfected populations. In this study, we used tetracycline curing followed by backcrossing to ensure that Wolbachia-infected and uninfected populations had matched mitochondrial and nuclear genomes. Cross-generational effects of tetracycline treatment on fitness have been speculated [42], with some studies accounting for potential disruptions to the microbiome by rearing treated mosquitoes in water from untreated mosquitoes (e.g. [43]). We found no evidence to suggest that tetracycline treatment influences fitness in the following generation, at least for a concentration commonly used to cure Wolbachia infections in insects [44]. While the effect of tetracycline treatment on the mosquito microbiome (aside from Wolbachia) remains to be tested, repeated backcrossing to a common background should help to minimize any differences between populations.
Our data, spanning up to a decade since the first releases of wMel in Australia, should be informative for Wolbachia release programs that have taken place more recently in other countries. Our results suggest that when the wMel infection is maintained in an Ae. aegypti population, the phenotypic effects associated with wMel invasion are likely to persist given that the infection

PLOS PATHOGENS
wMel stability in Aedes aegypti remains at a high density. However, more work is required to understand the extent of genetic changes in mosquito populations in response to Wolbachia releases. These may include host genomic changes in response to the presence of Wolbachia or indirect effects from the introduction of host genes from release stocks, such as the introduction of pesticide susceptibility genes in Tubiacanga, Brazil, which likely contributed to an unsuccessful Wolbachia release there [19]. Pesticide susceptibility is an issue in most countries where dengue is endemic and chemicals including pyrethroids, organophosphates, and insect growth regulators are applied to suppress Ae. aegypti populations [45]. However, this has not been an issue in the Cairns region where limited applications of pesticides have likely prevented the local evolution of resistance [46].
While wMel has retained complete cytoplasmic incompatibility and maternal transmission under laboratory conditions, these parameters are affected by environmental conditions [10] and require further evaluation under field conditions. The persistence of a high incidence of wMel in natural populations may be constrained by environmental conditions including high temperatures in some locations [25]. Furthermore, interactions with mosquito genetic background may influence the effects of Wolbachia infection on host fitness [41], which could help to explain variability in wMel establishment success in different countries. While several studies have measured wMel effects in local mosquito backgrounds (e.g. [47,48]), direct comparisons across multiple mosquito strains are needed to understand the contribution of host genetics to phenotype. Finally, ongoing monitoring remains important to identify any changes in Wolbachia infection frequency and inform the need for supplementary releases including those with different host genetic backgrounds and different Wolbachia strains.

Ethics statement
Blood feeding of female mosquitoes on human volunteers for this research was approved by the University of Melbourne Human Ethics Committee (approval 0723847). All adult subjects provided informed written consent (no children were involved).

Field collections and colony establishment
Aedes aegypti mosquitoes were collected as eggs from suburban Cairns in 2016, 2018, 2019 and 2020 from ovitraps. Felt strips from ovitraps were collected and processed identically to previously described methods . Subsets of Ae. aegypti larvae from central Cairns, Gordonvale and Yorkeys Knob were pooled from all traps within a suburb to establish laboratory populations for phenotypic comparisons. Thirty individuals from the F 1 and F 2 generations were screened for Wolbachia infection (see below) to confirm fixation of wMel within the laboratory populations. All populations were maintained at census size of~450 individuals per generation at 26˚C and a 12:12 light:dark cycle as described previously [50]. Female mosquitoes were fed on the forearm of a single human volunteer for egg production. collections, life stage was uncontrolled, with a mix of adults and larvae tested, thus data were only suitable for Wolbachia frequency estimates. Up to 10 individuals were screened per trap, with between 10 and 243 individuals screened per suburb. For the 2018 field collections, 30 4 th -instar larvae from 15 ovitraps (2 per trap) per suburb were screened for Wolbachia infection. Larvae were reared at a controlled density (50 larvae per tray) and stored at the same age (5 d post-hatching), allowing for a comparison of Wolbachia density across suburbs. In June 2019, we measured wMel density in adults (15 females and 15 males) from laboratory populations that were established from field collections at different times (Gordonvale at F 1 , Yorkeys Knob at F 1 and F 16 , central Cairns at F 9 and the wMel lab population). Aedes aegypti from field collections and laboratory populations were screened for Wolbachia infection status and density using a Roche Lightcycler 480 according to previously described methods [51]. Genomic DNA was extracted with 250 μL of 5% Chelex 100 Resin (Bio-Rad laboratories, Hercules CA) and 3 μL of Proteinase K (20 mg/mL) (Roche Diagnostics Australia Pty. Ltd., Castle Hill New South Wales, Australia). Tubes were incubated for 60 min at 65˚C then 10 min at 90˚C. Three primer sets were used to amplify markers specific to mosquitoes (mRpS6_F 5'AGTTGAACG-TATCGTTTCCCGCTAC3' and mRpS6_R 5' GAAGTGACGCAGCTTGTGGTCGTCC3'), Ae. aegypti (aRpS6_F 5'ATCAAGAAGCGCCGTGTCG3' and aRpS6_R 5'CAGGTGCAGGA TCTTCATGTATTCG3'), and wMel (w1_F 5'AAAATCTTTGTGAAGAGGTGATCTGC3' and w1_R 5' GCACTGGGATGACAGGAAAAGG3'). Relative Wolbachia densities were determined by subtracting the crossing point (Cp) value of the wMel-specific marker from the Cp value of the Ae. aegypti-specific marker. Differences in Cp were averaged across 3 consistent replicate runs, then transformed by 2 n .

Phenotypic comparisons
We performed two sets of experiments to evaluate the phenotypic effects of wMel derived from field and laboratory populations. In the first set, populations were cured with tetracycline to remove the wMel infection, maintaining similar genetic backgrounds between infected and uninfected lines. In the second set, we used backcrossing to introduce the wMel infection from different origins to a common genetic background, then compared populations to uninfected lines that had been crossed to the same background. Both sets of experiments involved the wMel Lab population which was collected from Cairns in 2014 and had spent at least 60 generations in the laboratory before the first set of experiments commenced. All experiments were performed at 26˚C and a 12:12 light:dark cycle. Experiment 1 was performed in September 2018 using wMel-infected populations collected from Yorkeys Knob in February 2018 (wMel YK) and Cairns in 2014 (wMel Lab). wMel YK and wMel Lab were divided into four population cages each. Two replicate populations from each line were treated for three consecutive generations with tetracycline hydrochloride (2 mg/mL) provided to adults in 10% sucrose solution to cure the wMel infection. Females were blood fed at 10 d old to ensure that they had fed on the antibiotic solution prior to blood feeding. The other two replicate populations from each line were left untreated but reared synchronously. All populations were reared in the absence of antibiotics for two generations before experiments commenced, with thirty adults from each population screened for Wolbachia infection to ensure complete removal of wMel in the treated lines (wMel YK.tet and wMel Lab. tet) and fixation of wMel in the untreated lines (wMel YK and wMel Lab). The wMel YK and wMel YK.tet populations were at F 8 in the laboratory when experiments commenced. Experiment 2 was performed in August 2020 using wMel-infected populations collected from Yorkeys Knob in February 2020 (wMel YK), Gordonvale in February 2020 (wMel GV), Cairns in 2014 (wMel Lab), as well as an uninfected population that had been cured of wMel PLOS PATHOGENS wMel stability in Aedes aegypti in the previous experiment (wMel Lab.tet). Two hundred females from each population were crossed to 200 males from a natively uninfected population collected from locations in Cairns prior to wMel releases. This process was repeated for two further generations to produce a similar nuclear background between populations. The backcrossed wMel YK and wMel GV populations were at F 5 in the laboratory when experiments commenced.
In both experiments, we measured larval development time, adult wing length, female fecundity and egg hatch proportions. Eggs (<1 week old) from each population were hatched in reverse osmosis (RO) water and 100 larvae (<1 d old) were counted into plastic trays filled with 500 mL RO water (with 6 replicate trays per population in experiment 1 and 12 replicate trays per population in experiment 2). Larvae were fed TetraMin Tropical Fish Food tablets (Tetra, Melle, Germany) ad libitum (with daily monitoring and removal of excess food to prevent overfeeding) until pupation. Development time was scored by counting and sexing pupae twice per day. Adults were pooled across replicate trays and released into cages. Wings of adults (20 males and 20 females per population) were dissected and measured for their length (from the alular notch to the wing tip). Females (5-7 d old, sugar-starved for 1 d) were blood fed and isolated in 70 mL specimen cups filled with 15 mL of larval rearing water, lined with a strip of sandpaper (Norton Master Painters P80; Saint-Gobain Abrasives Pty. Ltd., Thomastown, Victoria, Australia) and covered with a mesh lid. Twenty females were isolated per population in experiment 1. In experiment 2, fecundity and egg hatch proportions were tracked across four consecutive gonotrophic cycles by isolating 30 engorged females per population, returning females that laid eggs to cages, then blood feeding and isolating 30 engorged females again every 4-5 d. Eggs were collected 4 days after blood feeding, partially dried, then hatched 3 d after collection in trays filled with RO water and a few grains of yeast. Fecundity and egg hatch proportions were determined by counting the total number of eggs and the number of hatched eggs (with the egg cap clearly detached) under a dissecting microscope. In experiment 2, we also measured adult longevity (8 replicate cages of 25 males and 25 females per population) by maintaining adults in 3 L cages with cups of 10% sucrose and scoring mortality 3 times per week. Quiescent egg viability was measured in experiment 2 by storing eggs on sandpaper strips in a sealed container with a saturated solution of potassium chloride to maintain a constant humidity of~80%. Twelve replicate batches of eggs per population (median 68 eggs) were hatched on week 1 and 2, then every two weeks until week 22. Egg hatch proportions were scored as above for individual females.

Effects of tetracycline treatment on fitness
In the above experiments we used populations that had been cleared of Wolbachia infections through tetracycline treatment. Despite allowing for at least two generations of recovery, these treatments could potentially disrupt the microbiome (including the mitochondria), leading to fitness differences between lines that are independent of Wolbachia infection. To test whether antibiotic treatment has any effect on fitness, we fed wMel-infected (wMel Lab) and natively uninfected (F 29 in the laboratory) adults with 2 mg/mL tetracycline hydrochloride for 10 d before blood feeding, then measured fitness in the subsequent generation. We scored larval development time, survival to pupa and sex ratio (6 replicate trays per treatment) as well as female and male wing length, female fecundity and egg hatch proportions (20 individuals each). Experiments were performed identically to the phenotypic comparisons above.

Cytoplasmic incompatibility
We tested the ability of wMel-infected Ae. aegypti to induce cytoplasmic incompatibility in near-field and laboratory populations. Reciprocal crosses were performed in June 2019

PLOS PATHOGENS
wMel stability in Aedes aegypti between wMel-infected populations established from Gordonvale at F 1 (wMel GV F 1 ), Yorkeys Knob at F 1 (wMel YK F 1 ) or Cairns at F 60+ (wMel Lab), and a natively uninfected population (F 35 in the laboratory). Pupae from each population were sexed and released into separate 3 L cages to confirm accurate sex sorting. Groups of 40 females and 40 males (1 d old) were then aspirated into cages together and left for 5 d to mate. Females (starved of sugar for 24 hr) were blood-fed, isolated for oviposition and scored for fecundity and egg hatch proportions according to the methods described above for the phenotypic comparisons.

Wolbachia and mitochondrial whole genome sequencing
Previously, we sequenced the wMel and mitochondrial genomes of Ae. aegypti collected up to eight years after field releases [39]. Here, we extend these findings by sequencing the Wolbachia and mitochondrial genomes of wMel-infected Ae. aegypti sampled a decade apart from pre-release and post-release populations. The pre-release wMel-infected populations (wC45 F 9 and F 10 ) were sampled in 2010 and stored at -20˚C, while three wMel-infected populations were sampled in 2020: wMel GV F 2 collected from Gordonvale in 2020, wMel YK F 2 collected from Yorkeys Knob in 2020, and wMel Lab, collected from Cairns in 2014 and maintained under laboratory conditions until sampling. Genomic DNA was extracted from pooled samples containing five adult females. Sequencing libraries were then prepared as described previously [52].

Reference genome assembly
Sequencing reads were quality filtered using Trimmomatic [53]. The samples were trimmed in paired-end mode with the following parameter settings: leading = 20; trailing = 20; slidingwindow = 4:20; minlen = 70; adapter sequences were removed using the ILLUMINACLIP option, with maximum seed mismatches = 2 and the palindrome clip threshold = 30. Reads were aligned to a wMel reference genome (GenBank accession: NC_002978.6) and an Ae. aegypti mitochondrial reference genome (GenBank accession: MH348177.1) with the Burrows-Wheeler Aligner (BWA; [54]) using the bwa mem algorithm and default parameter settings. Quality filtering of alignments and variant calling was performed with SAMtools and BCFtools [55,56]. PCR duplicates were excluded from subsequent analyses by soft masking. Reads with a MAPQ score < 25 were removed from the alignment, except for reads with MAPQ = 0, which were permitted to allow for mapping to repetitive regions. Genotype likelihoods were calculated using a maximum of 2000 reads per position. For variant calling, ploidy was set to haploid. The variant call output was used to create a consensus nucleotide sequence, wherein genome positions with coverage < 5 were masked as 'N'. Loci were considered to be polymorphic within a sample if they had coverage � 30 and two or more alleles with a frequency of � 25%. Any such loci occurring within the 16S rRNA or 23S rRNA genes were excluded from analysis, as we observed a relatively high level of contaminant reads mapping to these regions. Genome sequences were inspected and aligned with Geneious v 9.1.8 (https://www. geneious.com).
Kraken2 [57] and the Standard-8 precompiled reference database (https://benlangmead. github.io/aws-indexes/k2; downloaded 17/9/21) were used to search for sequence contamination within the Wolbachia genomes. The sequencing reads mapped by bwa to the wMel reference genome were filtered to remove reads matching taxa other than Wolbachia, and genome assemblies were then repeated with the filtered datasets, using the above pipeline. The original pre-filtration genome sequences were edited to correct erroneous positions after comparison with the corresponding post-filtration genome sequences.

PLOS PATHOGENS
wMel stability in Aedes aegypti

Statistical analysis
Most experimental data (including Wolbachia density, development time, wing length, fecundity and egg hatch) were analyzed with general linear models (GLMs) while adult longevity data were analysed with log-rank tests. All analyses were carried out with IBM SPSS Statistics 26. Data were transformed where appropriate (with all proportional data being logit transformed). The first experiment had two replicate populations and data were initially analysed with replicate population (nested within population origin x Wolbachia infection status) included as a factor. Replicate populations were then pooled for a second analysis due to a lack of significant effect of replicate population (P > 0.1) for any trait. Wolbachia infection status and population origin were included as factors. In the second experiment, fecundity and egg hatch proportions were tracked across gonotrophic cycles with the same mosquitoes, so we ran separate analyses for each cycle. For quiescent egg viability, wMel-infected and uninfected populations were analysed both together and separately at the first (1 week) and last time (22 weeks) points. Field Wolbachia density data were analysed with release year and suburb (nested within release year) as factors. We performed Bonferroni corrections where multiple traits were evaluated in the same cohort of mosquitoes.  Table. Observed differences in Wolbachia and mitochondrial genomes. Positions are shown relative to the GenBank reference sequences used for read mapping. Some or all of the differences are likely to be sequencing artifacts (see main text). For comparison, the results of Huang et al. [39] and Dainty et al. [40] are also shown. ND = not determined. (XLSX) S2 Table. Genome positions found to show polymorphism within populations. In our study, loci were considered to be polymorphic within a sample if they had coverage � 30 and two or more alleles with a frequency of � 25%. Any such loci occurring within the 16S rRNA or 23S rRNA genes were excluded from analysis, as we observed a relatively high level of contaminant reads mapping to these regions. Positions are shown relative to the GenBank reference sequences used for read mapping. For comparison, the results of Huang et al. [39] and Dainty et al. [40] are also shown. For Dainty et al. [40], multiple samples per population were analysed, with each sample containing a single individual-allele frequencies are reported as number individuals/total individuals per population. NR = not reported; not scored as polymorphic by the authors of the original study; in most cases these values are likely to be at or close to 100%.