Transcriptomics and Metagenomics of Common Cutworm (Spodoptera Litura) and Fall Armyworm (Spodoptera Frugiperda) Reveal Their Differences on Detoxi cation and Development


 BackgroundSpodoptera litura is an important polyphagous pest that causes great damage to agriculture. We performed RNA-seq from 15 individuals of S. litura, including different larval (fifth and sixth instar larvae), chrysalis, and adult developmental stages. We also sequenced the Spodoptera frugiperda to compared with the S. litura. ResultsA total of 101,885 differentially expressed transcripts (DETs) were identified in the S. litura. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses indicated that the S. litura may undergo active development and basic metabolism, such as xenobiotic and detoxifying metabolism, during its larvae and adult stages, which explains the difficulty to control it. We also found that DETs of single-copy orthologous genes between S. litura and S. frugiperda were involved in basic metabolism and developmental, but S. litura were highly enriched in energy and metabolic processes, whereas S. frugiperda had stronger nervous and olfactory functions. Metagenomics analysis in larval S. litura and S. frugiperda revealed that microbiota can both participate in the detoxification and metabolism processes in them, but the relative abundance of detoxification-related microbiota was more abundant in S. frugiperda. Transcriptome results also confirmed the detoxification-related pathway of S. frugiperda was more abundant than that of S. litura. ConclusionsTherefore, we can conclude that the detoxification ability of S. frugiperda is stronger than that of in S. litura, which is related to the wide adaptation and strong resistance of S. frugiperda. Our study might provide clues to the development of novel pest control strategies.

As a destructive pest, numerous studies have been done to reveal the characteristic of S. litura. Previous studies primarily evaluated on the in uence of different host plant on the growth and development of S. litura [8-10] and its potential damage for different crops [2]. In addition, because of the poor e cacy of various insecticides against S. litura [3], researchers have also done a lot of investigations on its resistance to insecticides. Transcriptome analysis of S. litura uncover that P450s may be involved in the detoxi cation of uralaner in vivo [11]. Previous study found that the detoxion-related gene families (P450, GST, C0E, APN and ABC) were massively expanded in S. litura, which provided an explanation for the genetic basis for its high tolerance to pesticides [12]. Despite an enormous volume of research work done to understand the ecology and genetics of S. litura, the developmental characteristics of S. litura and the underlying changes about resistance to pesticides with development were still poorly understood. RNA sequencing has been widely used to obtain expression data throughout different developmental stages on agricultural pests and explore the key regulated genes related to development that can be targeted for the pest controlling. Through pairwise comparison of four developmental stages of Athetis lepigone, which was found that some differentially expressed genes were related to cuticle and wing formation as well as the growth and development [13]. Transcriptome analysis of Ostrinia furnacalis in four developmental stages revealed genes that associated with developmental pathways, cuticularization, wing formation and olfactory recognition [14]. To design e cient strategies to control this economic pest, it is necessary to carry out RNA-seq on S. litura to compare the differences in gene expression pro les at different stages.
The fall amyworm (Spodoptera frugiperda) is phylogenetically close to S. litura and belonged to the same family Noctuidae. In the larval stage, S. frugiperda experiences natural molting, extensive feeding, and phenotypic changes and during the chrysalis stage, it undergoes active energy and nutrient metabolism [15]. Despite the biological characteristics and morphological characteristics of each insect form were similar, for example, in morphology, high resistance to pesticides, strong reproducibility and migration [16], there exist a series of differences in feeding habit, reproductive behavior and damage degree between S. frugiperda and S. litura [17]. Spodoptera frugiperda prefers gramineous plants, whereas S. litura may prefer dicotyledons. Another difference is that S. frugiperda has more powerful locomotivity than S. litura. In addition, S. frugiperda shows higher mating frequency than that of S. litura. The genetic divergence between the two pests was revealed by comparative genomics [16], which found genes of S. frugiperda in the taste sensory perception and nervous system were enriched more abundantly than S. litura. Moreover, our team studied the gene expression differences between different developmental stages of S. frugiperda [15], which demonstrated that S. frugiperda experienced active metabolism, detoxifying and xenobiotic metabolism throughout its life, especially in the larval stage. Considering the observed differences in feeding habit and reproductive behavior between S. frugiperda and S. litura, it is possible that there are gene expression differences between them even during the same development stage. Therefore, comparative transcriptomics between the two pests at different development stages are helpful to better understand the genetic differences between them.
In addition, the insect gut microbes could contribute to the food digestion and protection against pathogens. Analysis of the metagenome of the larvae S. litura revealed that microbiota played a major role in digestion, detoxi cation, and nutrient supply [18,19]. Spodoptera frugiperda gut bacterial participated in modulate plant defense responses [20]. However, there is not yet available study about gut microbial differences between S. litura and S. frugiperda. Therefore, we aim to investigation the differences of gut micro ora composition and functional annotation of S. litura and S. frugiperda using metagenomic.
In this study, a total of 15 S. litura individuals, including larval stage ( fth instar larvae and sixth instar larvae), chrysalis stage and adult stage (female and male), were collected and were performed RNA-seq as well as metagenomic on larval stage. We identi ed differently expressed transcripts (DETs) between different developmental stages and identi ed metabolic pathways that changed with development. Furthermore, we also collected and re-analyzed the data of S. frugiperda from Wang et al., to identify the potential different genetic mechanisms with development between S. litura and S. frugiperda through comparative analysis [15].We aimed to reveal the gene expression changes between different developmental stages of S. litura, and the gene expression differences between S. litura and S. frugiperda in the same developmental stage. The metagenomics data showed there were signi cant differences in the microbial communities and function between S. litura and S. frugiperda. This study provided new insight into understanding and utilizing the genetic differences for formulating of personalized and phased control strategies of the two invasive pests.

Transcriptome sequencing and alignment
To explore the transcript expression patterns of S. litura at different developmental stages, we collected 15 S. litura samples from three developmental stages, larvae ( fth instar larvae and sixth instar larvae), chrysalis, and adults (females and males). Transcriptome sequencing of the cDNA library was performed on the Illumina HiSeq 2000 platform. A total of 446,991,085 raw reads and 134.1 Gb data were obtained from all samples. After quality control of the transcripts, 436,543,284 (97.66%) clean reads were obtained. Then we mapped the clean reads to the reference genome, the results showed that the mapping rate was 81.31% ~ 94.96% (Supplementary Table S1).

Transcript expression analysis
The transcripts which had 5 or more reads in any sample could be mapped to the reference genome were used for further analysis. We found that 23,020 transcripts were expressed in abundance, which accounted for 88.96% of the known transcripts. According to the developmental stage of S. litura, we combined the fth and sixth instar larvae as the larva group. Next, DEseq2 was used to standardize the raw expression matrix and then perform the differential transcripts analysis between groups. As shown in the principal component analysis (PCA) that the larvae, chrysalis, and adult groups of S. litura were clearly separated. At the same time, there was also a clear separation between the female and male adults ( Figure S1).
The developmental difference in the S. litura Through cross-comparison of DETs at three developmental stages, a total of 2,703 downregulated DETs and 4,270 upregulated DETs were identi ed in the LL group when comparing to LC. There were 3,145 downregulated DETs and 1,847 upregulated DETs in the LC group when comparing to LA. In addition, 2,979 downregulated DETs and 3,629 upregulated DETs were identi ed in the LL group when comparing to LA ( Figure 1A).
In order to understand the biological role of DETs, GO and KEGG function enrichment analysis was performed. The upregulated DETs in the LL group (i.e., LL vs. LC) were mainly enriched in biological process (BP) and molecular function (MF) GO terms, such as lipid catabolic process (GO:0016042), carbohydrate metabolic process (GO:0005975), glutathione metabolic process (GO:0006749) and trehalose transport (GO:0015771) in BP terms, which all related to basic metabolism processes. Besides, many enzyme-related MF terms were enriched by these DETs, such as glutathione transferase activity (GO:0004364) and oxidoreductase activity (GO:0016491) (Supplementary Table   S2). In KEGG enriched metabolic pathways, pathways related to detoxi cation metabolism were found, such as drug metabolism -cytochrome P450 (ko00982) and drug metabolism -other enzymes (ko00983) ( Figure 1B and Supplementary Table S2). Downregulated DETs in the LL group were mainly enriched in the development BP terms, such as locomotion (GO:0040011) and muscle contraction (GO:0006936). In MP terms, many cell activity processes were enriched, such as structural constituent of muscle (GO:0008307), actin lament binding (GO:0051015) and actin-dependent ATPase activity (GO:0030898) (Supplementary Table S2). Moreover, many signaling KEGG pathways, such as ErbB signaling pathway (ko04012), Rap1 signaling pathway (ko04015) and Hippo signaling pathway (ko04390) were enriched (Supplementary Table S2).
In the LC group (i.e., LC vs. LA), upregulated DETs in the LC group were mainly enriched in development BP terms, including muscle organ development (GO:0007517), tissue development (GO:0009888) and chitin metabolic process (GO:0006030). In MP terms, many terms were associated with chitin metabolism, such as chitin binding (GO:0008061), structural constituent of pupal chitin-based cuticle (GO:0008011) and chitinase activity (GO:0004568) (Supplementary Table S3). Some basic metabolisms in the KEGG pathways were enriched, such as citrate cycle (TCA cycle) (ko00020) and oxidative phosphorylation (ko00190) (Supplementary Table S3). The majority downregulated DETs were enriched in energy-related and behavior BP terms, including trehalose transport (GO:0015771), disaccharide transport (GO:0015766) and chorion-containing eggshell formation (GO:0007304). MP terms which were associated with enzyme activity were abundantly enriched, such as serine-type endopeptidase activity (GO:0004252) and serine hydrolase activity (GO:0017171) (Supplementary Table S3). For KEGG pathways analysis, many detoxi cation metabolism pathways were enriched, including metabolism of xenobiotics by cytochrome P450 (ko00980), drug metabolism -cytochrome P450 (ko00982) and drug metabolism -other enzymes (ko00983) (Supplementary Table S3). Upregulated DETs in the LL group (i.e., LL vs. LA) were mainly enriched in developmental BP terms, such as axoneme assembly (GO:0035082) and cytoplasmic translation (GO:0002181). In MF terms, many terms were associated with enzyme activity were enriched, including glutathione transferase activity (GO:0004364) and serine-type endopeptidase activity (GO:0004252) (Supplementary Table S4). In KEGG analysis, many pathways related to energy metabolism and detoxi cation metabolism were enriched, such as citrate cycle (TCA cycle) (ko00020) and metabolism of xenobiotics by cytochrome P450 (ko00980) (Supplementary Table S4). Downregulated DETs in the LL group were mainly enriched in BP terms which associated with nervous system, including axon guidance (GO:0007411), neuron differentiation (GO:0030182) and nervous system development (GO:0007399) (Supplementary Table  S4). There was only one KEGG pathway, lysosome (ko04142), which was enriched by these downregulated DETs (Supplementary Table S4).
In general, many BP terms which associated with development and basic metabolism and MF terms related to enzyme activity were enriched. We also found many terms were involved in chitin metabolism. Furthermore, some KEGG pathways related to detoxi cation were found.
The gene expression differences between the S. litura and S. frugiperda In order to analyze the gene expression differences between S. litura and S. frugiperda at the same age stage, we performed DETs analysis on the shared single-copy orthologous genes. The orthologous genes of two species were selected by the OrthoFinder, and 6,735 pairs of single-copy orthologs were identi ed after strict data screening. Among them, 6,728 were co-expressed in both S. litura and S. frugiperda. Gene expression between the two species can be quanti ed and compared after the length correction of corresponding gene CDS region. After that, we analyzed the differential expression of these single-copy orthologous genes in the same developmental stage between the two species, which identi ed 454 downregulated DETs and 440 upregulated DETs in the LL group when comparing to FL. A total of 493 downregulated DETs and 589 upregulated DETs were identi ed in the LC group when comparing to FC. In addition, 658 downregulated DETs and 611 upregulated DETs were identi ed in the LA group when comparing to FA.
DETs of single-copy orthologous genes between S. litura and S. frugiperda could reveal the differences in gene expression between these phylogenetic close related species ( Figure 1C). We further performed GO and KEGG enrichment analyzing on DETs of single-copy orthologous genes. The results showed that upregulated DETs in the LL group (i.e., LL vs. FL) were only enriched in one MF terms, structural constituent of cuticle (GO:0042302). There was no signi cant KEGG pathway enriched in upregulated DETs in the LL group. The GO analysis showed that downregulated DETs in the LL group were enriched in one BP term, olfactory learning (GO:0008355) and one CC term, synapse (GO:0045202), they were associated with sensory and nervous system development. There were four KEGG pathways enriched by these DETs, Rap1 signaling pathway (ko04015), axon regeneration (ko04361), adherens junction (ko04520) and tight junction (ko04530) (Supplementary Table S5).
The GO analysis showed that upregulated DETs in the LC group (i.e., LC vs. FC) were mainly enriched in, for example, extracellular space (GO:0005615), extracellular region (GO:0005576) and serine-type endopeptidase activity (GO:0004252). There was only one enriched KEGG pathway, fanconi anemia pathway (ko03460). Downregulated DETs in the LC group were enriched in nervous system developmental related CC terms, such as neuromuscular junction (GO:0031594). The nervous development and cell activity KEGG pathways were associated with these downregulated DETs included the axon regeneration (ko04361) and tight junction (ko04530) ( Supplementary Table S6).
Lastly, the GO annotation analysis of upregulated DETs in the LA group (i.e., LA vs. FA) were mainly enriched in energy related BP terms such as mitochondrial translation (GO:0032543), ATP synthesis coupled proton transport (GO:0015986) and aerobic respiration (GO:0009060), as well as binding-related MF terms, including chitin binding (GO:0008061) and heme binding (GO:0020037) ( Figure  1D). The metabolism-related KEGG pathways associated with these DETs included citrate cycle (TCA cycle) (ko00020), oxidative phosphorylation (ko00190) and 2-oxocarboxylic acid metabolism (ko01210) (Supplementary Table S7). The GO and KEGG enrichment analysis showed that there was no signi cant pathway enriched of downregulated DETs in the LA group.
In summary, the GO and KEGG enrichment analysis showed that the DETs of single-copy orthologous genes between S. litura and S. frugiperda were involved in basic metabolism and developmental. At the larvae stage, downregulated DETs in S. litura were enriched in a few pathways related to nervous system developmental and cell activity. At the chrysalis stage, there were more pathways associated with development in S. frugiperda than in S. litura. At the adult stage, the upregulated DETs of single-copy orthologous genes in S. litura were enriched with more terms and pathways related to energy metabolism than those of S. frugiperda.
Next, we focused on the detoxi cation gens. Certain gene families of encoding detoxi cation enzymes, such as cytochrome P450, carboxylesterase, and glutathione-s-transferase, were found expansion in S. litura, which would contribute to facilitate its adaptability to the hosts [21]. In order to further explore the difference of detoxi cation ability between S. litura and S. frugiperda, we compared the expression levels of their homologous genes involved in detoxi cation pathways related to P450, GST, and carboxylesterase, and the heatmap showed the expression of homologous genes involved in detoxi cation pathways shared by both species at different ages ( Figure 2). The results showed that the gene expression levels involved in detoxi cation pathways were different in three age groups of the same species. For example, the expression of genes involved in the glutathione transferase activity pathway was signi cantly higher at larval stage than that at other ages in S. litura ( Figure 2D). In addition, the gene expression levels involved in detoxi cation pathways were also different in different species at the same age group. For example, the expression of genes involved in the metabolic of xenobiotics by cytochrome P450 pathways was highest at the larval stage in S. litura, while the expression of genes involved in this pathway was high at almost all three age stages in S. frugiperda ( Figure 2D).
Microbiota composition differs between S. litura and S. frugiperda Metagenome sequencing analysis of larvae ( fth instar larvae and sixth instar larvae) stage of S. litura and S. frugiperda was carried out and we abtained 474,122,459 and 479,748,778 valid reads from 6 S. litura and 6 S. frugiperda samples, respectively (Supplementary Table S8). A total of 27 phyla, 49 classes, 114 orders, 152 families, 750 genera and 2,047 species were identi ed in S. litura, while 26 phyla, 52 classes,124 orders, 167 families, 881 genera and 2,631 species were found in S. frugiperda samples. Firmicutes, Proteobacteria and Actinobacteria were the dominant phyla in all samples. The proportion of Firmicutes was higher in S. litura compared with S. frugiperda (95.88% and 93.00%), whereas the relative percentage of Proteobacteria composition was higher in S. frugiperda (1.30% and 3.00%; Figure 3A). At the genus level, Enterococcus, Pediococcus and Weissella were the dominant genera of S. litura, while Enterococcus, Corynebacterium and Bacillus were the most abundant genera of S. frugiperda ( Figure 3B). We conducted a difference analysis for genus, the relative abundance of Pediococcus were signi cantly higher in S. litura than in S. frugiperda, however, the relative abundance of Corynebacterium, Clostridium, Glutamicibacter and so on were signi cantly higher in S. frugiperda compared with S. litura (LDA>3, p<0.05; Figure 3C). A heatmap was constructed based on the abundance of top 50 genera of each sample ( Figure 3D) also showed signi cant differences at the genus level.

Metagenomic functional analysis
To analyses the functional relationships of the composition of S. litura and S. frugiperda metagenome, genes were predicted from the CAZy, Humann3 and CARD database. A total of 302 annotated genes from six families, which contained glycoside hydrolases (GHs), glycosyltransferases (GTs), polysaccharide lyases (PLs), carbohydrate esterases (CEs), auxiliary activities (AAs), and carbohydrate binding modules (CBMs), were detected in the gut microbiota based on the Carbohydrate-Active Enzymes (CAZy) database (http://www.cazy.org). In S. litura, we obtained 170 GH families, 74 GT families, 48 CBM families, 17 CE families, 24 PL families, and 16 AA families. However, a search using S. frugiperda samples identi ed 134 GH families, 64 GT families, 40 CBM families, 17 CE families, 14 PL families and 11 AA families. The most abundant CAZyme families detected were represented by families GH38 and GH85 in GH families.
Moreover, LEfse identi ed GT31, GH13_25 and GT10 as the top three signi cantly higher CAZymes in the S. frugiperda than those in S. litura, while GH13_14, GH67 and GH94 as the top three signi cantly higher CAZymes in S. litura than those in S. frugiperda (LDA>2, P<0.05; Figure S2 and Supplementary  Table S9).
Metagenomic functional pro ling was performed using Humann3 and the pathway abundances were compared using LEfSe with the LDA>2. In the unstrati ed pathway analysis, 54 pathways differed between the S. litura and S. frugiperda groups. Among them, 48 pathways, such as mycothiol biosynthesis, were more abundant in the S. frugiperda group. In contrast, six pathways included the pentose phosphate pathway were more abundant in S. litura group ( Figure 3E). The relative abundance of all pathways were summarized in Supplementary data, Table S10.
Next, a total of 407 ARGs were detected in the gut microbiomes. We detected that S. frugiperda harbored 330 resistance gene types, while S. litura harbored less, only 117 gene types ( Figure 3F). Of all the AGRs, AbaF, Trimethoprim, and MexD predominated in S. frugiperda, and Trimethoprim, vanRC, and efrA genes prevailed in S. litura. Besides, 289 genes were only found in S. frugiperda and 76 genes were unique to S. litura (Supplementary Table S11). We then matched each resistance gene type to its corresponding antibiotic, found the cephalosporin and tmacrolide antibiotic dominated in all samples.

Discussion
Spodoptera litura is an important polyphagous pest agricultural pest, which is di cult to eliminate because it has evolved great resistance to many insecticides [22]. The different development stages transcriptomes of S. litura have been sequenced and characterized [23,24]. However, those studies only focus on the cuticular protein genes. In the present study, the transcriptome analyses of S. litura at three different developmental stages were carried out and revealed a total of 101,885 DETs. The most abundant DETs (6,973) were identi ed between LL and LC group, which might due to the physiological processes, such as metamorphosis and molting, during larval to chrysalis stage. Meanwhile, we found that most of the DETs were related to basic metabolism, detoxi cation, and development. In addition, the transcriptome comparative analysis of two polyphagous species, S. litura and S. frugiperda, revealed the DETs of single-copy orthologous genes of S. frugiperda were enriched in nerves system and olfactory functions, while those of S. litura were enriched in energy metabolism.
Spodoptera litura exhibits the typical developmental characteristics at different growth stages. In the larval stage, we discovered that involved in citric acid cycle and glycolysis, related to energy metabolic processes, were upregulated (Supplementary Table S2 and S4). These active metabolic processes were conducive to larval accumulate energy for its rapid growth and future metamorphosis and reproduction. The metamorphosis of insects is a complex process accompanied with dynamically morphological and physiological changes [25]. The energy decreases sharply after pupation and remains at a low level until the adult stage [26], this could explain the larval stage we found in our study had the most abundant energy metabolism function. In the chrysalis stage, the GO terms of organ and nervous system development were enriched in the upregulated transcripts compared to larvae and adults (Supplementary Table S2 and S3). Furthermore, some transcripts involved in chitin metabolism were enhanced signi cantly in the chrysalis stage. Chitin is involved in the formation of the epidermis, midgut and the peritrophic matrix of the insect, during the development of insects, the chitin content changes dramatically, which is closely related to the growth and development [23,26]. We supposed that the enriched terms of chitin metabolism not only associated with metamorphosis, but also could bene t to construct tougher cuticles and provide chemical and biological protection in the chrysalis stage of S. litura. This helps us to develop speci c insecticides that target genes related the metabolism of chitin [27]. In the adult stage, GO terms associated with reproduction were enriched, such as chorion-containing eggshell formation (GO:0007304), embryonic pattern speci cation (GO:0009880) and structural constituent of chorion (GO:0005213), which are essential to the behaviors like mating, oviposit, and incubation. Fatty acid derivatives were documented to be involved in the synthesis of sex pheromone [28]. Here, it was also supported by GO enrichment results, where the upregulated DEGs in the adult stage was enriched in unsaturated fatty acid biosynthetic process and fatty acid biosynthetic process terms. Spodoptera litura may have active ability of xenobiotic and detoxifying metabolism in the larval and adult stages. As an important member of the metabolic enzyme superfamily, P450 plays an important role in the metabolism of xenobiotics, such as pesticides and plant secondary metabolites [29]. In addition to detoxi cation reactions, P450s also are involved in various other processes, such as insect development and reproduction [30]. Glutathione (GSH) and glutathione s-transferase (GST) are important enzymes involved in many events of pesticide degradation [31]. In our study, the expression level of transcripts associated with P450 and GST in the larvae and adult stage was signi cantly higher than in the chrysalis stage, the larval stage was the most abundant, implying that the detoxi cation ability was the strongest in the larvae stage (Supplementary Table S2, S3, and S4). Previous studies have con rmed that fthinstar S. litura larvae began to respond to the insecticide of avermectin via the P450 and glutathione-stransferase enzymes metabolism mechanism [3]. Furthermore, it was demonstrated that the increased activity of P450 was an important factor in conferring indoxacarb resistance in the fourth-instar larva S. litura [32]. The detoxi cation pathways enriched in adults may be related to their strong migration ability, which is crucial for their adaptation to different hosts and survival. In general, this enlightens us that insecticide application in the chrysalis stage with relatively low expression of transcripts related to detoxi cation pathway of S. litura may assist in inhibiting the invasion of S. litura and yield satisfactory results. In addition, the plant xenobiotic azadirachtin [33], organophosphate insecticide chlorpyrifos [34] and newer chemistry insecticides indoxacarb [35], had been reported could cause the dysfunction of central nervous system of pests, so it was considered that these insecticides could be applied in the chrysalis stage with more active neurodevelopmental function. Moreover, the upregulated DETs in larval stage compared to chrysalis and adults were enriched in many basic metabolic processes (Supplementary Table S2 and S4), such as carbohydrate metabolic process (GO:0005975) and glutathione metabolic process (GO:0006749), indicating that larval stage has active metabolism ability. It is speculated that this may be related to the unique detoxi cation ability of the S. litura larvae. Larvae is the most important developmental stage that causes damage to host plants and thus the study of the high expression genes related to detoxi cation function at larval stage is of great signi cance for pest control.
In order to compare the differences in gene expression levels between S. litura and S. frugiperda, we rst screened single-copy orthologous genes of the two species and compared their differences. Then, the functional enrichment of the differential homologous genes was carried out. The results showed that S. frugiperda had stronger nervous and olfactory functions than S. litura in both larval and chrysalis stages, which were consistent with genome level studies that showed signi cant genetic richness in sensory perceptions, such as chemical stimulation and taste, as well as in the nervous system in S. frugiperda [16]. Sensory systems play important roles in insects' feeding and identi cation of host plants [21]. This may be one of the reasons why the invasive ability of S. frugiperda is stronger than that of S. litura [16]. In the present study, we found that the upregulated differential homologous in S. litura were highly enriched in energy and metabolic processes, which was the opposite of the genome-level results, that Cui et al., identi ed some S. frugiperda-speci c positive selective genes related to energy supply [16]. One possible explanation for this result is that we concentrated only on protein-coding genes in the transcriptome analysis, indicated a high proportion of functional genes related to energy for S. litura. Among the metabolic mechanisms of resistance that have been elucidated, the most common ones include enhanced detoxi cation of insecticide by upregulated P450, GST and the increased activities of carboxylesterase [35]. By comparing the differences in aforementioned detoxi cation pathways between the two species, we found that there existed differences in the expression levels of genes related to detoxi cation pathways in different age and species groups. This provides a theoretical basis for the rational application of insecticides in practice. For example, the overexpression of P450s is the main cause of pyrethroid resistance, while the contribution of GST and esterase to the resistance is small [36]. We found that the gene expression level involved in the metabolism of xenobiotics by cytochrome P450 pathways of S. litura was higher in the larval stage than in chrysalis and adult stages, while that of S. frugiperda was highly expressed in almost the whole life stage. Therefore, better effect may be achieved by using this insecticide in the chrysalis and adult stage of S. litura than that of S. frugiperda. These illustrate that the two species may develop their own unique reply mechanisms in terms of resistance to pesticides.
In this study we also investigated metagenomes of the gut microbiome in S. litura and S. frugiperda. A comparison of the relative abundance of the bacteria from S. litura and S. frugiperda revealed different microbiota compositions between the two pests. Firmicutes, Proteobacteria and Actinobacteria were the most dominant groups in all samples. Firmicutes and Proteobacteria are crucial to maintain the growth and development of insects during the metabolism of secondary metabolites in host plants [37]. The vast majority of Acinetobacter bacteria have strong drug resistance [38]. The nding underlined the importance of gut microbiota identity and composition in the adaptation of insects. At the genus level, the relative abundance of Clostridium, Glutamicibacter, Streptomyces, Escherichia, Enterobacter were signi cantly higher in S. frugiperda when compared to S. litura. These bacteria participated in food digestion, nutrition and the detoxi cation of plant defense compounds enriched. Clostridium has been shown to play a key role in the breakdown of celluloses and insecticide metabolism [39]. Streptomyces can produce antimicrobial secondary metabolites and provide chemical defenses for insects [40]. Enterobacter, participated in the encoding of carboxylesterase and glutathione S-transferase (GST), these enzymes involved in inducing xenobiotic detoxi cation, which reinforces the result that gut microbes can play an important role in insecticide resistance of S. frugiperda [41]. Additionally, the relative abundance of Pediococcus was signi cantly higher in S. litura when compared to S. frugiperda. Pediococcus have the ability to decompose organic materials and metabolize some plant secondary biomass to produce phenolic compounds, which are necessary for host insects [37]. These results demonstrate that the two pests have a signi cant difference at the genus level, which may depend on the difference of the diet, phylogeny of host and environmental adaptability between the two pests [42].
The S. litura and S. frugiperda gut metagenome exhibited gut bacteria may be involved in various metabolic pathways, such as the carbohydrate degradation and detoxi cation. CAZy database analysis showed that the gut microbes of S. litura and S. frugiperda were rich in carbohydrate degradation-related genes. GH family enzymes, including β-galactosidase, β-glucosidases, β-glucuronidase, and β-xylosidase, which are involved in the utilization of a variety of carbon sources. Pathway analysis revealed that some of the unstrati ed pathways differed between S. litura and S. frugiperda. Pentose phosphate pathway and glycogen biosynthesis pathways were more abundant in S. litura, which indicates that S. litura might have higher production of energy within the gut compared with S. frugiperda. In addition, the monothiol biosynthesis pathway was more abundant in S. frugiperda, monothiol was proved to have antioxidant and detoxi cation effects [43]. Overall, the data suggest S. litura has more abundant energy production pathways, while the pathway related to detoxi cation was more abundant in S. frugiperda, which is consistent with the ndings in our transcriptomics and the difference in energy pathways between the two species may re ect differences in resistance, as it has been reported that insecticide resistance consumes the energy reserved of the host and reduces the energy available for other metabolic functions [44].
Spodoptera litura and S. frugiperda also differ signi cantly in terms of ARGs pro les, S. frugiperda hit more ARGs than S. litura. The most abundant ARG in all samples was DfrA42, encoding dihydrofolate reductase that confers resistant to Trimethoprim. The second most abundant ARGs was Efra, was heterodimeric ABC transporter e ux pumps that conferred resistance to macrolide antibiotic, rifamycin antibiotic and uoroquinolone antibiotic [45]. Therefore, metagenomic analysis revealed that there were different at the species or pathway level when S. litura compared with S. frugiperda. The energy metabolism related pathway was more abundant in S. litura than in S. frugiperda. On the contrary, S. frugiperda had more abundant gut bacteria associated with detoxi cation, which contributed it stronger resistance ability than S. litura. The results allowed us to comprehensively analyze the potential functions of gut bacteria contributing to food digestion, nutrition, and metabolic detoxi cation in insects.

Conclusion
In the present study, the RNA-seq libraries was constructed to investigate the transcripts expression patterns at different developmental stages of S. litura. A large number of DETs were identi ed by pairwise comparisons among different developmental stages. Some of these DETs were related to detoxi cation and development. The comparative analysis between S. litura and S. frugiperda revealed differences in their host adaptation and energy metabolism. Besides, we identi ed the structure of the intestinal microbiota of S. litura and S. frugiperda based on metagenomics. Our results demonstrate that the relative abundance of gut microbiota involved in detoxi cation was signi cance higher in S. frugiperda and the pathway related to detoxi cation was also more abundant than S. litura. This could explain the resistance of S. frugiperda is stronger than S. litura. Our study not only provides an important basis for further understanding of the molecular mechanisms of pest development but also provides insight into the development of a novel pest control strategy. In the future, more stage-speci c pesticide products can be developed through further studies on genes related to detoxi cation metabolism, which is conducive to provide new strategies for effective S. litura and S. frugiperda control.

Samples collecting
In this experiment, the samples of S. litura were obtained from the corn eld in Xindu District, Chengdu, Sichuan Province, China. The main feed source of S. litura was powdery corn. The captive environment was as follows: keep the relative humidity at 70% ±10% and the temperature at 28°C ± 1°C in the arti cial box, the ratio of light to dark time was 18:6. In total of 15 individuals were used, including three developmental stages: fth instar larvae (N= 3), sixth instar larvae (N = 3), chrysalises (N = 3), female adults (N = 3), and male adults (N = 3). We collected intestinal tissues from the larva and whole body from the chrysalis and adult individuals respectively, then immediately frozen in liquid nitrogen and stored at − 80°C. Trizol reagent (Invitrogen, Carlsbad, CA) was used to extract the total RNA of the samples. The quantity and purity of total RNA were determined by 2% agarose gel electrophoresis. The total RNA purity with OD260/280 was detected by Nanodrop. The RNA density was determined by Qubit and the RNA integrity was analyzed by Agilent 2100.

Library preparation and sequencing
After the above operation, mRNA was collected by Oligo dT enriching beads. First, fragmentation buffer was added to the mRNA to break it into small fragments. Then, the rst strand of cDNA was synthesized using random hexamer primers and under the premise of regarding mRNA fragments as a template.
Second-strand was synthesized by adding buffer, dNTPs, DNA polymerase I and RNase H. AMPure XP beads was used to purify double-stranded cDNA. After performed puri cation, terminal repair, addition of A and ligation of sequencing adaptor, libraries were produced by performing PCR ampli cation. Finally, the PCR product was puri ed by AMpure XP beads to obtain the nal library. After the library was constructed, Qubit 2.0 was used for preliminary quanti cation, and the library was diluted to1.5ng/mL. Then, Agilent 2100 was used to detect the insert size of the library. After the insert met the expectation, Q-PCR method was used to accurately quantify the effective concentration of the library to ensure the quality of the library. The Illumina Hiseq 2000 platform at Novogene (Novogene, Beijing, China) was used to generate approximately 150 bp paired end (PE) raw reads.

Transcript data processing
The NGSQC Toolkit version 2.3.3 software [46] was employed to remove adapters and low quality reads (including reads with N bases > 10% and Q-value < 20) to obtain clean reads. The reference genome data and annotation les of S. litura were downloaded from the NCBI database (https://www.ncbi.nlm.nih.gov/genome/?term=Spodoptera+litura). The clean reads of each sample were mapped to the reference genome using HISAT2 v2.1.0 software [47]. Then Samtools v1.9 [48] was used to convert the output SAM les into BAM les and sort them by chromosome positions. Finally, we obtained a transcript GTF le with the expression information using StringTie v1.3.64 [49]. The raw expression matrix was generated using the prepDE.py script provided by StringTie.

Functional annotation
Because of lacking functional annotation information in the reference genome of S. litura, GO and KEGG functional annotations of mRNA sequences were performed. Based on the position annotation from the assembled GTF le of S. litura reference genome, we extracted the mRNA sequences. To obtain a comprehensive annotation information of S. litura transcript, the above-mentioned mRNA sequences were aligned to the NCBI (ftp://ftp.ncbi.nlm.nih.gov/blast/db/) nonredundant (NR) protein databases and the Swiss-Prot database using BLASTx with E-value cutoff of 1E-5 [50]. The Gene Ontology (GO) terms information were obtained from the Swiss-Prot results. To get a background data set of the Kyoto Encyclopedia of Genes and Genomes (KEGG) annotation information, nine species (Drosophila melanogaster, Bombyx Mori, Bombyx mandarina, apilio machaon, Pieris rapae, Danaus plexippus, Helicoverpa armigera, Trichoplusia ni, and Plutella xylostella) were considered as background sets to annotate the best alignment sequences using the KEGG database, which was accomplished by the KEGG Automatic Annotation Server (KAAS) [51]. Further enrichment analysis was performed based on the GO terms and KEGG pathways obtained from the mRNA sequences using ClusterPro ler v3.16.09 [52].
Differentially expressed transcript analysis and enrichment analysis All the samples were divided into different groups according to their developmental stages, including larvae group (hereafter LL, consisting of fth instar larvae and sixth instar larvae), chrysalis group (i.e., LC), and adult group (i.e., LA, three females and three males). Then, LL versus LC, LC versus LA, and LL versus LA were analyzed using DESeq2 in R packages that takes the raw read count matrix as the input (Love et al., 2014). Transcripts that in different compared groups with |log 2-fold change|≥ 2 and p value were corrected for multiple testing with the Benjamini-Hochberg false discovery rate (FDR ≤ 0.05) were considered as DETs.
Transcriptome data including three developmental stages of S. frugiperda were collected from a previous studies in our lab [15]. There were 15 samples of S. frugiperda, including larvae group (hereafter FL, consisting of 5th instar larvae and 6th instar larvae), chrysalis group (i.e., FC), and adult group (i.e., FA, three females and three males). We took the raw read count matrix as the input for subsequent analysis.
The OrthoFinder V2.3.117 [53] was used to identify the single-copy orthologous genes between S. litura and S. frugiperda at the same age stages. Then, GETMM 8 [54] (Gene Length Corrected TMM) method was applied to correct the single-copy orthologous genes that shared by S. litura and S. frugiperda in the expression matrix. DETs of the single-copy orthologous genes analysis between LL versus FL, LC versus FC, and LA versus FA were performed by DESeq2 after the expression matrix was normalized, set the FDR ≤ 0.05 and | log 2 FC ≥ 2 as the standard to select the DETs between the two species. The functional enrichment analysis of the DETs between developmental stages and the DETs between the S. litura and the fall armyworm were performed using the GO terms and KEGG pathways background data sets constructed according to S. litura mRNA sequences. FDR ≤ 0.05 were used as a judgment standard for the signi cantly enriched in GO terms and KEGG pathways by ClusterPro ler v3.16.0 [52].

Metagenomic sequencing and pro ling
For each sample from nine fecal samples of the S. litura and S. frugiperda ( fth instar larvae and sixth instar larvae), a total amount of 0.2 μg DNA per sample was used as input material for the DNA library preparations. After the genomic DNA of the sample passed the quality test, genomic DNA sample was fragmented by sonication to a size of 350 bp. Then DNA fragments were endpolished, A-tailed, and ligated with the full-length adapter for Illumina sequencing, followed by further PCR ampli cation. After PCR products were puri ed by AMPure XP system (Beckman Coulter, Beverly, USA), DNA concentration was measured by Qubit®3.0 Flurometer (Invitrogen, USA), libraries were analyzed for size distribution by NGS3K/Caliper and quanti ed by real-time PCR (3 nM [59] was applied to build non-redundant gene sets for all predicted genes with more than 95% identity and more than 90% coverage. The gene with the longest full length from each cluster was selected as the representative read of each gene set. The

Declarations
Ethics approval and consent to participate All the sampling processes in the study were approved by the Animal Ethics Committee of the College of Life Sciences of Sichuan University (Approval ID: 20210305001). All the experimental treatment of the animals used in this study was performed strictly in accordance with the current laws of animal welfare and research in China.

Not applicable
Availability of data and materials The raw sequencing reads from this study have been submitted to the CNGBdb with the project accession CNP0001314 (available at https://db.cngb.org/search/project/CNP0001314/).