Population Genomics of Infectious and Integrated Wolbachia pipientis Genomes in Drosophila ananassae

Abstract

Coevolution between Drosophila and its endosymbiont Wolbachia pipientis has many intriguing aspects. For example, Drosophila ananassae hosts two forms of W. pipientis genomes: One being the infectious bacterial genome and the other integrated into the host nuclear genome. Here, we characterize the infectious and integrated genomes of W. pipientis infecting D. ananassae (wAna), by genome sequencing 15 strains of D. ananassae that have either the infectious or integrated wAna genomes. Results indicate evolutionarily stable maternal transmission for the infectious wAna genome suggesting a relatively long-term coevolution with its host. In contrast, the integrated wAna genome showed pseudogene-like characteristics accumulating many variants that are predicted to have deleterious effects if present in an infectious bacterial genome. Phylogenomic analysis of sequence variation together with genotyping by polymerase chain reaction of large structural variations indicated several wAna variants among the eight infectious wAna genomes. In contrast, only a single wAna variant was found among the seven integrated wAna genomes examined in lines from Africa, south Asia, and south Pacific islands suggesting that the integration occurred once from a single infectious wAna genome and then spread geographically. Further analysis revealed that for all D. ananassae we examined with the integrated wAna genomes, the majority of the integrated wAna genomic regions is represented in at least two copies suggesting a double integration or single integration followed by an integrated genome duplication. The possible evolutionary mechanism underlying the widespread geographical presence of the duplicate integration of the wAna genome is an intriguing question remaining to be answered.

Keywords: endosymbiont; horizontal gene transfer; lateral gene transfer; whole-genome sequencing.

Fig. 1.—

Geographic origin of samples analyzed in this study. Drosophila ananassae strains with wAna^INF genomes are indicated with black colored stars and the localities are Cebu (Cebu, Philippines), GB1 (Mauritius), HNL0501 (Oahu, USA), KMJ1 (Kumejima, Japan), OGS-98K1 (Ogasawara, Japan), RC102 (Rwanda), TBU136 (Tonga), and VAV150 (Vava’u, Tonga). Drosophila ananassae strains with wAna^ITG genomes are indicated with red colored stars and the localities are BKK13 (Bangkok, Thailand), D38 (Coimbatore, India), EZ104 (Ethiopia), PNP1 (Phnom Pen, Cambodia), TB43 (Trinity Beach, Australia), TBU3 (Tonga), and TI8 (Thursday Island, Australia).

Fig. 2.—

Genome-wide read coverage for the eight wAna^INF genomes. Per-site read depth normalized against the D. ananassae nuclear genome (scaffold 13340) coverage is shown against the reference wRi genome coordinate. Note that regions of tRNA and rRNA, comprising 0.48% of the genome, were expected to have spurious reads from other endosymbiotic bacteria, leading to spikes of read coverage and were ignored from subsequent analysis (see text).

Fig. 3.—

Genome-wide read coverage for the seven wAna^ITG genomes. Per-site read depth normalized against the D. ananassae nuclear genome (scaffold 13340) coverage is shown against the reference wRi genome coordinate. Note that regions of tRNA and rRNA, comprising 0.48% of the genome, were expected to have spurious reads from other endosymbiotic bacteria, leading to spikes of read coverage and were ignored from subsequent analysis (see text).

Fig. 4.—

Computational and qPCR estimates of copy number for the LOW and HIGH regions of the wAna^ITG genome. Relative copy number of LOW and HIGH region is compared with D. ananassae rp49. The x axis represents read depth of LOW and HIGH region divided by read depth of rp49. The y axis represents relative concentration of LOW and HIGH region divided by relative concentration of rp49.

Fig. 5.—

PCR amplifications of large structural variations found in (A) wAna^INF and (B) wAna^ITG samples. PCR results for the nine large structural variations are shown according to the sample. KMJ1 TET, tetracycline-treated KMJ1 strain that is cured of W. pipientis infection and does not have evidence of the integration; Negative, negative control using water for the PCR reaction. Note SV8 and SV9 are large deletions that exist in all wAna^INF and wAna^ITG samples.

Fig. 6.—

Maximum-likelihood phylogeny of the (A) host D. ananassae mtDNA and (B) wAna genomes. Bootstrap values of greater than 95% are shown on the nodes of the phylogeny. Drosophila ananassae strains with wAna^INF are indicated with black colored branches, whereas those with wAna^ITG are indicated with red colored branches. Both trees are midpoint rooted.

Fig. 7.—

Analysis of sites called as heterozygotes in the wAna^ITG genomes. (A) For each wAna^ITG genome, a nonoverlapping sliding window analysis of comparing average read depth and the total number of single nucleotide variants. (B) Each data point represents a site called as a heterozygote within a wAna^ITG genome, with x axis representing the total read coverage and y axis representing the proportion of alternative alleles consisting the heterozygote variant. Spearman’s rho values are shown in top right corner of each figure and both correlations are highly significant (P < 0.0001).