Skip to main content
Download PDF

Integrated analysis of miRNA and mRNA expression profiles in tilapia gonads at an early stage of sex differentiation

BMC Genomics volume 17, Article number: 328 (2016) Cite this article

Abstract

Background

MicroRNAs (miRNAs) represent a second regulatory network that has important effects on gene expression and protein translation during biological process. However, the possible role of miRNAs in the early stages of fish sex differentiation is not well understood. In this study, we carried an integrated analysis of miRNA and mRNA expression profiles to explore their possibly regulatory patterns at the critical stage of sex differentiation in tilapia.

Results

We identified 279 pre-miRNA genes in tilapia genome, which were highly conserved in other fish species. Based on small RNA library sequencing, we identified 635 mature miRNAs in tilapia gonads, in which 62 and 49 miRNAs showed higher expression in XX and XY gonads, respectively. The predicted targets of these sex-biased miRNAs (e.g., miR-9, miR-21, miR-30a, miR-96, miR-200b, miR-212 and miR-7977) included genes encoding key enzymes in steroidogenic pathways (Cyp11a1, Hsd3b, Cyp19a1a, Hsd11b) and key molecules involved in vertebrate sex differentiation (Foxl2, Amh, Star1, Sf1, Dmrt1, and Gsdf). These genes also showed sex-biased expression in tilapia gonads at 5 dah. Some miRNAs (e.g., miR-96 and miR-737) targeted multiple genes involved in steroid synthesis, suggesting a complex miRNA regulatory network during early sex differentiation in this fish.

Conclusions

The sequence and expression patterns of most miRNAs in tilapia are conserved in fishes, indicating the basic functions of vertebrate miRNAs might share a common evolutionary origin. This comprehensive analysis of miRNA and mRNA at the early stage of molecular sex differentiation in tilapia XX and XY gonads lead to the discovery of differentially expressed miRNAs and their putative targets, which will facilitate studies of the regulatory network of molecular sex determination and differentiation in fishes.

Background

MicroRNAs (miRNAs), a class of small non-coding RNAs, have emerged as an important and abundant group of gene expression regulators in a wide range of species [1]. Mature miRNA sequences are approximately 20–24 nucleotides (nt) long, and are formed from longer primary miRNAs that contain hairpin stem-loop structures [2]. miRNAs recognize the complementary 3’-untranslated regions (3’-UTRs) of mRNAs to inhibit their expression [3–6]. Each miRNA may have multiple gene targets, and each gene target may also be regulated by multiple miRNAs [7, 8]. miRNAs are often expressed in a tissue-specific manner, and a large proportion of protein coding genes of animals are regulated by miRNAs [9]. miRNAs are key molecules in many biological processes, including apoptosis, cell proliferation and tissue development.

Among teleosts, miRNAs were first examined in zebrafish [10]. miRNA catalogs of various tissues in other fish species were characterized more recently [11–17]. Currently, 1623 miRNAs have been identified in 9 teleost species, compared to 2588 miRNAs in human (miRBase 21.0, http://www.mirbase.org/) [18], suggesting that more miRNAs await characterization in fishes. Next generation sequencing facilitates the profiling of both known and novel miRNAs, especially those expressed in low abundance. Data from such deep sequencing will help to characterize the structure of miRNA molecules, and will provide quantitative measurement of their expression.

Teleosts are the largest group of vertebrates and exhibit a remarkable variety of sexual systems [19, 20]. Genes, environmental factors, and even social interaction may contribute to sex determination and differentiation in a large number of fish species [21]. Compared to other vertebrates, the molecular mechanisms of sex determination and differentiation in teleosts are quite diverse [22]. Over the past few years, seven master sex-determining genes have been identified in teleosts, namely dmy/dmrt1b Y in Oryzias latipes and O. curvinotus [23, 24], gsdf Y in O. luzonensis [25], sox3 in O. dancena [26], amhy in Odontesthes hatchery [27] and Oreochromis niloticus [28], amhrII in Takifugu rubripes [29], gdf6Y in Nothobranchius furzeri [30] and sd Y in Oncorhynchus mykiss and several other salmonids [31, 32]. However, it is well-known that sex steroid hormones, especially estrogens, may influence early sex differentiation in fishes, in spite of variation in the master sex determining genes and the non-conserved subsequent genetic steps of sex differentiation [33]. Previous studies have been focused mainly on the function of Cyp19a1a, which encodes the key enzyme responsible for the conversion of androgens to estrogens [34–36]. The expression of Cyp19a1a in the gonads is mediated by both Foxl2 [37] and Dmrt1 [38, 39].

Several studies have explored the role of miRNAs in regulating sex differentiation. let-7 and miR-21 have been shown to regulate development of rainbow trout eggs [16]. miRNAs are also involved in oocyte development, hydration, and competence, indicating their importance in the regulation of oogenesis [40]. miR-141 and miR-429 have crucial functions in yellow catfish testis development and spermatogenesis [41]. Different expression pattern of miRNAs is observed between the embryos of females and males [42], as well as the ovaries and testes in tilapia gonads at later stages [17, 43], indicating their regulatory roles during fish reproduction. let-7a, miR-143, and miR-202 are upregulated to induce testes differentiation in halibut [44]. Therefore, miRNAs may represent novel regulators of gonadal development and sexual differentiation. However, in these studies the small RNA libraries were generated from whole embryos or from gonads at later developmental stages, which may mask important sex-biased miRNAs in gonads at the critical stage of sex differentiation. Importantly, the regulatory roles of miRNAs during early sex differentiation of fishes have not been investigated extensively. Functional prediction of miRNAs solely based on computational miRNA-mRNA interactions bears a high number of false positive predictions [45]. Thus, the joint investigation of miRNA and mRNA expression may help us to improve the quality of predicted interactions and understand the molecular mechanisms of post-transcriptional regulation.

Nile tilapia (O. niloticus) is one of the most important farmed fish with a production exceeding 4.5 million tons in 2014 (http://faostat3.fao.org/home/E). The growth rate of males is significantly higher than that of females in tilapia. All-male stocks are preferred because they improve the efficiency of commercial production. It is therefore important to understand the molecular mechanism of sex determination in this species. Thus, quantifying the expression of both miRNAs and mRNAs in the undifferentiated gonads at the critical stage of molecular sex determination [46], approximately 5 days after hatching (dah), might help to clarify the regulatory network during early sex differentiation and provide new information on the role of miRNAs in gonadal function.

Results

Identification of conserved miRNAs in the tilapia genome

The reference set of 1623 known teleost pre-miRNA sequences was used for a similarity search against the tilapia reference genome assembly. We manually checked the best hits to extract putative tilapia pre-miRNAs and confirmed that they were able to fold into the secondary stem-loop structure necessary for miRNA biogenesis. This resulted in the prediction of 279 distinct tilapia pre-miRNA genes (Additional file 1: Table S1). Pre-miRNAs were non-randomly distributed in the tilapia genome (X 2, p < 0.05). No pre-miRNAs were observed on LG3, while 29 were observed on LG7 (Fig. 1). Also, some pre-miRNAs were encoded at a single chromosomal location (e.g. miR-100 on LG10, and miR-106a on LG1), while others were encoded on multiple chromosomes (e.g. miR-103 on LG2, LG13 and LG19; miR-124 on LG9, LG15, LG19 and LG20). Comparisons of the number of pre-miRNAs identified in the tilapia genome with those reported in other fishes were shown in Fig. 2 and Additional file 2: Figure S1. The data suggested that most of the pre-miRNAs identified in tilapia have homologs in other fishes. To further study the conservation and evolution of the known pre-miRNAs in teleosts, sequence comparisons were also made with Astatotilapia burtoni, Neolamprologus brichardi, Pundamilia nyererei, Metriaclima zebra, Oryzias latipes, Gasterosteus aculeatus and Danio rerio. Based on these analyses, most of the pre-miRNAs showed a high level of conservation among the various fish species (Additional file 3: Figure S2).

Fig. 1
figure 1

Genomic distribution of identified pre-miRNAs in tilapia

Full size image
Fig. 2
figure 2

Comparisons of the number of pre-miRNAs of Oryzias latipes, Ictalurus punctatus, Danio rerio and Salmo salar in miRbase with the number of pre-miRNAs identified in O. niloticus

Full size image

High-throughput sequencing of miRNAs in tilapia gonads

Deep sequencing of miRNAs from tilapia XX and XY gonads at 5 dah was performed using the Illumina HiSeq 2000 platform. Approximately 12 million reads were obtained from each library. After filtering low-quality reads and empty adaptor sequences, 11,065,356 sequences from the XX gonad and 11,123,822 sequences from the XY gonad were mapped to the tilapia genome, amounting to 92.6 and 91.5 % of the total miRNA libraries, respectively (Table 1). A total of 8,203,487 (XX gonads) and 8,352,968 (XY gonads) reads were mapped to miRBase entries, representing a total of 635 miRNAs (Additional file 4: Table S2). In both XX and XY gonads, the majority of miRNAs were in the range of 21 to 24 nts, accounting for 70.18 and 64.05 % of the total mappable reads (Fig. 3). The average length of the mature miRNAs in the XX and XY gonad were 22.8 and 23.0 bp, respectively.

Table 1 Number of high-throughput reads generated from miRNA libraries of tilapia XX and XY gonads at 5 dah
Full size table
Fig. 3
figure 3

Length distribution of miRNA sequences in ovary and testes of tilapia by Illumina small RNA deep sequencing. Sequence length distribution of clean reads based on the abundance and distinct sequences; the most abundant size class was 22 nt, followed by 23 nt and 21 nt

Full size image

To profile the chromosomal distribution of active miRNAs detected, the read counts for each chromosome from both XX and XY gonads were consolidated and the percentage of transcribed miRNA reads per chromosome was computed. Approximately 40 % of the mature miRNAs were localized to LG 1, 2, 16, and 20 (Fig. 4). Only a few reads were mapped to LG3, which is large and known to be highly repetitive. The identified miRNAs in tilapia covered a large proportion of reported miRNAs in fishes, such as D. rerio, T. nigroviridis, T. rubripes, O. latipes (Additional file 5: Figure S3). Most of the previously reported cichlid miRNAs [47] exhibited orthology to the identified miRNAs in the present study (Additional file 6: Figure S4).

Fig. 4
figure 4

Circos circular visualization of the miRNA reads coverage on tilapia genome

Full size image

The top 7 most abundant miRNAs in tilapia gonads included miR-10a-5p, miR-10b, miR-10c-5p, miR-10d, miR-21, miR-100, and miR-181a-5p and accounted for 72.89 % (74.99 %) of the 8,203,487 (8,352,968) reads mapped to miRBase. Sequence analysis indicated that the relative abundance of miRNAs from one miRNA family varied greatly in the tilapia gonads, suggesting possible functional divergence among different members within the miRNA family. For example, abundance of different members of the miR-10 family varied from 9 reads (miR-10c-3p) to 2,765,827 and 2,848,186 (miR-10b) in tilapia XX and XY gonads, respectively.

Identification of novel sex-biased miRNAs in tilapia gonads

From the reads mapped to tilapia genome, we identified potentially novel miRNAs that were not matched to any sequences in the miRBase. Based on analyzing the surrounding sequences (50 nts on both directions) for the ability to form hairpin structures, 130 novel miRNAs displayed sex-biased expression (Additional file 7: Table S3). These putative miRNAs showed sequences variations at the 5’ or 3’ terminus compared with known miRNAs in miRBase. Of the 130 novel miRNAs, 49 and 45 miRNAs were identified uniquely in XX and XY gonads, respectively. A greater proportion of the novel miRNAs were sex-biased than were the set of known miRNAs.

Identification of miRNAs uniquely expressed in tilapia XX and XY gonads

Among the 635 mature miRNAs identified in the two libraries, 557 miRNAs were expressed in both libraries. We next directed our attention to the miRNAs that were identified exclusively from either XX or XY gonads. 39 miRNAs including miR-34c-5p, miR-153-5p, miR-749, and miR-2187-5p were expressed exclusively in XX gonads, while 49 miRNAs including miR-1306-5p, miR-132b and miR-18c were expressed exclusively in XY gonads. However, no miRNAs appeared exclusively in either sex with more than 15 reads.

Identification of miRNAs differentially expressed in tilapia gonads

We then compared the transcription level of differentially expressed miRNAs between XX and XY libraries to define miRNAs related to early ovary or testis differentiation. 62 miRNAs showed higher expression levels in XX gonads, whereas 49 were more highly expressed in XY gonads (Fig. 5). miR-9, miR-101, miR-192, miR-202, miR-212, miR-217a, and miR-737 showed greater expression in ovaries. The miRNAs upregulated in XY gonads included miR-27d, miR-29b, miR-92b, miR-140-3p, miR-144-5p, miR-375, miR-455 and miR-7977. Some abundantly expressed miRNAs also showed sex-biased expression with ~1.5 fold change, including miR-21, miR-96, and miR-200b.

Fig. 5
figure 5

Comparison of expression levels of miRNAs in XX and XY gonads. The X axis and Y axis show expression level of miRNAs in the two samples, respectively. Red points represent miRNAs with ratio > 2; green points represent miRNAs with 1/2 < ratio ≤ 2; blue points represent miRNAs with ratio ≤ 1/2. Ratio = normalized expression of XX gonad/normalized expression of XY gonad

Full size image

Identification of miRNAs involved in steroidogenic pathways and early sex differentiation

A total of 10,740 target genes were predicted for the 635 miRNAs by miRanda and Targetscan. The predicted modes of action for the predominant miRNAs in gonad included: 1) a single miRNA targeting a single or multiple position(s) in the 3’ UTR; 2) multiple miRNAs targeting different positions in the 3’ UTR of the same gene. The expression of target genes was quantified by mRNA sequencing at 5 dah. Sequencing of the mRNA libraries yielded ~800 million reads for female and male pool. In order to identify major miRNAs with possible role in early sex differentiation, we looked for target genes which showed inverse expression patterns compared to female and male upregulated miRNAs (Additional file 8: Table S4 and Additional file 9: Table S5). On average, miRanda and Targetscan predictions suggested about 4000 targets per miRNA compared to integrative analysis of miRNA and mRNA expression data which suggested about 90 target predictions per miRNA.

Subsequently, we examined the expression of key genes known to be involved in the steroid hormone biosynthesis pathway [48] and those that participate in vertebrate sex differentiation [46]. Cyp19a1a, Foxl2, Cyp11a1, Hsd3b and Amh showed female upregulated expression, while Dmrt1 and Hsd11b showed up-regulation in males at 5 dah. Nine microRNAs (miR-30a, miR-30a-3p, miR-30b-3p, miR-30d-3p, miR-30e-3p, miR-153--5p, miR-205b-3p, miR-737-5p, and miR-737-3p) were predicted to regulate the expression of Cyp19a1a. Among these 9 microRNAs, all members of miR-30 family showed down-regulation in females, though the difference was not a 2-fold decrease. Dmrt1 was predicted to be the target of miR-9-3p, miR-22-5p, miR-132, miR-212, miR-212-3p and miR-7641, and the expression of miR-212-3p was downregulated in males. miR-203a, 203-3p, miR-218, miR-429, miR-722, and miR-7977 were predicted to target Foxl2, and the expression of miR-7977 was downregulated in the XX gonad (Fig. 6 and Additional file 10: Figure S5). Interestingly, some miRNAs targeted two or more genes in the regulatory network of tilapia gonadal differentiation: miR-96 targeted both Amh and Hsd3b, while miR-737 targeted both Star1 and Hsd11b. Other miRNAs targeting key genes in the pathway for sex steroids synthesis and early sex differentiation, including Gsdf and Star1, were listed in Additional file 11: Table S6.

Fig. 6
figure 6

miRNA-gene network in the steroid hormone biosynthesis pathway [48] and regulatory network in early sex differentiation of Nile tilapia. The genes in boxes are differentially expressed in the mRNA-Seq and were validated by previous transcriptomic analyses [72]. Blue box nodes represent upregulated genes in XY gonads, red box nodes represent upregulated genes in XX gonads, and black box nodes represent genes with no significant difference between XX and XY. miRNAs in the ellipses potentially regulated the corresponding gene expression. Blue ellipses represent downregulated miRNAs in XX gonads; red ellipses represent downregulated miRNAs in XY gonads

Full size image

Discussion

miRNAs regulate a vast number of genes [9], which makes them likely to play a critical role in programming the ovary and testis. In the present study, both in silico prediction of conserved pre-miRNAs and an integrated analysis of deep sequencing of miRNA and mRNA from XX and XY gonads at the critical stage of sex determination were performed to identify conserved miRNAs in fish and sex-associated miRNAs in tilapia.

Previous bioinformatic studies identified 21 conserved miRNAs and 100 pre-miRNAs in cichlid fishes, using EST, GSS and partial genomic sequences [49, 50]. The present study identified a total of 279 conserved pre-miRNAs, by comparing the published tilapia genome to pre-miRNAs deposited in miRBase 21.0. The increase can be attributed to the submission of additional fish miRNAs to miRBase, and to an improved assembly of the tilapia genome.

The number of miRNA precursors varied significantly among chromosomes. The smallest number of miRNAs was identified on LG3 by both in silico prediction and miRNA reads mapping, even though it is the longest chromosome in the karyotype. LG7 carries a higher number of miRNA genes, suggesting an important role for this chromosome in regulating tilapia tissue development and differentiation as in human [51]. miRNAs are also reported to be non-randomly distributed in the human and mouse genomes [52–54]. In addition, the multiple locations of the same miRNA on different LGs, such as pre-miR-103 and pre-miR-124, suggest the possibility that transcription of the same miRNA from different sites may regulate the expression of target genes in different tissues or at different development stages [55].

In fish, the miRNA distribution in organs of the reproductive axis (brain, pituitary and gonad) has recently been reported in mature stages of the Atlantic halibut [44], yellow catfish [41], and Nile tilapia [17, 43]. However, to date there is little data on miRNAs that have a possible role in early sex differentiation, especially in fish gonads. Only one previous study identified 9 sexually differentially expressed miRNAs in early developing tilapia embryos [42]. Thus, we performed miRNA and mRNA Illumina sequencing on the same gonadal sample to identify the complete set of known and novel miRNAs, as well as explore their possible roles in early sex differentiation.

The miRNA libraries of XX and XY gonads displayed similar read length distribution after filtering high quality sequencing. The most abundant miRNAs in both libraries were 21 and 22 nt, followed by 23 nt, representing the typical size of Dicer-derived products, and is consistent with the length distribution of miRNA in mouse [56]. Ovary and testis at the mature stages in zebrafish and tilapia show significant piRNA expression with another peak at 26–28 nt due to the abundant expression of PIWI-interacting RNAs (piRNAs) [17, 57]. In the present study, there were more miRNAs with length of 21 and 22 nt in XX gonads, whereas in XY gonads there were more miRNAs with length of 23 nt, which has also been observed in a previous study of gonadal miRNAs in tilapia [17].

In tilapia gonads at 5 dah, the most abundant miRNAs were miR-10a-5p, miR-10b, miR-10c-5p, miR-10d, miR-21, miR-100, and miR-181a-5p, each with more than several thousand reads. As in human, members of miR-10 family dominate the miRNA population in both ovary and testis, consistent with a relationship of members of miR-10 with sex hormones [58]. Expression of miR-21 increases substantially following the induction of differentiation, suggesting its important roles in stem cell differentiation [59] and endocrine regulation [42, 60]. miR-100 is also highly expressed in mature ovaries and testes of yellow catfish [41], indicating its possible involvement in gonadal development. The miR-181a family is abundantly expressed in the gonads of tilapia [17], mouse [61], and human [62], which suggests that this family also may play a crucial role in sex determination and differentiation. The most abundantly expressed miRNAs in tilapia XX and XY gonads at 5 dah are functionally conserved with other vertebrates, indicating that these miRNAs may play conserved role in gonadal development and differentiation. However, we also identified 130 novel sex-biased miRNAs that were not similar to any known conserved miRNAs in miRBase. These novel miRNAs with overall low abundance may be newly evolved in tilapia, or may have been overlooked in other species. The preferential expression of one of the duplexes (−5p or -3p) of mature miRNA (e.g. miR-9b-5p, miR-9b-3p) has also been reported from different tissues of mouse [63]. These results indicate that different duplexes within one miRNA can have clearly different expression levels, mostly due to tissue or developmental stage-specific expression.

No sex-specific miRNAs with reads higher than 15 were identified in tilapia gonads, a result that is consistent with the previous study where the expression levels of the tissue-specific expressed miRNAs were low in Holstein cattle (1–23 reads) [64]. More miRNAs were significantly upregulated in ovaries, which is also observed in adult Nile tilapia [17], medaka [65] and Holstein cattle [64]. Among these upregulated miRNAs, some have been reported to be very important in ovarian development in various species, such as miR-101 in regulating ovary differentiation in embryonic chicken gonads [66], miR-212 in ovarian development with the putative target Wt1 and Gata2 [67, 68]. Our results (higher expression of miRNAs in XY gonads) are in agreement with previous studies. For example, miR-140-3p is also abundantly expressed in Sertoli cells of the developing mouse testis [69]. miR-375 has been suggested to play an important role in testicular differentiation in Xenopus [70]. As in tilapia, miR-499 is also a dominant miRNA in swine, indicating its involvement in testicular development [71]. miR-9 and miR-192 are also upregulated in tilapia ovary, while miR-27d, miR-29b, miR-92b, miR-144-5p and miR-455 show the same pattern of greater expression in testes at later stage [43], suggesting they might play crucial roles in tilapia sex differentiation. However, the function of most other sex-biased miRNAs was previously unknown. More studies are needed to uncover the regulatory roles of these miRNAs.

By applying integrated miRNA and mRNA analysis together with target prediction, we reduced the complexity of predicted miRNA-mRNA relations on average by more than 50-fold compared to pure in silico prediction. Thus, we were able to explore the functional details of sex-biased miRNAs in early sex differentiation with the combined data. The transcriptomic data collected in the current work and a previous study [72] further confirmed the critical role of estrogen in sex determination and differentiation. The differentially expressed miRNAs targeted multiple genes, including those genes encoding key enzymes in the steroid hormone biosynthesis pathway (Cyp11a1, Hsd3b, Cyp19a1a, Hsd11b) [48] and key factors involved in early sex differentiation (Foxl2, Amh, Star1, Sf1, Dmrt1, and Gsdf) [46]. Cyp19a1a encodes the key enzyme responsible for estrogen synthesis. XX gonads with higher expression of Cyp19a1a exhibited lower expression of miR-30, which suggests a critical role of the miR-30 family in early sex differentiation via ovarian cell steroidogenesis as in human [73]. However, how the ovary is differentiated through Cyp19a1a expression being regulated by the miR-30 family remains to be investigated. Furthermore, Foxl2 and Dmrt1 play antagonistic roles in sex differentiation in Nile tilapia via regulation of Cyp19a1a expression and estrogen production [39]. Elevated expression of Foxl2 and decreased expression of miR-7977 in tilapia XX gonads indicates the possible role of miR-7977 in early sex differentiation. The higher expression of miR-212, which may repress the expression of Dmrt1 in the XX gonad, suggests its possible role as a post-transcriptional regulator in the differentiation of granulosa cells, as demonstrated in mouse [74]. Interestingly, some miRNAs (e.g., miR-96 and miR-737) target multiple genes involved in steroid synthesis, suggesting the complex regulatory network in fish sex determination and differentiation. Hence, these results suggest a role for miRNAs in regulating the biosynthesis of steroid hormones during tilapia early sex differentiation.

Conclusions

The present study is the first to examine the miRNA and mRNA expression profile of tilapia XX and XY gonads in early sex differentiation. We identified 111 known miRNAs that were differentially expressed between the XX and XY gonads. The predicted targets of these sex-biased miRNAs included genes encoding key enzymes in steroidogenic pathways and key molecules involved in vertebrate sex differentiation. These genes also showed sex-biased expression in tilapia gonads at 5 dah based on transcriptomic analysis. Some miRNA (e.g., miR-96 and miR-737) targeted multiple genes involved in steroid synthesis, indicating the complex regulatory network in fish early sex differentiation. An additional 130 novel miRNAs with differential expression in the XX and XY gonad were also identified. These miRNAs are an interesting starting point for future research to understand miRNA-mRNA interaction in early sex differentiation of teleosts. The availability of mono-sex fish and genome editing in tilapia makes it an attractive model system for study of miRNA function during sex determination and differentiation.

Methods

Prediction of miRNA genes using in silico approaches

A database of 1623 known teleost precursor miRNA (pre-miRNA) sequences was downloaded from miRBase release 21.0 [18]. To detect miRNA genes in tilapia, we conducted a blastn similarity search of these pre-miRNAs against the tilapia genome sequences [47] with an E-value cutoff of 0.001. The blastn hits were then manually inspected and compared with their query sequences in order to extract adjacent nucleotides (~30 bp) in both directions that might form part of the pre-miRNA. RNA secondary structure of the tilapia putative miRNA sequences was predicted using Mfold [75] to ensure proper stem–loop folding. A reciprocal blastn of the tilapia pre-miRNAs against known teleost pre-miRNAs was performed to annotate these predicted pre-miRNAs and to assign orthology (Additional file 12: Figure S6). Specifically, we compared the high confidence set of tilapia pre-miRNA sequences with those of cichlids [47] and other fishes from miRBase release 21.0 [18].

Ethics

All fish experiments were conducted in accordance with the regulations of the Guide for Care and Use of Laboratory Animals and were approved by the Committee of Laboratory Animal Experimentation at Southwest University.

Fish materials

The founder strain of the Nile tilapia, which was first introduced from Egypt in Africa, was obtained from Prof. Nagahama (Laboratory of Reproductive Biology, National Institute for Basic Biology, Okazaki, Japan). Fishes were reared in re-circulating aerated freshwater tanks at 26 °C prior to use. The breeding of mono-sex fish (all-female or all-male progenies) facilitates sampling of the gonads before the morphological differentiation. All-XX progenies were obtained by crossing a pseudomale (XX male, producing sperm after hormonal sex reversal) with a wild type female (XX). All-XY progenies were obtained by crossing a supermale (YY) with a wild type female. XX and XY fish were dissected to obtain samples of gonads for RNA isolation at the critical stage of molecular sex determination [46] (5 dah, 150 gonads were pooled for each sex).

RNA isolation

Total RNA was extracted from gonads of XX and XY fishes frozen in liquid nitrogen using the Trizol Reagent (Invitrogen, Carlsbad, CA) according to the manufacturer’s instruction. The extracted RNA was further treated with DNaseI (RNase-free 5U/uL) to eliminate genomic DNA contamination. The integrity and concentration of RNA was determined using both agarose gel electrophoresis and Nanodrop spectrophotometer (Nanodrop Technologies, Wilmington, Delaware USA), respectively. RNA with OD260/280 = 1.8–2.1, OD260/230 ≥ 1.7, c ≥ 150 ng/μL, and RNA integrity number (RIN) ≥ 7.0 and 28S/18S > 1.5 was stored at −80 °C for subsequent library construction and sequencing.

Illumina paired-end cDNA library construction and sequencing (RNA-Seq)

Two small RNA libraries and two mRNA libraries were constructed from purified RNA of tilapia XX and XY gonads. sRNA-Seq (small RNA sequencing) and mRNA-Seq (mRNA sequencing) services were carried out using Illumina HiSeq 2000 platform (Illumina Inc., San Diego, CA, USA) provided by BGI-Shenzhen, China. For sRNA-Seq, small RNAs of 16–32 nt in length were isolated from the total RNA by size fractionation in a 15 % TBE urea polyacrylamide gel, small RNAs were ligated to 3’ adaptor and 5’ adaptor with T4 RNA ligase. Subsequently, PCR reaction was performed using primers complementary to the two adaptors. The amplification products were further purified on 6 % polyacrylamide TBE gel and used for sequence analysis. In brief, Oligo(dT) beads were used to isolate poly(A) mRNA from total RNA for mRNA-seq. The average insert size in an mRNA sequencing library is approximately 200 bp and both ends of the library were sequenced. The reads of sRNA-Seq and mRNA-Seq have been deposited in NCBI’s Sequence Read Archive (SRA) database with accession numbers SRS1212947 and SRS1212948.

miRNA annotation and expression profiling

For sRNA-Seq, Trimmomatic [76] and a custom PERL script was used to remove the low-quality reads and trim the adaptor (5’-GCCTTGGCACCCGAGAATTCCA-3’). After adaptor trimming, reads of 18–32 nt in length were kept for further bioinformatics analysis. These reads were mapped to Nile tilapia genome with a tolerance of one mismatch in the seed sequence using Bowtie2 [77]. Sequences with perfect matches or one mismatch were retained for further analysis. Subsequently, by blasting against the Rfam database (12.0, http://rfam.xfam.org/), the reads mapped to the Nile tilapia genome were analyzed to discard the rRNA, tRNA, snoRNA, and ncRNA sequences. The remaining sequences were identified as the conserved miRNAs in Nile tilapia by a blast search against miRBase 21.0 allowing no more than one mismatch. The detailed mapping process was performed as previously described [17]. The orthologous miRNA matches were named according to the original gene name except the species name. The remaining reads were used for novel miRNA prediction. To further analyze the RNA secondary structures for novel miRNA prediction, 50 nts of genomic sequence flanking each side were extracted, and the secondary structures with regard to miRNA precursor gene characteristics were predicted using Mfold [75] (Additional file 13: Figure S7). miRNA read counts were normalized to tags per million (TPM) counts in order to compare miRNA expression between XX and XY gonads. The TPM was calculated as follows: normalized expression, TPM = (actual miRNA count/number of total clean read) × 1,000,000. The fold-change |log2Ratio| ≥ 1 and P-value (<0.01) were calculated to identify differentially expressed miRNAs.

Target prediction of the sex-biased miRNAs and expression profile of mRNAs

To provide a better understanding of the roles of the sex-biased miRNAs and correlations between expression levels of miRNA and mRNA, the mRNAs potentially targeted by these miRNAs were predicted using both miRanda [78] and Targetscan [79]. The 3’ UTRs for all genes of tilapia were predicted using Transdecoder [80], and were used as the primary base-pairing region of miRNAs. The prediction was based on the degree of miRNA and target sequence complementation, as well as the free energy level of RNA-RNA duplexes. Expression of mRNAs was qualified using transcriptomic analyses. Reads from mRNA-seq were aligned to the Nile tilapia reference sequence with TopHat2 [81]. NCBI RefSeq annotations were used to guide the Cufflinks [82] assembly, and Cuffdiff was used to determine FPKM values for those gene models. The results were subsequently filtered to exclude gene models whose FPKM value was less than 0.5 in either females or males. Additionally, when comparison between FPKM of the two samples was carried out, if the FPKM value exceeded 0.5 in one sex and was zero in the other sex, it was considered an undefined bias favoring the sex as the similar criteria in previous study [83]. Female upregulated and male upregulated gene models were counted exclusively against female downregulated and male downregulated miRNAs.

Availability of supporting data

The data sets supporting the results of this article are included within the article and its additional files.

Abbreviations

Amh :

Anti-Müllerian hormone; steroidogenic acute regulatory protein 1

Cyp11a1 :

cytochrome P450 family 11 subfamily A member 1

Cyp19a1a :

cytochrome P450, family 19, subfamily A, polypeptide 1a

Dmrt1 :

doublesex- and mab-3-related transcription factor-1

Foxl2 :

forkhead box L2

Gsdf :

gonadal somatic cell derived factor

Hsd11b :

11-Beta hydroxysteroid dehydrogenase

Hsd3b :

3 beta- and steroid delta-isomerase 3

Sf1 :

steroidogenic factor-1

TPM:

tags per million counts

References

  1. Bartel DP. microRNAs: genomics, biogenesis, mechanism, and function. Cell. 2004;116(2):281–97.

    Article  CAS  PubMed  Google Scholar 

  2. Shabalina SA, Koonin EV. Origins and evolution of eukaryotic RNA interference. Trends Ecol Evol. 2008;23(10):578–87.

    Article  PubMed  PubMed Central  Google Scholar 

  3. Bartel DP. microRNAs: target recognition and regulatory functions. Cell. 2009;136(2):215–33.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  4. Bushati N, Cohen SM. MicroRNA functions. Annu Rev Cell Dev Bi. 2007;23:175–205.

    Article  CAS  Google Scholar 

  5. Lagos-Quintana M, Rauhut R, Yalcin A, Meyer J, Lendeckel W, Tuschl T. Identification of tissue-specific microRNAs from mouse. Curr Biol. 2002;12(9):735–9.

    Article  CAS  PubMed  Google Scholar 

  6. Stefani G, Slack FJ. Small non-coding RNAs in animal development. Nat Rev Mol Cell Bio. 2008;9(3):219–30.

    Article  CAS  Google Scholar 

  7. Krol J, Loedige I, Filipowicz W. The widespread regulation of microRNA biogenesis, function and decay. Nat Rev Genet. 2010;11(9):597–610.

    CAS  PubMed  Google Scholar 

  8. Ambros V, Chen XM. The regulation of genes and genomes by small RNAs. Development. 2007;134(9):1635–41.

    Article  CAS  PubMed  Google Scholar 

  9. Friedman RC, Farh KKH, Burge CB, Bartel DP. Most mammalian mRNAs are conserved targets of microRNAs. Genome Res. 2009;19(1):92–105.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Lim LP, Glasner ME, Yekta S, Burge CB, Bartel DP. Vertebrate microRNA genes. Science. 2003;299(5612):1540.

    Article  CAS  PubMed  Google Scholar 

  11. Xu ZQ, Chen JP, Li XG, Ge JC, Pan JL, Xu XF. Identification and characterization of microRNAs in channel catfish (Ictalurus punctatus) by using solexa sequencing technology. Plos One. 2013;8(1):e541747.

    Google Scholar 

  12. Bekaert M, Lowe NR, Bishop SC, Bron JE, Taggart JB, Houston RD. Sequencing and characterisation of an extensive Atlantic salmon (Salmo salar L.) microRNA repertoire. Plos One. 2013;8(7):e70136.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  13. Andreassen R, Worren MM, Hoyheim B. Discovery and characterization of miRNA genes in Atlantic salmon (Salmo salar) by use of a deep sequencing approach. BMC Genomics. 2013;14:482.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Vaz C, Wee CW, Lee GP, Ingham PW, Tanavde V, Mathavan S. Deep sequencing of small RNA facilitates tissue and sex associated microRNA discovery in zebrafish. BMC Genomics. 2015;16(1):950.

    Article  PubMed  PubMed Central  Google Scholar 

  15. Wongwarangkana C, Fujimori KE, Akiba M, Kinoshita S, Teruya M, Nezuo M, Masatoshi T, Watabe S, Asakawa S. Deep sequencing, profiling and detailed annotation of microRNAs in Takifugu rubripes. BMC Genomics. 2015;16:457.

    Article  PubMed  PubMed Central  Google Scholar 

  16. Ma H, Hostuttler M, Wei HR, Rexroad CE, Yao JB. Characterization of the rainbow trout egg microRNA transcriptome. PLoS One. 2012;7(6):e39649.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Xiao J, Zhong H, Zhou Y, Yu F, Gao Y, Luo YJ, Tang ZY, Guo ZB, Guo EY, Gan X, et al. Identification and characterization of microRNAs in ovary and testis of Nile tilapia (Oreochromis niloticus) by using solexa sequencing technology. PLoS One. 2014;9(1):e86821.

    Article  PubMed  PubMed Central  Google Scholar 

  18. Kozomara A, Griffiths-Jones S. miRBase: annotating high confidence microRNAs using deep sequencing data. Nucleic Acids Res. 2014;42:68–73.

    Article  Google Scholar 

  19. Godwin J, Luckenbach JA, Borski RJ. Ecology meets endocrinology: environmental sex determination in fishes. Evol Dev. 2003;5(1):40–9.

    Article  PubMed  Google Scholar 

  20. Martinez P, Vinas AM, Sanchez L, Diaz N, Ribas L, Piferrer F. Genetic architecture of sex determination in fish: applications to sex ratio control in aquaculture. Front Genet. 2014;5:340.

    PubMed  PubMed Central  Google Scholar 

  21. Baroiller JF, D'Cotta H. Environment and sex determination in farmed fish. Comp Biochem Physiol C Toxicol Pharmacol. 2001;130(4):399–409.

    Article  CAS  PubMed  Google Scholar 

  22. Cutting A, Chue J, Smith CA. Just how conserved is vertebrate sex determination? Dev Dynam. 2013;242(4):380–7.

    Article  CAS  Google Scholar 

  23. Matsuda M, Nagahama Y, Shinomiya A, Sato T, Matsuda C, Kobayashi T, Morrey CE, Shibata N, Asakawa S, Shimizu N, et al. DMY is a Y-specific DM-domain gene required for male development in the medaka fish. Nature. 2002;417(6888):559–63.

    Article  CAS  PubMed  Google Scholar 

  24. Nanda I, Kondo M, Hornung U, Asakawa S, Winkler C, Shimizu A, Shan ZH, Haaf T, Shimizu N, Shima A, et al. A duplicated copy of DMRT1 in the sex-determining region of the Y chromosome of the medaka, Oryzias latipes. P Natl Acad Sci USA. 2002;99(18):11778–83.

    Article  CAS  Google Scholar 

  25. Myosho T, Otake H, Masuyama H, Matsuda M, Kuroki Y, Fujiyama A, Naruse K, Hamaguchi S, Sakaizumi M. Tracing the Emergence of a Novel Sex-Determining Gene in Medaka, Oryzias luzonensis. Genetics. 2012;191(1):163–70.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Takehana Y, Matsuda M, Myosho T, Suster ML, Kawakami K, Shin-I T, Kohara Y, Kuroki Y, Toyoda A, Fujiyama A, et al. Co-option of Sox3 as the male-determining factor on the Y chromosome in the fish Oryzias dancena. Nat Commun. 2014;5:4157.

    Article  CAS  PubMed  Google Scholar 

  27. Hattori RS, Murai Y, Oura M, Masuda S, Majhi SK, Sakamoto T, Fernandino JI, Somoza GM, Yokota M, Strussmann CA. A Y-linked anti-Mullerian hormone duplication takes over a critical role in sex determination. Proc Natl Acad Sci U S A. 2012;109(8):2955–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Li M, Sun Y, Zhao J, Shi H, Zeng S, Ye K, Jiang D, Zhou L, Sun L, Tao W, et al. A tandem duplicate of anti-Mullerian hormone with a missense SNP on the Y chromosome is essential for male sex determination in Nile tilapia, Oreochromis niloticus. PLoS Genet. 2015;11(11):e1005678.

    Article  PubMed  PubMed Central  Google Scholar 

  29. Kamiya T, Kai W, Tasumi S, Oka A, Matsunaga T, Mizuno N, Fujita M, Suetake H, Suzuki S, Hosoya S, et al. A trans-species missense SNP in Amhr2 is associated with sex determination in the tiger pufferfish, Takifugu rubripes (Fugu). PLoS Genet. 2012;8(7):e1002798.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  30. Reichwald K, Petzold A, Koch P, Downie BR, Hartmann N, Pietsch S, Baumgart M, Chalopin D, Felder M, Bens M, et al. Insights into sex chromosome evolution and aging from the genome of a short-lived fish. Cell. 2015;163(6):1527–38.

    Article  CAS  PubMed  Google Scholar 

  31. Yano A, Guyomard R, Nicol B, Jouanno E, Quillet E, Klopp C, Cabau C, Bouchez O, Fostier A, Guiguen Y. An immune-related gene evolved into the master sex-determining gene in rainbow trout, Oncorhynchus mykiss. Curr Biol. 2012;22(15):1423–8.

    Article  CAS  PubMed  Google Scholar 

  32. Yano A, Nicol B, Jouanno E, Quillet E, Fostier A, Guyomard R, Guiguen Y. The sexually dimorphic on the Y-chromosome gene (sdY) is a conserved male-specific Y-chromosome sequence in many salmonids. Evol Appl. 2013;6(3):486–96.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Devlin RH, Nagahama Y. Sex determination and sex differentiation in fish: an overview of genetic, physiological, and environmental influences. Aquaculture. 2002;208(3–4):191–364.

    Article  CAS  Google Scholar 

  34. Trant JM, Gavasso S, Ackers J, Chung BC, Place AR. Developmental expression of cytochrome P450 aromatase genes (CYP19a and CYP19b) in zebrafish fry (Danio rerio). J Exp Zool. 2001;290(5):475–83.

    Article  CAS  PubMed  Google Scholar 

  35. Pannetier M, Fabre S, Batista F, Kocer A, Renault L, Jolivet G, Mandon-Pepin B, Cotinot C, Veltla R, Pailhoux E. FOXL2 activates P450 aromatase gene transcription: towards a better characterization of the early steps of mammalian ovarian development. J Mol Endocrinol. 2006;36(3):399–413.

    Article  CAS  PubMed  Google Scholar 

  36. Jorgensen A, Morthorst JE, Andersen O, Rasmussen LJ, Bjerregaard P. Expression profiles for six zebrafish genes during gonadal sex differentiation. Reprod Biol Endocrinol. 2008;6:25.

    Article  PubMed  PubMed Central  Google Scholar 

  37. Wang DS, Kobayashi T, Zhou LY, Paul-Prasanth B, Ijiri S, Sakai F, Okubo K, Morohashi KI, Nagahama Y. Foxl2 up-regulates aromatase gene transcription in a female-specific manner by binding to the promoter as well as interacting with Ad4 binding protein/steroidogenic factor 1. Mol Endocrinol. 2007;21(3):712–25.

    Article  CAS  PubMed  Google Scholar 

  38. Wang DS, Zhou LY, Kobayashi T, Matsuda M, Shibata Y, Sakai F, Nagahama Y. Doublesex- and Mab-3-related transcription factor-1 repression of aromatase transcription, a possible mechanism favoring the male pathway in tilapia. Endocrinology. 2010;151(3):1331–40.

    Article  CAS  PubMed  Google Scholar 

  39. Li MH, Yang HH, Li MR, Sun YL, Jiang XL, Xie QP, Wang TR, Shi HJ, Sun LN, Zhou LY, et al. Antagonistic roles of Dmrt1 and Foxl2 in sex differentiation via estrogen production in tilapia as demonstrated by TALENs. Endocrinology. 2013;154(12):4814–25.

    Article  CAS  PubMed  Google Scholar 

  40. Juanchich A, Le Cam A, Montfort J, Guiguen Y, Bobe J. Identification of differentially expressed miRNAs and their potential targets during fish ovarian development. Biol Reprod. 2013;88(5):128.

    Article  PubMed  Google Scholar 

  41. Jing J, Wu JJ, Liu W, Xiong ST, Ma WG, Zhang J, Wang WM, Gui JF, Mei J. Sex-biased miRNAs in gonad and their potential roles for testis development in yellow catfish. PLoS One. 2014;9(9):e107946.

    Article  PubMed  PubMed Central  Google Scholar 

  42. Eshel O, Shirak A, Dor L, Band M, Zak T, Markovich-Gordon M, Chalifa-Caspi V, Feldmesser E, Weller JI, Seroussi E, et al. Identification of male-specific amh duplication, sexually differentially expressed genes and microRNAs at early embryonic development of Nile tilapia (Oreochromis niloticus). BMC Genomics. 2014;15:774.

    Article  PubMed  PubMed Central  Google Scholar 

  43. Wang W, Liu W, Liu Q, Li B, An L, Hao R, Zhao J, Liu S, Song J. Coordinated microRNA and messenger RNA expression profiles for understanding sexual dimorphism of gonads and the potential roles of microRNA in the steroidogenesis pathway in Nile tilapia (Oreochromis niloticus). Theriogenology. 2016;85(5):970–8.

    Article  CAS  PubMed  Google Scholar 

  44. Bizuayehu TT, Babiak J, Norberg B, Fernandes JMO, Johansen SD, Babiak I. Sex-biased miRNA expression in Atlantic Halibut (Hippoglossus hippoglossus) brain and gonads. Sex Dev. 2012;6(5):257–66.

    Article  CAS  PubMed  Google Scholar 

  45. Muniategui A, Pey J, Planes FJ, Rubio A. Joint analysis of miRNA and mRNA expression data. Brief Bioinform. 2013;14(3):263–78.

    Article  CAS  PubMed  Google Scholar 

  46. Ijiri S, Kaneko H, Kobayashi T, Wang DS, Sakai F, Paul-Prasanth B, Nakamura M, Nagahama Y. Sexual dimorphic expression of genes in gonads during early differentiation of a teleost fish, the Nile tilapia Oreochromis niloticus. Biol Reprod. 2008;78(2):333–41.

    Article  CAS  PubMed  Google Scholar 

  47. Brawand D, Wagner CE, Li YI, Malinsky M, Keller I, Fan SH, Simakov O, Ng AY, Lim ZW, Bezault E, et al. The genomic substrate for adaptive radiation in African cichlid fish. Nature. 2014;513(7518):375–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Nagahama Y. 17 alpha,20 beta-dihydroxy-4-pregnen-3-one, a maturation-inducing hormone in fish oocytes: mechanisms of synthesis and action. Steroids. 1997;62(1):190–6.

    Article  CAS  PubMed  Google Scholar 

  49. Loh YHE, Yi SV, Streelman JT. Evolution of microRNAs and the diversification of species. Genome Biol Evol. 2011;3:55–65.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Li XH, Wu JS, Tang LH, Hu D. Identification of conserved microRNAs and their target genes in Nile tilapia (Oreochromis niloticus) by bioinformatic analysis. Genet Mol Res. 2015;14(1):2785–92.

    Article  CAS  PubMed  Google Scholar 

  51. Mathelier A, Carbone A. Large scale chromosomal mapping of human microRNA structural clusters. Nucleic Acids Res. 2013;41(8):4392–408.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Calin GA, Sevignani C, Dan Dumitru C, Hyslop T, Noch E, Yendamuri S, Shimizu M, Rattan S, Bullrich F, Negrini M, et al. Human microRNA genes are frequently located at fragile sites and genomic regions involved in cancers. Proc Natl Acad Sci U S A. 2004;101(9):2999–3004.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Sevignani C, Calin GA, Nnadi SC, Shimizu M, Davuluri RV, Hyslop T, Demant P, Croce CM, Siracusa LD. MicroRNA genes are frequently located near mouse cancer susceptibility loci. Proc Natl Acad Sci U S A. 2007;104(19):8017–22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Ghorai A, Ghosh U. miRNA gene counts in chromosomes vary widely in a species and biogenesis of miRNA largely depends on transcription or post-transcriptional processing of coding genes. Front Genet. 2014;5:100.

    Article  PubMed  PubMed Central  Google Scholar 

  55. Wu J, Wang DF, Liu YF, Wang L, Qiao X, Zhang SL. Identification of miRNAs involved in pear fruit development and quality. BMC Genomics. 2014;15:953.

    Article  PubMed  PubMed Central  Google Scholar 

  56. Ro S, Song R, Park C, Zheng H, Sanders KM, Yan W. Cloning and expression profiling of small RNAs expressed in the mouse ovary. RNA. 2007;13(12):2366–80.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  57. Houwing S, Kamminga LM, Berezikov E, Cronembold D, Girard A, van den Elst H, Filippov DV, Blaser H, Raz E, Moens CB, et al. A role for Piwi and piRNAs in germ cell maintenance and transposon silencing in zebrafish. Cell. 2007;129(1):69–82.

    Article  CAS  PubMed  Google Scholar 

  58. Liang SQ, Chen L, Huang HL, Zhi DS. The experimental study of miRNA in pituitary adenomas. Turk Neurosurg. 2013;23(6):721–7.

    PubMed  Google Scholar 

  59. Gangaraju VK, Lin HF. MicroRNAs: key regulators of stem cells. Nat Rev Mol Cell Bio. 2009;10(2):116–25.

    Article  CAS  Google Scholar 

  60. Li MZ, Liu YK, Wang T, Guan JQ, Luo ZG, Chen HS, Wang X, Chen L, Ma JD, Mu ZP, et al. Repertoire of porcine microRNAs in adult ovary and testis by deep sequencing. Int J Biol Sci. 2011;7(7):1045–55.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Saunders LR, Sharma AD, Tawney J, Nakagawa M, Okita K, Yamanaka S, Willenbring H, Verdin E. miRNAs regulate SIRT1 expression during mouse embryonic stem cell differentiation and in adult mouse tissues. Aging-Us. 2010;2(7):415–31.

    Article  CAS  Google Scholar 

  62. Sirotkin AV, Ovcharenko D, Grossmann R, Laukova M, Mlyncek M. Identification of microRNAs controlling human ovarian cell steroidogenesis via a genome-scale screen. J Cell Physiol. 2009;219(2):415–20.

    Article  CAS  PubMed  Google Scholar 

  63. Chiang HR, Schoenfeld LW, Ruby JG, Auyeung VC, Spies N, Baek D, Johnston WK, Russ C, Luo S, Babiarz JE, et al. Mammalian microRNAs: experimental evaluation of novel and previously annotated genes. Genes Dev. 2010;24(10):992–1009.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  64. Huang JM, Ju ZH, Li QL, Hou QL, Wang CF, Li JB, Li RL, Wang LL, Sun T, Hang SQ, et al. Solexa Sequencing of Novel and Differentially Expressed MicroRNAs in Testicular and Ovarian Tissues in Holstein Cattle. Int J Biol Sci. 2011;7(7):1016–26.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  65. Lau K, Lai KP, Bao JYJ, Zhang N, Tse A, Tong A, Li JW, Lok S, Kong RYC, Lui WY, et al. Identification and expression profiling of microRNAs in the brain, liver and gonads of marine medaka (Oryzias melastigma) and in response to hypoxia. PLoS One. 2014;9(10):e110698.

    Article  PubMed  PubMed Central  Google Scholar 

  66. Kang L, Cui X, Zhang Y, Yang C, Jiang Y. Identification of miRNAs associated with sexual maturity in chicken ovary by Illumina small RNA deep sequencing. BMC Genomics. 2013;14:352.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  67. Christenson LK. MicroRNA control of ovarian function. Anim Reprod. 2010;7(3):129–33.

    PubMed  PubMed Central  Google Scholar 

  68. Luesink M, Nigten J, Knops RHJN, de Witte TJM, van der Reijden BA, Jansen JH. MiR-132 and MiR-212 as regulators of WT1 and GATA2 in acute myeloid leukemia. Blood. 2009;114(22):526.

    Google Scholar 

  69. Rakoczy J, Fernandez-Valverde SL, Glazov EA, Wainwright EN, Sato T, Takada S, Combes AN, Korbie DJ, Miller D, Grimmond SM, et al. MicroRNAs-140-5p/140-3p modulate Leydig cell numbers in the developing mouse testis. Biol Reprod. 2013;88(6):143.

    Article  PubMed  Google Scholar 

  70. Kraggerud SM, Hoei-Hansen CE, Alagaratnam S, Skotheim RI, Abeler VM, Rajpert-De Meyts E, Lothe RA. Molecular characteristics of malignant ovarian germ cell tumors and comparison with testicular counterparts: implications for pathogenesis. Endocr Rev. 2013;34(3):339–76.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Zhang XJ, Li CM, Liu X, Lu CY, Bai CY, Zhao ZH, Sun BX. Differential expression of miR-499 and validation of predicted target genes in the testicular tissue of swine at different developmental stages. DNA Cell Biol. 2015;34(7):464–9.

    Article  CAS  PubMed  Google Scholar 

  72. Tao W, Yuan J, Zhou L, Sun L, Sun Y, Yang S, Li M, Zeng S, Huang B, Wang D. Characterization of gonadal transcriptomes from Nile tilapia (Oreochromis niloticus) reveals differentially expressed genes. PLoS One. 2013;8(5):e63604.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Sirotkin AV, Laukova M, Ovcharenko D, Brenaut P, Mlyncek M. Identification of MicroRNAs Controlling Human Ovarian Cell Proliferation and Apoptosis. J Cell Physiol. 2010; 223(1):49–56.

    CAS  PubMed  Google Scholar 

  74. Fiedler SD, Carletti MZ, Hong X, Christenson LK. Hormonal regulation of MicroRNA expression in periovulatory mouse mural granulosa cells. Biol Reprod. 2008;79(6):1030–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

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

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  76. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30(15):2114–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  77. Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012;9(4):357–U354.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

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

    Article  PubMed  PubMed Central  Google Scholar 

  79. 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.

    Article  CAS  PubMed  Google Scholar 

  80. Haas BJ, Papanicolaou A, Yassour M, Grabherr M, Blood PD, Bowden J, Couger MB, Eccles D, Li B, Lieber M, et al. De novo transcript sequence reconstruction from RNA-seq using the Trinity platform for reference generation and analysis. Nat Protoc. 2013;8(8):1494–512.

    Article  CAS  PubMed  Google Scholar 

  81. Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013;14(4):R36.

    Article  PubMed  PubMed Central  Google Scholar 

  82. Trapnell C, Hendrickson DG, Sauvageau M, Goff L, Rinn JL, Pachter L. Differential analysis of gene regulation at transcript resolution with RNA-seq. Nat Biotechnol. 2013;31(1):46–53.

    Article  CAS  PubMed  Google Scholar 

  83. Gammerdinger WJ, Conte MA, Acquah EA, Roberts RB, Kocher TD. Structure and decay of a proto-Y region in Tilapia, Oreochromis niloticus. BMC Genomics. 2014;15:975.

    Article  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

This work was supported by grants 91331119, 31502147, 31572609 and 31030063 from the National Natural Science Foundation of China; grant 2011AA100404 from the National High Technology Research and Development Program (863 program) of China; grant 2012CB723205 from the National Basic Research Program of China; grant 20130182130003 from the Specialized Research Fund for the Doctoral Program of Higher Education of China; grant cstc2013kjrc-tdjs80003 and cstc2014jcyjA80001 from the Natural Science Foundation Project of Chongqing, Chongqing Science and Technology Commission; grant XDJK2016B011, XDJK2016A003, XDJK2014B040 and XDJK2015A004 from the Fundamental Research Funds for the Central Universities.

Author information

Author notes

    Authors and Affiliations

    1. Key Laboratory of Freshwater Fish Reproduction and Development (Ministry of Education), Key Laboratory of Aquatic Science of Chongqing, School of Life Sciences, Southwest University, Chongqing, 400715, P. R. China

      Wenjing Tao, Lina Sun, Hongjuan Shi, Yunying Cheng, Dongneng Jiang & Deshou Wang

    2. State Key Laboratory of Freshwater Ecology and Biotechnology, Institute of Hydrobiology, Chinese Academy of Sciences, Wuhan, P. R. China

      Beide Fu

    3. Department of Biology, University of Maryland, College Park, MD, USA

      Matthew A. Conte, William J. Gammerdinger & Thomas D. Kocher

    Authors
    1. Wenjing Tao
      View author publications

      You can also search for this author in PubMed Google Scholar

    2. Lina Sun
      View author publications

      You can also search for this author in PubMed Google Scholar

    3. Hongjuan Shi
      View author publications

      You can also search for this author in PubMed Google Scholar

    4. Yunying Cheng
      View author publications

      You can also search for this author in PubMed Google Scholar

    5. Dongneng Jiang
      View author publications

      You can also search for this author in PubMed Google Scholar

    6. Beide Fu
      View author publications

      You can also search for this author in PubMed Google Scholar

    7. Matthew A. Conte
      View author publications

      You can also search for this author in PubMed Google Scholar

    8. William J. Gammerdinger
      View author publications

      You can also search for this author in PubMed Google Scholar

    9. Thomas D. Kocher
      View author publications

      You can also search for this author in PubMed Google Scholar

    10. Deshou Wang
      View author publications

      You can also search for this author in PubMed Google Scholar

    Corresponding authors

    Correspondence to Thomas D. Kocher or Deshou Wang.

    Additional information

    Competing interests

    The authors declare that they have no competing interests.

    Authors’ contributions

    This study was designed by DSW, TDK, and WJT and organized by DSW. DNJ managed the experimental fish, YYC, HJS, and LNS dissected the gonads at 5 dah, WJT and YYC carried out RNA extractions. WJT, BDF, MAC and WJG carried out the bioinformatics analyses. WJT, LNS, DSW, and TDK wrote the manuscript and all authors critiqued the manuscript for important intellectual content. All authors read and approved the final manuscript.

    Additional files

    Additional file 1: Table S1.

    Sequence of identified 279 conserved pre-miRNA genes in tilapia genome. (XLS 66 kb)

    Additional file 2: Figure S1.

    Comparisons of the number of pre-miRNAs of Paralichthys olivaceus, Hippoglossus hippoglossus, Fugu rubripe, Tetraodon nigroviridis, Cyprinus carpio, Oryzias latipes, Ictalurus punctatus, Danio rerio and Salmo salar in miRbase with the number of known miRNAs identified in Oreochromis niloticus. (TIF 486 kb)

    Additional file 3: Figure S2.

    Pre-miRNAs sequence alignments of O. niloticus, Astatotilapia burtoni, Neolamprologus brichardi, Pundamilia nyererei, Metriaclima zebra, O. latipes, G. aculeatus and D. rerio. A. pre-miR-140, B. pre-miR-499. (TIF 1897 kb)

    Additional file 4: Table S2.

    Sequence and expression profiles of known miRNAs in the ovary and testis of tilapia. (XLS 1613 kb)

    Additional file 5: Figure S3.

    Comparisons of the number of mature miRNAs in Danio rerio, Fugu rubripe, Tetraodon nigroviridis, and Oryzias latipe in miRbase with the number of known miRNAs identified in the gonads of O. niloticus. (TIF 434 kb)

    Additional file 6: Figure S4.

    Venn diagram of identified miRNAs in tilapia gonads in O. niloticus, A. burtoni, N. brichardi, P. nyererei and M. zebra published in previous study [47]. (TIF 369 kb)

    Additional file 7: Table S3.

    Sequence and expression profiles of novel miRNAs with in the ovary and testis of tilapia. (XLS 45 kb)

    Additional file 8: Table S4.

    Predicted targets of female upregulated miRNAs validated by mRNA sequencing. (XLS 1765 kb)

    Additional file 9: Table S5.

    Predicted targets of male upregulated miRNAs validated by mRNA sequencing. (XLS 2205 kb)

    Additional file 10: Figure S5.

    Target prediction for selected miRNA expressed in ovaries and testis in tilapia. (A) Cyp19a1a, (B) Dmrt1, and (C) Foxl2 are predicted as potential targets of miR-30a, miR-212 and miR-7977, respectively. (TIF 1109 kb)

    Additional file 11: Table S6.

    miRNAs and mRNAs in the steroid hormone biosynthesis pathway and involved in early sex differentiation. (XLS 29 kb)

    Additional file 12: Figure S6.

    Flowchart depicting the workflow for pre-miRNA prediction. A blastn similarity search of these pre-miRNAs against the tilapia genomic sequences [47] with an E-value cutoff of 0.001. The blastn hits were then manually inspected and compared with their query sequences in order to extract adjacent nucleotides (~30 bp) in both directions that might form part of the pre-miRNA. RNA secondary structure of the cichlid putative miRNA sequences was predicted using Mfold [75] to ensure proper stem–loop folding. A reciprocal blastn of the tilapia pre-miRNAs against known teleost pre-miRNAs was performed to identify the tilapia miRNA and to assign orthology. (TIF 1123 kb)

    Additional file 13: Figure S7.

    Flowchart depicting the workflow for miRNA profiling, annotation and novel miRNA prediction. The workflow comprises of five parts namely: filter of the reads, quantification of known miRNAs, the novel miRNA prediction pipeline, target prediction, quantification of mRNA. The quantification of known miRNAs and prediction of novel miRNAs was done using Mfold, and the quantification of mRNAs was performed by Tophat2 [81] and Cufflinks [82]. (TIF 1823 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

    Check for updates. Verify currency and authenticity via CrossMark

    Cite this article

    Tao, W., Sun, L., Shi, H. et al. Integrated analysis of miRNA and mRNA expression profiles in tilapia gonads at an early stage of sex differentiation. BMC Genomics 17, 328 (2016). https://doi.org/10.1186/s12864-016-2636-z

    Download citation

    • Received:

    • Accepted:

    • Published:

    • DOI: https://doi.org/10.1186/s12864-016-2636-z

    Keywords

    BMC Genomics

    ISSN: 1471-2164

    Contact us