Metastatic triple negative breast cancer adapts its metabolism to destination tissues while retaining key metabolic signatures.

Triple negative breast cancer (TNBC) metastases are assumed to exhibit similar functions in different organs as in the original primary tumor. However, studies of metastasis are often limited to a comparison of metastatic tumors with primary tumors of their origin, and little is known about the adaptation to the local environment of the metastatic sites. We therefore used transcriptomic data and metabolic network analyses to investigate whether metastatic tumors adapt their metabolism to the metastatic site and found that metastatic tumors adopt a metabolic signature with some similarity to primary tumors of their destinations. The extent of adaptation, however, varies across different organs, and metastatic tumors retain metabolic signatures associated with TNBC. Our findings suggest that a combination of anti-metastatic approaches and metabolic inhibitors selected specifically for different metastatic sites, rather than solely targeting TNBC primary tumors, may constitute a more effective treatment approach.

A reconstructed genome-scale metabolic model of Helicobacter pylori for predicting putative drug targets in clarithromycin and rifampicin resistance conditions.

BACKGROUND: Helicobacter pylori is considered a true human pathogen for which rising drug resistance constitutes a drastic concern globally. The present study aimed to reconstruct a genome-scale metabolic model (GSMM) to decipher the metabolic capability of H. pylori strains in response to clarithromycin and rifampicin along with identification of novel drug targets. MATERIALS AND METHODS: The iIT341 model of H. pylori was updated based on genome annotation data, and biochemical knowledge from literature and databases. Context-specific models were generated by integrating the transcriptomic data of clarithromycin and rifampicin resistance into the model. Flux balance analysis was employed for identifying essential genes in each strain, which were further prioritized upon being nonhomologs to humans, virulence factor analysis, druggability, and broad-spectrum analysis. Additionally, metabolic differences between sensitive and resistant strains were also investigated based on flux variability analysis and pathway enrichment analysis of transcriptomic data. RESULTS: The reconstructed GSMM was named as HpM485 model. Pathway enrichment and flux variability analyses demonstrated reduced activity in the ribosomal pathway in both clarithromycin- and rifampicin-resistant strains. Also, a significant decrease was detected in the activity of metabolic pathways of clarithromycin-resistant strain. Moreover, 23 and 16 essential genes were exclusively detected in clarithromycin- and rifampicin-resistant strains, respectively. Based on prioritization analysis, cyclopropane fatty acid synthase and phosphoenolpyruvate synthase were identified as putative drug targets in clarithromycin- and rifampicin-resistant strains, respectively. CONCLUSIONS: We present a robust and reliable metabolic model of H. pylori. This model can predict novel drug targets to combat drug resistance and explore the metabolic capability of H. pylori in various conditions.

Identification of Catechins' Binding Sites in Monomeric Aß42 through Ensemble Docking and MD Simulations.

The assembly of the amyloid-ß peptide (Aß) into toxic oligomers and fibrils is associated with Alzheimer's disease and dementia. Therefore, disrupting amyloid assembly by direct targeting of the Aß monomeric form with small molecules or antibodies is a promising therapeutic strategy. However, given the dynamic nature of Aß, standard computational tools cannot be easily applied for high-throughput structure-based virtual screening in drug discovery projects. In the current study, we propose a computational pipeline-in the framework of the ensemble docking strategy-to identify catechins' binding sites in monomeric Aß42. It is shown that both hydrophobic aromatic interactions and hydrogen bonding are crucial for the binding of catechins to Aß42. Additionally, it has been found that all the studied ligands, especially EGCG, can act as potent inhibitors against amyloid aggregation by blocking the central hydrophobic region of Aß. Our findings are evaluated and confirmed with multi-microsecond MD simulations. Finally, it is suggested that our proposed pipeline, with low computational cost in comparison with MD simulations, is a suitable approach for the virtual screening of ligand libraries against Aß.

Structural insights into the substrate-binding site of main protease for the structure-based COVID-19 drug discovery.

An attractive drug target to combat COVID-19 is the main protease (Mpro ) because of its key role in the viral life cycle by processing the polyproteins translated from the viral RNA. Studying the crystal structures of the protease is important to enhance our understanding of its mechanism of action at the atomic-level resolution, and consequently may provide crucial structural insights for structure-based drug discovery. In the current study, we report a comparative structural analysis of the Mpro substrate binding site for both apo and holo forms to identify key interacting residues and conserved water molecules during the ligand-binding process. It is shown that in addition to the catalytic dyad residues (His41 and Cys145), the oxyanion hole residues (Asn142-Ser144) and residues His164-Glu166 form essential parts of the substrate-binding pocket of the protease in the binding process. Furthermore, we address the issue of the substrate-binding pocket flexibility and show that two adjacent loops in the Mpro structures (residues Thr45-Met49 and Arg188-Ala191) with high flexibility can regulate the binding cavity' accessibility for different ligand sizes. Moreover, we discuss in detail the various structural and functional roles of several important conserved and mobile water molecules within and around the binding site in the proper enzymatic function of Mpro . We also present a new docking protocol in the framework of the ensemble docking strategy. The performance of the docking protocol has been evaluated in predicting ligand binding pose and affinity ranking for two popular docking programs; AutoDock4 and AutoDock Vina. Our docking results suggest that the top-ranked poses of the most populated clusters obtained by AutoDock Vina are the most important representative docking runs that show a very good performance in estimating experimental binding poses and affinity ranking.

Electrostatically induced pKa shifts in oligopeptides: the upshot of neighboring side chains.

pKa values of homorepeat hexapeptides with a 2,3-diazabicyclo[2.2.2]oct-2-ene (DBO) chromophore attached at the peptide C termini, through an asparagine derivative (Dbo), namely His6-Dbo (H6), Lys6-Dbo (K6), and Arg6-Dbo (R6), were determined by a novel fluorescence-based method. The fluorescence lifetime of Dbo in the peptides (τ) was measured as a function of pH. The side chains collide with Dbo intramolecularly and quench it efficiently only when they are deprotonated (i.e., pH ≥ side chain pKa). The pKa values of the H6, K6, and R6 peptides, attributable to side chain ionization, were found to be depressed compared to the pKa values of the His, Lys, and Arg residues in their free amino acid forms. We further looked into the structural changes of the peptides by molecular dynamics (MD) simulations; the peptides were structurally more expanded when their side chains are protonated. The structural expansion of the peptides reflects an electrostatic repulsion between the protonated side chain residues, which also accounts for the observed decrease in pKa values, which corresponds to a facilitated deprotonation, assisted by electrostatic repulsion.

Prediction of antimicrobial peptides toxicity based on their physico-chemical properties using machine learning techniques.

BACKGROUND: Antimicrobial peptides are promising tools to fight against ever-growing antibiotic resistance. However, despite many advantages, their toxicity to mammalian cells is a critical obstacle in clinical application and needs to be addressed. RESULTS: In this study, by using an up-to-date dataset, a machine learning model has been trained successfully to predict the toxicity of antimicrobial peptides. The comprehensive set of features of both physico-chemical and linguistic-based with local and global essences have undergone feature selection to identify key properties behind toxicity of antimicrobial peptides. After feature selection, the hybrid model showed the best performance with a recall of 0. 876 and a F1 score of 0. 849. CONCLUSIONS: The obtained model can be useful in extracting AMPs with low toxicity from AMP libraries in clinical applications. On the other hand, several properties with local nature including positions of strand forming and hydrophobic residues in final selected features show that these properties are critical definer of peptide properties and should be considered in developing models for activity prediction of peptides. The executable code is available at https://git.io/JRZaT .

Identification of a missense variant in CLDN2 in obstructive azoospermia.

Obstructive azoospermia (OA), defined as an obstruction in any region of the male genital tract, accounts for 40% of all azoospermia cases. Of all OA cases, ~30% are thought to have a genetic origin, however, hitherto, the underlying genetic etiology of the majority of these cases remain unknown. To address this, we took a family-based whole-exome sequencing approach to identify causal variants of OA in a multiplex family with epidydimal obstruction. A novel gain-of-function missense variant in CLDN2 (c.481G>C; p.Gly161Arg) was found to co-segregate with the phenotype, consistent with the X-linked inheritance pattern observed in the pedigree. To assess the pathogenicity of this variant, the wild and mutant protein structures were modeled and their potential for strand formation in multimeric form was assessed and compared. The results showed that dimeric and tetrameric arrangements of Claudin-2 were not only reduced, but were also significantly altered by this single residue change. We, therefore, envisage that this amino acid change likely forms a polymeric discontinuous strand, which may lead to the disruption of tight junctions among epithelial cells. This missense variant is thus likely to be responsible for the disruption of the blood-epididymis barrier, causing dislodged epithelial cells to clog the genital tract, hence causing OA. This study not only sheds light on the underlying pathobiology of OA, but also provides a basis for more efficient diagnosis in the clinical setting.

Hersintuzumab: A novel humanized anti-HER2 monoclonal antibody induces potent tumor growth inhibition.

Humanized monoclonal antibodies (mAbs) against HER2 including trastuzumab and pertuzumab are widely used to treat HER2 overexpressing metastatic breast cancers. These two mAbs recognize distinct epitopes on HER2 and their combination induces a more potent blockade of HER2 signaling than trastuzumab alone. Recently, we have reported characterization of a new chimeric mAb (c-1T0) which binds to an epitope different from that recognized by trastuzumab and significantly inhibits proliferation of HER2 overexpressing tumor cells. Here, we describe humanization of this mAb by grafting all six complementarity determining regions (CDRs) onto human variable germline genes. Humanized VH and VL sequences were synthesized and ligated to human γ1 and κ constant region genes using splice overlap extension (SOE) PCR. Subsequently, the humanized antibody designated hersintuzumab was expressed and characterized by ELISA, Western blot and flow cytometry. The purified humanized mAb binds to recombinant HER2 and HER2-overexpressing tumor cells with an affinity comparable with the chimeric and parental mouse mAbs. It recognizes an epitope distinct from those recognized by trastuzumab and pertuzumab. Binding of hersintuzumab to HER2 overexpressing tumor cells induces G1 cell cycle arrest, inhibition of ERK and AKT signaling pathways and growth inhibition. Moreover, hersintuzumab could induce antibody-dependent cell cytotoxicity (ADCC) on BT-474 cells. This new humanized mAb is a potentially valuable tool for single or combination breast cancer therapy.

Red/ox states of human protein disulfide isomerase regulate binding affinity of 17 beta-estradiol.

Human protein disulfide isomerase (hPDI) is a key redox-regulated thiol-containing protein operating as both oxidoreductase and molecular chaperone in the endoplasmic reticulum of cells. hPDI thiol-disulfide interchange reactions lead to the adoption of two distinct red/ox conformations with different substrate preferences. hPDI also displays high binding capacity for some endogenous steroid hormones including 17 beta-estradiol (E2) and thus contributes to the regulation of their intracellular concentration, storage and actions. The primary focus of this study was to investigate the impact of E2 binding on functional activity of recombinant hPDI. Then, we examined the effect of E2 binding on structural alteration of hPDI red/ox conformations and its influence on affinity and position of interaction using experimental and computational analysis. Our results revealed that interaction of one E2 per each hPDI molecule led to the inhibition of hPDI reductase activity and conformational changes in both oxidation states. Mutually, E2-binding position were also redox-regulated with higher affinity in oxidized hPDI compare to the reduced form. The importance of histidine-256 protonation states in distinct binding preferences of E2 were also demonstrated in hPDI red/ox conformations. These findings might pave the way for better understanding of the mechanisms behind the redox-dependent hormone-binding activity of hPDI.

Histidine substitution in the most flexible fragments of firefly luciferase modifies its thermal stability.

Molecular dynamics (MD) at two temperatures of 300 and 340 K identified two histidine residues, His461 and His489, in the most flexible regions of firefly luciferase, a light emitting enzyme. We therefore designed four protein mutants H461D, H489K, H489D and H489M to investigate their enzyme kinetic and thermodynamic stability changes. Substitution of His461 by aspartate (H461D) decreased ATP binding affinity, reduced the melting temperature of protein by around 25 °C and shifted its optimum temperature of activity to 10 °C. In line with the common feature of psychrophilic enzymes, the MD data showed that the overall flexibility of H461D was relatively high at low temperature, probably due to a decrease in the number of salt bridges around the mutation site. On the other hand, substitution of His489 by aspartate (H489D) introduced a new salt bridge between the C-terminal and N-terminal domains and increased protein rigidity but only slightly improved its thermal stability. Similar changes were observed for H489K and, to a lesser degree, H489M mutations. Based on our results we conclude that the MD simulation-based rational substitution of histidines by salt-bridge forming residues can modulate conformational dynamics in luciferase and shift its optimal temperature activity.

Insights into the molecular interaction between two polyoxygenated cinnamoylcoumarin derivatives and human serum albumin.

Ligand binding studies on human serum albumin (HSA) are crucial in determining the pharmacological properties of drug candidates. Here, two representatives of coumarin-chalcone hybrids were selected and their binding mechanism was identified via thermodynamics techniques, curve resolution analysis and computational methods at molecular levels. The binding parameters were derived using spectroscopic approaches and the results point to only one pocket located near the Trp214 residue in subdomain IIA of HSA. The protein tertiary structure was altered during ligand binding and formed an intermediate structure to create stronger ligand binding interactions. The best binding mode of the ligand was initially estimated by docking on an ensemble of HSA crystallographic structures and by molecular dynamics (MD) simulations. Per residue interaction energies were calculated over the MD trajectories as well. Reasonable agreement was found between experimental and theoretical results about the nature of binding, which was dominated by hydrogen bonding and van der Waals contributions.

Improving protein-ligand docking results using the Semiempirical quantum mechanics: testing on the PDBbind 2016 core set.

Molecular docking techniques are routinely employed for predicting ligand binding conformations and affinities in the in silico phase of the drug design and development process. In this study, a reliable semiempirical quantum mechanics (SQM) method, PM7, was employed for geometry optimization of top-ranked poses obtained from two widely used docking programs, AutoDock4 and AutoDock Vina. The PDBbind core set (version 2016), which contains high-quality crystal protein - ligand complexes with their corresponding experimental binding affinities, was used as an initial dataset in this research. It was shown that docking pose optimization improves the accuracy of pose predictions and is very useful for the refinement of docked complexes via removing clashes between ligands and proteins. It was also demonstrated that AutoDock Vina achieves a higher sampling power than AutoDock4 in generating accurate ligand poses (RMSD ≤ 2.0 Å), while AutoDock4 exhibits a better ranking power than AutoDock Vina. Finally, a new protocol based on a combination of the results obtained from the two docking programs was proposed for structure-based virtual screening studies, which benefits from the robust sampling abilities of AutoDock Vina and the reliable ranking performance of AutoDock4.Communicated by Ramaswamy H. Sarma.

Plasma acylcarnitines and amino acids in dyslipidemia: An integrated metabolomics and machine learning approach.

Purpose: The Discovery of underlying intermediates associated with the development of dyslipidemia results in a better understanding of pathophysiology of dyslipidemia and their modification will be a promising preventive and therapeutic strategy for the management of dyslipidemia. Methods: The entire dataset was selected from the Surveillance of Risk Factors of Noncommunicable Diseases (NCDs) in 30 provinces of Iran (STEPs 2016 Country report in Iran) that included 1200 subjects and was stratified into four binary classes with normal and abnormal cases based on their levels of triglyceride (TG), total cholesterol (TC), high-density lipoprotein cholesterol (HDL-C), and non-HDL-C.Plasma concentrations of 20 amino acids and 30 acylcarnitines in each class of dyslipidemia were evaluated using Tandem mass spectrometry. Then, these attributes, along with baseline characteristics data, were used to check whether machine learning (ML) algorithms could classify cases and controls. Results: Our ML framework accurately predicts TG binary classes. Among the models tested, the SVM model stood out, performing slightly better with an AUC of 0.81 and a standard deviation of test accuracy at 0.04. Consequently, it was chosen as the optimal model for TG classification. Moreover, the findings showed that alanine, phenylalanine, methionine, C3, C14:2, and C16 had great power in differentiating patients with high TG from normal TG controls. Conclusions: The comprehensive output of this work, along with sex-specific attributes, will improve our understanding of the underlying intermediates involved in dyslipidemia. Supplementary Information: The online version contains supplementary material available at 10.1007/s40200-024-01384-9.

Intermediate-aided allostery mechanism for α-glucosidase by Xanthene-11v as an inhibitor using residue interaction network analysis.

Exploring allosteric inhibition and the discovery of new inhibitor binding sites are important studies in protein regulation mechanisms and drug discovery. Structural and network-based analyses of trajectories resulting from molecular dynamics (MD) simulations have been developed to discover protein dynamics, landscape, functions, and allosteric regions. Here, an experimentally suggested non-competitive inhibitor, xanthene-11v, was considered to explore its allosteric inhibition mechanism in α-glucosidase MAL12. Comparative structural and network analyses were applied to eight 250 ns independent MD simulations, four of which were performed in the free state and four of which were performed in ligand-bound forms. Projected two-dimensional free energy landscapes (FEL) were constructed from the probabilistic distribution of conformations along the first two principal components. The post-simulation analyses of the coordinates, side-chain torsion angles, non-covalent interaction networks, network communities, and their centralities were performed on α-glucosidase conformations and the intermediate sub-states. Important communities of residues have been found that connect the allosteric site to the active site. Some of these residues like Thr307, Arg312, TYR344, ILE345, Phe357, Asp406, Val407, Asp408, and Leu436 are the key messengers in the transition pathway between allosteric and active sites. Evaluating the probability distribution of distances between gate residues including Val407 in one community and Phe158, and Pro65 in another community depicted the closure of this gate due to the inhibitor binding. Six macro states of protein were deduced from the topology of FEL and analysis of conformational preference of free and ligand-bound systems to these macro states shows a combination of lock-and-key, conformational selection, and induced fit mechanisms are effective in ligand binding. All these results reveal structural states, allosteric mechanisms, and key players in the inhibition pathway of α-glucosidase by xanthene-11v.

Design of a dual-function agent by fusing a designed anti-VEGF-A binder and CPG-2 enzyme.

Anti-VEGF therapies are common for the treatment of cancer. Carboxypeptidase G (CPG-2) enzyme is a zinc-dependent metalloenzyme that metabolizes non-toxic synthetic 'benzoic mustard prodrugs' to cytotoxic moieties in tumor cells. In this study, we designed a dual-activity agent by combining a designed anti-VEGF- and CPG-2 enzyme to convert methotrexate (MTX). VEGF-A was docked against a set of scaffolds, and suitable inverse rotamers were made. Rosetta design was used for the interface design. The top 1200 binders were chosen by flow cytometry and displayed in yeast. The activity of CPG-2 enzyme was analyzed at different temperature conditions and in the presence of the substrate, MTX. Optimal binders were selected and protein was eluted using immobilized metal affinity chromatography and size-exclusion chromatography. Both, native PAGE and on-yeast flow cytometry confirmed the binding of the binder to VEGF-A. The activity of truncated enzymes was slightly lower than that of full-length enzymes linked to VEGF-A. The method should be generally useful as a dual-activity agent for targeting VEGF-A and combination therapy with the enzyme CPG-2 for metabolizing non-toxic prodrugs to cytotoxic moieties.Communicated by Ramaswamy H. Sarma.

PANDORA v2.0: Benchmarking peptide-MHC II models and software improvements.

T-cell specificity to differentiate between self and non-self relies on T-cell receptor (TCR) recognition of peptides presented by the Major Histocompatibility Complex (MHC). Investigations into the three-dimensional (3D) structures of peptide:MHC (pMHC) complexes have provided valuable insights of MHC functions. Given the limited availability of experimental pMHC structures and considerable diversity of peptides and MHC alleles, it calls for the development of efficient and reliable computational approaches for modeling pMHC structures. Here we present an update of PANDORA and the systematic evaluation of its performance in modelling 3D structures of pMHC class II complexes (pMHC-II), which play a key role in the cancer immune response. PANDORA is a modelling software that can build low-energy models in a few minutes by restraining peptide residues inside the MHC-II binding groove. We benchmarked PANDORA on 136 experimentally determined pMHC-II structures covering 44 unique αß chain pairs. Our pipeline achieves a median backbone Ligand-Root Mean Squared Deviation (L-RMSD) of 0.42 Å on the binding core and 0.88 Å on the whole peptide for the benchmark dataset. We incorporated software improvements to make PANDORA a pan-allele framework and improved the user interface and software quality. Its computational efficiency allows enriching the wealth of pMHC binding affinity and mass spectrometry data with 3D models. These models can be used as a starting point for molecular dynamics simulations or structure-boosted deep learning algorithms to identify MHC-binding peptides. PANDORA is available as a Python package through Conda or as a source installation at https://github.com/X-lab-3D/PANDORA.

Ensemble learning from ensemble docking: revisiting the optimum ensemble size problem.

Despite considerable advances obtained by applying machine learning approaches in protein-ligand affinity predictions, the incorporation of receptor flexibility has remained an important bottleneck. While ensemble docking has been used widely as a solution to this problem, the optimum choice of receptor conformations is still an open question considering the issues related to the computational cost and false positive pose predictions. Here, a combination of ensemble learning and ensemble docking is suggested to rank different conformations of the target protein in light of their importance for the final accuracy of the model. Available X-ray structures of cyclin-dependent kinase 2 (CDK2) in complex with different ligands are used as an initial receptor ensemble, and its redundancy is removed through a graph-based redundancy removal, which is shown to be more efficient and less subjective than clustering-based representative selection methods. A set of ligands with available experimental affinity are docked to this nonredundant receptor ensemble, and the energetic features of the best scored poses are used in an ensemble learning procedure based on the random forest method. The importance of receptors is obtained through feature selection measures, and it is shown that a few of the most important conformations are sufficient to reach 1 kcal/mol accuracy in affinity prediction with considerable improvement of the early enrichment power of the models compared to the different ensemble docking without learning strategies. A clear strategy has been provided in which machine learning selects the most important experimental conformers of the receptor among a large set of protein-ligand complexes while simultaneously maintaining the final accuracy of affinity predictions at the highest level possible for available data. Our results could be informative for future attempts to design receptor-specific docking-rescoring strategies.

Binder design for targeting SARS-CoV-2 spike protein: An in silico perspective.

INTRODUCTION: The COVID-19 pandemic is now affecting all people around the world and getting worse. New antiviral medications are desperately needed other than the few approved medications that have shown no promising efficacy so far. METHODS: Here we report three blocking binders for targeting SARS-CoV-2 spike protein to block the interaction between the spike protein on the SARS-CoV-2 and the angiotensin-converting enzyme 2 (ACE2) receptors, responsible for viral homing into the alveolar epithelium type II cells (AECII). RESULTS: The design process is based on the collected natural scaffolds and using Rosetta interface for designing the binders. CONCLUSION: Based on the structural analysis, three binders were selected, and the results showed that they might be promising as new therapeutic targets for blocking COVID-19.

ET-score: Improving Protein-ligand Binding Affinity Prediction Based on Distance-weighted Interatomic Contact Features Using Extremely Randomized Trees Algorithm.

The molecular docking simulation is a key computational tool in modern drug discovery research that its predictive performance strongly depends on the employed scoring functions. Many recent studies have shown that the application of machine learning algorithms in the development of scoring functions has led to a significant improvement in docking performance. In this work, we introduce a new machine learning (ML) based scoring function called ET-Score, which employs the distance-weighted interatomic contacts between atom type pairs of the ligand and the protein for featurizing protein-ligand complexes and Extremely Randomized Trees algorithm for the training process. The performance of ET-Score is compared with some successful ML-based scoring functions and several popular classical scoring functions on the PDBbind 2016v core set. It is shown that our ET-Score model (with Pearson's correlation of 0.827 and RMSE of 1.332) achieves very good performance in comparison with most of the ML-based scoring functions and all classical scoring functions despite its extremely low computational cost. ET-Score's codes are freely available on the web at https://github.com/miladrayka/ET_Score.

Integration and gene co-expression network analysis of scRNA-seq transcriptomes reveal heterogeneity and key functional genes in human spermatogenesis.

Spermatogenesis is a complex process of cellular division and differentiation that begins with spermatogonia stem cells and leads to functional spermatozoa production. However, many of the molecular mechanisms underlying this process remain unclear. Single-cell RNA sequencing (scRNA-seq) is used to sequence the entire transcriptome at the single-cell level to assess cell-to-cell variability. In this study, more than 33,000 testicular cells from different scRNA-seq datasets with normal spermatogenesis were integrated to identify single-cell heterogeneity on a more comprehensive scale. Clustering, cell type assignments, differential expressed genes and pseudotime analysis characterized 5 spermatogonia, 4 spermatocyte, and 4 spermatid cell types during the spermatogenesis process. The UTF1 and ID4 genes were introduced as the most specific markers that can differentiate two undifferentiated spermatogonia stem cell sub-cellules. The C7orf61 and TNP can differentiate two round spermatid sub-cellules. The topological analysis of the weighted gene co-expression network along with the integrated scRNA-seq data revealed some bridge genes between spermatogenesis's main stages such as DNAJC5B, C1orf194, HSP90AB1, BST2, EEF1A1, CRISP2, PTMS, NFKBIA, CDKN3, and HLA-DRA. The importance of these key genes is confirmed by their role in male infertility in previous studies. It can be stated that, this integrated scRNA-seq of spermatogenic cells offers novel insights into cell-to-cell heterogeneity and suggests a list of key players with a pivotal role in male infertility from the fertile spermatogenesis datasets. These key functional genes can be introduced as candidates for filtering and prioritizing genotype-to-phenotype association in male infertility.