The complete prophage WOcauB3 (B3gp1–B3gp46), prophage WOcauB2 (B2gp1–B2gp47), the flanking-region genes from the prophage WOVitA1 (VA1gp52–VA1gp63), and two flanking-region genes from the prophage WORiB1 (WRi_005400–WRi_005900) were used as queries in a blastp search of the NCBI non-redundant protein database and a blastn search of the NCBI nucleotide collection (nr/nt) and whole-genome shotgun contigs (wgs) databases. The output E-value (<10⁻⁵) of the searches were used as criteria for data parsing. Sequences were aligned with ClustalW in BioEdit (Hall, 1999), and the Gblocks program (ver. 0.91b) (Castresana, 2000) was used to remove poorly aligned positions.

For phylogenetic analyses, ProtTest 3 (for amino acid sequences; Darriba et al., 2011) and jModelTest 2 (for nucleotide sequences; Darriba et al., 2012) were used to determine the best evolution model based on the corrected Akaike information criterion (AICc). PhyML 3.0 (Guindon et al., 2010) and Mrbayes 3.2 (Ronquist et al., 2012) were used to build phylogenetic trees with ML and BI methods respectively. The best models chosen by ProtTest 3, LG + I + G was used to generate the ML and BI tree for B3gp45. The best model chosen by jModelTest 2, GTR+G, was used to generate the ML tree for the Wolbachia MLST phylogeny.

To visualize the general compositional features of the putative horizontally transferred genes using GC-content, a cumulative GC profile was assembled (Gao and Zhang, 2006). The cumulative GC profile can identify genomic islands or HGTs through comparison of nucleotide compositional features (Gao and Zhang, 2006). The halting parameter was set to 7, and the minimum length to segment was set to 100.

MEGA6 was used to estimate the mean synonymous divergence for each group of sequences representing potential recent horizontal transfer of WO phages, other WO phages that seem not to results from recent horizontal transfer, and their corresponding Wolbachia hosts (Tamura et al., 2013). For each group of sequences, the Nei-Gojobori method was used to calculate the synonymous rate, and variance was computed using 1000 bootstrap replicates (Nei and Gojobori, 1986).

Musca domestica and Nasonia vitripennis were used in these experiments. The N. vitripennis populations were the Hangzhou strain (from the Gongyin Ye lab, ZheJiang University) (Zhang et al., 2005) infected with Wolbachia supergroup A (Liu et al., 2014). The housefly larvae were fed bran for 5–6 days until pupation. All wasps were reared on fresh house fly pupae at 25 ± 2°C under a 14 h light cycle in an atmosphere of 50–60% relative humidity, supplemented with a piece of cotton in a soft capsule shell of 10% honey water. The adult houseflies were kept at 25 ± 2°C, but supplied with a sugar/milk powder mixture (25/75%) and water instead. Adults of N. vitripennis were initially immersed in 95% ethanol at −20°C prior to DNA extraction.

Total N. vitripennis genomic DNA was extracted from a single wasp using the DNeasy Tissue Kit (Qiagen, Hilden, Germany) following the manufacturer's recommendations and resuspended in 20 μl double-distilled sterile water. DNA purity and concentration were determined with a NanoDrop 2000 Spectrophotometer (Thermo, Madison, WI, USA), and samples of poor quality were discarded. The identity of the DNA templates was confirmed by wsp 81f and 691r primers to amplify the Wolbachia surface protein gene (Zhou et al., 1998). The PCR reactions were performed using TransTaq DNA Polymerase HiFi Fidelity (TransGen Biotech, Beijing, China) with the recommended conditions and reagents. The resulting amplicons were electrophoresed on a 1% TBE agarose gel and photographed under UV illumination. The amplified PCR products were sequenced directly with an ABI3730 capillary autosequencer (Biosune, Beijing, China) after purification with the EasyPure PCR Purification Kit (TransGen Biotech, Beijing, China). If the products could not be sequenced directly, we cloned them into the pEASY-T5 vector (TransGen Biotech, Beijing, China), and a minimum of three positive clones were sequenced due to transformation-induced mutation. Sequence editing was performed with BioEdit (Hall, 1999).

Real-time qPCR was performed with a Stratagene Mx3000p qPCR System (Stratagene, La Jolla, CA, USA) (the primers are listed in Table S3). We used real-time qPCR to quantify the DNA copies of a putative transcriptional regulator gene (VA1gp53) and an Hsp20-family heat shock protein gene (VA1gp62) from the flanking region of WOVitA1; an ank gene (VA1gp3) from phage WOVitA1; and a heat-shock protein 60 gene (groEL) (Bordenstein et al., 2006) and cell division gene (ftsZ) from wVitA vs. prepared standard solutions. The amplified PCR products were sequenced directly to confirm the gene identity. A standard 10-fold dilution series from 10⁷ to 10³ copies were prepared and used to calculate the copy numbers of the genes. The genes's amplification efficiency in our experiments are 96.6–104.4%. Also, each melting curve showed that the primers amplify a single product.

The average copy number of the integrated phage was compared with the expected number and the difference was analyzed statistically with a two-tailed t-test (SAS Institute, Cary, NC, USA). With a single lysogenic copy of WOVitA1, the expected WOVitA1 number should always equal (no lytic activity) or exceed (with lytic activity producing multiple phage virions) the wVitA copy number. We normalized the small plate effects in real-time qPCR experiments as described previously (Wang et al., 2014). The compared percent nucleotide identity was analyzed by an Mann–Whitney U two-tailed test using Origin8.0.

De novo nucleotide sequences were deposited in GenBank under accession numbers KP966832–KP966840.

Previous searches of public databases suggested the transfer of bacterial spvB gene (B3gp45) between an unrelated bacterial genetic lineage and Wolbachia wCauB by WO (Tanaka et al., 2009). In the present study, we expand this finding by conducting a homology search for the complete phage WOcauB3 genes in all of the 32 reported Wolbachia genomes (information on all Wolbachia genomes in this study is listed in Table 1). Of the 32 genomes tested, we detected that only two, wCauB and wNo, have uniquely encoded spvB and nearby gene, B3gp46, which encodes a hypothetical protein that is packaged in phage WOcauB3 particles (Tanaka et al., 2009; Figure 2A). The Wolbachia wNo is a strain from supergroup B infecting D. simulans (Ellegaard et al., 2013).

To further trace the transfer trajectory of spvB and B3gp46, we compared the divergence between prophages WOcauB3 and WONo1–4 and between their Wolbachia hosts, wCauB and wNo. The Wolbachia strain wNo harbors four WO phages, WONo1–4 (Table 1).

Several lines of evidence support the possibility that WOcauB3 was transferred between wCauB and wNo and mediated the transfer of both genes.

In addition to WOcauB3, there is also evidence to support that WOcauB2 (from the same Wolbachia strain, wCauB) has experienced a transmission event. The transmission likely mediated the horizontal transfer of two associated ank genes.

In addition to WOcauB2 and WOcauB3, Wolbachia wVitA possesses another active phage, WOVitA1, which has been transferred to Wolbachia wVitB in N. vitripennis (the phage WOVitA1 in wVitB is WOVitB1; Kent et al., 2011b). It is interesting to note that both WOVitA1 and WOVitB not only have nearly identical phage regions (the WOVitA1 region is VA1gp1–VA1gp51) but also highly similar flanking regions with bacterial genes (VA1gp52–VA1gp63 in WOVitA1; Kent et al., 2011b). Here, by analyzing the attachment site (att) of WO phages (attP) and the density correlation of WO with its Wolbachia host, we demonstrate that the active phage WOVitA1 can mediate the transfer of its flanking bacterial region.

Inactive WO phages may also have been involved in gene transfer events. In eight Wolbachia genomes, we detected a conserved bacterial region extending over 20 kb that is highly homologous (>70% nt identity) to regions in the bacterial plasmids of Rickettsia buchneri sp. nov. and Rickettsia helvetica (Ishmael et al., 2009). These bacteria infect Ixodes scapularis (Kurtti et al., 2015) and Ixodes ricinus (Dong et al., 2012) ticks respectively. Interestingly, except in wNo and wAlb, this region in each of the other Wolbachia genomes is inserted in or near the phage WO, and some of the associated WO phages are degenerate (Figure 6A). We also detect homologs of some of the genes from this region in wBol1 (from Hypolimnas bolina; Duplouy et al., 2013), wVitB, wWil (from Drosophila willistoni) (Craig Venter Institute), wCoc (from Dactylopius coccus; Campana et al., 2015), and wRec (from Drosophila recens; Metcalf et al., 2014). However, all of the homologs are in scaffolds with small sizes, which prevents us from obtaining their flanking regions; therefore, we did not further analyze them. In public databases, there are no closely related orthologs of this bacterial region except in the Wolbachia strains and the two Rickettsia strains as mentioned above. There are at least three possible explanations for the distribution pattern of this region. First, this bacterial region is one of the modules of phage WO. Second, this conserved bacterial region has been frequently and independently inserted into Wolbachia at the same phage WO location, where there exists an active cloning location. Third, phage WO or plasmid mediates the transfer of this bacterial region among different Wolbachia strains, between Rickettsia from I. scapularis and I. ricinus, or between the two bacterial genera. However, the average nucleotide identity of the genes in this region (97.00%) is significantly higher than that of the other Wolbachia genes (93.28%; MWU, two-tailed, p < 0.01; Figure 6B and Table S4) and is also significantly higher than that of WORiB1, WORiB2, WOSol, WOMelB, WOSuz1, and WOAuB (88.67%; MWU, two-tailed, p < 0.01; Figure 6B). Furthermore, the synonymous mutation rate of the genes in the region from wRi, wCs, wMel, wSuz, and wAu (0.02) is significantly lower than the average synonymous mutation rate of the associated phage WO genes from WORiB1, WORiB2, WOSol, WOMelB, WOSuz1, and WOAuB (0.08; MWU, two-tailed, p < 0.01; Figure 6C). The average nucleotide identity and synonymous mutation rate analyses indicate this bacterial region is not phage WO module. Additionally, in the cumulative GC profile of phage WOMelB (Figure S5), the segmentation point, including the putative HGTs, also demonstrates phage WO and this bacterial region are from different DNA sources. Based on these analyses, we exclude the first explanation. However, we cannot exclude the second possibility that this bacterial region has been frequently and independently inserted into Wolbachia at the same phage WO location. Considering that phage and plasmid are two of the most common genetic vectors in nature (Syvanen, 1994; Canchaya et al., 2003), it is thus a more parsimonious explanation that phage WO or Rickettsia plasmid may have recently mediated the transfer of these genes.

Within the transfered region, there are 11 conserved genes (WRi_005730–WRi_005830; Figure 6A). These genes encode an NAD-dependent epimerase/dehydratase family protein, a glycosyltransferase, two putative L-allo-threonine aldolases, an ABC transporter permease, a GlpT/PgpT/UhpT transporter family protein, a UDP-glucose 6-dehydrogenase, and three conserved hypothetical proteins. Many of these proteins play a role in the synthesis and degradation of surface polysaccharides, which could alter the ability of different Wolbachia strains to interact with eukaryotic hosts (Ishmael et al., 2009).

Moreover, this transfered region is located adjacent to a conserved gene encoding a SNF2-family helicase, which was detected in a phage region via genome analysis (Ishmael et al., 2009). In eukaryotes, this gene may function in processes including transcriptional regulation, the maintenance of chromosome stability during mitosis and the processing of DNA damage (Eisen et al., 1995).

Transfer events mediated by phage WO can shape the genome composition of Wolbachia. For example, ANKs are rare in bacteria but common in eukaryotes and viruses (Bork, 1993; Li et al., 2006), while these genes are overrepresented in Wolbachia bacteria. For example, there are typically only 1–3 ank genes in the α-Proteobacteria (Andersson et al., 1998; Caturegli et al., 2000), but there are 60 ank genes in Wolbachia wPip from C. pipiens (the largest number of ank genes in any sequenced bacterial genome; Klasson et al., 2008), 35 in wRi (Klasson et al., 2009), and 23 in wMel (Wu et al., 2004). In this study, we demonstrate that WO phage particles can package and mediate the transfer of “extra” ank genes into Wolbachia genomes, e.g., B2gp46 and B2gp47 in WOcauB2, both of which are of non-Wolbachia origin. Moreover, WO can also mediate the transfer of “extra” ank genes between different Wolbachia strains. These results indicate that the WO-mediated transfer of “extra” ank genes may be a partial explanation for the abundance of ank genes in Wolbachia compared to other closely related bacteria.

The discovery of these horizontal transfer events raises the question of whether transferred genes play a role in Wolbachia or their eukaryotic hosts. Prophage-encoded virulence factors are important for a number of bacterial species, and these genes can increase pathogenicity or result in the emergence of new pathogens (Canchaya et al., 2003; Brüssow et al., 2004). This phenomenon has been recognized for the toxins of Vibrio cholerae (Waldor and Mekalanos, 1996), Streptococcus pyogenes (Broudy et al., 2002), and Hamiltonella defensa (Moran et al., 2005; Oliver et al., 2009), all of which are phage-encoded. Here, the presence of the transferred spvB motif gene in phage WO particles, and the role of this gene in WO infection of Wolbachia and the corresponding eukaryotic hosts needs further study.

Wolbachia mediated mosquito-borne disease control is a hot topic (Dobson et al., 2016; Loreto and Wallau, 2016; O'Neill, 2016; Waltz, 2016). Caged and open-field experiments showed that the wMel Wolbachia strain is able to block dengue transmission (Walker et al., 2011). However, there is also a potential risk that the Wolbachia strains, along with phage WO and other genes, may be transferred to other insects (Loreto and Wallau, 2016). Here we demonstate that phage WO can mediate HGT among different Wolbachia strains. Thus, future studies should also evaluate the biosafety of this phage vector when utilizing Wolbachia-infected mosquitos.

The phenomenon of eukaryotic host reproductive manipulation by Wolbachia is compelling, but the underlying mechanism still remains poorly characterized due to the lack of robust tools for transforming Wolbachia (Werren, 1997; LePage and Bordenstein, 2013). The phage WO has been proposed as the only potential transformation tool for Wolbachia (Fujii et al., 2004; Metcalf and Bordenstein, 2012). However, it remains unclear whether phage WO can be successfully used as such a tool, and there is little research concerning this issue (LePage and Bordenstein, 2013). Here, we show that phage WO can mediate gene transfer; the active phage WOVitA1 has typical characteristics of the core sequences. Additionally, the 3′ end of WO prophages sites might be used as multiple-cloning sites. All of these results further support that phage WO has the potential to be utilized as a genetic vector for the study of Wolbachia.

In the evolution of intracellular endosymbionts, genome reduction is the predominant trend differentiating endosymbionts from free-living bacteria. Additionally, intracellular endosymbionts are strictly constrained to living inside host-derived cells: their effective population size is reduced, which renders selection less efficient; they have limited opportunities to come into contact with other unrelated bacteria and have little chance to exchange genetic material; the stable and rich nutrients of the intracellular environment remove selection constraints on genes (like mobile DNA) that are no longer strictly required (Bordenstein and Reznikoff, 2005; Moya et al., 2008; Moran and Bennett, 2014). However, phage WO is widespread among Wolbachia genomes (present in about 89%; Bordenstein and Wernegreen, 2004) and even can comprise more than 20% of mobile DNA genes in Wolbachia (Chafee et al., 2010). What's more, the active mobile elements located within the genomes of endosymbionts can still mediate the deletion and insertion of genetic components at different locations in the genome. Based on several lines of evidence, the present study shows that the phage WO could mediate HGT between different Wolbachia strains of genes from Wolbachia and unrelated bacterial lineages, which showes Wolbachia genomes are not stable and might gain new genes by phage WO.

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.