- Research Article
- Open Access
De novo transcriptome and expression profile analyses of the Asian corn borer (Ostrinia furnacalis) reveals relevant flubendiamide response genes
BMC Genomics volume 18, Article number: 20 (2017)
The Asian corn borer (ACB), Ostrinia furnacalis (Guenée), has become the most damaging insect pest of corn in Asia. However, the lack of genome or transcriptome information heavily hinders our further understanding of ACB in every aspect at a molecular level and on a genome-wide scale. Here, we used the Ion Torrent Personal Genome Machine (PGM) Sequencer to explore the ACB transcriptome and to identify relevant genes in response to flubendiamide, showing high selective activity against ACB.
We obtained 35,430 unigenes, with an average length of 716 bp, representing a dramatic expansion of existing cDNA sequences available for ACB. These sequences were annotated with Non-redundant Protein (Nr), Gene Ontology (GO), Clusters of Orthologous Groups (COG) and Kyoto Encyclopedia of Genes and Genomes (KEGG) to better understand their functions. A total of 31 cytochrome P450 monooxygenases (P450s), 27 carboxyl/cholinesterases (CCEs) and 19 glutathione S-transferases (GSTs) were manually curated to construct phylogenetic trees, and 25 unigenes encoding target proteins (acetylcholinesterase, nicotinic acetylcholine receptor, gamma-aminobutyric acid receptor, glutamate-gated chloride channel, voltage-gated sodium channel and ryanodine receptor) were identified. In addition, we compared and validated the differentially expressed unigenes upon flubendiamide treatment, revealing that the genes for detoxification enzymes (P450s and esterase), calcium signaling pathways and muscle control pathways (twitchin and tropomyosin), immunoglobulin (hemolin), chemosensory protein and heat shock protein 70 were significantly overexpressed in response to flubendiamide, while the genes for cuticular protein, protease and oxidoreductase showed much lower expression levels.
The obtained transcriptome information provides large genomic resources available for further studies of ACB. The differentially expressed gene data will elucidate the molecular mechanisms of ACB in response to the novel diamide insecticide, flubendiamide. In particular, these findings will facilitate the identification of the genes involved in insecticide resistance and the development of new compounds to control the ACB.
The Asian corn borer (ACB), Ostrinia furnacalis (Guenée) (Lepidoptera: Crambidae) is the most destructive lepidopteran pest of corn in Asia, particularly in China and the Philippines. ACB larvae damage corn by attacking the new leaves, ears and cobs; these pests also bore into the stalks and cobs of corn . The damage manifested by ACB infestation constitutes a major constraint on agriculture in corn-producing countries worldwide. Yield losses caused by the ACB are estimated at 10 to 20% in an ordinary year but may be more than 30% or may even result in no harvest in an outbreak year . Moreover, the corn borer is considered an aggravating factor for the epidemiology of Fusarium ear rot in maize, and insect damage to ears can increase the fumonisin contamination of kernels [3, 4]. Novel insecticides with high selective activity are preferred to control this insect pest. Recently, two classes of diamide chemicals have emerged as novel insecticides targeting insect RyRs: the phthalic diamides, including flubendiamide [5, 6], and the anthranilic diamides, including chlorantraniliprole and cyantraniliprole . These diamides have potent species-specific insecticidal activity against a range of lepidopteran pest species and can elicit intracellular Ca2+ release from the sarcoplasmic reticulum or endoplasmic reticulum . Yang et al. reported that chlorantraniliprole provided excellent control efficacy against second-generation ACB and was relatively safe to the spider Xysticus ephippiatus. In addition, we previously demonstrated that flubendiamide also showed high efficacy against the ACB. Consequently, these diamide insecticides are suitable for integrated pest management programs for the regulation of the ACB . However, the molecular mechanisms of the ACB in response to diamide insecticides are not clearly understood.
Despite the insect’s importance in agricultural production, the sequence information available for the ACB lies in contrast to the economic importance of this pest. The current understanding of this pest has largely been hindered by the lack of thorough genetic information. Until recently, no complete picture has been achieved, even of specific gene families. For example, cytochrome P450s (CYPs) have been highly recognized as a supergene family, with 36 to 180 genes in insect genomes [9, 10], while only one CYP gene with a full-length coding region was accessed for the ACB in GenBank. This scarcity of genetic information for the ACB has also resulted in a paucity of genetic studies and integrated theories for understanding the basic biology of this pest.
The aim of the present study was to use PGM sequencing to produce a de novo transcriptome of the ACB as a resource for current and future studies of this pest species. This resource was subsequently used as a reference to provide insights into insecticide-related changes in gene expression in flubendiamide-treated ACB and to identify a broad range of genes encoding the target proteins and detoxification enzymes of insecticides currently in use for ACB control.
RNA-Seq and sequence assembly
When two ACB libraries (Control strain feed on a semi-artificial diet without insecticide treatment and flubendiamide-treated strain) were sequenced using the PGM platform, the output was 549.4 and 472.5 Mb. After trimming the adaptor sequences and removing the low-quality reads, 1,625,223 and 1,732,688 clean reads were generated from the control and flubendiamide-treated libraries, with average sequence lengths of 279 and 265 bp, respectively. The clean reads of the two libraries (control and flubendiamide-treated) were assembled into 31,059 and 31,478 unigene sequences, with average unigene lengths of 759 and 756 bp, respectively. A summary of the PGM sequencing and assembly is presented in Table 1, and the length distributions for the clean reads and unigenes are presented in Fig. 1.
Analysis of the transcripts from the Asian corn borer
BLASTx and BLASTn were used to compare each ACB transcript with a cut-off E-value of 1.0E−5 against GenBank entries. A total of 19,422 transcripts (55%) had BLASTx hits in Non-redundant Protein (Nr) databases, and 11,004 transcripts (31%) had BLASTn hits in Non-redundant Nucleotide Sequence (Nt) databases. Some of the ACB transcripts were homologous to those from more than one species, but in general, most ACB transcripts were homologous to Lepidoptera species, accounting for 4392 hits among the 11,004 BLASTn hits, including 1600 hits (14.5%) to Bombyx mori entries, 1363 hits (12.4%) to Ostrinia nubilalis, and 443 hits (4.0%) to Ostrinia nubilalis. The second highest hits were to Dipteran species, with 526 hits to Culex quinquefasciatus and 437 and 252 hits to the Anopheles gambiae and Drosophila ananassae, respectively. The top 15 insect species with significant BLASTn hits are shown in Fig. 2.
Unigene function annotation
To analyze the putative protein functions, the Clusters of Orthologous Groups (COG) database was used to predict and classify the potential functions of the identified unigenes. Using sequence homology, 15,188 unigenes (43%) were annotated and divided into 25 specific categories (Fig. 3). Among the COG classifications, the cluster of general function (21.57%) was the largest, followed by signal transduction mechanisms (11.42%) and posttranslational modification, protein turnover, and chaperones (8.00%). The categories nuclear structure (0.45%), defense mechanisms (0.38%) and cell motility (0.31%) represented the smallest groups.
We used the Gene Ontology (GO) database to obtain functional annotations of the ACB unigenes . A total of 11,349 transcripts (32%) were assigned to at least one GO term in three categories: biological processes, cellular component and molecular function (Fig. 4). The three categories were further classified into 50 functional subcategories (Additional file 1), among which cell and cell part in the cellular component, catalytic activity and binding in the molecular function, and cellular processes and metabolic processes in the biological processes represented the major subcategories, while fewer than 10 unigenes were observed in translation regulator activity (8), synapse (7), synapse part (3), virion (3), cell killing (2), and channel regulator activity (1). Some unigenes were assigned to multiple GO term categories, while others could not be assigned to a given GO term. The biological process terms were primarily associated with cellular processes, such as proteolysis, carbohydrate metabolic processes and oxidation-reduction utilization. Similar compositions and distributions of unigenes assigned to GO terms have been reported in transcriptomics descriptions from other insects [9, 12].
Unigene metabolic pathway analysis
The unigene metabolic pathway analysis was conducted using the Kyoto Encyclopedia of Genes and Genomes (KEGG) annotation system. A total of 8512 unigenes (24%) were mapped to 186 KEGG pathways (Additional file 2). The most enriched pathways included metabolic pathways, biosynthesis of secondary metabolites, purine metabolism, focal adhesion, phagosome, spliceosome, and protein processing in the endoplasmic reticulum and the ribosome. These annotations provide a valuable resource for the investigation of specific processes, functions and pathways in the ACB.
Transcripts encoding detoxification enzymes and insecticide targets
The known mechanisms underlying insecticide resistance in the ACB include decreased penetration, increased detoxification and target insensitivity. To screen genes that may evolve insecticide resistance, we mined the current transcriptomics data to identify unigenes encoding insecticide targets or detoxification enzymes. A number of sequences homologous to detoxification enzymes (P450s, CCEs and GSTs) and insecticide targets were identified.
A total of 174 P450-related unigenes were identified from the Nr annotation of the ACB transcriptome. After manually removing allelic variants of the same P450 gene and unigenes with short open reading frames (ORFs) <600 bp, the remaining 31 P450 unigene sequences (Additional file 3) were used to construct a phylogenetic tree. Based on the closest BLAST hits in the NCBI Nr database and the phylogenetic analysis with P450 genes from B. mori, the P450s were assigned to one of four CYP clans: CYP2, CYP3, CYP4 and the mitochondrial clan. CYP3 ranked as the largest clan, consisting of 11 members belonging to the CYP6 family and seven genes belonging to the CYP9 family. The CYP4 clan included five P450s from the CYP4 family, and the remaining P450 genes belonged to the mitochondrial (6) and CYP2 (1) clans, which might be involved in the ecdysteroid metabolism pathway (CYP301 family) and essential physiological functions (CYP303-305 and CYP15 families), respectively (Fig. 5).
We also identified 63 GST unigenes from the Nr annotation, 19 of which were unique and manually curated (Additional file 4). Based on the closest BLAST hits in the NCBI Nr database and the phylogenetic analysis, these unigenes were assigned to epsilon, delta, theta, zeta, omega, sigma, and microsomal classes (Fig. 6).
A total of 171 unigenes with CCE protein motifs were identified. Among these, 27 unigenes were manually curated (Additional file 5), as some of the original sequences were either allelic variants of the same CCE gene or contained short ORFs. Based on phylogenetic analyses with other known CCE genes from other insect species, these 27 sequences were assigned to the known classes of CCEs (Fig. 7).
In addition, a number of unigenes encoding insecticide target proteins, including acetylcholinesterase (AChE), the nicotinic acetylcholine receptor subunits (nAChRs), the gamma-aminobutyric acid (GABA) receptor, the glutamate-gated chloride channel (GluCl), the voltage-gated sodium channel (VGSC) and the ryanodine receptor (RyR), were also identified in the ACB transcriptome (Table 2 and Additional file 6).
Gene expression profile analysis of responses to flubendiamide
The ACB mRNA samples used for constructing the cDNA library were non-normalized and non-amplified by PCR; thus, the reads in the sequencing dataset most likely represent the relative abundance level of each assembled transcript. To determine whether flubendiamide treatment resulted in statistically significant changes in gene expression, the level of gene expression was determined after normalizing the number of genes in each library to reads per kb per million reads (RPKM). Compared with the control, a total of 1306 genes were affected in the flubendiamide-treated strain, including 415 up-regulated genes and 891 down-regulated genes (Additional file 7, Table 3). More genes were expressed in flubendiamide-treated ACB group (31,478) than in the control (31,059).
To confirm the quality of the transcriptome data and the differential expression results from sequencing and computational analyses, 15 up-regulated genes (three esterase genes: unigene 13173, unigene 29498, unigene 15941; two cytochrome P450 genes: unigene 17631, unigene 12321; three hemolin genes: unigene 12340, unigene 586, unigene 1173; two twitchin genes: unigene 12898, unigene 3568; two heat shock protein 70 genes: unigene 1911, unigene 3432; two tropomyosin genes: unigene 25476, unigene 4859; one chemosensory protein gene: unigene 3898) and 13 down-regulated genes (three cuticular protein genes: unigene 3900, unigene 4018, unigene 1448; three protease genes: unigene 2960, unigene 6059, unigene 7750; four oxidoreductase genes: unigene 5356, unigene 5215, unigene 6772, unigene 2729; one neuronal pentraxin: unigene 249; one glutathione S-transferase: unigene 10933; one cytochrome P450 gene: unigene 2468) were randomly selected for further analysis. These genes were amplified using qRT-PCR, and the qRT-PCR assessments confirmed the direction of change in gene expression based on differential gene expression (DGE) analysis (Fig. 8). Cuticular protein, protease and oxidoreductase had much lower expression levels in response to flubendiamide treatment, while detoxification enzyme (P450s and esterase), muscle contraction-related (twitchin and tropomyosin), immunoglobulin (hemolin), heat shock protein 70 and chemosensory protein genes were significantly overexpressed.
The ACB is the most destructive corn-stalk-boring pest in the countries of East and Southeast Asia, the Pacific and Australasia. In some cases, entire harvests are lost in an outbreak year . However, a lack of genetic information is still a barrier to our knowledge of this insect pest. NGS, developed for rapid sequencing and transcriptomics research, has enabled the use of RNA-Seq and DGE methods to analyze the ACB transcriptome and insecticide mechanisms [14, 15].
For the ACB, as a non-model insect without a reference genome sequence, the clean reads were assembled de novo using the short reads assembling programs Mira and CAP3. A total of 35,430 unigenes with a mean length of 716 bp were obtained from the transcriptome analysis. The average length of the unigenes for ACB in the present study was longer than those observed in Plutella xylostella (532 bp)  and Liposcelis bostrychophila (711)  but shorter than those in Blattella germanica (798 bp)  and Agrotis ipsilon (770 bp) . The length of the unigenes may be related to the sequencing technique and assembly tools used. BLAST searches against public databases, such as Nt, Nr, Swiss-Prot, GO, COG and KEGG, provided annotation data for the ACB. As expected, the ACB transcripts were homologous to those from more than one species, but most were homologous to Lepidoptera species. The silkworm, B. mori, returned the most BLAST hits against the ACB transcripts. This result could reflect the fact that B. mori is a Lepidoptera species with a completely sequenced genome .
Insecticide detoxification occurs in all insects, and a number of enzymes with insecticide detoxification mechanisms, encoded by members of the P450, GST and CCE multi-gene families, have been identified . Cytochrome P450s are an important enzyme class involved in the metabolism of endogenous substances, such as hormones, fatty acids and steroids, and xenobiotics, such as drugs, pesticides and plant secondary metabolites . The P450 genes represent a large and highly diverse gene family in different species. CYP3 and CYP4 P450 families in insect species have been implicated in the metabolism of plant secondary metabolites and synthetic insecticides . In ACB, the identified cytochrome P450s were assigned to four CYP clans: CYP2, CYP3, CYP4 and the mitochondrial clan; CYP3 was the largest clan. Together, the members of the CYP4, CYP6 and CYP9 families account for 74.2% of the total P450s, indicating that the ACB, similar to other insects, arms itself with a potent capacity for metabolizing various xenobiotics. It has been reported that some P450 genes belonging to the CYP 4, 6, 9 and 12 families are involved in insecticide resistance . For example, the cytochrome P450 gene CYP6BQ23 was recently shown to be the main factor causing pyrethroid resistance in M. aeneus across Europe through its hydroxylation of the alcohol moiety of pyrethroids to a less toxic hydroxy metabolite . In B. tabaci and M. persicae, the cytochrome P450 (CYP6CM1 and CYP6CY3, respectively) over expression contributes to neonicotinoid insecticide resistance because these enzymes can catalyze a more rapid conversion of imidacloprid to its less active form, 5-hydroxy-imidacloprid [23, 24]. Moreover, RNA interference-mediated gene silencing (RNAi) tests revealed that CYP6BG1 overexpression results in permethrin resistance . In the present study, flubendiamide treatment led to up-regulated expression levels of the unigenes 12321, 17631, 11949 and 5498, annotated to CYP337B3, CYP9G3, CYP4C1 and CYP9f2, respectively. CYP337B3 is a chimeric enzyme, capable of metabolizing fenvalerate into non-toxic 4′-hydroxyfenvalerate . Li et al. reported that CYP9G3 plays important roles in the responses to phoxim and cypermethrin metabolism in Bombyx mandarina . Methanol exposure significantly up-regulated CYP9f2 mRNA expression; thus, the resulting protein was suggested as a candidate enzyme for methanol metabolism in Drosophila melanogaster . In the present study, these overexpressed P450 genes might improve the metabolism and detoxification capacity of the ACB to flubendiamide, and unigenes 12321, 17631, 11949 and 5498 could be involved in the flubendiamide resistance of the ACB.
Known CCEs can be divided into 13 clades; clades A to C are involved in insecticide detoxification . Unigene 13173, localized with mitochondrial, cytosolic and secreted esterases, showed 14.5-fold up-regulation after flubendiamide treatment. This esterase gene might be related to flubendiamide detoxification in the ACB. Unigenes 29498 and 15941, annotated to antennal carboxylesterases, were also dramatically up-regulated. In insects, antennal carboxylesterase plays an important role in the rapid degradation of odorants.
A total of eight classes of GSTs are known; among these, delta and epsilon GSTs are unique to insects and appear to play important roles in xenobiotic detoxification and insecticide resistance [16, 29]. In the ACB transcriptome, five GSTs belonging to the delta and epsilon classes were identified; Hence, it is likely that these five GSTs may have similar roles in flubendiamide detoxification in ACB. Among the other classes of GSTs, such as the omega class, the number of genes identified in the ACB (6) was larger than those in B. mori, L. bostrychophila, A. gambiae and A. mellifera [16, 30]. For the sigma and theta classes of GSTs, the numbers of genes (2 and 1, respectively) in the ACB were consistent with those from B. mori, and three genes were shown to encode microsomal GSTs, members of which are known in eicosanoid and glutathione metabolism [30, 31].
Five types of genes encoding insecticide target proteins, including the AChE, nAChRs, GABA receptor, GluCl, VGSC and RyR, were identified in the ACB transcriptome. Although many of these genes were not full-length sequences, this information will nevertheless facilitate further characterization of these targets using PCR and/or random amplification of cDNA ends (RACE). Insect AChEs have received much interest, reflecting their role in neurotransmission and because these genes are the principal targets of organophosphate and carbamate insecticides. Three unigenes were identified as containing the ace gene in the present study. The nAChRs represent a diverse family of cys-loop ligand-gated ion channels. Ten to twelve nAChR-type receptor gene families have been reported in insects . In the present study, five alpha and two beta subunits were also identified as ACB nAChRs. In contrast, 12 nAChR gene families were identified in B. mori , while 10 nAChR gene families were identified in A. gambiae . The GABA receptors also belong to a superfamily of cys-loop neurotransmitter receptors. Mutation of the GABA-regulated chloride channel can lead to resistance to organochlorine and phenylpyrazole insecticides . Some insect pest populations have evolved modifications of the VGSC protein, which prevent the binding of insecticides (i.e., DDT, pyrethrins and pyrethroids), resulting in the insect developing resistance.
RyRs are members of a superfamily of intracellular Ca2+ channels. These channels regulate the release of calcium from the lumen of the sarcoplasmic/endoplasmic reticulum to the cytosol of muscle and non-muscle cells. Insect RyRs are the targets of two classes of diamide chemicals [8, 35]. The homomeric tetramer RyR is the largest known ion channel, with each monomer consisting of approximately 5000 amino acids. A total of six unigenes were identified as part of the O. furnacalis RyR (OfRyR) gene, with the longest unigene being 4420 bp. In the present study, the OfRyR gene was up-regulated after flubendiamide exposure compared with the control, but this difference was not statistically significant. Actually, the gene expression of RyR changed depending on insect’s life cycle. The expression level of RyR differed among different developmental stages. For example, the OfRyR showed the lowest expression level in egg and the highest in adult . Meanwhile, the PxRyR was highly expressed in the second-instar larvae and adults, and their relative expression levels were 12.6 and 11.8 times than that in the prepupae, which showed the lowest expression level of PxRyR . Furthermore, the genes involved in calcium signaling and muscle control pathways, such as the twitchin genes (unigenes 4859 and 12898), were also significantly overexpressed after flubendiamide treatment. These results suggested that pathways involving multiple genes, such as calcium signaling, muscle control and metabolic pathways, play dominant roles in flubendiamide exposure.
Studying the temporal transcriptome responses of insects to xenobiotics could lead to the discovery of novel molecular mechanisms contributing to insecticide detoxification and tolerance. Gene expression analysis after flubendiamide exposure will significantly complement and enrich the ACB gene expression data and facilitate the discovery of novel genes, gene functions and insecticidal targets of the ACB. In insects, the response to environmental stress triggers the expression of various proteins, including heat shock proteins, metallothioneins or p-glycoprotein . The data obtained in the present study revealed that the transcription levels of heat shock protein 70 genes were significantly up-regulated in response to flubendiamide exposure. However, some genes encoding enzymes involved in the production of cellular catabolism, such as trypsin and serine protease, were down-regulated in the ACB in response to flubendiamide, and the cuticular protein genes also showed much lower expression levels in flubendiamide treatment. The down-regulation of cuticular-related proteins indicated that flubendiamide might play a role in thinning the cuticles in the ACB. Indeed, thickened cuticles affect the rates of insecticide penetration . Thus, thinned cuticles will accelerate the transportation of flubendiamide to the target and improve its insecticidal efficacy.
In conclusion, the present study represents a functional transcriptome analysis of the ACB. The differential gene expression data elucidated the molecular mechanisms of the ACB in response to the novel diamide insecticide flubendiamide. In addition, this investigation will facilitate identification of the genes involved in insecticide resistance and may assist in designing new compounds or other strategies for the control of the ACB.
The ACB laboratory colony was originally obtained from a stalk of field corn at the Hengshui experimental field of the Chinese Academy of Agricultural Sciences, Hebei Province, China. The colony was maintained on an agar-free semi-artificial diet in the laboratory at 26 ± 1 °C, with 70% relative humidity and a 16:8 (L:D) photoperiod. ACBs for the generation of cDNA libraries were obtained from two different strains. One strain was used as a control strain and was fed on a semi-artificial diet without insecticide treatment, and the other strain was exposed to a diet containing flubendiamide at the 95% lethal concentration value (3.8 μg g−1). After 48 h, more than 50 surviving 2nd instar larvae of each strain were collected in two separate tubes. The ACB samples were immediately snap-frozen in liquid N2 and stored at −80 °C until further use.
RNA extraction, cDNA library construction, and PGM sequencing
Total RNA from each sample (flubendiamide-treated and control ACBs) was isolated using the RNeasy® Mini Kit (Qiagen, Germany), according to the manufacturer’s instructions. The integrity of the total RNA was confirmed using both 1.5% agarose gel electrophoresis and an Agilent 2100 Bioanalyzer (Palo Alto, USA). The quantity of the total RNA was determined using a NanoDrop 1000 spectrophotometer (Thermo, USA). Additionally, mRNA samples were isolated from each total RNA using PolyATtract® mRNA isolation systems (Promega, USA). According to the manufacturer’s instructions, the mRNA was purified using polyT oligo-attached magnetic beads. The mRNA was broken into short fragments (200–800 nucleotides) in the presence of fragmentation buffer at 70 °C for 1 min. These short fragments were used as templates for first-strand cDNA synthesis with random primers. Subsequently, second-strand cDNAs were synthesized using a reaction system containing buffer, dNTPs, RNaseH and DNA polymerase I. These cDNA fragments were purified using a QIAquick Gel Extraction Kit (Qiagen, Germany) and were subjected to end repair by the addition of an A base to the 3′ end. Subsequently, the cDNAs were connected with sequencing adapters. Suitable fragments (400–1000 bp) were extracted and amplified by emulsion PCR (emPCR) to construct a cDNA library. The cDNA libraries were sequenced on a PGM platform using the Ion PGM ™ sequencing 400 Kit (Life Technologies, USA).
De novo transcriptome assembly and analysis
The raw reads were produced by base calling, which transforms pyroluminescence intensity signals to nucleotide sequences. Raw reads were preprocessed using the in-house developed programs Dust , FastqMcf (http://code.google.com/p/ea-utils/wiki/FastqMcf) and seq_crumbs (http://bioinf.comav.upv.es/seq_crumbs) to trim the adapters and poly A/T tails and to remove short (<100 bp), low-quality contaminating sequences. The resulting clean reads from both treatments were assembled into unigene sequences using the Mira  and CAP3 (http://seq.cs.iastate.edu) assembler programs at default parameters.
All de novo assembled unigenes were BLASTx searched against the Nr protein and Swiss-Prot databases. The sequences with no BLASTx hits were searched against the NCBI Nt database using BLASTn . The best matches were used to identify coding regions and to determine the sequence direction. GO, COG and KEGG analyses were used for the functional classification of the annotated unigenes. For the GO analysis, both BLAST2GO and WEGO were employed. These programs were used to extract the GO terms associated with homologies identified with BLAST, returning a list of GO annotations represented as hierarchical categories of increasing specificity . COG and KEGG annotations were analyzed using the Blastall program against the COG (http://www.ncbi.nlm.gov/COG/) and KEGG (http://www.genome.jp/kegg/) databases . Potential ORFs were identified using the ORF finder (http://www.ncbi.nlm.nih.gov/gorf/gorf.html). Simple sequence repeats were identified using MISA (http://pgrc.ipk-gatersleben.de/misa/) .
Analysis of genes related to insecticide detoxification and target proteins
Sequences encoding insecticide target proteins, such as AChE, nAChRs, GABA, GluCl, VGSC and RyR were identified from the Nr database with an E-value cut-off of 10−5. The genes related to insecticide-detoxification enzymes (P450s, CCEs and GSTs) were also manually curated using the above methods. All of the confirmed protein sequences of P450s, CCEs and GSTs were used for alignment and phylogenetic analysis. MEGA5.0 software was used to perform multiple sequence alignment prior to phylogenetic analysis and to construct consensus phylogenetic trees using the neighbor-joining method with p-distance under default parameters. Bootstrap analysis of 1000 replications was performed to evaluate the branch strength of each tree.
Differentially expressed unigene detection
The gene expression level was calculated using the RPKM method . The criteria of P < 0.01 and the absolute value of log2 ratio ≥ 1 were used to judge the significance of differences in gene expression between the two samples (flubendiamide-treated ACB and control). More stringent criteria, with smaller P and larger fold-change values, were used to identify differentially expressed genes .
qRT-PCR was performed to confirm the differently expressed genes. Total RNA from two samples (flubendiamide-treated and control ACBs) was extracted as described above. Complementary DNA was synthesized using an equal amount of total RNA (1 μg) from the two samples according to the instructions of the PrimeScript™ RT Reagent Kit with gDNA Eraser (Perfect Real Time) (Takara, Dalian, China). The primers employed in qRT-PCR are listed in Additional file 8. Each 20 μL of PCR master mix contained 2 μL of diluted cDNA template, 10 μL of SYBR Premix Ex Taq™, 0.4 μL of forward primer (10 μM), 0.4 μL of reverse primer (10 μM), 0.4 μL of Rox Reference Dye II (50×) and 6.8 μL of RNase-free H2O, according to the instructions of SYBR® Premix Ex Taq™ (Perfect Real Time) (Takara, Dalian, China). Thermal cycling was completed on a 7500 Real-time PCR system (Applied Biosystems, USA) using the following program: 95 °C for 30 s, followed by 40 cycles of 95 °C for 5 s and 60 °C for 34 s. Subsequently, a melting curve analysis from 60 °C to 90 °C was applied to all reactions to verify a single PCR product. After which the average cycle threshold (Ct) was calculated per sample. The housekeeping gene β-actin from the ACB was used as a previously evaluated internal control and was applied as an appropriate reference for insecticide-induced gene expression profiling . The fold change of target genes was calculated using the relative quantitative method (2-△△Ct) . To check reproducibility, four technical replicates and three biological replicates were carried out for each treatment. Statistical analysis was performed using Duncan’s Multiple Range Test for significance (P < 0.05) in SPSS 17.0 (SPSS, Inc., Chicago, IL, USA).
Asian corn borer
Clusters of orthologous groups
Cytochrome P450 monooxygenase
Expressed sequenced tag
Gamma-aminobutyric acid receptor
Glutamate-gated chloride channel
Kyoto Encyclopedia of Genes and Genomes
nicotinic acetylcholine receptor
Non-redundant protein database
Open reading frame
Polymerase chain reaction
Personal genome machine
Voltage-gated sodium channel
He KL, Wang ZY, Bai SX, Zheng L, Wang YB, Cui HY. Efficacy of transgenic Bt cotton for resistance to the Asian corn borer (Lepidoptera: Crambidae). Crop Prot. 2006;25:167–73.
Jin TT, Chang X, Gatehouse AMR, Wang ZY, Edwards MG, He KL. Downregulation and mutation of a cadherin gene associated with Cry1Ac resistance in the Asian corn borer, Ostrinia furnacalis (Guenée). Toxins. 2014;6:2676–93.
Yang DB, Zhang LN, Yan XJ, Wang ZY, Yuan HZ. Effects of droplet distribution on insecticide toxicity to Asian corn borers (Ostrinia furnaealis) and spiders (Xysticus ephippiatus). J Integr Agric. 2014;13:122–31.
Folcher L, Jarry M, Weissenberger A, Gérault F, Eychenne N, Delos M, Regnault-Roger C. Comparative activity of agrochemical treatments on mycotoxin levels with regard to corn borers and Fusarium mycoflora in maize (Zea mays L.) fields. Crop Prot. 2009;28:302–8.
Kato K, Kiyonaka S, Sawaguchi Y, Tohnishi M, Masaki T, Yasokawa N, Mizuno Y, Mori E, Inoue K, Hamachi I, Takeshima H, Mori Y. Molecular characterization of flubendiamide sensitivity in lepidopterous ryanodine receptor Ca2+ release channel. Biochemistry. 2009;48:10342–52.
Ebbinghaus-Kintscher U, Luemmen P, Lobitz N, Schulte T, Funke C, Fischer R, Masaki T, Yasokawa N, Tohnishi M. Phthalic acid diamides activate ryanodine sensitive Ca2+ release channels in insects. Cell Calcium. 2006;39:21–33.
Cordova D, Benner EA, Sacher MD, Rauh JJ, Sopa JS, Lahm GP, Selby TP, Stevenson TM, Flexner L, Gutteridge S, Rhoades DF, Wu L, Smith RM, Tao Y. Anthranilic diamides: A new class of insecticides with a novel mode of action, ryanodine receptor activation. Pestic Biochem Physiol. 2006;84:196–214.
Cui L, Yang DB, Yan XJ, Rui CH, Wang ZY, Yuan HZ. Molecular cloning, characterization and expression profiling of a ryanodine receptor gene in Asian corn borer, Ostrinia furnacalis (Guenée). PLoS ONE. 2013;8:e75825.
Zhou XJ, Qian K, Tong Y, Zhu JJ, Qiu XH, Zeng XP. De novo transcriptome of the hemimetabolous German cockroach (Blattella germanica). PLoS ONE. 2014;9:e106932.
Feyereisen R. Insect CYP genes and P450 enzymes. In: Gilbert LI, editor. Insect Molecular Biology and Biochemistry. London: Academic; 2012. p. 236–345.
Harris MA, Clark J, Ireland A, Lomax J, Ashburner M, Foulger R, Eilbeck K, Lewis S, Marshall B, Mungall C, Richter J, Rubin GM, Blake JA, Bult C, Dolan M, Drabkin H, Eppig JT, Hill DP, Ni L, Ringwald M, Balakrishnan R, Cherry JM, Christie KR, Costanzo MC, Dwight SS, Engel S, Fisk DG, Hirschman JE, Hong EL, Nash RS, et al. The Gene Ontology (GO) database and informatics resource. Nucleic Acids Res. 2004;32:258–61.
Tang LD, Wang XM, Jin FL, Qiu BL, Wu JH, Ren SX. De novo sequencing-based transcriptome and digital gene expression analysis reveals insecticide resistance-relevant genes in Propylaea japonica (Thunberg) (Coleoptea: Coccinellidae). PLoS One. 2014;9:e100946.
Xu LN, Ling YH, Wang YQ, Wang ZY, Hu BJ, Zhou ZY, Hu F, He KL. Identification of differentially expressed microRNAs between Bacillus thuringiensis Cry1Ab resistant and -susceptible strains of Ostrinia furnacalis. Sci Rep. 2016;5:15461.
Xu Z, Zhu W, Liu Y, Liu X, Chen Q, Peng M, Wang XZ, Shen GM, He L. Analysis of insecticide resistance-related genes of the carmine spider mite Tetranychus cinnabarinus based on a de novo assembled transcriptome. PLoS One. 2014;9:e94779.
Lin Q, Jin F, Hu Z, Chen H, Yin F, Li Z, Dong XL, Zhang DY, Ren SX, Feng X. Transcriptome analysis of chlorantraniliprole resistance development in the Diamondback moth Plutella xylostella. PLoS One. 2013;8:e72314.
Dou W, Shen GM, Niu JZ, Ding TB, Wei DD, Wang JJ. Mining genes involved in insecticide resistance of Liposcelis bostrychophila Badonnel by transcriptome and expression profile analysis. PLoS ONE. 2013;8:e79878.
Gu SH, Wu KM, Guo YY, Pickett JA, Field LM, Zhou JJ, Zhang YJ. Identification of genes expressed in the sex pheromone gland of the black cutworm Agrotis ipsilon with putative roles in sex pheromone biosynthesis and transport. BMC Genomics. 2013;14:636.
Xia QY, Zhou ZY, Lu C, Cheng DJ, Dai FY, Li B, Zhao P, Zha XF, Cheng TC, Chai CL, Pan GQ, Xu JS, Liu C, Lin Y, Qian JF, Hou Y, Wu ZL, Li GR, Pan MH, Li CF, Shen YH, Lan XQ, Yuan LW, Li T, Xu HF, Yang GW, Wan YJ, Zhu Y, Yu MD, Shen WD. A draft sequence for the genome of the domesticated silkworm (Bombyx mori). Science. 2004;306:1937–40.
Li X, Schuler MA, Berenbaum MR. Molecular mechanisms of metabolic resistance to synthetic and natural xenobiotics. Annu Rev Entomol. 2007;52:231–53.
Zimmer CT, Maiwald F, Schorn C, Bass C, Ott M, Nauen R. A de novo transcriptome of European pollen beetle populations and its analysis, with special reference to insecticide action and resistance. Insect Mol Biol. 2014;23:511–26.
Karatolos N, Pauchet Y, Wilkinson P, Chauhan R, Denholm I, Gorman K, Nelson DR, Bass C, Hffrench-Constant R, Williamson MS. Pyrosequencing the transcriptome of the greenhouse whitefly, Trialeurodes vaporariorum reveals multiple transcripts encoding insecticide targets and detoxifying enzymes. BMC Genomics. 2011;12:56.
Zimmer CT, Bass C, Williamson MS, Kaussmann M, Wölfel K, Gutbrod O, Nauena R. Molecular and functional characterization of CYP6BQ23, a cytochrome P450 conferring resistance to pyrethroids in European populations of pollen beetle, Meligethes aeneus. Insect Biochem Mol Biol. 2014;45:18–29.
Karunker I, Morou E, Nikou D, Nauen R, Sertchook R, Stevenson BJ. Structural model and functional characterization of the Bemisia tabaci CYP6CM1vQ, a cytochrome P450 associated with high levels of imidacloprid resistance. Insect Biochem Mol Biol. 2009;39:697–706.
Puinean AM, Foster SP, Oliphant L, Denholm I, Field LM, Millar NS, Williamson MS, Bass C. Amplification of a cytochrome P450 gene is associated with resistance to neonicotinoid insecticides in the aphid Myzus persicae. PLoS Genet. 2010;6:e1000999.
Bautista MAM, Miyata T, Miura K, Tanaka T. RNA interference-mediated knockdown of a cytochrome P450, CYP6BG1, from the diamond back moth, Plutella xylostella, reduces larval resistance to permethrin. Insect Biochem Mol Biol. 2009;39:38–46.
Joußen N, Agnolet S, Lorenz S, Schöne SE, Ellinger R, Schneider B, Heckela DG. Resistance of Australian Helicoverpa armigera to fenvalerate is due to the chimeric P450 enzyme CYP337B3. Proc Natl Acad Sci U S A. 2012;109:15206–11.
Li F, Ni M, Zhang H, Wang BB, Xu KZ, Tian JH, Hu JS, Shen WD, Li B. Expression profile analysis of silkworm P450 family genes after phoxim induction. Pestic Biochem Physiol. 2015;122:103–9.
Wang SP, He GL, Chen RR, Li F, Li GQ. The involvement of cytochrome P450 monooxygenases in methanol elimination in Drosophila melanogaster larvae. Arch Insect Biochem Physiol. 2012;79:264–75.
Ranson H, Claudianos C, Ortelli F, Abgrall C, Hemingway J, Sharakhova MV, Unger MF, Collins FH, Feyereisen R. Evolution of supergene families associated with insecticide resistance. Science. 2002;298:179–81.
Claudianos C, Ranson H, Johnson RM, Biswas S, Schuler MA, Berenbaum MR, Feyereisen R, Oakeshott JG. A deficit of detoxification enzymes: pesticide sensitivity and environmental response in the honeybee. Insect Mol Biol. 2006;15:615–36.
Sheehan D, Meade G, Foley VM, Dowd CA. Structure, function and evolution of glutathione transferases: implications for classification of non-mammalian members of an ancient enzyme superfamily. Biochem J. 2001;360:1–16.
Shao YM, Dong K, Zhang CX. The nicotinic acetylcholine receptor gene family of the silkworm, Bombyx mori. BMC Genomics. 2007;8:324.
Jones AK, Grauso M, Sattelle DB. The nicotinic acetylcholine receptor gene family of the malaria mosquito, Anopheles gambiae. Genomics. 2005;85:176–87.
Hosie AM, Baylis HA, Buckingham SD, Sattelle DB. Actions of the insecticide fipronil, on dieldrin-sensitive and -resistant GAGB receptors of Drosophila melanogaster. Br J Pharmacol. 1995;115:909–12.
Sun LN, Qiu GS, Cui L, Ma CS, Yuan HZ. Molecular characterization of a ryanodine receptor gene from Spodoptera exigua and its upregulation by chlorantraniliprole. Pestic Biochem Physiol. 2015;123:56–63.
Guo L, Tang BZ, Dong W, Liang P, Gao XW. Cloning, characterisation and expression profiling of the cDNA encoding the ryanodine receptor in diamondback moth, Plutella xylostella (L.) (Lepidoptera: Plutellidae). Pest Manag Sci. 2012;68:1605–14.
Wood O, Hanrahan S, Coetzee M, Koekemoer L, Brooke B. Cuticle thickening associated with pyrethroid resistance in the major malaria vector Anopheles funestus. Parasit Vectors. 2010;3:1–7.
Morgμlis A, Gertz EM, Schaffer AA, Agarwala R. A fast and symmetric DUST implementation to mask low-complexity DNA sequences. J Comput Biol. 2006;13:1028–40.
Chevreux B, Pfisterer T, Drescher B, Driesel AJ, Muller WE, Wetter T, Suhai S. Using the miraEST assembler for reliable and automated mRNA transcript assembly and SNP detection in sequenced ESTs. Genome Res. 2004;14:1147–59.
Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Cherry JM, Davis AP, Dolinski K, Dwight SS, Eppig JT. Gene Ontology: tool for the unification of biology. Nat Genet. 2000;25:25–9.
Salomonis N, Hanspers K, Zambon AC, Vranizan K, Lawlor SC, Dahlquist KD, Doniger SW, Stuart J, Conklin BR, Pico AR. GenMAPP 2: new features and resources for pathway analysis. BMC Bioinformatics. 2007;8:217.
Thiel T, Michalek W, Varshney RK, Graner A. Exploiting EST databases for the development and characterization of gene-derived SSR-markers in barley (Hordeum vμlgare L.). Theor Appl Genet. 2003;106:411–22.
Mortazavi A, Williams BA, Mccue K, Schaeffer L, Wold B. Mapping and quantifying mammalian transcriptomes by RNA-Seq. Nat Methods. 2008;5:621–8.
Audic S, Claverie JM. The significance of digital gene expression profiles. Genome Res. 1997;7:986–95.
Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, Adiconis X, Fan L, Raychowdhury R, Zeng QD, Chen ZH, Mauceli E, Hacohen N, Gnirke A, Rhind N, Palma F, Birren BW, Nusbaum C, Lindblad-Toh K, Friedman N, Regev A. Full length transcriptome assembly from RNA-Seq data without a reference genome. Nat Biotechnol. 2011;29:644–52.
Pfaffl MW. A new mathematical model for relative quantification in real-time RT-PCR. Nucleic Acids Res. 2001;29:2002–7.
This work was financially supported by grants from the National Natural Science Foundation of China (31601667), the National Key Research and Development Program of China (2016YFD0200500) and the Public Welfare Project from the Ministry of Agriculture of the People’s Republic of China (201203038). The funding bodies had no role in the design of the study and collection, analysis and interpretation of the data and drafting of the manuscript.
Availability of data and materials
The raw nucleotide reads obtained using Ion Torrent PGM sequencing were submitted to the Sequence Read Archive (SRA) database at NCBI with accession number SRP073016 (http://trace.ncbi.nlm.nih.gov/Traces/sra_sub/sub.cgi?subid=664950&from=list&action=show:submission). All other data are included within the article and in the Additional files.
LC, CHR and HZY conceived and designed the study. LC, CHR, DBY and HZY drafted the manuscript. LC conducted the sample collection, library construction, data processing, bioinformatics analysis, RT-PCR and qRT-PCR. ZYW coordinated the study. All authors contributed to the preparation of the final version. All authors read and approved the final manuscript.
The authors declare that they have no competing interests.
Consent for publication
Ethics approval and consent to participate
GO classifications of the ACB unigenes. (DOCX 15 kb)
KEGG pathway annotation of the ACB unigenes. (DOCX 25 kb)
Cytochrome P450 nucleotide sequences of the ACB. (TXT 48 kb)
Glutathione S-transferase nucleotide sequences of the ACB. (TXT 17 kb)
Carboxyl/cholinesterase nucleotide sequences of the ACB. (TXT 43 kb)
Nucleotide sequences of the target sites of the most important insecticide classes. (TXT 25 kb)
Differentially expressed genes between flubendiamide-treated and control ACBs. (XLS 27928 kb)
Information for the primers used in the qRT-PCR analysis. (DOCX 16 kb)
About this article
Cite this article
Cui, L., Rui, C., Yang, D. et al. De novo transcriptome and expression profile analyses of the Asian corn borer (Ostrinia furnacalis) reveals relevant flubendiamide response genes. BMC Genomics 18, 20 (2017). https://doi.org/10.1186/s12864-016-3431-6
- Gene Ontology
- Insecticide Resistance
- Cuticular Protein