MicroRNAs in the oriental fruit fly, Bactrocera dorsalis: extending Drosophilid miRNA conservation to the Tephritidae
© Calla and Geib. 2015
Received: 23 April 2015
Accepted: 7 August 2015
Published: 5 October 2015
The oriental fruit fly, Bactrocera dorsalis, is an important plant pest species in the family Tephritidae. It is a phytophagous species with broad host range, and while not established in the mainland United States, is a species of great concern for introduction. Despite the vast amount of information available from the closely related model organism Drosophila melanogaster, information at the genome and transcriptome level is still very limited for this species. Small RNAs act as regulatory molecules capable of determining transcript levels in the cells. The most studied small RNAs are micro RNAs, which may impact as much as 30 % of all protein coding genes in animals.
We have sequenced small RNAs (sRNAs) from the Tephritid fruit fly, B. dorsalis (oriental fruit fly), specifically sRNAs corresponding to the 17 to 28 nucleotides long fraction of total RNA. Sequencing yielded more than 16 million reads in total. Seventy five miRNAs orthologous to known miRNAs were identified, as well as five additional novel miRNAs that might be specific to the genera, or to the Tephritid family. We constructed a gene expression profile for the identified miRNAs, and used comparative analysis with D. melanogaster to support our expression data. In addition, several miRNA clusters were identified in the genome that show conservancy with D. melanogaster. Potential targets for the identified miRNAs were also searched.
The data presented here adds to our growing pool of information concerning the genome structure and characteristics of true fruit flies. It provides a basis for comparative studies with other Dipteran and within Tephritid species, and can be used for applied research such as in the development of new control strategies based on gene silencing and transgenesis.
True fruit flies (Diptera: Tephritidae), constitute a family of mostly phytophagous species, many of which are considered to be serious pests of numerous plants. Within this family, species in the highly invasive genera Anastrepha, Bactrocera, Ceratitis and Rhagoletis are of worldwide economic importance, restricted by quarantine listings in Europe , and subject to constant eradication and establishment prevention in the United States. Research on this family of flies has been heavily concentrated on field strategies and quarantine methods . Despite the vast amount of information available from the closely related model organism Drosophila melanogaster, information at the genome and transcriptome level is still very limited in true fruit flies.
The regulation of gene expression in the cell is driven by fine-tuned molecular mechanisms that respond to developmental and environmental cues. Playing an important role in this system are the small RNAs (sRNA) which act as regulatory molecules capable of determining target mRNA expression levels [3, 4]. Several types of sRNAs have been documented to date. Among the most studied are micro RNAs (miRNAs) and small interfering RNAs (siRNAs); both types present in most species. Classes of sRNAs and their biogenesis pathways have been extensively described in the literature for many model species [3–6]. MiRNAs originate in genomic loci and are very often expressed in a tissue-specific manner. A large proportion, probably more than 30 %, of all protein coding genes of animals may be regulated by miRNAs.
We have sequenced sRNAs from the Tephritid fruit fly, Bactrocera dorsalis (Oriental fruit fly), specifically sRNAs corresponding to the 17 to 28 nucleotides long fraction of total RNA, looking at variation in composition and expression between different developmental stages, and between male and female sexes in the pupa stage. We were able to identify several miRNAs orthologous to known miRNAs and additional novel miRNAs that might be specific to the genera or to the Tephritid family. We constructed a profile of gene expression for the identified miRNAs, and used comparative analysis with D. melanogaster to support our expression data, identify conserved miRNA clusters in the genome, and mine for potential transcript targets for this miRNAs. The data presented here adds to the biological information concerning the genome structure and characteristics of true fruit flies. It provides a basis for comparative studies in other Dipteran species, and can be used for applied research such as in the development of new control strategies based on gene silencing and transgenesis.
Fly sample collection
A white pupal translocated strain (DTWP) was used, in which female pupae are white and male pupae are brown, allowing for separation of sex in the pupal stage . Flies were grown in liquid diet , (200 ml diet for approximately 3400 eggs) as previously described .
Three biological replications were carried out for sample collection, RNA extraction and sequencing. For each replicate, eggs were allowed to develop and samples collected at the following times and developmental stages: embryos (12 mg at 0–1 h after oviposition), young larvae (approximately 20 mg at 0–12 h after egg hatch), early male (brown) pupae (0–24 h after pre-pupal formation), and early female (white) pupae (0–24 h after pre-pupal formation). For embryos and larvae, samples were sieve washed and then blot dried. All samples were collected in 1.5 ml microcentrifuge tubes and flash frozen in liquid nitrogen immediately after collection. Samples were then stored at −80 °C until processed. For fertilized ovary collection, 7 day old females and males were left in a cage to mate. Thirty seven mating pairs were separated in cups and left for at least 90 min to make sure that the females were fertilized. Non-mating flies were removed and mating pairs were left until the next day. Mated females were sedated by exposing them to 4 °C for 10 min, and the ovaries were dissected. Ovaries from ten female flies were collected per replication in 1.5 ml microcentrifuge tubes and flash frozen in liquid nitrogen. Samples were preserved at −80 °C until RNA extraction.
RNA from each of the collected samples was extracted utilizing NucleoSpin® miRNA kit (Macherey-Nagel, Duren, Germany) following manufacturer’s protocol and recommendations. Initial tissue lysis was performed by grinding the frozen tissue in 1.5 ml tubes with plastic micro pestles, followed by addition of 300ul of lysis buffer. The NucleoSpin miRNA kit allows for separation of small RNA and large RNA fractions in silica membrane columns by differential ethanol concentrations. After purification, the quality and quantity of both the small and large RNA fractions for each sample was determined using a Quibit 2.0 fluorometer (Life Technologies, Carlsbad, California), and an Agilent 2100 Bioanalyzer with anAgilent small RNA kit (Ambion, Santa Clara, CA, USA).
Library preparation and sequencing
To prepare small RNA sequencing libraries, the Ion Total RNA-seq kit v2 for small RNA libraries was used following manufacturer protocols with some modifications. The small RNA fraction of each sample was ligated to adapters and reverse transcribed to cDNA. The cDNA was purified, size selected and each sample was differentially barcoded before amplification to allow subsequent sample identification. Amplified cDNA was checked for quality and size distribution using an Agilent 2100 Bioanalyzer with the High Sensitivity DNA assay kit. To further size select the cDNA, a BluePippin® (Sage Science, Beverly, MA) was used with a 3 % agarose gel cassette, enabling the enrichment of cDNAs between approximately 92 and 118 bp, (corresponding to small RNAs of 18–27 nucleotides length with the addition of sequencing adaptors and barcodes). Afterward, equimolar amounts of samples (identifiable by differential barcodes) were pooled, and the pooled library was quantified using a KAPA library quantification kit (KAPA Biosystems Woburn, MA) to assess the optimal amount for emulsion PCR. Emulsion PCR was performed on an Ion One Touch 2 System and the amplified beads were subjected to sequencing with the Ion Personal Genome Machine using an Ion PGM Sequencing 400 kit and Ion 318 Chip v2. To maximize reads per sample, sequencing runs were performed using 120 flows per run, followed by repeated sequencing of the library a total of 6 times across two initializations of the PGM instrument with the 400 kit. Post sequencing base calling, adapter trimming, and demultiplexing was performed using the Torrent Suite Software using default parameters for small RNA sequencing, and exported as fastq files.
Identification and classifications of small RNA sequences
The resulting 14 fastq files corresponding to a replicate sample from each life stage were entered separately in the mapper.pl module of the mirDeep2 package [10, 11], using a config.txt file to track the original sample in the final results. The mapper.pl module discarded reads smaller than 17 nt (option -l 17), collapsed identical reads and performed counting. Additionally, reads were mapped to the B. dorsalis genome (GeneBank Accession # JFBF00000000.1) using bowtie  with the following stringency parameters: 0 mismatches allowed in the seed, seed length of 18 nt, up to two mismatches after the seed, discarding reads mapping more than 5 times to the genome, and reporting only the best alignment for each read. The pipeline was designed to predict high-confidence miRNAs, and while discarding mature sequences that map to too many loci compromised detection, it also prevented false positive predictions. The collapsed reads obtained from the mapper.pl module were input into the miRdeep2 core module (miRdeep2.pl) with no reference miRNAs from a closely related species supplied, in this manner all the potential miRNAs and precursors from the data could be obtained. In the miRDeep2 module, read mappings were used to excise putative miRNA precursors according to stack height and all sequenced reads were aligned to these potential precursors. Additionally, the secondary hairpin structure and its stability were predicted for the excised precursors utilizing RNAfold [13–15] and Ranfold , and a score was assigned to each precursor. Scores are used to select precursors with highest probability of being genuine (the program kept precursors with scores higher than −50). To find other non-miRNA small RNAs, mapped reads (from mapper.pl module output) were scanned against the covariance models of the Rfam 11.0 release , using Infernal 1.1.1 cmscan .
The list of potential mature miRNAs and hairpin precursors obtained from the miRDeep2 core module were subjected to nucleotide BLAST search against the Sanger miRBase mature.fa and hairpin.fa database respectively (database release 21, http://www.mirbase.org) [19–21], this was done with the aim of discriminating known miRNAs and iso-miRNAs from potential novel miRNAs. Nucleotide BLAST for mature sequences was performed utilizing the blastn program with the blastn-short task default parameters, except for word-size which was changed from the default of 7 to 16; in this manner, only nearly identical sequences were identified. The resulting list was further filtered to keep only perfect and near-perfect matches to known miRNAS (near perfect defined as having one mismatch or up to a 2 nucleotide length difference between the query and subject sequences). BLAST for hairpin precursors was run with standard blastn with a word size 16 and cut-off e-value 1e-7, and the result was filtered to remove instances where the potential precursor was more than 8 bases longer than the aligned stretch.
To identify potential targets, the 3’UTR region was obtained from predicted gene annotations in the B. dorsalis genome assembly (GeneBank Accession # JFBF00000000.1), annotations were downloaded from USDA i5k web portal at https://i5k.nal.usda.gov/content/data-downloads. The 3’UTRs were compared against the identified miRNAs using miRANDA [22, 23]. A stringent threshold was applied for a conservative approach (pairing score: >155, energy score: <−7, gap opening penalty: −8 and Gap extension penalty: −8).
Differential expression of miRNAs
To assess the expression changes for the identified miRNAS across the life stages tested, the raw counts for the identified miRNAs obtained from the MirDeep2 quantifier.pl module were input into the EdgeR [24, 25], utilizing the pairwise exact-test modality to test each of the 11 possible comparisons on TMM (trimmed means of M-values) normalized counts. Correlations between miRNA expression in B. dorsalis and D. melanogaster were computed utilizing Spearman correlation (ρ). The method was chosen because we do not expect a perfect linear relationship between miRNA levels in both species due to developmental timing differences, Spearman correlation calculates coefficients based on rank. All figure generation and statistics were performed in R.
Results and discussion
Summary of sequencing results. Number of sequenced, processed and aligned reads
Total number of sequenced reads
Total usable reads equal to or longer than 17 nt
Total reads mapped
Number of unique sequences
Number of unique sequences mapping to at least one and up to 4 locations in the genome
Unique sequences matching an Rfam model (non-coding RNA)
46,142 (27.31 %)
3,471 (2.05 %)
114,655 (27.58 %)
6,723 (1.61 %)
97,218 (28.63 %)
5,117 (1.51 %)
55,738 (27.92 %
5,156 (2.58 %)
71,141 (23.47 %)
4,645 (1.53 %)
41,204 (23.88 %)
3,452 (2.00 %)
Male Pupae (rep1)
39,084 (15.10 %)
4,298 (1.66 %)
Male Pupae (rep2)
23,601 (14.32 %)
2,988 (1.81 %)
Male Pupae (rep3)
25,543 (14.46 %)
3,472 (1.96 %)
Female Pupae (rep1)
19,855 (11.40 %)
2,501 (1.43 %)
Female Pupae (rep2)
Female Pupae (rep3)
23,841 (15.52 %)
3,239 (2.1 %)
3,270 (17.92 %)
598 (3.27 %)
111,069 (27.85 %
6,114 (1.53 %)
145,910 (27.87 %)
8,025 (1.53 %)
(23.60 % of unique reads)
(27 % of unique reads)
Rfam classification. Number of unique sequences classified into each of the 12 main types of non-coding RNAs
Identification of known and novel miRNA in the dataset
MirDeep2 identified 149 potential miRNAs among the mapped reads in the pooled library, and the software assigned a provisional identification (provisional ID) to each of them. In five instances, two or three different provisional IDs were assigned to the exact same sequence because they aligned to more than one precursor excised from different genome scaffolds or from different regions in the same scaffold. The output was further filtered for miRNAs with a miRDeep2 score below or equal to 2, as a visual examination of the predicted hairpin structures and “stacks” built with the reads revealed that these hairpins had low signal to noise ratio (number of miRNA vs miRNA* and other reads in the region), and unusual secondary structures that did not resemble typical metazoan miRNA precursors [26, 27]. The remaining potential miRNA hairpins were screened to verify if they were contained within or overlapped predicted gene exons in the B. dorsalis genome. While the majority of the hairpin precursors were located outside of gene-coding regions and in introns, 13 of them were found either contained within an exon (10), contained a full exon (2), or overlapped an exon start (1). Whereas miRNAs are traditionally found in non-gene-coding-regions, we did not exclude these sequences from further analyses, as some more recently identified miRNAs in D. melanogaster were found within exonic regions, including 3’UTRs and coding sequences (CDSs) . In addition, because the gene set being used is largely computationally predicted, some errors in gene structure may be present.
Putative Novel miRNAs identified in Bactrocera dorsalis (orthologous to these miRNA sequences were not previously reported in other species). See also Fig. 1
miRBase assigned ID
Total read count
Consensus mature sequence
In animals, miRNAs regulate transcript abundance by complimentary base pairing to the 3’ UTR of the target RNA, with some exceptions . Potential targets for the 109 identified miRNAs were screened within the set of 3’UTR regions extracted from our most current B. dorsalis genome assembly (NCBI Assembly# ASM78921v2) utilizing the miRANDA software. The number of targets identified per each miRNA ranged between 5 and 469, even with the stringent parameters used for target detection, giving a total of 6506 predicted miRNA-target pairs, where many mRNA were targeted by multiple miRNAs. (Additional file 3: Table S3).
Differential expression of B. dorsalis miRNAs
The relative differential abundance of miRNAs between life stages was calculated using the negative binomial distribution with edgeR [25, 30]. In total, 47 of the 80 high confidence microRNAs identified showed differential expression (FDR corrected p-value <0.05) in at least one of the 10 possible pairwise comparisons between life stages. These significantly differentially regulated miRNAs included all 11 novel miRNA sequences. The comparison between female and male pupae yielded no differentially regulated miRNAs, and the most pronounced changes were observed between ovaries and pupae, and between eggs and pupae (Additional file 4: Table S4).
Visually inspecting each of the miRNA trend plots, there were 12 miRNAs with nearly identical expression patterns, and an additional miRNA (mir-276a), with a pattern varying only in the libraries derived from ovary tissue, which were the libraries with lowest correlation between the species. Among these twelve miRNAs with identical expression pattern, mir-282 was highly expressed in pupae of both D. melanogaster and B. dorsalis studies. Other studies have consistently found miR-282 to be expressed both in pupae and female adults of D. melanogaster [35, 38], and this miRNA has been shown to regulate viability and production of eggs through the targeting of the nervous-specific adenylate cyclase in pupae during metamorphosis . Mir-276a, had an equivalent trend in expression from embryo to pupa; however, the high expression of this miRNA in D. melanogaster ovary SSC was not observed in the fertilized B. dorsalis ovary library. This miRNA also had highest expression in pupae in both species, and is also implicated in neural development, specifically in olfactory response . Mir-137, mir-981, mir-87, and mir-927 had analogous expression patterns across both species, with highest expression observed in larval tissues, although their absolute expression was low compared to other significantly differentially regulated miRNAs. Even though these miRNAs were identified in genome-wide studies and computationally predicted in D. melanogaster [34, 36], we could only find one report of Dme-mir-87, showing its expression in relation to hormonal signaling , where mir-87 was found poorly expressed in early larva (in agreement with our data for B. dorsalis and the genome-wide studies in D. melanogaster), but then highly expressed in pupae. Finally, of the four conserved genomic miRNA clusters identified in B. dorsalis, only the mir-309-6 cluster (dme-mir-4, dme-mir-6-3, dme-mir-5, dme-mir-286, and dme-mir-309) had an analogous expression pattern across life stages.
Genomic clusters of miRNAs in B. dorsalis
miRNA clusters in B. dorsalis
Mature read count
The mir-309-6 cluster has been implicated in maternal transcript destabilization, the removal of transcripts maternally provided during oogenesis [44, 45]. Consistent with this observation, the expression of the mir-309-6 cluster was highest at the start of zygotic transcription in the early embryo in D. melanogaster. Transcripts for this in D. melanogaster were found subsequently depleted, but mature miRNAs could still be detected in larval stages [46, 47]. In our data for B. dorsalis, all miRNA members of this cluster were found to be very abundant in embryos, with counts per million reads comparable to those in D. melanogaster. While some miRNAs were still detected in the larval stage, no mature miRNAs were found in pupal samples, and a significant number of miRNAs were detected in fertilized ovary tissue (see next section). Taken together, this data indicates that this miRNA cluster is not specific to Drosophilids as previously believed, and that it has undergone extensive evolutionary divergence. Functionally, while the D. melanogaster mir-309-6 cluster was found to act as part of the zygotic machinery in the removal of maternal mRNAs, and is highly expressed in embryos, our B. dorsalis data indicates that the mir-309-6 cluster is functional before egg laying, either because zygotic transcription starts earlier than in D. melanogaster, or because these miRNAs are not only part of the zygotic machinery for maternal transcript decay, but also part of the maternal machinery, which plays the same role (removal of maternal mRNA) before the onset of zygotic transcription.
Deep sequencing of small RNAs has allowed the identification of miRNAs present at four different life stages of the Tephritid fruit fly Bactrocera dorsalis. Sixty-nine miRNAs homologs to miRNAs in other species were identified with high confidence, and sufficient evidence was gathered to identify eleven additional miRNAs that were not previously reported. The latter may include conserved miRNAs with relatively low expression, and/or miRNAs that have evolved independently and are specific to the Tephritid family. The three replications per library allowed for a robust differential expression analysis of the identified miRNAs and their classification into life stage specific groups; miRNAs falling in these categories could be considered of importance because they are likely involved during transitional stages in development. To complement the data, and to provide more biological insight, we attempted to provide a list of potential targets for the identified miRNAs. Given that in metazoans, perfect complementarity to only six nucleotides in the seed region of the small RNA and the target is sufficient to promote RNA silencing, like in Drosophila, the resulting list of candidate mRNA targets was very extensive [21, 48], even with the stringent parameters we set for the miRANDA output. Although the dataset of potential targets was narrowed down to a short list of miRNAs for which RNAseq information was available, there is still need for wet lab experiments for confirmation of the targets.
Taking advantage of the vast resources and data available for the model species D. melanogaster, comparative analysis across conserved orthologous miRNAs were utilized to further support our findings. High correlation was identified between datasets at the level of abundance across developmental stages. Moreover, groups of miRNAs that are physically clustered in genomic regions were also found to be conserved between both B. dorsalis and D. melanogaster. Although these miRNA clusters differed in genomic spacing between the two organisms, this difference was consistent for all the clusters, being B. dorsalis clusters arranged in regions up to three times longer than D. melanogaster miRNAs but with the same order of the miRNAs in the cluster. Only one of the miRNA clusters, the mir309-6 showed poor conservation including repeated miRNAs and a different arrangement. Although genome assembly errors were a possibility, the same cluster, with the same arrangement was found in the sequenced genomes of the two other Tephritids, namely B. cucurbitae, and C. capitata. With the data available for B. dorsalis, we hypothesize that this cluster, which has highly diverged from the D. melanogaster dme-mir-309-6 cluster, may also function in maternal transcript destabilization machinery as it does in Drosophilids, however, because it is also expressed in ovaries, it may not be specific for the zygotic machinery. Importantly, this cluster was previously reported as being specific for drosophilids, and we proved that is not, demonstrating that this dataset, and similar datasets from Tephritids can be used as comparative tools for flies and other insects, to draw more robust conclusions about evolutionary questions.
Knowledge on miRNAs in B. dorsalis could help in developing novel pest control strategies, for example, miRNAs that are specific for egg and larval stages, likely involved in key pathways for developmental transitions, can be further characterized and utilized in miRNA mimics feeding and plant expression [48, 49]. Because miRNAs are very important in controlling developmental states, miRNA mimics targeting female specific sex determination and development transcripts could be used to generate genetic sexing strains that can be utilized in Sterile Insect Technique (SIT). Finally, this dataset could be further explored to find other specific regulatory pathways of interest, and as an aid for functional characterization of genes.
We thank Steven Tam for assistance in colony rearing and fruit fly sample collections used in this study, and Brian Hall for assistance in computational analysis. Funding was provided by USDA-ARS, and the bioinformatic analysis was performed on computing resources at USDA-ARS Pacific Basin Agricultural Research Center (Moana cluster; Hilo, HI, http://moana.dnsalias.org) and the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number OCI-1053575XSEDE utilizing allocation TG-MCB140032 to SMG. Opinions, findings, conclusions, or recommendations expressed in this publication are those of the authors and do not necessarily reflect the views of the USDA. USDA is an equal opportunity provider and employer.
Open AccessThis 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.
- Smit J, Reijnen B, Stokvis F. Half of the European fruit fly species barcoded (Diptera, Tephritidae); a feasibility test for molecular identification. ZooKeys. 2013;365:279–305.View ArticlePubMedGoogle Scholar
- Clarke AR, Armstrong KF, Carmichael AE, Milne JR, Raghu S, Roderick GK, et al. Invasive phytophagous pests arising through a recent tropical evolutionary radiation: the Bactrocera dorsalis complex of fruit flies. Annu Rev Entomol. 2005;50(1):293–319.View ArticlePubMedGoogle Scholar
- Bartel DP. MicroRNAs: genomics, biogenesis, mechanism, and function. Cell. 2004;116(2):281–97.View ArticlePubMedGoogle Scholar
- Ghildiyal M, Zamore PD. Small silencing RNAs: an expanding universe. Nat Rev Genet. 2009;10(2):94–108.PubMed CentralView ArticlePubMedGoogle Scholar
- Kim VN, Han J, Siomi MC. Biogenesis of small RNAs in animals. Nat Rev Mol Cell Biol. 2009;10(2):126–39.View ArticlePubMedGoogle Scholar
- Bushati N. Cohen SM: microRNA Functions. Annu Rev Cell Dev Biol. 2007;23(1):175–205.View ArticlePubMedGoogle Scholar
- McCombs SD, Saul SH. Translocation-based genetic sexing system for the oriental fruit fly (Diptera: Tephritidae) based on pupal color dimorphism. Ann Entomol Soc Am. 1995;88(5):695–8.View ArticleGoogle Scholar
- Chang CL, Vargas RI. Wheat germ oil and its effects on a liquid larval rearing diet for oriental fruit flies (Diptera: Tephritidae). J Econ Entomol. 2007;100(2):322–6.View ArticlePubMedGoogle Scholar
- Geib S, Calla B, Hall B, Hou S, Manoukis N. Characterizing the developmental transcriptome of the oriental fruit fly, Bactrocera dorsalis (Diptera: Tephritidae) through comparative genomic analysis with Drosophila melanogaster utilizing modENCODE datasets. BMC Genomics. 2014;15(1):942.PubMed CentralView ArticlePubMedGoogle Scholar
- Mackowiak, SD. Identification of Novel and Known miRNAs in Deep-Sequencing Data with miRDeep2. Current Protocols in Bioinformatics. 2011;36:12.10:12.10.1–12.10.15.Google Scholar
- Friedländer 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.PubMed CentralView ArticlePubMedGoogle Scholar
- Langmead B, Trapnell C, Pop M, Salzberg S. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009;10(3):R25.PubMed CentralView ArticlePubMedGoogle Scholar
- Lorenz R, Bernhart SH, Höner zu Siederdissen C, Tafer H, Flamm C, Stadler PF, et al. ViennaRNA Package 2.0. Algorithms Mol Biol. 2011;6:26.PubMed CentralView ArticlePubMedGoogle Scholar
- Hofacker IL, Fontana W, Stadler PF, Bonhoeffer LS, Tacker M, Schuster P. Fast folding and comparison of RNA secondary structures. Monatshefte für Chemie / Chemical Monthly. 1994;125(2):167–88.View ArticleGoogle Scholar
- Hofacker IL, Stadler PF. Memory efficient folding algorithms for circular RNA secondary structures. Bioinformatics. 2006;22(10):1172–6.View ArticlePubMedGoogle Scholar
- Zuker M, Stiegler P. Optimal computer folding of large RNA sequences using thermodynamics and auxiliary information. Nucleic Acids Res. 1981;9(1):133–48.PubMed CentralView ArticlePubMedGoogle Scholar
- Nawrocki EP, Burge SW, Bateman A, Daub J, Eberhardt RY, Eddy SR, et al. Rfam 12.0: updates to the RNA families database. Nucleic Acids Res. 2015;43(Database issue):D130–7.Google Scholar
- Nawrocki EP, Eddy SR. Infernal 1.1: 100-fold faster RNA homology searches. Bioinformatics. 2013;29(22):2933–5.PubMed CentralView ArticlePubMedGoogle Scholar
- Griffiths-Jones, S. miRBase: microRNA Sequences and Annotation. Current Protocols in Bioinformatics. 2010;29:12.9:12.9.1–12.9.10.Google Scholar
- Kozomara A, Griffiths-Jones S. miRBase: annotating high confidence microRNAs using deep sequencing data. Nucleic Acids Res. 2014;42(D1):D68–73.PubMed CentralView ArticlePubMedGoogle Scholar
- Griffiths‐Jones S. The microRNA registry. Nucleic Acids Res. 2004;32 suppl 1:D109–11.PubMed CentralView ArticlePubMedGoogle Scholar
- Betel D, Koppal A, Agius P, Sander C, Leslie C. Comprehensive modeling of microRNA targets predicts functional non-conserved and non-canonical sites. Genome Biol. 2010;11(8):R90.PubMed CentralView ArticlePubMedGoogle Scholar
- Stokowy T, Eszlinger M, Swierniak M, Fujarewicz K, Jarzab B, Paschke R, et al. Analysis options for high-throughput sequencing in miRNA expression profiling. BMC Res Notes. 2014;7(1):144.PubMed CentralView ArticlePubMedGoogle Scholar
- McCarthy DJ, Chen Y, Smyth GK. Differential expression analysis of multifactor RNA-Seq experiments with respect to biological variation. Nucleic Acids Res. 2012;40(10):4288–97.PubMed CentralView ArticlePubMedGoogle Scholar
- Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26(1):139–40.PubMed CentralView ArticlePubMedGoogle Scholar
- Tarver JE, Donoghue PCJ, Peterson KJ. Do miRNAs have a deep evolutionary history? Bioessays. 2012;34(10):857–66.View ArticlePubMedGoogle Scholar
- Axtell M, Westholm J, Lai E. Vive la difference: biogenesis and evolution of microRNAs in plants and animals. Genome Biol. 2011;12(4):221.PubMed CentralView ArticlePubMedGoogle Scholar
- Berezikov E, Robine N, Samsonova A, Westholm JO, Naqvi A, Hung J-H, et al. Deep annotation of Drosophila melanogaster microRNAs yields insights into their processing, modification, and emergence. Genome Res. 2011;21(2):203–15.PubMed CentralView ArticlePubMedGoogle Scholar
- Ambros V, Bartel B, Bartel DP, Burge CB, Carrington JC, Chen X, et al. A uniform system for microRNA annotation. RNA. 2003;9(3):277–9.PubMed CentralView ArticlePubMedGoogle Scholar
- Robinson MD, Smyth GK. Small-sample estimation of negative binomial dispersion, with applications to SAGE data. Biostatistics. 2008;9(2):321–32.View ArticlePubMedGoogle Scholar
- Pasquinelli AE, Reinhart BJ, Slack F, Martindale MQ, Kuroda MI, Maller B, et al. Conservation of the sequence and temporal expression of let-7 heterochronic regulatory RNA. Nature. 2000;408(6808):86–9.View ArticlePubMedGoogle Scholar
- Sokol NS, Xu P, Jan Y-N, Ambros V. Drosophila let-7 microRNA is required for remodeling of the neuromusculature during metamorphosis. Genes Dev. 2008;22(12):1591–6.PubMed CentralView ArticlePubMedGoogle Scholar
- Sempere LF, Dubrovsky EB, Dubrovskaya VA, Berger EM, Ambros V. The expression of the let-7 small regulatory RNA is controlled by ecdysone during metamorphosis in Drosophila melanogaster. Dev Biol. 2002;244(1):170–9.View ArticlePubMedGoogle Scholar
- Stark A, Kheradpour P, Parts L, Brennecke J, Hodges E, Hannon G, et al. Systematic discovery and characterization of fly microRNAs using 12 Drosophila genomes. Genome Res. 2007;17:1865–79.PubMed CentralView ArticlePubMedGoogle Scholar
- Lai E, Tomancak P, Williams R, Rubin G. Computational identification of Drosophila microRNA genes. Genome Biol. 2003;4:R42.PubMed CentralView ArticlePubMedGoogle Scholar
- Ruby JG, Stark A, Johnston WK, Kellis M, Bartel DP, Lai EC. Evolution, biogenesis, expression, and target predictions of a substantially expanded set of Drosophila microRNAs. Genome Res. 2007;17(12):1850–64.PubMed CentralView ArticlePubMedGoogle Scholar
- Raff RA. The Shape of Life: Genes, Development, and the Evolution of Animal Form: University of Chicago Press; 1996.Google Scholar
- Fernandez-Valverde SL, Taft RJ, Mattick JS. Dynamic isomiR regulation in Drosophila development. RNA. 2010;16(10):1881–8.PubMed CentralView ArticlePubMedGoogle Scholar
- Vilmos P, Bujna Á, Szuperák M, Havelda Z, Várallyay É, Szabad J, et al. Viability, longevity, and egg production of Drosophila melanogaster are regulated by the miR-282 microRNA. Genetics. 2013;195(2):469–80.PubMed CentralView ArticlePubMedGoogle Scholar
- Li W, Cressy M, Qin H, Fulga T, Van Vactor D, Dubnau J. MicroRNA-276a functions in ellipsoid body and mushroom body neurons for naive and conditioned olfactory avoidance in drosophila. J Neurosci. 2013;33(13):5821–33.PubMed CentralView ArticlePubMedGoogle Scholar
- Sempere LF, Sokol NS, Dubrovsky EB, Berger EM, Ambros V. Temporal regulation of microRNA expression in Drosophila melanogaster mediated by hormonal signals and broad-Complex gene activity. Dev Biol. 2003;259(1):9–18.View ArticlePubMedGoogle Scholar
- Ninova M, Ronshaugen M, Griffiths-Jones S. Fast-evolving microRNAs are highly expressed in the early embryo of Drosophila virilis. RNA. 2014;20(3):360–72.PubMed CentralView ArticlePubMedGoogle Scholar
- Aravin AA, Lagos-Quintana M, Yalcin A, Zavolan M, Marks D, Snyder B, et al. The small RNA profile during Drosophila melanogaster development. Dev Cell. 2003;5(2):337–50.View ArticlePubMedGoogle Scholar
- Tadros W, Goldman AL, Babak T, Menzies F, Vardy L, Orr-Weaver T, et al. SMAUG is a major regulator of maternal mRNA destabilization in Drosophila and its translation is activated by the PAN GU kinase. Dev Cell. 2007;12(1):143–55.View ArticlePubMedGoogle Scholar
- Bushati N, Stark A, Brennecke J, Cohen SM. Temporal reciprocity of miRNAs and their targets during the maternal-to-zygotic transition in Drosophila. Curr Biol. 2008;18(7):501–6.View ArticlePubMedGoogle Scholar
- Biemar F, Zinzen R, Ronshaugen M, Sementchenko V, Manak JR, Levine MS. Spatial regulation of microRNA gene expression in the Drosophila embryo. Proc Natl Acad Sci U S A. 2005;102(44):15907–11.PubMed CentralView ArticlePubMedGoogle Scholar
- Thomsen S, Anders S, Janga SC, Huber W, Alonso C. Genome-wide analysis of mRNA decay patterns during early Drosophila development. Genome Biol. 2010;11(9):R93.PubMed CentralView ArticlePubMedGoogle Scholar
- Jayachandran B, Hussain M, Asgari S. An insect trypsin-like serine protease as a target of microRNA: utilization of microRNA mimics and inhibitors by oral feeding. Insect Biochem Mol Biol. 2013;43(4):398–406.View ArticlePubMedGoogle Scholar
- Agrawal N, Sachdev B, Rodrigues J, Sree KS, Bhatnagar RK. Development associated profiling of chitinase and microRNA of Helicoverpa armigera identified chitinase repressive microRNA. Sci Rep. 2013;3:2292.PubMed CentralPubMedGoogle Scholar