Data preprocessing and identification of DEGs
The gene expression profile with accession numbers GSE113079 was downloaded from GEO database. The results of before and after normalizing the microarray gene expression are shown in Fig. 1A and Fig. 1B. DEGs between peripheral blood mononuclear cells from CAD patients and peripheral blood mononuclear cells from healthy control were screened using limma package in R bioconductor (P-value <0.05, |logFC| > 0.97 for up regulated genes, and |logFC| < - 0.963 for down regulated genes). In this study, 1,036 total DEGs (525 up regulated genes and 511 down regulated genes, respectively) in GSE113079 was screened. The total number of DEGs collected for each subject in the differential gene expression analysis is given in Table 1. A volcano diagram was constructed for the DEGs and is shown in Fig. 2. The DEGs (up and down regulated genes) are presented by a cluster heatmap in Fig. 3A and Fig. 3B.
Pathway enrichment analysesof DEGs
Pathway enrichment analyses were performed using ToppGene, analyzing the pathway classification of DEGs (up and down regulated genes). Pathways of up regulated were mainly enriched in thyronamine and iodothyronamine metabolism, trehalose degradation, cytokine-cytokine receptor interaction, Olfactory transduction, FOXA1 transcription factor network, lissencephaly gene (LIS1) in neuronal migration and development, signaling by GPCR, GPCR downstream signaling, C21 seroid hormone metabolism, androgen and estrogen metabolism, ensemble of genes encoding extracellular matrix and extracellular matrix-associated proteins, ensemble of genes encoding ECM-associated proteins including ECM-affilaited proteins, ECM regulators and secreted factors, cortocotropin releasing factor receptor signalingpathway, 5HT4 type receptor mediated signaling pathway, corticotropin-releasing hormone signaling, G protein signaling via galphaq family, ornithine transcarbamylase deficiency (OTC deficiency and intracellular signaling through PGD2 receptor and prostaglandin D2 according to the BIOCYC, KEGG, PID, Reactome, MSigDB, GenMAPP, Pathway Ontology, PantherDB and SMPDB pathway analysis results (Table 2), whereas pathways of down regulated were mainly enriched in inosine-5'-phosphate biosynthesis, sulfate activation for sulfonation, antigen processing and presentation, graft-versus-host disease, Fc-epsilon receptor I signaling in mast cells, TGF-beta receptor signaling, signaling by interleukins, generic transcription pathway, glycosaminoglycan degradation, sterol biosynthesis, ras-independent pathway in NK cell-mediated cytotoxicity, hypoxia and p53 in the cardiovascular system, FAS signaling pathway, angiogenesis, pathway of folate cycle/metabolism, vascular endothelial growth factor signaling, sarcosinemia and purine metabolism according to the BIOCYC, KEGG, PID, Reactome, MSigDB, GenMAPP, Pathway Ontology, PantherDB and SMPDB pathway analysis results (Table 3).
Gene ontology (GO) enrichment analysis of DEGs
The Gene Ontology (GO) enrichment analyses were conducted using online tool ToppGene. GO terms of the up regulated and down regulated genes s were listed in Table 4 and Table 5, respectively. Gene Ontology (GO) enrichment analysis showed that the up regulated genes were mainly associated with nervous system process, G protein-coupled receptor signaling pathway, intrinsic component of plasma membrane, extracellular matrix, transmembranesignaling receptor activity and receptor regulator activity. The down regulated genes were mainly associated with cell cycle, regulation of immune system process, nuclear membrane, nuclear chromatin, DNA-binding transcription factor activity, RNA polymerase II-specific and signaling receptor binding.
PPI network construction and module analysis
The PPI network of the up and down regulated genes was analyzed by using online STRING database. A total of 3715 nodes with 6518 edges were reflected in PPI network of up regulated genes is shown in Fig 4A. CAPN13, EGFR, ACTBL2, ACTL8, ERBB3, PRMT5, GATA4, RHOV, CHD5, MAGEL2, THNSL2, SLC38A8, THPO and SPTSSB were the hub genes with high node degree distribution, betweenness centrality, stress centrality, closeness centrality and low clustering coefficient in the network are listed in Table 6. A total of 5135 nodes with 10628 edges were reflected in PPI network of down regulated genes is shown in Fig 4B. FYN, PAK2, CUL3, RPS6, NOTCH2, PDE4D, SPATA21, MYBL1, SMURF1, PDGFRB, DLG3, ADHFE1, NMB, SLC25A36, MLLT1 and RNF2 were the hub genes with high node degree distribution, betweenness centrality, stress centrality, closeness centrality and low clustring coefficient in the network are listed in Table 6.
Four significant modules were selected for each up and down regulated genes using the PEWCC1E plug-in. The top four modules for up regulated genes were selected (module 13, 105 nodes and 235 edges; module 20, 77 nodes and 97 edges; module 21, 73 nodes and 81 edges; module 34, 53 nodes and 58 edges) are shown in Fig. 5A. The results revealed that hub genes (ACTG2, GATA4, EGFR, TP73, ACTBL2, FOXJ1, BMP7 and CDK5R2) in these significant modules were mostly enriched in the muscle contraction, notch-mediated HES/HEY network, cytokine-cytokine receptor interaction, E2F transcription factor network, actin cytoskeleton, epithelial cell differentiation, biological adhesion and neuron projection. Similarly, top four modules for down regulated genes were selected (module 1, 92 nodes and 186 edges; module 2, 56 nodes and 187 edges; module 5, 49 nodes and 144 edges; module 11, 29 nodes and 57 edges) are shown in Fig. 5B. The results revealed that hub genes (RPS6, PAK2, PODN, LMNA, EIF1AX, RPS27, HSPA8, FYN and LMNB1) in these significant modules were mostly enriched in the mTOR signaling pathway, Fc-epsilon receptor I signaling in mast cells, ensemble of genes encoding extracellular matrix and extracellular matrix-associated proteins, caspase cascade in apoptosis, postsynapse, cell cycle, regulation of immune system process and positive regulation of signal transduction.
Construction of target gene - miRNA regulatory network
NetworkAnalyst was applied to screen the miRNAs of the up and down regulated genes. The miRNAs predicted by at least two databases (among the following databases: DIANA-TarBase and miRTarBase) were diagnosed as the miRNAs of the target genes. Then, Cytoscape software was used to draw the target gene - miRNA regulatory network. The target gene - miRNA regulatory network for up regulated genes included 1867 nodes and 3735 edges (Fig. 6A). As shown in Table 7, TRIM72 regulates 123 miRNAs (ex,hsa-mir-4537), TET3 regulates 105 miRNAs (ex,hsa-mir-3148), NFIB regulates 89 miRNAs (ex,hsa-mir-4517), SLC19A3 regulates 80 miRNAs (ex,hsa-mir-4500) and SMOC1 regulates 123 miRNAs (ex,hsa-mir-6133) were considered as target gene. We performed pathway and GO enrichment analysis of these predicted target genes, which mainly enriched in muscle contraction, FOXA1 transcription factor network, intrinsic component of plasma membrane and biological adhesion. The target gene - miRNA regulatory network for down regulated genes included 2529 nodes and 10243 edges (Fig. 6B). As shown in Table 7, PPP1R15B regulates 168 miRNAs (ex, hsa-mir-7150), WEE1 regulates 167 miRNAs (ex,hsa-mir-3926), RPRD2 regulates 152 miRNAs (ex,hsa-mir-4452), LCOR regulates 146 miRNAs (ex,hsa-mir-4310) and SAR1A regulates 145 miRNAs (ex,hsa-mir-5698) were considered as target gene. We performed pathway and GO enrichment analysis of these predicted target genes, which mainly enriched in regulation of hydrolase activity, cell cycle, gene expression, nuclear chromatin and protein processing in endoplasmic reticulum.
Construction of target gene - TF regulatory network
NetworkAnalyst was applied to screen the TFs of the up and down regulated genes. The TFs predicted by database (ChEA database) was diagnosed as the TFs of the target genes. Then, Cytoscape software was used to draw the target gene - TF regulatory network. The target gene - TF regulatory network for up regulated genes included 539 nodes and 5790 edges (Fig. 7A). As shown in Table 8, ACTL8 regulates 145 TFs (ex, EGR1), LHFPL3 regulates 132 TFs (ex,SOX2), CXCL12 regulates 119 TFs (ex,SUZ12), GLI2 regulates 117 TFs (ex, AR) and C7 regulates 114 TFs (ex, TP53) were considered as target gene. We performed pathway and GO enrichment analysis of these predicted target genes, which mainly enriched in epithelial cell differentiation, cytokine-cytokine receptor interaction, pathways in cancer and innate immune system. The target gene - TF regulatory network for down regulated genes included 608 nodes and 10262 edges (Fig. 7B). As shown in Table 8, PRIM2 regulates 218 TFs (ex, SOX2), regulates 211 TFs (ex, MYC), GMDS regulates 210 TFs (ex, SPI1), C5ORF58 regulates 190 TFs (ex, RUNX1) and C10orf88 regulates 180 TFs (ex, FLI1) were considered as target gene. We performed pathway and GO enrichment analysis of these predicted target genes, which mainly enriched in metabolic pathways, gene expression, asparagine N-linked glycosylation and cell cycle.
Validations of hub genes
The ten hub genes (up and down regulated) were selected for further validation of their potential prognostic value. Upon comparing the expression of hub genes in the human protein atlas database (Fig. 8), it showed that CAPN13, ACTBL2, ERBB3, GATA4 and GNB4 were highly expressed in peripheral blood mononuclear cells of bone marrow, while NOTCH2, EXOSC10, RNF2, PSMA1 and PRKAA1 were less expressed in peripheral blood mononuclear cells of bone marrow ROC analysis was performed from the 10 hub genes from GSE113079. The ROC curves of these ten hub genes all indicated favorable prognostic values for CAD. In addition, the area under ROC curve (AUC) of CAPN13, ACTBL2, ERBB3, GATA4, GNB4, NOTCH2, EXOSC10, RNF2, PSMA1 and PRKAA1 were 0.855 (p = 2.406664e-08), 0.923 (p = 5.413565e-10), 0.829 (p =2.857413e-08), 0.903 (p = 4.513268e-09), 0.918 (p = 6.358925e-09), 0.891 (p = 4.367291e-09), 0.927 (p = 1.344048e-10), 0.911 (p = 5.076899e-10), 0.892 (3.057148e-08) and 0.904 (p = 1.902203e-08), respectively (Fig. 9).
Molecular docking studies
The prevailing work aimsto discover the significant interactions responsible for complex stability with the receptor of the binding sites by docking studies. The docking studies was executed using Sybyl X 2.1 software on designed molecules which includes derivatives of dihydropyridine heterocyclic nucleus found in amlodipine a beta blocker normally used in in coronary artery disease. Beta-blockers suppress the heart's sympathetic activation, decreasing heart rate and contractility that lower the need for myocardial oxygen and thereby prevent or alleviate angina in CAD patients. Since beta-blockers suppress the heart rate during exercise, the initiation of angina or the ischemic threshold is postponed or stopped during exercise. In the treatment of exertional angina, all forms of beta-blockers tend to be equally successful. Based on the structure of amlodipine containing dihydropyridine heterocyclic nucleus the molecules containing dihydropyridine are designed to identify for docking studies in the present research. A total of 34 common dihydropyridine derivatives were developed and amlodipine was used as a standard for docking studies on over-expressed proteins, and the structures are shown in Fig.10, respectively. The one protein from each over expressed genes in coronary artery diseases such as ACTBL2 (Actin beta-like 2), CAPN13 (Calpain 13), ERBB3 (Erythroblasticleukemia viral oncogene homolog 3), GATA4 (GATA binding protein 4), GNB4 (Guanine nucleotide binding protein beta polypeptide 4)and their X-RAY crystallographic structure and co-crystallized PDB code 2FF3, 2I7A, 3LMG, 3DFV and 6UQ3 respectively were constructed for docking. To identify the potential molecule and its binding affinity to proteins, the docking was carried out on built molecules. It is said that the molecule with C-score greater than 5 are active and few molecules with particular proteins obtained greater than 8 respectively. Docking experiments were carried out on a total of 34 designed molecules, few obtained an outstanding C-score greater than 8 and few molecules obtained an optimum binding score of 4-4.9 then obtained less binding sore of 2.0-3.0 respectively. The molecule IM1, IM4, IM7, IM9, IM10,IM11, IM12, IM13, IM14, IM15, IM16, TZ19, TZ21, TZ25, TZ26, TZ28, TZ29 and IM1, IM9, IM10, IM11, IM13, TZ26 and IM12 with protein of PDB code 2FF3 and 3LMG and 3DFV respectively obtained excellent binding score of more than 7.Good binding score of 5 to 6.99 obtained from the molecules areIM2, IM3, IM5, IM6, IM8, IM17, TZ18, TZ20, TZ22, TZ23, TZ24, TZ27, TZ30, TZ31, TZ32, TZ33, TZ34 and IM7, IM11, IM12, TZ27 and IM2, IM3, IM4, IM5, IM6, IM7, IM8, IM12, IM14, IM16, IM17, TZ18, TZ19, TZ20, TZ21, TZ22, TZ23, TZ24, TZ25, TZ28, TZ29, TZ30, TZ31, TZ32, TZ33 and IM1, IM2, IM3, IM4, IM6, IM7, IM8, IM9, IM10, IM11, IM13, IM14, IM16, TZ18, TZ20, TZ23, TZ24, TZ26, TZ27, TZ28, TZ29 and IM7, IM8, IM10, IM11, IM12, IM13, IM13, IM16, TZ23, TZ25, TZ26 with PDB protein of 2FF3, 2I7A, 3LMG, 3DFV and 6UQ3 respectively. Molecules with optimum binding score are IM1, IM2, IM3, IM4, IM5, IM6, IM8, IM9, IM10, IM13, IM14, IM15, IM16, IM17, TZ18, TZ19, TZ20, TZ21, TZ22, TZ23, TZ24, TZ25, TZ26, TZ28, TZ29, TZ30, TZ31, TZ32, TZ33, TZ34 and TZ27, TZ34 and IM5, IM15, IM17, TZ19, TZ21, TZ22, TZ25, TZ30, TZ31, TZ32, TZ33, TZ34 and IM1, IM2, IM3, IM4, IM5, IM6, IM9, IM14, IM15, IM17, TZ18, TZ19, TZ20, TZ21, TZ22, TZ24, TZ27, TZ28, TZ29, TZ30, TZ31, TZ32, TZ33, TZ34with PDB code of 2I7A, 3LMG, 3DFV and 6UQ3 and the molecule IM7 obtained highest binding score of 9.00 greater than the standard amlodipine with PDB 2FF3 respectively the values are depicted in Table 9. The standard amlodipine obtained good binding score with 3LMG, 2FF3 and 6UQ3, and obtained optimum binding score with PDB 2I7A and 3DFV respectively. The Fig. 11 and Fig. 12 depicts 3D hydrogen bonding interactions of lignd with Protein, with aminoacids and other bonding interactions with amino acids and Fig. 13 depicts the 2D interactions with amino acids and their distance with protein code 2FF3 of molecule IM7 are depicted by 3D and 2D respectively.