Identification, Verification and Pathway Enrichment Analysis of Prognosis-Related Immune Genes in Patients With Hepatocellular Carcinoma

Abstract

Hepatocellular carcinoma is a common malignant tumor with poor prognosis, poor treatment effect, and lack of effective biomarkers. In this study, bioinformatics analysis of immune-related genes of hepatocellular carcinoma was used to construct a multi-gene combined marker that can predict the prognosis of patients. The RNA expression data of hepatocellular carcinoma were downloaded from The Cancer Genome Atlas (TCGA) database, and immune-related genes were obtained from the IMMPORT database. Differential analysis was performed by Wilcox test to obtain differentially expressed genes. Univariate Cox regression analysis, lasso regression analysis and multivariate Cox regression analysis were performed to establish a prognostic model of immune genes, a total of 5 genes (HDAC1, BIRC5, SPP1, STC2, NR6A1) were identified to construct the models. The expression levels of 5 genes in HCC tissues were significantly different from those in paracancerous tissues. The Kaplan-Meier survival curve showed that the risk score calculated according to the prognostic model was significantly related to the overall survival (OS) of HCC. The receiver operating characteristic (ROC) curve confirmed that the prognostic model had high accuracy. Independent prognostic analysis was performed to prove that the risk value can be used as an independent prognostic factor. Then, the gene expression data of hepatocellular carcinoma in the ICGC database was used as a validation data set for the verification of the above steps. In addition, we used the CIBERSORT software and TIMER database to conduct immune infiltration research, and the results showed that the five genes of the model and the risk score have a certain correlation with the content of immune cells. Moreover, through Gene Set Enrichment Analysis (GSEA) and the construction of protein interaction networks, we found that the p53-mediated signal transduction pathway is a potentially important signal pathway for hepatocellular carcinoma and is positively regulated by certain genes in the prognostic model. In conclusion, this study provides potential targets for predicting the prognosis and treatment of hepatocellular carcinoma patients, and also provides new ideas about the correlation between immune genes and potential pathways of hepatocellular carcinoma.

Keywords: bioinformatics analysis; hepatocellular carcinoma; immune genes; prognostic; signature.

Figure 1

Flow chart of this study.

Figure 2

Heat map (A) and volcano map (B) of differentially expressed genes; heat map (C) and volcano map (D) of immune differential genes; heat map (E) and volcano map (F) of differentially expressed genes related to the signal transduction pathway of P53 mediators. The abscissa of the heat map represents the sample: the blue area represents the adjacent tissue, the red area represents the hepatocellular carcinoma tissue; the ordinate represents the gene. The red dots on the volcano map represent genes whose expression levels are up-regulated, and the green dots represent genes whose expression levels are down-regulated.

Figure 3

(A) Forest map of 12 immune genes related to the prognosis of hepatocellular carcinoma, analyzed by univariate Cox regression, all 12 genes are high risk genes (HR >1, P < 0.0001). (B) LASSO coefficient spectrum of 12 immune genes, Generate a coefficient distribution map for a logarithmic (λ) sequence. (C) Selecting the best parameters for LIHC in the LASSO model (λ).

Figure 4

(A, B) Distribution of risk score in patients with hepatocellular carcinoma. The black dotted line serves as the dividing line between the high-risk group and the low-risk group. (C, D) Diagram of the relationship between risk score and patient survival time. (E, F) Heat map of five immune genes in prognostic model, the abscissa represents the sample: the red area is the low-risk group, and the blue area is the high-risk group. The result of (A, C, E) is based on TCGA data(training set), and the result of (B, D, F) is based on ICGC data (validation set).

Figure 5

(A, B) Kaplan-Meier survival curves of the training set and the validation set, the survival prognosis of the patients in the high-risk group was significantly worse than that of the patients in the low-risk group (P < 0.05). (C, D) The ROC curve of the prognostic model, the results of the training set (AUC = 0.762) and the validation set (AUC = 0.758) show that the predictive ability of the model is good. (The meaning of AUC value: 0.5-0.7 indicates acceptable predictive ability, 0.7-0.9 indicates good predictive ability, >0.9 indicates excellent predictive ability).

Figure 6

(A) Univariate independent prognostic analysis of the training set. Forest plot of the association between risk factors and overall survival of patients. (B) Multivariate independent prognostic analysis of the training set. The risk score based on the prognostic model can be used as an independent prognostic factor for hepatocellular carcinoma. (C, D) Independent prognostic analysis of univariate and multifactorial factors in the validation set. The forest map shows that the risk score can also be used as an independent prognostic factor for hepatocellular carcinoma in the validation set.

Figure 7

Correlation analysis between the expression of immune microenvironment cells and risk score. There are six types of immune cells involved in the correlation analysis, including B cells (A), CD4⁺ T cells (B), CD8⁺ T cells (C), dendritic cells (D), macrophages (E), and neutrophils (F).

Figure 8

(A) Histogram of immune cells. The abscissa is the sample selected from the TCGA-LIHC cohort (screening condition: P < 0.05), and the vertical row represents the composition of various immune cells in the sample. Each color represents a different cell type. (B–F) Violin diagrams of five genes for constructing the model, including HDAC1 (B), BIRC5 (C), SPP1 (D), STC2 (E) and NR6A1 (F). There are 22 kinds of immune cells on the abscissa, and the ordinate represents the content of immune cells. The red area represents the high expression group of genes, and the green area represents the low expression group of genes.

Figure 9

The five prognostic genes are upregulated in human HCC specimens: HDAC1 (A) , BIRC5 (B) , SPP1 (C) , STC2 (D) , NR6A1 (E). Training Set: TCGA data showing the expression profiles of the five prognostic genes in normal liver (n = 50) vs tumor tissue (n = 371). Validation Set: ICGC data showing the expression profiles of the five prognostic genes in normal liver (n = 202) vs tumor tissue (n = 243).

Figure 10

(A, B) The expression profiles of these five genes in normal liver tissues and hepatocellular carcinoma tissues. The results of immunohistochemistry were all from HPA database.

Figure 11

(A) Enrichment plot of 261 genes enriched on the signal transduction pathway of p53 class mediator. (B) Regulatory networks between prognosis-related immune genes and genes enriched in the pathway. There are three genes (BIRC5, HDAC1, NR6A1) in the model that participate in the hub of the network.