Back to Journals » Journal of Hepatocellular Carcinoma » Volume 7

Development and Verification of the Hypoxia-Related and Immune-Associated Prognosis Signature for Hepatocellular Carcinoma

Authors Hu B, Yang XB, Sang XT 

Received 10 August 2020

Accepted for publication 15 October 2020

Published 11 November 2020 Volume 2020:7 Pages 315—330

DOI https://doi.org/10.2147/JHC.S272109

Checked for plagiarism Yes

Review by Single anonymous peer review

Peer reviewer comments 4

Editor who approved publication: Dr Ahmed Kaseb



Bo Hu, Xiao-Bo Yang, Xin-Ting Sang

Department of Liver Surgery, Peking Union Medical College Hospital, Chinese Academy of Medical Sciences and Peking Union Medical College, Beijing 100730, People’s Republic of China

Correspondence: Xin-Ting Sang
Department of Liver Surgery, Peking Union Medical College Hospital, 1 Shuaifuyuan, Wangfujing, Beijing 100730, People’s Republic of China
Tel +86-10-69152836
Fax +86-10-69156043
Email [email protected]

Background: It has been widely suggested that the association of hypoxia with the immune status within the microenvironment of hepatocellular carcinoma (HCC) is of great clinical significance. The present work was carried out aiming to establish the hypoxia-related and immune-associated gene signature to stratify the risks in HCC.
Patients and Methods: The ssGSEA and t-SNE algorithms were utilized to estimate the immune and hypoxia statuses, respectively, using the TCGA database-derived cohort transcriptome profiles. Different immune groups are distinguished according to the ssGSEA scores, while the hypoxia-high and -low groups are inferred based on the distinct overall survival (OS) of the two groups of patients. Moreover, prognostic genes were identified using the Cox regression model in combination with the LASSO approach, which were later used to establish the hypoxia-related and immune-associated gene signature. At the same time, an ICGC cohort was used for external validation.
Results: A total of 13 genes, namely, HAVCR1, PSRC1, CCNJL, PDSS1, MEX3A, EID3, EPO, PLOD2, KPNA2, CDCA8, ADAMTS5, SLC1A7 and PIGZ, were discovered by the LASSO approach for constructing a gene signature to stratify the risk of HCC. Those low-risk cases showed superior prognosis (OS) to the high-risk counterparts (p< 0.05). Moreover, it was suggested by multivariate analysis that our constructed hypoxia-related and immune-associated prognosis signature might be used as the independent factor for prognosis prediction (p< 0.001). Patients in high-risk groups had severe hypoxia, higher immune checkpoint expression such as PD-L1, and different immunocyte infiltration states (eg, higher infiltration of regulatory T cells in the high-risk group) compared with those low-risk patients.
Conclusion: Our as-constructed hypoxia-related and immune-associated prognosis signature can be used as an approach to stratify the risk of HCC.

Keywords: hepatocellular carcinoma, microenvironment, hypoxia, immune, prognosis

Introduction

Liver cancer, of which hepatocellular carcinoma (HCC) is the dominant variety, represents the fourth leading cause of cancer deaths worldwide.1 Early-stage HCC is amenable to potentially curative treatment, which includes surgical treatment, radiofrequency ablation, and liver transplantation.2 However, given that the symptoms and physical characteristics of HCC cannot be easily detected, curative treatment is usually not possible at the time of diagnosis for >80% of patients.3 Exploring the HCC features is urgently needed, so as to exploit novel therapies.

Hypoxia is commonly seen in diverse solid tumors, which exerts a vital part in the genesis and progression of a tumor.4 Hypoxia can be detected in the context of HCC due to the destroyed balance between the supply and consumption of oxygen within the proliferating tumor.5 HCC and the normal liver show high vascularization; however, oxygen is greatly scavenged by the rapidly growing tumor cells in such nodules, finally generating the hypoxia microenvironment. As a matter of fact, HCC represents the malignancy experiencing the most severe hypoxia, with the median oxygen content of only 0.8%.6 The insufficient oxygen content in the tumor may induce numerous cellular and molecular responses, which may thereby affect the aggressiveness of the tumor and response to treatment. On the other hand, great attention has been paid to the development of immunotherapy to treat the solid tumors within the tumor immune microenvironment (TIME). Several immune-checkpoint inhibitors (ICIs), including anti-PD1 ligand (anti-PD-L1), anti-programmed death-1 (anti-PD1), as well as anti-cytotoxic T-lymphocyte-associated protein 4 (anti-CTLA-4) monoclonal antibodies,7 present favorable prognostic outcomes in selected advanced HCC cases. Besides, the interaction of tumor cells with the immune microenvironment components has been demonstrated as a key factor for HCC evolution as well as the possible immunotherapy response.8 It is interesting that more and more evidence discovers that hypoxia directly or indirectly interacts with the immune status in the HCC microenvironment,911 but the precise mechanisms have not been illustrated yet.

The present work conjectured that, the interaction between hypoxia and the immune status was of certain prognostic significance for HCC. The new prognosis signature was constructed by systemic analysis, in which the hypoxia and immune statuses were incorporated as the existing clinicopathological features and the staging system, with an aim of improving HCC prognosis.

Patients and Methods

Data Collection and Mining of mRNA Profiles

The level 3 messenger RNA (mRNA) expression patterns, together with the related clinical information, were obtained from 374 HCC and 50 healthy subjects from the TCGA database (https://cancergenome.nih.gov). In this study, cases with ≤30 days of survival or those with no survival data were eliminated, since they might die of the fetal complications (including hemorrhage, intracranial infection, and heart failure) rather than HCC. For further verification, the clinical data and transcriptional patterns were obtained from HCC cases in the International Cancer Genome Consortium (ICGC) database (https://icgc.org/). Among the mentioned cases, 203 HCC samples with complete mRNA expression data and corresponding clinical materials were selected for subsequent analyses.

Determination of Hypoxia Status Together with the Hypoxia-Associated Differentially Expressed Genes (DEGs)

The t-distributed Stochastic Neighbor Embedding (t-SNE) algorithm was used for identifying hypoxia status.12 t-SNE is a kind of unsupervised and non-parametric approach, which allows for dividing or condensing cases to multiple different clusters, according to the provided hallmarks or signatures. Besides, hallmark gene sets summarize and represent specific well-defined biological states or processes and display coherent expression. The set of hypoxia hallmark genes (n = 200), which is classic and has been used for hypoxia-related analysis of other tumors,13,14 was acquired based on the Molecular Signatures Database (MSigDB version 6.0), and the Euclidean distance was calculated between any two patients in the discovery cohort and condensed into two-dimensional points using a nonlinear dimensionality reduction algorithm t-SNE. In line with the clusters, two clusters (namely, “hypoxia-low” and “hypoxia-high” groups) were selected for assessing hypoxia status. Moreover, changes in target gene expression related to the HIF-1 signal transduction pathway were examined between both groups for exploring the relationship with hypoxia. Those target genes were obtained based on the Kyoto Encyclopedia of Genes and Genomes (KEGG) database (https://www.kegg.jp/; ID:04066), among which, 12 were related to “Increase oxygen delivery”, whereas 12 were associated with “Reduce oxygen consumption”. Meanwhile, DEGs between both groups were discovered using the “limma” algorithm15 of R software (Version 3.6.1; https://www.r-project.org/). Moreover, genes with the absolute log2 value (fold change, FC) >2 and the p-value < 0.001 after adjusting for false discovery rate (FDR) were deemed to be the hypoxia-associated DEGs.

Determination of the Immune Status Together with Immune-Associated DEGs

For the TCGA dataset, 29 immune signatures16 were first quantified for their enrichment degrees within respective HCC samples using single-sample gene-set enrichment analysis (ssGSEA) score.17,18 The 29 immune-associated gene sets represented diverse immune cell types, functions and pathways. The HCC samples were hierarchically clustered to high, medium or low immune group based on ssGSEA scores for those 29 immune signatures. In addition, the Estimation of Stromal and Immune cells in Malignant Tumor tissues using Expression data (ESTIMATE) algorithm,19 which has been developed recently, was adopted to measure stromal level (stromal score), immunocyte infiltration degree (immune score) and tumor purity in respective HCC samples as the validation of ssGSEA for immune status grouping. DEGs between low and high immune groups were studied by the “limma” package. Besides, genes with the absolute log2 value (FC) >1 and the p-value <0.05 after adjusting for FDR were deemed as the immune-associated DEGs. The gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) pathway analyses were conducted for these DEGs.

Construction and Verification of Our Prognosis Signature Associated with Hypoxia and Immune Statuses

Afterwards the intersection between immune- and hypoxia-associated DEGs was taken, and those overlapped genes were selected in later analyses. The prognosis immune-related model was constructed through univariate and multivariate Cox regression analysis, as well as the Least Absolute Shrinkage and Selection operator (LASSO), so as to predict overall survival (OS) for HCC cases. LASSO Cox regression was performed to avoid over-fitting while deleting those tightly correlated genes, and the vital genes were extracted from the genes screened out by univariate Cox regression. Then, multivariate Cox regression was also conducted to assess the contribution made by corresponding genes to prognosis. In addition, the risk scores were calculated by multiplying gene expression by the linearly combined regression coefficient acquired upon multivariate Cox regression. Based on the median risk score, all cases were classified as a high- or low-risk group.

Relationships of Prognostic Gene Signature with Hypoxia and Immunocyte Infiltration

Gene set enrichment analysis (GSEA) can be adopted for determining the statistical significance of the previously determined gene set as well as the consistent heterogeneities between two biological statuses.20 To assess the different hypoxia statuses in the high-risk group compared with the low-risk one, the present study carried out GESA by JAVA program https://www.broadinstitute.org/gsea). The FDR q < 0.25 and P < 0.05 indicated statistical significance. Besides, CIBERSORT algorithm was applied in analyzing the different immunocyte infiltration degrees in the low-risk group compared with the high-risk group with regard to 22 immunocyte subunits. CIBERSORT is a kind of deconvolution algorithm based on gene expression, which is developed to express cell composition in a tissue.21 Besides, the online database TIMER, which accounts for the web resource for the systemic evaluation of clinical significance of different immunocytes to different types of cancers, can be used to examine and visualize the levels of tumor-infiltrating immunocytes in addition to the CIBERSORT method.22 The TIMER database contains a total of 10,009 samples from 23 types of cancers derived from TCGA, which can assist in estimating the levels of six tumor-infiltrating immunocytes, namely, CD4 T cells, CD8 T cells, B cells, neutrophils, dendritic cells (DCs) and macrophages. In the present work, those immunocyte infiltration degrees were obtained from TCGA HCC cases, and the relationships between our constructed prognosis signature and immunocyte infiltration degrees were determined.

Statistical Methods

Statistical analysis was performed using R (v.3.6.1) software. Qualitative variables were explored by Fisher’s exact test or Pearson chi-square test. “Rtsne” of R package was employed to perform the t-SNE algorithm on the basis of nonlinear dimensionality reduction. In addition, “estimate” package was utilized to input immune scores. Meanwhile, “glmnet” package was used to model LASSO Cox regression. A difference of P < 0.05 indicated statistical significance unless specified otherwise.

Results

Hypoxia and the Hypoxia-Associated DEGs within HCC

The Euclidean distance between two cases was calculated in the discovery cohort using the expression matrix constituted by 200 hallmark genes of hypoxia based on MSigDB version 6.0, which was then concentrated to the 2D points by t-SNE, the nonlinear dimensionality reduction algorithm. Then, this study determined altogether three clusters and assigned every patient to the closest Euclidean distance (Figure 1). Altogether 169, 164 and 41 cases were collected into three clusters (Cluster 1–3), separately. Significant differences were detected across those three clusters upon survival comparison (p < 0.05 upon Log rank test). Previous studies have confirmed that hypoxia contributes to poor prognosis in HCC patients.23,24 Among the clusters, cases of Cluster3 had the best OS, whereas patients of Cluster1 displayed the poorest prognosis (Figure 1B), suggesting that Cluster1 and Cluster3 might be at the highest and lowest hypoxia states, separately. Meanwhile, the heatmap of the marker genes among the three clusters is shown in Figure 1C. As observed, the marker genes of hypoxia-high group (Cluster 1) were enriched into several hypoxia-related GO terms of biological process (BP), such as “Response to reduced oxygen level”, “Response to hypoxia”, “Glycolytic process” together with “Response to oxygen level” (Figure 1D; Supplementary Table 1). GO analysis concerning the marker genes of the hypoxia-low group (Cluster 3) was also conducted (Figure 1E; Supplementary Table 2), and “Glycosaminoglycan catabolic process”, “Glucose metabolic process” and “Aminoglycan catabolic process” were the top three GO terms of enrichment. Compared with Cluster 3, the marker genes in Cluster 1 were more enriched in several GO items directly related to hypoxia. Furthermore, changes in the target gene expression (hypoxia-high versus hypoxia-low) involved in the KEGG HIF-1 signal transduction pathway were explored. Genes participating in the “Increase oxygen delivery” (n = 12), together with those related to “Reduce oxygen consumption” (n = 12) were used. For those 12 genes involved in “Increase oxygen delivery”, 10 (83.33%) showed over-expression within the hypoxia-high group (Supplementary Figure 1), relative to that in the hypoxia-low group. Additionally, 8 out of 12 (66.67%) genes associated with “Reduce oxygen consumption” showed over-expression within the hypoxia-high group (Supplementary Figure 2). According to the above findings, those specific groups displayed significant association with hypoxia status. Therefore, cases of Cluster 1 and Cluster 3 were divided into “hypoxia-high” or “hypoxia-low” groups, separately. Moreover, “limma” package was used to identify a total of 1798 DEGs between hypoxia-high and hypoxia-low groups, among which, 1532 were up-regulated whereas 266 were down-regulated (Figure 1F and G).

Figure 1 Identification of hypoxia status and hypoxia-related differentially expressed genes (DEGs). (A) Dot plot for three distinct clusters identified by t-SNE algorithm based on 200 hypoxia hallmark genes. (B) Kaplan–Meier plot of overall survival for patients in three clusters. (C) Heatmap of the marker genes of the three clusters. The Gene Ontology (GO) of biological process (BP) analysis of the marker genes in Cluster 1 (D) and Cluster 3 (E) are shown. The thickness of the gray line represents the degree of association between GO terms. Heatmap (F) and volcano plot (G) demonstrating differentially expressed hypoxia-related genes between hepatocellular carcinoma (HCC) and non-tumor tissues. Red dots represent differentially up-regulated expressed genes, green dots represent differentially down-regulated expressed genes and black dots represent no differentially expressed genes.

Immune Status and Immune-associated DEGs in HCC

Afterwards, 29 immune-related genes representing different types of immunocytes were analyzed along with the corresponding pathways (Supplementary Table 3). Based on the ssGSEA scores that quantified the activities or abundances of immune signatures mentioned above in the cancer samples, HCC samples in the TCGA database were hierarchically clustered in immune-high, -medium or -low groups, which showed different activities and abundances of immunocytes, pathways or functions (Figure 2A and B). Besides, the ESTIMATE approach was used to calculate stromal scores and immune scores, and our results showed that the stromal scores and immune scores in the immune-high group remarkably increased compared with those in the immune-low group (Figure 2C and D). This study obtained the immune-associated DEGs through comparing gene expression in the immune-high group with that in the immune-low group, and altogether 1233 were identified (Figure 2E and F). The potential functions of these immune-related DEGs were then ascertained through GO annotation and KEGG enrichment pathway analysis (Figure 2G and H). Several immune-related GO terms of BP were identified, such as “Regulation of immune effector process”, “Regulation of B cell activation”, together with “Humoral immune response” (Supplementary Table 4). In addition, “Immunoglobulin complex”, “External side of plasma membrane” and “Immunoglobulin complex, circulating” were the most significantly enriched GO items of cellular component (CC) (Supplementary Figure 3A; Supplementary Table 5). Furthermore, the top 3 GO terms of molecular function (MF) were “Antigen binding”, “Immunoglobulin receptor binding” and “Extracellular matrix structural constituent” (Supplementary Figure 3B; Supplementary Table 6). As indicated from the KEGG pathway analysis, the DEGs were enriched via cancer and immune-associated pathways, such as the “T cell receptor signal transduction pathway”, “Th1 and Th2 cell differentiation”, “Th17 cell differentiation”, “PI3K−Akt signal transduction pathway” along with “TNF signal transduction pathway” (Supplementary Table 7).

Figure 2 Analysis of immune status of the Cancer Genome Atlas (TCGA) patients. (A) Hierarchical clustering of hepatocellular carcinoma (HCC) yields three subtypes (immune-high group, immune-medium group and immune-low group) in the TCGA database. (B) Different immune status of the three groups. Comparison of (C) immune score and (D) stromal score between immune-high and immune-low groups are shown. (E) Heatmap and (F) volcano plot demonstrating differentially expressed genes between immune-high and -low groups. Red dots represent differentially up-regulated expressed genes, green dots represent differentially down-regulated expressed genes and black dots represent no differentially expressed genes. (G) The Gene Ontology (GO) and (H) Kyoto encyclopedia of genes and genomes (KEGG) analysis of differentially expressed genes between immune-high and -low groups. ***P<0.001.

Development of the Risk Score Signature and Assessment of the Predicting Capacity

Those hypoxia-associated DEGs were intersected with the immune-associated ones, and altogether 320 overlapping genes were screened for subsequent analyses (Supplementary Table 8). Afterwards, this study chose 13 genes, hepatitis A virus cellular receptor 1 (HAVCR1), proline/serine-rich coiled-coil 1 (PSRC1), cyclin J-like (CCNJL), prenyl (decaprenyl) diphosphate synthase, subunit 1 (PDSS1), mex-3 RNA binding family member A (MEX3A), EP300 interacting inhibitor of differentiation 3 (EID3), erythropoietin (EPO), procollagen-lysine, 2-oxoglutarate 5-dioxygenase 2 (PLOD2), karyopherin alpha 2 (KPNA2), cell division cycle associated 8 (CDCA8), ADAM metallopeptidase with thrombospondin type 1 motif, 5 (ADAMTS5), solute carrier family 1 (glutamate transporter), member 7 (SLC1A7) and phosphatidylinositol glycan anchor biosynthesis, class Z (PIGZ) for constructing the prognosis signature via univariate and multivariate Cox regression analysis as well as LASSO (Figure 3A and B), aiming to categorize HCC patients into two groups with discrete OS, namely, high- or low-risk groups (Table 1). Thereafter, the risk score values were calculated according to the following formula: Risk score = [HAVCR1 expression* (0.0376)] + [PSRC1 expression * (0.0337)] + [CCNJL expression* (0.1417)] + [PDSS1 expression* (0.0530)] + [MEX3A expression* (0.0316)] + [EID3 expression* (0.2148)] + [EPO expression * (0.01476)] + [PLOD2 expression* (0.0081)] + [KPNA2 expression* (0.0296)] + [CDCA8 expression* (0.0381)] + [ADAMTS5 expression* (0.2877)] + [SLC1A7 expression* (0.0187)] + [PIGZ expression* (0.0065)]. Based on the median risk score, all cases were classified as high- or low-risk group. According to Kaplan-Meier (K-M) analysis (Figure 3C and D), high-risk patients had remarkably reduced OS relative to low-risk patients in different sets. Additionally, with regard to 1-year OS, the values of area under the curve (AUC) for TCGA HCC cohort and ICGC HCC cohort were 0.823 and 0.717, respectively (Figure 3E and F). In addition, the distributions of gene expression, risk score as well as survival status were plotted in the hypoxia- and immune-related signature of TCGA and ICGC HCC cohorts (Figure 3G and H). For better exploring the significance of our hypoxia- and immune-associated gene signature in independently predicting prognosis, univariate together with multivariate analysis was conducted, which revealed that the risk score value might serve as the independent factor for predicting the prognosis for TCGA HCC cohort (hazard ratio HR: 5.176, 95% confidence intervals CIs: 3.142−8.527, P<0.001) and ICGC HCC cohort (HR: 3.243, 95% CIs: 1.363−7.717, P<0.001) (Tables 1 and 2, respectively).

Table 1 Univariate and Multivariate Cox Regression Analyses of Clinicopathologic Characteristics Associated with Overall Survival in the Cancer Genome Atlas

Table 2 Univariate and Multivariate Cox Regression Analyses of Clinicopathologic Characteristics Associated with Overall Survival in the International Cancer Genome Consortium

Figure 3 LASSO coefficient profiles are depicted in (A). (B) shows the selection of the tuning parameter (lambda) in the LASSO model by tenfold cross-validation based on minimum criteria for OS; the lower X axis shows log (lambda), and the upper X axis shows the average number of OS-genes. The Y axis indicates partial likelihood deviance error. Red dots represent average partial likelihood deviances for every model with a given lambda, and vertical bars indicate the upper and lower values of the partial likelihood deviance errors. The vertical black dotted lines define the optimal values of lambda, which provides the best fit. Survival curves of patients in high risk group and low risk group of The Cancer Genome Atlas (TCGA) hepatocellular carcinoma (HCC) cohort (C) and the International Cancer Genome Consortium database (ICGC) HCC cohort (D) are exhibited. Patients in high-risk group suffered shorter overall survival. (E) and (F) show survival-dependent receiver operating characteristic (ROC) curves validation at 1 – year of prognostic value of the prognostic index in the two databases (TCGA and ICGC, respectively). Distribution of risk score, overall survival (OS), gene expression in (G) TCGA and (H) ICGC databases were also shown. Distribution of risk score, OS and heat map of the expression of eleven genes in low-risk and high-risk groups are listed in the picture from top to bottom.

Relationship Between the Prognosis Signature and the Clinicopathological Features

A total of 216 cases that had sufficient data on age, gender, clinical stage, T stage, tumor grade, platelet content, albumin content, alpha-fetoprotein (AFP) content, and vascular invasion were enrolled into the TCGA HCC cohort. The risk score values obtained based on the as-constructed signature showed significant correlation with tumor grade (Figure 4A) and gender (Figure 4B). In terms of the ICGC HCC cohort, a total of 203 patients with complete data on age, gender, clinical stage, portal vein invasion, hepatitis virus infection status, bile duct invasion, venous invasion, and hepatic fibrosis were collected. As a result, the risk score values were obviously related to clinical stage, portal vein invasion and venous invasion (Figure 4C–E).

Figure 4 Correlation of the prognostic signature with clinicopathological characteristics based on The Cancer Genome Atlas (TCGA) hepatocellular carcinoma (HCC) cohort and the International Cancer Genome Consortium database (ICGC) HCC cohort. The risk score was significantly correlated with tumor grade (A) and patient gender (B) in TCGA database, and clinical stage (C), vein invasion (D) and portal vein invasion (E).

Evaluation of Hypoxia Status within Both High- and Low-Risk TCGA Cohorts

To verify the different hypoxia statuses between the high- and low-risk TCGA cohorts, GSEA was implemented. Six hypoxia-related gene sets (including HARRIS_hypoxia M10508, MANALO_hypoxia_up M259, WINTER_hypoxia metagene M14072, REACTOME_cellular response to hypoxia M26925, WINTER_hypoxia_up M5466 and BIOCARTA_VEGF pathway M12975) were chosen for analysis. Remarkably, DEGs between the two cohorts showed significant enrichment into the above 6 gene sets (Figure 5A–F). The integrative diagrams that presents those above-mentioned items are exhibited in Figure 5G. Previous study has shown that hypoxia status is related to epithelial-mesenchymal transition (EMT) in HCC. Therefore, this study subsequently investigated the Twist family BHLH transcription factor 1 (TWIST-1), hypoxia-inducible factor-1α (HIF-1α), HIF-2α, snail family transcriptional repressor 1 (SNAI1), as well as EMT marker genes such as CDH1 (E-cadherin) and CDH2 (N-cadherin) for their expression quantities. Notably, the expression of all the genes mentioned above of the high-risk group remarkably increased compared with that of the low-risk group (Supplementary Figure 4).

Figure 5 Enrichment plots of hypoxia-related gene sets from gene set enrichment analysis (GSEA). GSEA results showing gene sets in (A) HARRIS_hypoxia (M10508), (B) MANALO_hypoxia_up (M259), (C) WINTER_hypoxia metagene (M14072), (D) REACTOME_cellular response to hypoxia (M26925), (E) WINTER_hypoxia_up (M5466) and (F) BIOCARTA_VEGF pathway (M12975) are differentially enriched in high risk phenotype. (G) summarizes the above six gene sets.

f

The Differences in Immunocyte Infiltration Degree Between High- and Low-Risk TCGA HCC Cohorts

In addition, the relationship of the hypoxia and immune-related prognosis model with the immunocyte infiltration degree of TCGA HCC cases was predicted for examining whether the risk score value partly reflected TIME status. Intriguingly, the findings in this study indicate that the macrophage (Cor = 0.501; p = 3.621e−23), neutrophil (Cor = 0.469; p = 3.518e−20) and DC (Cor = 0.389; p = 7.493e−14) contents were significantly positively correlated with the risk score (Figure 6A–C). Besides, CD8+ T cells (Cor = 0.276; p = 2.123e−07) (Figure 6D), CD4+ T cells (Cor = 0.250; p = 2.749e−06) (Figure 6E) and B cells (Cor = 0.238; p = 8.65e−06) (Figure 6F) also showed weak relationships with the risk score. Further, the CIBERSORT algorithm was also employed to evaluate different immunocyte infiltration degrees in the high-risk group compared with the low-risk group. As observed from Figure 7A, high-risk patients showed higher infiltration degrees of resting DCs, M0 macrophages and T regulatory cells (Tregs), while the low-risk group showed higher infiltration degrees of resting NK cells, CD8+ T cells and M1 macrophages.

Figure 6 Relationships between the hypoxia and immune-related prognostic model and infiltration abundances of six types of immune cells. The correlation was performed by using Pearson correlation analysis. (A) macrophages; (B) neutrophils; (C) dendritic cells; (D) CD8+T cells; (E) CD4+T cells and (F) B cells.

Figure 7 Analysis of different immune status in high- and low-risk groups of The Cancer Genome Atlas (TCGA) hepatocellular carcinoma (HCC) cohort. (A) Box plots visualizing significantly different immune cells between high-risk and low-risk patients. (BG) Box plots visualizing significantly different immune checkpoints between high-risk and low-risk cases. (HK) GSEA analysis revealing immune-related biological processes correlated with the signature. *P<0.05; **P<0.01; ***P<0.001.

Given than immunotherapy has become an established pillar of anti-cancer treatment, which improves the prognosis for cancer patients. Therefore, this study further detected the common immune checkpoint expression of both high- and low-risk groups. The results suggested that high-risk patients were associated with a significantly higher expression of T Cell Immunoreceptor with Ig and ITIM Domains (TIGIT), programmed death ligand 1 (PD-L1), program death-1 (PD-1), lymphocyte activation gene-3 (LAG3), T-cell immunoglobulin and mucin-domain containing-3 (TIM-3), and cytotoxic T-lymphocyte-associated protein 4 (CTLA-4) relative to that in the low-risk group (P<0.05) (Figure 7B–G). Moreover, we analyzed immune-associated biological processes correlated with the model. High risk HCC were significantly enriched in negative regulation of the immunity pathway, such as negative regulation of T cell differentiation, T cell signaling pathway, immune system process, and CD4+ positive αβ T cell activation (Figure 7H–K), which indicates an immunosuppressive microenvironment.

Discussion

Hypoxia has been reported to enhance the proliferation, metastasis, angiogenesis, radio- and chemo-resistance of HCC.25 Besides, HCC is the inflammation-related malignancy, in which the immune microenvironment has been recognized to promote cancer metastasis, invasion and development via constructing the symbiosis with tumor cells.26,27 Moreover, the hypoxia microenvironment further shapes the interaction between HCC and the stroma. It is found that immunomodulatory peptides, such as tumor necrosis factor-α (TNF-α) and interleukin-1 (IL-1), are able to trigger the HIF-1-dependent genes, suggesting the involvement of HIF-1 in immune responses.28

At present, few biomarkers have been developed based on mRNA expression profiles for estimating their status within HCC, due to the complexity of immune activity and hypoxia within the tumor microenvironment (TME).29 As a matter of fact, when hypoxia regions are developed in the tumor, the tumor is able to favorably respond to hypoxia and partially recover the supply of blood and nutrient within the tumor.30,31 Therefore, one single biomarker can hardly accurately examine the hypoxia status.32,33 Notably, t-SNE, the machine learning algorithm, has provided the refined and potent way of dimensionality reduction, thus assisting in exploring the possible subtypes within prostate cancer34 as well as breast cancer (BC).35 According to our results, t-SNE helped to identify different hypoxia TME patterns on the basis of 200 hallmark genes of hypoxia. In enrichment analysis, the marker genes of hypoxia-high group defined by us can be enriched in several GO terms directly related to hypoxia, such as “Response to reduced oxygen level” and “Response to hypoxia”. Besides, glycosaminoglycan catabolic process, which is significantly enriched in the hypoxia-low group, has been reported by Ströbel et al to occur in an environment of elevated oxygen content.36 In addition, changes in HIF-1 targeting gene expression levels were also examined for exploring the correlations with hypoxia. The above results confirm that the degree of hypoxia in the hypoxia-high group is more serious than that in the hypoxia-low group. Moreover, ssGSEA score was employed for quantifying the activities or abundances of the immune signatures mentioned above in the cancer samples. Besides, 3 groups were hierarchically clustered. ESTIMATE, the algorithm recently developed by utilizing specific cancer tissue transcriptional properties to deduce the cellularity and diverse infiltration degrees of normal cells,19 was used to examine the differences in immune status among different groups. Afterwards, the hypoxia- and immune-related risk score model was constructed based on the above-mentioned 13 genes. Notably, it was able to discriminate the high-risk population from the low-risk one, and the prognostic estimation was highly accurate. The low-risk HCC patients had prolonged OS compared with that in the high-risk HCC cohorts derived from TCGA and ICGC. With regard to clinical practicality, the prognostic model showed significant correlation with the gender and tumor grade of HCC patients from the TCGA dataset. Such a result indicated that the markedly increased risk score was determined in female patients and patients at advanced grade by our constructed model. For the ICGC HCC cohort, the risk score showed evident correlations with clinical stage, venous invasion and portal vein invasion, which highlighted that our model was related to vascular invasion. Moreover, high-risk patients selected by our model had a higher hypoxia status compared with their low-risk counterparts. Consequently, our constructed hypoxia- and immune-related signature might participate in HCC genesis and progression, rendering its potential as the valuable clinical biomarker.

Furthermore, findings in this study revealed that our signature was positively related to the infiltration degrees of 6 immunocytes, in particular for DCs, macrophages and neutrophils. Besides, CIBERSORT was also applied in calculating the percentages of 22 immunocyte subpopulations within HCC. The results suggested that, high proportions of resting DCs, Tregs and M0 macrophages were detected within high-risk patients compared with low-risk counterparts (P < 0.05). In addition, the TGF-β-miR-34a-CCL22 signaling-induced Tregs have been elaborated to promote the venous metastases of HBV-positive HCC, which is related to HCC progression and metastasis.37 Shen et al also demonstrated that, in TME of HCC, the enhanced prevalence as well as extensive functions of Tregs were related to cancer stage, which contributed to tumor prosperity and growth.38 It is worth mentioning that, the recruitment of Tregs is promoted by the hypoxia-induced Chemokine (C-C motif) ligand 28 (CCL28) in the HCC setting, which results in the enhanced angiogenesis and VEGF expression.39,40 Our model was also closely related to macrophage infiltration. Increasing evidence has revealed that DCs play a part as the negative predictor for HCC prognosis. For example, Zhou et al reported that plasmacytoid DCs infiltrated into the tumor, which might be used as a new biomarker for predicting adverse HCC prognosis, and this was attained through triggering immune tolerance as well as the inflammatory TME consisting of regulatory T cells along with IL-17-generating cells.41 Nevertheless, there are few studies on the effects of hypoxia on DCs for the time being. It is reported that, macrophages are recruited into hypoxia regions within cancer tissues, which thus activates the expression of pro-invasive, mitogenic, prometastatic and proangiogenic genes, such as TIE2 up-regulation.42,43 Such TIE2-positive macrophages (CD14+CD16+) promote angiogenesis, and their blood or tumor levels are tightly related to the micro-vessel density within HCC.44

In addition, hypoxia is suggested to promote HCC progression in murine through regulating IL-10 expression induced by HIF1α, and this contributes to the conversion of M1 macrophages within the tumor to M2 macrophages. The infiltration of M1 macrophages in high-risk patients calculated by our model significantly decreased compared with low-risk counterparts, and difference in M2 macrophage number between the two groups was not significant. M1 macrophages, which can activate the tumor-killing mechanisms and amplify the Th1 immunocytes responses, resist tumorigenesis.45 Tsuchiyama et al also found that the best monocyte chemoattractant protein-1 amount enhanced the anticancer activity of suicide gene therapy for the treatment of HCC through activating M1 macrophages.46 The potential antitumor effect of M1 macrophages might partially account for the superior prognosis for low-risk patients. Moreover, for low-risk patients, the CD8+ T cell infiltration degree increased relative to that of high-risk patients. As suggested by previous study, the strong tumor-associated antigen-specific CD8+ T cell responses suppressed HCC recurrence,47 demonstrating that CD8+ T cells exerted an antitumor role. However, HIF-1α together with the short isoform I.1 induced by its activation showed negative regulation on CD8+ T cells,48 which indicated that hypoxia might weaken CD8+ T cell function. Besides, hypoxia potently augmented the T-cell suppression in vitro mediated by macrophages depending on HIF-1α expression in macrophages.49 The infiltration degree of resting NK cells for low-risk patients remarkably increased relative to that of high-risk patients. Numerous immune recognition receptors are expressed in NK cells for recognizing the ligands on liver cells, Kupffer cells, stellate cells, and liver sinusoidal endothelial cells, thus maintaining the balance between immune tolerance and immune response in NK cells.50 Hypoxia has a profound effect on the NK cell infiltrate property, together with its influence on responses mediated by immune statuses in cancer tissues.51,52 Berchem et al illustrated that microvesicles derived from tumors under hypoxia had negative regulation on the functions of NK cells through the mechanism related to transferring miR23a and TGF-β,53 which suggested that the impaired cytotoxicity of NK cells might be induced by oxygen deficit. In brief, the prognostic gene signature developed in this study was able to accurately indicate the immunocyte infiltration degrees in high- and low-risk patients.

Furthermore, expression of immune checkpoints was detected in both high- and low-risk HCC cases, respectively. As a result, the expression of PD-L1, PD-1, TIM-3, LAG3 and CTLA-4 significantly increased in the high-risk HCC group compared with that in the low-risk group (P < 0.05), which meant that high-risk patients might benefit from the immune checkpoint inhibitors (ICIs) therapy. In view of the more severe hypoxia among high-risk patients, therapies targeting hypoxia may also be effective on high-risk patients; for instance, the bioreductive prodrugs that are activated by enzymatic reduction within the hypoxic tissues and the small molecular inhibitors related to hypoxic cell survival.54

Certain limitations should be noted in the present work. Firstly, the present work conducted independent external validation, yet all the changes across HCC cases derived from diverse geographical regions might not be covered, since all information and tissues were obtained from public databases in a retrospective way. Secondly, there are different microenvironment features among diverse tumor sites, like the invasion margin or the tumor core. Therefore, this study extracted all samples from the tumor core, making it impossible to assess the hypoxia and immune statuses among diverse tumor regions. As a result, our results should be further validated in the well-designed, multicenter, prospective studies.

Conclusions

To sum up, a prognostic signature is constructed and validated based on 13 hypoxia- and immune-associated genes, for the sake of predicting OS for HCC. This prognosis signature can help to select the individualized therapeutic strategy in clinical practice. In addition, our risk score model is connected with the hypoxia and immune statuses, which provides a comprehensive perspective for clarifying the underlying mechanisms that determine the prognosis for HCC.

Contribution to the field

Hypoxia and immunity are the hot spots of tumor related research, but for hepatocellular carcinoma, there is no research combining them. In this study, we combined hypoxia and immunity to construct a gene prognosis model. The model can well distinguish the high and low risk TCGA and ICGC HCC patients, so that the two groups of patients have different overall survival. In addition, the two groups of patients with hypoxia and immune cell infiltration state is also different. Our research is novel and valuable.

Data Sharing Statement

All data generated or analyzed during this study are included in this article and its supplementary information files.

Author Contributions

All authors made a significant contribution to the work reported, whether that is in the conception, study design, execution, acquisition of data, analysis and interpretation, or in all these areas; took part in drafting, revising or critically reviewing the article; gave final approval of the version to be published; have agreed on the journal to which the article has been submitted; and agree to be accountable for all aspects of the work.

Funding

There is no funding to report.

Disclosure

The authors declare that they have no competing interests.

References

1. El-Serag HB. Epidemiology of hepatocellular carcinoma. Liver. 2020:758–772.

2. Yang JD, Hainaut P, Gores GJ, Amadou A, Plymoth A, Roberts LR. A global view of hepatocellular carcinoma: trends, risk, prevention and management. Nat Rev Gastroenterol Hepatol. 2019:1–16.

3. Zongyi Y, Xiaowu L. Immunotherapy for hepatocellular carcinoma. Cancer Lett. 2020;470:8–17. doi:10.1016/j.canlet.2019.12.002

4. Harris AL. Hypoxia—a key regulatory factor in tumour growth. Nat Rev Cancer. 2002;2(1):38–47. doi:10.1038/nrc704

5. Menrad H, Werno C, Schmid T, et al. Roles of hypoxia-inducible factor-1α (HIF-1α) versus HIF-2α in the survival of hepatocellular tumor spheroids. Hepatology. 2010;51(6):2183–2192. doi:10.1002/hep.23597

6. Schlageter M, Terracciano LM, D’Angelo S, Sorrentino P. Histopathology of hepatocellular carcinoma. World J Gastroenterol. 2014;20(43):15955. doi:10.3748/wjg.v20.i43.15955

7. Hu B, Yang X-B, Sang X-T. Liver graft rejection following immune checkpoint inhibitors treatment: a review. Med Oncol. 2019;36(11):94. doi:10.1007/s12032-019-1316-7

8. Cariani E, Missale G. Immune landscape of hepatocellular carcinoma microenvironment: implications for prognosis and therapeutic applications. Liver Int. 2019;39(9):1608–1621. doi:10.1111/liv.14192

9. Chiu DK-C, Tse AP-W, Xu IM-J, et al. Hypoxia inducible factor HIF-1 promotes myeloid-derived suppressor cells accumulation through ENTPD2/CD39L1 in hepatocellular carcinoma. Nat Commun. 2017;8(1):1–12. doi:10.1038/s41467-017-00530-7

10. Ye L-Y, Chen W, Bai X-L, et al. Hypoxia-induced epithelial-to-mesenchymal transition in hepatocellular carcinoma induces an immunosuppressive tumor microenvironment to promote metastasis. Cancer Res. 2016;76(4):818–830. doi:10.1158/0008-5472.CAN-15-0977

11. Semaan A, Dietrich D, Bergheim D, et al. CXCL12 expression and PD-L1 expression serve as prognostic biomarkers in HCC and are induced by hypoxia. Virchows Archiv. 2017;470(2):185–196.

12. Cieslak MC, Castelfranco AM, Roncalli V, Lenz PH, Hartline DK. t-Distributed Stochastic Neighbor Embedding (t-SNE): a tool for eco-physiological transcriptomic analysis. Mar Genomics. 2019:100723.

13. Pant D, Narayanan SP, Vijay N, Shukla S. Hypoxia-induced changes in intragenic DNA methylation correlate with alternative splicing in breast cancer. J Biosci. 2020;45(1):1–14. doi:10.1007/s12038-019-9977-0

14. Lin W, Wu S, Chen X, et al. Characterization of hypoxia signature to evaluate the tumor immune microenvironment and predict prognosis in glioma groups. Front Oncol. 2020;10:796.

15. Phipson B, Lee S, Majewski IJ, Alexander WS, Smyth GK. Robust hyperparameter estimation protects against hypervariable genes and improves power to detect differential expression. Ann Appl Stat. 2016;10(2):946. doi:10.1214/16-AOAS920

16. He Y, Jiang Z, Chen C, Wang X. Classification of triple-negative breast cancers based on Immunogenomic profiling. J Exp Clin Cancer Res. 2018;37(1):327. doi:10.1186/s13046-018-1002-1

17. Barbie DA, Tamayo P, Boehm JS. Systematic RNA interference reveals that oncogenic KRAS-driven cancers require TBK1. Nature. 2009;462(7269):108–112. doi:10.1038/nature08460

18. Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinform. 2013;14(1):7. doi:10.1186/1471-2105-14-7

19. Yoshihara K, Shahmoradgoli M, Martínez E, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. 2013;4(1):1–11. doi:10.1038/ncomms3612

20. Subramanian A, Tamayo P, Mootha VK, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci. 2005;102(43):15545–15550. doi:10.1073/pnas.0506580102

21. Newman AM, Liu CL, Green MR, et al. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods. 2015;12(5):453–457.

22. Li T, Fan J, Wang B, et al. TIMER: a web server for comprehensive analysis of tumor-infiltrating immune cells. Cancer Res. 2017;77(21):e108–e110. doi:10.1158/0008-5472.CAN-17-0307

23. Liu L, Zhu X-D, Wang W-Q, et al. Activation of β-catenin by hypoxia in hepatocellular carcinoma contributes to enhanced metastatic potential and poor prognosis. Clin Cancer Res. 2010;16(10):2740–2750. doi:10.1158/1078-0432.CCR-09-2610

24. Zhang L, Huang G, Li X, et al. Hypoxia induces epithelial-mesenchymal transition via activation of SNAI1 by hypoxia-inducible factor-1α in hepatocellular carcinoma. BMC Cancer. 2013;13(1):1–9. doi:10.1186/1471-2407-13-108

25. Wu XZ, Xie GR, Chen D. Hypoxia and hepatocellular carcinoma: the therapeutic target for hepatocellular carcinoma. J Gastroenterol Hepatol. 2007;22(8):1178–1182. doi:10.1111/j.1440-1746.2007.04997.x

26. Fu Y, Liu S, Zeng S, Shen H. From bench to bed: the tumor immune microenvironment and current immunotherapeutic strategies for hepatocellular carcinoma. J Exp Clin Cancer Res. 2019;38(1):396. doi:10.1186/s13046-019-1396-4

27. Hu B, Yang X-B, Sang X-T. Development of an immune-related prognostic index associated with hepatocellular carcinoma. Aging. 2020;12(6):5010. doi:10.18632/aging.102926

28. Hellwig-Bürgel T, Stiehl DP, Wagner AE, Metzen E, Jelkmann W. hypoxia-inducible factor-1 (HIF-1): a novel transcription factor in immune reactions. J Interferon Cytokine Res. 2005;25(6):297–310. doi:10.1089/jir.2005.25.297

29. Staudacher JJ, Naarmann-de Vries IS, Ujvari SJ, et al. Hypoxia-induced gene expression results from selective mRNA partitioning to the endoplasmic reticulum. Nucleic Acids Res. 2015;43(6):3219–3236. doi:10.1093/nar/gkv167

30. Pugh CW. Modulation of the hypoxic response. Adv Exp Med Biol. 2016;903:259–271.

31. Fukumura D, Jain RK. Tumor microvasculature and microenvironment: targets for anti-angiogenesis and normalization. Microvasc Res. 2007;74(2–3):72–84. doi:10.1016/j.mvr.2007.05.003

32. van Malenstein H, Gevaert O, Libbrecht L, et al. A seven-gene set associated with chronic hypoxia of prognostic importance in hepatocellular carcinoma. Clin Cancer Res. 2010;16(16):4278–4288. doi:10.1158/1078-0432.CCR-09-3274

33. Pinato DJ, Pai M, Reccia I, et al. Preliminary qualification of a novel, hypoxic-based radiologic signature for trans-arterial chemoembolization in hepatocellular carcinoma. BMC Cancer. 2018;18(1):211. doi:10.1186/s12885-018-4120-4

34. Ahmed M, Lai TH, Zada S, et al. Functional linkage of RKIP to the epithelial to mesenchymal transition and autophagy during the development of prostate cancer. Cancers. 2018;10(8):273. doi:10.3390/cancers10080273

35. Guo L, Chen G, Zhang W, et al. A high-risk luminal A dominant breast cancer subtype with increased mobility. Breast Cancer Res Treat. 2019;175(2):459–472. doi:10.1007/s10549-019-05135-w

36. Ströbel S, Loparic M, Wendt D, et al. Anabolic and catabolic responses of human articular chondrocytes to varying oxygen percentages. Arthritis Res Ther. 2010;12(2):R34. doi:10.1186/ar2942

37. Yang P, Li Q-J, Feng Y, et al. TGF-β-miR-34a-CCL22 signaling-induced Treg cell recruitment promotes venous metastases of HBV-positive hepatocellular carcinoma. Cancer Cell. 2012;22(3):291–303. doi:10.1016/j.ccr.2012.07.023

38. Shen X, Li N, Li H, Zhang T, Wang F, Li Q. Increased prevalence of regulatory T cells in the tumor microenvironment and its correlation with TNM stage of hepatocellular carcinoma. J Cancer Res Clin Oncol. 2010;136(11):1745–1754. doi:10.1007/s00432-010-0833-8

39. Ren L, Yu Y, Wang L, Zhu Z, Lu R, Yao Z. Hypoxia-induced CCL28 promotes recruitment of regulatory T cells and tumor growth in liver cancer. Oncotarget. 2016;7(46):75763. doi:10.18632/oncotarget.12409

40. Facciabene A, Peng X, Hagemann IS, et al. Tumour hypoxia promotes tolerance and angiogenesis via CCL28 and T reg cells. Nature. 2011;475(7355):226–230.

41. Zhou Z-J, Xin H-Y, Li J, Hu Z-Q, Luo C-B, Zhou S-L. Intratumoral plasmacytoid dendritic cells as a poor prognostic factor for hepatocellular carcinoma following curative resection. Cancer Immunol Immunother. 2019:1–11.

42. Murdoch C, Tazzyman S, Webster S, Lewis CE. Expression of Tie-2 by human monocytes and their responses to angiopoietin-2. J Immunol. 2007;178(11):7405–7411. doi:10.4049/jimmunol.178.11.7405

43. Lewis C, Murdoch C. Macrophage responses to hypoxia: implications for tumor progression and anti-cancer therapies. Am J Pathol. 2005;167(3):627–635. doi:10.1016/S0002-9440(10)62038-X

44. He Y-F, Wang C-Q, Yu Y, et al. Tie2-expressing monocytes are associated with identification and prognoses of hepatitis B virus related hepatocellular carcinoma after resection. PLoS One. 2015;10(11). doi:10.1371/journal.pone.0143657.

45. Murray PJ, Wynn TA. Protective and pathogenic functions of macrophage subsets. Nat Rev Immunol. 2011;11(11):723–737. doi:10.1038/nri3073

46. Tsuchiyama T, Nakamoto Y, Sakai Y, Mukaida N, Kaneko S. Optimal amount of monocyte chemoattractant protein-1 enhances antitumor effects of suicide gene therapy against hepatocellular carcinoma by M1 macrophage activation. Cancer Sci. 2008;99(10):2075–2082. doi:10.1111/j.1349-7006.2008.00951.x

47. Hiroishi K, Eguchi J, Baba T, et al. Strong CD8+ T-cell responses against tumor-associated antigens prolong the recurrence-free interval after tumor treatment in patients with hepatocellular carcinoma. J Gastroenterol. 2010;45(4):451–458.

48. Lukashev D, Klebanov B, Kojima H, et al. Cutting edge: hypoxia-inducible factor 1α and its activation-inducible short isoform I. 1 negatively regulate functions of CD4+ and CD8+ T lymphocytes. J Immunol. 2006;177(8):4962–4965. doi:10.4049/jimmunol.177.8.4962

49. Doedens AL, Stockmann C, Rubinstein MP, et al. Macrophage expression of hypoxia-inducible factor-1α suppresses T-cell function and promotes tumor progression. Cancer Res. 2010;70(19):7465–7475. doi:10.1158/0008-5472.CAN-10-1439

50. Sun C, Sun H, Zhang C, Tian Z. NK cell receptor imbalance and NK cell dysfunction in HBV infection and hepatocellular carcinoma. Cell Mol Immunol. 2015;12(3):292–302. doi:10.1038/cmi.2014.91

51. Parodi M, Raggi F, Cangelosi D, et al. Hypoxia modifies the transcriptome of human NK cells, modulates their immunoregulatory profile, and influences NK cell subset migration. Front Immunol. 2018;9:2358. doi:10.3389/fimmu.2018.02358

52. Cui C, Fu K, Yang L, et al. Hypoxia-inducible gene 2 promotes the immune escape of hepatocellular carcinoma from nature killer cells through the interleukin-10-STAT3 signaling pathway. J Exp Clin Cancer Res. 2019;38(1):1–16. doi:10.1186/s13046-019-1233-9

53. Berchem G, Noman MZ, Bosseler M, et al. Hypoxic tumor-derived microvesicles negatively regulate NK cell function by a mechanism involving TGF-β and miR23a transfer. Oncoimmunology. 2016;5(4):e1062968. doi:10.1080/2162402X.2015.1062968

54. Wilson WR, Hay MP. Targeting hypoxia in cancer therapy. Nat Rev Cancer. 2011;11(6):393–410. doi:10.1038/nrc3064

Creative Commons License © 2020 The Author(s). This work is published and licensed by Dove Medical Press Limited. The full terms of this license are available at https://www.dovepress.com/terms.php and incorporate the Creative Commons Attribution - Non Commercial (unported, v3.0) License. By accessing the work you hereby accept the Terms. Non-commercial uses of the work are permitted without any further permission from Dove Medical Press Limited, provided the work is properly attributed. For permission for commercial use of this work, please see paragraphs 4.2 and 5 of our Terms.