Global Analysis in Nonobstructive Azoospermic Testis Identifies miRNAs Critical to Spermatogenesis

Introduction . The etiology of male infertility characterized by non-obstructive azoospermia is largely unknown, especially at the molecular level. Identifying dysregulated microRNAs (miRNAs) in male infertility would be useful to achieve a more profound understanding of its pathogenesis. Methods . Small RNA sequencing was performed on the testicular tissues of 10 nonobstructive azoospermic patients with the Sertoli cell only syndrome (SCOS) and 8 obstructive azoospermic individuals with normal spermatogenesis. The expressions of two dysregulated miRNAs were validated by quantitative real-time polymerase chain reaction, con ﬁ rming the results obtained by sequencing analysis. Bioinformatic analysis was undertaken to identify the main pathways impaired in complete spermatogenic failure. Results . A total of 136 miRNAs were detected to be di ﬀ erentially expressed in the Sertoli cell only syndrome group in comparison with the obstructive azoospermia group. Bioinformatic analysis suggested that the altered miRNAs were substantially involved in pathways related to spermatogenesis. Conclusions . Our study investigates the entire pro ﬁ le of miRNAs with emphasis on the crucial role of miRNAs in idiopathic Sertoli cell only syndrome, suggesting potential targets for employing molecular therapeutic strategies in the treatment of spermatogenic failure.


Introduction
Infertility affects nearly 15% of couples as a worldwide health problem, and half of the cases are attributable to male factors [1]. Around 30-40% of male infertility cases remain idiopathic, which also includes nearly 15% of azoospermic cases [2,3]. The etiology of azoospermia is largely unknown, especially at the molecular level. Spermatogenesis is an intricately organized and complex process, which is accomplished in three successive phases: consecutive mitotic divisions of spermatogonia, reductive meiotic division of spermatocytes, and differentiation of spermatids [4]. Several protein mediators are involved in the regulation of this complicated process; however, microRNAs (miRNAs) play critical roles at the posttranscriptional stage [5][6][7]. miRNAs are small (18)(19)(20)(21)(22) nucleotides long) noncoding RNAs, which can negatively and selectively modulate gene expression. miRNAs have been described as the molecular regulators of various cellular biological events such as growth, proliferation, differentiation, and programmed cell death [8][9][10]. Hence, alterations in the molecular signaling in the testis may result in impaired spermatogenesis and male infertility. Identification of the molecular players dysregulated in azoospermia/infertility would pave the way to the development of new therapeutic approaches in the future. Since miRNAs regulate a large number of genes, identification of dysregulated miRNAs in infertility would be useful in the diagnosis and treatment of infertility [11].
So far, a few studies have analyzed the miRNAs' expression levels in testicular biopsies obtained from azoospermic patients and normal individuals [12][13][14][15][16][17]. In the first study, Lian et al. performed microarray and identified altered expression of miRNAs in nonobstructive azoospermic (NOA) patients, with 19 and 154 up-and downregulated miRNAs, respectively [14]. Subsequently, another study analyzed 48 testicular biopsy samples obtained from patients with various degrees of spermatogenic impairment (germ cell arrest, mixed atrophy, and Sertoli cell only) and compared them with normozoospermic controls [12]. A number of dysregulated miRNAs were identified in the infertile group. Based on these results, Abu-Halima et al. selected a set of five miRNAs and analyzed them in a large set of infertile biopsy samples, emphasizing their value as diagnostic biomarkers of male infertility [13]. Munoz et al. screened 623 miRNAs in testicular samples of men with spermatogenic failure and identified a set of three miRNAs, including miR-122, miR-34c-5p, and miR-449a, with high predictive value for the presence of mature spermatozoa in the testis [15]. Noveski et al. analyzed testicular biopsies of infertile men with hypospermatogenesis, maturation arrest, hypospermatogenesis with azoospermia factor c (AZFc) deletion, and Sertoli cell only syndrome (SCOS), while men with normal spermatogenesis served as the control group [16]. They found some dysregulated miRNAs common among different histopathological groups. In addition, miRNA-mRNA network analysis revealed the role of several genes involved in the tumor necrosis factor-related apoptosis-inducing ligand (TRAIL) pathway. Recently, Zhang et al. used microarray to study the miRNA profile of testicular biopsies from NOA patients, finding that 129 miRNAs were aberrantly expressed in comparison to the controls [17]. In particular, this study reported that a panel composed of miR-10b-3p and miR-34b-5p may be used to predict the existence of testicular sperm in azoospermic patients.
Most of the above studies have laid the foundation for the identification of miRNA-based diagnosis of the pathogenesis that can unveil the path to therapy. Independent investigations in different populations would not only further strengthen these initial findings but also provide a way to finding consistent miRNAs that are critical to spermatogenesis and fertility. Therefore, in the current work, we performed small RNA sequencing in a selected set of patients with SCOS to uncover the global microRNA dysregulation in comparison to obstructive azoospermia (OA) individuals with normal spermatogenesis.

Material and Methods
2.1. Sample Collection. Azoospermic infertile men referred to the Royan Institute (Tehran, Iran) were recruited in this study after receiving their informed consent. Idiopathic NOA individuals (n = 10) characterized by SCOS constituted the case group while individuals with OA constituted the control group (n = 8). Testicular specimens of these individuals were obtained by testicular sperm extraction (TESE) procedure, performed as diagnostic biopsy for histopathological examination and sperm extraction for assisted reproductive treatment. Azoospermia was defined according to the WHO 2010 criteria [18]. No ejaculated spermatozoa were seen in two sequential semen analyses. SCOS was confirmed upon finding no spermatozoa within the testes and/ or epididymis upon microscopic examination of at least 50 tiny foci in biopsy specimens. OA was specified as containing a notable number of mature spermatozoa retrieved by TESE. The diagnosis was supported by the histological reports as well as hormonal findings for the two groups. All patients enrolled in this study had primary infertility. The NOA and OA groups were age-and body mass index (BMI) matched. They were nonsmokers and did not consume alcohol regularly.
Individuals with infertility related to defined medical conditions, such as pathologies involving vas deferens or epididymis, undescended testes (UDT), cryptorchidism, varicocele, a positive history of hormone therapy, retrograde ejaculation, current or previous infection with mumps, chromosome abnormalities, and AZF microdeletions were excluded from the study. After collecting, testicular samples were stored in RNALater® Solution (Behnogen, Iran) overnight at 4°C and placed at -80°C until further processing.
2.2. RNA Extraction. The total RNA was isolated from each sample using RiboEx reagent (GeneAll, Korea), following the manufacturer's recommendations. The elution was made in 30 μL of nuclease-free water. The quality and quantity of all RNA samples were evaluated by denaturing agarose gel electrophoresis and NanoDrop™ 1000 spectrophotometer, respectively.

Small RNA Library Preparation and Sequencing.
To undertaken sequencing, equal amounts of RNA samples were pooled to obtain 2 OA (4 samples each) and 3 NOA (4 : 3 : 3 samples each) groups. RNA integrity (RIN) of samples was tested using Agilent 2100 Bioanalyzer (Agilent RNA 6000 Nano Kit, USA), and samples with RIN values > 7 were considered acceptable for sequencing. BGISEQ-500 small RNA (sRNA) sample preparation kit was used to prepare the sRNA libraries. RNA (18-30 nt) molecules were ligated with adaptors; reverse transcribed and unique indices were incorporated during PCR amplification. After polyacrylamide gel electrophoresis (PAGE), the PCR amplicons were purified and dissolved in EB solution. The samples were loaded on the BGISEQ-500 platform at Yingbioteck (Shanghai, China), and sequencing by synthesis was conducted to generate single-end 50 base-pair reads [19].

Mapping Reads.
After trimming the adaptors and removing the contamination, the clean reads ≥ 18 nucleotides were mapped to themiRBase version 21 and Rfam databases. We used Bowtie2 [20] and cmsearch [21] to map filtered reads to the reference genome as well as other sRNA 2 Andrologia databases. Differentially expressed miRNA analysis was performed by DESeq2 algorithm, with a cutoff of log 2 FC > 1.

Bioinformatic
Analyses for Differentially Expressed miRNAs. In order to annotate and categorize the target genes attributed to cellular component (CC), biological process (BP), and molecular function (MF), Gene Ontology (GO) analysis (http://geneontology.org) was performed. Moreover, the KEGG (Kyoto Encyclopedia of Genes and Genomes) pathway enrichment analysis was used to categorize the genes targeted by the differentially expressed small RNAs (DESs) on the basis of the metabolic or signal transduction pathways they were involved in. Softwares including Miranda [22] and TargetScan (http://www.targetscan.org/ mmu_61/) [23] were used for predicting the targets. The default parameters were as follows: cutoff score ≥ 140 and energy ≤ −20 kcal/mol for Miranda and energy < −25 for TargetScan. The common targets identified by the two databases were used as the final miRNA target prediction results. We excluded the genes targeted by both low (<10 RPM) and highly (>100 RPM) abundant miRNAs to overcome the dominant contribution of these miRNA families leading to targeted gene downregulation as compared to those expressed at a moderate levels.
2.6. cDNA Synthesis and Quantitative Reverse Transcription PCR. To validate the expression level of selected miRNAs in the testis samples, cDNA synthesis for the same individual samples (n = 18) was performed with 1 μg of total RNA, utilizing TaKaRa PrimeScript II 1st strand cDNA synthesis kit (Takara Bio, Japan). The primers were designed according to the stem-loop RT method introduced by Chen et al. [24], and the details are described in Abbreviation: RT: reverse transcription; F: forward; R: universal reverse.  Color key    Table 2. The serum levels of follicle stimulating hormone (FSH) and luteinizing hormone (LH) in NOA patients were higher than those in OA patients, whereas the testosterone level in the NOA group was lower than that in the OA group.

Small RNA Sequencing and Mapping of Reads to Human
Reference Genome. Before carrying out data analysis, lowquality tags were removed, and after trimming, over 20 million reads were detected in each sample. The sequencing data for each sample, including the mapping ratio, are summarized in Table 3. The legnth distribution analysis showed that sRNAs had length between 18 and 30 nucleotides. The coverage percentage after cleaning was more than 90% in 4 samples.

Global miRNA Profiles between Azoospermia and
Normozoospermia. Using the pheatmap function, hierarchical clustering analysis of the miRNA's expression variation illustrated a set of aberrantly expressed miRNAs in the NOA patients (Figure 1). In total, 136 known miRNAs showed significant differential expression in the testis tissue between the NOA and OA individuals. Specifically, we observed that 120 miRNAs were downregulated in NOA patients in comparison to OA men, whereas 10 miRNAs were upregulated (adjusted P < 0:05, using fold change threshold values of ≥1 and ≤−1). The altered miRNAs were distributed across all chromosomes, except chromosomes 13,14,20,21,22, and Y. The maximum numbers of dysregulated miRNAs were mapped to chromosomes 11, 19, and X (Table 4).

Validation of Small RNA Sequencing Results by Real-Time PCR.
To validate the small RNA sequencing findings, qRT-PCR was conducted for two miRNAs; namely, miR-449a and miR-888, in individual tissue samples from NOA (n = 10) and OA (n = 8) patients. These miRNAs were chosen according to the fold change (FC). Both of these miRNAs were significantly downregulated in NOA samples (P < 0:05) (Figure 2). These findings were fully concordant with the differential expression observed by sequencing.

GO and Pathway Analyses for Targets of Differentially
Expressed miRNAs. To explore the functions of the genes regulated by the differentially expressed miRNAs, the target genes were subjected to Gene Ontology (GO) and pathway analyses. The GO analysis demonstrated the target genes to be enriched in GO terms of 25 biological processes, 18 molecular functions, and 16 cellular components (Figure 3). Cellular process, cell and cell part, and binding were the most enriched terms in the main categories of cellular components, biological processes, and molecular functions, respectively.
In addition, pathway enrichment analysis of the DES target genes was performed based on the KEGG database, and the results were reported as a bar plot for the statistics of KEGG term types (Figure 4). The most enriched pathway was signal transduction.

Discussion
In this study, we analyzed the miRNA expression profile in the testicular tissues of idiopathic azoospermic patients with SCOS and compared it with OA individuals with normal spermatogenesis. Using small RNA sequencing, 136 miR-NAs were recognized to be differentially expressed between the two groups. Most of these miRNAs belonged to the ascertained families, such as miR-515, miR-449, miR-34, miR-891, miR-743, miR-25, and miR-302, which are widely known and expressed as determined clusters. Specific or preferential expression of a number of miRNAs in the testis suggests their roles in spermatogenesis. In fact, some of the GO terms and KEGG pathways obtained in this study relate to sperm production, cell growth and death, signal transduction, replication and repair, energy metabolism, endocrine system, hormone secretion, reproductive process, cell motility, cell junction, and antioxidant activity (Figures 3 and 4).

Andrologia
CDKN1A, targeted by both miR-520h and miR-519c-3p, is found to be upregulated in undifferentiated spermatogonia due to DNA damage [28]. Similarly, the overexpression of CDKN1A has been reported in SCOS patients [29]. Down-regulations of its targetting miRNAs, miR-520h and miR-519c-3p, in our study is an interesting corroborating observation CCND1 is another gene, which is downregulated in SCOS patients, and is targeted by miR-520a-3p, miR-518a-5p, and miR-520c-3p. Strikingly, certain miR-515 family members also control the MAPK pathway genes ( Figure 5). Several studies have demonstrated that the MAPK pathway is involved in the differentiation and maturation of male germ cells. Specifically, cascade activation of MAPKs is suggested to be involved in regulating transcription and ectoplasmic specialization, flagellar motility of mature spermatozoa, acrosome reaction, and hyperactivation. For instance, p38δ (MAPK13) isoform of p38 mitogen-activated protein kinase is predominantly expressed and activated in XY germ cells [30]. Normal spermatogenesis is critically dependent on the integrity and maintenance of the bloodtestis barrier (BTB) as well as the adhesion between germ and Sertoli cells in the seminiferous tubules [31].
Accompanied by the ERK MAPK pathway activation, p38 may lead to the loss of germ cells from the epithelium by inducing the BTB disruption [32]. Moreover, p38 has been evidenced to play an essential role in sperm motility   Figure 3: The GO enrichment analysis for the targets of differentially expressed miRNAs. 9 Andrologia [33]. Homozygous deletion ofPRKACA (Prkaca -), another gene in this pathway, has been shown to result in substantial malformations in head morphology and forward motility of mature sperm in mice [34]. Two heat shock proteins (HSPs), HSPA8 and HSPA1B, targeted by the miR-515 family members, are involved in the MAPK pathway. HSPA8 is recognized to contribute to the maintenance of TSSK6 structure expressed in the late stages of spermiogenesis. Tssk6 disruption is linked with diminished sperm number, impaired condensation of DNA, and abnormal sperm shape and motility in male mice [35]. Rac1 in the Sertoli cells acts in the progression of spermatogenesis. Conditional knockout of Rac1 in adult mice results in spermatogenic arrest at the round spermatids stage and severe polarity disruption and increased apoptosis in the Sertoli cell [36]. The above information suggests that the miR-515 family regulates spermatogenesis through the MAPK pathway; however, the precise role of each miRNA is yet to be unveiled as the functions of many genes targeted by this family remain unknown. Dysregulation of miR-515, miR-518a-3p, miR-518d-5p, miR-518e-3p and -5p, miR-518b, miR-517-5p, miR-519d-3p, miR-520a-5p, miR-520f-5p, and miR-520h is consistent with their dysregulation in Klinefelter syndrome patients [37].
We observed that the miR-34 (miR-34b-3p, miR-34b-5p, miR-34c-3p, and miR-34c-5p) and miR-449 (miR-449a, miR-449b-3p, miR-449b-5p, and miR-449c-5p) family   10 Andrologia members located on chromosomes 5 and 11, respectively, were significantly dysregulated in SCOS patients. In a previous study, drastic reductions in the expressions of both the miR-34 and miR-449 families were observed in three groups of patients with impaired spermatogenesis, including SCOS, meiotic arrest, and mixed atrophy [12]. In another study, altered expressions of some members of the miR-449b and miR-34b/c families were observed even in patients who carried AZF microdeletions [16]. Both the above studies on testicular miRNAs were microarray based. In contrast to Abu-Halima et al. [13], our study was on a different population and used miRNA sequencing, which is an advanced and more sensitive method of gene expression analysis. In comparison to Noveski et al. [16], our samples were ruled out for Y-deletions. The miR449 family members (miR-449b and miR-449c) were not detectable in azoospermic patients with the Klinefelter syndrome [37]. Consistent dysregulation of these miRNAs across three studies signifies their role in spermatogenesis, making them excellent candidates for further investigations. Yuan et al. demonstrated that male miR-34b/c-or miR-449-knockout mice displayed normal spermatogenesis, resulting in normal fertilization and consequently normal preimplantation development, while miR-449 and miR-34b/c double knockout (miR-dKO) mice developed infertility in the form of severe spermatogenic impairment and oligoasthenoteratozoospermia (OAT) syndrome [38]. All these miRNAs, except miR-449b-3p, are involved in the regulation of hundreds of genes, which are parts of common biological processes and pathways, including mitotic cell cycle (GO:0000278), cell death (GO:0008219), response to stress (GO:0006950), cell cycle (GO:0007049), apoptosis (hsa04210), p53 signaling pathway (hsa04115), fatty acid metabolism (hsa01212), and fatty acid biosynthesis (hsa00061) [26]. Nonetheless, more investigations on these two families are needed since aberrations of single family did not lead to infertility. The miR-34 and miR-449 families are involved in ciliogenesis and cell cycle regulation through BMP and notch signaling, suggesting their importance in spermatogenesis [39] (Figure 5).
Two other families, including miR-891 (hsa-miR-891a-3p, hsa-miR-891a-5p, and hsa-miR-891b) and miR-743 (hsa-miR-888-5p, hsa-miR-890, hsa-miR-892a, and hsa-miR-892b), are located on the X chromosome as their members appear to express in the miR-888 cluster. No related GO terms or KEGG pathways have yet been identified for these miRNAs, although some roles in spermatogenesis have been suggested. miR-891b, miR-890, miR-892a, and miR-892b are speculated to be involved in sperm maturation insomuch that their expressions are predominately observed in the epididymis of human and other primates [40,41], indicating their contributionto sperm maturation ( Figure 5). miR-92a-2-5p and miR-25-5p, belonging to the miR-25 family, act at the mitotic G1/S transition checkpoint (GO:0044819). miR-92a-2-5p also plays a role in cell motility (GO:0048870) and cell proliferation (GO:0008283) [26]. miRNA-106a, located on chromosome X, is indispensable for self-renewal of the SSCs by targeting Ccnd1 and Stat3 in mouse [42]. The miR-25 and miR-302 families have been investigated for their roles in cell division and apoptosis in cancer, which suggests their significant contributions to the regulation of cell cycle and apoptosis in the testis ( Figure 5); however, further investigations are required to understand their exact roles in spermatogenesis.

Conclusions
In summary, this is the first study analyzing the expression profile of miRNAs through small RNA sequencing in the testicular tissue of idiopathic SCOS patients. We detected 136 differentially expressed miRNAs. While the differences in the germ cell composition in obstructive versus nonobstructive azoospermia cases could be behind the alterations identified in this study, irrespective of the reasons behind these differences, the identification of the missing/downregulated miRNAs would put forth the miRNAs that are critical miR-516b-5p miR-517b-3p miR-517b-5p miR-518a-3p miR-518b miR-519c-3p miR-520a-3p miR-520h miR-1283 11 Andrologia to infertility treatment in the future. These miRNAs may be used as prognostic markers for the treatment of azoospermia using stem cell therapy in the future. Further studies on these miRNAs may be required to investigate their biological functions in the testis. We admit the limitation of a limited sample size in this study, which was due to the difficulty in obtaining testicular biopsy samples.

Data Availability
The data from which the findings of this study have been extracted are available from the corresponding author upon reasonable request.

Ethical Approval
The current research was first approved by the Ethical Committee of Tarbiat Modares University (IR.TMU.REC.1395.382).

Consent
All individuals referring to the Royan Institute (Tehran, Iran) were recruited in this study after giving informed consent.

Disclosure
A preprint version of this manuscript has previously been published [43].