Endosymbiont Tremblaya phenacola influences the reproduction of cotton mealybugs by regulating the mechanistic target of rapamycin pathway

Abstract

The intricate evolutionary dynamics of endosymbiotic relationships result in unique characteristics among the genomes of symbionts, which profoundly influence host insect phenotypes. Here, we investigated an endosymbiotic system in Phenacoccus solenopsis, a notorious pest of the subfamily Phenacoccinae. The endosymbiont, “Candidatus Tremblaya phenacola” (T. phenacola PSOL), persisted throughout the complete life cycle of female hosts and was more active during oviposition, whereas there was a significant decline in abundance after pupation in males. Genome sequencing yielded an endosymbiont genome of 221.1 kb in size, comprising seven contigs and originating from a chimeric arrangement between betaproteobacteria and gammaproteobacteria. A comprehensive analysis of amino acid metabolic pathways demonstrated complementarity between the host and endosymbiont metabolism. Elimination of T. phenacola PSOL through antibiotic treatment significantly decreased P. solenopsis fecundity. Weighted gene coexpression network analysis demonstrated a correlation between genes associated with essential amino acid synthesis and those associated with host meiosis and oocyte maturation. Moreover, altering endosymbiont abundance activated the host mechanistic target of rapamycin pathway, suggesting that changes in the amino acid abundance affected the host reproductive capabilities via this signal pathway. Taken together, these findings demonstrate a mechanism by which the endosymbiont T. phenacola PSOL contributed to high fecundity in P. solenopsis and provide new insights into nutritional compensation and coevolution of the endosymbiotic system.

Introduction

Nutritional symbioses between animals and microorganisms contribute positively to the adaptability of both [1]. Some symbionts play vital roles in enabling the colonization of nutrient-deficient environments through the supplementation of the host’s metabolic requirements, which makes insects thrive on imbalanced carbohydrate-based diets [2, 3]. Long-term symbiotic relationships shape the bacterial genome and significantly affect the host phenotypes related to nutrition, metabolism, reproduction, immunity, and development [4–7]. For example, a genomic research on the pea aphid, Acyrthosiphon pisum, revealed an extensive metabolic exchange with its obligate nutritional endosymbiont, Buchnera aphidicola [8]. Such nutritional symbionts in aphids broaden the dietary range, influencing reproduction processes during aphids’ sexual phase [9]. Another interesting study of carpenter ant and its primary endosymbiont, Blochmannia floridanus, showed that the symbiont provides essential nitrogen and sulfur compounds to the ants, benefiting in return from the ants’ metabolic processes [10]. The intricate symbiotic partnerships between hosts and their endosymbionts lead to highly specialized interactions that can significantly affect the biology, ecology, and evolution of both partners.

Mealybugs (Hemiptera: Coccoidea: Pseudococcidae) are recognized as invasive pests in numerous regions across the world [11]. These diminutive sap-sucking insects feed on a broad range of vegetables, horticultural plants, and field crops, causing huge economic losses worldwide [12, 13]. Previous studies have identified two subfamilies of mealybugs, Pseudococcinae and Phenacoccinae [14], which have distinct lineages of bacterial endosymbionts [15] and patterns of nutritional complementation. A representative pest of the Phenacoccinae subfamily, cotton mealybug (Phenacoccus solenopsis Tinsley), infests over 200 plant species and has been observed in at least 35 geographic regions worldwide [16, 17]. The broad adaptability of this pest has been attributed in part to its beneficial mutualistic endosymbiotic system.

Mealybug primary endosymbionts (P-endosymbionts) are maternally inherited and localize to specialized host cells known as bacteriocytes. These cells accumulate to form larger structures, bacteriomes [18, 19]. Subfamily Pseudococcinae P-endosymbionts belong to the betaproteobacteria class “Candidatus Tremblaya princeps,” (T. phenacola PSOL), which in turn host co-obligate gammaproteobacterial endosymbionts [20]. This three-way symbiotic system has been extensively studied in terms of horizontal gene transfer, metabolic patchwork, and coevolution [18, 21–23]. By contrast, members of the subfamily Phenacoccinae contain only the betaproteobacterium “Candidatus Tremblaya phenacola” [24]. The sequencing and analysis of the genome of T. phenacola PPER (isolated from Phenacoccus peruvianus) revealed genome fusion between a betaproteobacterium and a gammaproteobacterium [25]. Furthermore, phylogenetic studies of T. phenacola and other insects within the subfamily Phenacoccinae have suggested independent symbiont coevolution in each insect host species [26, 27].

In addition to essential amino acids synthesis, mutualistic symbiosis in insects are capable of producing some important molecules that affect specific biological process. A recent study on aphid-Buchnera symbiosis highlights the differential expression of various eukaryotic cell signaling pathways in bacteriocytes of low-titer versus high-titer hosts. In low-titer genotypes, cell-growth pathways are upregulated, while in high-titer genotypes, pathways related to membrane trafficking, lysosomal processes, the mechanistic target of rapamycin (mTOR), and cytokine pathways are more active [28]. The mTOR signaling pathway is implicated in regulating both the reproductive capacity and lifespan of social insects, particularly termites [29]. The influence of essential amino acids produced via the mutualistic system on the reproductive processes of the cotton mealybug continues to be an area of uncertainty.

To gain a better understanding of coevolution and mutualistic symbiosis between the cotton mealybug and Tremblaya, we here examined the distribution of T. phenacola PSOL in male and female cotton mealybugs. We then assembled the T. phenacola PSOL genome to analyze the phylogenetic relationships among Tremblaya symbionts. Finally, we analyzed metabolic complementation within the symbiotic system and assessed the impacts of this system on cotton mealybug fecundity. These experiments directly demonstrate endosymbiont impacts on host insect reproductive phenotypes, providing deepened insights into the coevolution of an agriculturally important insect pest and its bacterial endosymbiont.

Materials and methods

Insect population

A laboratory colony of cotton mealybugs was used for these experiments. The progenitors were originally collected from the ornamental plant Rose of Sharon (Hibiscus syriacus L.) in Jinhua, China. Insects were reared on tomato plants (Solanum lycopersicum) in a climate-controlled chamber at 27 ± 1°C with 75% relative humidity.

Fluorescence in situ hybridization

To detect the distribution of the endosymbiont in the cotton mealybug, we used Fluorescence in Situ Hybridization (FISH), employed in a previous study [30]. In details, the second and third instar larvae, female adult, pupa, and male adult of cotton mealybugs were disinfected using 70% alcohol for 1 min, followed by 0.01% sodium hypochlorite for 1 min. They were then rinsed with sterile ddH₂O three times. The samples were fixed in Carnoy’s solution (chloroform: ethanol: acetic acid = 6:3:1) at 55°C for 30 min. After removing the fixatives, the samples were washed three times in sterile phosphate-buffered saline (PBS) buffer, then placed in 6% H₂O₂-ethanol solution to decolorize for 1 h. After the decolorizing agent was removed, hybridization solution (HYBA) Buffer (comprising 5 ml methanamide, 2.5 ml 20× saline sodium citrate, 100 μl heparin, 100 μl 10% Tween-20, and 105 μl 100 μg/ml salmon sperm DNA) was added to each sample for prehybridizing at 65°C for 2 h. The specific 16S rDNA probe (cy5–5’ATCTACGCATTTCACCGCTACTCCTGGAATTCTATCCCCCTCTTCCATACTCGAG-3′) was mixed with HYBA Buffer at a ratio of 1:99 and then incubated for 5 min at 80°C; 1 ml HYBA Buffer/probe mixture was then added to each sample and was incubated overnight at 65°C in the dark. Samples were washed with PBS + Tween, then incubated with 4′,6-diamidino-2-phenylindole for 30 min. Each sample was then observed and photographed on an LSM 800 laser confocal microscope (ZEISS, Oberkochen, Germany).

Bacteriome sequencing and assembly

About 150 female mealybugs were starved for 48 h, then the bacteriomes of these mealybugs were collected for T. phenacola PSOL DNA extraction. In detail, individuals were disinfected in 70% alcohol for 1 min followed by 0.01% sodium hypochlorite for 1 min; rinsed with sterile ddH₂O three times; then rinsed with sterile 0.2 M PBS (pH 6.8) three times. Bacteriomes were collected through insect dissection in a drop of sterile 0.2 M PBS in a sterile petri dish (1.5 × 9 cm) under an SMZ 645 stereo microscope (Nikon, Tokyo, Japan). The collected bacteriomes were gently washed twice with sterile 0.2 M PBS prior to DNA extraction with a DNA Isolation Mini Kit (Vazyme Biotech, Nanjing, China).

Extracted T. phenacola PSOL DNA was used to generate paired-end Illumina libraries with a 350-bp insert size. Libraries were sequenced on the HiSeq X-Ten platform (Illumina, San Diego, CA, USA). Duplicate reads were removed with FastUniq v1.1 [31] and the remaining reads were merged with FLASH v1.2.11 [32]. Merged reads were first mapped to the cotton mealybug genome assembly [33] with BWA v0.7.17 [34], followed by the filtration of host DNA contamination using Samtools v1.17 [35]. The whole-body PacBio reads of the cotton mealybug (PRJNA380754) were used as queries in BLAST searches against the cotton mealybug genome assembly using Megablast v2.14.0+ [36]. Unaligned sequences were collected with Seqkit v2.4.0 [37] and combined with the unmapped Illumina short reads; these sequences were regarded as candidate symbiont reads. A hybrid assembly was conducted in Unicycler using the “normal” mode [38]. Contigs of <500 bp in size were excluded from further analyses. The final assembly coverage was calculated with Samtools v1.17 [35].

T. phenacola PSOL genome annotation and taxonomic assignment

The T. phenacola PSOL genome was annotated to identify protein-coding, transfer RNA, and ribosomal RNA (rRNA) genes. These were predicted using the RAST server v1.1.0 [39] with default parameters and an E-value of 1e-5. The information on repetitive sequences was predicted using the RepeatMasker web server (https://www.repeatmasker.org/cgi-bin/WEBRepeatMasker). Conserved syntenic blocks between T. phenacola PSOL and T. phenacola PAVE or T. phenacola PPER were identified with TBLASTX and visualized in Circos [40].

Taxonomic assignments were conducted using two methods. First, all predicted protein-coding genes from T. phenacola PSOL were used as queries in BLAST searches against the NCBI nonredundant protein sequence (nr) database with the parameters “-evalue 1e-5 -outfmt 5.” The output xml file was then analyzed in Megan v6.12.2 [41]. Second, phylogenetic gene trees were constructed for each T. phenacola PSOL gene to determine the putative evolutionary position. A betaproteobacteria-gammaproteobacteria reference database was constructed from the protein sequences of all sequenced Tremblaya symbiovars, their obligatory gammaproteobacteria, and 24 489 betaproteobacterial and 122 778 gammaproteobacterial protein sequences downloaded from UniProtKB. Protein sequences of the alphaproteobacterium Rickettsia were downloaded from UniProtKB and used as an outgroup. For each T. phenacola PSOL protein, a phylogenetic tree was constructed from the 30 best BLASTP hits against the beta–gamma database, the two best BLASTP hits against the alpha database, and the target protein itself. Protein sequences were aligned by MAFFT v7.310 [42]. The program trimAL v1.2rev59 [43] with the “-automated1” parameter was used to exclude ambiguously aligned positions. The best substitution model was selected using Bayesian Information Criterion in IQ-TREE v1.5.5 [44] with the “test” mode. IQ-TREE was then used to construct a tree for each protein using the maximum likelihood (ML) method with 1000 bootstrap replicates. Single-gene trees were visualized by ggtree v2.4.0 [45]. Based on these steps, we finally obtained 92 beta-origin and 104 gamma-origin genes in T. phenacola PSOL.

Phylogenetic analyses of different mealybugs and their symbionts

Five mealybug genome assemblies were downloaded from the European Nucleotide Archive: Maconellicoccus hirsutus (PRJEB12066), Ferrisia virgata (PRJEB12067), Pseudococcus longispinus (PRJEB12068), Paracoccus marginatus (PRJEB12069), and Trionymus perrisii (PRJEB12071). These genomes were annotated using combined evidence from Augustus [46], BLAST, and HMMER [47]. The aphid A. pisum (GenBank accession GCA_000142985.2) was used as the outgroup. Proteins encoded by each gene family were aligned and trimmed as described above for single-gene tree construction. All of the trimmed single-copy protein sequences were then concatenated to generate one super-gene for each species. An ML tree was constructed based on the super-gene with the JTT + F + I + G4 method and 1000 bootstrap replicates.

For each class of proteobacteria, several mealybug symbionts and representative proteobacteria were selected for analysis. The T. phenacola PSOL and T. phenacola PPER proteins of betaproteobacterial and gammaproteobacterial origin (referred to as beta-origin and gamma-origin genes, respectively) were used in betaproteobacteria and gammaproteobacteria phylogenetic tree constructions, respectively. The alphaproteobacterium Rickettsia prowazekii (RefSeq accession GCF_000195735.1) was used as the outgroup. Phylogenetic trees were constructed as described above, with the betaproteobacterial tree based on 28 single-copy proteins (Table S1) using the LG + I + G4 model and the gammaproteobacterial tree based on 23 single-copy proteins (Table S2) using the Q.yeast+I + G4 model.

Substitution rate analysis

Single-copy beta-origin orthologs between T. phenacola PSOL and T. princeps PCIT and single-copy gamma-origin orthologs between T. phenacola PSOL and Moranella endobia PCIT were aligned with MAFFT v7.310. KaKs Calculator v2.0 [48] was employed to calculate the nonsynonymous substitution (dN) to synonymous substitution (dS) ratio between orthologous pairs with the parameters “-c 11 -m MA.” We also calculated the dN/dS ratio between orthologous pairs of T. phenacola PSOL and T. phenacola PPER using the same method. All statistical analyses and visualizations were performed in R v3.5.2.

Amino acid biosynthesis pathway construction

Genes involved in amino acid biosynthesis pathways were identified using an annotated cotton mealybug genome (http://v2.insect-genome.com/Organism/624), the T. phenacola PSOL genome annotated here, and the BlastKOALA tool on the Kyoto Encyclopedia of Genes and Genomes (KEGG) website (https://www.genome.jp/kegg/) [49]. Similarly, we analyzed the amino acid biosynthesis pathways of the other five mealybugs including F. virgata, M. hirsutus, P. longispinus, P. marginatus, T. perrisii. These genome data were obtained from a previous study in 2016 [18]. Complete metabolic pathways in these mealybug symbiotic systems were constructed using KEGG Mapper (https://www.genome.jp/kegg/mapper/) [50].

Mealybug phenotype analysis

Mealybugs were treated with antibiotics to reduce endosymbiont abundance. Tetracycline and gentamicin were obtained from Solarbio Technology Co., Ltd (Beijing, China). Antibiotic solutions (0.25 mg/ml and 0.75 mg/ml for both antibiotics) were prepared in sterile ddH₂O. Three-week-old tomato plants were positioned such that the branch roots were inside plastic cups (5 × 7 × 8 cm) filled with 20 ml of antibiotic solution or ddH₂O as a control. Tomato branches were vertically secured within the plastic cups using foam trays (0.4 × 5 cm thickness × diameter). Antibiotic solution was refreshed every 48 h until insects reached adulthood. The durations of preoviposition, oviposition, and postoviposition were recorded, as were the number of offspring and the offspring sex ratios. Each group contained 14–32 female adults, with 1 female representing one biological replicate.

To evaluate the physiological impact of antibiotics on mealybugs, we evaluated the effects of different antibiotic concentrations on mealybug survival rates, and each group contained 45 mealybug individuals. We also measured the impact of antibiotics on the ATPase content in cotton mealybugs using insect-ATPase ELISA kit (Kenuodibio, China). To process the insect samples for biochemical analysis, we homogenized 50 mg mealybugs for each group in 500 μl of ice-cold PBS. The homogenates were then centrifuged at 1000 × g for 20 min, and the supernatant was collected for further analysis. For the construction of a standard curve, we prepared dilutions of a standard substance to concentrations of 0.5, 1, 2, 4, 8, and 16 U. In a 96-well plate, we added 50 μl of the supernatant, followed by 100 μl of horseradish peroxidase-conjugated detection antibody, and it was incubated at 37°C for 60 min. After incubation, the liquid was discarded, and the plate was washed five times with washing buffer. Subsequently, 50 μl each of substrates A and B were added, and the plate was incubated at 37°C for 15 min. The reaction was stopped by adding 50 μl of stop solution, and the absorbance was measured at a wavelength of 450 nm.

Levels of free amino acids were tested in the control group and in the 0.25 mg/ml tetracycline group of mealybugs. Briefly, a collection of mealybugs from each treatment group weighing 50 mg in total were placed in separate centrifuge tubes and were homogenized in 500 μl of 2% sulfosalicylic acid. The homogenates were passed through a 0.22-μm filter membrane, then the free amino acid concentrations of each solution were measured with an L8900 amino acid analyzer (HITACHI, Tokyo, Japan). The data were visualized in GraphPad Prism 9.0.

RNA sequencing and gene expression analysis

First, second, third, pupa, male, and female of the cotton mealybugs were collected for RNA extraction. For each sample, the collection included 50 individuals for the group of first, second, and pupa; 100 individuals for male adults; 15 individuals for third and female adults. Total RNA was extracted from each replicate sample using TRIzol reagent (Thermo Fisher Scientific, Waltham, MA, USA). Sequencing libraries were constructed using an NEBNext Ultra RNA Library Prep Kit (New England Biolabs, Ipswich, MA, USA). All experiments were repeated in triplicate. In detail, mRNA enrichment was initially carried out using Oligo(dT) magnetic beads, followed by the fragmentation of mRNA. Single-stranded cDNA was then synthesized using random hexamer primers. Subsequently, double-stranded cDNA was produced by adding buffer, dNTPs, DNA polymerase I, and RNase H, then purified using AMPure XP beads. Finally, the double-stranded cDNA was processed for end repair, poly-A tailing, and adapter ligation, preparing it for sequencing. The resulting libraries were sequenced on the HiSeq4000 platform (Illumina) to generate paired-end 150-bp reads. Raw reads were filtered using fastp v0.23.4 [51] and then mapped to reference genome of cotton mealybug using HISAT v2.2.1 [52]. StringTie v2.2.1 [53] was used to assemble transcripts and estimate their abundance.

Gene expression levels were calculated in fragments per kilobase of transcript per million mapped reads. Based on these values, a gene coexpression network was constructed with weighted gene coexpression network analysis (WGCNA) [54, 55]. The TOMSimilarity module was used to calculate the coexpression similarity coefficients between genes with the following parameters: power = 6, TOMType = unsigned, deepSplit = 2, mergeCutHeight = 0.25, and numericLabels = TRUE. KEGG enrichment analysis of each cluster was performed with clusterProfiler v4.6.2 [56].

Quantitative reverse transcription

Total RNA of the cotton mealybugs was extracted using TRIzol reagent, and 1 μg of RNA was reverse-transcribed using HiScript III 1^st Strand cDNA Synthesis Kit (Vazyme Biotech). We designed gene-specific primers based on the sequence of 16S rRNA of T. phenacola PSOL for endosymbiont quantification, with the host’s Actin selected as the reference gene (Table S3). Quantitative reverse transcription was performed on a QuantStudio thermocycler (Thermo Fisher Scientific) in 20-μl reaction volumes containing 10 μl of HiScript III RT SuperMix (Vazyme), 1 μl of 100 ng/μl cDNA, and 1 μl of 50 ng/μl primers. The thermocycling program was as follows: initial denaturation, 95°C for 10 min; 40 cycles of denaturation (95°C for 15 s) and annealing (60°C for 30 s). The data were analyzed using the 2^−ΔΔC method [57]. These experiments were conducted with at least three independent biological replicates.

Statistical analyses

Differences between treatment groups in the preoviposition duration, postoviposition duration, and number of offspring were assessed with a Poisson regression or Poisson-like regression model. A logistic regression model was used to examine treatment-based differences in the sex ratio. Differences between treatment groups in other parameters were assessed with a one-way analysis of variance (ANOVA).

Results

Dynamic distribution of T. phenacola PSOL in cotton mealybugs

We first conducted FISH experiments to detect the endosymbiont distribution within the host using a T. phenacola PSOL-specific 16S rDNA probe. In female mealybugs, the endosymbiont was present in the bacteriome of the abdominal midsection throughout the entire lifecycle (Fig. 1A). By contrast, the bacteriome gradually disintegrated during pupation in males, resulting in a faint, dispersed signal in adult males (Fig. 1A). Subsequent qPCR experiments revealed significant differences in the abundance of T. phenacola PSOL between female and male mealybugs (Fig. 1B, P < .001). These results were comparable, with consistent levels throughout the lifecycle in the females and sharply decreased abundance in adult males. During the oviposition process in female mealybugs, we observed a significant increase in the abundance of T. phenacola PSOL (Fig. 1B, P < .01).

T. phenacola PSOL genome analysis reveals a fusion of beta and gamma proteobacterial genes

In a previous study, our lab isolated four scaffolds originating from T. phenacola PSOL using a PacBio library of DNA extracted from cotton mealybug [33]. To improve the endosymbiont genome assembly, we extracted bacteriome DNA to sequence T. phenacola PSOL directly. A total of 6.44 Gb clean data from T. phenacola PSOL bacteriomes were obtained, yielding a high-quality assembly comprising seven contigs with a genome size of 221.1 kb with 100 × coverage (Fig. 2A). The new genome assembly included 229 predicted protein-coding genes. Single-gene phylogenetic trees and analysis with Megan indicated that these genes had three distinct origins: betaproteobacteria (40.2%), gammaproteobacteria (45.4%), and unknown (14.4%) (Fig. 2A). Both beta-origin and gamma-origin genes of T. phenacola PSOL exhibited a high synteny with T. phenacola PPER (E-value <1e-5, identity >75%, length > 100 bp) (Fig. 2B), and most of the beta-origin genes also displayed high synteny with the T. phenacola PAVE genome (Fig. S1). Compared with other known Tremblaya genomes, T. phenacola PSOL has the largest genome size but the lowest GC content (34.8%) (Table S4). Compared with T. phenacola PPER, T. phenacola PSOL possesses a higher number of genes originating from gamma sources (Fig. S2). This difference might contribute to the lower GC content observed in T. phenacola PSOL. KEGG analysis indicated that beta-origin and gamma-origin genes generally had roles in distinct biological processes. For example, most genes predicted to be involved in environmental adaptation, transcription, and nucleotide metabolism are of beta-origin, whereas those involved in immune system and cell motility were gamma-origin genes (Fig. S3).

To elucidate the evolutionary relationships among mealybug endosymbionts, a phylogenetic tree was constructed. The T. phenacola PSOL was most closely related to T. phenacola PPER, whereas T. phenacola PAVE was more closely related to T. princeps symbiovars (Fig. 2C). The gammaproteobacteria present in T. phenacola PSOL and those isolated from T. phenacola PPER clustered together within the same monophyletic group, which was a sister taxon of the Sodalis lineage (Fig. 2D).

Genes originating from both beta and gammaproteobacteria are evolutionarily conserved

To evaluate the selective strength of genes in T. phenacola PSOL, we calculated the rates of synonymous (dS) and nonsynonymous (dN) substitutions between orthologous genes from the same lineage (beta or gamma) in T. phenacola PSOL and T. phenacola PPER. All T. phenacola PSOL genes had a dN/dS rate < 1 (mean = 0.20), but gamma-origin genes had significantly lower dN/dS rates than beta-origin genes (Wilcoxon rank sum test, P = 4 × 10⁻⁴) (Fig. 3A and B). For amino acid metabolism and translation processes, gamma-origin genes are under less selective pressure (Fig. 3C and D). Similar results were also obtained from substitution analysis in T. phenacola PSOL and the symbiotic system of T. princeps PCIT (beta) and M. endobia PCIT (gamma) (Fig. S4). Besides, gene expression analysis revealed that gamma-origin genes, especially those related to amino acid synthesis, are highly expressed during the female adult stage (Fig. S5B and Table S5).

Metabolic complementation for amino acid synthesis

Given high conservation of T. phenacola PSOL in the amino acid metabolism pathway and its nutritional complementation with the host, we constructed amino acid biosynthesis pathways using protein-coding genes from both cotton mealybug and the endosymbiont. The results indicated that the two species formed highly complementary “patchwork” metabolic pathways (Fig. 4), with varying components of each pathway present in only one partner. Endosymbiont genes participating in valine, leucine, and isoleucine biosynthesis were of betaproteobacterial origin (Fig. 4A), whereas genes involved in cysteine and histidine biosynthesis were of gammaproteobacterial origin (Fig. 4B and C). Other genes (such as those involved in tryptophan, tyrosine, and phenylalanine biosynthesis) were of both betaproteobacterial and gammaproteobacterial origin (Fig. 4D). Besides, upon comparing the amino acid synthesis pathways of various symbiotic systems, we observed that these pathways are largely conserved, yet exhibit species-specific characteristics across different mutualistic systems (Fig. S6). Compared to other systems, the genes for amino acid synthesis enzymes show less redundancy in T. phenacola PSOL (Fig. S6B).

T. phenacola PSOL promotes cotton mealybug fecundity

To gain deep insights into the effects of amino acids produced by this patchwork on insect physiology, the cotton mealybugs were exposed to antibiotics (Fig. 5A). Due to the inevitable physiological impact of antibiotics on insects, we assessed the effects of antibiotic toxicity on cotton mealybug survival rates and mitochondrial ATPase activity. Our results indicated that both gentamicin and tetracycline had no significant impact on these aspects (Fig. S7A and D). Treatment with tetracycline significantly inhibited the abundance of T. phenacola PSOL (Fig. S7B), leading to a decrease in partial free amino acid contents (Fig. S7C) and even impacting wax synthesis (Fig. S7E).

The preoviposition period was significantly extended in groups treated with 0.25 mg/ml or 0.75 mg/ml tetracycline (P < .05; Poisson regression) (Fig. 5B). Because oocysts persist in female adults upon death, the oviposition period and the postoviposition period were combined into a single period (postoviposition) for further analysis; this period was also extended in females treated with tetracycline (Fig. 5C). Female insects treated with either the low or high concentration of tetracycline showed a significant reduction in the number of offspring produced (Fig. 5D, P < .001; Poisson-like regression). Treatment with 0.25 mg/ml tetracycline led to a decline in the proportion of male offspring (P < .05; logistic regression), but no other antibiotic treatment significantly influenced the offspring sex ratio (Fig. 5E).

T. phenacola PSOL activated the host mechanistic target of rapamycin pathway

To further explore the mechanism by which the endosymbiont regulated host reproduction, we conducted a comprehensive gene coexpression analysis. Hierarchical clustering and pruning of adjacent dissimilarity yielded a gene-dependent dendrogram (Fig. S8A). The blue module exhibited a strong correlation with female developmental stage (Fig. S8B). Furthermore, genes in this module exhibited significant correlations between expression in T. phenacola PSOL and in several cotton mealybug processes, including the mTOR signaling pathway, oocyte meiosis, and oocyte maturation (Fig. 6A). KEGG enrichment analysis of endosymbiont genes in the blue module revealed a notable enrichment of genes associated with the amino acid synthesis pathway (Fig. 6B). To assess whether the endosymbiont mediates mTOR signaling through essential amino acid synthesis, thereby influencing host oocyte meiosis and maturation, we analyzed the expression levels of three genes (S6K, SKP2, and Tuberin) in this pathway enriched in the blue module. The results showed that antibiotic treatment significantly activated the mTOR pathway of the cotton mealybug (Fig. 6C–E).

Discussion

Many insects harbor symbionts that supply them with essential nutrients, providing benefits in development, immune metabolism, and various other biological processes [58–62]. In return, insects provide nutrient-rich intracellular niches for symbionts to colonize. Due to their extensive adaptation to the intracellular environment, nearly all endosymbionts are currently classified as unculturable. Thus, bacteriome genomic sequencing has emerged as a powerful tool to understand endosymbionts. Here, we conducted a comprehensive analysis revealing the vital role of the endosymbiont T. phenacola PSOL in the reproductive processes of its host, P. solenopsis.

We examined the distribution of T. phenacola PSOL across various developmental stages of the host. Our findings revealed a marked disparity in the presence of T. phenacola PSOL between female and male adults (Fig. 1). Similar results have also been reported in several Pseudococcinae [63]. Symbionts exhibit increased activity during oviposition (Fig. 1). This appears to be a mutual adaptation mechanism for both the symbiont and its host. In this process, symbionts enhance their chances of transmitting genetic material to offspring, whereas the host benefits from the nutrients, resulting in the production of a greater number of offspring (Figs 4 and 5). This dynamic relationship was also observed in the symbiotic system of the cereal weevil Sitophilus oryzae, and its midgut endosymbionts, Sodalis pierantonius. Following the adult metamorphosis of its host, the population of S. pierantonius significantly increases, but then rapidly declines and it is completely eliminated once the host finishes its cuticle synthesis [64].

In pursuit of comprehending this intriguing adaptation mechanism in the cotton mealybugs, we constructed a genome assembly of the endosymbiont T. phenacola PSOL (Fig. 2A). Due to limitations in sequencing methods and the repetitive sequences in the flanks (Fig. 2A), this assembly did not result in a single, complete contig for the endosymbiont. The synteny analysis with T. phenacola PPER indicated that this assembly is not the result of contamination or misassembly (Fig. 2B). Phylogenetic analysis indicated that T. phenacola PSOL was a chimeric endosymbiont resulting from genome fusion between betaproteobacteria and gammaproteobacteria. This is distinct from all other sequenced Tremblaya symbiovars, with the exception of T. phenacola PPER (Fig. 2C and D). Further analysis of the selective pressure on genes originating from betaproteobacteria and gammaproteobacteria indicates they have undergone purifying selection, with genes from gamma-proteobacteria showing a higher degree of conservation, particularly in amino acid metabolism and translation processes (Fig. 3). This is consistent with prior studies showing that a stable host environment, high nutritional supply, and redundancy with host genes lead to relaxed purifying selection within the symbiont genome [65].

Similar to other mealybug endosymbionts, the primary advantage of T. phenacola PSOL symbiosis to the host appeared to be provision of essential amino acids [23, 66]. Our findings here revealed a pronounced complementarity between host and endosymbiont amino acid synthesis genes, with limited functional redundancy identified between the symbiont and host genomes (Fig. 4 and Fig. S6). For example, cotton mealybug lost the cysK gene (a member of the cysteine biosynthesis pathway), but this loss was offset by the presence of a gamma-origin gene in T. phenacola PSOL (Fig. 4B). This stands in contrast to the M. hirsutus symbiotic system, which has lost the cysK gene entirely during the process of coevolution [18]. Additionally, beta-origin and gamma-origin symbiont genes in this pathway exhibited total complementarity (i.e. a complete lack of functional redundancy) (Fig. 4 and Fig. S6B). Previous studies have revealed similarly complementary patterns of gene loss and retention among the primary symbiont Tremblaya, the secondary symbiont Moranella, and the citrus mealybug host (Planococcus citri) [21, 67]. Overall, these data allowed us to observe the ways in which symbionts evolve in response to their own transmission, simultaneously assisting their host in overcoming adverse factors.

Mutualistic symbionts usually exert a positive effect on host’s reproduction. We here conducted the investigation of the reproductive impacts of a mutualistic endosymbiont, T. phenacola PSOL, on the cotton mealybug. Due to the inevitable physiological impact of antibiotics on insects, along with their suppressive effect on nonintracellular symbionts, we opted for an antibiotic with minimal intracellular action, gentamicin, and an antibiotic capable of killing intracellular bacteria, tetracycline [68, 69]. Unexpectedly, high concentrations of gentamicin still exhibited an inhibitory effect on T. phenacola PSOL (Fig. S7B). At the meanwhile, antibiotic interference led to reproductive disorders in the host, manifested as an extended preoviposition period and a notable decrease in offspring numbers (Fig. 5B–D). Groups treated with 0.25 mg/ml or 0.75 mg/ml gentamicin or with 0.75 mg/ml tetracycline exhibited nearly identical proportions of male offspring as untreated insects. This suggested that the presence of T. phenacola PSOL may not have been a significant factor in determining the sex ratio of cotton mealybug. However, several prior studies have demonstrated that nutritional bacteriocyte symbionts can manipulate the host sex ratio [61, 62, 70]. Our results underscore the considerable variation in reproductive impacts induced by different endosymbionts.

To further investigate the mechanisms through which amino acid compensation by an endosymbiont can positively impact host reproduction, we conducted a comprehensive transcriptomic analysis of the endosymbiont and the cotton mealybug at several developmental stages. WGCNA revealed one coexpression module containing T. phenacola PSOL genes that were coexpressed with host genes associated with amino acid biosynthesis, the mTOR pathway, and oocyte meiosis and maturation (Fig. 6A and B). This was consistent with prior studies showing that amino acids are crucial for animal reproductive processes [71–73]. Furthermore, branched-chain amino acids (namely leucine, L-isoleucine, and valine) are essential amino acids for mammals, and these metabolites activate the mTOR signal pathway through a variety of mechanisms [74]. In the present study, T. phenacola PSOL removal from cotton mealybug through antibiotic treatment upregulated several key genes in the mTOR pathway. One such gene was S6K, which is typically activated by growth factors, cytokines, and nutrients [75]. The stimulation of S6K in controlling cell growth has also been observed in Drosophila after amino acid deprivation [76]. Two other mTOR genes here found to be upregulated after antibiotic treatment, SKP2 and Tuberin, have multifaceted roles in cell cycle regulation, with primary roles in the regulation of protein synthesis [77, 78]. The expression patterns of these genes reflect a strong correlation between T. phenacola PSOL and the mTOR pathway of host, suggesting that the endosymbiont has essential roles in influencing the host reproductive phenotype. We reason that the cotton mealybug endosymbiont influenced host reproduction by regulating the expression of genes in the mTOR pathway.

In summary, we here generated a high-quality genomic resource that not only facilitated a comprehensive investigation of the cotton mealybug genomic landscape but will also serve as a valuable reference to dissect coevolution among members of an important symbiotic relationship. Investigation of the intricate dynamics by which the endosymbiont contributes to host reproduction sheds light on the interplay between their complementary modes of nutrient acquisition (Fig. 7). Our findings underscore broader trends in mutualistic symbiotic systems across diverse ecological niches. This study thus contributes to a deep comprehension of the remarkable diversity and intricacy inherent in mealybug symbiotic relationships, opening new avenues for the ecological and evolutionary understanding of endosymbiotic systems.

Acknowledgements

This research was supported by Key Project of Zhejiang Provincial Natural Science Foundation (LZ23C140002), National Key Research and Development Program (2019YFD1002105, 2022YFD1401600), National Natural Science Foundation of China (32102189), and Basic Public Welfare Research Program of Zhejiang Province (LY22C140005).

Author contributions

J.B., Z.Z., and H.D. performed the laboratory experiments and analyzed the data. M.L, H.T., Y.M., Y.X., K.H., and M.J. coordinated the study and validated the experimental results. S.W. and F.L. contributed to the conceptual design of the study, preparation of the manuscript, and securing funding support. All authors have read and agreed to the published version of the manuscript.

Conflicts of interest

The authors declare no financial conflicts of interest.

Data availability

Raw sequencing reads of mealybug bacteriomes have been deposited into the NCBI (Bio-Project: PRJNA550279). The genome assembly, genome annotation, and protein FASTA files have been deposited in http://v2.insect-genome.com/Organism/624.

References

Author notes