Identification and profiling of microRNAs during sheep’s testicular development

The normal development of the testis is essential for male reproduction, as it is the site of sperm production and a prerequisite for spermatogenesis. MiRNAs play crucial roles in various testicular biological processes, including cell proliferation, spermatogenesis, hormone secretion, metabolism, and reproductive regulation. In this study, we utilized deep sequencing data to analyze the expression patterns of small RNAs in testicular tissues of Southern × Hu sheep F1 hybrids at 0, 3, 6 months, and 1 year of age, thereby exploring the functions of miRNAs in testicular development and spermatogenesis. A total of 787 known miRNAs and 415 novel miRNAs were identified. We identified 217, 254, 405, 130, 305, and 138 DE miRNAs in the testes of M0 vs. M3, M0 vs. M6, M0 vs. Y1, M3 vs. M6, M3 vs. Y1, and M6 vs. Y1, respectively. GO annotation and KEGG pathway analysis of DE miRNA target genes revealed that target genes such as YAP1, ITGB1, DOT1L, SMAD4, and SOX9 may be involved in various biological processes, including reproductive pathways such as FOXO, Hippo, Wnt, cAMP, Rap1, and MAPK signaling pathways. The expression levels of 12 randomly selected miRNAs in testes at 0, 3, 6 months, and 1 year of age were detected by qRT-PCR, and the results were consistent with the sequencing data. This study characterized and investigated the differential expression of miRNAs in sheep testes at different developmental stages using deep sequencing technology. These findings will contribute to a deeper understanding of the functions of miRNAs in regulating testicular development and enhancing reproductive performance in male sheep.

Introduction

Hu sheep are mainly produced in Jiaxing and Taihu Lake regions of Zhejiang Province. As a unique breed of multi-parous sheep in our country, Hu sheep has many advantages such as perennial estrus, early sexual maturity, strong adaptability and suitability for crossing (1). The Southdown sheep, a short-wooled meat breed, is characterized by its dense, short, and light fleece, early maturity, ease of fattening, and tender meat, making it a preferred choice for crossbreeding sires (2). The enhancement of animal performance can be achieved through the introduction of breeds for crossbreeding improvement, a widely adopted strategy in animal production. Cross breeding can also combine the excellent characteristics of multiple varieties, create new characteristics that the original parents do not have, and enhance the vitality of the offspring (3). The Southdown × Hu F1 hybrids enhanced meat quality and increased lambing rates, retaining the high productivity of Hu sheep while inheriting the superior meat yield, quality, and rapid growth traits of Southdown sheep (4). As an important reproductive organ of male animals, the testis is mainly responsible for the production of sperm and androgens (5). Because normal testicular development is essential for breeding, studying testicular development in sheep is crucial for enhancing semen quality and increasing lamb production.

MicroRNA (miRNA) is a non-coding RNA between 18 and 24 nucleotides in length, which controls post-transcriptional gene silencing by interacting with the 3′-untranslated region (3′-UTR) of target mRNA, stimulating mRNA degradation or blocking translation (6). In the past few decades, many studies have investigated miRNA expression profiles using microarray technology, small RNA sequencing and reported that many miRNAs are expressed in mouse and human germ cells. In mammals, miRNAs are key controllers of cell differentiation and function, mediate a variety of cellular processes, and play an indispensable role in testicular development and spermatogenesis (7). Substantial evidence indicates that microRNAs (miRNAs) play pivotal roles in various facets of reproductive physiology, encompassing testicular development, spermatogenesis, and the regulation of functional effects (8–10). In the study of miRNA in Tibetan sheep, 1,118 immune-related miRNAs were found and speculated to be related to testicular development and meiosis in Tibetan sheep. Similarly, a recent study showed that miR-301b-3p and miR-3584-5p could promote the proliferation of rat immature testicular Sertoli cells by targeted inhibition of RASD1 gene expression (11). Yu et al. (12) found in mice that miR-34c could promote the differentiation and meiosis process of Spermatogonial stem cells (SSCs) by targeting Nanos2 gene and up-regulating the expression of meiosis-related genes Stra8 and Dazl. Smorag et al. (13) found that miR-34b-5p could regulate the meiotic process of mouse spermatocytes by targeting IGFBP2 gene. Overexpression of miR-10a in human and mouse testicular germ cells can target and inhibit the expression of Rad51 gene, leading to meiotic arrest and complete male infertility (14). Zhang et al. (15) found that miR-34b could target MAP2K1 (also known as MEK1) to induce apoptosis of bovine testicular Sertoli cells through MEK/ERK signaling pathway. However, miRNAs involved in testicular development and their underlying molecular mechanisms have not been identified in Southdown × Hu F1 hybrids.

In this study, we used deep sequencing technology to characterize and study the differential expression of miRNAs in the testis of sheep at different stages of development, understand the molecular regulation mechanism, and identify the key miRNA targets involved in sheep testicular development and spermatogenesis. This study will help to further understand the function of miRNA in sheep testicular development and to identify key miRNAs that enhance reproductive performance in male sheep in the future.

Materials and methods

Animals and sample preparation

Twelve healthy Southdown × Hu F1 sheep were castrated in Qinghuan Mutton Sheep Breeding Company in Gansu Province. The age of the sheep was obtained from sheep breeding records. There were 3 sheep at 0 months old (newborn, namely M0–1, M0–2, M0–3), 3 sheep at 3 months old (sexually immature, namely M3–1, M3–2, M3–3), 3 sheep at 6 months old (sexually mature, namely M6–1, M6–2, M6–3), and 3 sheep at 1 year old (adult, namely Y1–1, Y1–2, Y1–3). We removed testes from 12 sheep after anesthesia and then stored them in an RNA/DNA sample protector (Servicebio, Wuhan, China). The testis of each sheep was dissected longitudinally, and the right testicular tissue of each sheep was collected, part of which was immediately frozen in liquid nitrogen and stored at −80°C for total RNA and protein extraction. The remaining portion will be fixed in a 2.5% glutaraldehyde solution and subsequently processed for paraffin embedding and sectioning.

Small RNA library construction and sequencing

After total RNA was extracted by Trizol reagent kit (Invitrogen, Carlsbad, CA, United States), the RNA molecules in a size range of 18–30 nt were enriched by polyacrylamide gelelectrophoresis (PAGE). Then the 3′ adapters were added and the 36–48 nt RNAs were enriched. The 5′ adapters were then ligated to the RNAs as well. The ligation products were reverse transcribed by PCR amplification and the 140–160 bp size PCR products were enriched to generate a cDNA library and sequenced using Illumina HiSeq Xten by Gene DenovoBiotechnology Co. (Guangzhou, China)

Bioinformatics analysis

Basic reads are converted to sequence raw data by base calling. Low-quality reads were filtered to remove reads containing 5′ primer contaminants and poly (A). Reads without 3′ adapters and insert tags, as well as Reads shorter than 15 nt or longer than 41 nt in the raw data, were filtered out to obtain clean Reads. For preliminary analysis, the length distribution of clean sequences in the reference genome was determined. Non-coding RNAs are labeled as RNAs, tRNAs, small nuclear RNAs (snRNAs), miRNA, gene, rep, miRNA and unannotation. These RNAs were aligned and subsequently Bowtie (16) searches were performed against Rfam v.10.1¹ (17). Known miRNAs were identified by aligning sequencing reads to the miRBase v22 database² and the sheep reference genome (Qar Rambouillet v1.0). The known miRNA expression patterns in different samples were analyzed (18). Subsequently, the unannotated reads were analyzed by mirdeep2 to predict novel miRNAs (19). The hairpin structure of pre-miRNA and the miRBase database were used to identify the corresponding miRNA star sequence and miRNA mature sequence. The expression of known and novel miRNAs was analyzed using TPM (transcripts per million; miRNAs normalized by TPM) (20). Differential expression analysis of miRNAs between the two groups was performed using the DEG algorithm in the R package (21). The P- values were adjusted by the method of Benjamini and Hoch-berg to control the false discovery rate (FDR). Differentially expressed miRNAs were defined when the adjusted FDR was <0.05. In addition, we used the mean TPM value to calculate fold change (FC) between groups, defining up-regulated (log2FC ≥ 1) and down-regulated (log2FC ≤ 1) miRNAs, respectively. The target of DE miRNA was predicted using Miranda software with the following parameters: single residue pair matching score ≥ 150, ΔG ≤ −30 kcal/mol, and strict 5′ seed pairing required (22). Considering the hypergeometric distribution, R was used for Gene ontology (GO) enrichment and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of DE miRNA-target-get genes (23–25).

Validation of RNA sequence data

Twelve known sheep miRNAs were randomly selected for analysis, and qRT-PCR was used to validate the RNA-seq data. Real-time PCR was performed using a LightCycler 480 II real-time PCR apparatus (Roche, Swiss). Each sample was analyzed in triplicate. microRNA specific primer sequences were designed and synthesized in the laboratory by Qingdao Biotechnology Company based on microRNA sequences obtained from the miRBase database (Release 20.0). Primer information is given in Supplementary Table S1.

Results

Overview of sequencing results

We used Illumina Hiseq 2500 sequencing to analyze small RNA populations in 12 libraries obtained from testis of 0-month-old (M0), 3-month-old (M3), 6-month-old (M6), and 1-year-old (Y1) sheep. The results of miRNA sequence quality control are presented in Supplementary Table S2. A total of 157.27 Mb raw reads were obtained. After removing low-quality sequences, aptamers and discarding sequences shorter than 18 nt, 35.58 Mb, 32.67 Mb, 34.01 Mb and 41.01 Mb clean reads were obtained from M0, M3, M6 and Y1 libraries, respectively, for further analysis. 922,960, 1,780,035, 3,372,847, and 4,070,428 unique srRNAs were extracted from M0, M3, M6, and Y1 testes and mapped to the sheep reference genome (Supplementary Table S3). All clean reads were aligned to the miRBase database and recorded as one of the known RNA categories based on their biogenesis and annotation (Figure 1A). As shown in Figure 1A, Y1(7.21%) accounted for the largest proportion of known miRNAs, while M0(6.56%) and M3(5.85%) accounted for the second, and M6(3.3%) accounted for the smallest. However, the highest proportion of unannotated small RNAs was found in the M6(96.13%) and Y1(91.89%) libraries, which represent other types of small RNAs such as piRNAs. A principal component analysis (PCA) of all mapped genes showed that M0, M3, M6, and Y1 group could be distinguished by part along the axis of the first principal component (Figure 1B). Analysis of the read length distribution of all small RNA libraries showed that the dominant length of small RNAs was 22 nt, accounting for at least 37.72% (Figure 1C).

Figure 1

Figure 1. Identification and classification of miRNA in sheep testis. (A) Total number of unique sequences in libraries M0, M3, M6, and Y1. (B) PCA analysis of all mapped genes in M0, M3, M6, and Y1 groups. (C) Length distribution and abundance of sequences in the M0, M3, M6, and Y1 libraries.

Identification of known and novel miRNAs in sheep testes

To identify miRNAs in sheep testes, the sequences obtained after removal of small RNAs, including tRNA, snRNA, and rRNA, were compared with the miRBase database to obtain known miRNAs. Sequences that were not annotated into the miRBase database were predicted as novel miRNAs. A total of 787 known miRNAs and 415 novel miRNAs were identified across the 12 libraries (Supplementary Table S4). The frequency of expression of each miRNA in the 12 libraries varied widely, ranging from a few to hundreds of thousands of sequence reads. The expression of these novel miRNAs was very low, ranging from 18 to 25 nt in length with a distribution peak of 22 nt. Overall, the expression of most known miRNAs was highest in the testis of 3-month-old sheep, followed by 6-month-old and finally 0-month-old sheep (Figure 2A). The most highly expressed members in all libraries were members of the mir-199-x, mir-125-x, miR-99-z, miR-151-x, and miR-202-x families, each with more than 100,000 reads.

Figure 2

Figure 2. Comparative analysis of testicular miRNAs in sheep at different developmental stages. (A) miRNA expression levels in testis of 0, 3, 6 months and 1 year old sheep. (B) Barplot showing the number of up-and down-regulated DE miRNAs. (C) Scatter plot for differential comparison of miRNAs between pairs of groups; red, green, and blue dots indicate significantly up-regulated, down-regulated, and unchanged transcripts, respectively.

Differential expression of miRNAs in sheep testis tissue

Using the differential expression genes (DEGs) algorithm from an R package, the transcriptional changes of microRNAs (miRNAs) during sheep testis development were analyzed. Comparisons of miRNA expression profiles between different libraries were presented in Supplementary Tables S5–S10, visualized through scatter plots. Additionally, a scatter plot illustrated the pattern of DE miRNAs among all control groups, revealing that only 138 DE miRNAs were identified between M6 and Y1 testes, albeit with low expression levels (Figure 2B). Between M0 and M3, a total of 217 DE miRNAs were identified, with 139 upregulated and 78 downregulated. Between M0 and M6, 254 DE miRNAs were found, comprising 205 upregulated and 49 downregulated. From M0 to Y1, a remarkable change of 405 DE miRNAs was observed, with 319 upregulated and 86 downregulated. During the transition from M3 to M6, 130 DE miRNAs were identified, of which 97 were upregulated and 33 were downregulated. Lastly, between M3 and Y1, 305 DE miRNAs were detected, with 236 upregulated and 69 downregulated (Figure 2C).

Prediction of target genes for DE miRNAs

Using the Miranda software, the target genes of miRNAs were predicted. Subsequently, the potential functional roles of these differentially expressed miRNAs were investigated. Potential targets for the differentially expressed miRNAs were screened based on criteria of a minimum single-residue pair total match score of ≥150 and a total energy of ≤ − 30 kcal/mol. Among the 217, 254, 405, 130, 305, and 138 miRNAs identified between M0 vs. M3, M0 vs. M6, M0 vs. Y1, M3 vs. M6, M3 vs. Y1, and M6 vs. Y1, respectively, 5,233, 8,769, 23,588, 696, 8,021, and 180 targets were predicted in 885, 1,679, 2,679, 388, 2006, and 171 target genes, respectively. For the majority of differentially expressed miRNAs, multiple distinct target genes were present; however, for some differentially expressed miRNAs, only a single target gene was identified. Additionally, certain target genes were targeted by multiple differentially expressed miRNAs. For instance, oar-miR-431 was predicted to target 914 genes, miR-199-y was predicted to target 1,510 genes, while BMP4 was uniquely targeted by miR-142-y.

GO enrichment and KEGG pathway analysis of target genes

To gain a deeper understanding of the functions of DE miRNAs in sheep testis development, we performed Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses on the candidate target genes of all DE miRNAs. The GO analysis results are presented in Figure 3. Between M0 and M3, target genes were significantly enriched in categories such as animal organ development, cell morphogenesis, cell part morphogenesis, and Rab GTPase binding (Figure 3A). Between M0 and M6, target genes were significantly enriched in regulation of cellular biosynthetic process, dendrite morphogenesis, tube development, organ morphogenesis, and cell–cell adherens junction (Figure 3B). Between M0 and Y1, target genes were notably enriched in reproductive process, apical junction complex, apical plasma membrane, and ATP binding (Figure 3D). Between M3 and M6, target genes were significantly enriched in cell differentiation, cell morphogenesis involved in differentiation, Wnt signaling pathway, and apical junction complex (Figure 3C). Between M3 and Y1, target genes were enriched in tube morphogenesis, regulation of organelle organization, plasma membrane, Golgi complex, and microtubule binding (Figure 3F). Lastly, between M6 and Y1, target genes were significantly enriched in cell migration, dendrite morphogenesis, organ morphogenesis, and cell–cell junction organization (Figure 3E).

Figure 3

Figure 3. GO enrichment analysis of target genes of differentially expressed miRNAs. (A) M0 vs M3; (B) M0 vs M6; (C) M3 vs M6; (D) M0 vs Y1; (E) M6 vs Y1; (F) M3 vs Y1; M0 represents 0 months of age, M3 represents 3 months of age, M6 represents 6 months of age, and Y1 represents 1 year of age.

KEGG Pathway Annotation identified significant enrichment of target genes of DE miRNAs in various comparisons. In M0 vs. M3, 30 pathways were enriched, including the Ras, MAPK, Rap1, cAMP, and Notch signaling pathways (Figure 4A). Four key pathways related to testis development were cAMP, MAPK, PI3K-Akt, and Wnt signaling. In M0 vs. M6, target genes were enriched in 30 pathways, such as adherens junction, tight junction, spinocerebellar ataxia, phospholipase D signaling, and endocytosis (Figure 4B). Four testis development-related pathways were MAPK, cAMP, FOXO, and GnRH signaling pathways. Between M3 and M6, 30 pathways were enriched, including Notch, Ras, thyroid hormone, and ErbB signaling (Figure 4C). Five key pathways were Notch, ErbB, endocytosis, Ras, and thyroid hormone signaling.

Figure 4

Figure 4. KEGG analysis of target genes of differentially expressed miRNAs. (A) M0 vs M3; (B) M0 vs M6; (C) M3 vs M6; (D) M0 vs Y1; (E) M6 vs Y1; (F) M3 vs Y1; M0 represents 0 months of age, M3 represents 3 months of age, M6 represents 6 months of age, and Y1 represents 1 year of age.

In M0 vs. Y1, 49 pathways were enriched, including Notch, Ras, thyroid hormone, and ErbB signaling (Figure 4D). Five testis development-related pathways were FOXO, ErbB, endocytosis, Ras, and thyroid hormone signaling. Between M6 and Y1, 30 pathways were enriched, such as Hippo, adherens junction, cAMP, Notch, and Ras signaling (Figure 4E). Five key pathways were Notch, Hippo, cAMP, endocytosis, and FOXO signaling. Similarly, in M3 vs. Y1, 30 pathways were enriched, including tight junction, endocytosis, Ras, adherens junction, and spinocerebellar ataxia (Figure 4F). Four testis development-related pathways were Ras, MAPK, and cAMP signaling (Figure 5).

Figure 5

Figure 5. Cluster analysis of differential miRNA co-expression in pairwise comparison groups. (A) M0 vs M3; (B) M0 vs M6; (C) M0 vs Y1; (D) M3 vs M6; (E) M3 vs Y1; (F) M6 vs Y1; M0 represents 0 months of age, M3 represents 3 months of age, M6 represents 6 months of age, and Y1 represents 1 year of age.

Cluster analysis of miRNA co-expression

We performed co-expression cluster analysis of the differentially expressed miRNAs screened from each comparison group and screened the top 5 miRNAs with degree values. After comparison, the number of miRNA target genes in M0 vs. M3 group, M0 vs. Y1 group and M0 vs. M6 group was higher than that in other groups. Notably, oar-miR-29b, oar-miR-143, oar-miR-370-3p, oar-miR-19b, oar-miR-432, and oar-miR-29b were repeated for different times in the 6 comparison groups. TGF family genes related to testicular development were all regulated by differential miRNAs in the five comparison groups (Figure 5).

Quantitative RT-PCR validation

Randomly selected 12 known sheep miRNAs, including oar-miR-152, oar-miR-10b, oar-miR-21, oar-miR-22-3p, oar-miR-25, oar-miR-99a, oar-miR-379-3p, oar-miR-409-5p, oar-miR-1179, oar-miR-132, oar-miR-374a, and oar-miR-28-y, were used to validate the RNA-seq data. As shown in Figure 6, the qRT-PCR results exhibited similar expression trends compared to the small RNA-seq results.

Figure 6

Figure 6. Validation of the expression of differentially expressed miRNAs by qRT-PCR. The blue bar graph represents RT-qPCR, the red broken line represents RNA-seq, the left ordinate is PCR Expression, and the right ordinate is RNA-seq.

Discussion

The production of sperm and androgens is the primary function of the testis, both of which play a vital role in the entire reproductive system. The differentiation of somatic cells into Sertoli cells initiates male-specific development and guides the development of germ cells toward the spermatogenic lineage (26). Consequently, testicular development and spermatogenesis are crucial processes that affect the reproductive efficiency of male sheep. Additionally, several functional genes have been reported to regulate testicular development and spermatogenesis in sheep (27–29). Recently, miRNAs have been found to regulate biological processes by inhibiting the translation or degradation of target mRNAs (30–32). Since the discovery of miRNAs and their role in gene expression regulation, numerous studies have delved into their involvement in testicular hormone synthesis, cell proliferation, meiosis of haploid spermatocytes, and spermatogenesis (33–35). However, only a few studies have characterized the miRNAs involved in testicular development during early puberty in sheep. In this study, we employed Illumina Solexa technology to sequence small RNAs in testicular tissues from 0-month, 3-month, 6-month, and 1-year-old male sheep, analyzed DE miRNAs, predicted novel miRNAs, and conducted GO enrichment and KEGG pathway analyses of target genes across four miRNA libraries.

The dynamic testicular development of male sheep at 0-month, 3-month, 6-month, and 1-year-old corresponds to the neonatal, pre-pubertal, pubertal, and adult stages, respectively. By comparing the reference genome, we identified 217, 254, 405, 130, 205, and 138 known DE miRNAs in comparisons between M0 vs. M3, M0 vs. M6, M0 vs. Y1, M3 vs. M6, M3 vs. Y1, and M6 vs. Y1, respectively. Notably, more novel DE miRNAs were identified compared to known DE sheep miRNAs, contributing to the identification of key miRNAs in sheep testicular development and enriching the sheep miRBase database. Previous studies have reported that miR-34c, miR-21, and miR-499b regulate the reproductive performance of male mammals, including promoting spermatogonial meiosis, inhibiting Sertoli cell apoptosis, and maintaining the necessary germ cell population of spermatogonial stem cells (SSCs) (36–38). In our study, miR-34c, miR-21, and miR-499b showed significantly higher expression in the Y1 group compared to the M6 group, suggesting their potential involvement in sheep testicular development. Additionally, miR-31-5p downregulation has been detected in the miRNA profiles of seminal exosomes from patients with asthenozoospermia and azoospermia (39) and in non-obstructive azoospermia testes (40), indicating its potential as a predictive marker for azoospermia. Accordingly, miR-155 targets genes involved in MAPK and Wnt signaling pathways, affecting sperm motility (41). Overexpression of miR-10a in testicular germ cells of humans and mice can target and inhibit Rad51 gene expression, leading to meiotic arrest and complete male infertility (42). Similarly, these miRNAs may be associated with sheep testicular development, albeit their specific functions in sheep testes require further exploration and validation.

GO annotation and KEGG pathway analysis can provide a detailed and comprehensive understanding of the functions of DE miRNA target genes from the aspects of cellular components, biological processes, molecular functions and related pathways. GO and KEGG analyses showed that the target genes of DE miRNAs were involved in different biological processes, including growth and development, reproduction, cell proliferation, hormone secretion, and several other metabolic processes. Interestingly, the cAMP, Hippo, MAPK, Wnt, and FOXO signaling pathways are widely recognized to regulate reproduction.

FOXO1 is mainly expressed in Sertoli cells (43), while FOXO3a is expressed in other cells of the testis, such as Leydig cells (44). This cell type specificity allows FOXO family members to regulate the development and function of different cell types in the testis. FOXO further regulates testicular development by promoting gonadotropin secretion and testosterone synthesis. The most multifunctional miRNAs controlling target genes in FOXO signaling pathways included oar-miR-19b, oar-miR-200a, oar-miR-369-3p, oar-miR-200b, and oar-miR-493-5p. Each of them regulates key genes in the FOXO signaling pathway such as FOXO1 and FOXO6. The PI3K-Akt signaling pathway directly or indirectly maintains and promotes spermatogenesis by regulating the proliferation, survival, and anti-apoptosis of immature Sertoli cells, mature Sertoli cells, spermatogonial stem cells, and spermatogenic cells (45). In PI3K-Akt signaling pathway, oar-miR-1197-3p, oar-miR-136, oar-miR-199a-3p, oar-miR-410-5p, oar-miR-432 targeted WT1 gene regulation, WT1 can not only activate PI3K/Akt and other pathways to inhibit apoptosis, but also plays a crucial role in the assembly and maintenance of fetal testicular cords (46). The Wnt pathway affects the production and maturation of spermatogenic cells by regulating cell proliferation and differentiation. Studies have shown that abnormal activation or inhibition of Wnt may affect the normal proliferation and differentiation process of testicular cells, which in turn affects sperm production (47). As a key gene in this pathway, DOT1L is co-targeted by DE miRNAs such as oar-let-7f and oar-miR-22-3p. DOT1L promotes SSCs self-renewal by catalyzing Histone H3 lysine 79 methylation (H3K79me). It is essential for the maintenance of spermatogenesis (48). The expression of oar-miR-22-3p was found to be significantly changed from M0 to Y1, suggesting that oar-miR-22-3p plays an important role in spermatogenesis and testicular development in sheep. Evidence accumulated over the past few years has highlighted the major functions of MAPK, AMPK, and cAMP signaling pathways during spermatogenesis, with MAPK signaling regulating germ cell proliferation, meiosis, and Sertoli cell proliferation (49). In addition, the AMPK and FOXO signaling pathways influence the proliferation of Sertoli cells (50). Among the MAPK, cAMP, PI3K-Akt, FOXO, Wnt, and AMPK pathways, oar-miR-133, oar-miR-1197-3p, oar-miR-758-3pr, and oar-let-7f are the most potent miRNAs that regulate multiple genes within various testicular development and spermatogenesis-related pathways. Each of these miRNAs plays a crucial role in modulating diverse targets in these signaling cascades. In this study, the downregulation of oar-miR-26a, oar-miR-99a, and oar-miR-148a during the development of sheep testes from neonatal to adult stages suggests their significant contribution to promoting testicular development and spermatogenesis in sheep.

Conclusion

In summary, we have characterized the miRNA transcriptome across four testicular developmental stages, from birth to adulthood, in the F1 offspring of a cross between Southdown sheep and Hu sheep. A total of 788 known miRNAs and 415 novel miRNAs were identified across the 12 libraries. GO and KEGG pathway analyses of all DE miRNA targets in the libraries revealed significant enrichment in MAPK signaling pathway, Rap1 signaling pathway, Wnt signaling pathway, and FOXO signaling pathway. We believe these results will contribute to a deeper understanding of the roles of miRNAs in sheep testicular development and potentially facilitate the discovery of miRNAs to enhance the reproductive performance of male sheep in the future.

Data availability statement

The data presented in the study are deposited in the NCBI SRA repository, accession number PRJNA1195695.

Ethics statement

The animal studies were approved by Animal Management and Ethics Committee of Lanzhou Institute of Animal Husbandry and Medicine, Chinese Academy of Agricultural Sciences and met the requirements of the Institutional Animal Care and Use Committee (license number: SYXK-2019-010). The studies were conducted in accordance with the local legislation and institutional requirements. Written informed consent was obtained from the owners for the participation of their animals in this study.

Author contributions

BX: Formal analysis, Methodology, Validation, Visualization, Writing – original draft, Writing – review & editing. XA: Formal analysis, Investigation, Methodology, Validation, Visualization, Writing – original draft, Writing – review & editing. YYu: Validation, Visualization, Writing – original draft. HS: Methodology, Writing – original draft. GH: Investigation, Software, Writing – original draft. YYa: Conceptualization, Data curation, Funding acquisition, Project administration, Resources, Supervision, Writing – review & editing. SZ: Conceptualization, Data curation, Funding acquisition, Project administration, Resources, Supervision, Writing – review & editing.

Funding

The author(s) declare that financial support was received for the research and/or publication of this article. This research was funded by the Gansu Province Science and Technology Plan Project Technology Innovation Guidance Program-LuGan Science and Technology Collaboration topic of 23CXNA0020.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Generative AI statement

The authors declare that no Gen AI was used in the creation of this manuscript.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fvets.2025.1538990/full#supplementary-material

SUPPLEMENTARY TABLE S1 | Details of primer sequences of miRNAs used for RT-qPCR.

SUPPLEMENTARY TABLE S2 | miRNA Sequence Quality Control.

SUPPLEMENTARY TABLE S3 | Read mapping summary of sheep miRNA.

SUPPLEMENTARY TABLE S4 | The novel miRNAs and known miRNAs identified in all libraries.

SUPPLEMENTARY TABLE S5 | Differentially expressed miRNAs between M0 and M3.

SUPPLEMENTARY TABLE S6 | Differentially expressed miRNAs between M0 and M6.

SUPPLEMENTARY TABLE S7 | Differentially expressed miRNAs between M0 and Y1.

SUPPLEMENTARY TABLE S8 | Differentially expressed miRNAs between M3 and M6.

SUPPLEMENTARY TABLE S9 | Differentially expressed miRNAs between M3 and Y1.

SUPPLEMENTARY TABLE S10 | Differentially expressed miRNAs between M6 and Y1.

Footnotes

1. ^http://www.sanger.a.uk/softare/Rfam

2. ^http://www.mirbase.org/

References

1. Wang, J. Analysis of factors affecting farly growth traits of Hu sheep and estimation of genetic parameters (D). Xinjiang Agric Univ. (2022). doi: 10.27431/d.cnki.gxnyu.2022.000762

Crossref Full Text | Google Scholar

2. Liu, T, Lu, Z, Shi, H, Sun, W, Zhang, J, and Zhu, Z. Comparative study on mutton quality of F1 generation of different cross combinations. Acta Ecol Anim Domast. (2023) 44:46–52. doi: 10.3969/j.issn.1673-1182.2023.09.009

Crossref Full Text | Google Scholar

3. Zhang, R, An, X, Li, J, Lu, Z, Niu, C, and Xu, Z. Comparative analysis of growth performance, meat productivity, and meat quality in Hu sheep and its hybrids. Acta Pratacul Sin. (2024) 33:186–97. doi: 10.11686/cyxb2023157

Crossref Full Text | Google Scholar

4. Kong, L, Yue, Y, Zhang, C, and Lu, Z. Slaughtering performance and meat quality characteristics of Hu sheep and its hybrid offspring with Southdown sheep. Food Sci. (2023) 44:64–70. doi: 10.7506/spkx1002-6630-20221124-285

Crossref Full Text | Google Scholar

5. Marcon, E, Babak, T, Chua, G, Hughes, T, and Moens, PB. miRNA and piRNA localization in the male mammalian meiotic nucleus. Chromosom Res. (2008) 16:243–60. doi: 10.1007/s10577-007-1190-6

Crossref Full Text | Google Scholar

6. Afonso, J, Lima, AO, de Sousa, MAP, de Athayde, FRF, and Fortes, MRS. Transcription factors and miRNA act as contrary regulators of gene expression in the testis and epididymis of Bos indicus animals. Gene. (2024) 899:148133. doi: 10.1016/j.gene.2024.148133

PubMed Abstract | Crossref Full Text | Google Scholar

7. Hossain, MM, Sohel, MM, Schellander, K, and Tesfaye, D. Characterization and importance of microRNAs in mammalian gonadal functions. Cell Tissue Res. (2012) 349:679–90. doi: 10.1007/s00441-012-1469-6

Crossref Full Text | Google Scholar

8. McIver, SC, Roman, SD, Nixon, B, and McLaughlin, EA. miRNA and mammalian male germ cells. Hum Reprod Update. (2012) 18:44–59. doi: 10.1093/humupd/dmr041

Crossref Full Text | Google Scholar

9. Kotaja, N. MicroRNAs and spermatogenesis. Fertil Steril. (2014) 101:1552–62. doi: 10.1016/j.fertnstert.2014.04.025

PubMed Abstract | Crossref Full Text | Google Scholar

10. Chen, W, Cui, Y, Ning, M, Zhang, H, Yin, C, and He, Z. The mechanisms and functions of microRNAs in mediating the fate determinations of human spermatogonial stem cells and Sertoli cells. Semin Cell Dev Biol. (2022) 121:32–9. doi: 10.1016/j.semcdb.2021.05.003

PubMed Abstract | Crossref Full Text | Google Scholar

11. Yin, X, Ma, T, Han, R, Ding, J, Zhang, H, Han, X, et al. MiR-301b-3p/3584-5p enhances low-dose mono-n-butyl phthalate (MBP)-induced proliferation by targeting Rasd1 in Sertoli cells. Toxicol In Vitro. (2018) 47:79–88. doi: 10.1016/j.tiv.2017.11.009

Crossref Full Text | Google Scholar

12. Yu, M, Mu, H, Niu, Z, Zhu, Z, Zhu, H, and Jin, L. miR-34c enhances mouse spermatogonial stem cells differentiation by targeting Nanos2. J Cell Biochem. (2014) 115:232–42. doi: 10.1002/jcb.24655

PubMed Abstract | Crossref Full Text | Google Scholar

13. Smorag, L, Zheng, Y, Nolte, J, Zechner, U, Engel, W, and Pantakani, D. MicroRNA signature in various cell types of mouse spermatogenesis: evidence for stage-specifically expressed miRNA-221, −203 and -34b-5p mediated spermatogenesis regulation. Biol Cell. (2012) 104:677–92. doi: 10.1111/boc.201200014

PubMed Abstract | Crossref Full Text | Google Scholar

14. Wang, Z, Zhang, X, Li, Z, Abdalla, BA, Chen, Y, and Nie, Q. MiR-34b-5p mediates the proliferation and differentiation of myoblasts by targeting IGFBP2. Cells. (2019) 8:360. doi: 10.3390/cells8040360

Crossref Full Text | Google Scholar

15. Zhang, L, Ma, T, Tao, Q, Tan, W, Chen, H, Liu, W, et al. Bta-miR-34b inhibits proliferation and promotes apoptosis via the MEK/ERK pathway by targeting MAP2K1 in bovine primary Sertoli cells. J Anim Sci. (2020) 98:skaa313. doi: 10.1093/jas/skaa313

PubMed Abstract | Crossref Full Text | Google Scholar

16. Langmead, B, Trapnell, C, Pop, M, and Salzberg, SL. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. (2009) 10:R25. doi: 10.1186/gb-2009-10-3-r25

Crossref Full Text | Google Scholar

17. Griffiths-Jones, S, Bateman, A, Marshall, M, Khanna, A, and Eddy, SR. Rfam: an RNA family database. Nucleic Acids Res. (2003) 31:439–41. doi: 10.1093/nar/gkg006

PubMed Abstract | Crossref Full Text | Google Scholar

18. Griffiths-Jones, S, and Saini, HK. Van Dongen, S Enright AJ. miRBase: tools for microRNA genomics. Nucleic Acids Res. (2008) 36:D154–8. doi: 10.1093/nar/gkm952

PubMed Abstract | Crossref Full Text | Google Scholar

19. Friedlander, MR, Mackowiak, SD, Li, N, Chen, W, and Rajewsky, N. miRDeep2 accurately identifies known and hundreds of novel microRNA genes in seven animal clades. Nucleic Acids Res. (2012) 40:37–52. doi: 10.1093/nar/gkr688

PubMed Abstract | Crossref Full Text | Google Scholar

20. Sun, J, Wang, S, Li, C, Ren, Y, and Wang, J. Novel expression profiles of microRNAs suggest that specific miRNAs regulate gene expression for the sexual maturation of female Schistosoma japonicum after pairing. Parasite Vector. (2014) 7:177. doi: 10.1186/1756-3305-7-177

PubMed Abstract | Crossref Full Text | Google Scholar

21. Anders, S, and Huber, W. Differential expression analysis for sequence count data. Genome Biol. (2010) 11:R106. doi: 10.1186/gb-2010-11-10-r106

PubMed Abstract | Crossref Full Text | Google Scholar

22. Enright, AJ, John, B, Gaul, U, Tuschl, T, Sander, C, and Marks, DS. MicroRNA targets in Drosophila. Genome Biol. (2003) 5:R1. doi: 10.1186/gb-2003-5-1-r1

PubMed Abstract | Crossref Full Text | Google Scholar

23. Kanehisa, M, and Goto, S. KEGG: Kyoto encyclopedia of genes and genomes. Nucleic Acids Res. (2000) 27:29–34. doi: 10.1093/nar/27.1.29

PubMed Abstract | Crossref Full Text | Google Scholar

24. Kanehisa, M. Toward understanding the origin and evolution of cellular organisms. Protein Sci. (2019) 28:1947–51. doi: 10.1002/pro.3715

PubMed Abstract | Crossref Full Text | Google Scholar

25. Kanehisa, M, Furumichi, M, Sato, Y, Kawashima, M, and Ishiguro-Watanabe, M. KEGG for taxonomy-based analysis of pathways and genomes. Nucleic Acids Res. (2022) 51:D587–92. doi: 10.1093/nar/gkac963

PubMed Abstract | Crossref Full Text | Google Scholar

26. Mäkelä, JA, Koskenniemi, JJ, Virtanen, HE, and Toppari, J. Testis Development. Endocr Rev. (2019) 40:857–905. doi: 10.1210/er.2018-00140

PubMed Abstract | Crossref Full Text | Google Scholar

27. Li, S, Zhao, B, Yang, H, Dai, K, Cai, Y, Xu, H, et al. Comprehensive transcriptomic analysis revealing the regulatory dynamics and networks of the pituitary-testis axis in sheep across developmental stages. Front Vet Sci. (2024) 11:1367730. doi: 10.3389/fvets.2024.1367730

PubMed Abstract | Crossref Full Text | Google Scholar

28. Yang, H, Wang, F, Li, F, Ren, C, Pang, J, Wan, Y, et al. Comprehensive analysis of long noncoding RNA and mRNA expression patterns in sheep testicular maturation. Biol Reprod. (2018) 99:650–61. doi: 10.1093/biolre/ioy088

PubMed Abstract | Crossref Full Text | Google Scholar

29. Xu, H, Sun, W, Pei, S, Li, W, Li, F, and Yue, X. Identification of key genes related to postnatal testicular development based on transcriptomic data of testis in Hu sheep. Front Genet. (2022) 12:773695. doi: 10.3389/fgene.2021.773695

PubMed Abstract | Crossref Full Text | Google Scholar

30. Yun, Y, Wu, R, He, X, Qin, X, Chen, L, Sha, L, et al. Integrated transcriptome analysis of miRNAs and mRNAs in the skeletal muscle of Wuranke sheep. Genes (Basel). (2023) 14:2034. doi: 10.3390/genes14112034

PubMed Abstract | Crossref Full Text | Google Scholar

31. Wang, Q, Cao, H, Su, X, and Liu, W. Identification of key miRNAs regulating fat metabolism based on RNA-seq from fat-tailed sheep and F2 of wild argali. Gene. (2022) 834:146660. doi: 10.1016/j.gene.2022.146660

PubMed Abstract | Crossref Full Text | Google Scholar

32. Han, H, Chen, Q, Gao, Y, Li, J, Li, W, Dang, R, et al. Comparative transcriptomics analysis of testicular miRNA from cryptorchid and normal horses. Animals. (2020) 10:338. doi: 10.3390/ani10020338

PubMed Abstract | Crossref Full Text | Google Scholar

33. Pan, C, Fei, Q, Jin, J, Zheng, J, Wu, D, Li, H, et al. MicroRNA profiling reveals the role of miR-133b-3p in promoting apoptosis and inhibiting cell proliferation and testosterone synthesis in mouse TM3 cells. In Vitro Cell Dev Biol Anim. (2023) 59:63–75. doi: 10.1007/s11626-022-00745-z

PubMed Abstract | Crossref Full Text | Google Scholar

34. Yu, Z, Liu, Y, Li, Y, Zhang, J, Peng, J, Gong, J, et al. miRNA-338-3p inhibits glioma cell proliferation and progression by targeting MYT1L. Brain Res Bull. (2022) 179:1–12. doi: 10.1016/j.brainresbull.2021.11.016

PubMed Abstract | Crossref Full Text | Google Scholar

35. Walker, WH. Regulation of mammalian spermatogenesis by miRNAs. Semin Cell Dev Biol. (2022) 121:24–31. doi: 10.1016/j.semcdb.2021.05.009

PubMed Abstract | Crossref Full Text | Google Scholar

36. Niu, Z, Goodyear, SM, Rao, S, Wu, X, Tobias, JW, Avarbock, MR, et al. MicroRNA-21 regulates the self-renewal of mouse spermatogonial stem cells. Proc Natl Acad Sci USA. (2011) 108:12740–5. doi: 10.1073/pnas.1109987108

PubMed Abstract | Crossref Full Text | Google Scholar

37. Bouhallier, F, Allioli, N, Lavial, F, Chalmel, F, Perrard, MH, Durand, P, et al. Role of miR-34c microRNA in the late steps of spermatogenesis. RNA. (2010) 16:720–31. doi: 10.1261/rna.1963810

PubMed Abstract | Crossref Full Text | Google Scholar

38. Liu, WM, Pang, RT, Chiu, PC, Wong, BP, Lao, K, Lee, KF, et al. Sperm-borne microRNA-34c is required for the first cleavage division in mouse. Proc Natl Acad Sci USA. (2012) 109:490–4. doi: 10.1073/pnas.1110368109

Crossref Full Text | Google Scholar

39. Barcelo, M, Mata, A, Bassas, L, and Larriba, S. Exosomal microRNAs in seminal plasma are markers of the origin of azoospermia and can predict the presence of sperm in testicular tissue. Hum Reprod. (2018) 33:1087–98. doi: 10.1093/humrep/dey072

PubMed Abstract | Crossref Full Text | Google Scholar

40. Zhang, HT, Zhang, Z, Hong, K, Tang, WH, Liu, DF, Mao, JM, et al. Altered microRNA profiles of testicular biopsies from patients with nonobstructive azoospermia. Asian J Androl. (2020) 22:100–5. doi: 10.4103/aja.aja_35_19

PubMed Abstract | Crossref Full Text | Google Scholar

41. Liu, Y, Sun, Y, Li, Y, Bai, H, Xu, S, Xu, H, et al. Identification and differential expression of microRNAs in the testis of chicken with high and low sperm motility. Theriogenology. (2018) 122:94–101. doi: 10.1016/j.theriogenology.2018.09.010

PubMed Abstract | Crossref Full Text | Google Scholar

42. Gao, H, Wen, H, Cao, C, Dong, D, Yang, C, Xie, S, et al. Overexpression of MicroRNA-10a in germ cells causes male infertility by targeting Rad51 in mouse and human. Front Physiol. (2019) 10:765. doi: 10.3389/fphys.2019.00765

PubMed Abstract | Crossref Full Text | Google Scholar

43. Javdan, N, Ayatollahi, SA, Iqbal Choudhary, M, Al-Hasani, S, and Pazoki-Toroudi, H. FOXO1 targeting by capsaicin reduces tissue damage after testicular torsion. Andrologia. (2018) 50:e13032. doi: 10.1111/and.13032

PubMed Abstract | Crossref Full Text | Google Scholar

44. Zhao, X, Liu, F, Chen, Z, Gao, Y, Cui, K, and Shi, D. Phenotypic study of foxo3a knockout mice created by TALEN-mediated gene targeting. China Acad J Electr Publ House. (2018) 31:23–4. doi: 10.1038/nbt.2477

PubMed Abstract | Crossref Full Text | Google Scholar

45. Chen, K, Wei, B, Hao, S, and Yang, W. The PI3K/AKT signaling pathway: how does it regulate development of Sertoli cells and spermatogenic cells? Histol Histopathol. (2022) 37:621–36. doi: 10.14670/HH-18-457

PubMed Abstract | Crossref Full Text | Google Scholar

46. Wen, Q, Wang, Y, Tang, J, Cheng, C, and Liu, Y. Sertoli cell Wt1 regulates peritubular Myoid cell and fetal Leydig cell differentiation during fetal testis development. PLoS One. (2016) 11:e0167920. doi: 10.1371/journal.pone.0167920

PubMed Abstract | Crossref Full Text | Google Scholar

47. Lundgaard Riis, M, Delpouve, G, Nielsen, JE, Melau, C, Langhoff Thuesen, L, Juul Hare, K, et al. Inhibition of WNT/β-catenin signalling during sex-specific gonadal differentiation is essential for normal human fetal testis development. Cell Commun Signal. (2024) 22:330. doi: 10.1186/s12964-024-01704-9

PubMed Abstract | Crossref Full Text | Google Scholar

48. Lin, H, Cossu, IG, Leu, NA, Deshpande, AJ, Bernt, KM, Luo, M, et al. The DOT1L-MLLT10 complex regulates male fertility and promotes histone removal during spermiogenesis. Development. (2023) 150:dev201501. doi: 10.1242/dev.201501

PubMed Abstract | Crossref Full Text | Google Scholar

49. La, Y, Ma, X, Bao, P, Chu, M, Guo, X, Liang, C, et al. Identification and profiling of microRNAs during yak's testicular development. BMC Vet Res. (2023) 19:53. doi: 10.1186/s12917-023-03602-7

PubMed Abstract | Crossref Full Text | Google Scholar

50. Meroni, SB, Galardo, MN, Rindone, G, Gorga, A, Riera, MF, and Cigorraga, SB. Molecular mechanisms and signaling pathways involved in Sertoli cell proliferation. Front Endocrinol. (2019) 10:224. doi: 10.3389/fendo.2019.00224

PubMed Abstract | Crossref Full Text | Google Scholar