Integrative Bioinformatics Analyses of Potential Therapeutic Targets for Non-Obstructive Azoospermia

International Clinical and Medical Case Reports Journal Research Article (ISSN: 2832-5788) Int Clinc Med Case Rep Jour (ICMCRJ) 2023 | Volume 2 | Issue 6 Integrative Bioinformatics Analyses of Potential Therapeutic Targets for Non-Obstructive Azoospermia Yalin Ju†, Yuxin Ma†, Yaohui Jiang, Xiaoyu Zhou, Yujiao Dong, Shaobo Wu* Xi’an Jiaotong University Health Science Center, Xi’an 710061, Shaanxi Province, PR China †These authors contributed equally to this work and should be considered co-first authors. Citation: Yalin Ju, Yuxin Ma, Yaohui Jiang, Xiaoyu Zhou ,Yujiao Dong, Shaobo Wu. Integrative Bioinformatics Analyses of Potential Therapeutic Targets for Non-Obstructive Azoospermia. Int Clinc Med Case Rep Jour. 2023;2(6):1-9. Received Date: 07 February, 2023; Accepted Date: 10 February, 2023; Published Date: 12 February, 2023 *Corresponding author: Shaobo Wu. Xi’an Jiaotong University Health Science Center, Xi’an 710061, Shaanxi Province, PR China Copyright: © Shaobo Wu, Open Access 2023. This article, published in Int Clinc Med Case Rep Jour (ICMCRJ) (Attribution 4.0 International), as described by http:// creativecommons.org/licenses/by/4.0/. ABSTRACT Background: Non-obstructive azoospermia (NOA) is a vital form of male infertility, which is often suspected to be linked to undefined genetic abnormalities. It is essential to explore potential biomarkers and accurate molecular mechanisms that might be associated with NOA and move towards the potential application value for the diagnosis and treatment. Methods: We acquired the microarray expression profiles of NOA patients from the GEO database. Using the limma package in R software, differentially expressed genes (DEGs) in NOA were identified. Afterwards, we analyzed GO enrichment and KEGG pathway analysis, then used STRING to construct protein-protein interaction (PPI) network of DEGs. Cytoscape software was used to visualize the PPI network for DEGs and explore the hub genes that were contained in the PPI network. Finally, we tried to predict potential drugs or molecular compounds that interacted with the hub genes and visualized drug–gene interaction networks by the Cytoscape software. Results: We integrated co-expressed DEGs in two datasets and 136 up-regulated genes and 311 down-regulated genes were screened out. The results of GEO analysis displayed that CC and BP that DEGs were obviously enriched in were associated with spermatogenesis. KEGG indicated that DEGs were mainly enriched in Antigen Processing and Presentation and Staphylococcus Aureus Infection. We identified 10 central genes (CDCA8, CENPF, CCNB2, MND1, TYMS, RACGAP1, HMMR, RFC4, PTTG1 and KIF15) from the PPI network. 46 drugs or molecular compounds differentially regulated the expression of TYMS, five regulated HMMR, and 70 were found to interact with RACGAP1. Conclusion: These findings demonstrate that the identification of the above hub genes and potential drugs helps us get a deeper understanding of the mechanisms associated with NOA and offer potential biomarkers for its diagnosis and treatment. Keywords: Non-obstructive azoospermia (NOA); Biomarkers; Hub genes; Functional enrichment analysis International Clinical and Medical Case Reports Journal Research Article (ISSN: 2832-5788) Int Clinc Med Case Rep Jour (ICMCRJ) 2023 | Volume 2 | Issue 6 INTRODUCTION Male infertility has become an increasingly serious social problem in modern society, and nonobstructive azoospermia (NOA), which accounts for 10% of all infertile men, is an important cause for this.[1] Nonobstructive azoospermia patients, characterized by the lack of sperm in ejaculate,[2] mainly manifested as severe oligospermia, weak sperm, teratospermia or azoospermia.[3] Previously, patients with non-obstructive azoospermia were once considered to be infertile with few treatment options due to the absence of sperm in the ejaculate.[4] Clinical varicoceles have been reported in 4.3% to 13.3% of infertile men with nonobstructive azoospermia (NOA) and some studies have indicated that varicocelectomy, a cost-effective treatment, may improve sperm production.[5,6] Clinically, nonobstructive azoospermia can be treated only by micromanipulation-assisted fertilization using spermatogenic cells including spermatids.[7] However, as for nonobstructive azoospermia, only a few patients can benefit from intracytoplasmic sperm injection (ICSI), a technology that revolutionized the treatment of couples with male factor infertility,[8] and finally get their babies.[9] Thus, it is essential to explore potential biomarkers and molecular mechanisms that might be associated with nonobstructive azoospermia and move towards the potential application value for the diagnosis and treatment. The molecular basis of male infertility remains largely unknown.[10] non-obstructive azoospermia is often suspected to be linked to currently undefined genetic abnormalities.[11] Currently, many study have been devoted to finding novel targets for molecular therapy of male infertility. More and more studies have shown the importance of epigenetic modifications, including DNA methylation, in male infertility.[12] For example, Fu has demonstrated that miRNA-31-5p regulates the proliferation, DNA synthesis, and apoptosis of human SSCs by the PAK1-JAZF1-cyclin A2 pathway, which might provide novel targets for molecular therapy of male infertility.[13] Tu has suggested that Prmt6 is a non-obstructive azoospermia-susceptible locus.[14] In additions, accumulating studies have indicated the significance of noncoding RNAs (ncRNAs) in male infertility.[12] Recently with the advent of genomic techniques, a number of genes responsible for non-obstructive azoospermia were identified, NR5A1 (SF1), HSF2, SYCP3, SOHLH1, and TEX1110-14.[15] In our study, we identified differentially expressed genes (DEGs) for NOA by analysing two mRNA expression profiles that were downloaded from the GEO database. And then we found common DEGs between these two datasets via Venn diagrams. Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) were used to annotate the functions of these DEGs. Next, we constructed the PPI network with STRING and explored the hub genes with Cytohubba. Finally, we predicted the interaction between the hub genes and potential drugs or molecular compounds by the DGIdb database. Our realization is to explored promising biomarkers for potential conversion application value for the diagnosis or treatment of NOA. MATERIALS AND METHODS Gene expression microarray data We acquired the microarray expression profiles of non-obstructive azoospermia patients from the National Center for Biotechnology Information GEO database (http://www.ncbi.nlm.nih.gov/geo/) database. The GSE9210 dataset, contains 58 samples, including 47 non-obstructive azoospermia (NOA) and 11 obstructive International Clinical and Medical Case Reports Journal Research Article (ISSN: 2832-5788) Int Clinc Med Case Rep Jour (ICMCRJ) 2023 | Volume 2 | Issue 6 azoospermia (OA) patients that served as controls. The GSE145467 dataset, contains 20 samples, including 10 non-obstructive azoospermia (NOA) and 10 obstructive azoospermia (OA) patients. Identification of differentially expressed genes Firstly, the probe group is transformed into the gene symbol by the annotation information of the data set's respective platform. If multiple probe groups correspond to the same gene, the average value is calculated. The data were normalized and the differential genes (DEGs) were identified using the limma package (|logFC|>1 and p<0.05) in the R software (version 4.0.3).[16] Then, we used the pheatmap package to draw DEGs heat maps and the ggplots2 package to draw DEGs volcano maps for each dataset. An online tool (http://www.bioinformatics.com.cn) was used to plot Venn diagrams to detect common DEGs between two datasets. Functional and pathway enrichment analysis Go analysis is widely used to identify the characteristic biological properties of genes, including biological processes (BP), cellular components (CC) and molecular functions (MF). KEGG analysis provides a comprehensive set of genome sequences and biological interpretation of protein interaction network information. In our study, GO terms and KEGG pathway enrichment analysis of DEGs were completed by the clusterProfiler package17 and visualized by the ggplot2 package. The significance criterion to screen the significantly enriched GO terms and KEGG pathways was set at p <0.05. PPI network and identification of hub genes STRING (http://www.string-db.org/), which is an online tool to identify and predict interactions between genes or proteins. We used it to construct the PPI network of DEGs with the cut-off standard as a combined score >0.4.[18] Then, Cytoscape software was used to visualize the PPI network for DEGs. The MCODE[19] is a plug-in of Cytoscape, which was applied to identify significant modules. Moreover, CytoHubba[20] is another plug-in of Cytoscape, which was employed to study essential nodes in the network with 12 methods (The MCC shows satisfactory comparative performance) completed to explore the hub genes that were contained in the PPI network. Identification of the potential drugs The Drug Gene Interaction Database (DGIdb) version4.2.0 (https://www.dgidb.org)[21] is an available resource of the drug targeted and sensitive genome and drug-gene interaction. It includes 42 potentially druggable categories and least 49 interaction types as defined by source datasets. These include inhibitors, activators, cofactors, ligands, vaccines, and in many cases, interactions of unknown type. We used the DGIdb to predict potential drugs or molecular compounds that interacted with the hub genes and visualized drug–gene interaction networks by the Cytoscape software. RESULTS Identification of DEGs To identify DEGs from the non-obstructive azoospermia and obstructive azoospermia datasets, we extracted 794 and 4287 DEGs from GSE9210 and GSE145467, respectively, according to defined criteria. GSE9210 contains 313 up-regulated genes and 481 down-regulated genes, while GSE145467 contains 1997 up-regulated genes and International Clinical and Medical Case Reports Journal Research Article (ISSN: 2832-5788) Int Clinc Med Case Rep Jour (ICMCRJ) 2023 | Volume 2 | Issue 6 2290 down-regulated genes, which are shown by heat map and volcano map respectively (Figure 1). In addition, we used Venn map to integrate co-expressed DEGs, including 136 up-regulated genes and 311 down-regulated genes in NOA (Figure 2). GO and KEGG pathway enrichment analyses The GO function and KEGG enrichment were analyzed by using the clusterProfiler and ggplot2 package to further understand the function and mechanism of the identified DEGs. We selected the top five enriched terms in BP, CC, and KEGG are shown in Figure 3. The results of GO analysis displayed that DEGs were obviously enriched in cellular components containing MHC Protein Complex, Ciliary Part, Acrosomal Vesicle, Motile Cilium and Sperm Part. For biological processes, the DEGs were enriched in Fertilization, Germ Cell Development, Cellular Process Involved in Reproduction in Multicellular Organism, Spermatid Differentiation and Spermatid Development (Figure 3B, Table 1). In addition, the result of KEGG pathway enrichment indicated that DEGs were mainly enriched in Type I Diabetes Mellitus, Graft-Versus-Host Disease, Allograft Rejection, Antigen Processing and Presentation and Staphylococcus Aureus Infection (Table 2). PPI network analysis and hub gene selection The interactions between DEGs were downloaded from the STRING database and visualized using Cytoscape. Based on these standards, an important module (Figure 4A), including 18 nodes and 302 edges, was selected through cluster analysis of PPI network in MCODE. KEGG pathway enrichment analysis was then performed for the DEGs in the module, and the results are listed in Table 3. The DEGs in module was significantly enriched in Cell cycle, Oocyte meiosis, One carbon pool by folate, Human T-cell leukemia virus 1 infection and Mismatch repair. Moreover, hub genes were identified by CytoHubba in the present study; the 10 central genes (CDCA8, CENPF, CCNB2, MND1, TYMS, RACGAP1, HMMR, RFC4, PTTG1 and KIF15) were selected by the MCC method as shown in Figure 4B. Identification of the potential drugs DGIDB was used to determine drugs or molecular compounds that may reverse DEGS expression in nonobstructive azoospermia. As shown in the drug–gene interaction network (Figure 5), 46 drugs or molecular compounds included deoxy uridine monophosphate, trifluridine, thymidine monophosphate and raltitrexed, which differentially regulated the expression of TYMS. In addition, five drugs or molecular compounds that included hyaluronan, hyaluronate sodium, epirubicin, fluorouracil and cyclophosphamide regulated HMMR. It's worth noting that fluorouracil co-regulated the expression of HMMR and TYMS. Further, 70 drugs or molecular compounds, for instance, benzidine, pinoresinol, xamoterol fumarate and phenidone, were found to interact with RACGAP1 (Figure 5A). International Clinical and Medical Case Reports Journal Research Article (ISSN: 2832-5788) Int Clinc Med Case Rep Jour (ICMCRJ) 2023 | Volume 2 | Issue 6 Figure 1: Heatmap and Volcano plots of DEGs. A, Heatmap of GSE9210. B, heatmap of GSE145467. Red areas represent highly expressed genes and green areas represent lowly expressed genes in NOA. Darker color indicates the higher multiple of DEGs. C, Volcano plots of GSE9210. D, volcano plots of GSE145467. The xaxis label represents log2FoldChange and the y-axis label represents–log10 (adjusted p-value). Data points in red represent upregulated, and green represent downregulated genes. No signifificantly-changed genes are marked as gray points. DEGs: differentially expressed genes; NOA: non-obstructive azoospermia. International Clinical and Medical Case Reports Journal Research Article (ISSN: 2832-5788) Int Clinc Med Case Rep Jour (ICMCRJ) 2023 | Volume 2 | Issue 6 Figure 2: Venn diagram of co-expressed DEGs from GSE9210 and GSE145467 in NOA. A, 136 up-regulated genes. B, 311 down-regulated genes. Figure 3: GO enrichment and KEGG pathway analysis of common DEGs. A, the bubble chart of BP and CC mainly enriched in GO enrichment. B, the bar chart of BP and CC mainly enriched in GO enrichment. C, the bubble chart of top five enriched pathways in KEGG. D, the bar chart of top five enriched pathways in KEGG. BP, biological process; CC, cellular component; MF, molecular function; DEGs, differentially expressed genes. International Clinical and Medical Case Reports Journal Research Article (ISSN: 2832-5788) Int Clinc Med Case Rep Jour (ICMCRJ) 2023 | Volume 2 | Issue 6 Figure 4: PPI network construction and hub gene regulatory network. A, the PPI network of DEGs. B, the interactions between 10 central genes selected from the PPI network. PPI, protein–protein interaction; DEGs, differentially expressed genes. Figure 5: Identification of the potential drugs. A, 70 drugs or molecular compounds differentially regulating the expression of RACGAP1. B, 46 and 5 drugs or molecular compounds interacting with TYMS and HMMR respectively, among which fluorouracil co-regulated the expression of both HMMR and TYMS. Table 1: Go enrichment analysis of DEGs. Ontology ID Description P value BP GO:0007286 spermatid development 3.69E-13 BP GO:0048515 spermatid differentiation 8.00E-13 BP GO:0022412 cellular process involved in reproduction in multicellular organism 2.05E-12 BP GO:0007281 germ cell development 6.05E-11 BP GO:0009566 fertilization 9.72E-11 CC GO:0097223 sperm part 3.96E-11 CC GO:0031514 motile cilium 2.71E-08 CC GO:0001669 acrosomal vesicle 7.41E-08 CC GO:0044441 ciliary part 1.21E-05 CC GO:0042611 MHC protein complex 7.26E-05 Table 2: KEGG pathway enrichment analysis of DEGs. Ontology ID Description P value KEGG hsa05150 Staphylococcus aureus infection 2.70E-04 KEGG hsa04612 Antigen processing and presentation 4.18E-04 KEGG hsa05330 Allograft rejection 4.93E-04 KEGG hsa05332 Graft-versus-host disease 7.89E-04 KEGG hsa04940 Type I diabetes mellitus 8.81E-04 International Clinical and Medical Case Reports Journal Research Article (ISSN: 2832-5788) Int Clinc Med Case Rep Jour (ICMCRJ) 2023 | Volume 2 | Issue 6 Table 3: KEGG pathway enrichment analysis of module. Ontology ID Description P value KEGG hsa04110 Cell cycle 0.008 KEGG hsa04114 Oocyte meiosis 0.008 KEGG hsa00670 One carbon pool by folate 0.022 KEGG hsa05166 Human T-cell leukemia virus 1 infection 0.023 KEGG hsa03430 Mismatch repair 0.025 DISCUSSION Spermatogenesis is a complicated and subtle process including spermatogonial mitosis and meiosis and spermatid differentiation22. Du et al. believed that genetic regulation plays a vital role in NOA.[23] Therefore, in the present study, we aimed to find potential specific gene targets to provide more therapeutic thoughts. We integrated the GSE9210 and GSE145467 datasets and utilized Venn map to integrate co-expressed DEGs, including 136 up-regulated genes and 311 down-regulated genes in NOA. GO enrichment analysis of DEGs showed that the mainly enriched BP and CC are associated with spermatogenesis. Based on the module established through cluster analysis of PPI network in MCODE, we selected 10 central genes: CDCA8, CENPF, CCNB2, MND1, TYMS, RACGAP1, HMMR, RFC4, PTTG1 and KIF15, which we will describe and discuss in detail. During cytokinesis, the Chromosomal passenger complex (CPC) regulates furrow contractility and specifies the cleavage plane,[24,25] which is essential for transmission of the genome 26. In 2001, Walker found previouslyunrecognized CDCA8 as one of the components of CPC in mitosis[27]. A recent study shows that deletion of Haspin, a kinase to phosphorylate threonine 3 of histone H3, led to misalignment of chromosomes, apoptosis of metaphase spermatocytes and a decreased number of sperm by deregulation of CPC28.Cyclin B1 and B2 mRNAs were expressed in all stages of spermatogenic cells except for spermatozoa during the process of medaka spermatogenesis. Cyclin B1 controls the meiotic cell cycle, whereas cyclin B2 is involved in process other than meiosis.[29] Cyclin B2(CCNB2) is expressed in spermatogonia. It plays a vital role in the conversion of G2/M and in ensuring the smooth completion of the M phase[30]. Lin et al. found that the mRNA transcription levels of CCNB2 was significantly reduced in NOA patients with failed sperm extraction[31]. Therefore, CCNB2 might plays an important part in spermatogenesis. Centromere protein F (CENPF) is the principal sites of interaction between chromosomes with spindle microtubules.[32] CENPF regulates chromosome segregation during mitosis, and silencing of CENPF induces mitotic defects.[33] A study of embryogenesis in mice showed that depletion of CENPF in zygotes results in abnormalities and impairment of embryo development in specifically.[34] MND1, a spermatogenesis-related gene, is expressed before meiotic division and its corresponding protein contributes to the meiosis process.[35] TS, a folate-dependent protein encoded by TYMS, is responsible for the reductive methylation of dUMP to deoxythymidine monophosphate (dTMP) using the oxidation of 5,10-methylene tetrahydrofolate (5,10-MTHF) to dihydrofolat, playing an important role in DNA synthesis and repair.[36,37] RACGAP1 is a GTPase Activating Protein (GAP) and is prviously found to be highly produced in the human and mouse post-natal testis.[38] One study reported that the deletion of RACGAP1 led to azoospermia, a spermatogenesis arrest in the PLZF-positive population of undifferentiated spermatogonia.[39] Hyaluronan-mediated motility receptor (HMMR), a hyaluronan receptor, is localized along the tail, midpiece International Clinical and Medical Case Reports Journal Research Article (ISSN: 2832-5788) Int Clinc Med Case Rep Jour (ICMCRJ) 2023 | Volume 2 | Issue 6 and the head of sperm and is linked to sperm motility.[40] A recent research suggested that the abundance level of HMMR gene in normozoospermic men is higher than in oligoasthenozoospermic men.[41] RFC4 is a member of a pentameric complex consisting of RFC1-5 and functions to loads proliferating cell nuclear antigen (PCNA) onto DNA during DNA replication.[42] Gray et al. observed a strong staining of RFC4 in the nuclei of spermatogonia and primary spermatocytes in testis sections.[43] Two other central genes, namely PTTG1 and KIF15, have not been reported in the terms of spermatogenesis. In the end, we studied the relationship between the potential drug therapeutic targets and central genes and displayed three drug–gene interaction networks-TYMS, RACGAP1 and HMMR. As is mentioned above, the sperm count, motility, and morphology decreased because less HMMR binds to miR-23a/b-3p in the 3’UTRs in the sperm of men with NOA. Therefore, it is possible to utilize drugs to regulate relative gene expressions in order to reduce NOA. Potential drugs like fluorouracil, co-regulating the expression of TYMS and HMMR, hyaluronan, regulating HMMR and pinoresinol which interacts with RACGAP1. However, these drugs have not been studied in NOA and further research is required. DECLARATION OF COMPETING INTEREST The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. REFERENCES 1. Li X, Sun T, Wang X, Tang J, Liu Y. Restore natural fertility of Kit(w)/Kit(wv) mouse with nonobstructive azoospermia through gene editing on SSCs mediated by CRISPR-Cas9. Stem Cell Res Ther. 2019;10(1):271. 2. Hwang YS, Suzuki S, Seita Y, Ito J, Sakata Y, Aso H, et al. Reconstitution of prospermatogonial specification in vitro from human induced pluripotent stem cells. Nat Commun. 2020;11(1):5656-5656. 3. Lee JY, Dada R, Sabanegh E, Carpi A, Agarwal A. Role of genetics in azoospermia. Urology. 2011;77 (3):598-601. 4. Dabaja AA, Schlegel PN. Microdissection testicular sperm extraction: an update. Asian J Androl. 2013; 15(1):35-39. 5. Zhi EL, Liang GQ, Li P, Chen HX, Tian RH, Xu P, et al. Seminal plasma miR-192a: a biomarker predicting successful resolution of nonobstructive azoospermia following varicocele repair. Asian J Androl. 2018;20(4):396-399. 6. Pathak P, Chandrashekar A, Hakky TS, Pastuszak AW. Varicocele management in the era of in vitro fertilization/intracytoplasmic sperm injection. Asian J Androl. 2016;18(3):343-348. 7. Liu XY, Miao YL, Zhang J, Qiu JH, Cui XZ, Gao WQ, et al. Effects of activation on functional aster formation, microtubule assembly, and blastocyst development of goat oocytes injected with round spermatids. Cell Reprogram. 2012;14(5):436-447.