1. Introduction
Sex determination and sex differentiation are a classic and challenging topic that has fascinated researchers for generations [1]. Given that sexual dimorphism in growth, size, color, and other economic characteristics is prevalent, sex is also one of the most valuable traits in aquaculture animals [2,3]. For fish with sexual dimorphism in growth, sex-controlled breeding to produce monosexual populations can increase economic effectiveness [4]. Therefore, research shedding light on the molecular mechanisms of sex differentiation in teleosts has the potential to bring enormous value to the aquaculture industry.
Non-coding RNAs (ncRNAs), a class of transcribed RNA molecules that do not encode proteins, perform crucial roles in gene expression regulation and various biological processes [5]. Small RNAs (smaRNAs), long non-coding RNAs (lncRNAs), and circular RNAs (circRNAs) are the most common types of ncRNAs. The smaRNAs can be classified into several categories, including microRNAs (miRNAs), small interfering RNAs (siRNAs), piwi interacting RNAs (piRNAs), and others [6]. Of these, microRNAs (miRNAs) with a size of 21–25 nt are the most known and have been extensively studied, and they negatively regulate gene expression by pairing with the corresponding target mRNA sequences [7]. They play crucial roles in gonadal development and gametogenesis by inhibiting classically genes related to sex differentiation. For example, miR-124 inhibits the expression of sox9 in the ovaries [8], and miR-107, highly expressed in the testes, mediates the estrogen signaling pathway by directly targeting cyp19a1a [9]. The lncRNAs are long ncRNAs with a length over 200 nt [10]. The lncRNA TESCO interacts with the transcription factor SF1, which sustains the expression of Sox9 and testicular development [11]. The lncRNA DMR plays a role in sexual differentiation by negatively regulating the expression of Dmrt1 [12]. The circRNAs are characterized by a covalently closed-loop structure without terminal 5′ caps and 3′ polyadenylated tails [13,14]. They can act as miRNA sponges to regulate gene expression by preventing them from binding to their target mRNAs [15]. However, there are only a few studies on the role of circRNAs in sex differentiation of the teleost.
The common carp (Cyprinus carpio) is a XX/XY sex-determination system Cyprinid fish [16]. It is widely distributed worldwide and has become one of the most critical commercially cultivated freshwater fish species due to its adaptation to complex and changing environments and rapid growth rate [17]. The global production of common carp was about 4.41 million metric tons, accounting for roughly 5.2% of the total freshwater fishery production in 2021 [18]. More importantly, the common carp shows apparent sexual dimorphism in growth, with females growing faster than males [19]. Studying the molecular mechanisms of sexual differentiation will help conduct sex-controlled breeding. The chromosome-level whole-genome sequence of common carp is available [20,21]. Although many studies have been performed to analyze the differential expression (DE) genes between the testes and ovaries [18,22,23], the function of ncRNAs in sexual differentiation has not been reported. There is room for improvement in functional genomics research in common carp, and more omics data are required.
In the present study, we analyzed the expression profiles of ncRNAs in the ovaries and testes. We also analyzed ncRNA-mRNA interactions to investigate the potential role of competing endogenous RNA (ceRNA) crosstalk in sexual differentiation and maintenance. Our findings will contribute to the understanding of the genetic basis of sex differentiation in common carp.
2. Materials and Methods
2.1. Ethics Statement
All experiments complied with the principles of animal care and use for scientific purposes established by the Animal Care and Use Committee of the Chinese Academy of Fishery Sciences (Protocol code ACUC-CAFS-20220615). In addition, the animals were sedated with MS222 (40 mg/L) to minimize discomfort before surgery.
2.2. Sample Collection and RNA Extraction
The common carp used in this study were reared at the Fangshan Experimental Station in Beijing, China. Ovary and testis tissues, with three replications in each group, were dissected from sexually mature carp and then frozen in liquid nitrogen. All samples were stored in a −80 °C refrigerator for later use. Total RNA was extracted from ovaries and testes using TRIzol reagent (Invitrogen, Carlsbad, CA, USA). Briefly, the sample was fully homogenized in TRIzol and incubated at room temperature for 5 min. Next, we added chloroform, shook the mix thoroughly, incubated it for 3 min, and centrifuged it for 15 min at 12,000× g at 4 °C. The aqueous phase was transferred to a new tube and treated with DNase I to digest the residual DNA. Total RNA precipitate formed a white pellet after adding isopropanol and centrifugation. The concentration of total RNA was determined by using Nanodrop One spectrophotometer (Thermo Fisher, Waltham, MA, USA), and the quality was checked by using Agilent Bioanalyzer 2100 system.
2.3. Library Preparation and Sequencing
High-quality RNAs with the NanoDrop 260/280 ratios of approximately 2.0, RNA Integrity Number (RIN) ≥7.5 for lncRNA and circRNA, and RIN ≥8.5 for smaRNA were used for further library construction. After removing ribosomal RNAs, lncRNA libraries were generated using NEBNext® UltraTM RNA Library Prep Kit for Illumina® (NEB, Newton, MA, USA). Before constructing the circRNA libraries, the linear RNAs needed to be digested with RNase R (Epicentre, Verona, WI, USA), and the circRNA libraries were generated with the NEBNext® UltraTM Directional RNA Library Prep Kit for Illumina® (NEB). For smaRNA libraries, the 3′ and 5′ adaptors were ligated to the 3′ and 5′ ends of the smaRNAs, respectively, and then smaRNA libraries were constructed using NEB Next® Multiplex SmaRNA Library Prep Set for Illumina®. The Agilent 2100 Bioanalyzer system was used to assess the quality of libraries. Then, lncRNA and circRNA libraries were sequenced on the Illumina NovaSeq 6000 platform to produce 150 bp paired-end reads. The smaRNA libraries were sequenced on the Illumina Hiseq 2500 platform (Novogene, Beijing, China) to produce 50 bp single-end reads.
2.4. Data Quality Control and Read Alignment
Clean data were obtained by removing reads containing adapters, poly-N, and low-quality reads from the raw data. The common carp reference genome was downloaded from NCBI [21]. Paired-end clean reads generated from lncRNA and circRNA sequencing were aligned to the reference genome using Hisat2 (v2.0.5). For smaRNA analysis, clean reads of the length between 18 and 35 nt were selected for mapping the reference genome using Bowtie (v2.4.0) [24].
2.5. Bioinformatics Analysis of smaRNA Sequences
The reads mapped to the reference genome were aligned with the miRNA sequences in miRbase (v22.1) [25] to identify known miRNAs. The miREvo [26] and mir-deep2 [27] were used to predict novel miRNA based on the secondary hairpin structure of miRNA precursors. The target gene of miRNA was predicted using miRanda. The miRNA editing analysis was conducted using the MiRME method [28]. For miRNA family analysis, miFam.dat (
2.6. Processing of lncRNA Sequencing Data
StringTie (v1.3.1) was used to reassemble and quantify transcriptome based on the reads mapping to the reference genome to obtain better accuracy and integrity transcripts [29]. Three coding potential analysis software programs, including Coding-Non-Coding Index (CNCI) (v2.0) [30], Coding Potential Calculator (CPC) (v3.2.0) [31], and Pfam ScanC (PFAM) (v1.6) [32], were combined to screen non-protein-coding RNA candidates. The candidates that were not annotated in NCBI were predicted as lncRNAs. Target genes of lncRNAs were predicted by co-location and co-expression analysis between lncRNAs and protein-coding genes, and the threshold for co-location analysis was set as 100 kb upstream or downstream of the lncRNAs [33,34]. rMATS (v4.1.0) software was used for AS event analysis with default parameters.
2.7. Identification of circRNAs and Prediction of Coding Potential
The circRNAs were identified using find_circ (v1.0) [35] and CIRI (v2.0.5) [36] using reads mapping to the reference genome. The circRNAs with internal ribosome entry sites (IRES) have the potential to be translated into proteins. The IRES elements of circRNAs were first predicted using the IRESfinder software (v1.0), and then CNCI, CPC, and PFAM were combined to validate the coding potential of the circRNA sequences [37].
2.8. Differential Expression Analysis
Fragments Per Kilobase of transcript sequence per Millions of base pairs sequenced (FPKM) for lncRNAs and mRNAs and transcripts per millions (TPM) for miRNA and circRNA were used to normalize the expression levels. Differential expression analysis was performed using the edgeR R package (v3.22.5) based on FPKM values and DESeq2 R package (v1.20.0) based on TPM values [38,39]. Corrected p-values ≤ 0.05 were set as the threshold to mine expressed transcripts.
2.9. Enrichment Analysis of Differential Candidate Genes
The differential candidate genes were submitted for functional enrichment analysis using the ClusterProfiler R package (v3.9), which included Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway terms [40]. p-Values lower than 0.05 are considered significant, and the top terms were selected for visualization.
2.10. Quantitative PCR
The expression level of mRNA was detected by stem-loop quantitative PCR [41]. The complementary DNA (cDNA) for miRNA real-time PCR was reverse-transcribed using the GoScript™ Reverse Transcription System (Promega, Madison, WI, USA), and the cDNA for mRNA, lncRNA, and circRNA were synthesized using the Rever Tra-Ace M-MLV kit (TOYOBO, Japan). The real-time quantitative PCR was conducted on the Bio-Rad CFX96 Touch Detection System (Bio-Rad, Hercules, CA, USA). Each reaction system consisted of 10 μL of 2× SYBR Green mix (TOYOBO, Tokyo, Japan), 7 μL of ddH2O, 0.5 μL of forward and reverse primers, and 2 μL of cDNA. The reaction program was as follows: 95 °C for 3 min, followed by 40 cycles of 95 °C for 15 s, 60 °C for 15 s, and 72 °C for 15 s. For miRNA quantitative PCR, U6 was used as the reference gene, while 40 s was used as the reference gene for mRNA, lncRNA, and circRNA. The relative expression level was calculated using the 2−ΔΔCt method. Student t-test was used to analyze statistical significance. p ≤ 0.05 indicated a statistically significant difference.
2.11. ceRNA Analysis
The competing endogenous RNA (ceRNA) hypothesis, a core mechanism for the action of ncRNAs, suggests that lncRNAs or circRNAs can interact with miRNAs and keep them from binding with mRNAs. The miRanda software (v3.3a) [42] was used to predict the interactions of mRNAs-miRNAs, lncRNAs-miRNAs, and circRNAs-miRNAs. lncRNA-mRNA or circRNA-mRNA pairs with the same miRNA binding site were calculated (p-value ≤ 0.05). The ceRNA networks were visualized using Cytoscape (v3.10.0).
3. Results
3.1. Overview of Gonadal smaRNA Sequencing of Common Carp
The ovaries and testes of common carp were sampled and sequenced in the present study. smaRNA sequencing revealed a total of 70,471,545 raw reads from common carp gonads. After quality filtering, 34,805,668 clean readings were obtained from the ovaries and 34,658,888 from the testes. Reads ranging in length from 18 to 35 nt, making up 95.57% of all reads, were picked out for further analysis (Table S1). In both ovaries and testes, the most prevalent size of smaRNA is 27–28 nt, followed by 26 and 29 nt (Figure 1A). With 92.48% of filtered reads mapping the reference genome, smaRNA sequences were heavily mapped to chromosomes NC_056574.1, NC_056575.1, and NC_056599.1. The distribution density of reads on the top 10 chromosomes is shown in Figure 1B. Housekeeping ncRNAs such as ribosomal RNA (rRNA), transfer RNA (tRNA), small nuclear RNA (snRNA), and small nucleolar RNA (snoRNA) were filtered out and counted, as shown in Table S2.
3.2. Analysis of miRNA in Common Carp Gonads
In the present study, 2,765,443 reads were mapped to 141 known mature miRNAs by comparison with the miRBase database, and 57 novel mature miRNAs were identified by predicting the hairpin structure of their precursor sequences using the remaining reads without annotated information. The details of the miRNA sequence are shown in Table S3. The first base of a mature miRNA sequence is highly biased due to the cleavage of the mature miRNA from the precursor miRNA by the enzyme Dicer. A considerable proportion of known and novel miRNAs showed a base bias at the first position toward U, suggesting that the miRNAs identified in this study are reliable (Figure 1C). The family analysis was performed on the identified miRNAs to explore conservation and presence in other species (Table S4). The target genes of miRNAs may be different due to one or two base changes in the seed region [43], and the result of miRNA base editing analysis is shown in Table S5.
3.3. Differential Expression of miRNAs between Ovaries and Testes
We normalized the expression levels with TPM to detect the DE miRNAs between ovaries and testes, as shown in Figure 2A. A high correlation between the different samples indicates that the data met the requirements of the subsequent analysis (Figure 2B). Among all the miRNAs identified, 11 are only present in the ovaries, while 11 are only in the testes (Figure 2C). We identified 63 miRNAs up-regulated and 61 miRNAs down-regulated in the ovaries compared to the testes (Figure 2D; Table S6). Of all the DE miRNAs, there are 91 known and 33 newly identified miRNAs, including miR-125b and miR-153b, which are known to be involved in sex differentiation. The 124 DE miRNAs were grouped into six subclusters with various expression patterns. Subclusters 1, 4, and 6 were down-regulated in testes, while subclusters 2, 3, and 5 were up-regulated in testes (Figure 2E; Table S7).
3.4. Target Gene Prediction and Functional Enrichment of miRNAs
The miRNAs identified in this study were predicted to target 6776 genes (Table S8). miR-205, miR-214, and miR-460-5p may be involved in inhibiting estrogen receptors, novel_158 and novel_313 may target foxl3 and hsd17β3, respectively. The results of KEGG analysis showed that 73 terms were enriched, and the top 20 terms were shown in Figure 2F. Among these, the TGF-beta signaling pathway is a classical pathway involved in sex differentiation. GO analysis showed that the target genes of the DE miRNAs were involved in G-protein coupled receptor activity and membrane components and were heavily involved in hormone-related pathways such as hormone activity and cellular response to hormone stimuli (Figure 2G).
3.5. Overview of Gonadal lncRNA Sequencing
The lncRNA sequencing yielded 38.18 G raw data from the ovaries and 40.39 G from the testes. After removing the low-quality reads, we retained 36.91 G clean reads for the ovaries and 39.31 G for the testes. A total of 87.80% of clean reads in the ovary libraries and 83.38% in testis libraries were then efficiently mapped to the common carp reference genome (Table S9). We identified 13,911 lncRNA transcripts, including 7684 known and 6227 newly discovered (Table S10). A total of 56.7% of the newly discovered lncRNAs were mapped in the intergenic regions (Figure 3A). mRNA and lncRNA transcript characteristics, including length, exon number, and open reading frame length, suggest that our predictions are highly accurate and reliable (Figure 3B). In addition, 252 novel mRNAs were identified in this study (Table S10).
3.6. Differential Expression of lncRNAs
The FPKM box plot (Figure S1A) shows the overall expression of transcripts in each sample. To ensure the reliability of the data, we conducted a correlation analysis (Figure S1B) and PCA analysis (Figure S1C). A total of 7863 lncRNA transcripts were differentially expressed. Of DE lncRNAs, 1764 transcripts were up-regulated, while 6099 were down-regulated in ovaries compared to testes (Figure 4A,B; Table S11). The cluster analysis demonstrates consistent expressions within groups and differential expressions in different groups (Figure 4C). To identify the target genes of lncRNAs, we performed co-expression and co-location analysis of lncRNAs and protein-coding genes (Table S12). Target genes smad1 and esr1 were identified by co-expression analysis and zp3 by co-location analysis. Differentially expressed target genes of lncRNAs were then subjected to GO and KEGG enrichment analysis to investigate the potential function. For target genes predicted from co-expression, GO analysis revealed enrichment of transcription factor complexes, DNA conformation changes, and pathways related to the cytoskeleton (Figure 4D), while the GnRH signaling pathway and oocyte meiosis were enriched in the KEGG enrichment (Figure 4E). The KEGG enrichment of target genes predicted from co-location analysis showed that the steroid biosynthesis was enriched (Figure S2).
3.7. Differentially Expressed mRNAs
We also identified 34,405 sex-biased expression mRNAs in this study, including 15,653 highly expressed in the ovaries and 18,752 highly expressed in the testes (Figure S3A). The well-known genes dmrt1, amh, and gsdf were up-regulated in the testes, and the female differentiation star genes foxl2 and cyp19a1a and the ZP family genes were expressed highly in the ovaries (Table S13). The common carp is an allo-tetraploid species, and we also analyzed the subgenomes to which the DE mRNA belonged. There were 15,930 and 17,236 DE mRNAs belonging to subgenomes A and B, respectively. The remaining DE mRNAs were transcribed from mitochondria or unplaced genomic scaffolds (Table S13). GO analysis revealed that differentially expressed genes are involved in DNA metabolism, the regulation of the cell cycle, and transcription factor activity (Figure S3B). KEGG enrichment showed that DE genes are enriched mainly in the oocyte meiosis, cell cycle, and ErbB signaling pathways (Figure S3C). The ErbB signaling pathway can activate PI3K/AKT and MAPK signaling pathways to modulate gonadal differentiation.
Alternative splicing (AS), which can generate multiple transcripts from the same precursor mRNA, regulates various physiological processes, including testicular and ovarian development [44]. We analyzed five kinds of alternative splicing events: skipping exon, intron retention, alternative 3′ splice site, alternative 5′ splice site, and mutually exclusive exons. We identified 24,976 AS events. Skipping exon is the major AS event, accounting for 69.23%. Moreover, 3893 differential AS events exist between ovaries and testes (Table S14). A case of the AS pattern of akap13 is shown in Figure S4.
3.8. Identification of Gonadal circRNAs
In total, we obtained 55.93 G circRNA sequencing raw data, and after data trimming by removing low-quality reads, ~344 million clean reads were retained. On average, 93.53% of clean reads in the ovaries and 85.40% in the testes were mapped to the reference genome (Table S15). Then, 6611 circRNAs were identified using CIRI and find_circ software, and 6603 RNAs were retained for subsequent analysis after removing the duplicates (Table S16). A total of 90.29% of circRNAs were transcribed from exon regions, while only 4.63% and 5.08% of circRNAs were transcribed from intron and intergenic regions, respectively (Figure 5A). The average length of circRNAs in the ovaries was approximately 284 nt, while the average length in the testes was approximately 298 nt. More than 60% of circRNAs were between 200 nt and 400 nt in length, and only 3.37% were longer than 500 nt (Figure 5B).
3.9. Expression Profiles and Functional Analysis of circRNAs in Gonads
The expressions of circRNAs were adjusted using TPM, and the TPM box plots showed the expression profiles of each sample (Figure 6A). According to the volcano map, 2082 circRNAs showed differential expression between ovaries and testes, of which 1047 circRNAs were highly expressed in ovaries while 1035 were highly expressed in testes (Figure 6B; Table S17). Clustering analysis revealed that the expression patterns of circRNAs were consistent within each group and distinguishable between ovaries and testes (Figure 6C). We analyzed the function of genes transcribing DE circRNAs. GO enrichment revealed that they were associated with calcium channel activity, cytoskeleton protein binding, and transition metal ion binding (Figure 6D). KEGG analysis revealed enriched pathways involved in oocyte meiosis, WNT signaling, MAPK signaling, and calcium signaling (Figure 6E). Furthermore, certain circRNAs possess coding capacity, and we conducted a circRNA coding potential prediction, as presented in Table S18.
3.10. Validation of Expression Profiles by qPCR
The relative expression levels of four miRNAs, lncRNAs, mRNAs, and circRNAs were examined by qPCR, respectively. Primers used for qPCR are listed in Table S19. As shown in Figure 7, miR-122, miR-192, XR_002016692, XR_002020130, cyp19a1a, zp3, novel_circ_0006872, and novel_circ_0010650 showed ovary-biased expression, while miR-153b, miR-93, XR_006158258, XR_006156779, dmrt1, sox9a, novel_circ_0003495, and novel_circ_0011083 showed testis-biased expression. The expression trends of the 16 selected RNAs by qPCR were consistent with those of the transcriptomic profile analysis, demonstrating the reliability of the sequencing data.
3.11. Complex Interactions among mRNAs and ncRNAs
The crosstalk between the competing endogenous RNAs (ceRNAs) is a vital mechanism by which ncRNAs function. lncRNAs or circRNAs can sponge off miRNAs and then regulate gene expression. We analyzed the binding sites of miRNAs on lncRNAs and circRNAs and then established the lncRNA/circRNA-miRNA-mRNA regulatory relationship (Table S20). For example, 17-beta hydroxysteroid dehydrogenase 3, which catalyzes the conversion of androstenedione to testosterone, is encoded by the hsd17β3 gene. We predicted that multiple circRNAs and lncRNAs might compete with hsd17β3 mRNA for binding to a new miRNA named novel_313. These circRNAs and lncRNAs might act as spongers to regulate the expression of hsd17β3 (Figure 8; Table S21).
4. Discussion
In recent years, a growing body of research has suggested that RNA molecules exhibit more functions than simply serving as templates for protein synthesis [5]. The proportion of ncRNAs in all transcripts is much higher than mRNAs, and ncRNAs can interact with DNA, proteins, and other RNA molecules to orchestrate most biological processes [45,46]. The function of ncRNAs in regulating sex differentiation in teleosts requires further investigation. In the present study, we systematically analyzed the expression profiles of ncRNAs between ovaries and testes, revealing numerous sex-specific ncRNAs in common carp.
Many studies have shown that lncRNAs can regulate well-known genes or signaling pathways to perform biological functions [47,48]. We identified 13,911 lncRNA transcripts, of which 6227 were novel in the gonads of common carp. Comparison of ovarian and testicular differences revealed differential expression of 7863 lncRNAs whose target genes are associated with vital biological processes, including steroid biosynthesis, oocyte meiosis, and GnRH signaling pathways. It suggested that lncRNAs may play a crucial role in gonadal development by regulating hormone synthesis. Although sex determination genes vary considerably among different fish species, the downstream genes and pathways of sex differentiation tend to be more conserved [49,50]. We analyzed the differential expression of mRNA. The results showed that the well-studied genes foxl2 and cyp19a1a were highly expressed in the ovaries, while dmrt1 and amh were highly expressed in the testes, indicating the reliability of our results. Previous studies show that genome duplication can increase the complexity of gene regulatory systems and provide an opportunity for the gain or loss of gene function [51]. Our results showed 15,930 mRNAs transcribed from the A subgenome and 17,236 transcribed from the B subgenome, indicating that both A and B subgenomes play vital roles in regulating sex differentiation.
Most circRNAs contain multiple miRNA binding sites, allowing them to bind and sequester miRNAs, thereby preventing miRNAs from interacting with their target mRNAs [35]. In addition, some circRNAs can undergo translation and produce functional proteins [52]. In mice, a testis-specific circRNA, circSry, regulates the functional implementation of the sex-determination gene Sry by competitively binding to miRNA-138 [15]. We identified 6603 circRNAs in the gonads of common carp, of which 2064 showed differential expression between ovaries and testes. Enrichment analysis of DE circRNAs indicates that WNT signaling and MAPK signaling pathways were enriched. These pathways are known to play pivotal roles in preserving ovarian differentiation and regulating the proliferation of primordial germ cells, respectively [53,54]. circRNAs may potentially contribute to sex differentiation in common carp through the involvement of the signaling mentioned above pathways. In addition, we analyzed the coding potential of circRNAs in this study.
The miRNAs are some of the most widely studied smaRNAs, and they are conserved across different species [55,56]. The analysis of the length distribution of smaRNAs revealed that the size range 26–29 nt is the most prevalent size. This finding may be attributed to the higher expression of piRNAs, a class of smaRNAs with lengths of 25–30 nt, which are mainly expressed in the germline and play a crucial role in germline development [57]. Our research is consistent with the findings of previous studies conducted on the zebrafish (Danio rerio) and tilapia (Oreochromis niloticus) [58,59].
Mature miRNAs can modulate the expression of target genes by recognizing and binding to homologous mRNA sequences. Previous studies have demonstrated that miR-17-5p and miR-20a may inhibit the expression of dmrt1, while miR-138, miR-200a, and miR-338 may negatively regulate the expression of cyp17a2 in tilapias [60]. In the orange-spotted grouper (Epinephelus coioides), miR-26a-5p may act as a regulatory factor of cyp19a1a to modulate the production of estradiol-17β (E2), the primary functional estrogen [61]. miR-153b-3p, predominantly expressed in the testes of common carp, functions to fine-tune the complex process of male germ cell proliferation and differentiation by directly targeting amh [62]. In the current study, we identified 198 mature miRNAs, of which 57 are novel. There are 124 differentially expressed miRNAs, including miR-112, miR-125b, and miR-153b. The prediction of the miRNA target gene suggests that miR-205, miR-214, and miR-460-5p may participate in inhibiting several estrogen receptors. foxl3, which is related to influencing the sex fate decision in the medaka (Oryzias latipes) and tilapia [63,64], may be regulated by novel_158 in the common carp. Additionally, the differentially expressed miRNAs have a function in the TGF-beta signaling pathway, a traditional pathway that regulates sex differentiation and hormone-related pathways. These findings imply that miRNA plays a critical role in sex differentiation and maintenance in common carp.
The ceRNA crosstalk is recognized as a critical mechanism for ncRNA action. Many lncRNAs or circRNAs are known as miRNA sponges since they contain numerous high-affinity miRNA binding sites. They can bind and sequester particular miRNAs, preventing them from interacting with target mRNAs [65]. For example, a circRNA called circdmrt1 and a lncRNA called AMSDT were up-regulated in the Chinese tongue sole (Cynoglossus semilaevis) testes. These two ncRNAs can bind to cse-miR-196 and hence abolish its endogenous suppressive effect on the target gene gsdf [66]. In the current study, we analyzed ncRNA-mRNA interactions in common carp gonads to illuminate potential ceRNA crosstalk in sex determination and differentiation. We predicted that several circRNAs and lncRNAs might compete with hsd17β3 mRNA to regulate its expression.
5. Conclusions
We yielded a substantial pool of ncRNAs associated with sex differentiation in the present study. This valuable collection of scientific data will serve as a valuable resource not only for understanding the underlying genetic basis of gonad differentiation and development but also for sex-control breeding in common carp.
Conceptualization, J.L. and Y.Z.; methodology, S.Y., Z.M. and M.H.; software, J.Z. and Q.W.; validation, X.S., Y.Z. and M.H.; formal analysis, Y.C. (Yingjie Chen), K.W. and M.H.; investigation, R.Z. and Y.C. (Yiming Cao); resources, S.Y., Y.C. (Yingjie Chen), Z.M. and K.W.; data curation, Q.W.; writing—original draft preparation, M.H.; writing—review and editing, J.L. and M.H.; visualization, J.Z. and Q.W.; supervision, J.L.; project administration, J.L. and Y.Z.; funding acquisition, J.L. and M.H. All authors have read and agreed to the published version of the manuscript.
The animal study protocol was approved by the Animal Care and Use Committee of the Chinese Academy of Fishery Sciences (protocol code ACUC-CAFS-20220615, approved on 15 June 2022).
Not applicable.
The sequence data presented in the current study have been submitted to the National Center for Biotechnology Information (NCBI) BioProject database (Accession Numbers: PRJNA1020596).
We are grateful to Zhengsheng Zhao, Honglu chen, and Tianyang Zhang for their help with experimental materials and fish culture.
The authors declare no conflict of interest.
mRNA, messenger RNA; ncRNAs, non-coding RNAs; miRNAs, microRNAs; smaRNAs, small RNAs; lncRNAs, long non-coding RNAs; circRNAs, circular RNAs; ceRNAs, competing endogenous RNAs; siRNA, small interfering RNAs; piRNAs, piwi interacting RNAs; RIN, RNA integrity number; AS, alternative splicing; IRES, Internal Ribosome Entry Sites; FPKM, Fragments Per Kilobase of transcript sequence per Millions of base pairs sequenced; TPM, Transcripts Per Millions; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; CNCI, Coding-Non-Coding Index; CPC, Coding Potential Calculator; PFAM, Pfam ScanC; BP, biological process; CC: cell component; MF, molecular function; DE, differentially expressed; nt, nucleotides; hsd17β3, 17-beta hydroxysteroid dehydrogenase 3.
Footnotes
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
Enlarge this image.Figure 1. Characterization of gonadal smaRNAs and miRNAs in common carp. (A) Length distribution of smaRNA sequences in the ovaries and testes. (B) Distribution of gonadal miRNAs on chromosomes. The outermost circle is designated as the chromosome to be displayed. Orange and green represent reads mapped to positive and negative chains, respectively. (C) First base preference of known miRNAs ranging in length from 18 to 30 nt in the ovaries and testes. The horizontal axis represents the length of the miRNA, while the vertical axis represents the ratio of A/U/C/G in the first base.
Enlarge this image.Figure 2. Identification of differentially expressed miRNAs and functional analysis of their target genes in ovaries and testes. (A) TPM density distribution of different samples. The horizontal axis is the log10 (TPM + 1) value, and the vertical axis is the density distribution. (B) The correlation analysis of different samples. (C) Venn diagram of miRNAs expressed in ovaries and testes. (D) Volcano map of differentially expressed miRNAs. Blue dots represent the up-regulated miRNAs, and red dots represent the down-regulated miRNAs in ovaries compared to testes. (E) Clustering of the 124 DE miRNAs into six subclusters. Gray lines represent the relative expression levels of each gene; blue lines represent the mean value in each subcluster; red lines represent the reference. (F) The KEGG pathways with the lowest p-value by DE miRNA target genes. (G) The GO terms with the lowest p-value by DE miRNA target genes. BP: biological process; CC: cell component; MF: molecular function.
Enlarge this image.Figure 3. Characteristics of novel lncRNA transcripts identified in this study. (A) Classification and proportion of novel lncRNAs. antisense: overlapped with the antisense strand of one or more exons of a protein-coding gene; lincRNA: located in the intergenic region; sense overlapping: overlapped with sense strand of one or more exons of a protein-coding gene; sense intronic: located in the intronic region of a coding gene and does not intersect with any of exons on the same strand. (B) Comparison of the characteristics of all the lncRNAs and mRNAs.
Enlarge this image.Figure 4. Analysis of DE lncRNAs between ovaries and testes. (A) Volcano map of DE lncRNAs. The blue dots represent the up-regulated lncRNAs, and red dots represent the down-regulated lncRNAs in ovaries compared to testes. (B) The number of DE lncRNAs. (C) Cluster analysis of DE lncRNAs. (D) The GO terms with the lowest p-value by DE lncRNA target genes from co-expression analysis. BP: biological process; CC: cell component; MF: molecular function. (E) The KEGG pathways with the lowest p-value by DE lncRNA target genes from co-expression analysis.
Enlarge this image.Figure 5. Identification of novel circRNAs. (A) Proportion of regions from which circRNAs were transcribed. (B) Length distribution of circRNAs.
Enlarge this image.Figure 6. Expression and functional analysis of circRNAs in the gonads of common carp. (A) TPM box plot showing the expression profile of each sample. The horizontal axis is the value of log10 (TPM + 1). (B) The vocal map of DE circRNAs between ovaries and testes. Blue dots represent the up-regulated circRNAs, and red dots represent the down-regulated circRNAs in ovaries compared to testes. (C) Hierarchical clustering heatmap analysis of DE circRNAs between ovaries and testes. (D) The GO terms with the lowest p-value by genes transcribing DE circRNAs. (E) The KEGG pathways with the lowest p-value by genes transcribing DE circRNAs.
Enlarge this image.Figure 7. qPCR validation of miRNA, lncRNA, and circRNA sequencing data. The data shown are the mean ± SEM of three replicates, p ≤ 0.05.
Enlarge this image.Figure 8. A case of ceRNA crosstalk by Cytoscape. Blue and red nodes represented lncRNAs and circRNAs, respectively, and the details of lncRNAs and circRNA are shown in Table S21; the edge represents the interaction between the two nodes.
Supplementary Materials
The following supporting information can be downloaded at
References
1. Nagahama, Y.; Chakraborty, T.; Paul-Prasanth, B.; Ohta, K.; Nakamura, M. Sex determination, gonadal sex differentiation, and plasticity in vertebrate species. Physiol. Rev.; 2021; 101, pp. 1237-1308. [DOI: https://dx.doi.org/10.1152/physrev.00044.2019] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/33180655]
2. Mei, J.; Gui, J.F. Genetic basis and biotechnological manipulation of sexual dimorphism and sex determination in fish. Sci. China Life Sci.; 2015; 58, pp. 124-136. [DOI: https://dx.doi.org/10.1007/s11427-014-4797-9] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25563981]
3. Horne, C.R.; Hirst, A.G.; Atkinson, D. Selection for increased male size predicts variation in sexual size dimorphism among fish species. Proc. Biol. Sci.; 2020; 287, 20192640. [DOI: https://dx.doi.org/10.1098/rspb.2019.2640]
4. Li, X.Y.; Mei, J.; Ge, C.T.; Liu, X.L.; Gui, J.F. Sex determination mechanisms and sex control approaches in aquaculture animals. Sci. China Life Sci.; 2022; 65, pp. 1091-1122. [DOI: https://dx.doi.org/10.1007/s11427-021-2075-x] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35583710]
5. Beermann, J.; Piccoli, M.T.; Viereck, J.; Thum, T. Non-coding RNAs in development and disease: Background, mechanisms, and therapeutic approaches. Physiol. Rev.; 2016; 96, pp. 1297-1325. [DOI: https://dx.doi.org/10.1152/physrev.00041.2015] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/27535639]
6. Bartel, D.P. Metazoan microRNAs. Cell; 2018; 173, pp. 20-51. [DOI: https://dx.doi.org/10.1016/j.cell.2018.03.006]
7. Bracken, C.P.; Scott, H.-S.; Goodall, G.J. A network-biology perspective of microRNA function and dysfunction in cancer. Nat. Rev. Genet.; 2016; 17, pp. 719-732. [DOI: https://dx.doi.org/10.1038/nrg.2016.134]
8. Real, F.M.; Sekido, R.; Lupianez, D.G.; Lovell-Badge, R.; Jiménez, R.; Burgos, M. A microRNA (mmu-mir-124) prevents sox9 expression in developing mouse ovarian cells. Biol. Reprod.; 2013; 89, 78. [DOI: https://dx.doi.org/10.1095/biolreprod.113.110957]
9. Miao, N.; Wang, X.; Hou, Y.; Feng, Y.; Gong, Y. Identification of male-biased microRNA-107 as a direct regulator for nuclear receptor subfamily 5 group A member 1 based on sexually dimorphic microRNA expression profiling from chicken embryonic gonads. Mol. Cell Endocrinol.; 2016; 429, pp. 29-40. [DOI: https://dx.doi.org/10.1016/j.mce.2016.03.033]
10. Rinn, J.L.; Chang, H.Y. Genome regulation by long noncoding RNAs. Annu. Rev. Biochem.; 2012; 81, pp. 145-166. [DOI: https://dx.doi.org/10.1146/annurev-biochem-051410-092902]
11. Ryohei, S.; Robin, L.B. Sex determination involves synergistic action of SRY and SF1 on a specific Sox9 enhancer. Nature; 2008; 453, pp. 930-934. [DOI: https://dx.doi.org/10.1038/nature06944]
12. Zhang, L.; Lu, H.; Xin, D.; Cheng, H.; Zhou, R. A novel ncRNA gene from mouse chromosome 5 trans-splices with Dmrt1 on chromosome 19. Biochem. Biophys. Res. Commun.; 2010; 400, pp. 696-700. [DOI: https://dx.doi.org/10.1016/j.bbrc.2010.08.130]
13. Chen, L.L. The biogenesis and emerging roles of circular RNAs. Nat. Rev. Mol. Cell Biol.; 2016; 17, pp. 205-211. [DOI: https://dx.doi.org/10.1038/nrm.2015.32]
14. Shen, Y.; Guo, X.; Wang, W. Identification and characterization of circular RNAs in zebrafish. FEBS Lett.; 2017; 591, pp. 213-220. [DOI: https://dx.doi.org/10.1002/1873-3468.12500]
15. Hansen, T.B.; Jensen, T.I.; Clausen, B.H.; Bramsen, J.B.; Finsen, B.; Damgaard, C.K.; Kjems, J. Natural RNA circles function as efficient microRNA sponges. Nature; 2013; 495, pp. 384-388. [DOI: https://dx.doi.org/10.1038/nature11993]
16. Jiang, M.; Wu, X.; Chen, K.; Luo, H.; Yu, W.; Jia, S.; Li, Y.; Wang, Y.; Yang, P.; Zhu, Z. et al. Production of YY supermale and XY physiological female common carp for potential eradication of this invasive species. J. World Aquacult. Soc.; 2018; 49, pp. 315-327. [DOI: https://dx.doi.org/10.1111/jwas.12492]
17. Tessema, A.; Getahun, A.; Mengistou, S.; Fetahi, T.; Dejen, E. Reproductive biology of common carp (Cyprinus carpio Linnaeus, 1758) in Lake Hayq, Ethiopia. Fish. Aquat. Sci.; 2020; 23, 16. [DOI: https://dx.doi.org/10.1186/s41240-020-00162-x]
18. Wang, M.; Chen, L.; Zhou, Z.; Xiao, J.; Chen, B.; Huang, P.; Li, C.; Xue, Y.; Liu, R.; Bai, Y. et al. Comparative transcriptome analysis of early sexual differentiation in the male and female gonads of common carp (Cyprinus carpio). Aquaculture; 2023; 563, 738984. [DOI: https://dx.doi.org/10.1016/j.aquaculture.2022.738984]
19. Zhai, G.; Shu, T.; Chen, K.; Lou, Q.; Jia, J.; Huang, J.; Shi, C.; Jin, X.; He, J.; Jiang, D. et al. Successful production of an all-female common carp (Cyprinus carpio L.) population using cyp17a1-deficient neomale carp. Engineering; 2022; 8, pp. 181-189. [DOI: https://dx.doi.org/10.1016/j.eng.2021.03.026]
20. Xu, P.; Zhang, X.; Wang, X.; Li, J.; Liu, G.; Kuang, Y.; Xu, J.; Zheng, X.; Ren, L.; Wang, G. et al. Genome sequence and genetic diversity of the common carp, Cyprinus carpio. Nat. Genet.; 2014; 46, pp. 1212-1219. [DOI: https://dx.doi.org/10.1038/ng.3098] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25240282]
21. Li, J.T.; Wang, Q.; Huang Yang, M.D.; Li, Q.-S.; Cui, M.S.; Dong, Z.J.; Wang, H.W.; Yu, J.H.; Zhao, Y.J.; Yang, C.R. et al. Parallel subgenome structure and divergent expression evolution of allo-tetraploid common carp and goldfish. Nat. Genet.; 2021; 53, pp. 1493-1503. [DOI: https://dx.doi.org/10.1038/s41588-021-00933-9] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34594040]
22. Jia, Y.; Nan, P.; Zhang, W.; Wang, F.; Zhang, R.; Liang, T.; Ji, X.; Du, Q.; Chang, Z. Transcriptome analysis of three critical periods of ovarian development in Yellow River carp (Cyprinus carpio). Theriogenology; 2018; 105, pp. 15-26. [DOI: https://dx.doi.org/10.1016/j.theriogenology.2017.08.027] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/28923703]
23. Zhao, Y.; Chen, K.; Liu, F.; Jiang, M.; Chen, Z.; Chen, H.; Song, Y.; Tao, B.; Cui, X.; Li, Y. et al. Dynamic gene expression and alternative splicing events demonstrate co-regulation of testicular differentiation and maturation by the brain and gonad in common carp. Front. Endocrinol.; 2021; 12, 820463. [DOI: https://dx.doi.org/10.3389/fendo.2021.820463]
24. Langmead, B.; Trapnell, C.; Pop, M.; Salzberg, S.L. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol.; 2009; 10, R25. [DOI: https://dx.doi.org/10.1186/gb-2009-10-3-r25] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/19261174]
25. Sam, G.J.; Russell, J.G.; Stijn, V.D.; Alex, B.; Anton, J.E. miRBase: Tools for microRNA genomics. Nucleic Acids Res.; 2007; 36, pp. D154-D158. [DOI: https://dx.doi.org/10.1093/nar/gkm952]
26. Wen, M.; Shen, Y.; Shi, S.; Tang, T. miREvo: An integrative microRNA evolutionary analysis platform for next-generation sequencing experiments. BMC Bioinform.; 2012; 13, 140. [DOI: https://dx.doi.org/10.1186/1471-2105-13-140] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/22720726]
27. Friedlander, M.R.; Mackowiak, S.D.; Li, N.; Chen, W.; Rajewsky, N. miRDeep2 accurately identifies known and hundreds of novel microRNA genes in seven animal clades. Nucleic Acids Res.; 2011; 40, pp. 37-52. [DOI: https://dx.doi.org/10.1093/nar/gkr688]
28. Zheng, Y.; Ji, B.; Song, R.; Wang, S.; Li, T.; Zhang, X.; Chen, K.; Li, T.; Li, J. Accurate detection for a wide range of mutation and editing sites of microRNAs from small RNA high-throughput sequencing profiles. Nucleic Acids Res.; 2016; 44, e123. [DOI: https://dx.doi.org/10.1093/nar/gkw471]
29. Pertea, M.; Pertea, G.M.; Antonescu, C.M.; Chang, T.C.; Mendell, J.T.; Salzberg, S.L. StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat. Biotechnol.; 2015; 33, pp. 290-295. [DOI: https://dx.doi.org/10.1038/nbt.3122]
30. Sun, L.; Luo, H.; Bu, D.; Zhao, G.; Yu, K.; Zhang, C.; Liu, Y.; Chen, R.; Zhao, Y. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Res.; 2013; 41, e166. [DOI: https://dx.doi.org/10.1093/nar/gkt646]
31. Kang, Y.J.; Yang, D.C.; Kong, L.; Hou, M.; Meng, Y.Q.; Wei, L.; Gao, G. CPC2: A fast and accurate coding potential calculator based on sequence intrinsic features. Nucleic Acids Res.; 2017; 45, pp. W12-W16. [DOI: https://dx.doi.org/10.1093/nar/gkx428] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/28521017]
32. Mistry, J.; Bateman, A.; Finn, R.D. Predicting active site residue annotations in the Pfam database. BMC Bioinform.; 2007; 8, 298. [DOI: https://dx.doi.org/10.1186/1471-2105-8-298] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/17688688]
33. Kopp, F.; Mendell, J.T. Functional classification and experimental dissection of long noncoding RNAs. Cell; 2018; 172, pp. 393-407. [DOI: https://dx.doi.org/10.1016/j.cell.2018.01.011] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/29373828]
34. Bao, Z.; Yang, Z.; Huang, Z.; Zhou, Y.; Cui, Q.; Dong, D. LncRNADisease 2.0: An updated database of long non-coding RNA-associated diseases. Nucleic Acids Res.; 2019; 47, pp. D1034-D1037. [DOI: https://dx.doi.org/10.1093/nar/gky905] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30285109]
35. Memczak, S.; Jens, M.; Elefsinioti, A.; Torti, F.; Krueger, J.; Rybak, A.; Maier, L.; Mackowiak, S.D.; Gregersen, L.H.; Munschauer, M. et al. Circular RNAs are a large class of animal RNAs with regulatory potency. Nature; 2013; 495, pp. 333-338. [DOI: https://dx.doi.org/10.1038/nature11928] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/23446348]
36. Gao, Y.; Wang, J.; Zhao, F. CIRI: An efficient and unbiased algorithm for de novo circular RNA identification. Genome Biol.; 2015; 16, 4. [DOI: https://dx.doi.org/10.1186/s13059-014-0571-3] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25583365]
37. Zhao, J.; Wu, J.; Xu, T.; Yang, Q.; He, J.; Song, X. IRESfinder: Identifying RNA internal ribosome entry site in eukaryotic cell using framed k-mer features. J. Genet. Genom.; 2018; 45, pp. 403-406. [DOI: https://dx.doi.org/10.1016/j.jgg.2018.07.006]
38. Robinson, M.D.; McCarthy, D.J.; Smyth, G.K. edgeR: A Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics; 2010; 26, pp. 139-140. [DOI: https://dx.doi.org/10.1093/bioinformatics/btp616]
39. Love, M.I.; Huber, W.; Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol.; 2014; 15, 550. [DOI: https://dx.doi.org/10.1186/s13059-014-0550-8]
40. Yu, G.; Wang, L.-G.; Han, Y.; He, Q.-Y. ClusterProfiler: An R package for comparing biological themes among gene clusters. Omics; 2012; 16, pp. 284-287. [DOI: https://dx.doi.org/10.1089/omi.2011.0118]
41. Chen, C.; Ridzon, D.A.; Broomer, A.J.; Zhou, Z.; Lee, D.H.; Nguyen, J.T.; Barbisin, M.; Xu, N.L.; Mahuvakar, V.R.; Andersen, M.R. et al. Real-time quantification of microRNAs by stem-loop RT-PCR. Nucleic Acids Res.; 2005; 33, e179. [DOI: https://dx.doi.org/10.1093/nar/gni178] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/16314309]
42. Riffo-Campos, A.; Riquelme, I.; Brebi-Mieville, P. Tools for sequence-based miRNA target prediction: What to choose?. Int. J. Mol. Sci.; 2016; 17, 1987. [DOI: https://dx.doi.org/10.3390/ijms17121987] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/27941681]
43. Wei, Y.; Chen, S.; Yang, P.; Ma, Z.; Kang, L. Characterization and comparative profiling of the small RNA transcriptomes in two phases of locust. Genome Biol.; 2009; 10, R6. [DOI: https://dx.doi.org/10.1186/gb-2009-10-1-r6] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/19146710]
44. Nilsen, T.W.; Graveley, B.R. Expansion of the eukaryotic proteome by alternative splicing. Nature; 2010; 463, pp. 457-463. [DOI: https://dx.doi.org/10.1038/nature08909]
45. Miyawaki, S.; Tachibana, M. Role of epigenetic regulation in mammalian sex determination. Curr. Top Dev. Biol.; 2019; 134, pp. 195-221. [DOI: https://dx.doi.org/10.1016/bs.ctdb.2019.01.008]
46. Zhang, C.; Han, B.; Xu, T.; Li, D. The biological function and potential mechanism of long non-coding RNAs in cardiovascular disease. J. Cell Mol. Med.; 2020; 24, pp. 12900-12909. [DOI: https://dx.doi.org/10.1111/jcmm.15968]
47. Brown, C.; Hendrich, B.; Rupert, J.; Lafreniere, R.; Xing, Y.; Lawrence, J.; Willard, H. The human XIST gene: Analysis of a 17 kb inactive X-specific RNA that contains conserved repeats and is highly localized within the nucleus. Cell; 1992; 71, pp. 527-542. [DOI: https://dx.doi.org/10.1016/0092-8674(92)90520-M]
48. Yousefi, H.; Maheronnaghsh, M.; Molaei, F.; Mashouri, L.; Aref, A.R.; Momeny, M.; Alahari, S.K. Long noncoding RNAs and exosomal lncRNAs: Classification, and mechanisms in breast cancer metastasis and drug resistance. Oncogene; 2019; 39, pp. 953-974. [DOI: https://dx.doi.org/10.1038/s41388-019-1040-y]
49. Herpin, A.; Schartl, M. Plasticity of gene-regulatory networks controlling sex determination: Of masters, slaves, usual suspects, newcomers, and usurpators. EMBO Rep.; 2015; 16, pp. 1260-1274. [DOI: https://dx.doi.org/10.15252/embr.201540667]
50. Capel, B. Vertebrate sex determination: Evolutionary plasticity of a fundamental switch. Nat. Rev. Genet.; 2017; 18, pp. 675-689. [DOI: https://dx.doi.org/10.1038/nrg.2017.60]
51. Salzburger, W. Understanding explosive diversification through cichlid fish genomics. Nat. Rev. Genet.; 2018; 19, pp. 705-717. [DOI: https://dx.doi.org/10.1038/s41576-018-0043-9] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30111830]
52. Pamudurti, N.R.; Bartok, O.; Jens, M.; Ashwal-Fluss, R.; Stottmeister, C.; Ruhe, L.; Hanan, M.; Wyler, E.; Perez-Hernandez, D.; Ramberger, E. et al. Translation of circRNAs. Mol. Cell; 2017; 66, pp. 9-21. [DOI: https://dx.doi.org/10.1016/j.molcel.2017.02.021] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/28344080]
53. Lin, Y.T.; Capel, B. Cell fate commitment during mammalian sex determination. Curr. Opin. Genet Dev.; 2015; 32, pp. 144-152. [DOI: https://dx.doi.org/10.1016/j.gde.2015.03.003]
54. Feng, K.; Cui, X.; Song, Y.; Tao, B.; Chen, J.; Wang, J.; Liu, S.; Sun, Y.; Zhu, Z.; Trudeau, V. et al. Gnrh3 regulates PGC proliferation and sex differentiation in developing zebrafish. Endocrinology; 2020; 161, bqz024. [DOI: https://dx.doi.org/10.1210/endocr/bqz024]
55. Lee, R.C.; Feinbaum, R.L.; Ambros, V. The C. elegans heterochronic gene lin-4 encodes small RNAs with antisense complementarity to lin-14. Cell; 1993; 75, pp. 843-854. [DOI: https://dx.doi.org/10.1016/0092-8674(93)90529-Y] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/8252621]
56. Pasquinelli, A.; Reinhart, B.; Slack, F.; Martindale, M.; Kuroda, M.; Maller, B.; Hayward, D.; Ball, E.; Degnan, B.; Muller, P. et al. Conservation of the sequence and temporal expression of let-7 heterochronic regulatory RNA. Nature; 2000; 408, pp. 86-89. [DOI: https://dx.doi.org/10.1038/35040556] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/11081512]
57. Loubalova, Z.; Fulka, H.; Horvat, F.; Pasulka, J.; Malik, R.; Hirose, M.; Ogura, A.; Svoboda, P. Formation of spermatogonia and fertile oocytes in golden hamsters requires piRNAs. Nat. Cell Biol.; 2021; 23, pp. 992-1001. [DOI: https://dx.doi.org/10.1038/s41556-021-00746-2]
58. Xiao, J.; Zhong, H.; Zhou, Y.; Yu, F.; Gao, Y.; Luo, Y.; Tang, Z.; Guo, Z.; Guo, E.; 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, e86821. [DOI: https://dx.doi.org/10.1371/journal.pone.0086821]
59. Houwing, S.; Kamminga, L.M.; Berezikov, E.; Cronembold, D.; Girard, A.; van den Elst, H.; Filippov, D.V.; Blaser, H.; Raz, E.; Moens, C.B. et al. A role for Piwi and piRNAs in germ cell maintenance and transposon silencing in zebrafish. Cell; 2007; 129, pp. 69-82. [DOI: https://dx.doi.org/10.1016/j.cell.2007.03.026]
60. 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, pp. 970-978. [DOI: https://dx.doi.org/10.1016/j.theriogenology.2015.11.006]
61. Yu, Q.; Peng, C.; Ye, Z.; Tang, Z.; Li, S.; Xiao, L.; Liu, S.; Yang, Y.; Zhao, M.; Zhang, Y. et al. An estradiol-17β/miRNA-26a/cyp19a1a regulatory feedback loop in the protogynous hermaphroditic fish, Epinephelus coioides. Mol. Cell. Endocrinol.; 2020; 504, 110689. [DOI: https://dx.doi.org/10.1016/j.mce.2019.110689] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/31891771]
62. Ma, W.; Chen, K.; Jiang, M.; Jia, S.; Chen, J.; Tao, B.; Song, Y.; Li, Y.; Wang, Y.; Xiao, W. et al. MiR-153b-3p regulates the proliferation and differentiation of male germ cells by targeting amh in common carp (Cyprinus carpio). Aquaculture; 2021; 535, 736420. [DOI: https://dx.doi.org/10.1016/j.aquaculture.2021.736420]
63. Nishimura, T.; Sato, T.; Yamamoto, Y.; Watakabe, I.; Ohkawa, Y.; Suyama, M.; Kobayashi, S.; Tanaka, M. foxl3 is a germ cell–intrinsic factor involved in sperm-egg fate decision in medaka. Science; 2015; 349, pp. 328-331. [DOI: https://dx.doi.org/10.1126/science.aaa2657] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/26067255]
64. Dai, S.; Qi, S.; Wei, X.; Liu, X.; Li, Y.; Zhou, X.; Xiao, H.; Lu, B.; Wang, D.; Li, M. Germline sexual fate is determined by the antagonistic action of dmrt1 and foxl3/foxl2 in tilapia. Development; 2021; 148, dev199380. [DOI: https://dx.doi.org/10.1242/dev.199380] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/33741713]
65. Salmena, L.; Poliseno, L.; Tay, Y.; Kats, L.; Pandolfi, P.P. A ceRNA hypothesis: The Rosetta Stone of a hidden RNA language?. Cell; 2011; 146, pp. 353-358. [DOI: https://dx.doi.org/10.1016/j.cell.2011.07.014]
66. Tang, L.; Huang, F.; You, W.; Poetsch, A.; Nobrega, R.H.; Power, D.M.; Zhu, T.; Liu, K.; Wang, H.Y.; Wang, Q. et al. ceRNA crosstalk mediated by ncRNAs is a novel regulatory mechanism in fish sex determination and differentiation. Genome Res.; 2022; 32, pp. 1502-1515. [DOI: https://dx.doi.org/10.1101/gr.275962.121]
© 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/). Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.
Details
; Cao, Yiming 1 ; Yu, Shuangting 2 ; Wang, Kaikuo 3 ; Chen, Yingjie 3 ; Ma, Ziyao 3 ; Sun, Xiaoqing 1 ; Zhang, Yan 1 ; Li, Jiongtang 1 1 Key Laboratory of Aquatic Genomics, Ministry of Agriculture and Rural Affairs and Beijing Key Laboratory of Fishery Biotechnology, Chinese Academy of Fishery Sciences, Beijing 100141, China;
2 Key Laboratory of Aquatic Genomics, Ministry of Agriculture and Rural Affairs and Beijing Key Laboratory of Fishery Biotechnology, Chinese Academy of Fishery Sciences, Beijing 100141, China;
3 National Demonstration Center for Experimental Fisheries Science Education, Shanghai Ocean University, Shanghai 201306, China;