Prediction of marker genes associated with hypertension by bioinformatics analyses.

This study aimed to explore the underlying marker genes associated with hypertension by bioinformatics analyses. A gene expression profile (GSE54015) was downloaded. The differentially expressed genes (DEGs) between the normotensive female (NF) and hypertensive female (HF), and between the normotensive male (NM) and hypertensive male (HM) groups were analyzed. Gene Ontology (GO) and pathway enrichment analyses were performed, followed by protein-protein interaction (PPI) network construction. The transcription factors (TFs), and the common DEGs between the HF and HM groups were then analyzed. In total, 411 DEGs were identified between the HF and NF groups, and 418 DEGs were identified between the HM and NM groups. The upregulated DEGs in the HF and HM groups were enriched in 9 GO terms, including oxidation reduction, such as cytochrome P450, family 4, subfamily b, polypeptide 1 (Cyp4b1) and cytochrome P450, family 4, subfamily a, polypeptide 31 Cyp4a31). The downregulated DEGs were mainly enriched in GO terms related to hormone metabolic processes. In the PPI network, cytochrome P450, family 2, subfamily e, polypeptide 1 (Cyp2e1) had the highest degree in all 3 analysis methods in the HF group. Additionally, 4 TFs were indentified from the 2 groups of data, including sterol regulatory element binding transcription factor 1 (Srebf1), estrogen receptor 1 (Esr1), retinoid X receptor gamma (Rxrg) and peroxisome proliferator-activated receptor gamma (Pparg). The intersection genes were mainly enriched in GO terms related to the extracellular region. On the whole, our data indicate that the DEGs, Cyp4b1, Cyp4a31 and Loxl2, and the TFs, Esr1, Pparg and Rxrg, are associated with the progression of hypertension, and may thus serve as potential therapeutic targets in this disease.

Introduction

Hypertension, also known as arterial hypertension, is a chronic medical condition in which the blood pressure in the arteries is elevated (1). Hypertension is classified as either primary or secondary hypertension, and primary hypertension accounts for approximately 90–95% of cases (2). In primary hypertension, the cause of the hypertensioin is not clearly defined, whereas in secondary hyptertension, the cause of the high blood pressure is clearly determined. Hypertension can lead to cardiac hypertrophy, aortic dissection and heart failure, due to structural changes in blood vessel walls and in the heart muscle, which has become a major burden to healthcare systems (3). Thus, hypertension is the leading risk factor for cardiovascular disease (4). According to a report of the World Health Organization in 2011, cardiovascular diseases are the leading cause of morbidity and mortality worldwide (5).

As a multifactorial disease, hypertension results from a complex interaction of environmental factors and genetic factors (6). Numerous common genetic variants, as well as some rare genetic variants have been identified to affact blood pressure (7,8). For instance, bone morphogenetic protein receptor type 2, which belongs to the transforming growth factor (TGF)-β superfamily of receptors, has been found to be mutated in the majority of cases of hypertension (9). Interleukin (IL)-17 has also been thought to be involved in hypertension due to its role in increasing the production of some other immune system chemical signals, such as IL-1, IL-6 and tumor necrosis factor (TNF)-α (10). Additionally, Austin et al (11) suggested that mutations in caveolin-1 are associated with pulmonary arterial hypertension. Despite advances in understanding the etiology of hypertension, the genetic basis responsible for this disease remains poorly elucidated.

Presently, there are many species of animals used as models of hypertension. In our study, the microarray dataset, GSE54015, which contained data from Schlager BPH/2J and BPN/3J mouse models was downloaded for analysis. of note, Puig et al (12) also used the same mouse models to identify potential pharmacological targets of hypertension. In this study, we identified the differentially expressed genes (DEGs) between normotensive mice and hypertensive mice of both genders. Additionally, we performed functional enrichment analyses and protein-protein interaction (PPI) networks analysis to study and identify the pathways and DEGs associated with hypertension. Furthermore, transcription factors (TFs) and DEGs in both hypertensive female (HF) and male (HM) groups were analyzed. The results of this study may aid in the diagnosis and treatment of hypertension.

Data collection methods

Affymetrix microarray data

The microarray data of GSE54015 were downloaded, which was deposited in the National Center of Biotechnology Information (NCBI) Gene Expression Omnibus (GEO) (http://www.ncbi.nlm.nih.gov/geo/) database by Chiu et al (2). The platform is [MoGene-1_1-st] Affymetrix Mouse Gene 1.1 ST array [transcript (gene) version]. In the study by Chiu et al (2), 12- to 13-week-old Schlager hypertensive male and hypertensive female BPH/2J, and normotensive male and normotensive female BPN/3J mice (n=6/group) were used, and their kidney tissues were collected for the microarrays. In the present study, 4 types of mouse samples were used for analysis: normotensive female (NF) BPN/3J, hypertensive female (HF) BPH/2J, normotensive male (NM) BPN/3J and hypertensive male (HM) BPH/2J mice.

Data pre-processing and analysis of DEGs

The CEL raw data were downloaded and preprocessed using oligo (13) in R of Bioconductor (http://www.bioconductor.org/packages/release/bioc/html/) with the method of robust multiarray average (RMA) (14). The preprocessing processes included background correction, quartile data normalization and probe summarization.

The paired t-test based on the limma package (15) (no. 1406, Ritchie, 2015) in Bioconductor was used to identify the DEGs between the HF and NF groups, and between the HM and NM groups. The false discovery rate (FDR) was calculated to adjust the P-value with the Benjamini-Hochberg method (16). Only genes with |log2FC| ≥0.58 and an adjusted P-value <0.05 were regarded as DEGs.

Gene ontology (GO) and pathway enrichment analyses

The GO (http://www.geneontology.org) project is used to provide consistent descriptors for gene products in different databases, and to standardize classifications for sequences (17). Kyoto Encyclopedia of Genes and Genomes (KEGG) (http://www.genome.ad.jp/kegg/) is a knowledge base used to analyze gene functions, and to link genomic information with higher order functional information (18). In this study, the obtained DEGs were performed GO and KEGG pathway enrichment analyses to obtain the biological process (BP), molecular function (MF), cellular component (CC) terms and pathways. The analysis tool was the Database for Annotation, Visualization and Integrated Discovery (DAVID) (http://david.niaid.nih.gov) (19), which is an integrated biological knowledgebase used for the systematic and integrative analysis of large gene protein lists. The P-value <0.05 and enrichment count ≥2 were considered as thresholds.

PPI network construction

The Search Tool for the Retrieval of Interacting Genes (STRING) (http://string-db.org/) database (20) is a precomputed global resource used for evaluating PPI information. In the present study, the STRING online tool was applied to predict the PPI pairs of DEGs with a combined score of >0.4. The PPI network was constructed using cytoscape (21). Additionally, 3 methods including degree centrality, betweenness centrality and subgraph centrality were used to analyze the hub gene in the PPI network using the CytoNCA cytoscape plugin (http://apps.cytoscape.org/apps/cytonca) (22).

TF analysis

In order to predict the role of TFs in hypertension, we selected the TFs from the DEG corresponding proteins, and predicted the target associations between the TFs and the DEGs based on the TRANSFAC database in the iRegulon cytoscape plugin (http://iregulon.aertslab.org) (23), as previously described. FDR on motif similarity <0.05 was regarded as the threshold.

Analysis of DEGs in the HF and HM groups

In the present study, we screened the common DEGs between the HF and HM groups in order to analyze the similarity of hypertension in females and males. Subsequently, these selected DEGs were subjected to GO functional enrichment analysis using DAVID. Furthermore, the enriched GO terms were subjected to correlation analysis using the enrichmentmap cytoscape plugin (http://baderlab.org/Software/EnrichmentMap/) (24) (P<0.05).

Results

Identification of DEGs

In total, 201 upregulated and 210 downregulated DEGs were identified between the HF and NF groups. Moreover, 204 upregulated and 214 downregulated DEGs were identified between the HM and NM groups. The results are presented in the heatmaps in Fig. 1A and B.

Figure 1

Heatmap plot of differentially expressed genes accross all samples. (A) Hypertensive female (HF) vs. normotensive female (NF); (B) hypertensive male (HM) vs. normotensive male (NM). Red stands for high expression value and green stands for low expression value. The changes of color from green to red stands for the changes of expression value from low to high.

GO and pathway enrichment analyses of DEGs

The upregulated DEGs in the HF group were mainly enriched in GO terms related to ion transport and oxidation reduction, and in the pathway of histidine metabolism. The downregulated DEGs were mainly enriched in the GO terms related to steroid metabolic process and microsome, and in the pathway of arachidonic acid metabolism (Table I).

Table I

Functional enrichment analysis for the DEGs of the HF group.

		Term	Count	P-value
Up	BP	GO:0006811-ion transport	19	5.68E-05
		GO:0019751-polyol metabolic process	4	6.56E-03
		GO:0055114-oxidation reduction	13	1.66E-02
	CC	GO:0005576-extracellular region	32	4.98E-04
		GO:0042383-sarcolemma	4	1.17E-02
		GO:0005886-plasma membrane	41	1.75E-02
	MF	GO:0022836-gated channel activity	10	1.00E-03
		GO:0005216-ion channel activity	10	4.35E-03
		GO:0022838-substrate specific channel activity	10	5.32E-03
	KEGG	mmu00340:Histidine metabolism	3	1.96E-02
		mmu03320:PPAR signaling pathway	4	3.10E-02
Down	BP	GO:0008202-steroid metabolic process	9	7.92E-05
		GO:0055114-oxidation reduction	17	2.25E-04
		GO:0010565-regulation of cellular ketone metabolic process	4	9.96E-04
	CC	GO:0005792-microsome	9	2.16E-04
		GO:0042598-vesicular fraction	9	2.71E-04
		GO:0005783-endoplasmic reticulum	19	6.14E-04
	MF	GO:0030414-peptidase inhibitor activity	11	2.89E-06
		GO:0004866-endopeptidase inhibitor activity	10	1.06E-05
		GO:0004857-enzyme inhibitor activity	11	4.83E-05
	KEGG	mmu00590:Arachidonic acid metabolism	5	3.56E-03
		mmu00511:Other glycan degradation	3	6.44E-03
		mmu00830:Retinol metabolism	4	1.49E-02
		mmu04610:Complement and coagulation cascades	4	1.94E-02

DEGs, differentially expressed genes; HF, hypertensive female; up, upregulated DGEs; down, downregulated DEGs; BP, biological process; CC, cellular component; MF, molecular function; KEGG, Kyoto Encyclopedia of Genes and Genomes.

In the HM group, the upregulated DEGs were mainly enriched in GO terms related to co-factor metabolic process and oxidation reduction, and in the pathway of PPAR signaling. The downregulated DEGs were mainly enriched in GO terms related to hormone metabolic process and endoplasmic reticulum (Table II).

Table II

Functional enrichment analysis for the DEGs of the HM group.

		Term	Count	P-value
Up	BP	GO:0051186-cofactor metabolic process	8	1.04E-03
		GO:0055114-oxidation reduction	15	2.04E-03
		GO:0019748-secondary metabolic process	5	3.02E-03
	CC	GO:0005576-extracellular region	37	4.68E-06
		GO:0005777-peroxisome	5	2.03E-02
		GO:0042579-microbody	5	2.03E-02
	MF	GO:0019842-vitamin binding	10	1.64E-06
		GO:0016290-palmitoyl-CoA hydrolase activity	3	8.12E-03
		GO:0030170-pyridoxal phosphate binding	4	1.33E-02
	KEGG	mmu03320:PPAR signaling pathway	5	6.25E-03
		mmu00590:Arachidonic acid metabolism	5	7.43E-03
Down	BP	GO:0042445-hormone metabolic process	7	1.63E-04
		GO:0010817-regulation of hormone levels	7	1.07E-03
		GO:0042398-cellular amino acid derivative biosynthetic process	4	9.31E-03
	CC	GO:0005783-endoplasmic reticulum	15	1.29E-02
	MF	GO:0004866-endopeptidase inhibitor activity	7	3.52E-03
		GO:0030414-peptidase inhibitor activity	7	5.43E-03
		GO:0004867-serine-type endopeptidase inhibitor activity	5	2.09E-02

DEGs, differentially expressed genes; HM, hypertensive male; up, upregulated DGEs; down, downregulated DEGs; BP, biological process; CC, cellular component; MF, molecular function; KEGG, Kyoto Encyclopedia of Genes and Genomes.

The constructed PPI network based on the DEGs between the HF and NF groups is shown in Fig. 2A. The PPI network consisted of 193 nodes and 301 interaction pairs. The top 10 hub genes were selected based on the analysis of degree centrality, betweenness centrality and subgraph centrality. Among the 10 genes, cytochrome P450, family 2, subfamily E, polypeptide 1 (Cyp2e1) had the highest degree in all the 3 analysis methods (Table III).

Figure 2

The constructed protein-protein interaction (PPI) network for differentially expressed genes (DEGs). Pink stands for upregulated DEGs and green stands for downregulated DEGs. (A) Hypertensive female (HF) vs. normotensive female (NF); (B) hypertensive male (HM) vs. normotensive male (NM).

Table III

The top 10 network nodes in the HF group based on the analysis of degree centrality, betweenness centrality and subgraph centrality.

Gene	Degree		Betweenness		Subgragh
Cyp2e1	18	Cyp2e1	7915.4395	Cyp2e1	295.83835
Cyp4a12a	16	Aldh1a1	4667.3003	Cyp4a12a	279.15607
Cyp2b10	11	Pla2g5	4572.402	Cyp2b10	179.15228
Pla2g5	9	Srebf1	4370.598	Slc22a30	100.235176
Srebf1	9	Dkk1	4108.5645	Pla2g5	99.77931
Serpina1d	9	Cyp4a12a	3427.9216	Slco1a1	99.73588
Cyp7b1	9	Pld1	3171.1284	Cyp4f18	81.0327
Slc22a30	9	Serpina1d	3132.0713	Cyp4f14	81.03269
Anxa1	8	Prok1	3062.5642	Aldh1a1	73.60265
Gm266	8	Plcl2	2823.1924	Aldh1a7	69.66662

HF, hypertensive female.

The PPI network constructed by the DEGs between HM and NM included 182 nodes and 248 PPI interaction pairs (Fig. 2B). The top 10 hub genes respectively analyzed by 3 methods are shown in Table IV, namle UDP glucuronosyltransferase 2 family, polypeptide B37 (Ugt2b37), formimidoyltransferase cyclodeaminase (Ftcd) and cytochrome P450, family 4, subfamily a, polypeptide 14 (Cyp4a14).

Table IV

The top 10 network nodes in the HM group based on the analysis of degree centrality, betweenness centrality and subgraph centrality.

Gene	Degree		Betweenness		Subgragh
Ugt2b37	11	Ftcd	5603.764	Cyp4a14	57.048008
Cyp4a10	10	Lep	4363.4756	Cyp4a10	55.98465
Cyp4a14	10	Ugt2b37	4147.9277	Ppbp	52.882313
Lep	10	Serpina1d	4005.5542	Ugt2b37	50.838398
Ppbp	9	Cyp4a14	3500.8042	Aldh1a7	34.876118
Tgfb2	8	Fcgr2b	2844.4468	Aldh1a1	34.272568
Ftcd	8	Pld1	2784.8557	Aplnr	28.300505
Aldh1a7	7	Plcl2	2771.22	Grm7	28.300499
Aldh1a1	7	Il10rb	2728.78	Adra2b	28.300493
A2m	7	Synj2	2569.02	Oxgr1	28.300491

HM, hypertensive male.

Based on the TRANSFAC database in iRegulon, the TFs in the HF and HM samples were predicted and the transcriptional regulatory networks were constructed (Fig. 3A and B). In the HF group, the TFs were sterol regulatory element binding transcription factor 1 (Srebf1) and estrogen receptor 1 (Esr1). In the HM group, the TFs were retinoid X receptor gamma (Rxrg) and peroxisome proliferator-activated receptor gamma (Pparg).

Figure 3

The constructed transcriptional regulatory network. Octagon stands for transcription factor; circle stands for differentially expressed genes (DEGs) in hypertensive female (HF); square stands for DEGs in hypertensive male (HM). Pink stands for upregulated DEGs and green stands for downregulated DEGs. (A) HF group; (B) HM group.

There were 118 intersection genes among the upregulated DEGs in the HF and HM samples; they were mainly distributed in chromosomes 4 and 7. Among the downregulated DEGs, 112 were intersection genes and they were mainly distributed in chromosomes 11 and 14. Moreover, these intersection genes were mainly enriched in BP terms related to catecholamine biosynthetic process, CC terms related to extracellular region and MF terms related to endopeptidase inhibitor activity (Table V). For instance, lysyl oxidase-like 2 (Loxl2) was enriched in the CC term related to extracellular region. The interaction network for these terms are shown in Fig. 4.

Table V

The significant functions enriched by the intersection genes.

	Term	Count	P-value
BP	GO:0042423-catecholamine biosynthetic process	3	0.004975549
	GO:0006638-neutral lipid metabolic process	4	0.007615151
	GO:0042445-hormone metabolic process	5	0.008732351
CC	GO:0005576-extracellular region	31	4.56E-04
	GO:0000267-cell fraction	12	0.02798221
	GO:0005626-insoluble fraction	11	0.030455631
MF	GO:0004866-endopeptidase inhibitor activity	9	1.17E-04
	GO:0030414-peptidase inhibitor activity	9	2.16E-04
	GO:0004867-serine-type endopeptidase inhibitor activity	7	6.59E-04

BP, biological process; CC, cellular component; MF, molecular function.

Figure 4

The interaction network for the terms enriched by intersection genes. Node stands for Gene Ontology (GO) term, the more the number of gene that was enriched in term, the bigger the node. Edge stands for the interaction between terms, the more the number of intersection genes, the thicker the edge.

Discussion

In the present study, 411 DEGs were identified between the HF and NF samples, and 418 DEGs were identified between the HM and NM samples. The upregulated DEGs, Cyp4b1 and Cyp4a31, were mainly enriched in BP terms related to oxidation reduction in both the HF and HM samples. In addition, Esr1, Pparg and Rxrg were the key TFs involved in the transcriptional regulatory networks. Following the analysis of DEGs in the HF and HM samples, many intersection DEGs were identified between the 2 groups, such as Loxl2. These genes and TFs may play important roles in the progression of hypertension.

Previous studies have indicated that reactive oxygen species play an important pathophysiological role in the development of cardiovascular diseases, including hypertension (25,26). In human hypertension, biomarkers of systemic oxidative stress have been found to be elevated (27). In the present study, BP terms related to oxidation reduction were significant terms enriched by the DEGs in both the HF and HM samples. For instance, several cytochrome P450 (CYP)-associated genes, including Cyp4b1 and Cyp4a31 were enriched in this BP term. CYPs are important enzymes involved in the oxidative and reductive metabolism of numerous endogenous compounds in mammals (28). In the study by Chiu et al (2), Cyp4b1 was validated to be differentially expressed between hypertensive and normotensive mice by qPCR. CYP4A proteins have been found in vascular tissues (29) and have been suggested to be linked to hypertension in numerous experimental models (30). Sacerdoti et al (31) demonstrated that the inhibition or depletion of CYP4A activity can lower blood pressure in the spontaneously hypertensive rat. Of note, Cyp4a31 was found enriched in the PPAR signaling pathway in HF and HM samples. The PPARs are members of the nuclear receptor family of ligand-activated TFs, which play an important role in the regulation of cardiovascular circadian rhythms (32). Therefore, the upregulation of Cyp4b1 and Cyp4a31 in our study was in accordance with previous research. The 2 DEGs may be critical biomarkers associated with human hypertension.

In addition, Esr1 was a key TF in the transcriptional regulatory networks of the HF sample group. Esr1 is a member of the family of ESRs, belonging to superfamily of nuclear receptor TFs. This receptor is expressed in a wide range of tissues and cells, including vascular smooth muscle and vascular endothelial cells (33). Estrogen can increase vasodilatation and inhibit the development of atherosclerosis (33). Specifically, the incidence of cardiovascular disease, including hypertension differs between males and females (34). Accumulating evidence suggests that estrogen deficiency plays a major role in the pathogenesis of cardiovascular diseases in post-menopausal women (35). In our study, the TF, Esr1, was also identified in female hypertensive mice, which indicated that Esr1 may play an important role in the progression of hypertension in females.

Furthermore, the TFs, Pparg and Rxrg, were involved in the transcriptional regulatory networks in the HM samples. Pparg is a member of PPARs which, as mentioned above, play an important role in the regulation of cardiovascular circadian rhythms. Pparg has been implicated in atherosclerosis (36), coronary heart disease (37) and hypertension (38). RXRs are members of the NR2B family of nuclear receptors and are major heterodimerization partners of PPARs (39). RXRs are involved in mediating the anti-proliferative effects of retinoic acid, which can comprise a collection of molecules with retinoids to maintain homeostasis in adults (40). A growing number of studies have reported the role of RXRs and PPARs in controlling cellular metabolism in tissues, including the myocardium (40,41). Nohara et al (39) suggested that a variant of Rxrg contributes to human genetic dyslipidemia. Taken together, the TFs, Pparg and Rxrg, may have utility in the diagnose and treatment of hypertension.

Following the analysis of DEGs in the HF and HM samples, we found that there were many intersection DEGs between the 2 groups, such as Loxl2. This DEG was found enriched in the BP terms related to oxidation reduction and CC terms related to extracellular region. Loxl2 encodes a member of the lysyl oxidase gene family, which are extracellular copper-dependent amine oxidases initiating the formation of crosslinks in collagen and elastin (42). The collagen and elastin are important components of the extracellular matrix (43). It has been found that the activity of lysyl oxidase is essential to maintain the features of connective tissues, including cardiovascular systems (44). The deregulation of lysyl oxidase is associated with aneurysms, aortic rupture and different stages of the atherosclerotic process (45,46). Presently, there are few studies on the association between lysyl oxidase and hypertension. Thus, we speculated that Loxl2 may be a potential biomarker for the progression of hypertension.

In conclusion, our study provides a comprehensive bioinformatics analysis of biomarkers in hypertension, which may provide new insight into the mechanisms responsible for this disease. The DEGs, Cyp4b1, Cyp4a31 and Loxl2, and the TFs, Esr1, Pparg and Rxrg, may be associated with the progression of hypertension, and may serve as potential therapeutic targets in this disease.

However, there was a limitation to this study. The present study was implemented by means of bioinformatics; therefore, the screened genes and TFs have not been validated by actual experiments. As a result, further experimental studies with a larger sample size are warranted to confirm the results of our study.