Forward genetics in Wolbachia: Regulation of Wolbachia proliferation by the amplification and deletion of an addictive genomic island

Wolbachia is one of the most prevalent bacterial endosymbionts, infecting approximately 40% of terrestrial arthropod species. Wolbachia is often a reproductive parasite but can also provide fitness benefits to its host, as, for example, protection against viral pathogens. This protective effect is currently being applied to fight arboviruses transmission by releasing Wolbachia-transinfected mosquitoes. Titre regulation is a crucial aspect of Wolbachia biology. Higher titres can lead to stronger phenotypes and fidelity of transmission but can have a higher cost to the host. Since Wolbachia is maternally transmitted, its fitness depends on host fitness, and, therefore, its cost to the host may be under selection. Understanding how Wolbachia titres are regulated and other aspects of Wolbachia biology has been hampered by the lack of genetic tools. Here we developed a forward genetic screen to identify new Wolbachia over-proliferative mutant variants. We characterized in detail two new mutants, wMelPop2 and wMelOctoless, and show that the amplification or loss of the Octomom genomic region lead to over-proliferation. These results confirm previous data and expand on the complex role of this genomic region in the control of Wolbachia proliferation. Both new mutants shorten the host lifespan and increase antiviral protection. Moreover, we show that Wolbachia proliferation rate in Drosophila melanogaster depends on the interaction between Octomom copy number, the host developmental stage, and temperature. Our analysis also suggests that the life shortening and antiviral protection phenotypes of Wolbachia are dependent on different, but related, properties of the endosymbiont; the rate of proliferation and the titres near the time of infection, respectively. We also demonstrate the feasibility of a novel and unbiased experimental approach to study Wolbachia biology, which could be further adapted to characterize other genetically intractable bacterial endosymbionts.

We designed the screen to find new mutants of Wolbachia, and not host mutants, that lead 158 to the endosymbiont proliferation. However, EMS will most likely also induce mutations in 159 the host germline, in the nuclear or mitochondrial genomes, that can be transmitted. To 160 minimize the influence of new host nuclear mutations on our screen, we backcrossed the 161 EMS-treated females and their progeny, at every generation, with males from the control 162 isogenic line. To verify that new mutations in the host were not the cause of Wolbachia 163 over-proliferation, we replaced the first, second and third chromosomes of Drosophila 164 females carrying the over-proliferating Wolbachia variants, in line 1A and line 2A, with the 165 chromosomes of the control line, through the use of balancer chromosomes (S4 Fig.). We 166 repeated Wolbachia titres quantification and found that the over-proliferative phenotypes 167 were maintained (Fig 1C; lmm, p-value< 0.001 for both). 168 Since mitochondria are maternally transmitted, the experiments described above cannot 169 exclude the possibility that Wolbachia over-proliferation is mitochondria-determined 170 (directly or indirectly). Thus, the mitogenome of the lines showing higher Wolbachia titres 171 was Illumina sequenced and aligned to the mitochondrial reference genome release 6 172 (GenBank: KJ947872.2:1-14,000). A summary of the alignment features is given in S2 173 Table. We did not find SNPs or indels unique to the mitochondria of EMS-fed flies, which 174 shows that flies with over-proliferative Wolbachia did not inherit mutated mitochondria (S3 175   Table). Therefore, we concluded that the observed Wolbachia over-proliferative 176 phenotypes did not result from mutations in the host genome, neither nuclear nor 177 mitochondrial. 178 179

Identification of genetic basis of the new over proliferative variants 180
To identify the mutations associated with over-proliferation, we sequenced and assembled 181 the genomes of these over-proliferative Wolbachia. A summary of the sequencing results 182 is given in S2 table. We performed a hybrid assemble with short (Illumina) and long-reads 183 (Nanopore) and obtained single and circular Wolbachia chromosomes differing in length 184 (S4 Table). 185 To test our assembly pipeline we assembled a previously characterized Cluster III wMel 186 Wolbachia variant, named wMel [16], which derives from the line used for the original 187 wMel reference genome (GenBank: AE017196.1) [27]. The new wMel genome (GenBank: 188 CP046925.1) was also circular and comparable in size, structure and number of ORFs 189 with previously published wMel genomes [27,28], including the wMel reference genome 190 (S4 Table). We found, however, nine differences (2 SNPs and seven indels) relative to wMel reference genome, which we confirmed using Sanger sequencing (S5 Table). These 192 results validate our sequencing pipeline. 193 The only difference between the genome of the over-proliferative Wolbachia variant in Line 194 1A and wMelCS_b was an amplification of the Octomom region (Fig 2A and S1 Text). 195 There were three more copies of this region, giving a genome size difference of 62,814bp. 196 The Octomom region amplification, and lack of other differences, was also confirmed by   These results identify Octomom loss as the cause of its over-proliferative phenotype.  Mapping the Illumina sequence reads of this line identifies the loss of Octomom as the 236 only difference with wMelCS_b, and identifies no difference to wMelOctoless. We named 237 this line wMelOctoless2. These results further confirm that loss of the Octomom region 238 leads to an over-proliferative phenotype in Wolbachia. 239 We also sequenced by short and longs reads and made a hybrid assembly of wMelPop 240 genome and compared it to wMelPop2. We only detected and confirmed the two SNPs  (Table S3). We 245 confirmed this SNP using Sanger sequencing (S2 Text). 246 In summary, we were able to identify the genomic changes associated with the new over-247 proliferative variants and all map to loss or amplification of the Octomom region.

267
Wolbachia proliferation dynamics at 18ºC (A), 25ºC (B) and 29ºC (C). D. melanogaster males used in these 268 experiments developed at 25ºC, were collected on the day of eclosion and aged to specific time-point at a 269 given temperature (18ºC, 25ºC or 29ºC). Ten males were collected at each time-point for Wolbachia titre 270 measurement using qPCR. Wolbachia titres were normalized to that of 0-1 days-old wMelCS_b-infected 271 males. A replicate of the experiment is given in S8 Fig. Exponential

275
To analyse growth during adult life we fitted an exponential model to these data, from 276 which we also estimated doubling time of the Wolbachia variants, at the different 277 temperatures (Table 1). There is a wide range of doubling times for different Wolbachia 278  18ºC and 25ºC (p = 0.94), but increases at 29ºC (p < 0.001). 292 We also compared the effect of Octomom copy number in a model taking genotype 293 (wMelCS_b, for the control and new variants generated here, or wMelPop) and Octomom 294 copy number (zero, one, low and high) as factors. At 18ºC the amplification of the 295 Octomom copy number does not lead to a statistically significant increase of growth rate, 296 compared with wMelCS_b (lmm, p > 0.125 for both low and high Octomom copy number 297 compared with control). At this temperature only the deletion of Octomom leads to a 298 significantly higher growth, as seen above. However, variants with Octomom amplification grow more at 25ºC than at 18ºC, and grow more at 29ºC than at 25ºC (p < 0.001 for all 300 comparisons). At both 25ºC and 29ºC, variants with low Octomom copy numbers have a 301 higher growth rate than the variants with one copy or deletion of Octomom region (p < 302 0.001 for all comparisons). At these temperatures the variants with high Octomom copy 303 numbers grow faster than all the other variants (p < 0.003 for all comparisons). These 304 results confirm that the degree of amplification of the Octomom region controls the 305 intensity of the over-proliferation of these variants, as shown before [17]. The data also 306 demonstrate a strong interaction between temperature and the increased proliferation of 307 variants with amplification of the Octomom region. The effect of the amplification is not 308 significant at 18ºC and becomes increasingly stronger with higher temperature. On the 309 other hand, loss of Octomom leads to a smaller effect in growth, but similar at all 310 temperatures. Therefore, although both genomic mutations lead to an increase in 311 Wolbachia titres they have different impacts in the growth rates and their interaction with 312

temperature. 313
The statistical model taking into account the wMelCs_b or wMelPop genotype, which differ 314 only in two SNPs (see above), indicated a significant difference in growth between them, 315 at 25ºC (p < 0.001). This could indicate that these two SNPs also influence growth of 316 Wolbachia. However, this could also be due to the fact that the copy number of the 317 Octomom region was not equally controlled in wMelPop and wMelPop2 lines during these 318 experiments. wMelPop2 low copy number line had 2-3 copies of Octomom, while the 319 wMelPop line had 3 copies. To test if wMelPop and wMelPop2 indeed vary in proliferation 320 rate, we repeated this experiment with a more tightly controlled Octomom copy number in 321 these two lines, at 25ºC (S10 Fig A-D). Both wMelPop and wMelPop2 carrying 3 copies of 322 Octomom grow faster than wMelCS_b (lmm, p < 0.001 for both) and there is no difference 323 in growth between them (p = 0.32). This indicates that the genetic differences between 324 these lines do not affect their growth and that they are equally influenced by Octomom 325 copy number. 326 327

Rapid proliferation of Wolbachia during the host development 328
We also analysed the growth of wMelCS_b, wMelOctoless and wMelPop2 (8-9 Octomom 329 copies) during host development. D. melanogaster develops very fast from egg to adult, in 330 only 10 days at 25ºC. We predicted that Wolbachia would grow much faster during this 331 period than during adult life, in order to grow from the small population, present in the egg, 332 to the population spread throughout the adult tissues. We, therefore, estimated absolute numbers of Wolbachia genome copies in individuals at the different stages of  Table 2). Interestingly, at the end of the experiment, in newly eclosed adults, and as 339 also observed above, there are significant differences between the three variants (lm, p < 340 0.009 for all comparisons, Table 2, Fig 4). Also, males carry significantly less Wolbachia 341 than females (p = 0.029). Adults carried from approximately 600,000 to 4,500,000 342 Wolbachia genome copies.

346
Wolbachia genome copies per individual throughout Drosophila development. Dots represent either data 347 from a pool of 10 individuals (eggs and larvae) or from a single individual (pupae and adults). Wolbachia 348 proliferation in the first 120 hours were analysed using a exponential model. A summary is given in Table 1.

353
Mean (and standard error) of Wolbachia genome copies were estimated using qPCR.

355
Wolbachia growth seems to be restricted to the period between egg and white prepupae 356 rates is not statistically significant (p = 0.10 for interaction between Wolbachia variants and 363 growth). The growth rates of these variants are, therefore, very similar during this stage, 364 and much faster than in adults. At the same temperature, we estimated doubling times in 365 adults of wMelCS_b, wMelOctoless, and wMelPop2 (high-copy) to be, approximately, 13.9, 366 10.5, and 2.3 days, respectively. Therefore, Wolbachia growth at different stages of D. 367 melanogaster can vary dramatically, and the different variants respond differently to host 368 development stages. 369 We also asked if Wolbachia Octomom copy number changed in wMelCS_b and 370 wMelPop2, throughout development, as Wolbachia is proliferating so fast, and found no  The over-proliferation of wMelPop has been associated with a shortening of the host 378 lifespan [16,22]. We, therefore, tested if these new over-proliferative variants also shorten the lifespan of their host, at different temperatures, in males (Fig 5, S12 Fig). We also 380 performed this assay in females at 25ºC, with similar results to males at 25ºC (S12 Fig). There was a significant interaction between survival, Wolbachia variant, and temperature 382 (Cox proportional hazard model with mixed effects (CHR), p < 0.001). All lines, including 383 the Wolbachia-free line have a shorter lifespan at 25ºC than 18ºC, and even shorter at 384 29ºC (p < 0.001 for all these comparisons). wMelCS_b did not affect the host lifespan at 385 any temperature (p > 0.16 for all comparisons with the line without Wolbachia). 386 wMelOctoless strongly reduces host lifespan at all tested temperatures (p < 0.001, each 387 comparison with wMelCS_b) (Fig 5, S12 Fig). This deleterious effect is stronger at 18ºC, 388 where wMelOctoless is the tested variant with the highest impact on lifespan, although 389 very similar and not statistically different from wMelPop2 with high Octomom copy number 390   with high copy number of Octomom have the shortest lifespan of all tested lines (p < 0.001 416 for all comparisons). wMelPop2 and wMelPop with low copy number of Octomom always 417 have a weaker effect than high copy number variants (p < 0.001 for all these 418 comparisons), demonstrating the effect of the degree of amplification in this phenotype. As 419 observed with the high copy number variants their effect increases with temperature and is 420 stronger at 25ºC than 18ºC, and stronger at 29ºC than 25ºC (p < 0.05 for these 421 comparisons). In fact, they are not pathogenic at 18ºC and an effect on host lifespan is 422 only observed at 25ºC and 29ºC. These data confirm the association of Octomom region 423 amplification with host lifespan shortening, an increase in the severity of the phenotype 424 with increase in the Octomom copy number, and the increase in the severity of the 425 phenotype with temperature. 426 wMelPop2 and wMelPop, high and low copy number variants, had significantly different 427 effects in several comparisons at the three tested temperatures (Fig 5). This indicated that 428 there could be differences in this phenotype between these two lines. However, the 429 Octomom copy number was not exactly the same between these lines. Therefore, and as The life shortening phenotype of wMelPop has been associated with its over proliferation 448 since its discovery [22]. We tested if these phenotypes were correlated by taking advantage of the data on proliferation and shortening of lifespan, at different temperatures,   We found no correlation at either 18ºC (p = 0.21), the temperature where the flies were 487 kept after infection, or 25ºC (p = 0.35), the temperature in which flies developed and were 488 kept until being infected with DCV ( Fig S14). We also tested the correlation between the 489 antiviral strength of Wolbachia-induced anti-viral protection with Wolbachia titres in 0-1 490 day-old flies, as a proxy for Wolbachia titre at the day of infection, which is 3-5 day-old 491 flies. We found a significant correlation between these parameters (Fig 6F, p = 0

An unbiased approach for genetically intractable symbionts 513
A difficulty in studying obligate intracellular bacteria is that they are often genetically 514 intractable, mostly because of their dependency on the host cells to survive and difficulty in 515 isolating and growing them clonally. We decided to try to directly mutagenize and screen 516 Wolbachia in the host. The main difficulty of this approach is how to identify a new mutant 517 and link it with a phenotype in the Wolbachia population present within a host. We calculated here that newly emerged female adults, with wMelCS_b, carry approximately 519 840,000 Wolbachia genome copies, probably corresponding to the same number of 520 Wolbachia cells [30]. Since EMS induces random mutations, we expected mosaicism in 521 the Wolbachia population at the individual fly level. Each new mutant, when generated, 522 would be a unique cell within these 840,000 other cells. Since the Wolbachia phenotypes 523 are normally measured at the individual host levels (e.g. Wolbachia titres, antiviral 524 protection), the properties of individual or small numbers of mutant Wolbachia could be 525 diluted and unmeasurable. 526 We hypothesized, however, that over-proliferating Wolbachia cells could overtake the 527 population and that the resulting higher titres could be detectable. Indeed, fast proliferative 528 Wolbachia can be selected at the level of a single host [31]. To increase the probability of 529 isolating rare over-proliferating Wolbachia variants, we also relied on the bottleneck 530 imposed in the vertical transmission of Wolbachia. We calculated here that single embryos 531 carry approximately 4,000 Wolbachia genomes, which is consistent with previous 532 estimates [32]. Moreover, we treated flies with tetracycline to further reduce this 533 population, and increase drift. By screening at the immediate progeny (F1) of EMS-treated 534 females or three generations later (F4) we were able to select new over-proliferative 535 Wolbachia mutants. 536 537

Genetic bases of Wolbachia over-proliferation 538
After discarding the possibility that mutations in the host nuclear or mitochondrial genomes 539 were the cause of Wolbachia over-proliferation, we performed de novo assembly of the 540 ancestral, wMelCS_b, and the new mutant variants, wMelPop2 and wMelOctoless. The 541 assemblies generated complete full chromosomes of these Wolbachia. wMel and 542 wMelPop were also assembled. This allowed us to identify single nucleotide differences 543 and structural differences between these genomes [33]. To validate our pipeline we 544 compared our wMel genome assembly to the reference wMel genome. We identified only 545 seven indels and two SNPs, which we confirmed to be present in our line, by Sanger 546 sequencing. Our assembly results also allowed us to confirm two previously identified 547 SNPs between wMelCS_b (but also wMelOctoless and wMelPop2) and wMelPop. 548 Additionally, our assembly provides a significant improvement over the previous wMelPop 549 genome [26]. 550 The only differences between the new over-proliferative variants and wMelCS_b were 551 structural differences in the Octomom region. wMelPop2 has an amplification of this region. The assembly confirms that the Octomom region is amplified in tandem [17], and 553 that all copies are located in the Wolbachia genome. We had shown before that 554 amplification of the Octomom region and the degree of this amplification determined 555 wMelPop over-proliferative phenotype [16,17]. If the Octomom region copy number in 556 wMelPop reverted to one, the variant became phenotypically identical to wMelCS_b [17]. 557 We now show that if we start from wMelCS_b with one copy of Octomom, an increase in 558 the copy number of this genomic region leads to an over-proliferative phenotype. The bacterial community associated with the fly stocks was homogenized as in Pais et al. 682 [44], with minor modifications. Briefly, we collected eggs for 6 hours in fresh agar plates 683 supplemented with live yeast and sterilized the eggs surface by consecutive washes on 2.1% sodium hypochlorite (NaOCl) solution (10 minutes), 70% ethanol (5 minutes) and 685 sterile water (5 minutes). Next, we transferred axenic eggs to sterile fly food supplemented 686 with 40μL of 1:1 overnight culture of Acectobater OTU 2753 and Lactobacillus OTU 1865 687 [44]. We confirmed the presence of these bacterial species by squashing five females 688 aged 3-6 days in sterile 1x PBS, plating 30μL of the lysate in mannitol plates, incubate 689 them at 25ºC for 72h, and identify bacteria by colony morphology. Flies with wMelCS_b developed in fly food supplemented with tetracycline at the 702 concentrations 1.5625µg/ml, 3.125µg/ml, 6.25µg/ml, 12.5µg/ml, 25µg/ml, and 50µg/ml. 703 Three isofemale lines were established from each dose. In the F1, we randomly selected 704 four virgin females for egg-laying and Wolbachia titre measurement using qPCR. We set 705 up the next four generations using the progeny of a female with the median Wolbachia 706 titres. 707 708

Forward genetic screen 709
We attempted to mutagenize Wolbachia in vivo by feeding its host with the mutagen EMS. 710 DrosDel w 1118 isogenic flies carrying wMelCS_b were first raised in standard fly food or fly 711 food supplemented with tetracycline (from 1.5625 to 12.5µg/ml). Virgin females were 712 collected, starved for 6h, and then fed on a 1% sucrose solution with EMS concentrations 713 ranging from 10 to 8,000mM. Control flies fed on sucrose solution only. A dye was added 714 to the feeding solution to confirm intake and feeding proceeded for 13h (overnight). free Drosdel w 1118 isogenic males, egg-laying was allowed for 3-4 days, and parents 717 discarded. From the F1 progeny we collected virgin females, mated them individually with 718 2-3 Wolbachia-free Drosdel w 1118 isogenic males, and allowed egg laying for 3-4 days. 719 These females were collected when 10 days old, and Wolbachia titres determined by 720 qPCR. We followed the progeny of F1 females showing 50% or higher increase in 721 Wolbachia titres relative to control flies in the same conditions. We have also transferred 722 the progeny of these F1 for three more generations, without selection, and repeated the 723 determination of Wolbachia titres in F4 females. In the same batch of experiments we may 724 have tested more than one F1 or F4 progeny from each G0 female. Hence, over-725 proliferative Wolbachia variants isolated in the same batch of treated females may be a 726 result of a single event in the G0 female. 727 728

Real time quantitative PCR 729
DNA extraction for qPCR was performed as described before [17]. 730 The qPCR reactions were performed in the QuantStudioTM 7 Flex (Applied Biosystems). 731 The reaction mix and PCR program used were described before [16]. The specificity of the 732 amplicons were confirmed by analysing the melting curve profile.    We further refined our genome assemblies using a two-stepped pipeline. First, we mapped 815 Illumina reads to the corresponded draft genomes to identify mismatches, which were later 816 corrected via Sanger sequencing. Next, we improved genome assemblies by aligning 817 them using Mauve v2.4.0 [56], and the predicted differences were also confirmed by 818 Sanger sequencing. Only two predicted SNPs which are unique to wMelPop were real. 819 820

Read mapping and mutation calling 821
We identified mutations relative to a genome following a previously published pipeline [57]. To compare Wolbachia titres across multiple groups, we used log linear models (LMM) or 839 generalized linear mixed models assuming Gamma errors (GLMM). The effect of EMS on 840 Wolbachia titre was tested using non-linear regression. We estimated the doubling time of 841 Wolbachia variants using the equation log(2)/β, with β being the coefficients of a log linear 842

model. 843
The lifespan datasets and survival curves after challenge with DCV were analysed with 844 mixed effect Cox models [61]. 845 The significance of correlations were tested using Pearson's correlation coefficient. 846 In multiple comparisons were necessary, the p-values were adjusted as proposed by Holm 847

Development. 861
The fly work at the Fly Facility of Instituto Gulbenkian de Ciencia (Oeiras, Portugal), was 862 partially supported by the research infrastructure Congento, co-financed by Lisboa 863 Regional Operational Programme (Lisboa2020), under the PORTUGAL 2020 Partnership 864 Agreement, through the European Regional Development Fund (ERDF) and Fundação 865 para a Ciência e Tecnologia (Portugal) under the project LISBOA-01-0145-FEDER-866

1102
The deletion of individual Octomom genes was confirmed using qPCR. The copy number of three genes 1103 outside the Octomom region were also determined. Five females carrying wMelCS_b or wMelOctoless were 1104 used in the analysis. The copy number of wMelOctoless genes is relative to that of wMelCS_b.

1107
The relative copy number of genomic WD0513 in Wolbachia-carrying stocks throughout 30 fly generations.

1108
Each generation, 5-20 females were randomly collected for egg-laying for 3-4 days and used to determine 1109 the relative copy number of WD0513 as a proxy for the Octomom copy number. The progeny of a single 1110 female was used to set up the next generation. qPCR results were normalized to that of wMelCS_b, which 1111 has a single copy of Octomom per genome.

1172
Survival of males males after a DCV (A) or buffer (B) challenge. An experimental replicate is given in Fig 6.

1173
The experiment was set up and analysed as described in Fig 6. 1174 1175 S1 Table. Number of progeny of EMS-treated females screened for over-proliferating Wolbachia 1176 variants.