Identification of Hub Genes Associated with Hypertension and Their Interaction with miRNA Based on Weighted Gene Coexpression Network Analysis (WGCNA) Analysis.

BACKGROUND Hypertension is one of the most widespread health conditions in the world, and the molecular mechanism of it is still unclear. In this study, we identified the hub genes (hub miRNA genes) associated with hypertension and explored the relationship between hypertension miRNA-gene by constructing a mRNA co-expression network and a miRNA co-expression network, which can help to reveal the mechanism and predict the prognosis of hypertension progression. MATERIAL AND METHODS Based on gene expression profile data of hypertensive samples from the Gene Expression Omnibus database, WGCNA was used to detect hypertension-related biomarkers and key mRNA and miRNA modules. Then, DAVID was used to perform gene-annotation enrichment analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) and miRPath were used for pathway analysis of mRNA and miRNAs genes. RESULTS We identified 3 key modules relating to hypertension, 2 mRNA modules named Msaddlebrown and Mgreenyellow and 1 miRNA module named Msalmon. In addition, 12 hub genes (RPL21, RPS28, LOC442727/PTGAP10, LOC100129599/RPS29P14, TBXAS1, FCER1G, CFP, FURIN, PECAM1, IGSF6, NCF1C, and LOC285296/UNC93B3) and 7 hub miRNAs (hsa-miR-1268a/b, hsa-miR-513c-3p, hsa-miR-4799-5p, hsa-miR-296-3p, hsa-miR-5195-5p, hsa-miR-219-2-3p, and hsa-miR-548d-5p) relating to hypertension were identified. HIF-1 signaling pathway and insulin signaling pathway were closely related to the 3 key modules. We also discovered 4 miRNAs (hsa-miR-548am-3p, hsa-miR-513c-3p, hsa-miR-182-5p, and hsa-miR-548d-5p) and 6 genes (IGF1R, GSK3B, FOXO1, PRKAR2B, HIF1A, and PIK3R1) were the core nodes in the hypertension-related miRNA-gene network, and hsa-miR-548am-3p was at the center of the network. CONCLUSIONS These findings will help improve the understanding of the pathogenesis of hypertension, and the discovered genes can serve as signatures for early diagnosis of hypertension.

Background

Hypertension, also known as high blood pressure, is one of the most widespread health conditions in the world. It is one of the most dangerous factors affecting cardiovascular death. Every year, more than 17 million people die from cardiovascular disease. Hypertension is responsible for more than 45% of deaths due to heart disease, and 51% of deaths due to stroke [1]. The risk factors for hypertension include age, ethnicity, weight, diet, alcohol and tobacco use, gender, as well as existing health condition such as diabetes, high cholesterol levels, and chronic kidney disease. The pathogenic mechanism is even more complex, involving a variety of molecules and pathways [2].

Currently, hypertensive patients are treated with drugs such as diuretics, vasodilators, rapid-acting intravenous antihypertensive agents, and other drugs. Unfortunately, hypertension medications have many side effects such as cough, diarrhea, dizziness, drowsiness, depression, ulcers, and more. Since hypertension is a multiple-factor disease, further research is needed to reveal the molecular mechanism of hypertension and new biomarkers need to be discovered. The increased availability of high-throughput sequencing and microarrays along with the development of bioinformatics tools and algorithms has allowed for the discovery of several biomarkers that regulate blood pressure. For example, the anti-inflammatory cytokine interleukin-10 (IL-10) has a strong antihypertensive effect [3]; ENOS gene can increase plasma nitric oxide levels to reduce blood pressure [4]; and renin-angiotensin system (RAS) improves insulin resistance and prevents the development of renal hypertension [5]. Many miRNAs related to hypertension have also been discovered [6], such as miR-155 [7], miRNA-126 [8], miR-124 [9], and miR-150 [10]. Most of these studies only focused on one single gene or miRNA, and the identification of these gene targets are limited to just differential expression. Very few studies focused on expression profiles of multiple genes, and there is also insufficient attention to the high degree of interconnection between genes.

The network of interactions between biomolecules provides an important basis for systematic research on disease. WGCNA is an R package, which is based on the similarity between genes to construct a weighted correlation network [11]. It has unique advantages in handling complex data with multiple samples, which is a powerful method to uncover basic mechanism of gene–disease relationships [12]. Using WGCNA, Zhang et al. discovered ten hub genes that could be used as biomarkers for oral squamous cell carcinoma tumors [13]. In another study, six hub genes were found to regulate the signaling pathway of clear cell renal cell carcinoma (ccRCC) [14]. Wu et al. applied WGCNA to the identification of potential therapeutic targets for angiotensin II (Ang II) induced hypertension [7]. Here, we utilized the WGCNA method to construct gene and miRNA modules, identified new biomarkers related to hypertension, and explored the relationship between hypertension and marker genes.

Material and Methods

Data sources and preprocessing

GSE75360, GSE75670, and GSE117261 datasets were downloaded from Gene Expression Omnibus (GEO) (). GSE75360 contained 10 hypertensive and 11 normal human gene expression data, (Illumina HumanHT-12 v.4.0 Expression BeadChip). GSE75670 contained 6 hypertensive and 6 normal human miRNA expression data (Exiqon mercury™ LNA™ microRNA array, 7th generation [miRbase v18]). GSE117261 contained 58 pulmonary arterial hypertension and 25 normal human gene expression data ([HuGene-1_0-st] Affymetrix Human Gene 1.0 ST Array [transcript (gene) version]). Data preprocessing procedures such as correction of expression matrix, quality evaluation of expression data, and sample clustering were conducted in R version 3.5.2. Network analysis with WGCNA software package was also conducted in R version 3.5.2. The overall analysis workflow is showed in Figure 1.

Figure 1

The workflow of this study.

Module construction base on WGCNA algorithm

WGCNA explores the complex relationships between genes and phenotypes by constructing scale-free co-expression networks. It transforms gene expression data into co-expression modules and then identifies hub genes in the modules. In this research, we constructed 2 weighted co-expression networks based on mRNA and miRNA expression data separately. The construction processes were the same except that some parameter values were different, so we only introduced the construction process of gene co-expression network. First, the correlation of all gene pairs was calculated to construct a similarity matrix. Second, the soft threshold β was a weighted parameter of the adjacent function, which the optimal value was obtained by the pickSoftThreshold function in the R package WGCNA [15]. Third, the TOM similarity function was used to convert the adjacency value into a TOM matrix. Then, using dissimilarity matrix dissTOM=1-TOM, we clustered the genes into the hierarchy to get the system clustering tree [13]. Fourth, mRNAs with similar expression profile were divided to the same module. According to the number of the genes and miRNAs, the minModuleSize of the mRNA was set to 50 and the minModuleSize of the miRNA was set to 30 [16]. Finally, we calculated the differences of the modules eigengenes, and set an appropriate cutline for the modules dendrogram to merge highly similar modules.

Identification of key modules in co-expressed networks

Two methods were used to identify key modules related to hypertension. The first method was to calculate the Pearson correlation coefficient and significance P-value of module eigengenes (MEs) and hypertension trait. Here, MEs represents the overall expression level of the gene module [17]. The second method calculates the gene significance (GS) and the module significance (MS). Here, GS is the correlation between a gene and the clinical features. MS is the average GS of all genes in a module [13]. Generally, the higher the absolute value of MS and GS, the more relevant the gene module is to hypertension.

Enrichment analysis of key modules

To further understand the function of a key module and its biological significance, we used the online functional annotation database DAVID (). For the mRNA modules, we used Gene Ontology (GO) terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis to identify significantly enriched pathways (P-value <0.05). For miRNA modules, the online database miRPath V.3 () was used to predict miRNAs target genes for KEGG pathway enrichment [18].

Identify hub genes and hub miRNAs

This study used 2 methods to identify hub genes (miRNA) from key modules: 1) Importance threshold, and 2) MCC (maximal clique centrality function) algorithm which used the cytoHubba plugin. In the first method, hub genes are defined as genes with high GS value, high modular membership (MM) value, and low weighted P-value associated with genes and hypertension (P. weighted). MM was used to measure the importance of genes in modules [19]. The P-weighted value was calculated by the networkScreening function in WGCNA. A P. weighted value less than 0.05 was considered biologically significant [20]. In the second method, we used the MCC algorithm in cytoHubba plug-in from the Cytoscape software to identify hub genes [21].

Gene selection and hub gene verification

We used limma R package to analyze differently expressed genes (DEGs) between normal samples and hypertension samples in the dataset GSE75360 and set the cutoff value to log2FC< |0.182| and P-value <0.05. The volcanic map and hierarchical clustering analysis were performed by “ggplot2” and “pheatmap” package of R, respectively. Then, we used the jvenn () to draw Venn diagrams to overlap the genes in DEGs and hub genes [22]. We verified our hub genes using a different DEG dataset (GSE117261). We queried the role of hub genes through The Human Protein Atlas database () and used NCBI () to verify whether these hub genes could be used as biomarkers of hypertension.

Construction of miRNA-gene interactive network

DAVID and miRPath were used to identify enriched mRNA and miRNA pathways, respectively. Then, the overlapping pathways were chosen for network construction. Cytoscape was used to construct miRNA-gene interaction networks.

Results

Data preprocessing

After preprocessing, an expression matrix containing 29 595 genes was obtained. We calculated standard deviation (SD) of all genes, then, ranked the SD from large to small, and selected the top 6000 genes as the input data for the construction of the gene expression network. We performed the same preprocessing procedures on the miRNA expression matrix. A total of 1916 miRNAs were selected for co-expression network construction of miRNAs.

Construction of weighted co-expression network

To obtain a network that meets the scale-free topology criterion, we calculated network structures by using different soft-thresholding power range from 1 to 20. The scale-free topological fitting index of mRNA co-expression network reaches 0.8, when the soft-thresholding power was 14 (Supplementary Figure 1A, 1B), which met the scale-free network criterion. Then, the genes were clustered into modules by hierarchical clustering according to expression value, and the most similar modules were merged by setting the MEDissThres cutting line to 0.2 (Supplementary Figure 1C). Finally, 19 mRNA gene modules were identified (Figure 2). For miRNA co-expression networks, the soft-thresholding power was set at 9 and the scale-free topological fitting index at 0.88 (Supplementary Figure 2A, 2B). The MEDissThres cutting line was set to 0.1 (Supplementary Figure 2C), and 14 miRNA gene modules were identified (Figure 3). The relationship for each mRNA module was analyzed by drawing a network heat map (Supplementary Figure 3). The network heat map was produced with the modules in the miRNA co-expression network (Figure 4). The different colors in the vertical and horizontal axes stand for different miRNA modules. The yellow color in the middle area indicates a degree of connection for each miRNA module. The figure shows no significantly differences in the interaction between the modules, indicating that the miRNA modules had relatively high independence.

Figure 2

The cluster dendrogram of mRNA in mRNA expression data, each branch represents a gene, and each color below represents a co-expression module. The first ribbon represents the module detected by dynamic tree cutting, and the second ribbon represents the module after merging the similar module.

Figure 3

The cluster dendrogram of miRNAs in miRNA expression data.

Figure 4

The interaction relationship of co-expressed miRNAs. The different colors in the vertical and horizontal axes stand for different miRNAs modules. The yellow color in the middle area indicated a degree of connection for each miRNA module.

Identifying of key modules

We used 2 methods to identify key modules. The first method calculate the Pearson correlation and significance P-value between MEs and hypertension, and the other method calculates GS. We found 2 mRNA modules with significant correlations to hypertension, named Msaddlebrown and Mgreenyellow. Compared to other mRNA modules, the correlation coefficient of theese two modules were the largest (Figure 5A). In the miRNA co-expression network, the Msalmon module was highly related to hypertension (Figure 5B). In order to ensure the identified modules were significantly associated with hypertension, we calculated the GS in the modules and verify the key modules again (Figure 6A, 6B). For the mRNA co-expression network, the Msaddlebrown and Mgreenyellow modules had the largest GS scores at 0.451 and 0.410, respectively. For miRNA co-expression networks, the absolute GS value of Msalmon was 0.398. These GS values all indicate that these modules were significantly associated with hypertension.

Figure 5

(A) mRNA module-trait relationship. The MESaddlebrown module was most significantly related to hypertension, and the MEGreenyellow module was the second one. (B) miRNA module-trait relationship. The MESalmon module was most significantly related to hypertension.

Figure 6

(A) Module significance of mRNA. Distribution of average mRNA significance and errors in the modules related to hypertension. Salmon module is associated with high blood pressure. (B) Module significance of mRNA.

We performed GO analysis and KEGG analysis on the genes in the Msaddlebrown and Mgreenyellow modules to explore their biological significance (Supplementary Tables 1, 2). The Msaddlebrown module was related to RNA transcriptional translation, such as SRP-dependent co-translational protein targeting to membrane, translational initiation, and RNA binding. The Mgreenyellow module was involved in a variety of immune and metabolic processes, such as positive regulation of IL-6 production, positive regulation of transcription from RNA polymerase II promoter, glycoside catabolic process, and so on. Furthermore, KEGG pathway analysis showed enrichment in ribosome, HIF-1 signaling pathway, and osteoclast differentiation pathways.

The significant pathways in Msalmon module were Hippo signaling pathway, adherens junction, and proteoglycans in cancer (Supplementary Table 3). Among them, the HIF-1 signaling pathway and the insulin signaling pathway were shared by both mRNA and miRNA modules.

Identification of hub genes and miRNAs

According to the definition of module connectivity, we calculated the MM and GS of the genes (miRNAs) in each of the key modules to select the hub genes (miRNAs). Then, we used the networkScreening function to obtain the P. weighted of each gene (miRNA). |MM| >0.8, |GS| >0.2 and P. weighted <0.05 were used as the identification criteria. Finally, 26 genes were obtained in the Msaddlebrown module (Supplementary Table 4), 53 genes were obtained in the Mgreenyellow module (Supplementary Table 5), and 22 miRNAs were obtained in the Msalmon module (Supplementary Table 6). We also import the files of these 3 key modules into Cytoscape software and used the plugin cytoHubba to identify and visualize pivotal genes and miRNAs (Supplementary Figure 4A–4C). We finally determine that there were 4 hub genes in Msaddlebrown module, 8 hub genes in Mgreenyellow module, and 7 hub miRNAs in Msalmon module (Table 1 and Supplementary Figure 5A–5C).

Table 1

Hub genes (MiRNAs) in key modules.

Modules	Name	DEGs	GS	MM	Log2FC
Msaddlebrown	RPS28	Down	−0.65	0.884	−0.282
Msaddlebrown	RPL21	Down	−0.555	0.863	−0.290
Msaddlebrown	LOC442727(PTMAP10)	Down	−0.527	0.895	−0.345
Msaddlebrown	LOC100129599(RPS29P14)	Down	−0.441	0.848	−0.303
Mgreenyellow	TBXAS1	Down	−0.517	0.926	−0.316
Mgreenyellow	FCER1G	Down	−0.597	0.922	−0.478
Mgreenyellow	CFP	Down	−0.550	0.867	−0.347
Mgreenyellow	FURIN	Down	−0.544	0.872	−0.316
Mgreenyellow	PECAM1	Down	−0.487	0.844	−0.243
Mgreenyellow	IGSF6	Down	−0.485	0.860	−0.412
Mgreenyellow	NCF1C	Down	−0.450	0.879	−0.354
Mgreenyellow	LOC285296(UNC93B3)	Down	−0.451	0.852	−0.174
Msaddlebrown	hsa-miR-1268a/hsa-miR-1268b	Up	0.416	0.987	0.550
Msaddlebrown	hsa-miR-513c-3p	Up	0.354	0.984	0.641
Msaddlebrown	hsa-miR-4799-5p	Up	0.392	0.975	1.054
Msaddlebrown	hsa-miR-296-3p	Up	0.417	0.974	0.408
Msaddlebrown	hsa-miR-5195-5p	Up	0.429	0.972	0.915
Msaddlebrown	hsa-miR-219-2-3p	Up	0.451	0.972	0.746
Msaddlebrown	hsa-miR-548d-5p	Up	0.399	0.968	0.901

The verification and functional analysis of hub genes

Since Msaddlebrown and Mgreenyellow modules were negatively correlated with hypertension status, we wondered if the hub genes in these 2 modules were also negatively correlated with hypertension. We verified this by correlating hub genes expression values with hypertension status. There was a significant difference in hub genes (miRNAs) expression levels between normal and hypertensive group Supplementary Figure 6A–6C. We detected the expression of 12 hub genes by DEGs on original dataset GSE75360. Figure 7 shows the DEGs in a volcano map. The hierarchical clustering heat map of DEGs is shown in Supplementary Figure 7. We used the jvenn tool to overlap the hub genes in the Msaddlebrown and Mgreenyellow module and DEGs, respectively, and found 11 hub genes were both in the DEGs gene list and the hub genes lists (Supplementary Figure 8A, 8B). We used the same method to verify the hub gene in another dataset GSE117261 from the GEO database, and found TBXAS1, FCER1G, and IGSF6 were in the DEGs gene list and hub genes lists (Supplementary Figure 9). We compared the expression status of these three hub genes in normal and hypertensive patients from the GSE117261 dataset, and the results were consistent with those from the GSE75360 dataset (Supplementary Figure 10).

Figure 7

Volcano map of DEGs for GSE75360 dataset. The red dots on the right corresponds to 2-fold up changes with P-value less than 0.05, and the blue dots on the left means 2-fold down changes with P-value less than 0.05.

Analysis of miRNA-gene interaction networks

We discovered that the HIF-1 signaling pathway and the insulin signaling pathway were also found in the miRNA networks. In order to better understand the regulatory relationship between genes and miRNAs, the miRNA-gene interaction network was constructed, which was based on genes and miRNAs involving the same pathway. There was a total of 46 nodes (21 genes and 25 miRNAs) and 112 pairs of interactions in the miRNA-gene interaction network (Figure 8). We used the MCC algorithm in cytoHubba plugin to screen the top 10 miRNA-genes in the network. Four miRNAs (hsa-miR-548am-3p, hsa-miR-513c-3p, hsa-miR-182-5p, and hsa-miR-548d-5p) and 6 genes (IGF1R, GSK3B, FOXO1, PRKAR2B, HIF1A, and PIK3R1) were the core nodes of the network (Table 2, Figure 9), and hsa-miR-548am-3p was considered as the core regulator because it targeted these 6 genes.

Figure 8

The regulatory network of miRNAs and target genes. Among them, the green prisms are mRNA, the purple triangles are miRNA, pink solid lines represent the HIF-1 pathway, and the green dotted lines represent the insulin pathway.

Table 2

Top 10 in network MiRNA-gene ranked by MCC method.

Rank	Name	Score	DEGs
1	hsa-miR-548am-3p	14	Up
2	hsa-miR-513c-3p	11	Up
2	hsa-miR-182-5p	11	Up
4	hsa-miR-548d-5p	8	Up
4	IGF1R	8	Down
6	GSK3B	7	Up
7	FOXO1	6	Up
7	PRKAR2B	6	Up
7	HIF1A	6	Down
7	PIK3R1	6	Up

miRNA – microRNA; MCC – Maximal Clique Centrality function.

Figure 9

Top 10 nodes in miRNA-gene network.

Discussion

In this study, mRNA and miRNA co-expression networks of hypertensive patients samples were constructed by using the WGCNA method. Out of the 19 identified mRNA modules, Msaddlebrown and Mgreenyellow had the most significant correlations to hypertension. For miRNAs, we identified 14 modules, of which Msalmon was the most significantly associated with hypertension. Finally, we identified 4 hub genes in the Msaddlebrown module, 8 hub genes in the Mgreenyellow module, and 7 hub miRNAs in the Msalmon module that were correlate with hypertension. Four miRNAs and 6 genes were also associated with the genetic susceptibility to hypertension. Our findings help us better understand the pathogenesis of hypertension, which in turn will provide us candidate biomarkers for clinical decision-making, potential therapeutic targets for accurate diagnosis, and treatment targets of hypertension.

In the Msaddlebrown module, the GO analysis indicates that Msaddlebrown mainly refered to membrane formation, translation, and ribosome and rRNA processing. For the Mgreenyellow module, the insulin signaling pathway and the HIF-1 signaling pathway were enriched. The HIF-1 signaling pathway is reported to be important for development of pulmonary hypertension in chronic hypoxia [23]. The HIF-1 pathway is related to proliferation of pulmonary arterial smooth muscle cells (PASMCs), which is also a central pathological component for a kind of hypertension [24]. We found 12 hub genes associated with hypertension (RPL21, RPS28, LOC442727/PTGAP10, LOC100129599/RPS29P14, TBXAS1, FCER1G, CFP, FURIN, PECAM1, IGSF6, NCF1C, and LOC285296/UNC93B3). Other studies have shown that some of these hub genes are related to hypertension, such as TBXAS1, FCER1G, FURIN, and PECAM1. The enzyme encoded by the TBXAS1 gene can catalyze many reactions involving cholesterol, steroids, drug metabolism and other lipid synthesis; and a study has reported that TBXAS1 is a potent vasoconstrictor [25]. FCER1G is the Fc fragment of FcepsilonRI (IgE) receptor Ig and IgE expression on the cell surface of human platelets. It has been reported that FCER1G has a regulatory effect on diabetic kidneys, and as hypertension is closely related to this disease, it can be inferred that it has a certain effect on hypertension [26]. FURIN encodes proteases, and Li et al. reported that FURIN may be a candidate gene for human hypertension [27]. PECAM1 encoded immunoglobulin are involved in processes of leukocyte migration and angiogenesis, and studies have shown that PECAM-1 expression is upregulated in pregnancy-induced hypertensive patients [28].

Some hypertension-related pathways were found in the Msalmon module, which indicates that miRNA modules may also play a role in hypertension. The pathways include the adipocytokine signaling pathway, the thyroid hormone signaling pathway, the insulin signaling pathway, and the HIF-1 signaling pathway. The 3 key modules (Msaddlebrown, Mgreenyellow and Msalmon) have significant relations to hypertension, so we speculate that there may be regulatory relationships of the miRNAs and genes in the same pathway that affect blood pressure. In order to find out their relationship, we constructed a miRNA-gene regulatory network (Figure 8). We found 10 core nodes of the network (hsa-miR-548am-3p, hsa-miR-513c-3p, hsa-miR-182-5p, hsa-miR-548d-5p, IGF1R, GSK3B, FOXO1, PRKAR2B, HIF1A, and PIK3R1), and hsa-miR-548am-3p was considered the core regulator. Although these 4 miRNAs have not been reported to regulate hypertension, their target genes have blood pressure regulation effect. It is found that IGF1R, the target gene of hsa-miR-548am-3p, hsa-miR-513c-3p, hsa-miR-182-5p, and hsa-miR-548d-5p, may be associated to the genetic susceptibility to hypertension [29]. Due to its function in lung development, IGF-1R may play a role in the pathogenesis of pulmonary hypertension [30]. GSK3B, a target gene of hsa-miR-548am-3p, hsa-miR-513c-3p and hsa-miR-182-5p, encodes for the enzyme glycogen synthase kinase 3 beta (GSK3B). Platelet-derived factor (PDGF) can affect the abnormal growth of pulmonary arterial smooth muscle cells (PASMC) in pulmonary hypertension by inhibiting β-catenin (βC) activation of GSK3B [31]. FOXO1, a target gene of hsa-miR-548am-3p, hsa-miR-513c-3p, hsa-miR-182-5p, and hsa-miR-548d-5p, can control angiotensinogen (AGT) and Ang II levels and regulate blood pressure by regulating liver gene expression [32]. PRKAR2B, a target gene of hsa-miR-548am-3p and hsa-miR-548d-5p, may be a candidate gene for spontaneously hypertensive rats [33]. HIF1A, a target gene of hsa-miR-548am-3p, hsa-miR-513c-3p, and hsa-miR-182-5p, plays a vital role in the pathophysiology of embryonic angiogenesis and ischemic diseases. Sheng et al. found that gene polymorphism of HIF1A was associated with left ventricular hypertrophy in essential hypertension [34]. PIK3R1, a target gene of hsa-miR-548am-3p, hsa-miR-182-5p, and hsa-miR-548d-5p, plays a vital role in the metabolism of insulin, and PIK3CA/PIK3R1 may be involved in chronic thromboembolic pulmonary hypertension pathophysiology [35]. These results indicate that these 4 miRNAs may influence the occurrence of hypertension by regulating related target genes.

This study had several limitations. First, our sample size was small and may not fully represent hypertension patients. Second, there is no biological experimental verification of the hub genes (miRNAs) and relationship between genes and miRNA. In a follow-up study, the molecular verification experiment will be conducted to uncover the molecular level mechanisms of miRNA-gene interactions and their relations to hypertension. Verified molecular markers can be used as new diagnostic indexes of hypertension in the future.

Conclusions

This study established a WGCNA-based gene expression data process workflow, identified 2 mRNA modules and 1 miRNA module related to hypertension, and provided potential candidate biomarkers for hypertension treatment. Our analysis revealed novel miRNA-gene interactions as well as central miRNAs and genes that play critical roles in hypertension.

(A) Determining the soft-thresholding power in WGCNA. The left graph is a scale-free fitting index for analyzing various soft-thresholding powers (β). The graph on the right analyzes the average connectivity of various soft-thresholding powers. (B) The histogram verifies the selected β value is approached without scale. When the logarithm (log(k)) of the number of k-nodes is negatively correlated with the logarithm of the probability of occurrence of the node (log(p(k))), and the correlation coefficient is greater than 0.8, it is shown that the selected soft threshold accords with the standard. The left graph is the histogram of connectivity distribution when β is 14. The right image shows the scale-free topology checked when β=14. (C). The cluster dendrogram of gene modules eigengenes.

(A) The average connection between the scale-free fitting index of various soft threshold power and all kinds of soft threshold power in miRNA expression data. The left graph is a scale-free fitting index for analyzing various soft-thresholding powers (β). The graph on the right analyzes the average connectivity of various soft-thresholding powers. (B) Verify whether the selected β value is close to scale-free. When β is 9, the scale-free topological fitting index reaches 0.89, which meets the scale-free network standard. (C) The cluster dendrogram of miRNA modules eigengenes.

Analysis of the interaction relationship of co-expressed genes. The different colors of the horizontal axis and vertical axis represent different modules. The yellow brightness in the middle indicates the degree of connection of different modules. There is no significant difference in the interaction between the different modules, indicating a high degree of independence between these modules.

(A) Diamonds represent genes. From red to yellow, the top 10 hub genes in the saddlebrown module are ordered in descending order. (B) Diamonds represent genes. From red to yellow, the top 10 hub genes in the greenyellow module are ordered in descending order. (C) Diamonds represent miRNAs. From red to yellow, the top 10 hub miRNAs in the salmon module are ordered in descending order.

(A) Venn diagram of hub genes in saddlebrown module. (B) Venn diagram of hub genes in green-yellow module. (C) Venn diagram of hub miRNAs in salmon module.

(A) The expression status of the hub genes in saddlebrown is negatively correlated with hypertension, indicating that it plays an important role in inhibiting the occurrence of hypertension, and the results shown in the figure. The expression status of the hub genes in greenyellow is negatively correlated with hypertension, indicating that it plays an important role in inhibiting the occurrence of hypertension, and the results shown in the figure are in accordance with the results of WGCNA. (B) The expression status of the hub genes in greenyellow is negatively correlated with hypertension, indicating that it plays an important role in inhibiting the occurrence of hypertension, and the results shown in the figure are in accordance with the results of WGCNA. (C) The expression status of the hub miRNA in salmon is negatively correlated with hypertension, indicating that it plays an important role in promote the occurrence of hypertension, and the results shown in the figure are in accordance with the results of WGCNA.

Heat map hierarchical clustering reveals the comparison between high blood pressure samples and normal samples in DEGs of the GSE75360 dataset.

(A) The common genes between DEGs and saddlebrown module were screened by Venn. It was found that 4 hub genes in saddlebrown were also in DEGs. (B) The common genes between the DEGs and greenyellow module were screened by Venn. Seven of the 8 hub genes in the greenyellow module are in DEGs.

The common genes between the DEGs of the GSE117261 dataset and hub genes were screened by Venn hub genes in the DEGs.

Expression of 3 hub genes in the GSE117261 datasets.

Supplementary Table 3

The results of salmon module analysis in mirPath.

KEGG pathway	p-value	#genes	#miRNAs
Hippo signaling pathway	2.18E-07	94	28
Adherens junction	3.19E-07	55	26
Proteoglycans in cancer	3.19E-07	123	29
Pathways in cancer	1.90E-06	229	33
Prostate cancer	3.25E-06	64	27
TGF-beta signaling pathway	4.56E-06	51	26
Pancreatic cancer	1.62E-05	48	24
Transcriptional misregulation in cancer	1.62E-05	105	29
Protein processing in endoplasmic reticulum	2.64E-05	106	27
Chronic myeloid leukemia	3.12E-05	50	26
FoxO signaling pathway	8.58E-05	83	29
Cell cycle	9.30E-05	76	28
Ubiquitin mediated proteolysis	0.000137428	89	28
ErbB signaling pathway	0.000220758	56	23
AMPK signaling pathway	0.000220758	77	28
Hepatitis B	0.000220758	81	28
Prion diseases	0.000234754	14	12
Renal cell carcinoma	0.000234754	45	24
Fatty acid biosynthesis	0.000328055	8	11
Wnt signaling pathway	0.000328055	84	26
Endocytosis	0.000328055	120	31
Colorectal cancer	0.000538812	43	24
Focal adhesion	0.00055548	122	27
Thyroid hormone signaling pathway	0.000591287	71	25
Neurotrophin signaling pathway	0.000591287	72	27
Bacterial invasion of epithelial cells	0.000592348	48	24
Shigellosis	0.001159846	42	25
SNARE interactions in vesicular transport	0.001507172	22	14
Acute myeloid leukemia	0.001507172	37	22
Non-small cell lung cancer	0.001507172	35	22
Viral carcinogenesis	0.001664854	104	31
Glycosaminoglycan biosynthesis – chondroitin sulfate/dermatan sulfate	0.001674489	11	11
Glioma	0.00182823	40	23
Axon guidance	0.00228799	73	27
Adipocytokine signaling pathway	0.002582229	46	20
mRNA surveillance pathway	0.002582229	56	24
Choline metabolism in cancer	0.003052312	63	26
Endometrial cancer	0.00337732	34	24
Arrhythmogenic right ventricular cardiomyopathy (ARVC)	0.004005275	38	21
Small cell lung cancer	0.00695403	53	22
p53 signaling pathway	0.00695403	42	26
Melanoma	0.010615755	42	26
PI3K-Akt signaling pathway	0.010615755	182	29
Alanine, aspartate and glutamate metabolism	0.013556846	23	18
Sphingolipid signaling pathway	0.013556846	67	25
Signaling pathways regulating pluripotency of stem cells	0.019376413	78	29
Insulin signaling pathway	0.020843228	80	24
Rap1 signaling pathway	0.022783789	113	28
Fc gamma R-mediated phagocytosis	0.032659471	53	20
MAPK signaling pathway	0.034939635	133	31
RNA degradation	0.039468488	49	20
RNA transport	0.039695937	94	28
HIF-1 signaling pathway	0.045306902	60	27

Supplementary Table 4

26 hub genes identified in the saddlebrown.

Gene name	Module color	GS	MM	p.MMsaddlebrown
LOC346950	Saddlebrown	−0.601197709	0.957602673	1.01E-11
LOC730288	Saddlebrown	−0.604978411	0.948872979	5.79E-11
LOC651453	Saddlebrown	−0.459749264	0.942082904	1.84E-10
LOC729255	Saddlebrown	−0.587440957	0.922040399	2.86E-09
LOC441743	Saddlebrown	−0.478152959	0.912936642	7.89E-09
LOC644790	Saddlebrown	−0.579349333	0.912306001	8.43E-09
LOC730382	Saddlebrown	−0.565669291	0.900406106	2.69E-08
LOC442727	Saddlebrown	−0.526685302	0.89531027	4.24E-08
LOC648343	Saddlebrown	−0.487547164	0.889430022	6.95E-08
LOC442454	Saddlebrown	−0.615757481	0.888692889	7.38E-08
RPS28	Saddlebrown	−0.646553637	0.884044065	1.07E-07
LOC645630	Saddlebrown	−0.63719721	0.874165775	2.23E-07
RPL21	Saddlebrown	−0.55460734	0.86291539	4.81E-07
LOC100134273	Saddlebrown	−0.461686834	0.861381943	5.31E-07
LOC653156	Saddlebrown	−0.531697938	0.855435749	7.73E-07
LOC100129599	Saddlebrown	−0.440650986	0.847964345	1.21E-06
LOC100130154	Saddlebrown	−0.648033957	0.846382684	1.33E-06
LOC389156	Saddlebrown	−0.576581109	0.837823455	2.14E-06
IL27RA	Saddlebrown	−0.566595063	0.837574429	2.17E-06
LOC645968	Saddlebrown	−0.574222298	0.833838557	2.66E-06
RPL12P6	Saddlebrown	−0.370854211	0.819547605	5.48E-06
LOC647673	Saddlebrown	−0.550047745	0.811269744	8.11E-06
HS.24119	Saddlebrown	0.422593803	−0.80822296	9.32E-06
LOC641849	Saddlebrown	−0.561382353	0.807938149	9.44E-06
LOC643997	Saddlebrown	−0.62418576	0.807894354	9.46E-06
LOC402644	Saddlebrown	−0.602398481	0.802312918	1.21E-05

Supplementary Table 5

53 hub genes identified in the greenyellow.

Gene name	Module color	GS	MM	p.MM greenyellow
PGD	Greenyellow	−0.434129475	0.946648977	8.60E-11
APLP2	Greenyellow	−0.613094499	0.935412596	5.05E-10
TYROBP	Greenyellow	−0.501678945	0.930806437	9.55E-10
TBXAS1	Greenyellow	−0.517207792	0.926143763	1.74E-09
FCER1G	Greenyellow	−0.596814877	0.921776439	2.95E-09
FKBP1A	Greenyellow	−0.57900758	0.914880662	6.41E-09
LOC25845	Greenyellow	0.439117962	−0.907524776	1.37E-08
NUP214	Greenyellow	−0.584252139	0.90255565	2.21E-08
SKAP1	Greenyellow	0.44444412	−0.900542702	2.66E-08
LOC642489	Greenyellow	−0.547397605	0.896794749	3.72E-08
CST3	Greenyellow	−0.479727023	0.89467232	4.48E-08
LTBR	Greenyellow	−0.50621654	0.892023335	5.61E-08
FCGRT	Greenyellow	−0.548226646	0.889646993	6.83E-08
TTYH3	Greenyellow	−0.453836688	0.880705335	1.38E-07
NCF1C	Greenyellow	−0.450320897	0.879107249	1.56E-07
PYCARD	Greenyellow	−0.454380899	0.875035746	2.10E-07
ESYT1	Greenyellow	0.514042042	−0.873855037	2.28E-07
RXRA	Greenyellow	−0.671014834	0.872732002	2.47E-07
ZNF792	Greenyellow	0.463899499	−0.87219588	2.57E-07
FURIN	Greenyellow	−0.54355238	0.871954617	2.61E-07
SYK	Greenyellow	−0.460231959	0.870122146	2.97E-07
FUCA2	Greenyellow	−0.453218553	0.868278067	3.37E-07
CFP	Greenyellow	−0.550302824	0.866886875	3.70E-07
RAB5C	Greenyellow	−0.554564791	0.864131533	4.44E-07
C9ORF167	Greenyellow	−0.365162049	0.864061526	4.46E-07
SH3TC1	Greenyellow	−0.418683416	0.863783292	4.55E-07
IGSF6	Greenyellow	−0.484564612	0.859764188	5.89E-07
ACP6	Greenyellow	0.610300842	−0.85351514	8.69E-07
CYBB	Greenyellow	−0.567450199	0.852636639	9.17E-07
LOC100133163	Greenyellow	0.516405892	−0.852177584	9.43E-07
LOC285296	Greenyellow	−0.45061083	0.851785134	9.65E-07
LOC730278	Greenyellow	−0.702260803	0.851209949	9.99E-07
RGS19	Greenyellow	−0.555367503	0.847795031	1.22E-06
PECAM1	Greenyellow	−0.486962729	0.84365182	1.55E-06
ARID3A	Greenyellow	−0.462869481	0.843573809	1.56E-06
ZNF385A	Greenyellow	−0.50600348	0.841357803	1.76E-06
PPPDE2	Greenyellow	−0.542019307	0.84101991	1.80E-06
LRP1	Greenyellow	−0.623124649	0.838515852	2.06E-06
IFI30	Greenyellow	−0.595101224	0.83487417	2.51E-06
ARF3	Greenyellow	−0.389712065	0.833677878	2.68E-06
CTSD	Greenyellow	−0.528731344	0.830185491	3.22E-06
SIGLEC9	Greenyellow	−0.485469865	0.824853637	4.22E-06
PHF19	Greenyellow	0.447166744	−0.822797953	4.67E-06
C14ORF139	Greenyellow	0.410259281	−0.822707045	4.69E-06
C15ORF39	Greenyellow	−0.428338711	0.822289038	4.79E-06
UBTD1	Greenyellow	−0.416538289	0.82137832	5.01E-06
DAPK1	Greenyellow	−0.655158887	0.819826079	5.41E-06
CYFIP2	Greenyellow	0.475140293	−0.819355168	5.53E-06
LOC653888	Greenyellow	−0.36054457	0.818955519	5.64E-06
LOC644086	Greenyellow	−0.451322436	0.817870278	5.94E-06
ZNF827	Greenyellow	0.53343424	−0.814712324	6.90E-06
HMOX1	Greenyellow	−0.478808329	0.806051764	1.03E-05
LOC100129201	Greenyellow	−0.531458889	0.801345505	1.27E-05

Supplementary Table 6

22 hub genes identified in the salmon.

miRNA name	Module color	GS	MM	p.MMsalmon
hsa-miR-1268a/hsa-miR-1268b	Salmon	0.416064753	0.986610639	3.31E-09
hsa-miR-513c-3p	Salmon	0.353693204	0.983548774	9.23E-09
hsa-miR-4676-5p	Salmon	0.427676991	0.977111887	4.76E-08
hsa-miR-3065-5p	Salmon	0.372723442	0.97539399	6.82E-08
hsa-miR-4799-5p	Salmon	0.391832663	0.975301669	6.94E-08
hsa-miR-296-3p	Salmon	0.417125099	0.973985636	8.98E-08
hsa-miR-4452	Salmon	0.456458843	0.972033422	1.29E-07
hsa-miR-5195-5p	Salmon	0.428595201	0.971735808	1.35E-07
hsa-miR-219-2-3p	Salmon	0.450679649	0.971591899	1.39E-07
hsa-miR-548d-5p	Salmon	0.39856976	0.96823523	2.41E-07
hsa-miR-3120-5p	Salmon	0.429196273	0.966398946	3.19E-07
hsa-miR-4478	Salmon	0.505013687	0.964646031	4.10E-07
hsa-miR-4789-3p	Salmon	0.459045717	0.955529254	1.27E-06
hsa-miR-548am-3p	Salmon	0.45331186	0.955059255	1.34E-06
hsa-miR-500a-5p	Salmon	0.393282326	0.952195235	1.81E-06
hsa-miR-4712-3p	Salmon	0.427439039	0.950729799	2.10E-06
hsa-miR-3119	Salmon	0.450893848	0.950205284	2.22E-06
hsa-miR-1264	Salmon	0.566131244	0.947371832	2.91E-06
hsa-miR-212-5p	Salmon	0.483676354	0.932827254	9.62E-06
hsa-miR-548am-5p/hsa-miR-548au-5p/hsa-miR-548c-5p/hsa-miR-548o-5p	Salmon	0.469740795	0.904479189	5.32E-05
hsa-miR-4710	Salmon	0.5188767	0.888272911	0.000113346
hsa-miR-4768-5p	Salmon	0.529008968	0.870622705	0.000228833