Development and Validation of a Combined Hypoxia and Immune Prognostic Classifier for Lung Adenocarcinoma

Abstract

Lung cancer is the leading cause of cancer-related deaths worldwide. Hypoxia is a crucial microenvironmental factor in lung adenocarcinoma (LUAD). However, the prognostic value based on hypoxia and immune in LUAD remains to be further clarified. The hypoxia-related genes (HRGs) and immune-related genes (IRGs) were downloaded from the public database. The RNA-seq expression and matched complete clinical data for LUAD were retrieved from The Cancer Genome Atlas (TCGA) and Gene Expression Omnibus (GEO) database. The least absolute shrinkage and selection operator (LASSO) Cox regression analysis was applied to model construction. Hypoxia expression profiles, immune cell infiltration, functional enrichment analysis, Tumor Immune Dysfunction and Exclusion (TIDE) score and the somatic mutation status were analyzed and compared based on the model. Moreover, immunofluorescence (IF) staining in human LUAD cases to explore the expression of hypoxia marker and immune checkpoint. A prognostic model of 9 genes was established, which can divide patients into two subgroups. There were obvious differences in hypoxia and immune characteristics in the two groups, the group with high-risk score value showed significantly high expression of hypoxia genes and programmed death ligand-1 (PD-L1), and maybe more sensitive to immunotherapy. Patients in the high-risk group had shorter overall survival (OS). This model has a good predictive value for the prognosis of LUAD. We constructed a new HRGs and IRGs model for prognostic prediction of LUAD. This model may benefit future immunotherapy for LUAD.

Keywords: hypoxia; immune; immunotherapy; lung adenocarcinoma; prognosis.

Figure 1

Flowchart of this study.

Figure 2

Construction of a predictive model of LUAD. (A) Venn diagram showing 1423 genes identified in the TCGA and GSE31210 cohorts. (B) Venn diagram showing 54 prognostic associated genes identified in the TCGA and GSE31210 cohorts. (C, D) Screening the optimal parameter (lambda), which is represented by the vertical black line in the plot. (E) Correlation between the risk score value and the nine genes.

Figure 3

Evaluation of the nine-genes signature in the TCGA cohort. (A) The distributions of the risk score. (B, C) t-SNE and PCA analysis showed significant differences between groups of patients. (D, E) The distributions of OS status and OS of patients between high-risk and low-risk groups, patients in the high-risk group had higher score values and mortality. (F) Time-independent ROC analysis of the risk score for prediction of the OS, the area under the curve for 1, 2, and 3 years reached 0.726, 0.736, and 0.710, respectively.

Figure 4

Validation of the nine-genes signature in the GSE31210 cohort. (A) The distributions of the risk score. (B, C) t-SNE and PCA analysis showed significant differences between groups of patients. (D, E) The distributions of OS status and OS of patients between high-risk and low-risk groups, patients in the high-risk group had higher score values and mortality. (F) Time-independent ROC analysis of the risk score for prediction of the OS, the area under the curve for 1, 2, and 3 years reached 0.904, 0.808, and 0.733, respectively.

Figure 5

Hypoxia profiles in the risk score stratified groups. (A, B) Box plot showed that the expression of key hypoxia genes was significantly higher in the high-risk group of the TCGA (A) and GSE31210 (B) cohorts. (C, D) GSEA demonstrated that hypoxia-related biological processes enriched in the high-risk group of the TCGA (C) and GSE31210 (D) cohorts. *P < 0.05, **P < 0.01, ***P < 0.001, and ns P > 0.05.

Figure 6

Immune profiles in the risk score stratified groups. (A-H) Correlation analysis of the risk scores value and microenvironment scores (A, E), stromal scores (B, F), hematopoietic stem cells (C, G), T cell CD4+ Th2 (D, H) in the cohorts. Spearman's rank correlation analysis was used for data analysis. (I, J) The comparison of immune cell fractions between groups of the TCGA (I) and GSE31210 (J) cohorts via the CIBERSORT method. The Wilcoxon signed-rank test was used for data analysis.

Figure 7

The expression level of PD-L1, TMB, IF, and TIDE score in the risk score stratified groups. (A, B, C) PD-L1 is more expressed in high-risk groups in the TCGA (A), GSE31210 (B) and GSE72094 (C) cohorts. (D) TMB was higher in the high-risk group in the TCGA cohort (P = 0.084). (E, F) Representative images of using IF to detect CA-IX (red), PD-L1 (green) and DAPI (blue) in the TJMUGH cohort. (G, H) Statistical analysis shows the differences in the mean fluorescence intensity (MFI) of PD-L1 (G) and CA-IX (H) between groups, the expression of PD-L1 and CA-IX was higher in the high-risk patients. (I, J) The differences of TIDE score (I) and immunotherapy sensitivity (J) between groups, high-risk group have lower TIDE scores and higher response rates in the TCGA cohort. *P < 0.05, ***P < 0.001.

Figure 8

Forrest plot of the univariate and multivariate Cox regression analyses regarding OS in the TCGA (A, B) and the GSE31210 (C, D) cohorts.