Bioinformatic gene analysis for possible biomarkers and therapeutic targets of hypertension-related renal cell carcinoma.

BACKGROUND: Renal cell carcinoma (RCC) is one of the most prevalent malignant tumors of the urinary system. Hypertension can cause hypertensive nephropathy (HN). Meanwhile, Hypertension is considered to be related to kidney cancer. We analyzed co-expressed genes to find out the relationship between hypertension and RCC and show possible biomarkers and novel therapeutic targets of hypertension-related RCC.
METHODS: We identified the differentially expressed genes (DEGs) of HN and RCC through analyzing Gene Expression Omnibus (GEO) datasets GSE99339, GSE99325, GSE53757 and GSE15641 by means of bioinformatics analysis, respectively. Then we evaluated these genes with protein-protein interaction (PPI) networks, Gene Ontology (GO) analysis, Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis and CTD database. Subsequently, we verified co-expressed DEGs with Gene Expression Profiling Interactive Analysis (GEPIA) database. Finally, corresponding predicted miRNAs of co-expressed DEGs were identified and verified via mirDIP database and Starbase, respectively.
RESULTS: We identified 9 co-expressed DEGs, including BCAT1, CORO1A, CRIP1, ESRRG, FN1, LYZ, PYCARD, SAP30, and PTRF. CRIP1 and ESRRG and their corresponding predicted miRNAs, especially hsa-miR-221-5p, hsa-miR-205-5p, hsa-miR-152-3p and hsa-miR-137 may be notably related to hypertension-related RCC.
CONCLUSIONS: CRIP1 and ESRRG genes have great potential to become novel biomarkers and therapeutic targets concerning hypertension-related RCC. 2020 Translational Andrology and Urology. All rights reserved.

Introduction

Renal cell carcinoma (RCC) is a common subtype of malignant tumor in the urinary system. Its incidence rate was 11.281 per 100,000 person-years. The incidence of RCC in men is higher than that in women (1). Around the world, the incidence of RCC is higher in developed countries than in developing countries (2). According to cancer statistics 2020, in 2020 about 73,750 (4.1%) newly diagnosed kidney cancer cases are expected and around 14,830 (2.4%) will die of this cancer (3). In China, the incidence of RCC is the second highest in malignant tumors of the urinary system. Its incidence and mortality vary greatly from region to region and are higher in cities than in rural areas (4). Clear cell renal cell carcinoma (ccRCC) is the most important type in RCC. According to the statistics, it accounts for around 80–90% of RCC (5). Although the early detection of RCC has been improved in recent years, more than one third of patients have local advanced tumor or metastatic disease at the time of diagnosis (6). The treatment remains a challenge and reliable biomarkers are essential to prevent metastasis and improve the quality of patients’ life (7). Although some possible biomarkers have been found such as Ki-67, p53, VEGF, there is currently no reliable biomarker for RCC prediction (8). RCCs are usually resistant to conventional chemotherapy and almost all chemotherapeutic agents are ineffective against metastatic RCC (9). The use of targeted therapy is contributed to surgery in patients with locally advanced or metastatic RCC. Although some therapeutic targets have been applied, such as VEGFR and HIF-1α, our exploration of therapeutic targets is far from enough (7).

Hypertension is a widespread chronic disease with an increased incidence worldwide (10). Hypertension is an important risk factor for kidney cancer (11). Hypertensive nephropathy (HN) is one of the main complications of hypertension. Over the past few years, the relationship between hypertension and the risk of kidney cancer have been explored by many prospective studies (12-14). In the VITAL study, hypertension was independently associated with RCC risk (HR 1.70; 95% CI, 1.30–2.22) (15). A meta-analysis of prospective studies shows a strong positive correlation between hypertension and kidney cancer. A history of hypertension was associated with 67% increased risk of RCC. Considering heterogeneity and publication bias, each 10 mmHg increase in blood pressure was associated with 10–22% increased risk of RCC (16). However, the molecular mechanism of the relationship between kidney cancer and hypertension is not clear. For these reasons, the identification of key molecular involved in hypertension-related kidney cancer is urgent and highly demanded to improve the clinical outcome. In our research, HN-related DEGs (HN-DEGs) and RCC-related DEGs (RCC-DEGs) were identified by bioinformatic analysis. Then co-expressed DEGs (co-DEGs) of HN and RCC were found. Furthermore, further analysis and verification of DEGs and predicted targeted miRNAs were conducted for HN patients inclined to RCC to find possible molecular mechanisms. We present the following article in accordance with the MDAR reporting checklist (available at http://dx.doi.org/10.21037/tau-20-817).

Methods

Gene expression profiles data

We downloaded GSE99339, GSE99325, GSE53757 and GSE15641 datasets from GEO (http://www.ncbi.nlm.nih.gov/geo/). The expression profiling arrays of GSE99339 were generated applying GPL19109 [HG-U133_Plus_2] Affymetrix Human Genome U133 Plus 2.0 Array [CDF: Brainarray HGU133Plus2_Hs_ENTREZG_v18], including 15 Hypertensive Nephropathy specimens and 14 Tumor Nephrectomy specimens. The expression profiling arrays of GSE99325 were generated using GPL19184 [HG-U133A] Affymetrix Human Genome U133A Array [Custom Brainarray v18 ENTREZG CDF], including 20 Hypertensive Nephropathy specimens and 4 Cadaveric Donor specimens. Next, we used GPL570 [HG-U133_Plus_2] Affymetrix Human Genome U133 Plus 2.0 Array to generate the expression profiling arrays of GSE53757, including 72 ccRCC specimens and 72 normal specimens. Moreover, we also used GPL96 [HG-U133A] Affymetrix Human Genome U133A Array to generate the expression profiling arrays of GSE15641, including 49 RCC tissues and 23 normal tissues. We used the two gene expression profiles of GSE99339 and GSE99325 to filter differentially expressed genes (DEGs) of HN. Similarly, the two gene expression profiles of GSE53757 and GSE15641 were applied to filter DEGs of RCC. All procedures performed in this study were in accordance with the ethical standards of the institutional and/or national research committee and with the Helsinki Declaration (as revised in 2013).

Data processing

We used R packages of “affy”, “affyPLM”, and “limma” (http://www.bioconductor.org/packages/release/bioc/html/affy.html) to assess GSE99339, GSE99325, GSE53757 and GSE15641 RAW datasets. Background correction, probe summarization, quantile normalization, and log2-transformation were used to generate a robust multi-array average (RMA), a log-transformed mismatch, and perfect match probe methods. We applied the Benjamini-Hochberg method to adjust original P values and the false discovery rate (FDR) procedure to calculate fold-changes (FC). We used Genes expression values of the|log2 FC| >1.5 or <0.667 and adjusted P<0.05 to screen HN-DEGs. Moreover, we used the |log2 FC| >2 or <0.5 and adjusted P<0.05 to filter RCC-DEGs. Furthermore, volcano plots and Venn diagrams were made for co-DEGs of HN- and RCC-DEGs.

Analysis of protein-protein interaction (PPI) networks

We used the search tool for the retrieval of interacting genes (STRING) database (V11; http://string-db.org/) to analyze PPI networks of HN- and RCC-DEGs. Analytic data of the STRING database were downloaded with a combined score >0.4. Then we applied Cytoscape software (V3.5.1; http://cytoscape.org/) to analyze and visualize node degrees and biological networks.

Functional enrichment analysis

We conducted Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of HN- and RCC-DEGs with the database for annotation, visualization and integrated discovery (DAVID) bioinformatics resources (http://david.abcc.ncifcrf.gov/).We identified the significantly enriched GO terms and KEGG pathway maps related to biofunctions based on a P<0.05. Furthermore, we identified enriched functions of HN- and RCC-DEGs in molecular functions, biological processes, and cellular components, respectively.

Afterward, we used the AmiGO database (v2.0; http://amigo.geneontology.org/amigo/) to analyze the GO consortium of chosen co-DEGs in order to check the annotate and accuracy biofunctions of confirmed co-DEGs.

Identification of co-DEGs related to hypertension or renal cancer

We used the comparative toxicogenomics database (CTD, http://ctdbase.org/) to describe chemical-gene/protein interactions and chemical-disease and gene-disease relationships in order to know their associations. These data were used to identify the genes related to hypertension or renal cancer and scored genes indirectly through the CTD.

Co-DEGs validation

We used the Gene Expression Profiling Interactive Analysis (GEPIA) database (http://gepia.cancer-pku.cn/) to verify the differential expression of these genes in ccRCC. Overall survival curve was used to analyze survival differences.

Functional and pathway enrichment related to predicted miRNAs and Co-DEGs and validation of predicted miRNAs

We applied microRNA Data Integration Portal (mirDIP) (http://ophid.utoronto.ca/mirDIP) to predict miRNAs that may be regulated by 9 genes. Five top candidate miRNAs were determined according to predicted scores for each Co-DEG. We used Diana-miRPath (v3.0; http://www.microrna.gr/miRPathv3) to conduct GO terms and pathway enrichment analysis of these predicted miRNAs. Subsequently, the online tool Starbase (http://starbase.sysu.edu.cn/) was used to identify the relationship between predicted miRNAs and corresponding genes.

Statistical analysis

Data of DEGs were analyzed using R 4.0.0. Student’s t-test was performed for comparisons between two groups, whereas ANOVA was performed for repeated measures. The χ2 test was used to analyze the statistical significance of GO and pathway. Differences with P<0.05 were considered statistically significant.

Data of Co-DEGs validation were based on GEPIA. We chose the dataset of KIRC. One-way ANOVA was used for differential analysis. The expression data were log2(TPM+1) transformed for differential analysis and the log2FC was defined as median (Tumor) – median (Normal). The screening criterion was P value <0.05. In survival analysis, we used Log-rank test for hypothesis test (P<0.05). We calculated the hazards ratio (HR) based on Cox PH Mode and added the 95% CI as dotted line.

We used Pearson correlation analysis to identify the relationship between predicted miRNAs and Corresponding genes based on Starbase. The expression data of ccRCC were from TCGA. P<0.01 was considered statistically significant.

Results

Identification of DEGs

21,592 probes in GSE99339 and GSE99325 datasets and 76,937 probes in GSE53757 and GSE15641 were identified. Then we confirmed HN- and RCC-DEGs. In GSE99325, we identified 1,147 DEGs in HN tissue specimens compared with normal tissue specimens, including 585 upregulated genes and 144 downregulated genes (). Similarly, we also described the DEGs of the other three datasets in the form of volcano plots (). The common 60 genes of GSE99339 and GSE99325 were confirmed, including 38 upregulated and 22 downregulated genes, which are the DEGs of HN (). In the same way, the common 679 genes of GSE53757 and GSE15641 are the DEGs of RCC (). Heatmaps of HN-DEGs in relation to cellular response to tumor necrosis factor, response to lipopolysaccharide, multicellular organismal homeostasis and neutrophil degranulation were conducted for genes expression (). The value of RCC-DEGs expression concerning extracellular structure organization, leukocyte migration, neutrophil activation and response to oxygen levels has been shown in .

Figure 1

Volcano plots and Venn diagrams. (A-D) The volcano plots image the differentially expressed genes (DEGs) in GSE99325, GSE99339, GSE15641, and GSE53757. The blue and red dots represent downregulated and upregulated genes, respectively. (E) Venn diagrams of hypertensive nephrology (HN) related DEGs in GSE99325 and GSE99339. (F) Venn diagrams of renal cell carcinoma (RCC) related DEGs in GSE15641 and GSE53757.The purple and blue graphics represent downregulated genes. The grey and orange graphics represent upregulated genes.

Figure 2

Heatmaps of hypertensive nephrology (HN) related differentially expressed genes (DEGs). (A-D) Results of DEGs expression concerning cellular response to tumor necrosis factor, response to lipopolysaccharide, multicellular organismal homeostasis, and neutrophil degranulation. Red, greater expression. Blue, less expression.

Figure 3

Heatmaps of renal cell carcinoma (RCC) related differentially expressed genes (DEGs). (A-D) Results of DEGs expression concerning extracellular structure organization, leukocyte migration, neutrophil activation and response to oxygen levels. Red, greater expression. Blue, less expression.

Functional enrichment in Co-DEGs

illustrates 9 Co-DEGs of HN- and RCC-DEGs, including Sin3A associated protein 30 (SAP30), polymerase I and transcript release factor (PTRF), lysozyme (LYZ), PYD and CARD domain containing (PYCARD), estrogen related receptor gamma (ESRRG), coronin 1A (CORO1A), fibronectin 1 (FN1), branched chain amino acid transaminase 1 (BCAT1), cysteine rich protein 1 (CRIP1). We confirmed GO term enrichment concerning biological processes, molecular functions, and cellular components with AmiGO database. We found that Co-DEGs were related to many processes (Table S1).

Figure 4

Venn diagrams and protein-protein interaction (PPI) network. (A) Common differentially expressed genes (DEGs) between renal cell carcinoma (RCC) and hypertensive nephrology (HN); (B) PPI network of DEGs of HN; (C) PPI network of DEGs of RCC. Purple, greater degree. Blue, lesser degree.

PPI network analysis, GO analysis and KEGG pathway enrichment analysis

Fifty-eight and 618 nodes have been identified from PPI network of HN- and RCC-DEGs, respectively (). The hub nodes including fibronectin 1 (FN1; degree =21), vascular cell adhesion molecule 1 (VCAM1; degree =14), alpha 2-HS glycoprotein (AHSG; degree =12), nuclear receptor subfamily 1 group H member 4 (NR1H4; degree =10), PPARG coactivator 1 alpha (PPARGC1A; degree =10), coronin 1A (CORO1A; degree =10) and activating transcription factor 3 (ATF3; degree =10) are hub genes of HN. Similarly, the hub genes including kininogen 1 (KNG1; degree =58), complement C3 (C3; degree =55), fibronectin 1 (FN1; degree =38), TIMP metallopeptidase inhibitor 1 (TIMP1; degree =36), albumin (ALB; degree =36), prolyl 4-hydroxylase subunit beta (P4HB; degree =34) are illustrated in RCC-DEGs with a high degree.

We identified the GO terms involved in biological processes among HN-DEGs by using the DAVID database. They were mainly associated with cellular response to tumor necrosis factor (Fold Enrichment: 7.88; P value: 2.82E-05), response to tumor necrosis factor (Fold Enrichment: 7.35; P value: 4.4E-05), response to lipopolysaccharide (Fold Enrichment: 6.95; P value: 6.26E-05), response to molecule of bacterial origin (Fold Enrichment: 6.68; P value: 7.98E-05) and cellular response to lipopolysaccharide (Fold Enrichment: 6.68; P value: 7.99E-05). There is evident correlation in secretory granule lumen (Fold Enrichment: 6.47; P value: 3.21E-04), cytoplasmic vesicle lumen (Fold Enrichment: 6.14; P value: 4.22E-04), vesicle lumen (Fold Enrichment: 6.12; P value: 4.29E-04) and cytosolic part (Fold Enrichment: 7.00; P value: 7.30E-04) in connection with cellular components. Furthermore, the terms concerning molecular functions were primarily referred to primary amine oxidase activity (Fold Enrichment: 107.25; P value: 1.41E-04), transcription coregulator activity (Fold Enrichment: 4.58; P value: 3.22E-04) and oxidoreductase activity, acting on the CH-NH2 group of donors, oxygen as acceptor (Fold Enrichment: 40.22; P value: 0.001). RCC was also analyzed. The biological processes terms including leukocyte migration (Fold Enrichment: 4.19; P value: 2.13E-24), extracellular structure organization (Fold Enrichment: 4.19; P value: 1.28E-20), neutrophil activation (Fold Enrichment: 3.56; P value: 4.95E-18) and response to oxygen levels (Fold Enrichment: 4.16; P value: 8.77E-20) were significantly enriched. Similarly, the terms of collagen-containing extracellular matrix (Fold Enrichment: 5.60; P value: 5.73E-29), extracellular matrix (Fold Enrichment: 4.41; P value: 9.82E-27) and endoplasmic reticulum lumen (Fold Enrichment: 4.69; P value: 7.46E-19) concerning cellular components were primarily enriched. In addition, the molecular functions terms of antigen binding (Fold Enrichment: 4.66; P value: 1.38E-13), peptide antigen binding (Fold Enrichment: 13.10; P value: 3.94E-13) and glycosaminoglycan binding (Fold Enrichment: 4.17; P value: 3.20E-12) were notably enriched ().

Figure 5

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis. (A,B) Hypertensive nephrology and renal cell carcinoma related GO analysis for differentially expressed genes (DEGs). The ordinates represent the function name of the DEGs, and the abscissas represent negative Lg-P values. (C,D) KEGG pathway of hypertensive nephrology and renal cell carcinoma related DEGs. Rectangular length and rectangular colors represent the number of DEGs contained in the pathway and P. adjust, respectively.

KEGG pathway analysis data is shown in . It indicated that the HN-DEGs were greatly enriched in pathways of tyrosine metabolism (P value: 0.006) (). However, KEGG terms including phagosome (P value: 6.96E-13), rheumatoid arthritis (P value: 1.48E-10) and viral myocarditis (P value: 5.03E-10) were enriched in RCC-DEGs ().

Co-DEGs validation and detection

The CTD showed that Co-DEGs BCAT1, CORO1A, CRIP1, ESRRG, FN1, LYZ, PYCARD, and SAP30 targeted hypertension and renal cancer and these data appear in . Compared with para cancerous normal tissues, the expression levels of 7 genes except BCAT1 were obviously higher in ccRCC in view of the GEPIA database (). Among the 8 key genes, we found that CRIP1, ESRRG, LYZ, and PYCARD were obviously associated with the overall survival of ccRCC patients ().

Figure 6

Co-expressed genes related to hypertension and kidney cancer based on the comparative toxicogenomics database (CTD). (A-H) Inference Score of BCAT1, CORO1A, CRIP1, ESRRG, FN1, LYZ, PYCARD, SAP30.

Figure 7

Gene expression and overall survival based on GEPIA database. (A-H) The gene expression levels of BCAT1, CORO1A, CRIP1, ESRRG, FN1, LYZ, PYCARD and SAP30 in normal kidney and renal clear cell carcinoma (ccRCC) tissues. (I-L) Overall survival analysis of 4 genes in ccRCC. Expression levels of CRIP1, ESRRG, LYZ and PYCARD are associated with the overall survival of patients with ccRCC based on P<0.05.

We used mirDIP and DIANA-MirPath database to predict the top 5 targeted miRNAs of each Co-DEG associated with HN-related RCC (Table S2). In , we identified the relationship between predicted miRNAs and Corresponding genes in ccRCC. Hsa-miR-429 and hsa-miR-200b-3p are negatively correlated with FN1, respectively. Hsa-miR-30e-5p and hsa-miR-30b-5p are negatively correlated with SAP30, respectively. Hsa-miR-221-5p is negatively correlated with CRIP1. Hsa-miR-205-5p, hsa-miR-152-3p and hsa-miR-137 are negatively correlated with ESRRG, respectively.

Figure 8

The relationship between predicted miRNAs and Corresponding genes. The relative expression levels between predicted miRNAs and Corresponding genes in renal clear cell carcinoma in starBase. (A,B) Hsa-miR-429 and hsa-miR-200b-3p are negatively correlated with FN1, respectively (P<0.01). (C,D) Hsa-miR-30e-5p and hsa-miR-30b-5p are negatively correlated with SAP30, respectively (P<0.01). (E) Hsa-miR-221-5p is negatively correlated with CRIP1 (P<0.01). (F-H) Hsa-miR-205-5p, hsa-miR-152-3p and hsa-miR-137 are negatively correlated with ESRRG, respectively (P<0.01).

Discussion

Numerous observational studies have systematically reported an increased risk of RCC in patients with hypertension (17,18). The biological mechanism of the relationship between hypertension and RCC is indefinite. However, it is speculated that it is related to chronic kidney hypoxia and lipid peroxidation with the generation of reactive oxygen species (ROS) (19,20). Patients with hypertension may cause chronic kidney hypoxia due to the transcription of hypoxia inducible factors, which can promote tumor cell angiogenesis and proliferation (21). we identified 8 key genes of hypertension-related RCC, including BCAT1, CORO1A, CRIP1, ESRRG, FN1, LYZ, PYCARD, and SAP30. These genes will help us further explore the mechanism of hypertension-related RCC and may become biomarkers of hypertension-related RCC. And they are also important for exploring new therapeutic targets of hypertension-related RCC.

One study showed that CRIP1 has a strong connection with pressure at the population level. The effect of CRIP1 on blood pressure may be achieved by mediating SH2B3 (22). Meanwhile, Macrophages that lack SH2B3 expression are easily activated, producing more ROS (23), which can cause kidney injury. We generally believe that SH2B3 negatively inhibits proinflammatory cell signaling within the kidney in normal and pathological states (24). Long-term effects of these mechanisms may lead to the development of hypertension-related RCC. The detection of CRIP1 helps us to detect RCC early in patients with hypertension, so as to improve the survival rate of RCC. Recent studies have revealed that ESRRG is a key transcriptional regulator of mitochondrial oxidative phosphorylation (OxPhos) and fatty acid oxidation (FAO) (25). Through the analysis of epigenomic elements of the promoter, it is found that ESRRG is a new obesity-susceptibility gene (26). It is well known that obesity and hypertension are highly correlated. Based on this, we speculate that ESRRG is associated with hypertension. ESRRG is highly expressed in the kidney and plays a major role in normal embryonic kidney development (27). The lack of ESRRG in renal epithelial cells (RECs) causes serious renal energy and absorption dysfunction and renal failure. The expression of ESRRG is positively correlated with renal function and decreases in patients with chronic kidney disease (CKD) (28). CKD is one of the causes of RCC (29). ESRRA and ESRRG drive changes in HIF1A and HIF2A, which are important in developing ccRCC molecular phenotypes (30). In addition, ESRRG also plays a role in other tumors, including breast cancer, endometrial cancer, gastric cancer, liver cancer, and prostate cancer. Further exploration of the role of ESRRG in the occurrence and development of hypertension-related RCC will help us find new therapeutic targets for RCC, thereby alleviating the situation of drug resistance in the treatment of renal cancer. FN1 is a gene that has been widely explored. The expression of FN1 is increased in patients with hypertension (31). Several studies show a strong correlation between FN1 and RCC (32,33). FN1 expression in RCC is related to a higher disease-related mortality rate, indicating a probable role in RCC progression (34). Apoptosis genes such as PYCARD have shown potential to improve the prognosis of other cancers and may be demonstrated by further research to have the same potential in RCC (35). Transcription factor SAP30 is important to activation of expression of NETO2 gene in ccRCC. Meanwhile, mRNA level of SAP30 increased significantly and was positively correlated with NETO2 gene expression (36).

In our data, miR-429, miR-200b-3p, miR-30e-5p, miR-30b-5p, miR-221-5p, miR-205-5p, miR-152-3p and miR-137 may sponge with corresponding genes and act on hypertension-related RCC. Has-miR-429 is one of the most studied miRNAs in RCC. It is reported that miR-429 suppresses tumor cell proliferation, metastasis and epithelial-mesenchymal transition by direct targeting of BMI1, E2F3 and VEGF in RCC (37,38). MiR-137 inhibits the growth and invasion of ccRCC cells, and induces apoptosis, acting as a tumor suppressor gene. MiR-137 works by targeting 3'-UTR of RLIP76 which is an oncogene identified in ccRCC (39).

Conclusions

According to existing literature reports, several genes such as CRIP1, ESRRG, FN1, PYCARD, and SAP30 are related to hypertension or RCC. Through the verification of CTD, GEPIA and Starbase databases, we found that the two genes CRIP1, ESRRG and their corresponding predicted miRNAs are most likely instructive for further exploration of hypertension-related RCC and contribute to finding new therapeutic targets.

However, some limitations also exist. primarily, this is a microarray analysis study and our data were acquired from a publicly available database. Given that gene expression is not necessarily equivalent to protein expression, we need further in vivo and in vitro experiments to clarify molecular mechanisms of key genes for clinical applications. Secondly, prospective clinical studies may validate our point of view better.

The article’s supplementary files as