Detoxification Gene Families at the Genome-Wide Level of Rhus Gall Aphid Schlechtendalia chinensis

The Rhus gall aphid Schlechtendalia chinensis uses the species Rhus chinensis as its primary host plant, on which galls are produced. The galls have medicinal properties and can be used in various situations due to their high tannin content. Detoxification enzymes play significant roles in the insect lifecycle. In this study, we focused on five detoxification gene families, i.e., glutathione-S-transferase (GST), ABC transporter (ABC), Carboxylesterase (CCE), cyto-chrome P450 (CYP), and UDP-glycosyltransferase (UDP), and manually annotated 144 detoxification genes of S. chinensis using genome-wide techniques. The detoxification genes appeared mostly on chromosome 1, where a total of two pair genes were identified to show tandem duplications. There were 38 gene pairs between genomes of S. chinensis and Acyrthosiphon pisum in the detoxification gene families by collinear comparison. Ka/Ks ratios showed that detoxification genes of S. chinensis were mainly affected by purification selection during evolution. The gene expression numbers of P450s and ABCs by transcriptome sequencing data were greater, while gene expression of CCEs was the highest, suggesting they might be important in the detoxification process. Our study has firstly identified the genes of the different detoxification gene families in the S. chinensis genome, and then analyzed their general features and expression, demonstrating the importance of the detoxification genes in the aphid and providing new information for further research.


Introduction
The Rhus gall (or sumac-gall) aphids switch host plants between Rhus species and mosses to finish their life cycles, and form galls on their primary host plants Rhus (Anacardiaceae) species [1][2][3]. The galls are often known as the Chinese galls and they are rich in tannins and economically important in Asia because galls have medicinal properties and represent sources of industrial tannin [4,5]. This aphid group belonged to the subtribe Melaphidina of tribe Fordini (Aphididae: Eriosomatinae) [6][7][8], including six genera and 12 species [3], among which S. chinensis is the most common and wide-spread species, with R. chinensis as its unique primary host plant and Mniaceae species as its secondary hosts, as well as a life cycle including both sexual and asexual reproduction stages [9,10].
Mature and fresh Rhus galls were collected from the garden in Shanxi University, Shanxi, China. The gall was cut open and live aphids were used for third-Generation Sequencing. In addition, we selected the aphid samples from three mature galls, numbered A4601, A4603, and A4621 for transcription sequencing to characterize the gene expression pattern of detoxification genes in S. chinensis. The specimens were stored at the School of Life Science, Shanxi University, China.

Identification of Detoxification Genes from S. chinensis
We sequenced the whole genome of S. chinensis by third-generation sequencing on the Pacbio platform. The protein-coding genes in the genome were annotated by integrating three approaches, namely de novo prediction, homology search, and transcript-based assembly. The protein sequences of the detoxification genes were obtained by searching in the annotation table of S. chinensis using key words of detoxification genes, and then confirmed in the NCBI Conserved Domain Database (NCBI-CDD) (e-value = 1 × 10 −3 ), excluding those lacking conserved domains from the analysis [23].

Phylogenetic Tree, Motif Pattern, Domain, Gene Structure of Detoxification Genes
The protein coding sequences of detoxification genes of the three aphid species, A. pisum, Cinara cedri and Myzus persicae, were downloaded from the Insect BASE website (http://v2.insect-genome.com/, accessed on 12 May 2022). We performed the multiplesequence alignment of detoxification genes using ClustalW software (created by Kumar et al.; Philadelphia, PA, USA) [24]. The alignment results were exported to fasta format, and then opened using tbtools to trim using Trimmer [25]. The protein sequences with large differences were filtered out. We constructed a phylogenetic tree using the Neighbor-Joining (NJ) method with the parameters of Poisson model, complete deletion and 1000 bootstrap replicates, and visualized and improved the tree using the program Evolview (http://www.evolgenius.info/evolview/, accessed on 18 March 2022) [26]. The relative frequency of the corresponding amino acid at each position was calculated on the Web Logo online website (http://weblogo.threeplusone.com/, accessed on 18 March 2022). We conducted motif analysis of the detoxification genes on the MEME online server (http://meme-suite.org/tools/meme, accessed on 20 March 2022) with parameters "minimum width = 6, maximum width = 50, number of motif to find = 10" [27], and the analysis of conserve domain by the NCBI Conserved Domain Database (NCBI-CDD) (e-value = 1 × 10 −3 ) [23]. The exon and intron structures were displayed in all protein sequences using the Gene Structure Display Server (GSDS) (http://gsds.cbi.pku.edu.cn/, accessed on 22 March 2022) [28]. Finally, the corresponding results were visualized by the program TBtools [29].

Chromosomal Locations, Collinearity and Selection Pressure
The detoxification genes' positions on chromosomes were displayed using TBtools (version 1.09876, created by Chen et al.; Wuhan, China) [29]. We conducted and visualized collinearity analysis by MCScanX (created by Wang et al.; Wuhan, China) [30] and Circos (version 2.50, created by Krzywinski et al., Vancouver, Canada) [31]. In addition, we used Ka Ks_Calculator 2.0 (version 2.0, created by Zhang et al.; Beijing, China) to calculate the ratio between the non-synonymous replacement rate (Ka) and the synonymous replacement rate (Ks) of two protein-coding genes, which is an indicator of molecular evolution to determine whether there is selective pressure on a protein-coding gene [32].

Expression Profile of Detoxification Genes
We extracted total mixed RNA from the three S. chinensis samples by the Trizol method [33], and constructed a sequence library using the Illumina TruseqTM RNA sample prep Kit [34]. On the Illumina HiSeq2500 platform, we carried out next generation sequencing for the sequence library. We measured the original raw reads and removed the low-quality, repetitive data with adapters to obtain clean reads. Finally, we assembled the transcriptome data from scratch to obtain the unigenes sequence by the Trinity assembly software [35]. Meanwhile, we compared unigenes sequence with six major databases, i.e., NR, Swiss-Prot, Pfam, COG, GO, and KEGG, to complete functional annotation and classification analysis [36][37][38][39][40][41][42][43]. We screened detoxification-related genes with our gene annotation lists. We compared the CDS sequences of the genomic detoxification genes with the transcriptome detoxification gene by MAFFT Alignment to determine the gene expression [44].

Identification of Detoxification Genes of S. chinensis
We identified nine genes in GSTs, 55 genes in ABCs, 18 genes in CCEs, 48 genes in CYPs and 14 genes in UDPs in the S. chinensis genome (Table 1). These genes were divided into three subfamilies for GSTs, i.e., delta, theta, and sigma, six subfamilies for ABCs, i.e., A, B, C, D, E, F, and G, four subfamilies for P450s, i.e., CYP2, CYP6, CYP4, and mitochondrion, two subfamilies UGT4 and UGT5 for UDPs. Detoxification gene numbers displayed heterogeneity in the annotated Aphidinae genomes. For example, A. pisum displayed more genes of GSTs than other aphid species due to having more genes in delta subfamily. S. chinensis had the greatest numbers in ABCs due to more genes in the C subfamily, and fewer in CCEs and CYPs due to fewer genes in Esterase and CYP6 subfamilies, respectively. These gene families, i.e., GSTs, ABCs, and CCEs, clearly exhibited expansion and contraction.

Characteristic of the Five Detoxification Genes of S. chinensis
We performed a characteristic analysis of the detoxification genes, including motif, domain, and the number of exons, and constructed the phylogenetic tree of the protein sequences of five detoxification gene families from the four aphid species S. chinensis, A. pisum, C. cedri, and M. persicae .

CCEs
The phylogenetic tree of the CCEs divided the detoxification genes into three subfamilies, i.e., Esterase, Acetylcholinesterase and Neuroligin (Figure 3a). The motif order of the three classes was different, but most of them had the common order 2-6-4-1-3-8-7-10-9 ( Figure 3b). Ten conserved motifs of CCEs varied from 15 to 45 amino acids in length ( Figure 2b). The conserved domain of CCEs included coesterase, Abhydrolase and esterase _lipase ( Figure 3c). The analysis on the coding sequence (CDS) and untranslated regions (UTRs) of the CCEs showed that exon numbers ranged from two to 18, and a total of 41 members exhibited 5 and 3 UTRs, eight members presented no UTR and the remaining seven members had either a 5 or 3 UTR (Figure 3d). Differences in the motif pattern and gene structure of the different classes might be the reason for the differences in their physiological functions.

CCEs
The phylogenetic tree of the CCEs divided the detoxification genes into three subfamilies, i.e., Esterase, Acetylcholinesterase and Neuroligin ( Figure 3a). The motif order of the three classes was different, but most of them had the common order 2-6-4-1-3-8-7-10-9 ( Figure 3b). Ten conserved motifs of CCEs varied from 15 to 45 amino acids in length ( Figure 2b). The conserved domain of CCEs included coesterase, Abhydrolase and esterase _lipase ( Figure 3c). The analysis on the coding sequence (CDS) and untranslated regions (UTRs) of the CCEs showed that exon numbers ranged from two to 18, and a total of 41 members exhibited 5′and 3′ UTRs, eight members presented no UTR and the remaining seven members had either a 5′or 3′ UTR ( Figure 3d). Differences in the motif pattern and gene structure of the different classes might be the reason for the differences in their physiological functions.

UDPs
The phylogenetic tree divided the protein sequence of UDPs into two subfamilies, i.e., UGT4 and UGT5 ( Figure 4a). All genes of UDPs showed the same motif order 3-7-10-4-6-5-8-1-2-9, only few genes lacked a 1-2 motif, i.e., Schi02G001900, Schi01G030560 ( Figure 4b). Ten conserved motifs of UDPs varied from 15 to 50 amino acids in length ( Figure 2c). The conserved domain of UDPs included Glycosyltransferase _GTB, GT1_Gtf and MCS (Figure 4c). Analysis of the coding sequence (CDS) and untranslated regions (UTRs) of the UDPs showed that the numbers of exons ranged from three to 11. Only two genes had 11 exons, 27 UDP genes (17.3%) had five exons, and 25 UDP genes (81.7%) had four exons. A total of 43 members exhibited 5 and 3 UTRs, 13 members presented no UTR, and the remaining five members had either 5 or 3 UTR. The result shows that the same class has the same number of exons in the UDP gene family (Figure 4d).

GSTs
The phylogenetic tree divided the protein sequence of GSTs into three class, i.e., Delta, Theta and Sigma (Figure 5a). The delta subfamily included four genes, i.e., Schi05G007220, Schi05G002040, Schi05G002060, and Schi05G007230, which shared the same motif order (1-6) as the nine genes in the Sigma subfamily. Schi01G005590, Schi01G013290, Schi01G013300, Schi03G000690 and Schi03G000710 belonged to theta subfamily and shared the same motif order (3-1-4-2-5) (Figure 5b). Details of six putative motifs are outlined in Figure 2d. These conserved motifs ranged from 15 to 50 amino acids in length. The conserved domain of GSTs included GST _ Delta _ Epsilon, Gst A, GST _Theta, GST _Sigma and Thioredoxin (Figure 5c). In the delta and theta subfamily, the numbers of most exons were four, a few exons were just one. In the sigma subfamily, the numbers of most exons were 11, and other members had five exons (Figure 5d). The genes in the same groups had similar motif patterns and numbers of exons, indicating that they were highly conserved and the inferred functions were similar.

GSTs
The phylogenetic tree divided the protein sequence of GSTs into three class, i.e., Delta, Theta and Sigma (Figure 5a). The delta subfamily included four genes, i.e., Schi05G007220, Schi05G002040, Schi05G002060, and Schi05G007230, which shared the same motif order (1-6) as the nine genes in the Sigma subfamily. Schi01G005590, Schi01G013290, Schi01G013300, Schi03G000690 and Schi03G000710 belonged to theta subfamily and shared the same motif order (3-1-4-2-5) (Figure 5b). Details of six putative motifs are outlined in Figure 2d. These conserved motifs ranged from 15 to 50 amino acids in length. The conserved domain of GSTs included GST _ Delta _ Epsilon, Gst A, GST _Theta, GST _Sigma and Thioredoxin (Figure 5c). In the delta and theta subfamily, the numbers of most exons were four, a few exons were just one. In the sigma subfamily, the numbers of most exons were 11, and other members had five exons (Figure 5d). The genes in the same groups had similar motif patterns and numbers of exons, indicating that they were highly conserved and the inferred functions were similar.

ABCs
The phylogenetic tree of the protein sequences of ABCs of S. chinensis were distributed to G, C, A, E, C, B, F and D subfamilies (Figure 6a). The motif order for G class was 2-5-7-1-3-4-10-6-8-9. The C class motif order was 2-5-7-1. A class motif order for most members was 2-5-7-1, with a few members (Schi03G001540) having special motif 8. The E, C, B, F and D class motif order was 2-7-1 (Figure 6b). The lengths of ABC conserved motifs ranged from 28 to 41 amino acids (Figure 2e). The conserved domain of ABCs included 3a1204, CcmA, YadH, Rli1, Uup, ATM1, MTABC and ABC2_membrane_3 ( Figure  6c). The numbers of exon ranged from six to 28 by the predictions of the gene structure. The exon number of most members was 14 -28. Analysis of the coding sequence (CDS) and untranslated regions (UTRs) of the ABC gene family found that just six members did not have a UTR (Figure 6d).

ABCs
The phylogenetic tree of the protein sequences of ABCs of S. chinensis were distributed to G, C, A, E, C, B, F and D subfamilies (Figure 6a). The motif order for G class was 2-5-7-1-3-4-10-6-8-9. The C class motif order was 2-5-7-1. A class motif order for most members was 2-5-7-1, with a few members (Schi03G001540) having special motif 8. The E, C, B, F and D class motif order was 2-7-1 (Figure 6b). The lengths of ABC conserved motifs ranged from 28 to 41 amino acids (Figure 2e). The conserved domain of ABCs included 3a1204, CcmA, YadH, Rli1, Uup, ATM1, MTABC and ABC2_membrane_3 (Figure 6c). The numbers of exon ranged from six to 28 by the predictions of the gene structure. The exon number of most members was 14-28. Analysis of the coding sequence (CDS) and untranslated regions (UTRs) of the ABC gene family found that just six members did not have a UTR (Figure 6d).

Chromosomal Location and Collinearity of Detoxification Genes of S. chinensis
The location and collinearity analysis of all detoxification genes showed that 144 genes were unevenly distributed on chromosomes 1-9 (Figure 7a). Chromosome 1 had the most members of detoxification genes with 33 genes, among which there were eight genes in ABCs, four genes in CCEs, three genes in GSTs, six genes in P450s and 12 genes in UDPs. A total of two pair genes showed tandem duplications on Chromosome 1. The two genes Schi01G040920 and Schi01G044340 of the P450 family, and Schi01G004180 and Schi01G00620 of the ABCs showed tandem duplications, respectively. Only two pairs of genes in the detoxification genes had collinearity (Figure 7b). Chromosome 7 had the fewest detoxification genes with only one member. The distribution of detoxification genes on chromosomes with no bias to the 5 or 3 ends may be related to their function.

Chromosomal Location and Collinearity of Detoxification Genes of S. chinensis
The location and collinearity analysis of all detoxification genes showed that 144 genes were unevenly distributed on chromosomes 1-9 ( Figure 7a). Chromosome 1 had the most members of detoxification genes with 33 genes, among which there were eight genes in ABCs, four genes in CCEs, three genes in GSTs, six genes in P450s and 12 genes in UDPs. A total of two pair genes showed tandem duplications on Chromosome 1. The two genes Schi01G040920 and Schi01G044340 of the P450 family, and Schi01G004180 and Schi01G00620 of the ABCs showed tandem duplications, respectively. Only two pairs of genes in the detoxification genes had collinearity (Figure 7b). Chromosome 7 had the fewest detoxification genes with only one member. The distribution of detoxification genes on chromosomes with no bias to the 5′ or 3′ ends may be related to their function. The collinear comparison map of the detoxification gene family was established by MC Scan X between S. chinensis and A. pisum (Figure 7c). There were 38 pairs of collinearity (homologous gene pairs) in A. pisum, including nine in P450s, 16 in ABCs, three in CCEs, three in GSTs, and seven in UDPs. There were more homologous gene pairs of ABCs and P450s in the S. chinensis and A. pisum, which may be related to the large number of these two gene families. It is inferred that these two gene families are relatively conservative in the evolutionary process and have relatively stable functions.
The Ka/Ks analysis of 38 pairs of homologous genes existing in S. chinensis and A. pisum was carried out, and the results are shown in Table 2. The ratios of Ka/Ks between gene pairs were all <1, which indicated that the detoxification genes in S. chinensis were mainly affected by purification selection during evolution.

Expression Profiles of Detoxification Genes in S. chinensis
We examined 61 detoxification genes in the transcriptome data of S. chinensis, among which there were 13 genes in ABCs, nine genes in CCEs, seven genes in GSTs, 19 genes in P450s and 11 genes in UDPs, respectively. The gene number expressed in the P450s was the highest and the gene number expressed in GSTs was the least (Table 3). In terms of gene expression, the overall expression of CCEs was the highest, while that of P450s was the lowest. Individual gene expression was particularly high in some gene families. For example, the expression level of Schi01G005590 (GSTs) was 251.68, Schi01G003300 (P450s) was 166.92., and Schi05G003270(CCEs) was 127.19, from which we infer that these genes play an important role in the detoxification process of S. chinensis. The collinear comparison map of the detoxification gene family was established by MC Scan X between S. chinensis and A. pisum (Figure 7c). There were 38 pairs of collinearity (homologous gene pairs) in A. pisum, including nine in P450s, 16 in ABCs, three in CCEs, three in GSTs, and seven in UDPs. There were more homologous gene pairs of ABCs and P450s in the S. chinensis and A. pisum, which may be related to the large number of these two gene families. It is inferred that these two gene families are relatively conservative in the evolutionary process and have relatively stable functions.
The Ka/Ks analysis of 38 pairs of homologous genes existing in S. chinensis and A. pisum was carried out, and the results are shown in Table 2. The ratios of Ka/Ks between gene pairs were all < 1, which indicated that the detoxification genes in S. chinensis were mainly affected by purification selection during evolution.

Physicochemical Properties of Detoxification Gene Products in S. chinensis
The molecular weight of the detoxification genes products ranged from 20,632.58 Da to 199,951.12 Da (Table 4), the number of amino acids from 178 to 1764, and the aliphatic index from 74.78 to 112.86. The isoelectric points (pI) ranged from 4.93 to 9.58. The instability index of which 51% of the detoxification gene products were more than 40, indicated that those genes were unstable and easily degraded. The aliphatic amino acid index of detoxification gene proteins was 74.78-112.86, and the average of the hydropathicity indices ranged from −0.81 to 0.306. The proteins with hydropathicity value less than 0 were hydrophilic proteins.

Prediction of Protein Multi-Level Structures of Detoxification Gene Products in S. chinensis
Since the members of detoxification gene family products have similar protein structures, we selected the protein members of GSTs to predict secondary and tertiary structures. The predicted secondary structures of the nine GST proteins are shown in Table 5. The GST protein was composed of four parts: α-helix, extended chain, β-turn and random coil, among which the α helix ratio of the GST protein was the highest, followed by random coil, and the ratio of β-turn was the smallest. Our Signal P prediction showed that there were no signal peptides in all the GST proteins of S. chinensis. The tertiary structures of GST proteins are shown in Figure 8. Schi01G005590, Schi01G013290, Schi01G013300, Schi03G000690 and Schi03G000710, Schi05G002040 and Schi05G002060, Schi05G007220, and Schi05G007230 have similar protein structures, respectively. These genes belonged to sigma, theta, and delta subfamilies, respectively. The result indicated that gene members belonging to the same subfamily have similar protein structures, indicating they have similar biological functions.

Expansion and Contraction of Detoxification Genes in S. chinensis
The sequenced genomes of A. pisum, M. persicae, and C. cedri with manually annotated detoxification gene families served as comparisons for our studies [50,51]. The number of detoxification genes in C. cedri was derived from the annotation table in NCBI with Accession No. GCA_902439185.1, and those of A. pisum and M. persicae were derived from published articles [16,[52][53][54].
The number of detoxification genes of S. chinensis was 144, less than that of A. pisum, C. cedri and M. persicae, which might be related to its unique host plant. It was predicted that polyphagous insects require a greater complement of detoxification-related enzymes for they were usually exposed to a higher diversity of plant secondary metabolites than oligophagous ones [54,55]. The gene number of GSTs and UDPs in S. chinensis was nine and 14, respectively, which were relatively less and conserved by comparison with the other three aphids. The fewer detoxification-related genes may be due to the unobvious duplication events of detoxification-related genes. Genes with conserved roles usually have relatively stable copies, while those with diversified functions have higher rates of gain-and-loss with random changing degrees [56]. The numbers of ABCs of S. chinensis was 55, that was more than other aphids and showed significant genetic expansion. It is inferred that the ABCs play an important role in the degradation of secondary metabolites of S. chinensis. The number of P450s was 48 in S. chinensis, which was clearly contracted by comparison with the other three aphids. CYPs in many insects are associated with the

Expansion and Contraction of Detoxification Genes in S. chinensis
The sequenced genomes of A. pisum, M. persicae, and C. cedri with manually annotated detoxification gene families served as comparisons for our studies [50,51]. The number of detoxification genes in C. cedri was derived from the annotation table in NCBI with Accession No. GCA_902439185.1, and those of A. pisum and M. persicae were derived from published articles [16,[52][53][54].
The number of detoxification genes of S. chinensis was 144, less than that of A. pisum, C. cedri and M. persicae, which might be related to its unique host plant. It was predicted that polyphagous insects require a greater complement of detoxification-related enzymes for they were usually exposed to a higher diversity of plant secondary metabolites than oligophagous ones [54,55]. The gene number of GSTs and UDPs in S. chinensis was nine and 14, respectively, which were relatively less and conserved by comparison with the other three aphids. The fewer detoxification-related genes may be due to the unobvious duplication events of detoxification-related genes. Genes with conserved roles usually have relatively stable copies, while those with diversified functions have higher rates of gain-and-loss with random changing degrees [56]. The numbers of ABCs of S. chinensis was 55, that was more than other aphids and showed significant genetic expansion. It is inferred that the ABCs play an important role in the degradation of secondary metabolites of S. chinensis. The number of P450s was 48 in S. chinensis, which was clearly contracted by comparison with the other three aphids. CYPs in many insects are associated with the metabolism or detoxification of key endogenous substrates and xenobiotics, such as steroid hormones and lipids, plant natural products, and pesticides, which are key components for the successful adaptation of insects to their host plants [57].

Characteristics and Expression of the Detoxification Genes of S. chinensis
The current study compared characteristics of the five detoxification gene families of S. chinensis with other three aphids, A. pisum, M. persicae, and C. cedri. Ten conserved motifs were found in four gene families, except for GSTs with six conserved motifs. The conserved motifs in detoxification genes are very important in the functional domains, and the highest one was its key structure, while the motif patterns can finely tune the function of detoxification genes [58]. Structural variation affects gene evolution [59]. The detoxification gene proteins encoded by the subfamily members usually had the same motif orders and exons, which indicated that genes in closely related groups were highly conserved and might have similar functions.

Collinearity, Chromosome Position and Evolutionary Rate of Detoxification Genes of S. chinensis
We performed analysis on collinearity relationships to further investigate the gene duplication events within detoxification genes. There were two pairs of gene duplications in the detoxification genes of S. chinensis, including P450s and ABCs. The differential expansion events during the species evolution might result in the phenomena that the number of family members was not correlated with genome size. Tandem duplications generate a large number of genes, which is considered as the most effective mechanism for producing and maintaining gene copies [57]. It was reported that gene duplications were critical for the evolution of new genes and novel functions and were major forces driving gene family expansion [60]. For example, CYP genes, often clustering in genomes, were considered as a result of gene duplication events [56]. We calculated the Ka and Ks values to estimate the evolutionary trend and revealed a functional selection pressure between duplicated gene pairs. The ratio of Ka to Ks in protein-coding genes can determine whether there is selection pressure acting in the process of gene evolution. Detoxification genes of S. chinensis showed a strong purifying selection during evolution, which suggested that their functions may be evolutionarily conserved.

Conclusions
Here, we performed a comprehensive genome-wide analysis of the detoxification gene family in S. chinensis and compared the results with the genomes of A. pisum, M. persicae, and C.cedri. We manually annotated 144 genes of S. chinensis, including nine in GSTs, 55 in ABCs, 18 in CCEs, 48 in P450s, and 14 in UDPs. We constructed the phylogenetic trees, motif patterns, domains, and gene structures of detoxification genes from these four aphids and further analyzed the chromosomal location, collinearity, evolution rates, and their expression. Finally, we predicted characteristics, physicochemical properties, and protein multi-level structures of detoxification gene products of S. chinensis. Our results provide comprehensive information, molecular data, and gene candidates for further analyses. S. chinensis can survive in galls with the tannin content up to 70%, and we infer that it has a strong ability to reduce secondary metabolites. The phenomenon may be related not only to its own detoxification genes, but also to the existence of endosymbionts or a long-term obligate parasitic relationship.
Author Contributions: H.H. and Z.R. designed the study, interpreted all the data and findings and wrote the manuscript, and made equal contributions as major authors, while M.J.C.C. validated, revised, and edited the final manuscript. All authors have read and agreed to the published version of the manuscript.