Abstract Numerous variables that regulate the metabolism of Sertoli cells and sperm have been identified, one of which is sex steroid hormones. These hormones play a vital role in maintaining energy homeostasis, influencing the overall metabolic balance of the human body. The proper functioning of the reproductive system is closely linked to energy status, as the reproductive axis responds to metabolic signals. The aim of this study was to investigate the gene expression patterns of metabolite interconversion enzymes in testicular cells (Sertoli cells and spermatogonia) of non-obstructive azoospermia (NOA) patients, as compared to normal controls, to understand the molecular mechanisms contributing to NOA. We used microarray and bioinformatics techniques to analyze 2912 genes encoding metabolite interconversion enzymes, including methyltransferase, monooxygenase, transmembrane reductase, and phosphohydrolase, in both testicular cells and normal samples. In sperm, the upregulation of MOXD1, ACAD10, PCYT1A, ARG1, METTL6, GPLD1, MAOA, and CYP46A1 was observed, while ENTPD2, CPT1C, ADC, and CYB5B were downregulated. Similarly, in the Sertoli cells of three NOA patients, RPIA, PIK3C3, LYPLA2, CA11, MBOAT7, and HDHD2 were upregulated, while NAA25, MAN2A1, CYB561, PNPLA5, RRM2, and other genes were downregulated. Using STRING and Cytoscape, we predicted the functional and molecular interactions of these proteins and identified key hub genes. Pathway enrichment analysis highlighted significant roles for G1/S-specific transcription, pyruvate metabolism, and citric acid metabolism in sperm, and the p53 signaling pathway and folate metabolism in Sertoli cells. Additionally, Weighted Gene Co-expression Network Analysis (WGCNA) and single-cell RNA sequencing (scRNA-seq) were performed to validate these findings, revealing significant alterations in gene expression and cellular distribution in NOA patients. Together, these results provide new insights into the molecular mechanisms underlying NOA and identify potential therapeutic targets. Supplementary Information The online version contains supplementary material available at 10.1038/s41598-024-77875-9. Keywords: Non-obstructive azoospermia, Gene expression, Sperm, Sertoli cell, Metabolite Subject terms: Cell biology, Genetics, Endocrinology, Urology Introduction Male factor infertility is increasing in developed countries due to lifestyle factors affecting spermatogenesis^[36]1. Sertoli cells and spermatogonia stem cells (SSCs) play crucial roles in supporting germ cell development, primarily metabolizing glucose into lactate, which serves as the main energy source for germ cells^[37]2. The lactate metabolism is regulated by various mechanisms, including glucose transporters, lactate dehydrogenase, and monocarboxylate transporters^[38]3. Hormonal factors, particularly sex steroid hormones, influence Sertoli cell metabolism, impacting overall energy balance. Normal reproductive function is closely linked to metabolic status, as the reproductive axis responds to metabolic cues^[39]4. During spermatogenesis, Sertoli cells face challenges from both external and internal factors like environmental influences, hormonal imbalances, diseases, and oxidative stress^[40]5,[41]6. These factors can disrupt the ideal environment needed for germ cell development. Germ cells have unique nutritional needs during different stages of sperm formation, adjusting their metabolism accordingly^[42]7. Additionally, the testes naturally have low oxygen levels, which may contribute to germ cells using diverse metabolic pathways for energy at different developmental stages. Although glycolysis is a conserved process across species, many glycolytic enzymes have specialized forms found primarily in testicular cells, particularly during the postmeiotic phase of spermatogenesis^[43]8. Germ cells rely heavily on carbohydrate metabolism, utilizing both aerobic and anaerobic pathways. Lactate plays a crucial role in stimulating RNA and protein synthesis in spermatids while also protecting germ cells from apoptosis^[44]9,[45]10. Although the glycolytic potential in spermatids is lower compared to earlier developmental stages, both glycolytic and gluconeogenic pathways may be active in round spermatids due to lactate recycling^[46]11. Spermatozoa have the highest glycolytic activity among germ cells and primarily use glucose or fructose for energy, with limited tricarboxylic acid (TCA) cycle activity. Energy production in spermatozoa is compartmentalized, with mitochondria and oxidative phosphorylation occurring in the midpiece, while glycolysis occurs in the principal piece^[47]12,[48]13. The role of energy in the acrosome reaction is uncertain. The pentose phosphate pathway is active in germ cells, particularly in spermatocytes, supporting RNA biosynthesis and producing nicotinamide adenine dinucleotide phosphate (NADPH) and ribose 5-phosphate^[49]14. Lipid metabolism involves transport of long polyenes from Sertoli cells to germ cells, with spermatozoa undergoing lipid remodeling in the epididymis, and pachytene spermatocytes metabolizing fatty acids. The energy production capacity of germ cells relies on functional Sertoli cells, as their failure leads to unsuccessful spermatogenesis despite the germ cells’ ability to utilize various energy sources at different stages^[50]14,[51]15. Spermatogenesis is a highly conserved process across vertebrates, ensuring cyclic waves of sperm production and fertility by balancing self-renewal and differentiation of spermatogonial stem cells in the adult testis. This cyclic activity involves high rates of mitochondrial oxygen consumption and reactive oxygen species (ROS) generation^[52]16–[53]18. Enzyme complexes in the mitochondrial respiratory chain, such as xanthines, nicotinamide adenine dinucleotide phosphate (NADPH) oxidase, and cytochrome P450, contribute to ROS production. ROS encompass free radicals and oxygen derivatives like superoxide anion, hydrogen peroxide, hydroxyl radical, lipid hydroperoxides, peroxyl radicals, and peroxynitrite^[54]18. Their effects can be either beneficial or harmful depending on their type, concentration, location, and duration of exposure within biological systems^[55]19,[56]20. Metabolite interconversion enzymes, such as methyltransferases, monooxygenases, transmembrane reductases, and phosphohydrolases, play crucial roles in maintaining cellular energy homeostasis, which is vital for spermatogenesis. These enzymes are involved in processes such as DNA methylation (methyltransferases), which regulates gene expression; oxidative stress regulation (monooxygenases); membrane transport and signaling (transmembrane reductases); and energy metabolism (phosphohydrolases). Disruptions in these metabolic pathways could potentially contribute to the development of NOA by impairing the normal function of Sertoli cells and spermatogonia, leading to defective spermatogenesis. On the other hand, high doses or insufficient removal of reactive oxygen species (ROS) can pose significant dangers, leading to modifications of vulnerable molecules like DNA, lipids, and proteins^[57]21,[58]22. Various mechanisms such as ionizing radiation, xenobiotic metabolism, inflammation, increased cellular metabolism, and activation of oxidases and oxygenases can contribute to ROS generation^[59]23. The testis, rich in highly unsaturated fatty acids, is particularly susceptible to ROS damage. However, the low oxygen levels in testicular tissue serve as a defense mechanism against free radical-mediated damage during spermatogenesis and Leydig cell steroidogenesis^[60]24,[61]25. Additionally, an elaborate array of antioxidant enzymes and scavengers help maintain the functions of spermatogenesis and Leydig cells despite ROS exposure. To understand how ROS are neutralized or limited in the testis by antioxidant systems, we outline the major pathways of ROS generation and the mechanisms of antioxidative defense^[62]26. Specialized enzymes like xanthine oxidases and NADPH oxidases, as well as cellular metabolism, particularly the mitochondrial electron transport chain, generate superoxide radicals, which are then converted to hydrogen peroxide by superoxide dismutase (SOD)^[63]26,[64]27. Hydrogen peroxide, along with superoxide radicals and iron, forms more reactive species, ultimately leading to lipid peroxidation. This lipid peroxide is detoxified to water by enzymes like glutathione peroxidase (GPx) or glutathione-S-transferase (GST). The defense against SOD is typically facilitated by catalase or peroxidases like GPxs, which utilize reduced glutathione (GSH) as an electron donor^[65]28,[66]29. Glutathione maintains cellular redox balance by acting as an electron donor for various antioxidative enzymes. It also forms conjugates with harmful endogenous and xenobiotic compounds through the catalysis of glutathione-S-transferase (GST)^[67]30. The levels of reduced glutathione (GSH) are regulated through ATP-consuming steps involving γ-glutamylcysteine synthetase (cGCS) and glutathione synthetase. Alternatively, glutathione reductase (GR) recycles oxidized glutathione (GSSG) back to GSH in an NADPH-dependent manner^[68]31. In the presence of reactive oxygen species (ROS), GSH acts as an electron donor, leading to the formation of oxidized glutathione (GSSG), which is either recycled by GR or removed from the cells. Thus, GR indirectly contributes to cellular protection against oxidative stress^[69]32,[70]33. Additionally, the testis contains small molecular weight antioxidant substances such as zinc, various free radical scavengers, vitamins C and E, melatonin, and cytochrome C, which further protect against oxidative damage^[71]34. The role of metabolite interconversion enzyme genes, such as hydroxysteroid 11-beta dehydrogenase 1, glucosamine (UDP-N-acetyl)-2-epimerase/N-acetylmannosamine kinase, hydroxysteroid 11-beta dehydrogenase 2, and serine hydroxymethyltransferase 2, in sperm function and fertility remains uncertain, particularly in vitro settings. Some studies have investigated their association with specific stages of the human seminiferous epithelium. However, there is limited research on the expression of these genes in male germ cells. In this experimental study, we conducted microarray and in-silico analyses to demonstrate and analyze the expression of metabolite interconversion enzyme genes in sperm and Sertoli cell tubules. Material and methods Patient and control selection According to our previous studies, three patients with NOA induced by abnormal spermatogenesis were compared to three control cases. We were scheduled for a testicular biopsy for genetic testing and then an intracytoplasmic injection of the oocyte. The testicular biopsies and sperm extraction were not scheduled to coincide with the oocyte retrieval and intracytoplasmic sperm injection (ICSI) cycles, which were performed 4–5 months later^[72]35–[73]38. Histological sections The testes were horizontally divided into two distinct parts. From each testis, three testicular slices were obtained, encompassing the top, lower, and middle regions. A total of 48 h were allocated for the fixation of all testicular slices in Bouin’s solution, which consisted of 7.5 ml of saturated picric acid, 2.65 ml of glacial acetic acid, and 2.5 ml of 7% formaldehyde. The samples underwent the removal of Bouin’s fixative using a solution of 70% ethyl alcohol. Subsequently, the fixed cross sections of the testis were embedded in paraffin wax using established methodologies. The paraffin blocks were stained into 5 m thick sections (PAS) using hematoxylin and eosin (H&E), Masson’s trichrome, and periodic acidSchiff staining. Selection of human sertoli cells The dissociated tubules from cryopreserved testicular tissue were enzymatically digested for 30 min at 37 °C with 750 U/mL collagenase Type IV (Sigma, St Louis, MO, USA), 0.25 g/mL dispase II (Roche, Ludwigsburg, Germany), and 5 g/mL DNase in HBSS buffer with Ca^2+ and Mg^2+ (PAA) to create the digestion was halted using 10% ES cell-qualified FBS. The cell suspension was centrifuged for 15 min at 1000 rpm after passing through a 100-m cell strainer^[74]36. Our prior study isolated and identified human Sertoli cells. After three washes, testicular tissues were sliced into 0.2 cm pieces and treated with Enzyme I (10 mL of DMEM containing 2 mg/mL type IV collagenase and 10 mg/mL DNase I) for 15 min at 34 °C. According to our recent studies, these cells were chosen for microarray investigations. Highly pure Vimentin+, Sox9+, Gata4+, Vasa−, and UTF1-human Sertoli cell cultures were used in our investigations^[75]37–[76]42. Labelling probes, RNA extraction, preparing RNA and hybridization of cDNA To extract total RNA for Western blotting, the RNeasy Mini Kit (Qiagen, Germantown, MD, USA) was used. Amplification using the Ambion MessageAmp aRNA Kit followed. The University Clinic in Tübingen, Germany received samples for microarray analysis. The Human U133 + 2.0 Genome oligonucleotide array (Affymetrix, Santa Clara, CA, USA) was used to study gene expression. Raw data was normalized and biostatistically analyzed by Berlin-based MicroDiscovery GmbH. Total RNA was extracted from lectin-selected Sertoli cells for microarray analysis. RNA was extracted using the RNeasy Mini Kit (Qiagen) and amplified using the MessageAmp aRNA Kit (Ambion, Austin, TX, USA). The micromanipulation apparatus also captured 100 cells per probe every sample. These cells were quickly placed in a 10 L RNA direct lysis solution and kept at − 80 °C. RNA was collected from cell samples, and microarray datasets were produced by hybridizing human mRNA-derived cDNA after Super-AmpTM amplification. Bioinformatics was used to analyze Agilent Whole Human Genome Oligo Microarrays 8 × 60 K v2, made by Miltenyi in Germany^[77]13,[78]43,[79]44. Trizol RNA isolation reagent (Gibco BRL, Grand Island, USA) was used to extract sperm RNA following manufacturer’s instructions. Ejaculated sperm from 10 healthy adults was compared to normal and affected. Gel electrophoresis and UV spectrophotometry measured RNA quality and amount. A reverse transcription solution of 50 mg total spermatozoa RNA, oligo(dT), M-MLV reverse transcriptase, and 200 mCi alpha-33PdATP generated cDNA probes. After 3 h at 68 °C with pre-hybridization solution (6SSC, 0.5 percent SDS, 5Denhardt’s, 100 mg/ml denatured salmon sperm DNA), 33Plabeled human spermatozoa cDNA probes hybridized cDNA microarrays for 20 h. After washing with wash solution (0.1SSC, 0.5 percent SDS) for 1 h at 65 degrees Celsius, microarray membranes were phosphor screened overnight. The membrane was then examined with a Fuji Photo Film FLA-3000A plate/fluorescent image analyzer. Array Gauge (Fuji Photo Film) linearly scanned each point’s radioactive intensity using 65,536 50-mm grayscale pixels. After removing background from areas lacking PCR results, the program identified positive clones with signal intensities over 10. Two spots in each DNA sample had valid signal intensities with 0.3 standard deviation. Processing of data and analysis of differential expressions The expression of 50,000 transcripts in non-obstructive azoospermia was assessed using a microarray. Subsequently, these genes were juxtaposed with regular cells. The initial configuration consisted of an adj. p-value of 0.05 and a log2 fold change (FC) greater than 1.5. DEGs were graphically described using the volcanic map and heat map. The R program was used to generate the volcano map and heat map. Constructing a PPI network and identifying hub genes A testis tissue-specific protein–protein interaction (PPI) network comprising differentially expressed genes (DEGs) was constructed using the NetworkAnalyst website ([80]https://www.networkanalyst.ca/). Following this, the data obtained from the study of the NetworkAnalyst website were imported into cytoscape software (version 3.10.0) in order to form the PPI network. The cytohubba plug-in was utilized to discover hub genes in the PPI network using four calculation methods: degree, density of maximum neighborhood component (DMNC), maximal clique centrality (MCC), and maximum neighborhood component (MNC). The Molecular Complex Detection (MCODE) plug-in in Cytoscape was employed to identify the modules with the highest level of interaction within the aforementioned PPI network. The established criteria were as follows: a maximum depth of 100, a k-core of 4, a node score cutoff of 0.2, and a degree cutoff of 4. The gene clusters that exhibited scores greater than 5 were given priority. CytoHubba was used to identify the hub genes in the network. The top 10 core genes were identified using four algorithms: MCC, EPC, MNC, and Degree. The final hub genes were obtained by crossing the results of the four algorithms. Functional enrichment analysis of shared DEGs The functional enrichment analysis of shared DEGs was conducted using the clusterProfiler R package (version 4.10.1). The Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis were employed for this purpose. The gene sets were limited to a minimum of 10 and a maximum of 4000. A P-value of less than 0.01 was deemed statistically significant. The Gene Ontology (GO) words, which enumerate genes, can be grouped according to biological process (BP), molecular function (MF), and cellular component (CC). These clusters were plotted using the “ggplot” (version 3.5.1)^[81]45 R tools. The “GOplot” (version 1.0.2)^[82]46 R tools were used to visualize the KEGG pathways. RNA sequencing data The GEO database (accession numbers [83]GSE216907^[84]47 and [85]GSE149512^[86]48) was used to gather single-cell RNA sequencing (ScRNA-seq) data from 8 NOA samples and 2 normal samples devoid of significant diseases. The ScRNA-seq data underwent processing using the limma R package^[87]49, employing the log-normalization technique. This was done by ensuring that the log2(fold change (FC)) value was greater than 1 and the p-value was less than 0.01. Cell filtering, normalization, dimensionality reduction, and clustering were performed using the Seurat R package methodology for scRNA-seq data. A canonical-correlation analysis (CCA) approach was used to conduct the IntegrateData function in order to reduce batch effects. Presumably, due to the quality control process used to the downloaded data, 4200 cells out of a total of 5600 cells met the specified criteria: minGene = 300, maxGene = 7000, and mitochondrial content < 15%. The dimensionality reduction and clustering were performed using Uniform Manifold Approximation and Projection (UMAP) with a resolution of 0.8. Cell type annotation was ultimately accomplished by the utilization of marker genes identified in prior studies. The Seurat tool was utilized to visualize the spatial transcriptome data. Unsupervised clustering and enrichment analysis The process of unsupervised clustering was performed using the ConsensusClusterPlus software^[88]50. The Kuhn-Munkres (KM) method was executed, with the measurement distance being expressed in “euclidean”. A sampling procedure was performed using 1000 bootstraps, with each bootstrap consisting of 75% of the total samples. This resulted in the identification of four prominent clusters. The identification of differential genes within these four clusters was conducted by comparing each cluster to the remaining clusters using the following criteria: a log2(FC) value greater than 1.5 and a p-value less than 0.05. Next, the differential genes were analyzed using the clusterProfiler software through GO enrichment analysis. The Benjamini–Hochberg method was utilized to adjust the p-values, with a significance level of p < 0.01. Weighted gene co-expression network analysis The R package WGCNA^[89]51 is utilized to construct gene network modules associated with diseases through the application of Weighted Gene Co-expression Network Analysis (WGCNA), which is a robust bioinformatics approach. The generation of a weighted adjacency matrix involved the calculation of Pearson’s correlation coefficient for each pair of genes. A soft threshold parameter of 8 was chosen to classify genes into distinct gene modules, with a minimum size requirement of 40 genes. For future study, the top 150 genes with the highest node centrality, determined using principal component analysis within each module, were retrieved. One of the primary goals of the WGCNA study is to ascertain the gene modules that have a strong correlation with clinical characteristics. With the objective of examining the relationship between phenotype and each gene module, the Bicor function from the WGCNA package was utilized. This analysis adopted the biweight midcorrelation (bicor) approach. Integration of single‐cell sequencing data and cell clustering and annotation The genome sequencing (scRNA) datasets were subjected to analysis using the seurat (version 5.0) and harmony package (version 1.6). The FindAllMarkers function was utilized to acquire marker genes for cell clustering. The singler package (version 1.6.2) was then automatically annotated, along with manual querying of the Human Cell Landscape (HCL) database ([90]http://bis.zju.edu.cn/HCL/), for the purpose of cell clustering^[91]48,[92]52. Constructing trajectories of hub genes in single‐cell datasets Monocle3^[93]53 (version 4.1) is a r package specifically designed for analyzing scRNA-seq data. It offers a range of functions including clustering, pseudotime, differential analysis, and more. This study use various functions, including reduce_dimension, learn_graph, and plot_genes_in_pseudotime, to construct the developmental trajectory of sperm cells and hub genes. Identification of hub miRNAs in NOA The degree of module connection can be determined by calculating the absolute value of Pearson’s correlation. Moreover, the association between clinical traits and genes can be quantified using the absolute value of Pearson’s correlation. In order to identify hub miRNAs, two web-based technologies, namely microT-CDS3 and TargetScan, were utilized to predict potential miRNAs for hub genes. A threshold was set at a score of 0.9 for microT-CDS and a context++ score of 0.1 for TargetScan. Subsequently, the miRNAs that were identified as frequent candidates in the hub module and predicted using microT-CDS and TargetScan were classified as genuine hub miRNAs. Informed consent Informed consent was obtained from all subjects and/or their legal guardian(s). The authors affirm that informed permission was gained, indicating that the participants were fully informed about the study’s objectives, potential hazards, and advantages. The execution of all procedures adhered strictly to the applicable norms and regulations. All patients have completed the process of reading and providing information, and have been given the chance to ask any questions they may have. They comprehend that their involvement is optional and that they have the freedom to exit at any moment, without providing an explanation and without incurring any expenses. They comprehend that I will get a duplicate of this permission form. They willingly consent to participate in this research. The techniques were conducted in compliance with applicable standards and regulations (See Informed Consent and Ethical Committee approved in the related file). Ethics approval and consent to participate The local ethics committees granted permission for all experiments involving human material conducted at the Amol University of Special Modern Technologies (Approval code: Ir.ausmt.rec.1402.07) as well as the University Hospitals of Tubingen and Heidelberg. All procedures were performed in accordance with the ethical standards of the responsible committee on human experimentation and with the Helsinki Declaration of 1975, as revised in 2008. Any additional information that could potentially identify individual participants was anonymized to ensure privacy and confidentiality. Results Isolation of Sertoli cells, spermatogonia and histological sections The gene expression profiles of microarrays from different cell groups were originally assessed for study in comparison to the normal cases. Supplementary [94]1 provides a concise overview of the various groups, laboratory experiments, and computational analysis. Colonies of spermatogonia were observed after digestion in the CD49f MACS and matrix fractions (collagen nonbinding, laminin binding) of enriched spermatogonia (Fig. [95]1A,B). However, these colonies were not observed in the negative selected fraction of cells or in patients without spermatogonia. The proliferation of hFibs in primary cell cultures was reduced through the use of MACS and matrix selection, resulting in their persistence in certain cell populations. Instances of testicular histology and the distribution of testicular cell types can be detected in males diagnosed with non-obstructive aortic aneurysm (NOA) (Fig. [96]1A) as well as in individuals with normal cell (Fig. [97]1B and Supplementary [98]2). Fig. 1. [99]Fig. 1 [100]Open in a new tab Image showing nonobstructive azoospermia and normal cells in a bright field and gene expression. Here, a high-resolution photograph of azoospermia is shown. (A) The presentation features a bright-field picture of normal cells, along by a scale bar measuring 100 μm. (B) The presentation features a bright-field picture of NOA cells, along by a scale bar measuring 100 μm. (C) The DEGs study yields data in both Sperm. (D) The DEGs study yields data in both Sertoli cells via volcano plot. The process of displaying genes that are expressed differently on the microarrays of Sertoli cells. (E,F) The red data points represent genes that have had a substantial increase in expression, whereas the blue data points represent genes that have experienced a major decrease in expression. The data points in blue indicate genes that do not demonstrate a statistically significant difference. The specifications of the DEG exhibit a log2FC value over 1.5 and an adj. P-value below 0.05. The examination of the heat map provides insights into the gene expression levels of metabolite interconversion enzymes in Sperm and Sertoli cells. Up and down gene expression between sperm and Sertoli cell of NOA and normal case by microarray Microarray analysis was used to evaluate a comprehensive set of 2519 genes responsible for encoding metabolite interconversion enzymes. This work has provided significant knowledge on the correlation between genes responsible for converting metabolites, such as MOXD1, ACAD10, PCYT1A, ARG1, METTL6, GPLD1, MAOA, CYP46A1, ENTPD2, CPT1C, ADC, and CYB5B, and their presence in sperm. The findings of the microarray study conducted on three human subjects with distinct NOA demonstrated that MOXD1, ACAD10, PCYT1A, ARG1, METTL6, GPLD1, MAOA, and CYP46A1 displayed increased expression, whereas the expression of ENTPD2, CPT1C, ADC, and CYB5B exhibited decreased expression in the normal case of sperm (Fig. [101]1C,E). Moreover, this research has provided significant findings on the changes in gene expression inside Sertoli cells with genes involved in metabolite interconversion enzymes. The genes present in Sertoli cells include NAA25, MAN2A1, CYB561, PNPLA5, RRM2, MTHFD1L, FTO, L2HGDH, GALNTL2, OXCT1, GSTT2, HSD17B7, PDPR, SESN1, ESCO2, SYNJ2, EBPL, DHFR, SORD, CES1, RPIA, PIK3C3, LYPLA2, CA11, MBOAT7, and HDHD2. The microarray analysis of three human cases with different NOA showed that the expression of RPIA, PIK3C3, LYPLA2, CA11, MBOAT7, and HDHD2 was increased, while the expression of NAA25, MAN2A1, CYB561, PNPLA5, RRM2, MTHFD1L, FTO, L2HGDH, GALNTL2, OXCT1, GSTT2, HSD17B7, PDPR, SESN1, ESCO2, SYNJ2, EBPL, DHFR, SORD, and CES1 was decreased in the normal case in Sertoli cell (Fig. [102]1D,F) (Supplementary [103]3). Identification of gene clusters and hub genes The STRING website received a total of 38 DEGs that met the criteria of a confidence level above 0.9. The PPI network, consisting of 186 nodes and 366 edges, was generated and shown using Cytoscape, as shown in Fig. [104]2. The size and color of the nodes were indicative of the level of neighborhood connectivity, while the magnitude of the nodes signified the degree value. The color of each edge between the two proteins serves as an indicator of the confidence level derived by STRING analysis. The “MCODE” plug-in was used to examine the modules exhibiting a significant degree of connection, with a special focus on clusters with scores above 5. Three working subnet modules were selected from the PPI network. The top-scoring MCODE cluster in Fig. [105]2A represents the cluster of gene expression variations in sperm. This cluster has 12 nodes and 156 edges, with an MCODE score of 20. The cluster in question is included throughout the whole of the STRING network. The gene expression variations observed in Sertoli cells were organized into a cluster of 26 nodes and 232 edges. The genes that were up-regulated are shown by the red node, while the genes that were down-regulated are represented by the green node. The hub genes were found by the use of CytoHubba in Cytoscape, wherein the results obtained from four algorithms (MCC, EPC, MNC, and Degree) were utilized for comparison. The data pertaining to CytoHubba are shown in Supplementary [106]3. After the calculation was finished, the gene that showed the most agreement among the four methods was chosen. Following this, six genes that exhibited down-regulation and were positioned in the top ten were determined to be the most likely hub genes. According to Fig. [107]2A, the sperm has five hub genes, namely CYP46A1, ENTPD2, CPT1C, ADC, and CYB5B. The Sertoli cells are characterized by the presence of five hub genes, namely RPIA, PIK3C3, LYPLA2, CA11, and MBOAT7 (Fig. [108]2B). Fig. 2. [109]Fig. 2 [110]Open in a new tab The PPI network was analyzed by cytoHubba and MCODE. (A) The PPI network in Sperm is encoded by hared DEGs and five cluster modules are recovered using MCODE. It consists of 12 nodes and 156 edges. (B) In Sertoli cells, the PPI network consists of 26 nodes and 232 edges. The gene degree is quantified by the magnitude of the node. Functional and signaling pathway enrichment of DEGs The functions of DEGs were investigated by uploading 260 common DEGs for enrichment analysis using GO and KEGG. Figure [111]3 illustrates that the top five enriched Gene Ontology (GO) keywords in the BP of DEGs for sperm up and down gene expression were ribonucleoprotein complex biogenesis, ribosome biogenesis, DNA replication, chromosomal segregation, and mitotic cell cycle phase transition (Fig. [112]3A1) (Supplementary [113]4 sheet 1). The top five enriched GO terms in the BP of DEG’s Sertoli cell up and down gene expression were the regulation of cytokinesis, lipid catabolic process, positive regulation of chromosome separation, mitotic DNA damage checkpoint signaling, and establishment of RNA localization (Fig. [114]3A2) (Supplementary [115]4 sheet 2). Fig. 3. [116]Fig. 3 [117]Open in a new tab GO pathway enrichment analysis of common DEGs. (A1) The top 5 enriched biological process pathways are presented in Sperm, (A2) The top 5 enriched biological process pathways are presented in Sertoli cell, (B1) the top 5 GO terms associated with the molecular process in Sperm are listed, (B2) the top 5 GO terms associated with the molecular process in Sertoli cell are listed, (C1) The top 4 enriched cellular process pathways are presented and (C2) the top 4 enriched cellular process pathways are presented. The MF of DEG’s sperm up and down gene expression displayed five enriched GO terms: ATP hydrolysis activity, helicase activity, ribonucleoprotein complex binding, DNA replication origin binding, and damaged DNA binding (Fig. [118]3B1). Similarly, the MF of DEG’s Sertoli cell up and down gene expression showed five enriched GO terms: RNA methyltransferase activity, histone binding, DNA-directed 5′–3′ RNA polymerase activity, protein kinase regulator activity, and ATP-dependent DNA damage sensor activity (Fig. [119]3B2) (Supplementary [120]4 sheet 3). In the MF of DEG’s sperm up and down gene expression, the top five enriched GO terms were condensed chromosome, kinetochore, and lysosomal membrane (Fig. [121]3C1). In the MF of DEG’s Sertoli cell up and down gene expression, the top five enriched GO terms were protein-DNA complex, DNA replication preinitiation complex, condensed nuclear chromosome, and DNA packaging complex (Fig. [122]3C2) (Supplementary [123]4 sheet 5). Pathway enrichment analysis using KEGG and Reactome The DAVID database examined the KEGG, Reactome, and WikiPathway datasets to identify the frequently altered DEGs. The DEGs in NOA sperm primarily participate in G1/S-specific transcription, pyruvate metabolism, citric acid metabolism, lipid metabolism, phospholipid metabolism, reversal of alkylation damage caused by DNA dioxygenases, formation of xylulose-5-Phosphate, folate and pterine metabolism, and fructose metabolism. Sertoli cells exhibited a predominant presence of DEGs that were involved in various biological processes, including oxidation by cytochrome P450, amino acid metabolism, acetylcholine synthesis, serotonin transporter activity, biosynthesis and regeneration of tetrahydrobiopterin, catabolism of phenylalanine, the p53 signaling pathway, and base excision repair. The data shown in Fig. [124]4 are visually represented via graphs, as well as in Supplementary [125]5. Fig. 4. [126]Fig. 4 [127]Open in a new tab Several crucial signaling pathways are associated with the increase or decrease in expression of genes. (A) The Venn diagram was constructed by integrating three signaling pathway databases, namely KEGG, Reactome, Bioplanet, and WikiPathways. Additionally, it was observed that some signaling pathways exhibited a correlation with the upregulation or downregulation of certain genes. Furthermore, the genes implicated in major signaling pathways in sperm were also identified. A Venn diagram was created by integrating three signaling pathway databases (KEGG, Reactome, BioPlanet, and WikiPathways). (B) Certain signaling pathways were found to be associated with genes that were either up-regulated or downregulated. Genes implicated in key signaling pathways were identified in Sertoli cells. Landscape of scRNA-seq of NOA in sperm In addition, we conducted a thorough examination of the up/down gene expression in NOA at the level of individual cells. The scRNA-seq data of 29 samples was acquired, consisting of 7 NOA samples and 22 normal samples. Upon doing quality control, normalization, and eliminating the batch impact, we identified a grand total of 61,622 cells that exhibited distinct colors based on their respective sample origins. Based on the findings of prior studies, we have successfully classified seven distinct cell types through the utilization of the UMAP method (Fig. [128]6). These cell types encompass Leydig cells, which express IGF1, peritubular cells, which express WFDC1, Sertoli cells, which express SOX9, germ cells, which express TNP1, endotheliocytes, which express VWF, smooth muscle cells, which express NOTCH3, and macrophages, which express CD163. The cells derived from both NOA and normal samples exhibited a uniform distribution, suggesting the elimination of the batch effect. There is a significant disparity in the distribution of Leydig cells, peritubular cells, and germ cells between the NOA and normal samples, as seen in Fig. [129]2F. The subpopulations of germ cells, namely spermatogonia (expressing SMS, ALDH1A1, KIT, and TKTL1), spermatocyte (expressing PIWIL1 and OVOL2), and spermatid (expressing PRM3, ACR, TEX29, and ACRV1), were further identified in our study due to their significant impact on the etiology of NOA. The cell type proportion plot demonstrates a significant abundance of spermatogonia, whereas the spermatocyte and spermatid proportions are relatively low in NOA. In order to verify, we used single-cell RNA (scRNA) and demonstrated that METTL6, PCYT1A, MAOA, and GPLD1 exhibit upregulation, whereas CPT1C and CYB5B are elevated in Sperm. The UMAP analysis demonstrates the expression of this gene in both spermatocytes and spermatids, providing validation for our microarray data analysis (Fig. [130]5). Figure [131]6 displays the cells and their respective cell types, such as Leydig cells, peritubular cells, Sertoli cells, spermatogonia, spermatocyte, spermatid, and macrophages. The accuracy of cell type annotation was verified by the genome browser views of recognized marker genes, including SOX9 for Sertoli cells, IGF1 for Leydig cells, SMS for spermatogonia, PIWIL1 for spermatocyte, ACR for spermatid, CD163 for macrophages, and WFDC1 for peritubular cells (Fig. [132]5A). The UMAP projections reveal the distribution of the hub gene across several cell types in NOA compared to normal individuals. Fig. 6. [133]Fig. 6 [134]Open in a new tab Expression of METTL6, PCYT1A, MAOA, GPLD1, CPT1C and CYB5B in scRNA-seq in Sperm. Gene expression of sc-Seq for the marker genes in scRNA-seq in (A) METTL6, (B) PCYT1A, (C) MAOA, (D) GPLD1, (E) CPT1C and (F) CYB5B. Fig. 5. [135]Fig. 5 [136]Open in a new tab The identification of a hub gene module unique to Sperm of NOA. (A) The consensus clustering matrix was used to divide the 29 samples into four distinct groups. The standard samples were classified into several groups. (B) UMAP plots depicting the genetic markers associated with each major cell type. A total of 29 NOA and normal samples, germ cells, and spermatogonia were subjected to UMAP plots. UMAP maps of all cells, stained according to tissue source. The distribution of cell clusters in NOA and normal samples may be determined. (C) This study examines the relative proportions of germ cell clusters in NOA and normal samples. Gene expression of sc-seq for the marker genes in scRNA-seq in METTL6, PCYT1A, MAOA, GPLD1, CPT1C and CYB5B NOA versus normal. Weighted gene co-expression network analysis For this investigation, we used WGCNA to group together genes that are strongly connected with NOA. Following the screening process, all samples were included in the analysis. In the present investigation, a soft threshold of 8 was used to establish the scale-free network (Fig. [137]7A). Subsequently, 10 co-expression modules were evaluated by combining them based on the predetermined cut-off values (Fig. [138]7B). In addition, an examination of the module correlations was conducted, leading to the identification of the most significant link between the blue module and NOA (Fig. [139]7C). In conclusion, the 552 genes identified in the MEblue module were selected as target genes based on their strong association with the prevalence of NOA. This study examines the gene expression patterns of sc-seq in the Sertoli cell using scRNA-seq (Fig. [140]8A-F). The expression levels of SESN1, HDHD2, MBOAT7, PIK3C3, MTHFD1L, and MAN2A1 NOA were compared to the normal expression levels, as shown in Fig. [141]7D. Figure [142]9 displays the genome browser views of single-cell sequencing (scSeq) and single-cell RNA sequencing (scRNA-seq) for the marker genes METTL6, PCYT1A, MAOA, GPLD1, CPT1C, and CYB5B. Fig. 7. [143]Fig. 7 [144]Open in a new tab The present study aims to identify a hub gene module that exhibits a correlation with NOA in Sertoli cells. (A) Branches of the cluster dendrogram of the most connected genes gave rise to 11 gene co-expression modules. (B) Heat map of the correlation between module eigengenes and phenotype. (C) Intergenic connectivity of genes in turquoise module. (D) Gene expression of sc-seq for the Sertoli cell in scRNA-seq in SESN1, HDHD2, MBOAT7, PIK3C3, MTHFD1L, and MAN2A1 NOA versus normal. Fig. 8. [145]Fig. 8 [146]Open in a new tab This study aims to validate hub genes in a single-cell transcriptome dataset obtained from both normal and NOA Sperm. The screen displayed many cell kinds, each represented by a distinct hue. The UMAP projections reveal the distribution of the hub gene across several cell types in NOA compared to normal individuals. Genome browser views of scSeq for the marker genes in scRNA-seq in (A) METTL6, (B) PCYT1A, (C) MAOA, (D) GPLD1, (E) CPT1C and (F) CYB5B. Fig. 9. [147]Fig. 9 [148]Open in a new tab This study investigates the integration of miRNA-mRNA interaction networks for hub genes. (A1) The blue and brown circles correspond to hub genes that have been down-regulated and up-regulated in testicle cells, respectively. (A2) Box plot miRNA involved networks for hub genes in testicular cells. (B1) The blue and brown circles correspond to hub genes that have been down-regulated and up-regulated in Sertoli cells, respectively. (B2) Box plot miRNA involved networks for hub genes in Sertoli cells. MiRNA is shown by blue triangles, indicating its linkage with hub genes. The degree’s value is denoted by the node’s size. Construction of miRNA–mRNA networks The possible miRNAs of the hub genes were predicted using the miRNet database. The hub genes’ 36 mRNAs, acquired via the PPI network, were uploaded and a total of 122 targeted miRNAs were identified. The miRNA and mRNA co-expression regulatory network was created using Cytoscape, based on the prediction findings. The network comprised of 140 nodes and 255 edges, as shown in Fig. [149]10. Based on the findings from the enrichment functional analysis of miRNA, it can be inferred that the pathway in question may be subject to regulation by hsa-mir-2964-5p, has-mir-412-5p, hsa-mir-4536, mir-704, and has-mir-187 in the context of sperm. Furthermore, hsa-mir-2964-5p and has-mir-412-5p were shown to affect two out of the six hub genes, namely MOXD1, ACAD10, PCYT1A, CYB5B, PCYT1A, and METTL6. In respect to Sertoli cells, there was a significant association seen between the genes has-mir-4764, has-mir-4636, has-mir-4637, and mir-3074 and the up and down expression of certain genes, namely CES1, SESN1, HSD17B7, GALNTL2, and EBPL. Fig. 10. [150]Fig. 10 [151]Open in a new tab The Reactome analysis reveals that the upregulation or downregulation of genes may have a role in the expression of many kinds of metabolite interconversion enzymes. Discussion Our study shows that some gene was changed in sperm and Sertoli cells that related to the converter enzyme and they are so important for spermatogenesis and metabolism. Sestrins, previously discovered as p53 target genes that diminish with age, are involved in nutritional stress response, metabolic control, and autophagy^[152]54. Sestrins are strong antioxidants that activate the Nrf2-Keap1 pathway to promote p62-dependent autophagic breakdown. Previous research showed that Sestrin2 and Sestrin3 protect spermatogenesis against heat stress^[153]55. We analysis single-cell RNA sequencing (scRNA-seq) on 29 samples (7 NOA and 22 normal) to explore differences in cell type distribution. Significant alterations were found in NOA samples, with a marked increase in spermatogonia and a decrease in spermatocytes and spermatids, indicating a blockage in spermatogenesis. These findings highlight a disrupted germ cell maturation process in NOA, providing a deeper understanding of its pathology, which involves both germ cell defects and alterations in the supporting somatic cell environment. Our findings indicate significant alterations in the expression of key metabolite interconversion enzymes in both Sertoli cells and spermatogonia from NOA patients compared to normal controls. Notably, we observed upregulation of enzymes such as MOXD1 (monooxygenase), which is involved in oxidative stress regulation and neurotransmitter metabolism, and PCYT1A (phosphatidylcholine cytidylyltransferase), which plays a role in membrane lipid biosynthesis, suggesting that alterations in membrane dynamics may influence cellular communication and spermatogenesis. The downregulation of ENTPD2 (ecto-nucleoside triphosphate diphosphohydrolase 2) in sperm and NAA25 (N-alpha-acetyltransferase) in Sertoli cells, both involved in ATP and acetyl-CoA metabolism, highlights potential disruptions in energy regulation within testicular cells, which could contribute to impaired spermatogenesis in NOA patients. The enzymes identified in this study have been linked to essential processes such as DNA repair (methyltransferase activity), redox balance (monooxygenase activity), and lipid metabolism (phosphohydrolase activity), all of which are crucial for maintaining the metabolic health of testicular cells. Man2a2 global knockout tests were smaller, had less spermatids, and generated immature sperm. Man2a2 mutant germ cells had decreased adherence to Sertoli cells, suggesting an N-glycan ending in GlcNAc is needed^[154]56. Global deletion of the other alpha mannosidase II encoded by the Man2a1 gene does not impair viability or male fertility, but whole-body knockout of both genes kills neonates. E15.5 Man2a1/Man2a2 double knockout embryo tissues had unusual hybrid N-glycans, whereas germ cells missing MAN2A2 did not^[155]57. During prenatal life, Man1a2 gene disruption, which encodes an alpha mannosidase that prunes oligomannose N-glycans, is fatal^[156]58. The methyltransferase complex includes METTL5, tRNA methyltransferase activator subunit 11-2 (TRMT112), Cbl proto-oncogene like 1 (CBLL1), zinc finger protein 217 (ZFP217), RNA binding motif protein 15 (RBM15), RBM15B, and vir-like m6A methyltransferases^[157]59,[158]60. Dysfunction or expression of these writers disrupted m6A regulation, mRNA metabolism, and gene expression. METTL3, an S-adenosylmethionine-dependent methyltransferase, detects the conserved RRACH motif on RNAs^[159]61. METTL14 and METTL3 create a stable heterodimer that boosts catalytic activity and RNA methylation site recognition^[160]62. Through sperm count and motility variation, neurotransmission gene polymorphisms HTR2A, MAOA, and SLC18A1 are separately linked with asthenozoospermia without allelic or genotype interaction^[161]63. Cell membrane-anchored proteins’ inositol phosphate connections are hydrolyzed by GPLD1. GPLD1-overexpressing mice had improved cognition and neurogenesis. Age-dependent cognitive decline is delayed in Snell dwarf (DW) and GKO mice^[162]64. A complex network of endocrine, paracrine, and autocrine signals controls vertebrate spermatogenesis. They may be separated into general transcription factors, nuclear receptor superfamily, testicular transcription factors, testis-specific gene transcription, and cell junction dynamics transcriptional regulators^[163]65. Recent studies show that in response to hypothalamic gonadotropin hormone release (GnRH), the pituitary gland releases LH and FSH, which regulate spermatogenesis and other transcription factors. LH controls testosterone release by somatic Leydig cells in the interstitium between seminiferous tubules; FSH stimulates Sertoli cell signaling, gene expression, and peptide and signaling molecule secretion^[164]66. cAMP Response Element-Binding Protein (CREB) transcription factor, a key FSH signal transducer in Sertoli cells. ROS overexpression activates many signaling pathways that regulate gene expression via CREB transcription factors^[165]67. Alternatively, the spliced variant CREB mRNAs in the rat testis are spermatogenic, cycle dependent, and expressed throughout germ and Sertoli cell development, suggesting that the CREB isoforms may be key participants in spermatogenesis. The transcription factor CREM is widely expressed in male germ cells and controls post-meiotic genes including transition proteins and protamines, making it the probable gene expression regulator during spermatogenesis^[166]68,[167]69. Targeted CREM gene disruption stops the initial phase of spermiogenesis differentiation. These results suggest that hormonal cues directly regulate spermatogenesis genes via CREM in post-meiotic germ cell differentiation^[168]70–[169]72. It is now apparent that germ cells and testicular somatic cells (Leydig and Sertoli cells) contain P450arom mRNA, which is converted into an estrogen-producing enzyme. Thus, germ cells have a modified androgen/estrogen ratio, and while testosterone regulates testicular functions, estrogens also regulate gonadotropin secretion, Leydig cell development, steroidogenesis, and spermatogenesis and spermiogenesis in some mammals^[170]73. Human cells generate energy by oxidizing diet molecules and synthesizing essential molecules not in the diet and inactivating and eliminating toxic ones from the endogenous or extracellular environment^[171]74. Synthetic processes produce nucleotides, amino acids, and porphyrins, which are easily categorized by chemical nature. Chemical toxin properties also classify detoxification processes (biological oxidations). All these processes are intimately connected. Energy-generating intermediates are starting materials for amino acid and other compound biosynthesis, broad-specificity oxidoreductase enzymes can be used in detoxification and biosynthesis, and hormone-mediated signaling processes coordinate energy-generating and energy-storing reactions and couple them to other biosynthetic processes. Microarray and bioinformatics show that change in metabolite interconversion enzymes, including methyltransferase, monooxygenase, transmembrane reductase, and phosphohydrolase can cause NOA and disrupt spermatogenesis metabolism (Fig. [172]10). Conclusion Finally, the expression of MOXD1, ACAD10, PCYT1A, ARG1, METTL6, GPLD1, MAOA, and CYP46A1 was upregulated, while expression of MOXD1, ACAD10, PCYT1A, ARG1, METTL6, GPLD1, MAOA, and CYP46A1 was downregulated in sperm. There was an increase in RPIA, PIK3C3, LYPLA2, CA11, MBOAT7, and HDHD2 expression in the Sertoli cells of three human cases with NOA, whereas expression of NAA25, MAN2A1, CYB561, PNPLA5, RRM2, MTHFD1L, FTO, L2HGDH, GALNTL2, OXCT1, GSTT2, HSD17B7, PDPR, SESN1, ESCO2, SYNJ2, EBPL, DHFR, SORD, and CES1 were downregulated. In our study, we suggested that protein interactions, pathway signaling, single cell sequencing and miRNAs can be affected by these genes. Collectively, these findings provide a new mechanism for NOA and new therapeutic targets for NOA patients. Electronic Supplementary Material Below is the link to the electronic supplementary material. [173]Supplementary Material 1^ (16.3KB, xlsx) [174]Supplementary Material 2^ (30.2KB, xlsx) [175]Supplementary Material 3^ (264.2KB, pdf) [176]Supplementary Material 4^ (221.9KB, xlsx) [177]Supplementary Material 5^ (1MB, xlsx) Acknowledgements