Chilo suppressalis is a widespread rice pest that poses a major threat to food security in China. This pest can develop resistance to Cry toxins from Bacillus thuringiensis (Bt), threatening the sustainable use of insect-resistant transgenic Bt rice. However, the molecular basis for the resistance mechanisms of C. suppressalis to Cry1C toxin remains unknown. This study aimed to identify genes associated with the mechanism of Cry1C resistance in C. suppressalis by comparing the midgut transcriptomic responses of resistant and susceptible C. suppressalis strains to Cry1C toxin and to provide information for insect resistance management.
A C. suppressalis midgut transcriptome of 139,206 unigenes was de novo assembled from 373 million Illumina HiSeq and Roche 454 clean reads. Comparative analysis identified 5328 significantly differentially expressed unigenes (DEGs) between C. suppressalis Cry1C-resistant and -susceptible strains. DEGs encoding Bt Cry toxin receptors, aminopeptidase-P like protein, the ABC subfamily and alkaline phosphatase were downregulated, suggesting an association with C. suppressalis Cry1C resistance. Additionally, Cry1C resistance in C. suppressalis may be related to changes in the transcription levels of enzymes involved in hydrolysis, digestive, catalytic and detoxification processes.
Our study identified genes potentially involved in Cry1C resistance in C. suppressalis by comparative transcriptome analysis. The assembled and annotated transcriptome data provide valuable genomic resources for further study of the molecular mechanisms of C. suppressalis resistance to Cry toxins.
The striped stem borer Chilo suppressalis Walker (Lepidoptera: Pyralidae) is a widespread rice pest that poses a major threat to food security throughout Asia, including China, India, Sri Lanka, Japan and Malaysia [1,2,3]. Excessive application of chemical pesticides has led to the rapid evolution of insect resistance . Consequently, insect-resistant genetically modified rice expressing Bacillus thuringiensis (Bt) has been used as an alternative method to control rice stem borers. Indeed, insect-resistant transgenic rice lines expressing Cry1A or Cry1C insecticidal proteins have been effective in controlling C. suppressalis [5,6,7,8,9,10]. In particular, the cry1C rice line has robust prospects for commercial use in China because it shows high efficiency in suppressing the damage caused by the rice pest complex C. suppressalis and Sesamia inferens. Nonetheless, the future application of Bt rice is severely threatened by the evolution of insect resistance. Therefore, the development and implementation of precautionary insect resistance management (IRM) strategies are critical for the sustainable use of Bt rice.
Understanding the mode of Cry toxin action and the mechanisms that confer resistance of insects to Cry toxins not only enhances the control efficacy of pests but also facilitates the development of IRM strategies to delay insect resistance . There are two hypotheses that have been proposed regarding the mode of Cry toxin action: the pore formation model and the signal transduction model [12,13,14]. Based on the pore formation model, protoxin is solubilized in the alkaline gut, digested by gut proteinase, and converted to the activated toxin, which then binds to the cadherin receptor located in the cell membrane of the insect midgut. The interaction between the activated toxin and cadherin receptor results in proteolytic cleavage and activated toxin oligomerization. These oligomers bind to glycosylphosphatidylinositol (GPI)-anchoring receptors, aminopeptidase N (APN) and alkaline phosphatase (ALP), among others; the oligomers are inserted into lipid raft membranes and form pores that affect the ionic balance across the cell membrane, causing cell death due to osmotic lysis [12, 13]. In contrast, the signal transduction model considers that the binding of the activated Bt toxin to specific receptors stimulates the G-protein coupled signalling pathway, leading to the activation of protein kinases . Regardless of the model, conversion of the protoxin to the activated toxin by insect midgut proteases and the binding of the activated toxins to receptors on the surface of midgut cells are recognized as essential steps for toxicity . Thus, decreased conversion of the protoxin to activated toxin and reduced toxin binding of Cry toxin to receptors are considered the most common mechanisms resulting in insect resistance to Cry toxins due to downregulation or mutation of the midgut proteinase and receptors [15, 16]. For example, decreased expression of the trypsin gene leading to incomplete activation of protoxin results in Cry1A resistance in Helicoverpa armigera and Ostrinia nubilalis [17,18,19,20,21]. Currently, more than four types of functional receptors, including cadherin , aminopeptidase N (APN) , alkaline phosphatase (ALP) , ATP-binding cassette (ABC) transporters [25, 26] and others, have been identified and verified in lepidopteran species. Additionally, reduced binding ability of Cry toxin to midgut receptors via a decrease in the activity and transcription of ALP or APN as well as mutations of APN, cadherin, and ABCC2 [14, 22, 25,26,27,28] lead to Cry1A resistance in H. armigera, Heliothis virescens, Plutella xylostella and Bombyx. mori.
In contrast to the Cry1A toxin, the Cry1C mode of action has not been described in detail, though it is also a pore-forming toxin. For lepidopteran species, APN and cadherin have been identified as Cry1C receptors in Spodoptera littoralis and Spodoptera exigua [29,30,31]. Moreover, reduced expression of SeAPN1 has been found in an S. exigua Cry1C-resistant strain . Previous studies have also shown that Cry1C receptors exhibit low or no competition with Cry1A toxins [29, 33]. Moreover, Cry1C and Cry1A toxins specifically bind to distinct isoforms of APN in the brush border membrane of Manduca sexta . In our previous study, we observed that Cry1A and Cry1C were able to recognize different binding proteins in the midgut of C. suppressalis . Furthermore, we found that cadherin from the midgut of C. suppressalis played differential functional roles with regard to Cry1Ab and Cry1C intoxication . Based on these findings, Cry1A and Cry1C toxins appear to have different modes of action in insects. Thus, understanding the molecular mechanism involved in Cry1C action is very important to better implement Cry1C alone or pyramided with Cry1A.
Relative receptors of Cry toxin, APN [36, 37], cadherin [35, 38] and ALP , have been identified in C. suppressalis, though some new receptors that may be involved in Cry toxin action have not been identified. In particular, the resistance mechanism of Cry toxin in C. suppressalis remains unknown. Furthermore, Cry1C has been shown to be particularly effective in controlling the rice pest complex C. suppressalis and S. inferens. Therefore, clarification of the molecular mechanism of Cry1C resistance in C. suppressalis appears to be particularly important.
RNA-sequencing (RNA-Seq) technology is widely applied for studying changes in gene expression in animals, plants and other organisms. In particular, Solexa/Illumina RNA-Seq and digital gene expression-based next-generation sequencing technology have been employed to identify Cry toxin resistance-related genes in Ostrinia furnacalis, H. armigera and P. xylostella [40,41,42].
In this study, we examined transcriptional differences between laboratory selected Cry1C-resistant (FZ1C) and -susceptible strains (FZS) of C. suppressalis by bioinformatics and RNA-Seq technologies. The differentially expressed transcripts were identified and further validated by quantitative real-time (qRT-PCR) analysis. The results of this study provide prospective genes and pathways contributing to resistance evolution and reveal the molecular mechanism of insect resistance to Cry toxin.
Illumina sequencing analysis and de novo assembly
The number of raw RNA sequencing reads from susceptible (FZS) and resistant (FZ1C) strains of C. suppressalis ranged from 55,728,490 to 84,319,118, respectively (Table 1). The corresponding number of clean reads ranged from 54,571,128 to 83,168,312. After the clean reads were assembled by Trinity software, 139,206 unigenes of both susceptible and resistant strains of C. suppressalis with a mean length of 1290 bp and N50 value of 2343 bp (Table 2) were obtained. Among these, 86,912 unigenes were greater than 500 bp, accounting for 62.43% of the assembled unigenes.
Annotation of assembled unigenes
Annotation of the 139,206 unigene sequences from the midgut transcriptome of the susceptible (FZS) and resistant (FZ1C) strains of C. suppressalis was performed by Blast searches with a cut-off E-value of 1e-5 in the following databases: NR (NCBI non-redundant protein sequences), NT (NCBI non-redundant nucleotide sequences), PFAM (Protein family), KOG/COG (Clusters of Orthologous Groups of proteins), Swiss-Prot (a manually annotated and reviewed protein sequence database), KO (KEGG Orthology) and GO (Gene Ontology). According to the results, 40.18% of the unigenes were annotated in NR, 27.72% in NT, 31.41% in PFAM, 15.17% in KOG, 25.9% in Swiss-Prot, 14% in KO and 31.45% in GO (Additional file 1: Figure S1,). In Blastx homology searches with a cut-off E-value of 1e-5, 33.8, 18.0, 17.3, 1.9, 1.7, and 27.5% of the unigenes matched with B. mori, P. xylostella, Danaus plexippus, Lasius. niger, Tribolium castaneum and other insects, respectively (Additional file 2: Figure S2,). The two best matches were the lepidopteran species B. mori and P. xylostella.
Gene ontology (GO) classification
GO annotation was conducted using Blast2GO software with a cut-off E-value 1e-6 to classify the unigenes into three GO classes, i.e., biological processes, cellular components and molecular functions. The unigenes were assigned into 56 functional groups (Fig. 1). Most were categorized as cellular process, metabolic process, single-organism process, biological regulation, cell, cell part, binding, and catalytic activity. Additionally, a high percentage of unigenes were assigned to the categories organelle, regulation of biological process, macromolecular complex, membrane, localization, membrane part, response to stimulus, organelle part, cellular component organization or biogenesis, signalling, transporter activity and multiorganism process (Fig. 1).
KOG is a database for the classification of orthologous genes. To further test the integrity and effectiveness of the annotation process, the unigenes annotated in the NR database were classified in KOG. In total, 21,127 unigenes were categorized into 25 KOG categories (A–W, Y and Z; Fig. 2), with “General function prediction only” accounting for the highest proportion among all the categories, followed by “Signal transduction mechanisms”, “Posttranslational modification, protein turnover, chaperones”, and “Transcription”. “Nuclear structure” and “Cell motility” occupied the lowest proportions of all categories (Fig. 2).
Functional classification by KEGG
To gain further insight into the global network underlying the experimental treatments, KEGG pathway analysis was performed, whereby 19,498 annotated unigenes were mapped to a reference pathway. A total of 230 pathways were obtained, assigned into five major groups: cellular components (3400 unigenes), environmental information processing (5266 unigenes), genetic information processing (3920 unigenes), metabolism (8947 unigenes) and organismal systems (8783 unigenes) (Fig. 3).
Differentially expressed unigenes (DEGs) among midgut transcripts
Before the analysis of differentially expressed unigenes between the resistant and susceptible strains, the hierarchical clustering analysis on the midgut samples collected for transcriptome sequencing were conducted based on the filtered and normalized count matrix using the R package PVClust v.2.0. The hierarchical clustering analysis showed that three sample replications from C. suppressalis Cry1C-resistant or Cry1C-susceptible strain grouped in an individually respective cluster (Additional file 3: Figure S3). Moreover, the two clusters grouped were completely the same (Additional file 3: Figure S3), indicating no difference for the samples between Cry1C-resistant and Cry1C-susceptible strains. These suggest that the data from C. suppressalis midgut transcriptomes are repeatable and feasible. Then, the differentially expressed unigenes (DEGs) were analysed by comparing the midgut transcriptome of the FZ1C and FZS strains. The results revealed 5328 unigenes to be significantly differentially expressed in the midgut transcriptomes of the FZ1C and FZS strains (padj< 0.05, |log2fold-change| > 1). Among these, 2908 unigenes were upregulated and 2420 unigenes downregulated. Of all DEGs, 3 unigenes were downregulated and 1 upregulated more than 10 times, 96 unigenes were downregulated and 242 upregulated between 5 and 10 times, and 2321 unigenes were downregulated and 2665 upregulated between 1 and 5 times (Fig. 4).
The DEGs encoding potential Cry toxin receptors or those involved in the Cry toxin action pathway were further investigated (Table 3), including aminopeptidase-N, the ABC transporter family, alkaline phosphatase-like, trypsin, cadherin, V-type proton ATPase catalytic subunit A and heat shock proteins. Among them, nine unigenes encode APNs, but they were all upregulated in the FZ1C strain. Two unigenes annotated as an aminopeptidase P-like protein (APP-like) were both downregulated in the FZ1C strain; moreover, one of them was downregulated 10.741-fold compared with that in the FZS strain. Seventeen unigenes were associated with the ABC transporter family, 11 of which were upregulated and 6 downregulated in the FZ1C strain compared with the FZS strain. Ten unigenes encode ALPs: seven of them were downregulated and 3 of them upregulated in the FZ1C strain. Thirteen unigenes were annotated as trypsin or trypsin/chymotrypsinogen-like protein, with 12 of them being over-transcribed and 1 under-transcribed in the FZ1C strain. Three unigenes encode cadherin. One of them was downregulated, and 2 of them were upregulated. Two unigenes encoding a V-type proton ATPase catalytic subunit A were both upregulated. Four unigenes associated with heat shock proteins were all upregulated. The heatmap of DEGs encoding potential Cry toxin receptors or those involved in pathways of the mode of Cry toxin action are displayed (Additional file 4: Figure S4).
In this study, we not only investigated DEGs encoding potential Cry toxin receptors but also examined DEGs encoding detoxification enzymes, such as cytochrome P450 (P450), carboxylesterase (CaE) and glutathione S-transferase (GST) (Table 3). There were 15 unigenes encoding CaE; 12 of them were upregulated, and 3 of them were downregulated. Four unigenes encoding P450 were upregulated, and the other 2 unigenes were downregulated. The 2 DEGs encoding GST were both downregulated. Moreover, some putative glutamate receptors targeting insecticides were identified.
Immune genes involved in haemolymph melanization, such as serpin protease and serpin protease inhibitor, were also examined (Table 3). The results showed 17 unigenes encoding serpin proteases. Seven of them were upregulated and 10 downregulated in the resistant strain compared with the susceptible strain. Twenty-two unigenes encode serpin protease inhibitors: two of them were upregulated and the other 20 downregulated in the Cry1C-resistant strain compared to the susceptible strain (Table 3). A heatmap of DEGs encoding serpin proteases and serpin protease inhibitors is exhibited (Additional file 4: Figure S4, F and G).
To further analyse the biological functions of the significant DEGs between the resistant (FZ1C) and susceptible (FZS) strains of C. suppressalis, GO functional and pathway enrichment analyses were carried out. The upregulated unigenes were annotated to 1727, 378 and 761 GO terms in the biological process, cellular component and molecular function categories, respectively. The top 30 GO terms (according to corrected P-values) are displayed in Fig. 5, including 7 GO terms for biological processes, 4 GO terms for cellular components and 19 GO terms for molecular functions. Among the 7 biological processes, the amide biosynthetic process (GO: 0043604, corrected P-value = 3.87e− 07) was the most enriched GO term. The cellular amide metabolic process (GO: 0043603, corrected P-value = 1.89e− 06) was the largest, with 132 upregulated unigenes. Among the 4 cellular component GO terms, the proteasome core complex term (GO: 0005839, corrected P-value = 1.61e− 06) was the most represented, and the ribosome (GO: 0005840, corrected P-value = 4.89e− 06) was the largest category, with 88 upregulated unigenes. For the significantly enriched GO terms in “molecular function”, small molecule binding (GO: 0036094, corrected P-value = 1.01e− 10) was the most enriched and largest category, with 296 upregulated unigenes.
The downregulated unigenes were annotated to 1540, 331 and 603 GO terms for biological process, cellular component and molecular function, respectively. Among the top 30 GO terms (according to corrected P-value) 26 terms were for biological process, 1 for cellular component and 3 for molecular function (Fig. 6). With regard to the 26 biological process GO terms, inositol catabolic process was the most enriched category (GO: 0019310, corrected P-value = 1.35e− 05), and the oxidation-reduction process (GO: 0055114, corrected P-value = 7.24e− 05) was the largest, which included 129 downregulated unigenes. In terms of cellular components, the proteinaceous extracellular matrix was the most enriched term (GO: 0005578, corrected P-value = 0.017839), with 15 downregulated unigenes. Among the 3 molecular function GO terms, inositol oxygenase activity (GO: 0050113, corrected P-value = 1.35e− 05) was the most represented, and catalytic activity (GO: 0003824, corrected P-value = 0.0053335) was the largest, with 565 downregulated unigenes.
Validation of differentially expressed genes by qRT-PCR
To validate the expression patterns of the unigenes, 12 differentially expressed unigenes (including Cluster-21,897.60875, 21,897.37634, 25,303.0, 21,897.94254, 21,897.2225, 21,897.56483, 21,897.57101, 21,897.55566, 21,897.52385, 21,897.54886, 21,897.24953 and 21,897.56465) were selected for qRT-PCR according to their fold change values in the expression profiles from the dataset. The elongation factor was set as the candidate reference gene for RT-qPCR normalization. The results showed that the expression patterns presented by the comparative transcriptome were consistent with the qRT-PCR results (Fig. 7). Overall, relative expression of the top 8 unigenes was significantly upregulated in the resistant FZ1C strain of C. suppressalis, whereas that of the last 4 unigenes was remarkably reduced in the resistant strain (Fig. 7).
The mode of Cry toxin action is a very complex process involving the toxin structure, many enzymes and receptors, among others, and a change in any step of the toxicology process inevitably leads to insect resistance. In this study, comparative transcriptome analysis was conducted between susceptible and resistant strains of C. suppressalis. A total of 139,206 unigenes were de novo assembled from 373 million clean reads in the C. suppressalis transcriptome. Among them, 5328 significantly DEGs were identified between C. suppressalis Cry1C-resistant and Cry1C-susceptible strains. Among these DEGs, unigenes encoding Bt Cry toxin receptors, aminopeptidase-P like protein, the ABC subfamily and alkaline phosphatase were downregulated and suggested to be associated with C. suppressalis Cry1C resistance. Additionally, the transcription level of unigenes encoding enzymes involved in hydrolysis, digestive, catalytic and detoxification processes were significant increased, which appeared to be related to Cry1C resistance in C. suppressalis. Presumably, multiple genes and different pathways are involved in Cry1C toxin resistance.
Among these DEGs, unigenes associated with the Cry toxin receptor were identified, such as aminopeptidase N, ABC transporter family members, alkaline phosphatase, cadherin, trypsin and V-type proton ATPase catalytic subunit A. Many previous studies on glycosylphophatidylinositol (GPI)-anchored APN1 from several insects have consistently demonstrated that APN1 is an important midgut receptor for Cry toxins and is associated with insect resistance [25, 43,44,45,46]. Our previous study found that knockdown of aminopeptidase-N isoform (APNs) results in decreased susceptibility of C. suppressalis larvae to Cry1A toxins, indicating that APNs are important functional receptors of Cry1A toxins in C. suppressalis . Qiu et al.  demonstrated that downregulation of APNs by RNAi reduced the susceptibility of C. suppressalis larvae to Cry1C and Cry1Ab/Cry1Ac fusion proteins expressed by transgenic Bt rice, suggesting that APNs are involved in the mode of Cry1Ab/Cry1Ac and Cry1C action. In this study, five unigenes were annotated to APN1 (ABC69855.3), which included Cluster-21,897.55689, Cluster-21,897.55577, Cluster-21,897.55566, Cluster-21,897.49282 and Cluster-21,897.54390. One unigene annotated to APN2 (AGG36452.1) was Cluster-21,897.57790. Three unigenes were annotated to APN5 (AFU51581.1), which contained Cluster-21,897.54666, Cluster-21,897.56678 and Cluster-21,897.57101. Nine unigenes encoding APNs in the current study were all upregulated in the resistant strain (FZ1C) compared with the susceptible strain (FZS) of C. suppressalis. These results are consistent with those of previous studies [41, 47, 48]. Presumably, these upregulated APNs may compensate for the fitness costs of resistance . In addition, we observed significant downregulation of two unigenes annotated as aminopeptidase P-like proteins in the resistant strain (FZ1C) of C. suppressalis; in particular, expression of one unigene was decreased 10.74-fold in the FZ1C strain compared with the FZS strain. Khajuria et al.  reported that an aminopeptidase P-like protein is possibly involved in Cry resistance in O. nubilalis. These results indicate that an aminopeptidase P-like protein is likely involved in the mode of Cry1C toxin action and accounts for C. suppressalis resistance to Cry1C.
Another important group of receptors involved in the mode of Cry toxin action is ABC transporters. Many studies have proven that the resistance of three lepidopteran species, H. virescens, P. xylostella and Trichoplusia ni, to Cry1Ac is linked to ABC transporters [27, 50]. ABC transporters, such as ABCG1, ABCC2, and ABCC3, are reported to be Cry toxin receptors, and their mutation or downregulation is proven to be related to insect resistance [25, 51, 52]. In P. xylostella, mutations in ABCC2 disrupt the binding of Cry1Ac to membrane vesicles and lead to Cry toxin resistance. Transcriptome analysis of P. xylostella showed that eight unigenes annotated as ABCC2 were detected in the Cry1Ac-resistant strain, with the majority being downregulated . In this study, we identified 17 unigenes encoding ABC transporter family members, including ABCA1, ABCB1, ABCB7, ABCC1, ABCC3, ABCC4, ABCD2, ABCF2, ABCG1 and ABCG4, but no differentially expressed unigenes were annotated as ABCC2. Furthermore, six downregulated unigenes annotated as ABC transporters (Cluster-21,897.52385, Cluster-21,897.87849, Cluster-21,897.77904, Cluster-21,897.854, Cluster-21,897.54886 and Cluster-21,897.58674) were identified in the C. suppressalis resistant (FZ1C) strain and were considered to be associated with the resistance of C. suppressalis to Cry1C toxin. Nevertheless, further study should be conducted to reveal the exact function of ABC transporters in the resistant strain of C. suppressalis.
Previous studies have also demonstrated that ALP is an important receptor that interacts with Cry toxins [53,54,55,56]. In this study, we identified 10 unigenes matching ALP1, ALP3 and ALP-like that were significantly differentially expressed between Cry1C-resistant and Cry1C-susceptible strains of C. suppressalis. Among them, 7 unigenes (Cluster-21,897.109470, Cluster-21,897.58470, Cluster-21,897.75922, Cluster-21,897.56465, Cluster-21,897.58789, Cluster-21,897.40486 and Cluster-21,897.24953) were downregulated. We speculate that these unigenes with decreased transcript levels may account for C. suppressalis resistance to Cry1C toxin, in accordance with a previous study .
Three unigenes were annotated as cadherin (Cluster-21,897.53425, Cluster-21,897.66616 and Cluster-21,897.72833), but they were all upregulated in the resistant strain. This result is consistent with our previous study, which found that cadherin did not play a major role in the mode of Cry 1C action in C. suppressalis .
The conversion of the Cry protoxin to the cytotoxic form by proteases is an inevitable process for Cry toxicity. However, overexpression of midgut proteinases may enhance the possibility of degrading activated Cry toxin and cause a decrease in toxicity . In this study, 12 unigenes annotated as trypsin or trypsin/chymotrypsinogen-like proteinase were overexpressed in the Cry1C-resistant strain compared with the susceptible strain. Similar results were also reported in O. furnacalis Cry1Ab-resistant strains  and Bt-resistant strains of mosquitoes . Over-transcription of these proteinases in the resistant strain of C. suppressalis possibly reflects higher activity of midgut proteinases for degrading toxins. Additionally, some unigenes encoding serpin protease and serpin protease inhibitors were found to be associated with immune function for the regulation of haemolymph melanization. In addition, unigenes matching the detoxification enzymes cytochrome P450 (P450), carboxylesterase (CaE) and glutathione S-transferase (GST) were identified. In particular, P450 family members play an important role in degrading chemical insecticides. They have also been observed to respond to Cry toxin in lepidopterous insects [40, 57]. In the midgut transcriptome of O. nubilalis larvae , six unigenes encoding P450 were found to be downregulated and eight unigenes upregulated after exposure to Cry toxins. However, in O. furnacalis , eight unigenes matching P450 were all downregulated in a Cry toxin-resistant strain. Thus, the exact role of P450 in the mode of Cry toxin action is uncertain, as they are generally thought to be involved in the degradation of xenobiotics . In the present study, a large proportion of unigenes (12/21) annotated as detoxification enzymes were upregulated in the Cry1C-resistant strain. These enzymes for detoxification may act as a general response to environmental stress (such as Cry toxin) or may help insects avoid the damage of Cry toxins and compensate for fitness costs [40, 57].
GO and KEGG analyses provided important clues to reveal the potential mechanisms involved in the development of resistance. In our study, a large number of unigenes showed significant enrichment in digestive enzymes, hydrolase enzymes, detoxification enzymes, and catalytic enzyme-related GO categories (Additional file 5). Previous studies reported that these enzymes are considered to be the most important factors in reducing fitness costs [17, 19, 41, 59]. These findings suggest that C. suppressalis resistance to Cry1C toxin may be associated with increased digestive activity, catalytic activity, detoxification activity and hydrolase activity in the larval midgut. Moreover, KEGG category analyses revealed significantly enrichment of a series of unigenes in xenobiotic stimulation and immunity-related KEGG pathways (Additional file 6). Similar findings were also found for O. furnacalis [40, 57] and H. armigera . These results also indicated that Cry1C resistance may be associated with the immune response of C. suppressalis.
In conclusion, this study is the first to report the mechanism of C. suppressalis resistance to Cry toxin based on genetic information from the sequenced transcriptome, revealing a large number of unigenes with greatly enriched sequence information for C. suppressalis. Based on our findings, several factors are related to Cry1C toxin resistance in this pest. The candidate receptors aminopeptidase-P-like protein, the ABC subfamily and alkaline phosphatase, which were downregulated, might be associated with Cry1C resistance in C. suppressalis. Additionally, changes in enzymes related to digestive activity, catalytic activity, detoxification activity, hydrolase activity and related genes may account for C. suppressalis resistance to Cry1C toxin. We also identified multiple pathways that might be involved in Cry1C resistance in C. suppressalis. Therefore, we speculate that the mechanism of C. suppressalis resistance to Cry1C toxin is a complicated process involving a series of genetic and metabolic factors. With these important genetic resources, further studies should be conducted to validate the gene functions associated with Bt resistance in C. suppressalis using RNA interference or genome-editing technology.
Insect rearing and selection of the resistant strain
The Fuzhou susceptible strain (FZS) of C. suppressalis was initially collected from paddy fields in Fuzhou (26.1°N, 119.3°E), Fujian Province, China, in 2012. The insect colony was continually reared on an artificial diet as described previously  without exposure to any Cry toxins. Based on the FZS strain, the Cry1C-resistant strain (FZ1C) was obtained from the fifth generations colony of the susceptible strain selected with an artificial diet with trypsin-activated Cry1C toxin (Meiyan (Beijing) Agricultural Technology Co., Ltd.). The Cry1C-selected strain (FZ1C) was initially exposed to a sublethal dose of the Cry1C diet throughout larval development. The toxin concentration was steadily increased in succeeding generations to achieve 40–60% mortality of the exposed larvae. After 41 generations of selection, the strain (FZ1C) developed 42.6-fold resistance . In this study, the Cry1C-resistant strain (FZ1C) reared on an artificial diet containing 80 μg/mL activated Cry1C toxin, which had been selected for 57 generations with 120-fold resistance, was used to detect C. suppressalis Cry1C resistance-related genes. In parallel, the susceptible strain (FZS) reared in the absence of Cry1C toxin was used as the negative control.
All cultures were maintained under constant conditions (27 ± 1 °C, 70–80% RH, and a 16:8-h light/dark photoperiod).
Dissection of the midgut and extraction of RNA
Before dissection, the FZS strain was reared on the diet without any Cry toxins or other chemical insecticides, and the FZ1C strain selected for 57 generations was reared on the artificial diet with 80 μg/mL activated Cry1C toxin. Larval midguts of both strains of C. suppressalis were dissected when the larvae reached the 4th instar. The respective midgut tissues were washed with 0.7% NaCl (w/v) to remove debris and haemolymph, and the samples were stored in RNA holding solution (TransGen Biotech Co. Ltd., Beijing, China). Thirty individual midguts from 4th-instar larvae of FZS and FZ1C were collected in an RNase-free 1.5 mL tube as one biological replicate and each strain contained six biological replicates. Three replicates were used to Illumina RNA-Seq and the gene expression profile analysis, the others were used for RT-qPCR analysis to validate the results of differentially expressed genes. TRIzol reagent was used to extract total RNA from each sample according to the manufacturer’s instructions, in which RNase-free DNase I was used to treat all samples to remove genomic DNA contamination. The RNA was suspended in 100 μL DEPC-treated water. After assessment of the quantity and purity of the total RNA, the RNA samples were delivered to Beijing Novogene Technology Company for RNA sequencing.
RNA-Seq library preparation and Illumina sequencing
The RNA-Seq libraries were generated from a total amount of 1.5 μg RNA per biological replicate using NEBNext® Ultra™ RNA Library Prep Kit for Illumina® (NEB, USA) according to the manufacturer’s instructions, and index codes were added to identify the sequences of each sample. Briefly, mRNA was purified from total RNA using poly-T oligo-attached magnetic beads and fragmented by using divalent cations under elevated temperature. Synthesis of first-strand cDNA was performed using random hexamer primers and RNase H. Second-strand cDNA synthesis was subsequently performed using DNA polymerase I and RNase H. The remaining overhangs of the cDNA was converting into blunt ends by exonuclease/polymerase activities of DNA polymerase I, the 3′ ends of the DNA fragments were adenylated and NEBNext adaptors with hairpin loop structures were ligated. In order to obtain cDNA fragments of 250 ~ 300 bp in length, the cDNA library fragments were purified with the AMPure XP system (Beckman Coulter, Beverly, USA). Subsequently, the cDNA fragments were treated with adding 3 μL USER Enzyme (NEB, USA) and incubated at 37 °C for 15 min followed by 5 min at 95 °C before Polymerase Chain Reaction (PCR) reaction. The cDNA fragments were amplificated by PCR using Phusion High-Fidelity DNA polymerase, universal PCR primers and index Primer. After purified by AMPure XP system, the library quality was evaluated using Agilent Bioanalyzer 2100 system.
The cBot Cluster Generation System was applied to cluster the index-coded samples using a TruSeq PE Cluster Kit v3-cBot-HS (Illumina) according to the manufacturer’s protocols. After cluster generation, the library preparations were sequenced using the Illumina HiSeq 4000 platform to generate 150 bp paired-end raw reads. The clean reads were obtained by removing reads containing adapters, reads containing poly-N and low-quality reads from the raw reads by first processing through in-house Perl scripts. At the same time, the Q20, Q30, GC content and sequence duplication level of the clean reads were calculated in this step. All of the downstream analyses were conducted with clean reads of high quality.
Identification of differentially expressed genes (DEGs)
Gene expression levels were estimated by RSEM . The DESeq R package (1.10.1) was applied to analyse differentially expressed genes (DEGs) between the susceptible and resistant strains basing on the negative binomial distribution model. The resulting P values were adjusted using Benjamini and Hochberg’s approach for controlling the false discovery rate. Genes with an threshold of |log2(fold- change)| and adjusted P value (P < 0.05) found by DESeq were considered significantly differentially expressed.
Validation of differentially expressed genes by real-time quantitative PCR
Transcriptome results were verified using quantitative real-time PCR (qRT-PCR). Total RNA was extracted from each sample using TRIzol reagent according to the product instructions, with RNase-free DNase I to digest all samples to avoid genomic DNA contamination. The RNA was suspended in 100 μL DEPC-treated water. First-strand cDNA was synthesized by a Fast Quant RT Kit (With gDNase, KR106) (TIANGEN Biotech Co. Ltd., Beijing, China). Twelve genes (Additional file 7) identified as differentially expressed were selected based on their fold changes in the expression profile. Primer pairs (Additional file 8) were designed using Primer3plus (http://www.bioinformatics.nl/cgi-bin/primer3plus/primer3plus.cgi). qPCR was performed using the ABI 7500 system (ABI, USA) with a reaction volume of 20 μL containing 1 μL of 1:5 diluted cDNA in RNase-Free ddH2O, 10 μL 2 × SuperReal PreMix Plus (SYBR Green, FP205, TIANGEN Biotech Co. Ltd., Beijing, China), 0.6 μL of 10 μM each primer, 0.4 μL of 50× ROX reference dye, and 7.4 μL of RNase-Free ddH2O. The qPCR conditions were 95 °C for 15 min, followed by 40 cycles of 95 °C for 10 s for denaturation and 60 °C for 32 s for annealing and extension. The melting curve was produced at 95 °C for 15 s, 60 °C for 1 min, 95 °C for 15 s, and 60 °C for 15 s. Elongation factor 1 (EF-1) was used as an internal reference gene. Each of the three biological replicates was measured with four technical replicates, and the expression levels of candidate genes were normalized with the EF-1 gene.
Availability of data and materials
The raw reads of the two transcriptomes in this study have been deposited in the NCBI SRA database under accession numbers SRR12170219 (midgut transcriptome of susceptible strain of Chilo suppressalis) and SRR12170218 (midgut transcriptome of Cry1C resistant strain of Chilo suppressalis).
Differentially expressed unigenes (DEGs)
Insect resistance management
Aminopeptidase P-like protein
Cytochrome P450 monooxygenase
Han LZ, Li SB, Liu PL, Peng YF, Hou ML. New artificial diet for continuous rearing of Chilo suppressalis (Lepidoptera: Crambidae). Ann Entomol Soc Am. 2012;105(2):253–8.
Ma WH, Zhao XX, Yin CL, Fan J, Du XY, Chen TY, et al. A chromosome-level genome assembly reveals the genetic basis of cold tolerance in a notorious rice insect pest, Chilo suppressalis. Mol Ecol Resour. 2019;00:1–15.
Deka S, Barthakur S. Overview on current status of biotechnological interventions on yellow stem borer, Scirpophaga incertulas (Lepidoptera: Crambidae) resistance in rice. Biotechnol Adv. 2010;28:70–81.
Gao YL, Fu Q, Wang F, Lai FX, Luo J, Peng YF, et al. Effects of transgenic rice harboring cry1Ac and CpTI genes on survival of Chilo suppressalis and Sesamia inferens and field composition of rice stem borers. Chin J Rice Sci. 2006;20:543–8.
Han LZ, Hou ML, Wu KM, Peng YF, Wang F. Lethal and sub-lethal effects of transgenic rice containing cry1Ac and CpTI genes on the pink stem borer, Sesamia inferens (Walker). Sci Agri Sin. 2011;10:384–93.
Han LZ, Han C, Liu ZW, Chen FJ, Jurat-Fuentes JL, Hou ML, et al. Binding site concentration explains the differential susceptibility of Chilo suppressalis and Sesamia inferens to Cry1A-producing rice. Appl Environ Microbiol. 2014;80:5134–40.
Li ZY, Sui H, Xu YB, Han LZ, Chen FJ. Effects of insect-resistant transgenic Bt rice with a fused cry1Ab+cry1Ac gene on population dynamic of the stem borers, Chilo suppressalis and Sesamia inferens, occurring in paddy field. Acta Ecol Sin. 2012;32:1783–9.
Wang YN, Zhang L, Li YH, Liu YM, Han LZ, Zhu Z, et al. Expression of Cry1Ab protein in a marker-free transgenic Bt rice line and its efficacy in controlling a target pest, Chilo suppressalis (Lepidoptera: Crambidae). Environ Entomol. 2014;43:528–36.
Zhang XB, Candas M, Griko NB, Taussig R, Bulla LA. A mechanism of cell death involving an adenylyl cyclase/PKA signaling pathway is induced by the Cry1Ab toxin of Bacillus thuringiensis. Proc Natl Acad Sci U S A. 2006;103(26):9897–902.
Rajagopal R, Arora N, Sivakumar S, Rao NGV, Nimbalkar SA, Bhatnagar RK. Resistance of Helicoverpa armigera to Cry1Ac toxin from Bacillus thuringiensis is due to improper processing of the protoxin. Biochem J. 2009;419:309–16.
Atsumi S, Miyamoto K, Yamamoto K, Narukawa J, Kawai S, Sezutsu H, et al. Single amino acid mutation in an ATP-binding cassette transporter gene causes resistance to Bt toxin Cry1Ab in the silkworm, Bombyx mori. Proc Natl Acad Sci U S A. 2012;109:E1591–8.
Jurat-Fuentes JL, Karumbaiah L, Jakka SRK, Ning C, Liu C, Wu KM, et al. Reduced levels of membrane-bound alkaline phosphatase are common to lepidopteran strains resistant to cry toxins from Bacillus thuringiensis. PLoS One. 2011;6:e17606.
Ren XL, Chen RR, Zhang Y, et al. A Spodoptera exigua cadherin serves as a putative receptor for Bacillus thuringiensis Cry1Ca toxin and shows differential enhancement of Cry1Ca and Cry1Ac toxicity. Appl Environ Microb. 2013;79(18):5576–83.
Wang XY, Du LX, Liu CX, Gong L, Han LZ, Peng YF. RNAi in the striped stem borer, Chilo suppressalis, establishes a functional role for aminopeptidase N in Cry1Ab intoxication. J Invertebr Pathol. 2017;143:1–10.
Qiu L, Fan JX, Zhang BY, Liu L, Wang XP, Lei CL, et al. RNA interference knockdown of aminopeptidase N genes decrease the susceptibility of Chilo suppressalis larvae to Cry1Ab/Cry1Ac and Cry1C-expressing transgenic rice. J Invertebr Pathol. 2017;145:9–12.
Qiu L, Wang P, Wu T, Li B, Wang X, Lei C. Down regulation of Chilo suppressalis alkaline phosphatase genes associated with resistance to three transgenic Bacillus thuringiensis rice lines. Insect Mol Biol. 2018;27(1):83–9.
Zhang SP, Cheng HM, Gao YL, Wang GR, Liang GM, Wu KM. Mutation of an aminopeptidase N gene is associated with Helicoverpa armigera resistance to Bacillus thuringiensis Cry1Ac toxin. Insect Biochem Mol Biol. 2009;39:421–9.
Flores-Escobar B, Rodríguez-Magadan H, Bravo A, Soberón M, Gómez I. Differential role of Manduca sexta aminopeptidase-N and alkaline phosphatase in the mode of action of Cry1Aa, Cry1Ab, and Cry1Ac toxins from Bacillus thuringiensis. Appl Environ Microbiol. 2013;79:4543.
Després L, Stalinski R, Tetreau G, Paris M, Bonin A, Navratil V, et al. Gene expression patterns and sequence polymorphisms associated with mosquito resistance to Bacillus thuringiensis israelensis toxins. BMC Genomics. 2014;15:926.
Khajuria C, Buschman LL, Chen MS, Siegfried BD, Zhu KY. Identification of a novel Aminopeptidase P-like gene (OnAPP) possibly involved in Bt toxicity and resistance in a major corn Pest (Ostrinia nubilalis). PLoS One. 2011;6(8):e23983.
Guo ZJ, Kang S, Zhu X, Xia JX, Wu QJ, Wang SL, et al. Down-regulation of a novel ABC transporter gene (Pxwhite) is associated with Cry1Ac resistance in the diamondback moth, Plutella xylostella (L.). Insect Biochem Mol Biol. 2015;59:30–40.
Tanaka S, Endo H, Adegawa S, Iizuka A, Imamura K, Kikuta S, et al. Bombyx mori abc transporter c2 structures responsible for the receptor function of Bacillus thuringiensis Cry1Aa toxin. Insect Biochem Mol Biol. 2017;91:44–54.
Moonsom S, Chaisri U, Kasinrerk W, Angsuthanasombat C. Binding characteristic to mosquito-larval midgut proteins of the cloned domain II-III fragment from the Bacillus thuringiensis Cry4Ba toxin. J Biochem Mol Biol. 2007;40:783–90.
Vellichirammal NN, Wang HC, Eyun S, Moriyama EN, Coates BS, Miller NJ, et al. Transcriptional analysis of susceptible and resistant European corn borer strains and their response to Cry1F protoxin. BMC Genomics. 2015;16:558.
Cao GC, Zhang LL, Liang GM, Li XC, Wu KM. Involvement of nonbinding site proteinases in the development of resistance of Helicoverpa armigera (Lepidoptera: Noctuidae) to Cry1Ac. J Econ Entomol. 2013;106:2514–21.
We thank the Beijing Novogene Technology Company for the assistance in original data processing and related bioinformatics analysis. This research was funded by the National Genetically Modified Organisms Key Breeding Projects of China (2016ZX08011-001, 2018ZX0801101B and 2016ZX08012-005), the National Natural Science Foundation of China (31572336) and the Scientific and Technological Innovation Project of the Chinese Academy of Agricultural Sciences.
This research was funded by the National Genetically Modified Organisms Key Breeding Projects of China (2016ZX08011–001, 2018ZX0801101B and 2016ZX08012–005), the National Natural Science Foundation of China (31572336) and the Scientific and Technological Innovation Project of the Chinese Academy of Agricultural Sciences.
Authors and Affiliations
Department of Entomology, College of Plant Protection, Nanjing Agricultural University, Nanjing, 210095, China
Geng Chen, Yanhui Wang, Yanmin Liu & Fajun Chen
The State Key Laboratory for Biology of Plant Diseases and Insect Pests, Institute of Plant Protection, Chinese Academy of Agricultural Sciences, Beijing, 100193, China
GC carried out the experiments and drafted the manuscript. YHW and YML performed the experiments and helped with the data analysis. LZH and FJC conceived of the study, participated in its design and coordination and helped to revise the manuscript. All authors read and approved the final manuscript.
We claimed that the following statement was accurate. Chilo suppressalis is a common agricultural pest, which is not included in the “List of Endangered and Protected Animals in China”. The original collecting site is not private and the field collection of Chilo suppressalis was permitted. All animal procedure in this study was performed according to guidelines developed by the ethics committee of the State Key Laboratory for Biology of Plant Diseases and Insect Pests, Institute of Plant Protection, Chinese Academy of Agricultural Sciences and the Experimental Animal Welfare Ethics Committee of Nanjing Agricultural University. In order to make all efforts to minimize suffering, all animal procedure was performed under anesthesia, and wounds were cleaned that before they got infected.
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Hierarchical cluster analysis of C. suppressalis midgut samples from Cry1C-resistant (FZ1C) and Cry1C-susceptible (FZS) strains. FZ1Ca_1, FZ1Ca_2, FZ1Ca_3: three sample replications from the FZ1C strain; FZS_1, FZS_2, FZS_3: three sample replications from the FZS strain.
Primer sequences for qRT-PCR of twelve selected unigenes annotated as potential Cry toxin receptors.
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.
Chen, G., Wang, Y., Liu, Y. et al. Differences in midgut transcriptomes between resistant and susceptible strains of Chilo suppressalis to Cry1C toxin.
BMC Genomics21, 634 (2020). https://doi.org/10.1186/s12864-020-07051-6