Differentially Expressed MicroRNAs in Testes of Dominant and Subordinate Nile Tilapia Males and Identification of Oni-miR-499 as Regulator of amh Gene Expression Open Access

Non-coding RNAs (ncRNAs) affect the processes of transcription and translation and thus are central regulators of gene expression [1]. MicroRNAs (miRNAs) are a class of endogenous ncRNAs that control gene expression post-transcriptionally in their mature form. Mostly starting from long microRNA genes, so called pri-miRNAs transcribed in the nucleus. They are further processed by nuclear enzyme Drosha into pre-miRNAs. After export by Exportin-5 into the cytoplasm, the Dicer complex produces a very short mature miRNA. With a final length of about 22 nucleotides, they belong to the small ncRNAs but play a major role as epigenetic modulators. When they interact with complementary mRNA sequences, these mRNA transcripts are destabilized or their translation is repressed [2, 3]. In animals, miRNAs usually bind within the area of the 3′ untranslated region (3′-UTR) of a target mRNA. They play a central role in cellular mechanisms like differentiation, development, proliferation, apoptosis, and tumorigenesis [4, 5].

miRNAs have already been associated with mammalian reproduction [6]. In male mammals, there is some evidence that they influence gonadal and germ cell development as they regulate the self-renewal, maintenance and differentiation of spermatogonial stem cells and the process of spermatogenesis [7, 8]. However, only limited information is available about the role of miRNAs in reproductive processes in teleost fish outside of the mammalian lineage [9, 10].

In fish, as in all vertebrates, reproduction is regulated by the endocrine system. The BPG axis is of central importance because it regulates pubertal development as well as the maintenance of adult reproductive ability [11]. The gonadotropin-releasing hormone neurons located in the hypothalamic region of the brain combine internal (nutritional and hormonal) and external (environmental and social) information to control reproduction. They can send projections directly to the pituitary gland, where gonadotropin-releasing hormone stimulates the synthesis of the follicle-stimulating hormone (Fsh) and luteinizing hormone (Lh) [12‒14]. In males, these gonadotropins regulate steroidogenesis and via the sex hormones also spermatogenesis. Gonadal feedback on the gonadotropin secretion imparted by sex steroids in fish has already been shown some decades ago [15]. Since then, various factors influencing the reproductive axis in teleost fish have been discovered [16‒20].

In this regard, the anti-Müllerian hormone (Amh) became a focus of research. It was described in tetrapods, in which it induces the regression of the Müllerian ducts in male embryos and drives the progression of the male reproductive tract [21]. In contrast, teleost fish do not develop Müllerian ducts but exhibit amh orthologous genes [22]. In teleosts, Amh plays important roles in sex determination, gametogenesis, and gonad development. Recent works report the Amh or its type II receptor (Amhr2) genes as male master sex determining genes in certain teleost species, mostly via linage-specific amh duplications, as is the case with some Nile tilapia strains. There, the task of sex determination is taken over by a tandem duplication of amh (amh∆Y and amhY). Such male-specific amh copies are expressed early before the onset of sexual differentiation of the gonads, possibly by regulating early germ cell development (survival, differentiation, and proliferation) similarly to one of the known functions of autosomal amh copies [23]. Besides inhibiting spermatogenesis and gonia proliferation in fish [24‒26], experiments in zebrafish provided evidence that it also negatively regulates steroidogenesis and inhibits Fsh-induced 11-KT release from Leydig cells [26, 27]. Furthermore, in vitro experiments with testis explants from Japanese eel and zebrafish have shown that 11-KT-stimulated proliferation of type A spermatogonia is suppressed by recombinant Amh [24, 26]. The gonadotropin Fsh is a regulator of gametogenesis and steroidogenesis [26, 28‒34] and might in turn inhibit the expression of amh in fish [26, 27, 33, 34]. This could be necessary for the realization of the gonadotropic effects.

Oreochromis niloticus is an African cichlid fish that has not only become one of the most important species in tropical and subtropical aquaculture [35] but also an important model for reproduction and gonadal development studies [36‒40]. In recent years, there have been several reports on the differences in miRNAs expression between males and females [41‒44]. In contrast, we analyzed miRNA expression in the same sex but with different social rank. The Nile tilapia shows territorial behavior in nature and the establishment of stable long-term social hierarchies in aquaculture [45‒47], which makes it a suitable model for our investigations. In such a hierarchy, social status can strongly influence physiology and reproductive behavior [48] because only the dominant male can hold its territory and gets the opportunity to reproduce. Besides the obvious distinct coloration of dominant and subordinate males, they show significant differences in gonadal development and hormone levels, which appoint toward a reduced steroidogenic and spermatogenic potential in subordinate individuals [46]. By analyzing the testicular and pituitary gene expression of dominant and subordinate Nile tilapia, we have already been able to identify differentially expressed (DE) factors along the BPG axis [49]. With that study, the evidence that steroidogenesis-related genes, growth factors, and gonadotropins could be involved in the establishment and maintenance of social dominance was strengthened. A short graphical summary of the results from this previous study is presented in Figure 1. The Nile tilapia strain we used does not carry the tandem duplication (amh∆Y and amhY) [49].

In the study presented here, we used the same set of samples to investigate the expression of miRNAs in the testes of dominant and subordinate individuals with the aim of finding miRNAs regulating reproduction and verifying their effects on differentially expressed genes. miRNA target prediction and correlation to differentially expressed genes from our previous study identified cellular processes involved in cell motility and communication, as well as metabolic pathways involved in the endocrine system. We also identified miR-499 as upregulated in dominant males, and we examined how miR-499 affects the BPG axis through its regulation of the central reproductive factor amh.

The miRNA expression analysis and validation presented in this study was performed with the same testis samples as in our previous study on transcriptome analysis in dominant and subordinate male Nile tilapia [49]. The details of the origin of our Nile tilapia brood stock as well as fish maintenance and sampling have already been extensively described there. In brief: O. niloticus used in this study were reared in 550 or 820 L tanks (Aqua Schwarz GmbH, Göttingen, Germany) in the institute’s fish facility in family groups of 20–30 siblings. The fish were housed at 26°C in a recirculating system (water exchange 5 percent per day) and a photoperiod of 14 h light and 10 h darkness. The electric conductivity was around 600 µS/cm. The fish were fed ad libitum with F-3P Optiline BC (6 mm, Skretting, Norway) once a day. Twice a week, the food was supplemented with vitamins (Multivit, Wiegandt, Germany) and once a month, trace elements (Tracevit, Wiegandt, Germany) were added to the water. This study used 3 dominant and 5 subordinate adult males for RNA sequencing and one additional dominant and subordinate adult male for validation (online suppl. Table S1; for all online suppl. material, see https://doi.org/10.1159/000546304) to increase statistical power. Fish were selected to ensure a comparable age distribution between groups, with ages ranging from 9 to 42 months. For the primary cell culture (PCC) experiment, the testes of two dominant males were used after morphological examination to be well developed (online suppl. Table S1). The rearing density for the adult fishes was 8.3–13 kg/m³ for the younger fishes (9–12.5 months) and 18–21 kg/m³ for the older families (28–42 months) at the time of sampling (online suppl. Table S1).

For sample collection, the animals were anaesthetized with 0.2 g/L benzocaine and euthanized by severance of the spinal cord. The testes were removed and rinsed in sterile PBS. After that, slices were cut from the middle sections of the organs, which were snap frozen in liquid nitrogen and stored at −80°C until the RNA extraction was performed. For the PCC experiment, testes were processed directly after PBS rinse without freezing.

The RNA extraction is described in Thönnes et al. [49]. In short: RNA integrity was verified by 2100 Bioanalyser (Agilent, USA) at the Max-Planck Institute of Molecular Cell Biology and Genetics Dresden (Germany). RIN values ranged between 6.5 and 9.1. RNA sequencing was performed at the Deep Sequencing Facility of the CMCB technology platform (Dresden, Germany) using total RNA (50 ng) that was isolated from the testes of the subordinate fish 7S, 9S, 17S, 30S, and 31S, as well as the dominant fish 8D, 18D, and 32D (online suppl. Table S1). Total RNA was converted into small RNA libraries using the BioScientific NEXTFlex Small RNA-Seq Kit (v3). Enrichment of miRNA molecules from amplified cDNA was done by size selection using an 8% PAGE 10-well gel in the respective range (150 bp band). Libraries were barcoded, pooled and sequenced for 10 million reads on an Illumina NextSeq machine in 75 bp SE mode. From each library, approximately 9.3–15.9 million reads were obtained. For miRNA analysis, the NextFlex small RNA adapter was trimmed with cutadapt (v1.8.1) [50]. Reads without an adapter sequence were discarded. The unique molecule identifier introduced by the NextFlex protocol was removed as well. Trimmed reads were mapped against tRNA sequences as predicted by tRNA scan [51] as well as known rRNA sequences and the whole genome (Ensembl) with bowtie [52]. Reads mapping to non-miRNA loci were removed and filtered to contain only read lengths of 18–25 base pairs. miRDeep_star [53] was used to predict expressed miRNA loci. miRNA annotation was downloaded from miRBase (Release 22) and converted to Ensembl reference to identify known miRNAs. For further analysis, loci of known miRNAs or with a miRDeep_star score of at least 7 were used. Using Infernal 1.1 [54], precursor sequences of putative novel miRNAs were searched for known miRNA secondary structures in Rfam [55] and those with hits were included in the analysis. Counts were normalized by the counts-per-million (cpm) method. miRNA with multiple dominant or subordinate samples having a cpm below 1 were discarded to remove very low-expressed miRNAs. Raw counts were normalized and tested for differential expression using DESeq2 v1.20.0 [56] with Benjamini-Hochberg multiple testing correction. miRNAs with a false discovery rate below 10% (p_adj <0.1) were considered significantly DE. Fastq files are available at the NCBI Sequence Read Archive (SRA) under the BioProject ID PRJNA817031.

The RNA was isolated from the testis slices using peqGOLD TriFast (Peqlab, Germany) according to the manufacturer’s protocol. The extracted RNA was quantified using a NanoDrop 1000 Spectrophotometer (Peqlab, Germany) and the integrity was determined on a 1.5% agarose gel. For miRNA sequencing validation RNAs from the subordinate individuals 7S, 9S, 17S, 30S, 31S, and 194S and dominant individuals 8D, 18D, 32D, and 195D (online suppl. Table S1) were used.

For poly(A) reverse transcription, total RNA (1 µg) was reverse transcribed using the qScript microRNA cDNA Synthesis Kit (Quantabio, USA) according to the manufacturer’s protocol. RNA isolated from LNCaP cells (human prostate cancer cell line) was used for performing the internal positive control of the kit by proving successful cDNA synthesis with the provided primers for the detection of the human small RNA SNORD44.

For stem-loop reverse transcription, the RNA samples were treated with DNaseI (recombinant, RNase-free, Roche, Switzerland) and reverse transcribed by using pretreated stem-loop primers, which were designed after Tong et al. [57]. The approach of pre-treatment was based on the protocol created by Kramer [58] without the usage of mineral oil and performed in a Biometra T-Gradient Cycler (Biometra GmbH, Germany). Each stem-loop primer (100 µm) was incubated at 95°C for 10 min before the temperature was stepwise reduced by 2°C and maintained for 1 min until 75°C was reached. The temperature was maintained at 75°C, 68°C, 65°C, and 62°C for 1 h each and at 60°C for 3.5 h. The pretreated primers were stored at −20°C. For reverse transcription, 1 µL of the respective stem-loop primer (1 µm) and 1 µg of the DNaseI-treated RNA were incubated at 16°C for 30 min in a total volume of 9 µL. After the addition of reaction buffer (1x) and reverse transcriptase (200 units) (both: RevertAid Reverse Transcriptase kit, Thermo Fisher Scientific, USA) as well as dNTP Mix in a final concentration of 1 mm (Thermo Fisher Scientific, USA), the mixture with a total volume of 20 µL was incubated at 42°C for 1 h. The reaction was terminated by heat inactivation of the enzyme according to the manufacturer’s protocol.

Each real-time quantitative PCR (RT-qPCR) reaction was carried out in a total volume of 10 µL containing 5 µL PowerUp SYBR Green Master Mix (Thermo Fisher Scientific, USA). For miRNA amplification, final primer concentrations of 350 nm were used. 18S-rRNA and gapdh served as reference genes for normalization and were amplified using final primer concentrations of 250 nm. All primers and their efficiency are given in online supplementary Table S2. The reactions contained different amounts of cDNA, which are listed in online supplementary Table S3. RT-qPCR was performed in triplicate for each sample using a QuantStudio 5 Real-Time PCR System (Thermo Fisher Scientific, USA) and the following parameters: 50°C for 2 min, 95°C for 5 min, followed by 40 cycles at 95°C for 15 s and 60°C for 1 min. Melting curve analysis of the amplification products was performed in the range of 60°C–95°C. All analyzed products showed one distinct peak (results not shown). No-template controls were always included. Serially diluted cDNA was applied for standard curves to determine primer efficiencies, which lay between 97% and 107%. Ct values with a difference of more than 0.5 in relation to the average of the other two values of the respective triplicate were excluded. Relative miRNA expression levels were calculated after Taylor et al. [59].

Previously performed RNA-Seq of the samples (BioProject ID PRJNA817031, [49]) was recomputed using the updated O_niloticus_UMD_NMBU Assembly (GCF_001858045.2, NCBI Annotation Release 104). Raw reads were trimmed with Trimmomatic v0.36 [60] using the following parameters: ILLUMINACLIP:TruSeq3-SE:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36. Mapping of trimmed reads to the reference and abundance estimation were performed by RSEM v1.3.1 [61] using the integrated Bowtie 2 v2.3.5.1 mapping algorithm [62]. Differential gene expression was calculated using the DESeq2 R package v1.38.3 [56]. Normalized counts of samples from dominant and subordinate animals were tested using the Wald test, including age group as a covariant, and corrected for the multiple testing problem with Benjamini-Hochberg. Genes with a false discovery rate below 5% (p_adj <0.05) were considered significantly DE (see online suppl. Table S4).

Possible targets of the DE miRNAs were predicted with miRmap 1.2.0 [63] using the REST API. Mature sequences of significant upregulated miRNAs were aligned against the 3′-UTR sequences (O_niloticus_UMD_NMBU Assembly, GCF_001858045.2, NCBI Annotation Release 104) of significant downregulated mRNAs and vice versa. To examine miRNA regulation of amh, all known mature miRNA sequences (miRBase, Release 22.1) were aligned to its 3′-UTR sequence. Furthermore, 3′-UTR sequences of amh from selected species (NCBI Annotation Release 104) were analyzed against their respective homologous miR-499 sequence (miRBase, Release 22.1), to examine if the effect prediction is identified across different vertebrates. Genes were classified as miRNA targets if a miRmap score was computed.

Gene Ontology and KEGG enrichment analysis were performed using the STRING database v12.0 [64] by uploading the NCBI identifiers of all predicted miRNA targets. Of the 1,858 miRNA targets input, 1,440 were identified and compared to the whole genome Nile tilapia dataset. Gene Ontology term and KEGG pathway enrichments with a false discovery rate below 5% (p_adj <0.05) were considered significantly enriched.

HeLa cells were maintained in DMEM/Ham’s F-12 with stable glutamine (Biochrom, Germany) supplemented with 10% fetal bovine serum (Gibco via Thermo Fisher Scientific, USA), 1% non-essential amino acids (Biochrom, Germany) and 1% Penicillin/Streptomycin (10,000 U/mL, 10,000 μg/mL; Biochrom, Germany). Cells were cultivated at 37°C and 5% CO₂.

The UTR of amh (Ensembl Release 103: ENSONIG00000004781) was amplified using the primers amhUTR-F and amhUTR-R (online suppl. Table S2) with 0.5 units GoTaq^® G2 Flexi DNA Polymerase (Promega, USA) from 50 ng Nile tilapia genomic DNA. The product comprised the artificially added restriction sites (DraI and XbaI) and 601 bp of the amh-3′-UTR downstream of the stop codon and was ligated into a DraI-XbaI opened pmirGLO Dual-Luciferase miRNA Target Expression Vector (Promega, USA). After transformation in E. coli and plasmid isolation, the identity of the recombinant plasmid pmirGLO_amhUTR was validated by sequencing using the primers T7-Terminator and pmirGLO-F (online suppl. Table S2). Sequencing results showed 100% base pair identity with ENSONIG00000004781 (O_niloticus_UMD_NMBU, version 103).

Transfections were carried out in white 96-well cell culture treated plates (Nunc, Thermo Fisher Scientific, USA) with Attractene Transfection Reagent (QIAGEN, Germany) according to the manufacturer’s fast-forward protocol. Per well, 40,000 HeLa cells in a volume of 100 µL were co-transfected with 100 ng of the respective plasmid (pmirGLO or pmirGLO_amhUTR, plasmid DNA prepared with ZymoPURE II Plasmid Maxiprep Kit, Zymo Research Europe GmbH, Germany) in 50 µL of DMEM/Ham’s F-12 with stable glutamine (Biochrom, Germany) and 25 nm, 50 nm or 100 nm of the Agomir or negative control Agomir MAH00000 (both: Applied Biological Materials, Canada) by using 0.75 µL Attractene Transfection Reagent. The transfected cells were incubated at 37°C and 5% CO₂ for 24 h. All transfections were performed in technical triplicates and at least three biological replicates.

For luminescence measuring, the Dual-Glo Luciferase Assay System (Promega, USA) and a plate reader (Tecan, Austria) were used. The transfected cells were washed twice with Hanks′ Balanced Salt Solution (Biochrom, Germany). Afterward, 50 µL fresh cell culture medium, which had been brought to room temperature, was added to the cells. The Dual-Glo Luciferase Assay System was used according to the manufacturer’s protocol, but the volumes of Dual-Glo Reagent and Dual-Glo Stop & Glo Reagent were reduced to 50 µL. To remove air bubbles, the plates were centrifuged shortly. After 20 min, the firefly or Renilla luminescence was measured with an integration time of 1,000 ms.

The testes from two fish were processed according to Tokalov and Gutzeit [65]. In brief, slices from the middle part of the removed testes were cut into 2 mm³ pieces using scissors and washed three times with PBS. PBS was then replaced with Leibowitz’s L-15 medium (Pan Biotech, Germany) and the suspension was centrifuged at 200 g for 3 min. The supernatant was discarded, and 15 mL of Leibowitz’s L-15 medium was added. The suspension was centrifuged again, and 10 mL Leibowitz’s L-15 medium containing 10 mg collagenase D (11088858001, Roche, Switzerland) and 2 mg DNaseI (10104159001, Roche, Switzerland) were added. The suspension was incubated for 2 h at 26°C at 25 rpm on an orbital shaker and the tissue was manually pipetted up and down every 20 min. After incubation, 10 mL of Leibowitz’s L-15 medium was added, and the cell suspension was successively filtered through CellTrics^® meshes with 150 µm, 50 µm, and 30 µm pore sizes (Sysmex Partec, Germany), respectively. To reduce the proportion of spermatozoa and haploid spermatids significantly, a Percoll™ (Cytiva, Sweden) gradient centrifugation step as described by [65] was performed. The cells obtained were washed with and resuspended in 10 mL complete culture medium (Leibowitz’s L-15, 10% FBS [Thermo Fisher Scientific, USA], 1% Penicillin/Streptomycin [10,000 U/mL/10,000 μg/mL, Biochrom, Germany]). The cells were then examined under a light microscope and the cell density was determined using a counting chamber. Spermatozoa or haploid spermatids were no longer visible. 2 mL cell suspension each were seeded into 5 wells of a 6-well plate (92006, TPP AG, Switzerland) with a cell density of 1.72 × 10⁶ cells/mL. After an initial incubation at 26°C for 24 h, the supernatants were discarded, the adherent cells were washed with sterile PBS and then overlaid with 1 mL of fresh culture medium. In total, 410,000 cells per well were determined after trypsination of a representative well. For treatment with oni-miR-499 agomir and the negative control agomir MAH00000 (random sequence; both: Applied Biological Materials, Canada), a transfection mix was made as follows: 100 µL of Leibowitz’s L-15 medium (without serum, 1% Penicillin/Streptomycin) containing 1.1 µm agomir or control agomir and 4 µL Attractene Transfection Reagent (QIAGEN, Germany). The transfection mix was vortexed for 20 s, pre-incubated at RT for 10 min, and added drop by drop to the 24-hour-old adherent cells. The cells were then incubated for 24 h at 26°C in technical duplicates.

RNA extraction was performed as described by [49]. The total RNA was quantified using a Nanophotometer NP80 (Implen, Germany), and the integrity was determined on a 1.5% agarose gel. The extracted RNA was digested with DNaseI (recombinant, RNase-free, Roche, Switzerland), and the absence of genomic DNA was ensured by an intron-spanning PCR for the housekeeping gene actin (see online suppl. Table S2). For poly(A) reverse transcription, 1 µg of total RNA was processed using the High Capacity cDNA Reverse Transcription Kit (Thermo Fisher Scientific, USA) according to the manufacturer’s protocol.

In addition to a morphological observation of the primary cells by phase-contrast microscopy (Axiovert 135, Zeiss, Germany), the presence of various testis cell types was analyzed with marker genes for somatic and germline cells using qualitative PCR and gene specific primers (see online suppl. Table S2). Standard PCR conditions using GoTaq^® Polymerase (Promega, USA) according to the manufacturer‘s instructions were applied in a total volume of 15 µL (12–45 ng template, 0.5 U GoTaq^® Polymerase, primers 500 nm f.c., 1.5–2 mm MgCl₂ f.c.; annealing temperature of 58–62°C and 36 cycles). PCR products were confirmed by gel electrophoreses using a 1% agarose gel.

RT-qPCR for amh transcript detection was performed as described above with the following changes: 18S-rRNA served as reference gene for normalization, and the final primer concentration was 125 nm. Detailed information on the primers used is given in online supplementary Table S2. The amount of cDNA used per reaction was 12 ng.

For statistical analysis, the software OriginPro 2019 (OriginLab Corporation, USA) was used. Normal distribution was checked using the Shapiro-Wilk test. Regarding luciferase assay experiments, oni-miR-499 transfected cells were set in proportion to the agomir negative control transfected cells. Statistical differences of the relative luciferase activity were calculated by a one-sample t test. A p value below 0.05 was considered significant.

Significant differences between the log2-transformed normalized miRNA expression data of the dominant and subordinate Nile tilapia males as obtained by RT-qPCR was checked by a two-sample t test. p values below 0.05 were considered significant. The same parameters were applied for the evaluation of the RT-qPCR analysis of amh expression in the PCCs.

Currently, 695 mature miRNA sequences from O. niloticus are documented in miRBase (Release 22.1). Using RNA sequencing, we identified 262 mature miRNAs in testis samples of Nile tilapia individuals living in stable hierarchies (Fig. 2a). Among them, we found 158 known miRNAs that were documented in miRBase (Release 22) and 104 putative novel miRNAs (Fig. 2b). Samples of dominant and subordinate males clustered during principal component analysis (Fig. 2d) and further analysis revealed that 8.8% (23 miRNAs) of all miRNAs were DE between the gonads of dominant and subordinate individuals (p < 0.1) (Fig. 2c). 17 of those 23 DE miRNAs were identified as being upregulated in dominant males. All of them could be found in miRBase (Release 22), whereas only 3 out of the 6 miRNAs that were found to be upregulated in subordinate individuals were already known. Among the 17 miRNAs upregulated in dominant males, 7 show a log2 fold change value above 1. Only the oni-miR-499 exhibited a log2 fold change exceeding the value of 2 (Fig. 3b). Based on these data, this miRNA could be detected as the most upregulated one in the testes of dominant tilapia males. Complete sequencing data can be found in online supplementary Table S5.

For the validation of our RNA-Seq results, we have chosen three miRNAs that may play a role within the BPG axis in fish at the level of the gonads. One of them is oni-let-7d, for which the regulation of the androgen receptor alpha was predicted in rainbow trout [67]. The next one, oni-miR-22a, could be involved in male sex determination and differentiation as well as the process of spermatogenesis because there is some evidence for its targeting of gonadal somatic cell-derived factor in carp [68]. The third, oni-miR-499, is interesting because of its considerable upregulation in dominant males in our study (Fig. 3b) and the prediction that the amh-3′-UTR could be a target of that miRNA (online suppl. Table S6).

Due to the short length of the miRNAs, reverse transcription requires elongation of the miRNA template. We decided to apply two different reverse transcription methods to detect mature miRNAs (Fig. 4). In the first one, extension was carried out using the poly(A) method and the second one involved the use of a stem-loop primer for cDNA synthesis, which has the advantage of high specificity and can discriminate even between small differences in nucleotide sequence such as in the let-7 family [69]. In brief, poly(A) reverse transcription results in the addition of a poly(A) tail to all miRNAs, which they naturally lack. An oligo-dT primer is used for cDNA synthesis and this primer has additional nucleotides including the sequence for a universal primer. Now, RT-qPCR can be performed with a miRNA-specific forward primer and a universal reverse primer (Fig. 4a). Stem-loop reverse transcription is based on the use of stem-loop primers for cDNA synthesis. The specific structure is achieved by two complementary sections that allow the formation of the stem, which includes a miRNA-specific sequence at its 3′-end. The loop further contains the sequence of the universal reverse primer that is used for RT-qPCR along with the miRNA-specific forward primer (Fig. 4b).

With these two methods, the expression of the three selected miRNAs (oni-let-7d, oni-miR-22a, and oni-miR-499), all of them upregulated in the testes of dominant males (Fig. 3a), were re-examined. In comparison to the RNA-Seq experiment, one individual of each rank was added for heightening of statistical power (online suppl. Table S1). Statistical analysis was performed on log2-transformed normalized data (online suppl. Table S7). The significant upregulation of oni-miR-499 and oni-miR-22a in dominant males was confirmed by both RT-qPCR methods (Fig. 4c). The significant differential expression of oni-let-7d was only verifiable by using the poly(A) method (Fig. 4c). This may be attributed to the occurrence of other oni-let-7 family members which were not amplified by using the more specific stem-loop primers [69].

For further insight into the miRNA regulated differences in dominant Nile tilapia, targets of the DE miRNAs were predicted (see online suppl. Table S6). Gene ontology and KEGG pathway enrichment analysis of the miRNA targets (Fig. 5) revealed the enrichment of cellular processes, such as focal adhesion (GO:0005925/onl04510) and regulation of actin cytoskeleton organization (GO:0032956/onl04810), as well as environmental information processing, such as integrin-mediated signaling pathway (GO:0032956), ECM-receptor interaction (onl04512) and neuroactive ligand-receptor interaction (onl04080). Additionally the enrichment of biological processes concerning organismal systems was observed as, for example, multicellular organism development (GO:0007275), anatomical structure development (GO:0048856), and system development (GO:0048731). Furthermore, pathways involved in the metabolism of carbohydrates, lipids, amino acids, and nucleotides were enriched in dominant males, as well as processes of the endocrine system, such as the PPAR signaling pathway (onl03320), steroid biosynthetic process (GO:0006694), and C21-steroid hormone metabolic process (GO:0008207). Complete results of the gene ontology and KEGG enrichment analyses are given in the online supplementary Table S8.

The previous study showed that dominant and subordinate males have different steroid hormone plasma levels [49]. To examine the role of miRNAs in their production, major metabolic pathways were analyzed, describing the metabolism of cholesterol to testosterone, 11-ketotestosterone, and estradiol among others. Multiple genes coding for enzymes involved in the steroid hormone biosynthesis were found to be upregulated in dominant males, such as cyp11a1/b1, cyp17a2, cyp21a2, hsd3b, hsd17b3, and hsd11b2, some of which were identified as possible miRNA targets, and we have shown these miRNAs downregulated in dominant males in Figure 6. Although not all proteins involved in the production of testosterone and 11-ketotestosterone were identified in the data or were not DE, our results indicate that the upregulation of key players on a transcriptional level via the downregulation of specific miRNAs, as shown in Figure 6, may result in higher hormone plasma levels. The gene (cyp19a1 alias aromatase) coding for Cyp19a1, which is responsible for the conversion of testosterone to estradiol, is not significantly DE between dominant and subordinate males at the transcript level, even though higher estradiol plasma levels in dominant males occurred [49].

Since miRNAs are known for their ability to affect various biological processes, including reproduction, we wanted to analyze their potential role in targeting Amh as a central factor within the BPG axis in teleost fish. Putative miRNAs targeting the 3′-UTR of amh were identified using miRmap and raw values were computed to quantify the repression strength. The lower the raw miRmap score, the higher the predicted repression. Overall, 183 miRNAs were identified (Fig. 7c, complete results in online suppl. Table S9), among them six miRNAs were found to be DE between dominant and subordinate males: oni-miR-184a, oni-miR-184b, oni-miR-499, oni-miR-135b, oni-miR-7, and oni-miR-139 with miRmap scores of 0.05, 0.07, 0.10, 0.10, 0.10, and 0.12, respectively. Though oni-miR-139 has the lowest miRmap score (0.05) and therefore the highest amh repression as predicted in silico, oni-miR-499 with miRmap score 0.10 was picked for further analysis because it showed the highest expression of the DE miRNAs in dominant males (Fig. 3b).

A putative binding site for miR-499 in the amh-3′-UTR could also be identified in the related Cichliformes Neolamprologus brichardi and Pundamilia nyererei, with miRmap scores of 0.12 and 0.09, respectively. In the amh-3′-UTR of Danio rerio, a member of the Teleostei, binding motifs for dre-miR-499 could be identified (miRmap score: 0.33); however, this was not the case in other teleosts we analyzed, such as Oryzias latipes, Gadus morhua, and Ictalurus punctatus. Other commonly analyzed species such as Homo sapiens, Mus musculus, Rattus norvegicus, and Xenopus tropicalis also did not show a putative binding site but Ornithorhynchus anatinus (Platypus) did (miRmap score: 0.19).

To prove the predicted interaction between oni-miR-499 and the Nile tilapia amh-3′-UTR, we established a Dual-Luciferase Assay in HeLa cells. This assay makes use of two luciferase reporter genes that are provided by the pmirGLO vector. Firefly luciferase acts as an experimental reporter, while Renilla luciferase serves as a control reporter for normalizing measured values. If a miRNA binds to the introduced target sequence, the firefly luciferase expression is usually reduced, lowering its bioluminescence signal. We performed transfection experiments with the pmirGLO_amhUTR plasmid to analyze the effect on the luciferase activity when co-transfected with three different oni-miR-499 agomir or agomir negative control concentrations and the pmirGLO plasmid to check whether miRNA binding occurs without the presence of the amh-3′-UTR. We found a significant downregulation of the relative luciferase activity with all agomir concentrations when co-transfected with pmirGLO_amhUTR, while the cotransfection with pmirGLO showed no significant change in luciferase activity (Fig. 8, complete data in online suppl. Table S10), indicating the interaction of oni-miR-499 with the amh-3′-UTR in vitro.

In order to investigate the interaction between oni-miR-499 and the amh-3′-UTR in Nile tilapia testis, primary testis cells were isolated and cultured (Fig. 9a). After 24 h, the cells start to attach and irregularly shaped cells as well asfloating aggregates, preferably consisting of germ line cells from spermatogonia to spermatids [65], occurred (Fig. 9b). After removal of the supernatant, the remaining attached cells and the associated aggregates were exposed to oni-miR-499 agomir for another 24 h (Fig. 9b). The oni-miR-499 agomir/ATR transfection reagent mix resulted in a significant decrease in amh expression (p = 0.034) in comparison to the control agomir (Fig. 9d). Marker gene analysis from a parallel experiment without the addition of ATR transfection reagent confirmed the identity of the adherent cells in the PCC after 48 h of cultivation (Fig. 9c).

Expression patterns and functions of miRNAs are a young but rapidly growing field of research. In particular, next generation sequencing and bioinformatics analysis foster our understanding of the mechanisms of miRNA action. miRNAs represent an additional level of regulation of gene expression and therefore could be involved in all aspects of reproduction and at all levels of the BPG axis [4, 6]. With regard to the gonads, the data available for fish are still far behind that of mammals [7‒10, 70], but both, conserved miRNAs as well as teleost- or species-specific miRNAs, have been identified in Nile tilapia [71].

Current work on miRNA in fish gonads preferentially refer to the ovary, and functional data are rare [72]. Knowledge on miRNAs in Nile tilapia gonads so far has focused either on a general comparison of ovaries and testes in adults [41, 44] or on immature gonads of embryos and larvae during the sex determination period of genetic males and females [42, 43]. In contrast to these studies, we used an experimental system in which we compare the same organ (testis) in a different activation state caused by the different social status of the males. The selected males are stable in their position for at least 4 weeks until months. Potentially ascending or unclearly identifiable males were not taken into account, so that stable, different hormonal states exist for the two compared groups of dominant and subordinate males [49]. This system significantly differs from analyses using gene knock outs because spermatogenesis is still active in subordinate Nile tilapia males. Using genetically non-interrupted signaling pathways is a promising approach to identify regulatory pathways involved in the control of spermatogenesis and testis development. In the dominant males used in the presented study, the gsi was increased roughly 3–5-fold compared to the subordinates. The same applies to hormone levels; 11-KT and testosterone were significantly elevated in dominants, at least 3-fold up to more than 10-fold. Lh was always at least twice as high in dominant males [49]. This is in concordance with the histological observations and previous gene expression analyses. Huge islets of Leydig cells and a thickened tunica albuginea with many Cyp11b (key enzyme for 11-KT biosynthesis)-positive myoid cells characterize the testis from long-term dominant males [46]. Elevated levels of typical Leydig cells genes (e.g., sf1, star2, hsd17b3, hsd3b1, hsd11b2, cyp11a2, cyp11b1, cyp11c1, cyp17a2) and reduced expression of amh and amhrII indicate strongly activated spermatogenesis in dominant males, too [46, 49]. Other spermatogenesis-related genes upregulated in dominant males are those coding for androgen, estrogen, prostaglandin, gonadotropin, and retinoic acid receptors, as well as a number of Igf-binding proteins, Wnt pathway-related components, Hedgehog receptor smoothened and insl3. The dataset further indicates a higher metabolic activity in testes of dominant males and that estrogen signaling also plays a role in the fish testis [49]. All of this emphasizes the suitability of our experimental system to detect regulative differences in testis development and the entire BPG axis.

Using this socially determined system in a comparative approach, we identified 262 mature miRNAs (158 known and 104 novel) in the testes of Nile tilapia, of which 23 were proven to be DE in dominant and subordinate males. As we had already observed for differential mRNA expression [49], we see a higher proportion (75%) of upregulated miRNAs in the testis of dominant males than downregulated miRNAs (Fig. 3). All 17 upregulated miRNAs are registered in miRBase (Release 22), whereas 3 out of the 6 downregulated miRNAs are novel (novel-mir-510, novel-mir-352, and novel-mir-1576).

We assigned the novel-miRNAs to a miRNA family by comparing their sequence with the Rfam database. Novel-mir-510 matches miR-210, novel-mir-352 matches miR-1388, and novel-mir-1576 fits to miR-268. While the first two match with a very high probability, the assignment of novel-mir-1576 (miR-268) still remains subject to uncertainty (Rfam_evalues in online suppl. Table S5). Literature search did not reveal any relation of miR-268 to gonadal development or any other tissue, whereas the other two had already been previously identified in the gonads of fishes and also various branches of the animal kingdom, from protostomes (Chinese mitten crab) and basal deuterostomes (Sea urchin) to mammals (online suppl. Table S11 for details). In the testis of Japanese flounder, miR-1388 affects nectin2l expression and therefore spermatid maturation as in mammals [73]. More examples for miR-1388 (novel-mir-352) in online supplementary Table S11 are based on the expression information available from the supplementary material of gene expression studies used testicular tissue. In contrast to miR-1388, miR-210 appears to be a more widespread gonadal miRNA, frequently found in literature (online suppl. Table S11). miR-210 is related to hypoxia-coupled cellular responses and thus also to male infertility in mammals. For example, miR-210 is upregulated when spermatogenesis is impaired in a varicocele rat model [74]. Two studies in fish using ovarian cell cultures showed connection of miR-210 to apoptosis regulation in marine medaka and diminished oocyte meiosis-related genes in zebrafish [75, 76], respectively. In this zebrafish, ovarian cell culture an increase of igf2b was observed after miR-210 addition, which is the opposite of what has been reported for human spermatogenesis, where upregulation of miR-210 in testis correlates with IGF2 repression and male infertility [77]. The latter is consistent with the results from our study, where we see reduced miR-210 levels in correlation with elevated igf2 expression and activated spermatogenesis in the testis of the dominants. Maybe, miR-210-igf2 interaction could therefore represent a conserved mechanism in the regulation of spermatogenesis in vertebrates. In addition, miR-210 (novel-mir-510) is the only miRNA of the 6 downregulated miRNAs in dominants that does not appear to interfere with the regulation of steroid biosynthesis genes as shown by Figure 6, while the other two novel miRNAs identified in our study could influence spermatogenesis in this way.

For the other three miRNAs that were downregulated in the testes from dominant individuals (miR-462, miR-25, miR-135b), there are already examples for testis expression (online suppl. Table S11). We only found fish examples for miR-462. This seems consistent because miR-462 has been described as a teleost-specific miRNA [78]. Expression was reported in Tiger pufferfish and Zebrafish testes as well as in Zebrafish spermatozoa [79, 80], respectively. Furthermore, Presslauer et al. [81] detected sex-specific expression of miR-462-5p in testes from 9- and 12-week-old Zebrafish, and a testis-biased expression in the gonads of the Amur sturgeon has also been described [82]. This shows the occurrence of miR-462 also in Chondrostei outside the Teleostei. Testis expression for miR-25 and miR-135b can also be found in the large phylogenetic superclass of bony fish, although the examples for miR-25 are still limited. miR-25 plays a role in testis maturation in Atlantic salmon [83] and Rainbow trout [67] and was also found in the spermatozoa of Zebrafish [80]. For mammals, miR-25 was detected in murine Sertoli cell cultures and in the testes of pigs and humans (online suppl. Table S11). From experiments from humans or non-gonadal tissue from fish, one could deduce, for example, that miR-25 has a possible influence on apoptosis induction and cell cycle re-entry inhibition [84, 85]. Similar to miR-25 expression in Atlantic salmon and Rainbow trout, miR-135b correlates with the pubertal stage of testis development in these two species [67, 83]). In the Nile tilapia testis, miR-135b occurs in higher abundance than in the ovary [41]. Recently, an increased expression of miR-135b in the testis of the protogynous hermaphroditic ricefield eel was reported, and an involvement of miR-135b in gonadal transformation to the testis was deduced [86]. Due to its wide distribution in the testes of vertebrates and their occurrence in the testis of a crab, miR-135b (as well as its homologs and ancestors) could represent a well-conserved testis-related miRNA (online suppl. Table S11 for further examples).

We found several reports for testis expression in fish and mammals (online suppl. Table S11) for all 17 upregulated miRNAs in our study, but functional studies are still limited. For 14 from this pool of upregulated miRNAs we found examples for testis expression also outside the vertebrate lineage in several branches of the Bilateria, indicating considerable conservation of these miRNAs with respect to testis expression. Only 3 (miR-499, miR-148, and miR-24) of the 17 upregulated miRNAs are still without precedent outside of vertebrata in our literature search.

The attention to the role of miRNAs in reproduction and testis function is based on their involvement in various fertility disorders in humans as well as in testicular germ cell tumors and aims were finding infertility markers and to be able to fight tumors. The majority of literature on the involvement and function of miRNA in the testis comes from this field and results in a series of male fertility miRNAs from mammals [7, 87, 88]. Among them, there are several miRNAs that we were able to identify as upregulated in our functional approach with dominant and subordinate fish, which are miR-10a, miR-125a, miR-184, miR-22, miR-24, miR-30a, and let-7 family members. A proposed function for the miR-125a, miR-30a, miR-184, and let-7 family in mammalian testis is a role in tumor suppression and the regulation of cell differentiation and proliferation can be assumed [87]. Others, like miR-22 in sheep, for example, could interfere with estrogen signaling in fetal mouse Sertoli cells [89] and miR-22 is regulated by BDNF in mouse Leydig cell cultures [90]. Further functional experiments from mouse revealed that miR-10a inhibits Rad51, which is necessary for double strand break repair and correct proceeding of meiosis [91]. For miR-24, an elevated expression level in a mole was reported during breeding season and an interaction with the polyubiquitin-B gene, a factor that marks proteins for degradation [92]. If available, further information on the function of the identified miRNAs is also listed in the Remarks column in online supplementary Table S11.

The miRNA most DE in the testes of dominant and subordinate Nile tilapia males is miR-499, which has so far been primarily associated with skeletal and cardiac muscle development in vertebrates ([93] for review; [94] – Nile tilapia; [95] – human; [96] – chicken). miRNA-499 is part of the so called myomiR regulatory network which encompasses muscle tissue enriched miRNAs and conserved regulators of muscle development, and which promotes myoblast proliferation and differentiation of slow-twitch muscle fibers [97, 98]. miR-499 is an intronic miRNA in myosin heavy chain 14 genes (myh14/MYH7b) reflecting a close coupling of the gene and its miRNA with muscle fiber development. During teleost evolution, the duplication of the myh14/miR-499 locus was followed by several events of gene loss and an additional duplication in zebrafish [93], explaining different genomic constellations in ray-finned fish. As examples, there are two myh14 copies in Nile tilapia, one of them with the miR-499 intronic integration; and there are three copies in zebrafish, each with an intronic miR-499 insertion; and finally, both myh14 genes were lost in medaka, but a miR-499 coding sequence inclusive conserved regulatory elements remained [93]. This example for divergent resolution of myh14 and miR-499 genes underlines the importance of this miRNA in vertebrates, and it is a common phenomenon in teleost genome evolution [99, 100].

For testicular miR-499 expression, we found more than 20 examples within the vertebrates (online suppl. Table S11). A recent work on a mole (Eospalax baileyi) showed increased miR-499 expression during the breeding season, which corresponds to an activation of spermatogenesis [92]. Expression in the ovary has also been reported ([44] – Nile tilapia; [101] – salamander; [102] – human). However, a direct interaction of miR-499 with human AMH expression is not expected because inside AMH-3′-UTR of humans and other mammals we did not find binding motifs for miR-499 by our exemplary analysis (online suppl. Table S12). In contrast, we predicted miR-499 binding motifs for the amh-3′-UTRs in more ancient organisms like platypus and some teleosts including Nile tilapia. Our performed in vitro reporter gene assay and the detection of altered gene expression in primary testis cells after miR-499 application confirmed the connection between miR-499 and amh expression level in Nile tilapia. In addition, a total of 183 miRNAs are putative interaction partners with the amh-3′-UTR (Fig. 7c, online suppl. Table S9), and six of them are DE between dominant and subordinate males (upregulated: mir-499, miR-139, miR-184a, miR-184b, mir-7 and downregulated miR-135b in dominants). In vivo and in vitro experiments in common carp showed that miR-153b-3p suppresses amh expression [103]. We also found miR-153d in our predictions for the Tilapia amh-3′-UTR but did not find a differential expression, thus not supporting a testicular function in dominant and subordinate Nile tilapia males.

Interestingly, miR-499 should also have the potential to affect the amh copies amh∆Y and amhY, which are not present in the strain used in our study but in some other Nile tilapia strains [104, 105] because the 3′-UTRs are identical in this region for all three genes. miR-499 is actually expressed in the early Nile tilapia gonad at 5 dah [43], which would allow such an influence on that copies in the early gonads.

The interaction of miR-499 with the amh-3′-UTR that we have shown experimentally for Nile tilapia could be a regulatory mechanism limited to certain fish groups and others branches of the vertebrates like Monotremata to boost spermatogenesis via reduction of amh expression. Amh is essential for the balance of germ cell differentiation and inhibits steroid biosynthesis and action in fish [23, 25, 26]. In sum, miR-499 is emerging as a conserved testis/gonadal miRNA in vertebrates acting via different mechanisms in different animal groups.

Without functional experiments as described for miR-499 and miR-153b, bioinformatic target prediction is one option to generate hypotheses about the possible functions of identified miRNAs. In our analysis, we only considered the DE mRNAs and 23 DE miRNAs from dominants and subordinates, and among them, only the significant upregulated mRNAs to the significant downregulated miRNAs and vice versa for target identifications. However, miRNA-mRNA target prediction usually results in a plethora of hits, and in our case, mRNA target prediction revealed several hundred target mRNAs per identified miRNA despite the restriction to the same tissue type and previous filtering of data (online suppl. Table S6, S8). The reason for that complexity of miRNA networks lies in the fact that one miRNA can influence many different genes and that a single mRNA can be influenced by many different miRNAs at the same time [106].

With our pre-filtered data, we carried out the presented KEGG and GO term enrichment analyses of the identified miRNA targets. These data indicate molecular mechanisms associated with activated testis development and spermatogenesis in dominant individuals (Fig. 5; online suppl. Table S8). Enriched pathways of general metabolism and metabolism of carbohydrates, lipids, amino acids and nucleotides, as well as the circulatory system, match the highly developed gonads with high gsi in the dominants. The same holds true to the other enrichments of KEGG terms like focal adhesion (cellular communication), regulation of actin cytoskeleton (cell motility), ECM-receptor interaction (signal molecules and interaction), and calcium signaling pathway (signal transduction). From the pathways and terms presented in Figure 5, we selected genes that have already been linked to gonadal development and function in the literature. This is intended to indicate the potential of genes identified in our study and possibly regulated by miRNA in the testis, as well as to allow the derivation of hypotheses regarding their possible functions. We have included the results of our literature search as an extended discussion of this chapter here at the supplementary material (Text file S13).

Due to the activated BPG axis in dominant males, we expected terms and processes correlated to steroidogenesis in our KEGG and Gene ontology analyses. We did not get the KEGG pathway “steroid hormone biosynthesis” as an enriched pathway with the applied thresholds FDR ≤0.05 for significance, but we still achieved for this one a 2.9 times more frequent pathway annotation than expected in a dataset of this size. However, in the GO term analysis, the biological processes “steroid biosynthetic process,” “steroid metabolic process” and “C21-steroid hormone metabolic process” appear (online suppl. Table S8). GO terms can capture, in some cases, more narrowly defined processes than KEGG analyses [107]. In the gene lists within the enriched GO terms, typical genes for biosynthesis steps of steroid hormones are included. In order to visualize the connections between genes of individual steps along the biosynthetic pathway and the associated miRNAs found, we combined a scheme derived from the onl00140 KEGG pathway “steroid hormone biosynthesis” with information about the differential expression of mRNAs and miRNAs (Fig. 6). It is shown there that 5 out of 6 of the downregulated miRNAs found in dominants (miR-462, miR-25, miR-135b, novel-mir-352, novel-mir-1576) could be involved in the activation of genes for enzymes of steroid biosynthesis. This includes the C21 pathway from cholesterol via pregnenolone to cortisol, as well as the branch from the C19 steroids to testosterone (Fig. 6), [108]. Both hormones (testosterone and cortisol) are considerably involved in the regulation of spermatogenesis in the fish testis [109, 110]. It shows us which biosynthetic pathways of steroids are affected by the regulatory influence of miRNAs and thereby assigns a possible function to the downregulated miRNAs found in our study in the context of maintaining high androgen levels in stable dominant Nile tilapia males [46, 49, 108].

Numerous genes discussed here and in the supplemental material can be found in different pathways at the same time. This lie in the general complexity of miRNA/mRNA regulation as already explained above, as well as in the involvement of essential signaling pathways in the activation of spermatogenesis and their intrinsic overlap. Overall, we see a higher weighting of the upregulated mRNAs and thus also of the corresponding downregulated miRNAs in the in silico data (online suppl. Table S8). However, the true impact of each single miRNA has to be experimental proven in the future.

In this study, we used an experimental system with stable dominant and subordinate Nile tilapia males to identify miRNAs involved in the maintenance of active spermatogenesis and the implementation of the socially activated BPG axis in the testis of dominants. We have identified 23 DE miRNAs. Except for two, these miRNAs show a conserved testis expression in vertebrates, 17 also in basal branches of the animal kingdom. Within the two exceptions, miR-462 could be specific for fish and novel-mir-1576 could be a real novelty. The latter is one of three new predicted miRNAs in Nile tilapia alongside novel-mir-510 (miR-210) and novel-mir-352 (miR-1388). To identify the regulatory networks and effectors linked to differently developed testes, we used a targeted data filtering strategy (GO term and KEGG enrichment analyses with predicted miRNA target mRNAs based on differentially up-/downregulated miRNA/mRNA pools and vice versa). We have demonstrated a number of biological processes and provided examples from the literature for their involvement in controlling testis development and function. These include steroidogenesis, integrins and focal adhesion signaling, extracellular matrix interaction, various signaling pathways connected to differentiation of spermatogonia and meiosis, G-coupled receptors for Leydig and Sertoli cell functions, apoptosis inhibiting effects as well as possible negative feedback mechanism via PPAR signaling. There is a large overlap between the biological processes involved. Although our approach already achieves a certain functionality of the statements due to the similarity of the tissues, the identified biological processes and pathways are based on miRNA target predictions and there is a need for experimental evidence of the postulated interactions. A luciferase reporter assay with the amh-3′-UTR and a miR-499 agomir application on primary testis cells proved the interaction of miR-499 and Nile tilapia amh, most probably a mechanism to stimulate spermatogenesis via Amh reduction in dominant males. Our study is in line with others about seasonal, pubertal, social, and stress-induced changes in gonad development, and future studies and experimental models are needed to elucidate the importance of miRNAs in the control of male fertility.

We would like to thank Jannette Wober, Antje Beyer, and Franziska Hoffmann for cell culture and reporter gene assays support, Georg Kretzschmar for providing plasmids and helpful discussions, Christin Froschauer and Klaus Reinhardt for support with RT-qPCR, Nadja Zimmermann for fish care, and Alexander Froschauer for his continuous support in all matters. Many thanks to Andreas Petzold, Andreas Dahl, and the facility members of the DRESDEN-concept Genome Center at CMCB Technology Platform for RNA sequencing and initial data analysis.

All procedures regarding the animals and euthanasia were performed according to the EU directive 2010/63/EU guidelines as well as to the German national regulations and animal welfare (Tierschutzgesetz§11, Abs. 1, No. 1). Fish husbandry and procedures were approved by administrative regulations (Landesdirektion Sachsen, File number 25-5131/346/6). This approval process by the Landesdirektion Sachsen required prior examination and approval by the Technische Universität Dresden Animal Care and Use Committee (authorization number 25-5131/525/2).

The authors have no conflicts of interest to declare.

This work was supported by the Deutsche Forschungsgemeinschaft (PF683/5-1) to F.P., M.T., and J.E. and by funds from the “Forschungspool der DFG-Programmpauschale/BMBF-Projektpauschale” of the TU Dresden to R.P. The article processing charge (APC) was funded by the joint publication funds of the TU Dresden, including Carl Gustav Carus Faculty of Medicine, and the SLUB Dresden as well as the Open Access Publication Funding of the Deutsche Forschungsgemeinschaft. The funders had no role in the design, data collection, data analysis, and reporting of this study.

R.P.: data acquisition, bioinformatic computation and interpretation, manuscript writing, and creation of illustrations; J.E.: miRNA expression validation, luciferase assays, manuscript writing, and creation of illustrations; J.S.: primary cell cultures and expression analysis and creation of illustrations; M.T.: study design, sampling and dissection, and RNA extraction; F.P.: conceptualization and funding, sampling and dissection, data analysis, and manuscript writing. All authors approved the final manuscript.

Rebecca Prause and Josephin Eckart contributed equally to this work.

Fastq files we evaluated in our study are openly available at the NCBI Sequence Read Archive (SRA) under the BioProject ID PRJNA817031.