Skip to main content
Download PDF

Characterization and expression profiling of microRNAs in response to plant feeding in two host-plant strains of the lepidopteran pest Spodoptera frugiperda

BMC Genomics volume 19, Article number: 804 (2018) Cite this article

Abstract

Background

A change in the environment may impair development or survival of living organisms leading them to adapt to the change. The resulting adaptation trait may reverse, or become fixed in the population leading to evolution of species. Deciphering the molecular basis of adaptive traits can thus give evolutionary clues. In phytophagous insects, a change in host-plant range can lead to emergence of new species. Among them, Spodoptera frugiperda is a major agricultural lepidopteran pest consisting of two host-plant strains having diverged 3 MA, based on mitochondrial markers. In this paper, we address the role of microRNAs, important gene expression regulators, in response to host-plant change and in adaptive evolution.

Results

Using small RNA sequencing, we characterized miRNA repertoires of the corn (C) and rice (R) strains of S. frugiperda, expressed during larval development on two different host-plants, corn and rice, in the frame of reciprocal transplant experiments. We provide evidence for 76 and 68 known miRNAs in C and R strains and 139 and 171 novel miRNAs. Based on read counts analysis, 34 of the microRNAs were differentially expressed in the C strain larvae fed on rice as compared to the C strain larvae fed on corn. Twenty one were differentially expressed on rice compared to corn in R strain. Nine were differentially expressed in the R strain compared to C strain when reared on corn. A similar ratio of microRNAs was differentially expressed between strains on rice. We could validate experimentally by QPCR, variation in expression of the most differentially expressed candidates. We used bioinformatics methods to determine the target mRNAs of known microRNAs. Comparison with the mRNA expression profile during similar reciprocal transplant experiment revealed potential mRNA targets of these host-plant regulated miRNAs.

Conclusions

In the current study, we performed the first systematic analysis of miRNAs in Lepidopteran pests feeding on host-plants. We identified a set of the differentially expressed miRNAs that respond to the plant diet, or differ constitutively between the two host plant strains. Among the latter, the ones that are also deregulated in response to host-plant are molecular candidates underlying a complex adaptive trait.

Background

Successful adaptation to the host-plant is of fundamental importance to herbivorous insects. It requires that the adult female accepts the host-plant for oviposition and that the host allows feeding and proper development of larvae. Understanding of the underlying genetic mechanisms used by insects in response to their host-plants ([1] for review) recently progressed thanks to availability of new reference genomes and transcriptomes of major polyphagous herbivorous insect pests. Genomic sequences analyses highlighted expansion of chemosensory and/or detoxification genes in generalist herbivores compared to specialist ones reflecting their larger diets [2,3,4,5,6]. RNA-seq analyses revealed that generalist herbivores use transcriptional plasticity of various categories of genes in response to their diet [3, 6,7,8,9]: detoxification, digestion, cuticular and ribosomal genes. Transcriptional plasticity of specific sets of genes has also been shown in an oliphagous lepidopteran species, Manduca sexta, when it is reared on host -or acceptable non-host plants [10]. These data show that a change of host-plant in herbivorous insect requires large scale transcriptional changes involving combinations of various gene family members.

While the nature and role of protein coding genes involved in adaptation to the host-plant begin to emerge, the putative role of non-coding genes has yet to be explored. Among them, microRNAs (miRNAs) are a class of small non coding RNAs (sncRNAs) of approximately 22 nt in length, which act as post-transcriptional regulators of gene expression and are known to help fine-tune complex genetic networks ([11], for review). The mode of action of miRNAs results in relatively weak modulation of less than twofold both at the RNA and protein levels [12]. Two other classes of small non coding RNAs combat the invasion and the expansion of transposable elements (TE), the short interfering RNA (siRNA) pathway suppresses TEs in all tissues of plant and animals, whereas the Piwi-interacting RNA pathway (piRNA) seems more specific to the gonads of metazoans (see [13] for review).

Since miRNAs have been shown to be involved in many physiological or cellular processes such as differentiation, proliferation, apoptosis and development [14], we hypothesized that they may play a role in adaptation of a phytophagous insect to its host-plants and that they may show different expression patterns in different host-plant races of the insect. To test this hypothesis, we used the noctuid moth Spodoptera frugiperda, which consists of two host-plant strains, one mostly associated to corn (C strain or SfC), the other to rice (R strain or SfR), and whose genomes are recently available [4]. We performed reciprocal transplant (RT) experiments of the two strains on the two host-plants and isolated and sequenced sncRNAs from feeding larvae. We present the differential expression patterns of miRNAs and their putative coding genes targets involved i) in phenotypic plasticity (the ability of a single genotype to produce multiple phenotypes in response to variation in the environment) of each strain in response to corn or rice ii) in adaptive evolution or genetic drift, by additional comparison of the two strains on the same host-plant.

Results

Deep sequencing of S. frugiperda small RNA

To characterize S. frugiperda miRNA, small RNA libraries were constructed from whole body of corn-fed and rice-fed larvae of S. frugiperda C and R strains. Two independent libraries (biological replicates) were prepared and independently sequenced using Illumina technology. Between 31.5 to 57 millions high quality sequence reads were obtained after adapter trimming in each library (See Additional file 1: Table S1). Their size distribution shows an over-abundance of sequences at 22 nt, typical of miRNAs (Fig. 1a-b). On average, we detected 11.12+/− 3.67% of sequences of size 22 nt on corn and 5.49+/− 1.06% on rice. We detected a second peak between 25 and 33 nt that could correspond to expression of piRNAs. After collapsing the reads, we identified between 2 to 5 million unique sequences, 4% of which (4.2+/− 0.2% on corn and 3.71+/− 0.52% on rice) on average corresponding to miRNAs of 22 nt (Fig. 1c, d). However, in term of diversity, piRNAs seem more abundant than miRNAs (Fig. 1c, d) which may reflect abundance and diversity of transposable elements (TE) in the genome of S. frugiperda (29.14% and 29.10% of genome coverage by TE in C and R strain, respectively [4]).

Fig. 1
figure1

Size profiling of small non coding RNAs and their homology to different RNA classes or to Transposable Elements (TE). The percentage of sncRNA reads is plotted as a function of their size (between 15 nt to 40 nt corresponding to the size range that has been selected from the gel for library construction), a and c SfC, c and d SfR, in green on corn, in red on rice. CC: SfC on corn, CR: SfC on rice, RC: SfR on rice, RR: SfR on rice. a and b total reads, c and d unique reads. e and f Pie charts representing the average % of reads (total counts from 2 replicates on corn for SfC (e) or SfR (f)) mapping either to SfC or plant miR precursors, or TE (SfC TE copies) as expected for putative endo-siRNA or piRNAs, or mRNA (SfC OGS2.2), or SfC tRNA, SfC rRNA (18S and 28S RNAs)

Full size image

In eukaryotic cells, the vast majority of cellular RNA consists of ribosomal RNA (rRNA) - ~ 80 to 90% of total RNA for most cells, followed by transfer RNA (tRNA) - 10 to 15%, messenger RNA (mRNA) - 3 to 7%, miRNA - 0.003 to 0.02% [15]. The small ncRNA sequences of SfC or SfR were aligned against these different references, in addition to the sequences of TE copies, which should highlight putative piRNAs or endo-siRNAs (Fig. 1e, f), The most abundant hits (29–29.1%) were found matching to Sf miRNAs precursors (precursors of known or novel miRNAs sequences) or Sf TE copies (15.5–18.4%) corresponding to putative piRNAs or endo-siRNAs. 16–16.7% of them matched to rRNA, and 0.2–2.4% to tRNA, the most abundant RNA classes expected in the cells. These data show that the sRNA sequences were enriched in functional miRNAs or TE-interacting RNAs compared to those resulting from degradation of ribosomal or transfer RNAs. We found also sncRNA sequences matching to miRNA precursors of plants (Zea or Oryza). It is expected since sncRNAs have been extracted from whole larvae feeding on plants. The lack of sequence homology between miRNA families in plants and animals (with exception of one family [16]) as well as differences in their biogenesis and mode of action suggested that miRNAs have evolved independently in both kingdoms from an ancient siRNA mechanism already present in the last common ancestor of all eukaryotes [12]. The presence of plant miRNA reads in larval samples is thus not expected to interfere with the study of animal miRNAs performed in this paper.

Annotation of S. frugiperda miRNA genes

miRDeep2 analysis

To detect miRNA genes, raw sequencing data were analyzed with miRDeep2 software. miRDeep2 [17, 18] maps the sRNA reads to the genome and excises potential miRNA precursors sequences from the genome. The secondary structures of the miRNA precursors are predicted and their stability is estimated by RNAfold. mirDeep2 uses a probabilistic model of miRNA biogenesis by the Dicer protein to score frequency and compatibility of mapping of the small RNA sequence reads (the signature) on the secondary structure of the miRNA genomic precursors (the structure) as compared to a non-miRNA precursor hairpin [17, 18]. Read stacks correspond to mature miRNA sequences. The score reflects the likelihood of each precursor to be a genuine miRNA. Furthermore, since the algorithm may generate hairpins with read stacks that have no connection with miRNA biology, corresponding to false positives, miRDeep2 estimates the rate of false positive by shuffling the observed combinations of structures and signatures and submit them to the core algorithm. The difference in score distribution between the genuine combinations and the control ones is used to estimate the number of false positive novel miRNAs for varying score cut-offs. The sequence of mature predicted miRNAs are compared to mature miRNA sequences contained in miRBase (release 21) which allows to sort them in two classes, known or novel depending if they existed or not in miRBase. The genomes of the C strain (v3.1) and the R strain (v1.0) were used as references using a pool of sequence reads from either the C or the R strain, respectively. As shown on Table 1, we obtained 76 and 68 known miRNAs predicted genes (70 and 64 unique ones) and 139 and 171 novel ones (126 and 158 unique ones) in the C and R strain.

Table 1 Number of miRNAs genes predicted in the Spodoptera frugiperda genome
Full size table

We calculated for the 139 novel miRs predicted in SfC, the processing precision frequency, defined as the ratio of reads corresponding exactly to the mature miRNA and miRNA* sequences, divided by the total amount of reads mapping to the hairpin [19] (See Additional file 2: Excel file S1, tabs “SfC Novel”). A value close to one indicates high precision and a value close to zero indicates the production of very few miRNA/miRNA* duplexes, with a cut-off of less than 0.1 corresponding to low processing precision [20]. One hundred thirty six novel miRs (97.8%) have an efficiency of more than 0.1, and 111 (79.8%) have an efficiency of > 0.5 thereby showing a medium to high processing precision in S. frugiperda and*/or that the mirDEEP2 prediction generated few false positives.

After filtering for the prediction showing a score > 4 (the lowest score cut-off corresponding to a prediction signal-to-noise ratio r > 10, r = total miRNA hairpins reported/ mean estimated false positive miRNA hairpin over 100 rounds of permutated control) and a significant randfold p-value, we obtained 66 and 59 genes for known miRNAs and 78 and 102 for novel ones in C and in R strain. All miRNA genes prediction can be found in Supplementary excel file 1 (See Additional file 2). Compared to the miRNAs identified from ovary cell lines of S. frugiperda (Sf21) [21] whose sequences are not in the miRBase release 21, we found 40 to 36 additional predictions of known miRNAs genes and 103 to 155 novel ones in C and R strain.

Other miRs involved in response to Baculovirus infection have been described from another ovary cell lines of S. frugiperda, Sf9, with precursor sequences mapped on Bombyx mori genome [22].

The specific nucleotide occurrence was analyzed in the obtained miRNA sequences, S. frugiperda showed a dominant bias for uracil (U) at the first nucleotide (Fig. 2). The dominance of U at first position towards 5’end is a conserved feature of miRNAs [23].

Fig. 2
figure2

Base composition of known or novel mature miRNAs

Full size image

Orthology

To identify orthologous miRNAs between C and R strain a reciprocal blastn was performed between mature sequences of known or novel miRNAs at an e-value < 0.001. We identified 57 orthologs among known miRNAs and 75 among novel ones. The list of orthologs can be found in Additional file 2: Excel file S1. This step allowed comparison of expression levels of miRNAs according to the genetic background (C or R strain).

miRNA expression variation

DESeq2 analysis according to host-plant

Expression values generated by miRDeep2 in the RT experiments were used to analyze differential expression of S. frugiperda miRNAs according to the host-plants in each strain. In Additional file 3: Figure S1, for SfC in the upper part and SfR in the lower part, the MA-plots show the log2 fold change on rice versus corn over the mean of normalized counts, i.e. the average of counts normalized by size factors. The points with FDR less than 0.05 are colored in red. We also visualized samples (treatments, replicates) by Principal Component Analysis (PCA) as shown in Additional file 4: Figure S2. Most variation was linked to treatment, i.e., change of host-plant, the smallest variation was found between biological replicates.

Based on read counts, 34 known or novel miRNAs (out of 144, 23.6%) were differentially expressed (FDR < 0.05) in the C strain larvae fed on rice as compared to the C strain larvae fed on corn (Details of DESeq2 analyses can be found in Additional file 5: Excel file S2.Twenty one known or novel miRNAs (out of 161, 13%) were differentially expressed on rice compared to corn in R strain. For example, in the C strain, known miR-34-5p and miR-190-5p were overexpressed on rice (in red on Fig. 3a), while novel miR-375-5p was overexpressed on corn (in green). In the R strain (Fig. 3b), novel tca-miR-375-5p and known miR-190-5p were overexpressed on rice as compared to corn (in red). Novel mmu-miR-155-3p was overexpressed on corn compared to rice (in green). We used TaqMan RT-qPCR miRNA assays on total RNA extracted from the RT experiments to validate expression variation based on read counts, for some of the most differentially expressed miRNAs (Table 2). We could confirm that known miR-34-5p and miR-190-5p were significantly overexpressed on rice compared to corn in the C strain and that known miR-190-5p and novel tca-miR-375-5p were significantly overexpressed on rice as compared to corn in the R strain (Table 2).

Fig. 3
figure3

Differential expression of miRNAs genes on rice compared to corn in Spodoptera frugiperda larvae (L4 instar), after rearing for 3 generations on whole plants. The miRNAs showing a significant differential expression after DESeq2 analysis (log2foldchange > 1 or < 1 and FDR < 0.05) are shown. a In SfC b In SfR. In red, miRNAs up-regulated on rice, in green miRNAs up-regulated on corn

Full size image
Table 2 Experimental validation of variation in miRNA expression
Full size table

miRNA expression variation between strains

We compared expression values between the two strains when they were reared on the same host plant to detect differences linked to the genetic background. In Additional file 6: Figure S3, on corn (top panel) and on rice (bottom panel), the MA-plots show the log2 fold change in SfR versus SfC over the mean of normalized counts. Again by PCA analysis of samples, more variation was found between strains than between biological replicates (See Additional file 7: Figure S4).

Based on reads counts, nine known or novel miRNAs (out of 129, 6.97%) were differentially expressed (FDR < 0.05) in the R strain compared to C strain when reared on corn (See Additional file 8: Excel file S3). A similar number of miRNAs (10 out of 129, 7.75%) were differentially expressed between strains on rice.

For example, known miR-34-5p and novel dme-miR-275-3p were overexpressed in SfR compared to SfC on corn (in red on top of Fig. 4), the latter being also overexpressed in SfR compared to SfC on rice (in red on bottom of Fig. 4). Using TaqMan RT-qPCR miRNA assays, we could confirm upregulation of miR-34-5p in SfR compared to SfC on corn, and that of the novel dme-miR-275-3p in SfR compared to SfC on rice (Table 2), however not on corn, maybe due to the presence of a RNA molecule competing with the microRNA for the TaqMan probe in this condition.

Fig. 4
figure4

Differential expression of miRNAs genes according to the genetic background. The relative expression of miRNAs in SfR compared to SfC is shown, either on corn (a) or on rice (b). In red, miRNAs up-regulated in SfR, in green miRNAs up-regulated in SfC

Full size image

Expression differences both between strains and between plants

The expression differences between host-strains may result from genetic drift and also possibly from divergent selection by the environment, in this case the different host-plants. To identify the ones putatively involved in adaptation to the host-plant, we looked for the miRNA genes of the known class that were differentially expressed both between strains and between plants (FDR < 0.05). As shown on Fig. 5a, four miRs (miR-10-5p, miR-34-5p, miR-263a-5p, miR278-5p) were differentially expressed both constitutively between strains on corn, and within SfC when reared on different plants. The two miRNAs genes that were differentially expressed between strains on rice (miR-279b-3p and miR-308-3p) were also differentially expressed in SfC on different plants. On Fig. 5b, we found that among the 4 miRNAs that were differentially expressed between SfR and SfC on corn, only one, miR-278-5p, was also differentially expressed when SfR was reared on rice compared to corn. None of the two miRNAs that were differentially expressed between strains on rice was differentially expressed in SfR according to the host-plant. Most constitutive differences between strains are involved in interaction with the plant in SfC, less of them in SfR, suggesting that adaptation to the plant played a more pronounced role in SfC evolution than in SfR.

Fig. 5
figure5

Are the constitutive expression differences between strains involved in phenotypic plasticity within strains? This Venn diagram highlights the miRNAs that are differentially expressed (FDR < 0.05) both between strains on the same plant and within strain (SfC: (a), SfR: (b)) on rice compared to corn

Full size image

Potential target genes regulated by miRNA

mRNA targets of known miRNAs were searched using TargetScan, miRanda, Rna22 and miRmap (See Methods) against the 3’ UTR of the C strain gene set, OGS2.2. The target gene list of the seven DE miRNAs identified in this study (miR-34, miR-190, miR-1a-5p, miR-998, miR-278, miR-263a-5p, miR-10-5p) can be found in Additional file 9: Excel file S4. Among the gene targets predicted by TargetScan, we filtered out those showing variation in their expression in the RT experiments (according to [7]). We also checked whether they were predicted as miR targets by at least one of the three other softwares. The list of target genes that are differentially expressed in SfC on different host plants can be found in Additional file 10: Excel file S5. The target genes that are differentially expressed in SfR compared to SfC on corn are listed in Additional file 11: Excel file S6.

Potential regulated target genes

To reduce the number of putative false positive coding gene targets among those predicted by TargetScan, we assumed that true targets should be expressed in the same experimental conditions as miRNA genes, and downregulated when miRNA genes were overexpressed. To identify these candidate genes, we used RNA-Seq results obtained in the same reciprocal transplant experiments after two generations on plants [7]. We limited the analysis of target genes to the Official Gene Set of the C strain, which is the reference and has been manually curated [4], and to RT experimental conditions in which the RNA-Seq data of putative target genes were available in duplicates in [7]. We thus provide a detailed analysis of targets genes in the following conditions: 1) SfC reared on corn and on rice 2) SfC and SfR when reared on corn. We focused on the most differentially expressed miRNAs of “known” class (FDR < 0.05 and absolute value of log2 fold change > 1.3) in these conditions. The complete list of differentially expressed (FDR < 0.05) coding gene targets of these known miRNA genes, with their log2 fold change and annotation can be found in Additional file 10: Excel file S5 and Additional file 11: Excel file S6. MiR-34 and miR-190 are up-regulated in SfC when reared on rice, a non-host plant. Among down-regulated targets of miR-34, we found members of different gene families (Table 3), some of which are also targets of miR-190: First, a representative of the takeout gene family, then three genes encoding cuticle proteins. A gene encoding a subunit of the 26S proteasome was specifically targeted by miR34. As miR-190 specific targets, were found a cuticle protein gene, an acyl-CoA desaturase gene, and a member of the Osiris gene family, osi9a [7].

Table 3 List of target genes of known miRs and their differential expression
Full size table

Among miRNAs that are overexpressed in SfC on corn compared to rice, three, miR-998, miR-263a-5p and miR-10-5p share two targets that encode transporters (facilitated trehalose transporter Tret1 and a sodium- and of chloride-dependent glycine transporter 1) that are downregulated on corn. In addition, a gene encoding a spermine oxidase enzyme is the target of both miR-998 and miR-263a-5p.

When we compared the two strains on the same host-plant corn (See Additional file 11: Excel file S6), we found that miR-34 was up-regulated in SfR and that miR-263a-5p, miR10-5p and miR-278-5p were up-regulated in SfC. Among the down-regulated targets of miR-34, we found a representative of the transcription factor family AP-2 (Table 3). Among the down-regulated targets of miR-263a-5p, were found a protein of unknown function and a cuticular protein. Among the targets of miR10-5p, two genes coding for transporters (proton-coupled folate transporter, facilitated trehalose transporter Tret1) were found.

Because of the energy required for freeing base-pairing interactions within mRNA in order to allow miRNA binding, secondary structures have been shown to contribute to target recognition by miRNAs [24]. Using mfold [25] to predict secondary structures of 3’UTR of predicted targets, we could show that miR-34 and miR-190 map to loops rather than stems in the secondary structures of their targets, takeout (GSSPFG00021718001-RA) and acyl-CoA desaturase (GSSPFG00006314001-RA), their most differentially expressed candidates respectively (Fig. 6a and b), which may facilitate the interaction.

Fig. 6
figure6

Examples of complementarity between miRNAs seed sequences and the secondary structure of their putative targets. Using mfold [25] to predict and draw secondary structures of 3’UTR of the predicted targets, we could show that miR-34 and miR-190 map to loops rather than stems in the secondary structures of their targets, takeout (GSSPFG00021718001-RA) in (a) and acyl-CoA desaturase (GSSPFG00006314001-RA) in (b), their most differentially expressed candidates, respectively

Full size image

Discussion

In another paper [7], we have shown that the two host-plant strains of Spodoptera frugiperda show phenotypic and transcriptional plasticity when reared on their preferred versus alternative host-plants, corn or rice. In this paper, using similar reciprocal transplant experiments, we have shown variation in expression profile of miRNAs in the two lineages according to the host plant or between lineages on the same host-plant, either corn or rice. We identified putative targets of differentially expressed miRNAs in the 3’UTR of coding genes, and present a detailed analysis of those that are expressed i) in the same experimental condition as miRNAs ii) in the opposite direction (downregulated when miRNAs are upregulated and vice versa), that we consider the most reliable candidates.

Among down-regulated targets of miR-34, we found (See Additional file 10: Excel file S5, Table 3) members of different gene families, some of which are also targets of miR-190: First, a representative of the takeout gene family, putatively involved in feeding behavior and response to starvation as in D. melanogaster [26]. This representative is targeted both by miR-34 and miR-190. Second, three genes encoding cuticle proteins whose downregulation may reflect slower development of SfC on rice compared to corn [7]. Another cuticle component is also targeted by miR-190. Among miR-34 specific targets, we found a gene encoding a subunit of the 26S proteasome, which can be involved in protein degradation in response to oxidative stress. In the case of phytophagous insects, oxidative stress can be generated by prooxidant allelochemicals produced by host-plants. Among miR-190 specific targets, an acyl-CoA desaturase gene, necessary for fatty acid biosynthesis, may be repressed because rice is a poor food for SfC, and a member (Osiris 9) of the Osiris gene family, putatively involved in response to plant toxins as in Drosophila sechellia [27]. A recent analysis of the conserved patterns of Osiris gene expression in different insect species, suggests that Osiris genes may play a central role in insect adaptive evolution [28].

The three miRNAs miR-998, miR-263a-5p and miR-10-5p that are overexpressed in SfC on corn share two targets that encode transporters (facilitated trehalose transporter Tret1 and a sodium- and of chloride-dependent glycine transporter 1) that are downregulated on corn (up-regulated on rice). In most insects, trehalose (a-D-glucopyranosyl-(1,1)-a-D-glucopyranoside) is the main haemolymph sugar. In Drosophila, Tret 1 is necessary for the transport of trehalose produced in the fat body and its uptake into other tissues that require a carbon source, and thereby regulates trehalose levels in the hemolymph [29]. The glycine transporter GLYT1, by controlling the reuptake of glycine at synapses [30], regulates neurotransmission, where glycine plays the role of inhibitory neurotransmitter. In addition a gene encoding a spermine oxidase enzyme is the target of both miR-998 and miR-263a-5p. Polyamines (PA), comprising spermine (Spm), spermidine (Spd) and putrescine (Put), are ubiquitous polybasic molecules, with many important biological functions, like cell growth, differentiation, and apoptosis. They interact reversibly with nucleic acids, regulating chromatin status and gene expression, and modulating ion-channels’ function and stability (reviewed in [31]). Polyamine oxidases (PAO) include spermine oxidase. These enzymes, containing a flavin adenine dinucleotide (FAD), catalyze the oxidation of polyamines, and lead to formation of hydrogen peroxide. They regulate cellular polyamine concentration.

When we compared the two strains on the same host-plant corn (See Additional file 11: Excel file S6), we found that miR-34 was up-regulated in SfR and that miR-263a-5p, miR10-5p and miR-278-5p were up-regulated in SfC. Among the down-regulated targets of miR-34, we found a representative of the transcription factor family AP-2. From [32], AP-2 is expressed during Drosophila embryogenesis in the maxillary segment and neural structures, whereas during larval development, it is expressed in the central nervous system (CNS) and the leg, antennal and labial imaginal disks [33, 34]. In a Drosophila AP-2 mutant, defects in proboscis development and leg-joint formation have been described [35, 36]. We found a homolog of the graves disease carrier protein, a protein of as yet uncharacterized function that belongs to the mitochondrial metabolite carrier family (which includes the ADP/ATP translocator, the phosphate carrier and the hydrogen ion uncoupling protein). Among the down-regulated targets of miR-263a-5p, were found a protein of unknown function and a cuticular protein. Among the down-regulated targets of miR-10-5p, two genes coding for transporters (proton-coupled folate transporter, facilitated trehalose transporter Tret1) were identified. Folates are a family of B9 vitamins found in nature primarily as 5-methyltetrahydrofolate (5-methylTHF). 5-MethylTHF provides the methyl group for the synthesis of methionine from homocysteine and therefore is necessary for the formation of S-adenosylmethionine, which is required for a variety of methylation reactions [37]. In addition, folates are the sources of a methylene moiety for the de novo synthesis of thymidylate from deoxyuridylate and two formate moieties for the de novo synthesis of the purine ring. PCFT has been extensively studied in humans due to the importance of folates in cancer progression, however it has not been studied in insects. In humans, PCFT is expressed at the acidic microenvironment of the apical brush-border membrane of the proximal small intestine and allows the intestinal absorption of folates. In humans and mice, loss of function mutations in PCFT, results in severe systemic folate deficiency with anemia, sometimes pancytopenia, hypo-immunoglobulinemia and gastrointestinal defects [38].

This analysis of putative targets of DE miRNA reveals that they control important pathways necessary for cellular or organismal homeostasis during insect-plant interaction. Although we provide transcriptomic evidence of downregulation of putative targets, their experimental validation will require additional efforts like the use of biotin tagged miRNAs to capture them and/or translation profiling [39, 40], which will be the focus of future work.

For expression analysis of the predicted target genes, insect samples were collected after two generations on plants whereas for the study of miRNA expression they were collected after three generations. This was to ensure dilution of miRNAs putatively maternally inherited (it is the case for miR-34 in D. melanogaster for instance [41]) and related to the artificial diet on which the insect fed before the reciprocal transplant experiment on plants started. We considered that in controlled experimental conditions of growth, gene expression during each successive insect developmental cycle is reproducible and comparable from one generation to another (Since S. frugiperda is a quarantine organism in France, we performed the RT experiment in large incubators with controlled hygrometry, temperature and light conditions). If the plant exerted a selection pressure, gene allele frequencies may have changed but not significantly between two generations since the laboratory strain that we used has a limited genetic polymorphism.

For differential expression analysis of miR genes and their targets, two biological replicates of the reciprocal transplant experiment have been performed. Since the first use of RNA-Seq to analyze transcriptomic data [42], we know more on the parameters necessary to optimize detection of differentially expressed genes. Lamarre et al. recently showed that the read depth in part compensates for the number of replicates to increase the ratio of differentially expressed genes detected [43]: with 20 million reads (for 20,000 DE genes - corresponding to a coverage of 1000 reads per gene) and 2 replicates, 85% of DE genes are found, and 8 replicates are needed to reach a ratio of 100%, while with 2.5 million reads (~ coverage 125 reads/gene) and 2 replicates, only 15% of DE genes are detected and the use of 8 replicates increases the ratio of DE genes to only 60%. In our miR differential expression analysis, we performed only two biological replicates. However, the read depth in each SfC library (Additional file 1: Table S1) was 42 million on average, with more than 20% of the sequences having the expected size for miRs (21 to 23 nt), this corresponds to an average coverage of > 39,688 reads per miR gene, a higher coverage than the maximal one used in the study of Lamarre et al.. Lamarre et al. also showed that the rate of false positives obtained with DESeq2 is minimal with 2 replicates and increases with the number or replicates. According to them, 70% of true positives can be detected with two replicates and this number increases with the number of replicates. We conclude that our experimental design enabled detection of a reasonable number of reliable DE miR candidates although more repetitions may be necessary to deepen the study.

In the same line, our differential expression analysis was done from RNA extracted from whole body of the larvae feeding on plants. We are aware that this experimental design may underestimate the number of miRs showing variation in their expression due to the fact that overexpression in one tissue may be masked by down regulation in another one. It will be the subject of future effort to look for expression variation tissue by tissue.

The miRNAs expression differences that we uncovered between the host-strains may result from genetic drift and also possibly from divergent selection by the environment, in this case the different host-plants. To identify the latter, we searched for the miRNAs that were differentially expressed both between strains on the same plant and within strains in response to the host-plant. Interestingly, we found that all the miRNAs expression difference between strains on plants (6/6, miR-10-5p, miR-34-5p, miR-263a-5p, miR-278-5p, miR279b-3p, miR-308-3p) were also involved in response to the host-plant in SfC. By contrast, one miRNA only (miR-278-5p) out of the six that showed constitutive difference between strains, showed also variation in response to the host-plant in SfR. By measuring the fitness of the insects on plants, we had found better survival of the C strain on corn compared to rice [7] suggesting that SfC was adapted to corn. The miRNAs that are both deregulated between strains and in response to plant may be involved in this adaptive evolution of SfC.

MiRNAs have been shown to be involved in response to various environmental changes, like response to starvation in C. elegans [44], to freezing and anoxia stress in the freeze tolerant fly Eurosta solidaginis [45], to thermal plasticity of the Senegalese sole [46], to drought in Tobacco [47]. In the case of Spodoptera frugiperda, by regulating phenotypic plasticity, miRNAs may have also played a role in evolutionary adaptation as has been discussed in the case of human miRNAs [48]. Further work is needed to show that the miR candidates identified in this study are directly involved in fitness of the insect on its host-plant, however the possibility of transgenerational inheritance of these molecules in male or female gametes [41, 49, 50] suggests that they could facilitate transmission of complex adaptive traits from parents to offspring.

Conclusions

In the current study, we performed the first systematic analysis of miRNAs in Lepidopteran pests reared on whole host and non-host plants. We identified a set of the differentially expressed miRNAs that respond to the plant diet, or differ constitutively between the two host plant strains. The analysis of the putative targets of these DE miRNAs revealed that they control important pathways necessary for cellular or organismal homeostasis during insect-plant interaction. Since two classes of non coding RNAs, siRNAs and miRNAs have been used to control insect development [51, 52], this study of regulatory molecules of insect - plant interaction can bring clues on novel environment friendly biological control of crop pests.

Methods

Biological material

Two laboratory strains of S. frugiperda were used in this study: the corn-strain, originated from French Guadeloupe and the rice-strain, originated from Florida, USA. Each strain was reared on its principal and alternative host plant (Zea mays L. cv B73 or Oryza sativa L. japonica cv Arelate) under controlled conditions (temperature: 24 °C, photoperiod 16:8 light:dark, relative humidity: 65%). Eggs from S. frugiperda were deposited on plant leaves and fed ad libitum. After three generations on plants, larvae of fourth instars (L4) were collected. Therefore, larvae selected for small RNA extraction and sequencing comprised four groups: C strain reared on corn (CC), R strain reared on corn (RC), C strain reared on rice (CR) and R strain reared on rice (RR). Each experiment comprised two independent replicates.

Small RNA extraction and sequencing

Small RNA was extracted from 12 L4 larvae using the mirVana miRNA Isolation Kit (Ambion) according to the manufacturer’s instructions. Extracted small RNA was quantified using a NanoDrop 2000 spectrophotometer (Thermo Scientific) and analyzed by PAGE and silver-staining. Small RNA samples from each replicate of CC, RC, CR and RR larvae were used for cDNA library preparation using the TruSeq ®Small RNA Sample Prep Kit (Illumina) and sequencing on Illumina HiSeq 2500. A size selection targeting snRNAs in the range of 15 to 40 nt was done. Library preparation and sequencing were performed by the platform MGX-Montpellier GenomiX (Montpellier, France).

Computational analysis of small RNA sequencing data and miRNA identification.

High throughput sequencing generates small RNA reads of 50 nucleotides in length (single reads). Raw reads were trimmed to remove adapter sequences using cutadapt software (version 1.4.1) [53]. For downstream analysis, only reads having length between 15 and 40 bases were considered.

miRNA genes annotation

miRDeep2 algorithm [18] was used to detect miRNA from small RNA deep sequencing data. This algorithm uses a probabilistic model to score the fit of sequenced RNA to the biological model of miRNA biogenesis [17]. Briefly, reads are aligned to the S. frugiperda genome (Corn variant assembly version 3.1 or Rice variant assembly version 1.0), and only reads that do not map more than five times to the genome were used for miRNA detection. Then, using the read mappings as guidelines, potential miRNA precursors are excised from the genome and the miRDeep2 core algorithm scores their likelihood to be a real miRNA precursor. The output is a scored list of known and novel miRNA in the deep sequenced sample. Known miRNAs were identified by similarity to miRNA sequences from miRBase database (release 21).

Orthology analysis

To explore the conservation of miRNAs between C and R strains, a reciprocal blastn was performed to search for orthologs [54, 55]. BLAST for miRNA mature sequences was run with blastn with default parameters except for a lower e-value threshold of 1e-3.

Analysis of differential miRNA expression

miRDeep2 algorithm also provides read counts for the detected miRNA. To assess changes in miRNA expression between S. frugiperda variants and host-plant conditions, the read counts data for known and novel miRNA were used as input for the R package DESeq2 [56]. DESeq2 uses negative binomial generalized linear models to test for differential expression. An adjusted p-value for multiple testing was computed with the Benjamini-Hochberg procedure to control false discovery rate (FDR). Results with a FDR < 0.05 were considered statistically significant.

Experimental validation of differential expression

The miRNA expression levels were quantified using TaqMan small RNA assay system from Life Technologies. Briefly, total RNA from samples was isolated using Trizol. One to 10 ng of total RNA was used for reverse transcription using a specific RT primer with the following conditions, i.e., 16 °C: 30 min, 42 °C: 30 min, 85 °C: 5 min. Subsequently, the cDNA was used for qRT analysis with TaqMan probes according to the manufacturer’s instructions. For all qRT based TaqMan assays, the qRT-PCR quantifications were performed on two biological replicates (pools of 12 larvae) of the reciprocal transplant experiments, with 3 technical replicates. During the qRT analysis, the 2-ΔΔCT method was employed and each Ct value of the test miR was normalized to that of an endogenous miRNA (Sf mir-31 5p) whose expression remained stable in the different experimental conditions that were tested in C and R strain. To check that the miR expression level differed significantly between two experimental conditions, we used a pair wise fixed reallocation randomization statistical test [57] (2000 iterations, p-value< 0.001) which avoids making any assumptions about distributions compared to standard parametric tests such as analysis of variance or t-tests.

Detection and functional annotation of potential target genes regulated by miRNA

TargetScan [58–60] predicts biological targets of miRNAs by searching for the presence of conserved 8-mer and 7-mer sites that match seed region of each miRNA by calculating thermodynamic free energy using the RNAFold package [61]. Predictions are ranked using the site number, site type, and site context. TargetScan (version 5.0) was run with default parameters using mature known miRNAs and the 3’UTR of predicted coding genes set OGS2.2_UTR3 of the C strain available on the webportal [62]. Predictions obtained by TargetScan were confirmed by at least one of the three following softwares: MiRanda, Rna22 and miRmap. MiRanda (version v3.3a) [63] allows one wobble pairing in the seed region when it is compensated by matches in the 3′ end of the miRNA, it calculates the binding energy of the duplex structure and its position within the 3’UTR, it was used with the same parameters as in [64]. Rna22 (version v2) [65] is a tool based on a search for patterns that are statistically significant miRNA motifs created after a sequence analysis of known mature miRNAs, it was used with default parameters. MiRmap [66] is a web-based application [67] that combines many thermodynamic, evolutionary, probabilistic and sequence-based features.

References

  1. 1.

    Simon JC, d'Alencon E, Guy E, Jacquin-Joly E, Jaquiery J, Nouhaud P, Peccoud J, Sugio A, Streiff R. Genomics of adaptation to host-plants in herbivorous insects. Brief Funct Genomics. 2015;14(6):413–23.

    CAS  Article  Google Scholar 

  2. 2.

    Xu W, Papanicolaou A, Zhang HJ, Anderson A. Expansion of a bitter taste receptor family in a polyphagous insect herbivore. Sci Rep. 2016;6:23666.

    CAS  Article  Google Scholar 

  3. 3.

    Pearce SL, Clarke DF, East PD, Elfekih S, Gordon KHJ, Jermiin LS, McGaughran A, Oakeshott JG, Papanikolaou A, Perera OP, et al. Genomic innovations, transcriptional plasticity and gene loss underlying the evolution and divergence of two highly polyphagous and invasive Helicoverpa pest species. BMC Biol. 2017;15(1):63.

    CAS  Article  Google Scholar 

  4. 4.

    Gouin A, Bretaudeau A, Nam K, Gimenez S, Aury JM, Duvic B, Hilliou F, Durand N, Montagné N, Darboux I, et al. Two genomes of highly polyphagous lepidopteran pests (Spodoptera frugiperda, Noctuidae) with different host-plant ranges. Sci Rep. 2017;7(1):11816.

  5. 5.

    Cheng T, Wu J, Wu Y, Chilukuri RV, Huang L, Yamamoto K, Feng L, Li W, Chen Z, Guo H, et al. Genomic adaptation to polyphagy and insecticides in a major east Asian noctuid pest. Nat Ecol Evol. 2017;1(11):1747–56.

    Article  Google Scholar 

  6. 6.

    Mathers TC, Chen Y, Kaithakottil G, Legeai F, Mugford ST, Baa-Puyoulet P, Bretaudeau A, Clavijo B, Colella S, Collin O, et al. Rapid transcriptional plasticity of duplicated gene clusters enables a clonally reproducing aphid to colonise diverse plant species. Genome Biol. 2017;18(1):27.

    Article  Google Scholar 

  7. 7.

    Orsucci M, Moné Y, Audiot P, Gimenez S, Nhim S, Nait-Saidi R, Frayssinet M, Dumont G, Pommier A, Boudon JP, et al. Transcriptional plasticity evolution in two strains of Spodoptera frugiperda (Lepidoptera: Noctuidae) feeding on alternative host-plants: BioRxiv; 2018. https://doi.org/10.1101/263186.

  8. 8.

    Celorio-Mancera MDLP, Heckel DG, Vogel H. Transcriptional analysis of physiological pathways in a generalist herbivore: responses to different host plants and plant structures by the cotton bollworm, Helicoverpa armigera. Entomol Exp Appl. 2012;144:123–33.

    Article  Google Scholar 

  9. 9.

    de la Paz Celorio-Mancera M, Wheat CW, Vogel H, Söderlind L, Janz N, Nylin S. Mechanisms of macroevolution: polyphagous plasticity in butterfly larvae revealed by RNA-Seq. Mol Ecol. 2013;22:4884–95.

    CAS  Article  Google Scholar 

  10. 10.

    Koenig C, Bretschneider A, Heckel DG, Grosse-Wilde E, Hansson BS, Vogel H. The plastic response of Manduca sexta to host and non-host plants. Insect Biochem Mol Biol. 2015;63:72–85.

    CAS  Article  Google Scholar 

  11. 11.

    Cora D, Re A, Caselle M, Bussolino F. MicroRNA-mediated regulatory circuits: outlook and perspectives. Phys Biol. 2017;14(4):045001.

    Article  Google Scholar 

  12. 12.

    Moran Y, Agron M, Praher D, Technau U. The evolutionary origin of plant and animal microRNAs. Nat Ecol Evol. 2017;1(3):27.

    Article  Google Scholar 

  13. 13.

    Ghildiyal M, Zamore PD. Small silencing RNAs: an expanding universe. Nat Rev Genet. 2009;10(2):94–108.

    CAS  Article  Google Scholar 

  14. 14.

    Catalanotto C, Cogoni C, Zardo G. MicroRNA in Control of Gene Expression: An Overview of Nuclear Functions. Int J Mol Sci. 2016;17(10):1712.

    Article  Google Scholar 

  15. 15.

    Palazzo AF, Lee ES. Non-coding RNA: what is functional and what is junk? Front Genet. 2015;6:2.

    Article  Google Scholar 

  16. 16.

    Arteaga-Vazquez M, Caballero-Perez J, Vielle-Calzada JP. A family of microRNAs present in plants and animals. Plant Cell. 2006;18(12):3355–69.

    CAS  Article  Google Scholar 

  17. 17.

    Friedlander MR, Chen W, Adamidi C, Maaskola J, Einspanier R, Knespel S, Rajewsky N. Discovering microRNAs from deep sequencing data using miRDeep. Nat Biotechnol. 2008;26(4):407–15.

    Article  Google Scholar 

  18. 18.

    Friedlander MR, Mackowiak SD, Li N, Chen W, Rajewsky N. miRDeep2 accurately identifies known and hundreds of novel microRNA genes in seven animal clades. Nucleic Acids Res. 2012;40(1):37–52.

    Article  Google Scholar 

  19. 19.

    Ma Z, Coruh C, Axtell MJ. Arabidopsis lyrata small RNAs: transient MIRNA and small interfering RNA loci within the Arabidopsis genus. Plant Cell. 2010;22(4):1090–103.

    CAS  Article  Google Scholar 

  20. 20.

    Pradhan B, Naqvi AR, Saraf S, Mukherjee SK, Dey N. Prediction and characterization of tomato leaf curl New Delhi virus (ToLCNDV) responsive novel microRNAs in Solanum lycopersicum. Virus Res. 2015;195:183–95.

    CAS  Article  Google Scholar 

  21. 21.

    Kakumani PK, Chinnappan M, Singh AK, Malhotra P, Mukherjee SK, Bhatnagar RK. Identification and characteristics of microRNAs from army worm, Spodoptera frugiperda cell line Sf21. PLoS One. 2015;10(2):e0116988.

    Article  Google Scholar 

  22. 22.

    Mehrabadi M, Hussain M, Asgari S. MicroRNAome of Spodoptera frugiperda cells (Sf9) and its alteration following baculovirus infection. J Gen Virol. 2013;94(Pt 6):1385–97.

    CAS  Article  Google Scholar 

  23. 23.

    Lau NC, Lim LP, Weinstein EG, Bartel DP. An abundant class of tiny RNAs with probable regulatory roles in Caenorhabditis elegans. Science. 2001;294(5543):858–62.

    CAS  Article  Google Scholar 

  24. 24.

    Kertesz M, Iovino N, Unnerstall U, Gaul U, Segal E. The role of site accessibility in microRNA target recognition. Nat Genet. 2007;39(10):1278–84.

    CAS  Article  Google Scholar 

  25. 25.

    Zuker M. Mfold web server for nucleic acid folding and hybridization prediction. Nucleic Acids Res. 2003;31(13):3406–15.

    CAS  Article  Google Scholar 

  26. 26.

    Sarov-Blat L, So WV, Liu L, Rosbash M. The Drosophila takeout gene is a novel molecular link between circadian rhythms and feeding behavior. Cell. 2000;101(6):647–56.

    CAS  Article  Google Scholar 

  27. 27.

    Andrade Lopez JM, Lanno SM, Auerbach JM, Moskowitz EC, Sligar LA, Wittkopp PJ, Coolon JD. Genetic basis of octanoic acid resistance in Drosophila sechellia: functional analysis of a fine-mapped region. Mol Ecol. 2017;26(4):1148–60.

    CAS  Article  Google Scholar 

  28. 28.

    Smith CR, Morandin C, Noureddine M, Pant S. Conserved roles of Osiris genes in insect development, polymorphism, and protection. J Evol Biol. 2018;31(4):516–29.

    CAS  Article  Google Scholar 

  29. 29.

    Kanamori Y, Saito A, Hagiwara-Komoda Y, Tanaka D, Mitsumasu K, Kikuta S, Watanabe M, Cornette R, Kikawada T, Okuda T. The trehalose transporter 1 gene sequence is conserved in insects and encodes proteins with different kinetic properties involved in trehalose import into peripheral tissues. Insect Biochem Mol Biol. 2010;40(1):30–7.

    CAS  Article  Google Scholar 

  30. 30.

    Fernandez-Sanchez E, Diez-Guerra FJ, Cubelos B, Gimenez C, Zafra F. Mechanisms of endoplasmic-reticulum export of glycine transporter-1 (GLYT1). Biochem J. 2008;409(3):669–81.

    CAS  Article  Google Scholar 

  31. 31.

    Miller-Fleming L, Olin-Sandoval V, Campbell K, Ralser M. Remaining mysteries of molecular biology: the role of polyamines in the cell. J Mol Biol. 2015;427(21):3389–406.

    CAS  Article  Google Scholar 

  32. 32.

    Eckert D, Buhl S, Weber S, Jager R, Schorle H. The AP-2 family of transcription factors. Genome Biol. 2005;6(13):246.

    Article  Google Scholar 

  33. 33.

    Bauer R, McGuffin ME, Mattox W, Tainsky MA. Cloning and characterization of the Drosophila homologue of the AP-2 transcription factor. Oncogene. 1998;17(15):1911–22.

    CAS  Article  Google Scholar 

  34. 34.

    Monge I, Mitchell PJ. DAP-2, the Drosophila homolog of transcription factor AP-2. Mech Dev. 1998;76(1–2):191–5.

    CAS  Article  Google Scholar 

  35. 35.

    Kerber B, Monge I, Mueller M, Mitchell PJ, Cohen SM. The AP-2 transcription factor is required for joint formation and cell survival in Drosophila leg development. Development. 2001;128(8):1231–8.

    CAS  PubMed  Google Scholar 

  36. 36.

    Monge I, Krishnamurthy R, Sims D, Hirth F, Spengler M, Kammermeier L, Reichert H, Mitchell PJ. Drosophila transcription factor AP-2 in proboscis, leg and brain central complex development. Development. 2001;128(8):1239–52.

    CAS  PubMed  Google Scholar 

  37. 37.

    Visentin M, Diop-Bove N, Zhao R, Goldman ID. The intestinal absorption of folates. Annu Rev Physiol. 2014;76:251–74.

    CAS  Article  Google Scholar 

  38. 38.

    Zhao R, Goldman ID. The proton-coupled folate transporter: physiological and pharmacological roles. Curr Opin Pharmacol. 2013;13(6):875–80.

    CAS  Article  Google Scholar 

  39. 39.

    Tan SM, Lieberman J. Capture and identification of miRNA targets by biotin Pulldown and RNA-seq. Methods Mol Biol. 2016;1358:211–28.

    CAS  Article  Google Scholar 

  40. 40.

    Thomson DW, Bracken CP, Goodall GJ. Experimental strategies for microRNA target identification. Nucleic Acids Res. 2011;39(16):6845–53.

    CAS  Article  Google Scholar 

  41. 41.

    Soni K, Choudhary A, Patowary A, Singh AR, Bhatia S, Sivasubbu S, Chandrasekaran S, Pillai B. miR-34 is maternally inherited in Drosophila melanogaster and Danio rerio. Nucleic Acids Res. 2013;41(8):4470–80.

    CAS  Article  Google Scholar 

  42. 42.

    Wang Z, Gerstein M, Snyder M. RNA-Seq: a revolutionary tool for transcriptomics. Nat Rev Genet. 2009;10:57–63.

    CAS  Article  Google Scholar 

  43. 43.

    Lamarre S, Frasse P, Zouine M, Labourdette D, Sainderichin E, Hu G, Le Berre-Anton V, Bouzayen M, Maza E. Optimization of an RNA-Seq differential gene expression analysis depending on biological replicate number and library size. Front Plant Sci. 2018;9:108.

    Article  Google Scholar 

  44. 44.

    Garcia-Segura L, Abreu-Goodger C, Hernandez-Mendoza A, Dimitrova Dinkova TD, Padilla-Noriega L, Perez-Andrade ME, Miranda-Rios J. High-throughput profiling of Caenorhabditis elegans starvation-responsive microRNAs. PLoS One. 2015;10(11):e0142262.

    Article  Google Scholar 

  45. 45.

    Lyons PJ, Storey KB, Morin P Jr. Expression of miRNAs in response to freezing and anoxia stresses in the freeze tolerant fly Eurosta solidaginis. Cryobiology. 2015;71(1):97–102.

    CAS  Article  Google Scholar 

  46. 46.

    Campos C, Sundaram AY, Valente LM, Conceicao LE, Engrola S, Fernandes JM. Thermal plasticity of the miRNA transcriptome during Senegalese sole development. BMC Genomics. 2014;15:525.

    Article  Google Scholar 

  47. 47.

    Yin F, Qin C, Gao J, Liu M, Luo X, Zhang W, Liu H, Liao X, Shen Y, Mao L, et al. Genome-wide identification and analysis of drought-responsive genes and microRNAs in tobacco. Int J Mol Sci. 2015;16(3):5714–40.

    CAS  Article  Google Scholar 

  48. 48.

    Voskarides K. Plasticity vs mutation. The role of microRNAs in human adaptation. Mech Ageing Dev. 2017;163:36–9.

    CAS  Article  Google Scholar 

  49. 49.

    Rodgers AB, Morgan CP, Leu NA, Bale TL. Transgenerational epigenetic programming via sperm microRNA recapitulates effects of paternal stress. Proc Natl Acad Sci U S A. 2015;112(44):13699–704.

    CAS  Article  Google Scholar 

  50. 50.

    Rassoulzadegan M, Grandjean V, Gounon P, Vincent S, Gillot I, Cuzin F. RNA-mediated non-mendelian inheritance of an epigenetic change in the mouse. Nature. 2006;441:469–74.

    CAS  Article  Google Scholar 

  51. 51.

    Zhang YL, Huang QX, Yin GH, Lee S, Jia RZ, Liu ZX, Yu NT, Pennerman KK, Chen X, Guo AP. Identification of microRNAs by small RNA deep sequencing for synthetic microRNA mimics to control Spodoptera exigua. Gene. 2015;557(2):215–21.

    CAS  Article  Google Scholar 

  52. 52.

    Nandety RS, Kuo YW, Nouri S, Falk BW. Emerging strategies for RNA interference (RNAi) applications in insects. Bioengineered. 2015;6(1):8–19.

    CAS  Article  Google Scholar 

  53. 53.

    Martin M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnetjournal. 2011;17(1):10–2.

    Google Scholar 

  54. 54.

    Altschul SF, Madden TL, Schaffer AA, Zhang J, Zhang Z, Miller W, Lipman DJ. Gapped BLAST and PSI-BLAST: a new generation of protein database search programs. Nucleic Acids Res. 1997;25(17):3389–402.

    CAS  Article  Google Scholar 

  55. 55.

    Camacho C, Coulouris G, Avagyan V, Ma N, Papadopoulos J, Bealer K, Madden TL. BLAST+: architecture and applications. BMC Bioinformatics. 2009;10:421.

    Article  Google Scholar 

  56. 56.

    Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol. 2010;11:R106.

    CAS  Article  Google Scholar 

  57. 57.

    Pfaffl MW, Horgan GW, Dempfle L. Relative expression software tool (REST) for group-wise comparison and statistical analysis of relative expression results in real-time PCR. Nucleic Acids Res. 2002;30(9):e36.

    Article  Google Scholar 

  58. 58.

    Lewis BP, Shih IH, Jones-Rhoades MW, Bartel DP, Burge CB. Prediction of mammalian microRNA targets. Cell. 2003;115(7):787–98.

    CAS  Article  Google Scholar 

  59. 59.

    Lewis BP, Burge CB, Bartel DP. Conserved seed pairing, often flanked by adenosines, indicates that thousands of human genes are microRNA targets. Cell. 2005;120(1):15–20.

    CAS  Article  Google Scholar 

  60. 60.

    TargetScanHuman Prediction of microRNA targets. https://genesmitedu/targetscan.

  61. 61.

    Hofacker IL. Vienna RNA secondary structure server. Nucleic Acids Res. 2003;31(13):3429–31.

    CAS  Article  Google Scholar 

  62. 62.

    BIPAA Bioinformatics Platform for Agroecosystem Arthropods. In: https://bipaagenouestorg/is/lepidodb/spodoptera_frugiperda/

  63. 63.

    Enright AJ, John B, Gaul U, Tuschl T, Sander C, Marks DS. MicroRNA targets in Drosophila. Genome Biol. 2003;5(1):R1.

    Article  Google Scholar 

  64. 64.

    John B, Enright AJ, Aravin A, Tuschl T, Sander C, Marks DS. Human MicroRNA targets. PLoS Biol. 2004;2(11):e363.

    Article  Google Scholar 

  65. 65.

    Miranda KC, Huynh T, Tay Y, Ang YS, Tam WL, Thomson AM, Lim B, Rigoutsos I. A pattern-based method for the identification of MicroRNA binding sites and their corresponding heteroduplexes. Cell. 2006;126(6):1203–17.

    CAS  Article  Google Scholar 

  66. 66.

    Vejnar CE, Blum M, Zdobnov EM. miRmap web: comprehensive microRNA target prediction online. Nucleic Acids Res. 2013;41(Web Server issue):W165–8.

    Article  Google Scholar 

  67. 67.

    MiRmap. In. https://mirmap.ezlab.org.

Download references

Acknowledgments

We are grateful to Clotilde Gibard and Gaëtan Clabots for maintaining the insect collections of the DGIMI laboratory in Montpellier. We thank Clotilde Gibard for her help to rear Spodoptera frugiperda on whole plants during reciprocal transplant experiments.

Funding

This work was supported by a grant from the French National Research Agency (ANR-12-BSV7–0004-01; http://www.agence-nationale-recherche.fr/) for E. A. including a post-doctoral fellowship for Y.M. H.P. acknowledges financial support from the “France Génomique” National infrastructure, funded as part of “Investissement d’Avenir”(ANR-10-INBS-09).The funding bodies had no role in the design of the study and collection, analysis, and interpretation of data and in writing the manuscript.

Availability of data and materials

The datasets generated during this study are available from NCBI SRA (SRP132165).

https://www.ncbi.nlm.nih.gov/Traces/study/?acc=SRP132165

https://www.ncbi.nlm.nih.gov/sra/?term=SRP132165

The Genbank accession numbers of microRNA precursors identified in this study are MH340548 - MH341001.

Author information

Author notes

  1. Yves Moné, Sandra Nhim and Sylvie Gimenez contributed equally to this work.

Affiliations

  1. DGIMI, Univ Montpellier, INRA, Montpellier, France

    Yves Moné, Sandra Nhim, Sylvie Gimenez, Imène Seninet, Nicolas Nègre & Emmanuelle d’Alençon

  2. INRA, UMR Institut de Génétique, Environnement et Protection des Plantes (IGEPP), BioInformatics Platform for Agroecosystems Arthropods (BIPAA), Campus Beaulieu, Rennes, France

    Fabrice Legeai

  3. INRIA, IRISA, Genscale, Campus Beaulieu, Rennes, France

    Fabrice Legeai

  4. MGX, CNRS, INSERM, Univ Montpellier, Montpellier, France

    Hugues Parrinello

Authors
  1. Yves Moné
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Sandra Nhim
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Sylvie Gimenez
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Fabrice Legeai
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Imène Seninet
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Hugues Parrinello
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Nicolas Nègre
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Emmanuelle d’Alençon
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

Conception and design: YM, NN, EA Acquisition of data, analysis and interpretation: YM, SN, SG, FL, IS, HP Drafting and revision of the manuscript: YM, NN, EA. All authors have read and approved the manuscript.

Corresponding author

Correspondence to Emmanuelle d’Alençon.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1:

Table S1. Number of sequence reads in each small non coding RNAs library. (DOCX 14 kb)

Additional file 2:

Excel file S1. Predictions of miR genes by mirDEEP2 in SfC and SfR and orthology table. (XLSX 92 kb)

Additional file 3:

Figure S1. MA-plots showing the relative expression of known or novel miR according to the host-plant (Rice compared to corn) in each strain. Top panel, in SfC, bottom panel In SfR. (PPTX 5002 kb)

Additional file 4:

Figure S2. Variation between samples (treatments, replicates) of larvae exposed to different plants displayed by Principal Component Analysis (PCA). Top panel: SfC, bottom panel: SfR. (PPTX 4001 kb)

Additional file 5:

Excel file S2. Relative expression resulting from DESEQ2 analysis within strains according to the host-plant. In red novel miR, in black known miRs. (XLSX 53 kb)

Additional file 6:

Figure S3. MA-plots showing the relative expression of known or novel miR according to the genetic background. Top panel, relative expression analyzed by DESEQ2 in SfR compared to SfC on corn, bottom panel, relative expression in SfR compared to SfC on rice. (PPTX 4997 kb)

Additional file 7:

Figure S4. Variation between samples (treatments, replicates) of larvae exposed to different plants displayed by Principal Component Analysis (PCA). Top panel: SfR compared to SfC on corn, bottom panel: SfR compared to SfC on rice. (PPTX 3008 kb)

Additional file 8:

Excel file S3. Relative expression resulting from DESEQ2 analysis according to the genetic background on the same host-plant. In red novel miR, in black known miRs. (XLSX 39 kb)

Additional file 9:

Excel file S4. List of gene target predictions of known differentially expressed miRs. (XLSX 1578 kb)

Additional file 10:

Excel file S5. Target genes of known miRs and their relative expression in SfC on rice compared to corn. (XLSX 35 kb)

Additional file 11:

Excel file S6. Target genes of known miRs and their relative expression in SfR compared to SfC on corn. (XLSX 19 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. 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.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Moné, Y., Nhim, S., Gimenez, S. et al. Characterization and expression profiling of microRNAs in response to plant feeding in two host-plant strains of the lepidopteran pest Spodoptera frugiperda. BMC Genomics 19, 804 (2018). https://doi.org/10.1186/s12864-018-5119-6

Download citation

Keywords

  • Adaptation
  • Phenotypic plasticity
  • Regulation of gene expression
  • Insect plant interaction
  • microRNAs

BMC Genomics

ISSN: 1471-2164

Contact us