RNA virome diversity and Wolbachia infection in individual Drosophila simulans flies

The endosymbiont bacteria of the genus Wolbachia are associated with multiple mutualistic effects on insect biology, including nutritional and antiviral properties. Members of the genus Wolbachia naturally occur in fly species of the genus Drosophila, providing an operational model host for studying how virome composition may be affected by its presence. Drosophila simulans populations can carry a variety of strains of members of the genus Wolbachia , with the wAu strain associated with strong antiviral protection under experimental conditions. We used D. simulans sampled from the Perth Hills, Western Australia, to investigate the potential virus protective effect of the wAu strain of Wolbachia on individual wild-caught flies. Our data revealed no appreciable variation in virus composition and abundance between individuals infected or uninfected with Wolbachia associated with the presence or absence of wAu. However, it remains unclear whether wAu might affect viral infection and host survival by increasing tolerance rather than inducing complete resistance. These data also provide new insights into the natural virome diversity of D. simulans. Despite the small number of individuals sampled, we identified a repertoire of RNA viruses, including nora virus, galbut virus, thika virus and La Jolla virus, that have been identified in other species of the genus Drosophila. Chaq virus-like sequences associated with galbut virus were also detected. In addition, we identified five novel viruses from the families Reoviridae, Tombusviridae, Mitoviridae and Bunyaviridae. Overall, this study highlights the complex interaction between Wolbachia and RNA virus infections and provides a baseline description of the natural virome of D. simulans.


INTRODUCTION
The members of the alpha-proteobacterium genus Wolbachia (order Rickettsiales) are widespread endosymbionts of arthropods and nematodes (i.e. filarial and plant-parasitic nematodes) that can establish interactions with their hosts ranging from parasitic to mutualistic [1,2]. The genetic diversity of the genus Wolbachia is substantial and currently represented by 11 distinctive supergroups (denoted A-J), although the majority of strains of members of the genus Wolbachia belong to supergroups A and B [3] that are estimated to have diverged around 50 million years ago [4]. Although these bacteria are commonly found in reproductive tissues and the germline of their hosts, they have also been found in somatic tissues, such as the brain, salivary glands and gut [5][6][7][8][9], such that understanding infection dynamics in detail is not a trivial matter [7]. Members of the genus Wolbachia primarily spread OPEN ACCESS by vertical inheritance through transovarian transmission. However, the presence of members of the genus Wolbachia in a diverse range of host species indicates that horizontal transmission, probably through antagonistic interactions (i.e. herbivory, parasitism and predation), also contributes to the dissemination of the bacteria in nature [4,10].
The occurrence of bacteria of the genus Wolbachia in insects is often associated with their ability to manipulate host reproductive mechanisms and induce a range of alterations, including parthenogenesis, feminization, cytoplasmic incompatibility and sex-ratio distortion [11]. Among these, cytoplasmic incompatibility is the most common phenotypic effect, and as such represents an appealing approach for vector population control. In this case, embryonic lethality is contingent on the infection status and the strain type harboured by males and females [2]. In addition, the study of Wolbachia-host interactions has revealed a variety of mutualistic effects on host biology [1,12]. For instance, in filarial nematodes and the parasitoid wasp Asobara tabida, the presence of some strains of members of the genus Wolbachia has been positively associated with developmental processes, fertility and host viability [12][13][14]. Furthermore, nutritional mutualism between members of the genus Wolbachia and the bedbug Cimex lectularius as well as Wolbachia-infected planthoppers, has been suggested as a means to explain B vitamin supplementation [15][16][17].
Arguably the most important outcome of Wolbachia infection in insects is its potential for virus-blocking, which also provides a basis for intervention strategies based on the control of arbovirus transmission. This seemingly antiviral effect of members of the genus Wolbachia has been well documented in some species of insects, including flies and mosquitoes. A striking example involves the transinfection of Aedes aegypti mosquitoes with a strain of the genus Wolbachia infecting Drosophila melanogaster (wMel). A. aegypti is the primary vector of a number of important arboviruses, including members of the species Dengue virus, Zika virus and Chikungunya virus, and the establishment of the wMel strain in wild mosquito populations represents a powerful and promising approach to decrease virus transmission [18,19]. Although the underlying mechanisms remain to be fully determined, it has been suggested that members of the genus Wolbachia can modify the host environment or boost basal immunity to viruses by pre-stimulating the immune response of their hosts [20]. Potential antiviral mechanisms affected by members of the genus Wolbachia include gene expression of the Toll pathway, RNA interference, and modification of the host oxidative environment that probably trigger an antiviral immune response and hence limit infection [20][21][22].
Unlike A. aegypti mosquitoes, members of the genus Wolbachia naturally occur in species of the genus Drosophila, providing a valuable model system to study Wolbachia-related virus protection [23,24]. Natural populations of species of the genus Drosophila can carry a diverse array of insectspecific viruses belonging to the families Picornaviridae, Dicistroviridae, Bunyaviridae, Reoviridae and Iflaviridae amongst others [25]. The co-occurrence of members of the genus Wolbachia in D. melanogaster has been associated with increased survival and different levels of resistance to laboratory viral infections in fly stocks under experimental conditions [23,26]. For example, Wolbachia-infected flies containing the dicistrovirus Drosophila C virus (DCV) showed a delay in mortality compared with Wolbachia-free flies [26]. In contrast, other studies found no or limited effect of members of the genus Wolbachia on viral protection, as well as on virus prevalence and abundance in field-collected flies [25,27]. Such contrasting data emphasize the need for further research efforts to characterize the effect of strains of the genus Wolbachia on virus composition in members of the genus Drosophila in nature.
Although the origin of D. simulans is thought to have been in East Africa or Madagascar, this species now has a cosmopolitan distribution [28]. In Australia, D. simulans has been recorded along both the east and west coasts as well as Tasmania, with the earliest record dating to 1956 [29]. Human mobility and human-mediated activities have been associated with the introduction and spread of both D. simulans and Wolbachia into Australia, where wild fly populations occur near human settlements, feeding and breeding on a variety of horticultural crops [30,31]. Several strains of members of the genus Wolbachia from supergroups A and B can naturally occur in populations of D. simulans (e.g. wAu, wRi, wHa, wMa and wNo) [32,33]. Of these, wAu is associated with strong antiviral protection against Flock House virus (FHV) (Nodaviridae) and DCV (Dicistroviridae) under experimental conditions [32]. The wAu infection in Australia was one of the first infections by members of the genus Wolbachia identified as showing no cytoplasmic incompatibility, despite being widespread at a low to intermediate frequency [34]. wAu increased in frequency along the east coast of Australia until it was replaced by wRi that exhibits cytoplasmic incompatibility. However, unlike wAu, wRi has not yet reached the Australian west coast [30]. In this study, we used a meta-transcriptomic (i.e. RNA shotgun sequencing) approach to determine the virome diversity of individual field-collected D. simulans flies from Western Australia, and investigated how this virome diversity might be affected by the presence of the wAu strain of Wolbachia.

D. simulans collection and taxonomic identification
Flies used for the virus work performed here were collected at Raeburn Orchards in the Perth Hills in Western Australia (latitude −32.1036°, longitude 116.0695°), in July 2018 using banana bait. The frequency of members of the genus Wolbachia at two other locations in the area (Roleystone, latitude −32.1396°, longitude 116.0701°; Cannington, latitude −32.0243°, longitude 115.9363°) was also established with additional samples. Taxonomic identification to the species level was conducted on the basis of the morphology of reproductive traits of males and via DNA barcoding (cox1 gene marker). Field-collected flies were maintained at 19 °C under standard laboratory conditions until F1 offspring were raised. Parental and F1 generations were then stored at −80 °C until molecular processing.

Wolbachia detection
Wolbachia infection of field females was determined using F1 offspring from each field female. Note that wAu is transmitted at 100 % from field females to the F1 laboratory generation [34]. DNA extraction from heads was performed using Chelex 100 Resin (Bio-Rad Laboratories,) [35] as adapted y Shi et al. [27]. Screening of natural infection with members of the genus Wolbachia was conducted using a real-time PCR/ high-resolution melt assay (RT/HRM) and strain-specific primers targeting a 340 bp region of the surface protein of Wolbachia (wsp) gene for wRi and wAu strains. The assay was run following the protocol of Kriesner et al. [30]. In addition, reads were mapped to reference Wolbachia wsp gene sequences for wRi (CP001391.1) and wAu (LK055284.1) with BBMap v.37.98 (minid=0.95) (available at https:// sourceforge. net/ projects/ bbmap/).

RNA extraction and meta-transcriptome sequencing
We screened a total of 16 individual flies to assess the effect of infection with members of the genus Wolbachia on virome composition in D. simulans. Specimens were rinsed three times in RNA and DNA-free PBS solution (GIBCO). Total RNA from individual flies was extracted using the RNeasy Plus Mini Kit (Qiagen) following the manufacturer's instructions. RNA-seq libraries were constructed using a TruSeq total RNA Library Preparation Kit (Illumina). Host ribosomal depletion was performed using a Ribo-Zero Gold rRNA Removal Kit (Human/Mouse/Rat) (Illumina) and paired-end transcriptome sequencing was performed on the HiSeq2500 platform (Illumina). Libraries from Wolbachia-negative and -positive infected flies were run in two separate lanes.

De novo meta-transcriptome assembly and viral genome annotation
The overall quality assessment of reads was conducted in FastQC and Trimmomatic [36]. A de novo assembly of RNA-Seq data was performed using MEGAHIT v.1.1.3, with default parameters [37]. Assembled contigs were then annotated through comparisons against the NCBI nonredundant (NCBI-nr) database using DIAMOND v2.0.4 [38], with a cutoff e-value <1e −05 . To identify protein-encoding sequences, open reading frames (ORFs) were predicted in positive and reverse-complement strands, with a minimum length of 600 nt between two stop codons using the GetOrf program (EMBOSS) [39]. Functional annotation was carried out using InterProScan v5.39-77.0 [40], and the HMMer software (http:// hmmer. org/) was used to perform sequence-profile searches against the Pfam HMM database. To expand the de novo assembled contigs of known viruses, the reads were mapped against reference genomic sequences. Provisional virus names were derived from geographic locations in the Perth Hills, Western Australia.

Estimates of viral abundance
Viral abundance was assessed using the number of reads per million (RPM). This metric quantifies the number of reads per million mapped to a given contig assembly over the total number of reads. RPM values lower than 0.1 % of the highest count for each virus across samples were presumed to be index-hopping artefacts and excluded from the remaining analyses [41]. To compare abundance levels, reads were mapped to reference ribosomal and mitochondrial genes from Wolbachia (16S and cox1) and D. simulans (rpl32 and cox1), as well as against all the RNA viruses identified in the annotation analyses. Mapping was performed using BBMap v.37.98 (minid=0.95) (available at https:// sourceforge. net/ projects/ bbmap/).

Sequence alignment and phylogenetic analysis
RNA viral sequences identified in D. simulans were compared with homologous reference sequences retrieved from the NCBI GenBank database and aligned with MAFF v7.450 (E-INS-I algorithm) [42]. Phylogenetic trees for these data were then inferred using sequences of the conserved RNAdependent RNA polymerase (RdRp) gene. To this end, both the best-fit model of amino acid substitution and phylogenetic relationships were estimated using the maximum likelihood (ML) [43] approach implemented in IQ-TREE v1.6.12 [44]. Nodal support was estimated combining the SH-like approximate likelihood ratio test (SH-aLRT) and the Ultrafast Bootstrap Approximation (UFboot) [45]. Redundant contigs with over 99 % amino acid similarity were excluded. For those libraries containing viruses that were unlikely to be associated with Drosophila, taxonomic profiling and read mapping to components of the fly microbiome and diet were conducted using the CCMetagen software (default settings) [46] and BBMap v.37.98 (minid=0.95).

Statistical analysis
The assumption of data normality was assessed by visual inspection and using Kolmogorov-Smirnov (K-S) and Shapiro-Wilk's tests. As the data was not normally distributed, a Mann-Whitney-Wilcoxon test was used to compare the RNA virome composition with respect to the presence/ absence of members of the genus Wolbachia. Comparisons were made using raw and log-transformed data corresponding to RPM values (i.e. viral abundance) for each library. All analyses were performed using R software package rstatix (available at https:// rpkgs. datanovia. com/ rstatix/).

RESULTS
A total of 272 female flies were wild-caught in the Perth Hills, Western Australia and tested for Wolbachia infection through their F1s. The overall prevalence of members of the genus Wolbachia was 63.6 % (173 out of 272), with frequencies at the three sampled locations varying from 54.8 % (Raeburn Orchard, n=73) to 63.8 % (Roleystone, n=130) and 72.5 % (Cannington, n=69). We randomly selected 16 flies from the Raeburn Orchard field females for individual sequencing and RNA virus screening, representing eight Wolbachia-positive and eight Wolbachia-negative specimens.
We identified the Wolbachia strain in D. simulans using sequence-specific primers targeting the wsp gene. We further confirmed the occurrence of Wolbachia by mapping the reads back to the wRi and wAu wsp genes. Most of the Wolbachiainfected flies showed a median coverage >100 reads, number of mapping reads >40, and coverage percentage >90 % to the reference wAu strain, confirming that infected flies harbour wAu rather than wRi. No reads mapped to the wsp gene for library RAPP88 (Table S1, available in the online version of this article) despite the positive infection status determined using a Wolbachia-specific qPCR assay.
For comparison of virus diversity among libraries we mapped the reads of each library to stably expressed genes: 16S and cox1 in Wolbachia and rpl32 and cox1 in D. simulans. This provided an internal control to identify any effect on viral abundance due to potential biases introduced during RNA extraction or library preparation. Although, as expected, there was moderate variation in the abundance values, expression levels of reference maker genes were relatively stable across libraries in both Wolbachia and D. simulans (Fig. 1).
Some of the newly discovered RNA viruses identified here were probably infecting hosts other than D. simulans, and hence might be associated with the fly diet or microbiome. Specifically, these viruses were closely related to Phytomonas sp. TCC231 leishbunyavirus 1 (in the case of Raeburn bunyalike virus), Leptomonas pyrrhocoris RNA virus (Cannin tombus-like virus) and two mito-like viruses (Araluen  Table  S3), that are associated with trypanosomatid protozoans and fungal hosts, respectively. In addition, taxonomic composition analyses as well as read mapping to common components of Drosophila diet or microbiome revealed that 0.1 and 0.05 % of all non-rRNA reads mapped to fungi and trypanosomatids (cox1 gene marker), respectively. Hence, multiple microorganisms were present within individual fly libraries which may explain the occurrence of viruses not directly associated with Drosophila (files available at https:// doi. org/ 10. 6084/ m9. figshare. c. 5466690). In contrast, Lesley reo-like virus is likely to be a bona fide arthropod virus since it grouped with viruses previously detected in odonates and mosquitoes. In addition, it exhibited only approximately 24 % nucleotide similarity with those reoviruses previously reported to be contaminants in Drosophila cell culture [25]. This indicates that Lesley reo-like virus is not a component of known contaminants and more is likely to be part of the natural D. simulans virome. The five newly identified viruses in this study corresponded to full or nearly complete genomes (see below). However, for the majority of the known viruses of Drosophila we only were able to identify ORFs encoding the RdRp: the exceptions were La Jolla virus and thika virus, for which we also predicted structural components corresponding to coat and capsid proteins.
We next characterized the virome profile present in D. simulans in relation to the wAu infection status (Fig. 2, Tables 1 and S4). Accordingly, we identified a slightly higher number (n=9) of viruses in Wolbachia-negative flies compared with Wolbachia-positive flies (n=6). Among these, galbut virus, nora virus, thika virus, as well as three novel viruses identified in this study -Raeburn bunya-like virus, Araluen mito-like virus and Cannin tombus-like virus -were present in D. simulans regardless of Wolbachia infection. Likewise, 'chaq virus-like' sequences were observed co-occurring with galbut virus in the two groups of D. simulans. In contrast, La Jolla virus, as well as the novel Carmel mito-like virus and Lesley reo-like virus, were only found in wAu-negative flies. Overall, assembled viral contigs displayed high sequence similarity at the nucleotide and amino acid level within and between libraries and regardless of the presence or absence of members of the genus Wolbachia (Table S3).
We also assessed the potential effect of infection with members of the genus Wolbachia on the abundance of RNA viruses present in wAu-infected and wAu-uninfected flies. Overall, the number of non-rRNA reads represented approximately 50 % of the total number of reads (n=743 389 696 pairend reads) (Fig. S1). Furthermore, the RPM values among viruses infecting Wolbachia-negative and -positive infected flies was highly heterogeneous, ranging from 47 to 232 346 and 5 to 37 688 virus RPM, respectively. With the exception of thika virus, viruses present in both wAu-positive and wAu-negative flies were 1.87-40.17-fold more abundant in  the wAu-negative individuals than wAu-positive D. simulans.
In contrast, the abundance of thika virus was 0.39-fold higher in the Wolbachia-positive flies (Fig. 3, Table S2). However, despite this variation in virus abundance levels between groups, there was a non-significant difference between wAunegative and wAu-positive D. simulans (Mann-Whitney-Wilcoxon test; Fig. 5). In the case of the viruses only detected in the wAu-negative flies, La Jolla virus was present in a single library in moderate abundance (RPM=378), whilst the newly discovered Lesley reo-like virus was detected in four out of eight libraries (RPM=3360-8749) (Table S2). Although an interesting result, the limited sample size (n=16) means that these observations should be taken with caution and that larger sample sizes are needed for corroboration.

DISCUSSION
The occurrence and spread of Wolbachia infection has been widely documented in natural populations of Drosophila [10,30,47]. Indeed, D. simulans is commonly used as an experimental model to investigate the interactions within the tripartite Drosophila-Wolbachia-virus system. In Australia, D. simulans can be naturally infected with two strains of Wolbachia from supergroup A wAu and wRi. While wRi has been gradually displacing wAu in eastern Australia, reflected in the changing infection frequencies in surveyed populations since 2004, D. simulans from the west coast of Australia only harbour wAu [30]. A simple and plausible explanation for this difference is the geographic separation of D. simulans populations inhabiting the east and west coasts of Australia and the challenging environmental conditions posed by the intervening desert [30].
We corroborated the presence of Wolbachia infection across samples by identifying the wsp, 16S and cox1 marker genes.
The lack of reads mapping to the library RAPP88 might reflect either low levels of wsp RNA molecules present in the input for library preparation or high variability compared with the reference sequence. Although Wolbachia density was not experimentally assessed, the similar levels of 16S and cox-1 abundance across libraries indicate no appreciable biases in the library preparation and RNA sequencing steps.
Estimates from previous surveys indicated that the frequency of the wAu strain in Western Australia exceeded 50 % in D. simulans [30]. This is consistent with the data provided here and indicates that Wolbachia might be present in a significant proportion of the natural fly population, at least around Perth.
Although wAu does not cause cytoplasmic incompatibility, its spread is suggested to confer fitness advantages (increased survival and/or reproduction) to the host, including antiviral protection [48,49], that might favour its spread and prevent the bacteria from being eliminated from D. simulans populations [30,50].  with wMel in eastern Australia [27], in which virus protection was not observed regardless of the Wolbachia infection status and Wolbachia density. Even so, the absence of a significant association between wAu infection and virus diversity does not necessarily translate into a homogeneous effect of wAu on the different viruses identified here. For example, it is plausible that the restricted presence of La Jolla virus and the newly identified Lesley reo-like virus in Wolbachia-free flies could reflect some effect of antiviral protection in D. simulans [27,51]. Indeed, contrasting results were observed in D. melanogaster, where La Jolla virus was widely distributed across different libraries [27]. Although this might provide insights into wAu-virus interactions, studies based on larger sample sizes are clearly needed to determine whether the apparent association between La Jolla virus and Wolbachia-uninfected flies observed here is an artefact due to small sample sizes. Indeed, it is notable that La Jolla virus was so rarely detected in the D. simulans flies studied here.
It has previously been shown that the wAu strain of Wolbachia has a protective role against virus infection in D. simulans when flies are challenged with FHV and DCV in a laboratory setting [24,32]. Moreover, the wAu strain is protective against the dengue (DENV) and zika (ZIKV) viruses in Aedes aegypti mosquitoes [52]. Although our observation of an apparent lack of Wolbachia-mediated antiviral protection contrasts with those obtained previously, it is likely that differences may depend on Wolbachiahost species combinations and natural or artificial viral infections, which may also explain the contrasting results for La Jolla virus. Indeed, most of the available studies have documented the antiviral effect in transinfected insect hosts with non-natural strains of Wolbachia and viruses under laboratory conditions, as opposed to the study of the natural virome undertaken here.
It is noteworthy that ecological variables such as temperature might affect Wolbachia-virus-host interactions. Here, we collected flies during the Western Australian winter (mean temperature of 21 °C daylight time) and the specimens were maintained at 19 °C under laboratory conditions. Lower temperatures have been associated with an increase in viral resistance against DCV in D. melanogaster infected with wMel and wMelCS [53]. Similarly, variations in host developmental temperature have been associated with differences in Wolbachia-mediated virus blocking in natural populations [53]. In this context, flies developed at lower temperature (18 °C) exhibited a reduction in Wolbachia-conferred antiviral protection. On the other hand, the presence of Wolbachia has been suggested to influence host temperature preferences. For instance, wRi and wHo strains seem to manipulate D. simulans flies to seek cooler temperatures [54]. Although the effect of temperature on wAu and D. simulans need to be tested, this indicates that the results observed here as well as a protective scenario might be temperature-dependent. This highlights the importance of careful future studies of the interactions within the host-virus-Wolbachia system along with environmental factors in natural populations [55][56][57]. As well as the small sample size, an important caveat of our work is that we explored the Wolbachia-mediated virus protection in terms of virus abundance levels reflected in RPM values. This provides insights into virus resistance, but not into tolerance or host survival. Thus, it is still possible that Wolbachia is increasing tolerance to virus infection as has been documented for DCV [32]. In addition, although we were not able to assess Wolbachia density, results of previous studies have indicated that wAu is maintained at high-density in D. simulans and has a role in virus blocking [58]. Further research is clearly needed to assess these features in natural populations to determine any link with antiviral protection.
Collectively, comparisons of the virome composition in wAu infected and uninfected D. simulans revealed the presence of natural and relatively highly abundant Drosophila-associated viruses in both groups [25,27,59]. Consistently with results from previous studies, we noted the co-occurrence of chaq virus-like sequences and galbut virus, supporting the idea that chaq virus might be part of a satellite-helper virus system or an additional segment associated with galbut virus [27,60]. In addition to insectassociated viruses, we identified viruses that are likely to infect other hosts and hence were likely to be associated with components of the diet or microbiome of D. simulans [61]. For instance, novel viruses from the families Tombusviridae and Bunyaviridae were related to viruses in trypanosomatid protozoa (Leptomonas and Leishmania). Similarly, given their normal host range distribution, the novel viruses from the family Narnaviridae might be associated with fungal hosts. Evidence of trypanosomatids and fungi in the gut of several species of Drosophila, with effects on larvae eclosion and pupation times [61,62], has been reported. This, in turn, highlights the extent to which Australian D. simulans can be parasitized in nature [62][63][64][65][66].
In sum, we provide a preliminary framework for assessing the effect of the wAu strain on the virome of D. simulans, using a meta-transcriptomic analysis of individual wAu-infected and uninfected flies. In doing so we identified Drosophila-associated viruses along with five novel viruses likely to be associated with fly diet or microbiome. Although our sample size is small, we saw no detectable Wolbachia-associated antiviral effect on virus composition and abundance, although the approach taken prevented us from drawing conclusions on virus tolerance. Further research employing larger sample sizes over broad spatial scales, including different Wolbachia-Drosophila combinations, will enable a more nuanced understanding of Wolbachia-virus dynamics in wild Drosophila populations.