MicroRNAs in metamorphic and non-metamorphic transitions in hemimetabolan insect metamorphosis
© Rubio et al.; licensee BioMed Central Ltd. 2012
Received: 6 June 2012
Accepted: 26 July 2012
Published: 10 August 2012
Skip to main content
© Rubio et al.; licensee BioMed Central Ltd. 2012
Received: 6 June 2012
Accepted: 26 July 2012
Published: 10 August 2012
Previous work showed that miRNAs play key roles in the regulation of metamorphosis in the hemimetabolan species Blattella germanica. To gain insight about which miRNAs might be important, we have constructed two miRNA libraries, one of the penultimate, pre-metamorphic nymphal instar (N5) and the other of the last, metamorphic nymphal instar (N6).
High throughput sequencing gave 61 canonical miRNAs present in the N5 and N6 libraries, although at different proportions in each. Comparison of both libraries led to the identification of three and 37 miRNAs significantly more expressed in N5 and N6 respectively. Twelve of these 40 miRNAs were then investigated further by qRT-PCR and results indicated that miR-252-3p was well expressed in N5 but not in N6, whereas let-7-5p, miR-100-5p and miR-125-5p showed the reverse pattern. 20-Hydroxyecdysone (20E) tended to stimulate miRNA expression, whereas juvenile hormone (JH) inhibited the 20E stimulatory effect. Expression of let-7, miR-100 and miR-125 was increased by 20E, which has also been observed in D. melanogaster. The only miRNA that was inhibited by 20E was miR-252-3p. The involvement of let-7, miR-100 and miR-125 in metamorphosis has been demonstrated in other insects. Depletion of miR-252-3p caused growth and developmental delays, which suggests that this miRNA is involved in regulating these processes prior to metamorphosis.
The comparative analysis of miRNA libraries from pre-metamorphic (N5) and metamorphic stages (N6) of B. germanica proved to be a useful tool to identify miRNAs with roles in hemimetabolan metamorphosis. Three miRNAs emerged as important factors in the metamorphic stage (N6): let-7-5p, miR-100-5p and miR-125-5p, whereas miR-252-3p appears to be important in the pre-metamorphic stage (N5).
Metamorphosis has been a key innovation in insect evolution, which had a decisive influence on the dramatic diversification of this animal class, with around one million presently described species . There are two types of metamorphosis: hemimetaboly, or gradual metamorphosis, and holometaboly or abrupt metamorphosis , the latter being evolutionarily the most successful . In all cases, insect metamorphosis is hormonally regulated, and the most important hormones are the ecdysteroid 20-hydroxyecdysone (20E), which determines the successive moults, and the terpenoid juvenile hormone (JH), which represses the expression of adult features. Each hormone exerts its respective action through a cascade of transcription factors that transduce the hormonal signal to the effector genes [3, 4].
Recently, a number of reports have pointed to microRNAs (miRNAs) as important players in the metamorphic transition, in hemimetabolan as well as in holometabolan species . miRNAs are small non-coding RNAs of about 21–22 nucleotides that modulate gene expression at a post-transcriptional level, often in the context of developmental and morphogenetic processes [6–8]. They undergo molecular processing before becoming mature miRNAs [9, 10]. First, miRNAs are transcribed as part of a primary transcript (pri-miRNA), which contains one or more miRNA precursors (pre-miRNAs). In the nucleus, the pri-miRNAs are processed into hairpin-structured pre-miRNAs by the ribonuclease Drosha, and exported to the cytoplasm. Once in the cytoplasm, the pre-miRNAs are cleaved by the ribonuclease Dicer-1 into an imperfectly-paired duplex, whose 5’- and 3’-strands can either give two respective miRNAs or only one. Generally, only one of the strands give rise to a mature miRNA, whereas the other, called miRNA* or the star strand, is degraded . However, given that both strands can give mature miRNAs in some cases , the miRNAs are presently labelled as miRNA-5p or miRNA-3p according to the position of the original strand in the pre-miRNA. This is the convention used in the present paper.
One of the more dramatic demonstrations that miRNAs are involved in insect metamorphosis was reported by Gomez-Orte & Belles , who silenced dicer-1 expression by RNAi in the last nymphal instar of the cockroach Blattella germanica, and obtained supernumerary nymphs after the following moult, instead of adults . These results showed that Dicer-1 and miRNAs are important factors in the mechanisms regulating metamorphosis. However, which particular miRNAs were involved remained unclear.
B. germanica is a polyneopteran exopterigote insect that shows a gradual morphological transformation during the life cycle, and is a good representative of the less modified hemimetaboly. This, and the fact that it is especially sensitive to gene silencing by RNAi , have made B. germanica a favourite hemimetabolan model to study metamorphosis. However, the genome of B. germanica is yet to be sequenced. The first step to approach the study of miRNAs involved in metamorphosis was therefore to establish a baseline miRNA catalogue using high-throughput sequencing . To gain more information about which miRNAs may be important for B. germanica metamorphosis, we constructed two miRNA libraries, one with RNA extracted around the peak of 20E of the penultimate nymphal instar (N5), and the other constructed equivalently but during the last instar (N6). The endocrine context differs between these two libraries with the presence of JH in N5 and its absence in N6. We therefore aimed to identify miRNAs differentially expressed in N5 and N6. Those that were identified with this approach were further studied by establishing the expression pattern in N5 and N6, and by investigating the effect of 20E alone and that of 20E plus JH in their expression. Finally, we selected miR-252-3p, the most interesting miRNA emerging from these studies, to carry out a functional analysis using a specific anti-miRNA.
Specimens of B. germanica were obtained from a colony reared in the dark at 30 ± 1°C and 60–70% relative humidity. Freshly ecdysed nymphs or adult females were selected and used at the appropriate ages. All dissections, treatments and tissue sampling were carried out on carbon dioxide-anaesthetized specimens. Samples for RNA extraction and quantification were frozen immediately after dissection, and stored at −80°C until use.
To construct the two miRNA libraries, we carried out total RNA extraction from the whole body (excluding the head and the digestive tube, to avoid ocular pigments and intestine parasites) of B. germanica specimens, using the miRNeasy extraction kit (QIAGEN). Nine individuals taken in days 3, 4 and 5 (three each day) of fifth nymphal instar were pooled for the N5 library. The equivalent sampling was carried out on days 5, 6 and 7 of sixth (last) nymphal instar for the N6 library. The days chosen correspond to the 20E respective moulting peaks in N5 and N6 of B. germanica.
Solexa technology (using an Illumina HiSeq 2000 sequencer) was used for small RNA sequencing in B. germanica. Sequences were received in FASTQ format, with the adaptor sequences trimmed. The first filter applied was based on quality values; considering that Illumina quality values range from 0 to 40 (using ASCII 64–104 in fastq) , we eliminated those reads having 80% of the sequence with values lower than 20. Then we applied a trimming process to delete regions with repetitive nucleotides or short motifs, and low quality regions, either at the beginning or at the end of the sequence. Sequences with the same length were grouped, then we obtained the counts collapsing the sequences using the suite of programs Fastx-toolkit (http://hannonlab.cshl.edu/fastx_toolkit), and those sequences with 10–28 nucleotides were kept. To eliminate artefacts, we mapped the sequences against Pfam database (excluding miRNAs), and against the genome sequences of Blattabacterium cuenoti, a bacteroid endosimbiont of B. germanica, B. germanica densovirus  and Escherichia coli. The remaining sequences were considered putative miRNAs. We identified the canonical miRNAs by Bowtie alignment  against mature miRNAs of the mirBASE database (November 2011 version) . We followed the -3p or -5p nomenclature, instead of mature and star strand, assuming that the miRNA position in the precursor of B. germanica was conserved with respect to D. melanogaster. The seed was defined as the first 17 nucleotides and not mismatches were allowed . Sequences shorter than 17 nucleotides were labelled as a given miRNA when they had 100% matching with the corresponding region of the miRNA seed. For expression analysis, different miRNAs length variants, or isoMIRs  of the same miRNA were consider together. Finally, the poorest represented miRNAs, i.e., the pool of those representing only 0.01% of the total reads, were eliminated from the differential expression analysis.
In order to compare the N5 and N6 libraries, the counts were normalized by transforming them into RPKM (Reads Per Kilobase Mapped)  and were analysed by the NOISeq R package . The counts of the different miRNA length variants were merged. The method of NOISeq applied was NOISeqsim, which assumes that there are not replicates for the experimental conditions and it simulates them. We used the default values to define those miRNAs differentially expressed (probability > 0.8). The evaluation of the results was based on the NOISeq scores: probability of being differentially expressed, M-value and D-value. The results were represented by a Volcano plot, as described by Cui and Churchill .
miRNA expression patterns were established at 48 h intervals in female specimens of N5 and N6 and in freshly emerged adults. Three biological replicates of each age were studied. Whole body sampling and total RNA extraction was carried out as in the library construction. Quantification of miRNA levels was carried out by quantitative real-time PCR (qRT-PCR). An amount of 400 ng of total RNA was reverse transcribed with the NCode miRNA first-strand synthesis and qRT-PCR kit (Invitrogen), following the manufacturer protocol. Amplification reactions were carried out using IQTM SYBR Green Supermix (BioRad) and the following protocol: 95°C for 2 min, and 40 cycles of 95°C for 15 s and 60°C for 30 s, in a MyIQ Real-Time PCR Detection System (BioRad). A dissociation curve was obtained to ensure that there was only one product amplified after the amplification phase. All reactions were run in triplicate. The endogenous reference gene used was U6, as in ; therefore, results are given as copies of RNA per 1000 copies of U6. Primer sequences for each miRNA are available on request.
In a first set of experiments, specimens were treated as freshly emerged sixth instar nymph (N6), by injecting 1 μL of 20E (Sigma) at a concentration of 1 μg/μL dissolved in water with 10% ethanol; controls where treated equivalently with 1 μL of water with 10% ethanol. Injection was performed through the membrane between fourth and fifth ventrites. In another experimental set, specimens were treated at the same stage, by applying topically 1 μL of JH III (Sigma) dissolved in acetone at a concentration of 2 μg/μL, and with 20E in the same conditions as in the previous set of experiments; controls were treated equivalently with 1 μL of acetone and 1 μL of water with 10% ethanol. JH III is the native JH of B. germanica, and given that the commercial JH III is a mixture of isomers containing ca. 50% of the biologically active (10R)-JH III, thus the active dose applied was around 1 μg per specimen. Total RNA extraction and quantification of dicer-1 expression and miRNAs were performed 24 h after the treatment following the procedures described above, and using three biological replicates. Statistical analysis of relative expression measurements was carried out with the REST software tool .
To deplete miR-252-3p expression we used miRCURY LNA™ microRNA Power Inhibitor (Exiqon). The sequence of LNA anti-miR-252-3p (GATAAGCACTTGAGCAGCAGG) is quite similar to the sequence of miR-252-5p (CTAAGTACTAGTGCCGCAGGA), which suggests that it might act as miR-252-5p mimic. However, the seed of LNA-antimiR-252-3p (ATAAGC) is very different to that of miR-252-5p (TAAGTA), which makes very improbable this possibility. A dose of 1 μL of LNA at 50 μM was injected in the abdominal cavity 24 h after the moult to the fifth instar nymph (N5), and a second injection at the same concentration was carried out 48 h later. Preliminary experiments administering a single dose 24 h after the moult to N5 depleted miR-252-3p less efficiently. Controls were injected with miRCURY LNA™ microRNA Inhibitor Negative Control A (Exiqon) at the same stages and concentrations. Specimens for miRNA quantification were collected 24 h after the second injection. RNA extraction, quantification of miRNA levels and statistical analysis were performed as described above. Specimens to study the phenotype were left alive until the adult moult.
Counts and RPKM of the canonical miRNAs found in the N5 and N6 libraries
miRNAs differentially expressed in N5 sorted by the absolute M-value
miRNAs differentially expressed in N6, sorted by the M-value
Notably, both strands of the miRNA precursor miR-252 are well present in the two libraries, suggesting that both of them are biologically active miRNAs. The 5’ strand (miR-252-5p) is about 12-fold more abundant than the 3’ strand (miR-252-3p) in the N5 library and ca. 800-fold times more abundant in N6. Furthermore, from N5 to N6 there is an increase of miR-252-5p and a significant decrease of miR-252-3p (Table 1). miR-276 represents a similar case, with both strands represented in both libraries. In this case, however, the 3’ strand (miR-276-3p) is the predominant one in both libraries, being clearly more highly expressed in N6 with respect to N5, whereas the 5’ strand (miR-276-5p) is more highly expressed in the N5 with respect to N6.
For expression pattern studies and hormonal treatments, we selected the three miRNAs more highly expressed (M-value between 1.6 and 5.3) in the N5 library (miR-252-3p, miR-276-5p and miR-190-5p) and the four more differentially expressed (M-value between −1.4 and −2.0) in the N6 library (bantam-3p, miR-100-5p, miR-125-5p and let-7-5p). We also added three miRNAs that are similarly represented in N5 and N6: miR-252-5p and miR-276-3p, whose respective partner strands were chosen due to their higher expression in the N5 library, and miR-1-3p, which was the most abundant miRNA (accounting for more than 50% of total reads). Finally, we included two additional miRNAs (miR-14-3p and miR-34-5p) because of their potential interest in the context of moulting and metamorphosis. miR-14-3p has been described targeting the ecdysone receptor (EcR) in D. melanogaster, and miR-34-5p is inducible by the JH analogue methoprene in Drosophila S2 cells .
The patterns of the remaining miRNAs examined show two bursts of expression. In the cases of miR-276-5p, miR-190-5p, miR-252-5p and bantam-3p, there were expression bursts seen both in N5 and in N6, approximately corresponding to the peaks of 20E (Figure 2). These data do not correlate with quantitative sequencing data, which predicted that miR-276-5p and miR-190-5p should have shown much higher expression levels around the 20E peak of N5, whereas those of miR-252-5p and bantam-3p should have been much lower at this stage. Two expression bursts of miR-1-3p, miR-34-5p and mir-276-3p occurred in N6, the first just after emergence and the other coinciding or close to the peak of 20E of this stage (Figure 2). Comparing the levels of expression around the peaks of 20E in N5 and N6, these approximately patterns correspond to what would be expected based on the sequencing data. Minor discrepancies, like in the case of miR-276-5p, miR-276-3p, miR-190-5p, can be due to individual variability, given that sequencing data come from a pool of nine specimens, whereas data of expression patterns is based on three particular specimens.
On the other hand, the endocrine situation of a nymph-nymph transition at the peak of 20E (like in N5) was mimicked by treating equivalently staged specimens with 2 μg of JH and 1 μg of 20E. Results (Figure 3B) indicate that all the miRNAs studied, except miR-1-3p and miR-34-5p, had expression levels significantly lower than those of the controls, although the expression of the latter miRNA showed a tendency to be inhibited by the double hormonal treatment. The expression of dicer-1 was not significantly affected by the treatment (Figure 3B).
High throughput sequencing of whole body miRNAs around the peak of 20E in the N5 and N6 libraries of B. germanica gave a total of 61 canonical miRNAs that were present in both libraries, although in different proportions in each. In a previous whole body miRNA library prepared during N6, but representing all eight days of the stage, we obtained a total of 49 canonical miRNAs . Therefore, the number of miRNAs recovered in the present study, using animals collected around the peak of 20E, is higher. This suggests that most miRNAs might be more efficiently expressed around the peak of 20E, and that there was a dilution effect in the former library, which had all the days of the instar represented. Practically all canonical miRNAs reported by Cristino et al.  are present (except miR-1-5p and miR-184-5p), and are the most abundant in the present libraries. Some of the miRNAs (miR-8, miR-9a, miR-10, miR-71, miR-252, miR-276, miR-281) are represented by both strands, -3p and -5p, as was also observed in the previous N6 library .
The comparative analysis of N5 and N6 libraries allowed the identification of three and 37 miRNAs significantly more expressed in N5 and N6 respectively (Figure 1). Subsequently, 12 of these 40 miRNAs were investigated further by qRT-PCR to establish their expression patterns during N5 and N6. Results indicate that miR-252-3p is the only clear representative of miRNA that is well expressed in N5 and with decreased expression in N6, whereas let-7-5p, miR-100-5p and miR-125-5p are clear representatives of the reverse pattern (Figure 2). Whereas the expression pattern of miR-252-3p had yet to be studied in other insects, those patterns of let-7, miR-100 and miR-125 have been thoroughly examined. In D. melanogaster, these three miRNAs cluster in the same primary transcript, as first revealed by Bashirullah et al.  and Sempere et al. , and later confirmed in this fly and other insects [5, 32]. Thus, the members of this miRNA cluster are expressed simultaneously, and studies in the holometabolan species, like D. melanogaster and B. mori have shown that expression concentrates in pre-metamorphic stages. In D. melanogaster, expression starts in late third (last) instar larvae, around the peak of 20E that triggers puparium formation, and continues until the imaginal moult [30, 33]. In B. mori, let-7 expression starts in the moult leading to the penultimate larval instar, and continues until the emergence of the adult . In the hemimetabolan B. germanica, the present study shows that expression of let-7, miR-100 and miR-125 concentrates in the last nymphal instar, having an expression burst around the peak of 20E that triggers metamorphosis (Figure 2).
Of note, among the miRNAs selected to study the expression pattern, we included miR-252 and miR-276, which were represented by both the -3p and -5p strands. In principle, the two respective strands originate from the same precursor and both should presumably be expressed at similar rates through parallel patterns. But this was not found to be the case. In the case of miR-252, the -3p strand is expressed at much lower levels than those of the -5p strand, whereas in the case of miR-276, the situation is the reverse; therefore, in each case the pattern of expression of one strand is clearly different to the expression of the other (Figure 2). These observations suggest that the expression is regulated differently in both strands, and that they play distinct miRNA functions . In support of this notion, Liu et al.  reported that both strands of miR-276 are differentially expressed in different tissues and stages of B. mori, although miR-276-3p is always the most highly expressed strand.
Among the patterns examined, most of the expression peaks, either in N5 or N6, occur close to a peak of 20E (Figure 2), which suggests cause-effect relationships. This possible hormonal regulation of miRNA expression was therefore tested experimentally. In general, the experiments of hormonal treatments suggested that 20E alone tends to stimulate miRNA expression, whereas JH attenuates, or even inhibits the 20E stimulatory effect (Figure 3). Among those miRNAs that increased in expression after 20E treatment, we found let-7, miR-100 and miR-125. The increase in expression of miR-100-5p after 20E treatment is statistically significant, whereas let-7-5p and miR-125-5p had a tendency to increase (Figure 3A). This is consistent with previous studies in D. melanogaster that revealed that 20E increases the levels of these miRNAs [30, 36] by inducing the expression of the polycistronic primary transcript that contains the let-7, miR-100 and miR-125 precursors . In D. melanogaster, the burst of expression of let-7-5p in prepupae is important to ensure appropriate remodelling of the abdominal neuromusculature during metamorphosis . Moreover, an independent work of Caygill and Johnston  showed that during D. melanogaster metamorphosis, the absence of let-7-5p and miR-125-5p results in temporal delays in the terminal cell-cycle exit in the wing, and in the maturation of neuromuscular junctions of imaginal abdominal muscles. The authors focused on the latter process by identifying the abrupt (ab) gene as a let-7-5p target, and by providing evidence that showed that let-7-5p regulates the maturation rate of abdominal neuromuscular junctions during metamorphosis by modulating ab expression .
Conversely, the significant stimulatory effect of 20E on miR-1-3p and the tendency of 20E to stimulate miR-34-5p that we observed in B. germanica contrasts with the results obtained in D. melanogaster, where 20E did not affect miR-1-3p and inhibited miR-34-5p . We interpret these differences as reflections of different modes of regulation of these miRNAs in the two model species. The differences are also illustrated by the results obtained after JH treatment of B. germanica, which abolished the stimulatory effects of 20E on practically all miRNAs (Figure 3B), which was also the case for miR-34-5p. This contrasts again with the observations of Sempere et al. , which suggested that JH has a stimulatory effect on this miRNA in D. melanogaster.
Possibly the most interesting result of the experiments of hormonal treatments is that of miR-252-3p, which was the only miRNA that was inhibited by 20E alone (Figure 3A). Moreover, JH did not seem to counteract the 20E effect, as compared to other miRNAs (Figure 3B). Inhibition by 20E suggests that the decrease of miR-252-3p in N6, and the increase interrupted by a plateau around the 20E peak in N5 (Figure 2), might be due to the action of 20E. The other strand, classically considered the “mature” strand of miR-252, miR-252-5p, has been found in a number of insects , and its expression has been carefully studied in B. mori, where it shows a continuous high expression from the spinning larvae to pupal and adult stages . Conversely, miR-252-3p has never been studied, nor has it been related with insect moulting and metamorphosis. It was for these reasons that we studied it functionally.
Depletion of miR-252-3p levels in the penultimate (N5) nymphal instar caused retarded growth and developmental delays within the instar (Figure 4). At the end, however, the treated insects were able to moult to the last (N6) nymphal instar and to the adult normally. Observed delays in N5 reached as long as 14 days in some cases, and levels of miR-252-3p measured in specimens on N5D13, which were presumably going to moult within the next one or two days, had somewhat recovered compared with levels measured on N5D6, although they were still significantly lower than levels of the miRNA of the controls (Figure 4). The data suggest that miR-252-3p regulates transcripts that are important for growth and development in N5, in the transition to the developmentally important N6 which precedes the imaginal moult. The fact that the treated specimens were finally able to moult to N6 and to the adult stage even with reduced levels of miR-252-3p, suggests that these levels, although still quite low, are above the operative threshold for this miRNA.
1. The comparative analysis of a miRNA library of the penultimate nymphal instar (N5) with an equivalent library from the last instar nymph (N6) of B. germanica revealed 40 canonical miRNAs that are differentially expressed in each library.
2. In terms of counts, three of these miRNAs were found to be significantly more expressed in N5, whereas 37 were found to be more expressed in N6. qRT-PCR studies confirmed that miR-252-3p is well expressed in N5 but not in N6, whereas let-7-5p, miR-100-5p and miR-125-5p showed the reverse pattern.
3. Hormonal treatments indicated that 20E increases the expression of let-7, miR-100 and miR-125, but inhibits that of miR-252-3p.
4. The involvement of let-7, miR-100 and miR-125 in metamorphosis has been shown in other insects. Depletion of miR-252-3p in B. germanica caused growth and developmental delays, which suggests that this miRNA is involved in regulating these processes prior to metamorphosis.
Support for this research was provided by the Spanish MICINN (grant CGL2008-03517/BOS to X.B.) and by the CSIC (grant 2010TW0019, from the Formosa program, to X.B. and a predoctoral fellowship to M.R. from the JAE program) and from the Spanish National Institute for Bioinformatics (http://www.inab.org).
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.