Identification of candidate biomarkers and therapeutic agents for heart failure by bioinformatics analysis

Abstract

Introduction

Heart failure (HF) is a heterogeneous clinical syndrome and affects millions of people all over the world. HF occurs when the cardiac overload and injury, which is a worldwide complaint. The aim of this study was to screen and verify hub genes involved in developmental HF as well as to explore active drug molecules.

Methods

The expression profiling by high throughput sequencing of GSE141910 dataset was downloaded from the Gene Expression Omnibus (GEO) database, which contained 366 samples, including 200 heart failure samples and 166 non heart failure samples. The raw data was integrated to find differentially expressed genes (DEGs) and were further analyzed with bioinformatics analysis. Gene ontology (GO) and REACTOME enrichment analyses were performed via ToppGene; protein–protein interaction (PPI) networks of the DEGs was constructed based on data from the HiPPIE interactome database; modules analysis was performed; target gene—miRNA regulatory network and target gene—TF regulatory network were constructed and analyzed; hub genes were validated; molecular docking studies was performed.

Results

A total of 881 DEGs, including 442 up regulated genes and 439 down regulated genes were observed. Most of the DEGs were significantly enriched in biological adhesion, extracellular matrix, signaling receptor binding, secretion, intrinsic component of plasma membrane, signaling receptor activity, extracellular matrix organization and neutrophil degranulation. The top hub genes ESR1, PYHIN1, PPP2R2B, LCK, TP63, PCLAF, CFTR, TK1, ECT2 and FKBP5 were identified from the PPI network. Module analysis revealed that HF was associated with adaptive immune system and neutrophil degranulation. The target genes, miRNAs and TFs were identified from the target gene—miRNA regulatory network and target gene—TF regulatory network. Furthermore, receiver operating characteristic (ROC) curve analysis and RT-PCR analysis revealed that ESR1, PYHIN1, PPP2R2B, LCK, TP63, PCLAF, CFTR, TK1, ECT2 and FKBP5 might serve as prognostic, diagnostic biomarkers and therapeutic target for HF. The predicted targets of these active molecules were then confirmed.

Conclusion

The current investigation identified a series of key genes and pathways that might be involved in the progression of HF, providing a new understanding of the underlying molecular mechanisms of HF.

Introduction

Heart failure (HF) is a cardiovascular disease characterized by tachycardia, tachypnoea, pulmonary rales, pleural effusion, raised jugular venous pressure, peripheral oedema and hepatomegaly [1]. Morbidity and mortality linked with HF is a prevalent worldwide health problem holding a universal position as the leading cause of death [2]. The numbers of cases of HF are rising globally and it has become a key health issue. According to a survey, the prevalence HF is expected to exceed 50% of the global population [3]. Research suggests that modification in multiple genes and signaling pathways are associated in controlling the advancement of HF. However, a lack of investigation on the precise molecular mechanisms of HF development limits the treatment efficacy of the disease at present.

Previous study showed that HF was related to the expression of MECP2 [4] RBM20 [5], CaMKII [6], troponin I [7] and SERCA2a [8]. Toll-Like receptor signaling pathway [9], activin type II receptor signaling pathway [10], CaMKII signaling pathways [11], Drp1 signaling pathways [12] and JAK-STAT signaling pathway [13] were liable for progression of HF. More investigations are required to focus on treatments that enhance the outcome of patients with HF, to strictly make the diagnosis of the disease based on screening of biomarkers. These investigations can upgrade prognosis of patients by lowering the risk of advancement of HF and related complications. So it is essential to recognize the mechanism and find biomarkers with a good specificity and sensitivity.

The recent high-throughput RNA sequencing data has been widely employed to screen the differentially expressed genes (DEGs) between normal samples and HF samples in human beings, which makes it accessible for us to further explore the entire molecular alterations in HF at multiple levels involving DNA, RNA, proteins, epigenetic alterations, and metabolism [14]. However, there still exist obstacles to put these RNA seq data in application in clinic for the reason that the number of DEGs found by expression profiling by high throughput sequencing were massive and the statistical analyses were also too sophisticated [15,16,17,18,19]

In this study, first, we had chosen dataset GSE141910 from Gene Expression Omnibus (GEO) (http:// www.ncbi.nlm.nih.gov/geo/) [20]. Second, we applied for limma tool in R software to obtain the differentially expressed genes (DEGs) in this dataset. Third, the ToppGene was used to analyze these DEGs including biological process (BP), cellular component (CC) and molecular function (MF) REACTOME pathways. Fourth, we established protein–protein interaction (PPI) network and then applied Cytotype PEWCC1 for module analysis of the DEGs which would identify some hub genes. Fifth, we established target gene—miRNA regulatory network and target gene—TF regulatory network. In addition, we further validated the hub genes by receiver operating characteristic (ROC) curve analysis and RT-PCR analysis. Finally, we performed molecular docking studies for over expressed hub genes. Results from the present investigation might provide new vision into potential prognostic and therapeutic targets for HF.

Materials and methods

Data resource

Expression profiling by high throughput sequencing with series number GSE141910 based on platform GPL16791 was downloaded from the GEO database. The dataset of GSE141910 contained 200 heart failure samples and 166 non heart failure samples. It was downloaded from the GEO database in NCBI based on the platform of GPL16791 Illumina HiSeq 2500 (Homo sapiens).

Identification of DEGs in HF

DEGs of dataset GSE141910 between HF groups and non heart failure groups were respectively analyzed using the limma package in R [21]. Fold changes (FCs) in the expression of individual genes were calculated and DEGs with P < 0.05, |log FC|> 1.158 for up regulated genes and |log FC|< − 0.83 for down regulated genes were considered to be significant. Hierarchical clustering and visualization were used by Heat-map package of R.

Functional enrichment analysis

Gene Ontology (GO) analysis and REACTOME pathway analysis were performed to determine the functions of DEGs using the ToppGene (ToppFun) (https://toppgene.cchmc.org/enrichment.jsp) [22] GO terms (http://geneontology.org/) [23] included biological processes (BP), cellular components (CC) and molecular functions (MF) of genomic products. REACTOME (https://reactome.org/) [24] analyzes pathways of important gene products. ToppGene is a bioinformatics database for analyzing the functional interpretation of lists of proteins and genes. The cutoff value was set to P < 0.05.

Protein–protein interaction network construction and module screening

PPI networks are used to establish all protein coding genes into a massive biological network that serves an advance compassionate of the functional system of the proteome [25]. The HiPPIE interactome (https://cbdm.uni-mainz.de/hippie/) [26] database furnish information regarding predicted and experimental interactions of proteins. In the current investigation, the DEGs were mapped into the HiPPIE interactome database to find significant protein pairs with a combined score of > 0.4. The PPI network was subsequently constructed using Cytoscape software, version 3.8.2 (www.cytoscape.org) [27]. The nodes with a higher node degree [28], higher betweenness centrality [29], higher stress centrality [30] and higher closeness centrality [31] were considered as hub genes. Additionally, cluster analysis for identifying significant function modules with a degree cutoff > 2 in the PPI network was performed using the PEWCC1 (http://apps.cytoscape.org/apps/PEWCC1) [32] in Cytoscape.

Target gene—miRNA regulatory network construction

The miRNet database (https://www.mirnet.ca/) [33] contains information on miRNA and the regulated genes. Using information collected from the miRNet database, hub genes were matched with their associated miRNA. The target gene—miRNA regulatory network then was constructed using Cytoscape software. MiRNAs and target are selected based on highest node degree.

Target gene—TF regulatory network construction

The NetworkAnalyst database (https://www.networkanalyst.ca/) [34] contains information on TF and the regulated genes. Using information collected from the NetworkAnalyst database, hub genes were matched with their associated TF. The target gene—TF regulatory network then was constructed using Cytoscape software. TFs and target genes are selected based on highest node degree.

Receiver operating characteristic (ROC) curve analysis

Then ROC curve analysis was implementing to classify the sensitivity and specificity of the hub genes for HF diagnosis and we investigated how large the area under the curve (AUC) was by using the statistical package pROC in R software [35].

RT-PCR analysis

H9C2 cells (ATCC) were cultured in Dulbecco’s minimal essential medium (DMEM) (Sigma-Aldrich) supplemented with 10% fetal calf serum (Sigma-Aldrich) and 1% streptomycin (Sigma-Aldrich) at 37 °C in 5% CO₂. HL-1 cells (ATCC) was culture in Claycomb medium (Sigma-Aldrich) supplemented with 10% fetal bovine serum (Sigma-Aldrich), 1% streptomycin (Sigma-Aldrich), 1% glutamax (Sigma-Aldrich) and 0.1 mM norepinephrine (Sigma-Aldrich) at 37 °C in 5% CO₂. Total RNA was isolated from cell culture of H9C2 for HF and HL-1 for normal control using the TRI Reagent (Sigma, USA). cDNA was synthesized using 2.0 μg of total RNA with the Reverse transcription cDNA kit (Thermo Fisher Scientific, Waltham, MA, USA). The 7 Flex real-time PCR system (Thermo Fisher Scientific, Waltham, MA, USA) was employed to detect the relative mRNA expression. The relative expression levels were determined by the 2-ΔΔCt method and normalized to internal control beta-actin [36]. All RT-PCR reactions were performed in triplicate. The primers used to explore mRNA expression of ten hub genes were listed in Table 1.

Table 1 The sequences of primers for quantitative RT-PCR

Identification of candidate small molecules

SYBYL-X 2.0 perpetual drug design software has been used for surflex-docking studies of the designed novel molecules and the standard on over expressed genes of PDB protein. Using ChemDraw Software, all designed molecules and standards were sketched, imported and saved using open babel free software in sdf. template. The protein of over expressed genes of ESR1, LCK, PPP2R2B, TP63 and their co-crystallised protein of PDB code 4PXM, 1KSW, 2HV7, 3VD8 and 6RU6 were extracted from Protein Data Bank [37,38,39,40]. Optimizations of the designed molecules were performed by standard process by applying Gasteiger Huckel (GH) charges together with the TRIPOS force field. In addition, energy minimization was achieved using MMFF94s and MMFF94 algorithm methods. The preparation of the protein was done after protein incorporation. The co-crystallized ligand and all water molecules have been eliminated from the crystal structure; more hydrogen’s were added and the side chain was set, TRIPOS force field was used for the minimization of structure. The interaction efficiency of the compounds with the receptor was expressed in kcal/mol units by the Surflex-Dock score. The best location was integrated into the molecular region by the interaction between the protein and the ligand. Using Discovery Studio Visualizer, the visualisation of ligand interaction with receptor is performed.

Results

Identification of DEGs in HF

We identified 881 DEGs in the GSE141910 dataset using the limma package in R. Based on the limma analysis, using the adj P val < 0.05, |log FC|> 1.158 for up regulated genes and |log FC|< − 0.83 for down regulated genes, a total of 881 DEGs were identified, consisting of 442 genes were up regulated and 439 genes were down regulated. The DEGs are listed in Additional file 1: Table S1. The volcano plot for DEGs is illustrated in Fig. 1. Figure 2 is the hierarchical clustering heat-map.

Fig. 1

Volcano plot of differentially expressed genes. Genes with a significant change of more than two-fold were selected. Green dot represented up regulated significant genes and red dot represented down regulated significant genes

Fig. 2

Heat map of differentially expressed genes. Legend on the top left indicate log fold change of genes. (A1 – A200 = heart failure samples; B1 – B166 = non heart failure samples)

Functional enrichment analysis

Results of GO analysis showed that the up regulated genes were significantly enriched in BP, CC, and MF, including biological adhesion, regulation of immune system process, extracellular matrix, cell surface, signaling receptor binding and molecular function regulator (Table 2); the down regulated genes were significantly enriched in BP, CC, and MF, including secretion, defense response, intrinsic component of plasma membrane, whole membrane, signaling receptor activity and molecular transducer activity (Table 2). Pathway analysis showed that the up regulated genes were significantly enriched in extracellular matrix organization and immunoregulatory interactions between a lymphoid and a non-lymphoid cell (Table 3); the down regulated genes were significantly enriched in neutrophil degranulation and SLC-mediated transmembrane transport (Table 3).

Table 2 The enriched GO terms of the up and down regulated differentially expressed genes

Table 3 The enriched pathway terms of the up and down regulated differentially expressed genes

Protein–protein interaction (PPI) network and module analysis

Based on the HiPPIE interactome database, the PPI network for the DEGs (including 6541 nodes and 13,909 edges) was constructed (Fig. 3A). Up regulated gene with higher node degree, higher betweenness centrality, higher stress centrality and higher closeness centrality were as follows: ESR1, PYHIN1, PPP2R2B, LCK, TP63 and so on. Down regulated genes had higher node degree, higher betweenness centrality, higher stress centrality and higher closeness centrality were as follows PCLAF, CFTR, TK1, ECT2, FKBP5 and so on. The node degree, betweenness centrality, stress centrality and closeness centrality are listed in Table 4.

Fig. 3

PPI network and the most significant modules of DEGs. A The PPI network of DEGs was constructed using Cytoscape. B The most significant module was obtained from PPI network with 4 nodes and 6 edges for up regulated genes. C The most significant module was obtained from PPI network with 6 nodes and 10 edges for down regulated genes. Up regulated genes are marked in green; down regulated genes are marked in red

Table 4 Topology table for up and down regulated genes

Additionally, two significant modules, including module 1 (10 nodes and 24 edges) and module 2 (5 nodes and 10 edges) (Fig. 3B) and module 3 (55 nodes and 115 edges), were acquired by PEWCC1 plug-in (Fig. 3C). Furthermore, GO terms and REACTOME pathways were significantly enriched by module 1, including adaptive immune system, immunoregulatory interactions between a lymphoid and a non-lymphoid cell, hemostasis, biological adhesion and regulation of immune system process. Meanwhile, the nodes in module 2 were significantly enriched in GO terms and REACTOME pathways, including neutrophil degranulation and secretion.

Target gene—miRNA regulatory network construction

Associations between 2063 miRNAs and their 319 target genes were collected from the target gene—miRNA regulatory network (Fig. 4). MiRNAs of hsa-mir-4533, hsa-mir-548ac, hsa-mir-548i, hsa-mir-5585-3p, hsa-mir-6750-3p, hsa-mir-200c-3p, hsa-mir-1273 g-3p, hsa-mir-1244, hsa-mir-4789-3p and hsa-mir-766-3p, which exhibited a high degree of interaction, were selected from this network. Furthermore, the results also showed that FSCN1 was the target of hsa-mir-4533, ESR1 was the target of hsa-mir-548ac, TMEM30B was the target of hsa-mir-548i, SCN2B was the target of hsa-mir-5585-3p, CENPA was the target of hsa-mir-6750-3p, FKBP5 was the target of hsa-mir-200c-3p, PCLAF was the target of hsa-mir-1273g-3p, CEP55 was the target of hsa-mir-1244, ATP2A2 was the target of hsa-mir-4789-3p and TK1 was the target of hsa-mir-766-3p, and are listed in Table 5.

Fig. 4

Target gene—miRNA regulatory network between target genes. The light orange color diamond nodes represent the key miRNAs; up regulated genes are marked in green; down regulated genes are marked in red

Table 5 miRNA—target gene and TF—target gene interaction

Target gene—TF regulatory network construction

Associations between 330 TFs and their 247 target genes were collected from the target gene—TF regulatory network (Fig. 5). TFs of ESRRA, RERE, HMG20B, THRAP3, ATF1, MXD3, ARID4B, CBFB, TAF7 and CREM, which exhibited a high degree of interaction, were selected from this network. Furthermore, the results also showed that FSCN1 was the target of ESRRA, APOA1 was the target of RERE, COL1A1 was the target of HMG20B, HBB was the target of THRAP3, LCK was the target of ATF1, SOCS3 was the target of MXD3, BCL6 was the target of ARID4B, FKBP5 was the target of CBFB, ANLN was the target of TAF7 and ATP2A2 was the target of CREM, and are listed in Table 5.

Fig. 5

Target gene—TF regulatory network between target genes. The sky blue color triangle nodes represent the key TFs; up regulated genes are marked in green; down regulated genes are marked in red

Receiver operating characteristic (ROC) curve analysis

First of all, we performed the ROC curve analysis among 10 hub genes based on the GSE141910. The results showed that ESR1, PYHIN1, PPP2R2B, LCK, TP63, PCLAF, CFTR, TK1, ECT2 and FKBP5 achieved an AUC value of > 0.7, demonstrating that these ten genes have high sensitivity and specificity for HF, suggesting they can be served as biomarkers for the diagnosis of HF (Fig. 6).

Fig. 6

ROC curve analyses of hub genes. A ESR1, B PYHIN1, C PPP2R2B, D LCK, E TP63, F PCLAF, G CFTR, H TK1, I ECT2, J FKBP5

RT-PCR analysis

RT-PCR was used to validate the hub genes between normal and HF cell lines. The results suggested that the mRNA expression level of ESR1, PYHIN1, PPP2R2B, LCK and TP63 were significantly increased in HF compared with that in normal, while PCLAF, CFTR, TK1, ECT2 and FKBP5 were significantly decreased in HF compared with that in normal and are shown in Fig. 7.

Fig. 7

RT-PCR analyses of hub genes. A ESR1, B PYHIN1, C PPP2R2B, D LCK, E TP63, F PCLAF, G) CFTR H TK1, I ECT2, J FKBP5

Identification of candidate small molecules

In the present study docking simulations are performed to spot the active site and foremost interactions accountable for complex stability with the receptor binding sites. In heart failure recognized over expressed genes and their proteins of x-ray crystallographic structure are chosen from PDB for docking studies. Most generally, medications containing benzothiadiazine ring hydrochlorothiazide are used in heart failure either alone or in conjunction with other drugs, based on this the molecules containing heterocyclic ring of benzothiadiazine are designed and hydrochlorothiazide is uses as a reference standard. Docking experiments using Sybyl-X 2.1.1. drug design perpetual software were used on the designed molecules. Docking studies were performed in order to understand the biding interaction of standard hydrochlorothiazide and designed molecules on over expressed protein. The X- RAY crystallographic structure of one proteins from each over expressed genes of ESR1, LCK, PPP2R2B, PYHIN1, TP63 and their co-crystallised protein of PDB code 4PXM, 1KSW, 2HV7, 3VD8 and 6RU6 respectively were selected for the docking studies to identify and predict the potential molecule based on the binding score with the protein and successful in heart failure. For the docking tests, a total of 34 molecules were built and the molecule with binding score greater than 5 is believed to be good. The designed molecules obtained docking score of 5 to 7 were HIM10, HTZ5, HIM6, HTZ31, HIM3, HIM14, HIM1, HIM7 and HIM11, HIM16, HTZ9, HIM17, HIM12, HTZ12, HIM6, HTZ7, HIM10, HTZ3 and HIM8, HTZ9, HIM6, HIM4, HIM13, HTZ16, HIM9, HIM7, HTZ5, HIM16, HTZ7, HIM10, HIM5, HIM12, HIM15, HTZ12, HIM3, HIM14 and HIM14, HIM6, HIM17, HTZ7, HIM10, HIM1, HTZ9, HIM3, HIM16, HIM15, HIM8, HIM9, HIM7, HTZ10, HTZ3, HTZ5, HTZ1, HIM13, HTZ4, HIM11, HTZ12, HTZ14, HIM2 and HIM7, HTZ13, HTZ5, HIM15, HIM12, HIM6, HTZ11, HIM14, HTZ9, HIM11, HIM13, HIM9, HIM8, HIM10, HIM1, HIM5, HIM4, HTZ12, HIM2, HIM17, HIM3, HTZ1, HTZ8, HIM3, HTZ14, HTZ3 with proteins 4PXM and, 1KSW and 2HV7 and 3VD8 and 6RUR respectively (Fig. 8). The molecules obtained binding score of less than 5 were HTZ13, HTZ12, HTZ10, HIM3, HIM15, HIM16, HIM13, HIM8, HTZ16, HIM2, HIM4, HIM17, HTZ17, HIM11, HTZ5, HTZ3, HIM9, HTZ15, HTZ5, HTZ9, HTZ11, HIM5, HTZ8 and HTZ14, HIM14, HTZ13, HIM13, HTZ16, HIM2, HIM3, HTZ10, HIM7, HIM1, HTZ1, HTZ4, HIM8, HIM5, HTZ2, HIM9, HTZ5, HTZ15, HTZ3, HIM4, HIM15, HTZ17, HTZ8, HTZ11 and HTZ14, HIM2, HIM1, HTZ11, HIM17, HTZ13, HTZ4, HTZ2, HIM3, HTZ15, HTZ8, HTZ17, HTZ1, HTZ3 and HTZ8, HIM4, HTZ16, HTZ15, HIM5, HTZ11, HTZ13, HIM3, HTZ17, HTZ2 and HTZ7, HTZ4, HTZ2, HTZ17, HTZ15 with proteins 4PXM and, 1KSW and 2HV7 and 3VD8 and 6RUR respectively. The molecules obtained very less binding score are HTZ1, HIM12, HTZ2, HTZ4 with protein 4PXM and the standard hydrochlorothiazide (HTZ) obtained less binding score with all proteins, the values are depicted in Table 6.

Fig. 8

Structures of designed molecules

Table 6 Docking results of Designed Molecules on Over Expressed Proteins

Discussion

HF is the most prevalent form of cardiovascular disease among the elderly. A complete studies of HF, comprising pathogenic factors, pathological processes, clinical manifestations, early clinical diagnosis, clinical prevention, and drug therapy targets urgency to be consistently analyzed. In the present investigation, bioinformatics analysis was engaged to explore HF biomarkers and the pathological processes in myocardial tissues, acquired from HF groups and non heart failure groups. We analyzed GSE141910 expression profiling by high throughput sequencing obtained 881 different genes between HF groups and non heart failure groups, 442 up regulated and 439 down regulated genes. HBA2 and HBA1 have a key role in hypertension [41], but these genes might be linked with development HF. SFRP4 was linked with progression of myocardial ischemia [42]. Emmens et al. [43] and Broch et al. [44] found that PENK (proenkephalin) and IL1RL1 were up regulated in HF. ALOX15B has lipid accumulation and inflammation activity and is highly expressed in atherosclerosis [45]. Studies have shown that expression of MYH6 was associated with hypertrophic cardiomyopathy [46].

In functional enrichment analysis, some genes involved with regulation of cardiovascular system processes were enriched in HF. Liu et al. [47], Kosugi et al. [48], McMacken et al. [49], Pan and Zhang [50], Li et al. [51] and Jiang et al. [52] presented that expression of HLA-DQA1, KDM5D, UCHL1, SAA1, ARG1 and LYVE1 were associated with progression of cardiomyopathy. Hou et al. [53] and Olesen et al. [54] demonstrated that DACT2 and KCND3 were found to be substantially related to atrial fibrillation. Ge and Concannon [55], Ferjeni et al. [56], Anquetil et al. [57], Glawe et al. [58], Kawabata et al. [59], Li et al. [60], Buraczynska et al. [61], Amini et al. [62], Yang et al. [63], Du Toit et al. [64], Hirose et al. [65], Zhang et al. [66], Griffin et al. [67], Zouidi et al. [68], Trombetta et al. [69], Alharbi et al. [70], Ikarashi et al. [71], Dharmadhikari et al. [72], Sutton et al. [73] and Deng et al. [74] reported that UBASH3A, ZAP70, IDO1, ITGAL (integrin subunit alpha L). ITGB7, RASGRP1, CNR1, SLC2A1, SLC11A1, GPR84, SSTR5, KCNB1, GLUL (glutamate-ammonia ligase), BANK1, CACNA1E, LGR5, AQP3, SIGLEC7, SSTR2 and DNER (delta/notch like EGF repeat containing) could be an index for diabetes, but these genes might be responsible for progression of HF. Experiments show that expression of FAP (fibroblast activation protein alpha) [75], THBS4 [76], CD27 [77], LEF1 [78], CTHRC1 [79], ESR1 [80], CXCL9 [81], SERPINA3 [82], TRPC4 [83], F13A1 [84], PIK3C2A [85], KCNIP2 [86] and GPR4 [87] contributed to myocardial infarction. MFAP4 [88], ALOX15 [89], COL1A1 [90], APOA1 [91], PDE5A [92], CX3CR1 [93], THY1 [94], GREM1 [95], FMOD (fibromodulin) [96], NPPA (natriuretic peptide A) [97], LTBP2 [98], LUM (lumican) [99], IL34 [100], NRG1 [101], CXCL14 [102], CXCL10 [103], ACE (angiotensin I converting enzyme) [104], CFTR (ystic fibrosis transmembrane conductance regulator) [105], S100A8 [106], S100A9 [106], HP (haptoglobin) [107], AGTR1 [108], ATP2A2 [109], IL10 [110], EDN1 [111], TLR2 [112], MCEMP1 [113], TPO (thyroid peroxidase) [114], CD163 [115], IL18R1 [116], KCNA7 [117] and CALCRL (calcitonin receptor like receptor) [118] have an important role in HF. Li et al. [119], Deckx et al. [120], Ichihara et al. [121] and Paik et al. [122] showed that the SERPINE2, OGN (osteoglycin), AGTR2 and WNT10B promoted cardiac interstitial fibrosis. Cai et al. [123], Mo et al. [124], Sun et al. [125], Martinelli et al. [126], Zhao et al. [127], Assimes et al. [128] and Piechota et al. [129] showed that CCR7, FCN1, ESM1, F8 (coagulation factor VIII), C1QTNF1, ALOX5 and MSR1 were an important target gene for coronary artery disease. STAB2 have been suggested to be associated with venous thromboembolic disease [130]. Genes such as COMP (cartilage oligomeric matrix protein) [131], CHI3L1 [132], PLA2G2A [133], P2RY12 [134], CR1 [135], HPSE (heparanase) [136], PTX3 [137] and SERPINE1 [138] were related to atherosclerosis. CCDC80 [139], CMA1 [140], MDK (midkine) [141], GNA14 [142], SCG2 [143], NPPB (natriuretic peptide B) [144], FGF10 [145], ARNTL (aryl hydrocarbon receptor nuclear translocator like) [146], WNK3 [147], EDNRB (endothelin receptor type B) [148], THBS1 [149], SELE (selectin E) [150], SLC4A7 [151], AQP4 [152] and KCNK3 [153] are thought to be responsible for progression of hypertension, but these genes might to be associated with progression of HF. CNTNAP2 [154], GLI2 [155], DPT (dermatopontin) [156], AEBP1 [157], ITIH5 [158], CXCL11 [159], GDNF (glial cell derived neurotrophic factor) [160], MCHR1 [161], FLT3 [162], ELANE (elastase, neutrophil expressed) [163], OSMR (oncostatin M receptor) [164] and IL15RA [165] are involved in development of obesity, but these genes might be key for progression of HF. CTSG (cathepsin G) is a protein coding gene plays important roles in aortic aneurysms [166]. Evidence from Safa et al. [167], Chen et al. [168], Zhou et al. [169], Hu et al. [170], Lou et al. [171], Zhang et al. [172] and Chen et al. [173] study indicated that the expression of CCL22, CCR1, FPR1, KNG1, CRISPLD2, CD38 and GPRC5A were linked with progression of ischemic heart disease. Li et al. [174] showed that STEAP3 expression can be associated with cardiac hypertrophy progression.

The HiPPIE interactome database was used to construct the PPI network, and modules analysis was performed. We finally screened out up regulated hub genes and down regulated hub genes, including ESR1, PYHIN1, PPP2R2B, LCK, TP63, PCLAF, CD247, CD2, CD5, CD48, CFTR, TK1, ECT2, FKBP5, S100A9 and S100A8 from the PPI network and its modules. TP63 might serve as a potential prognostic factor in cardiomyopathy [175]. The expression of FKBP5 is related to the progression of coronary artery disease [176]. CD247 plays a central role in hypertension [177], but this gene might be involved in the HF. PYHIN1, PPP2R2B, LCK (LCK proto-oncogene, Src family tyrosine kinase), PCLAF (PCNA clamp associated factor), TK1, ECT2, CD2, CD5 and CD48 might be the novel biomarker for HF.

The miRNet database and NetworkAnalyst database were used to construct the target gene—miRNA regulatory network and target gene—TF regulatory network. We finally screened out target genes, miRNA, TFs, including FSCN1, ESR1, TMEM30B, SCN2B, CENPA, FKBP5, PCLAF, CEP55, ATP2A2, TK1, APOA1, COL1A1, HBB, LCK, SOCS3, BCL6, ANLN, hsa-mir-4533, hsa-mir-548ac, hsa-mir-548i, hsa-mir-5585-3p, hsa-mir-6750-3p, hsa-mir-200c-3p, hsa-mir-1273 g-3p, hsa-mir-1244, hsa-mir-4789-3p, hsa-mir-766-3p, ESRRA, RERE, HMG20B, THRAP3, ATF1, MXD3, ARID4B, CBFB, TAF7 and CREM from the target gene—miRNA regulatory network and target gene—TF regulatory network. SCN2B [178] and SOCS3 [179] are considered as a markers for HF and might be a new therapeutic target. BCL6 levels are correlated with disease severity in patients with atherosclerosis [180]. A previous study showed that hsa-mir-1273 g-3p [181], hsa-mir-4789-3p [182] and ATF1 [183] could involved in hypertension, but these markers might be responsible for progression of HF. hsa-miR-518f, was demonstrated to be associated with cardiomyopathy [184]. An evidence demonstrating a role for ESRRA (estrogen related receptor alpha) [185] and THRAP3 [186] in diabetes, but these genes might be liable for development of HF. FSCN1, TMEM30B, CENPA (centromere protein A), CEP55, HBB (hemoglobin subunit beta), ANLN (anillin actin binding protein), hsa-mir-4533, hsa-mir-548ac, hsa-mir-548i, hsa-mir-5585-3p, hsa-mir-6750-3p, hsa-mir-200c-3p, hsa-mir-1244, RERE(arginine-glutamic acid dipeptide repeats), HMG20B, MXD3, ARID4B, CBFB (core-binding factor subunit beta), TAF7 and CREM (cAMP response element modulator) might be the novel biomarker for HF.

The molecules HIM6, HIM10 obtained good binding score of more 5 to 6.999 with all proteins and the molecules HIM11, HIM12, HIM14, HTZ9, HTZ10 and HTZ12 obtained binding score above 5 and less than 9 with PDB protein code of 2HV7, 3VD8 and 6RUR respectively. The molecule HIM11 obtained highest binding score of 8.678 with 2HV7 and its interaction with amino acids are molecule HIM11 (Fig. 9) has obtained with a high binding score with PDB protein 2HV7, the interactions of molecule is the C6 side chin acyl carbonyl C=O formed hydrogen bond interaction with amino acid GLN-207 with bond length 1.92 A° and 3’ N–H group of imidazole ring formed hydrogen bond interaction with VAL-305 with bond length 2.36 A° respectively. It also formed other interactions of carbon hydrogen bond of –CH₃ group of carboxylate at C6 with PRO-304 and amide-pi stacked and pi–pi stacked interaction of electrons of aromatic ring A with ALA-204 and ring C with HIS-155 and HIS-308. Molecule formed pi-alkyl interaction of ring B with PRO-304 and all interactions with amino acids and bond length are depicted by 3D and 2D figures (Fig. 10 and Fig. 11).

Fig. 9

Structure of active designed molecule of HIM11

Fig.10

3D binding of molecule HIM11 with 2HV7

Fig.11

2D binding of molecule HIM11 with 2HV7

Conclusions

The present investigation aimed at characterizing the expression profiling by high throughput sequencing of the HF patients. Our bioinformatics analyses revealed key gene signatures as candidate biomarkers in HF. Hub genes (ESR1, PYHIN1, PPP2R2B, LCK, TP63, PCLAF, CFTR, TK1, ECT2 and FKBP5) were diagnosed as an essential genetic factors in HF. In general, DEGs linked with HF genes, including already known markers of HF and other HF related diseases, and novel biomarkers, were diagnosed. Potential implicated miRNAs and TFs were also diagnosed. The diagnosed hub genes might represent candidate diagnostic and prognostic biomarkers, and therapeutic targets. The current investigation reported novel genes and signaling pathways in HF, and further investigation is required.