Stratified analysis identifies HIF-2 α as a therapeutic target for highly immune-infiltrated melanomas

Abstract

While immune-checkpoint blockade (ICB) has revolutionized treatment of metastatic melanoma over the last decade, the identification of broadly applicable robust biomarkers has been challenging, driven in large part by the heterogeneity of ICB regimens and patient and tumor characteristics. To disentangle these features, we performed a standardized meta-analysis of eight cohorts of patients treated with anti-PD-1 (n=290), anti-CTLA-4 (n=175), and combination anti-PD-1/anti-CTLA-4 (n=51) with RNA sequencing of pre-treatment tumor and clinical annotations. Stratifying by immune-high vs -low tumors, we found that surprisingly, high immune infiltrate was a biomarker for response to combination ICB, but not anti-PD-1 alone. Additionally, hypoxia-related signatures were associated with non-response to anti-PD-1, but only amongst immune infiltrate-high melanomas. In a cohort of scRNA-seq of patients with metastatic melanoma, hypoxia also correlated with immunosuppression and changes in tumor-stromal communication in the tumor microenvironment (TME). Clinically actionable targets of hypoxia signaling were also uniquely expressed across different cell types. We focused on one such target, HIF-2α, which was specifically upregulated in endothelial cells and fibroblasts but not in immune cells or tumor cells. HIF-2α inhibition, in combination with anti-PD-1, enhanced tumor growth control in pre-clinical models, but only in a more immune-infiltrated melanoma model. Our work demonstrates how careful stratification by clinical and molecular characteristics can be leveraged to derive meaningful biological insights and lead to the rational discovery of novel clinical targets for combination therapy.

Keywords: HIF-2α; anti-PD-1; hypoxia; immunotherapy; melanoma.

Figure 1.. Study overview, cohort selection, and clinical characteristics.

a. Diagram depicting study overview. The study has four major components: (1) aggregate bulk RNA-sequencing data from 8 cohorts for a meta-analysis, (2) identify features of the tumor microenvironment (TME) associated with response and resistance, (3) analyze cell types associated with such features in orthogonal scRNA-sequencing cohorts to evaluate the potential for targeting such signatures, and (4) testing in pre-clinical models. b. Eight cohorts included in bulk transcriptomic analysis (n=694). c. Response rates by treatment group. Only samples which passed clinical inclusion criteria with response data were included; 33 samples were missing response annotations. d. Progression-free survival rates stratified by treatment group. Samples which passed clinical exclusion criteria with survival data were considered. 149 samples were missing survival data or were treated sequentially; these were excluded from the analysis. The survival rates between the four treatment groups were significantly different (Log-rank p < 0.0001). Significance of pairwise log-rank test also shown, where ∗ ∗ ∗p < 0.001, ∗ ∗ p < 0.01, ∗p < 0.05.

Figure 2.. Immune stratification reveals increased response rates to combination ICB in patients with immune-high tumors.

a. Heatmap showing the scaled ssGSEA scores of 19 immune cell type signatures across all bulk RNA-seq samples. Across all samples, hierarchical clustering is used to define two clusters: immune-high and immune-low. b. Co-correlation matrix between the immune score and other immune infiltrate metrics showing high levels of correlation between the immune score and other metrics (between 0.56 and 0.94). From left to right and top to bottom, the columns of the heatmap are: (1) expanded IFNγ signature [12], (*) the immune score (outlined in red), (2) Hallmark IFNγ signature [26], (3) TLS marker genes [27], (4) CIBERSORTx Tirosh immune proportion [70], (5) CIBERSORTx LM22 absolute immune score [71], (6) non-expanded IFNγ signature [12], (7) CIBERSORTx LM22 CD8+ T cell proportion [71], (8) CIBERSORTx Tirosh CD8+ T cell proportion [70], (9) Cabrita TLS signature [27]. c. In samples with both whole exome sequencing and bulk RNA-seq (n=162), tumor purity (estimated by FACETS [77]) anti-correlates with the immune score (Spearman’s ρ = −0.64, p < 0.0001). d. The immune score by responders and non-responders in each treatment group. e-f. Frequency of response by treatment group in (e) immune-low and (f) immune-high. g. Odds ratio and significance by one-sided Fisher’s Exact Test for response in immune-high versus immune-low populations for each treatment group. h. Odds ratio and significance by one-sided Fisher’s Exact Test for response in combination versus anti-PD-1 in immune-high and immune-low. For d-f., 20 samples with RNA-seq but no response annotations are not shown. ∗ ∗ ∗p < 0.001, ∗ ∗ p < 0.01, ∗p < 0.05; n.s. denotes “not significant.”

Figure 3.. Hypoxia and related pathways enrich in non-responders to anti-PD-1 in immune-high tumor samples.

a. Enrichment of 254 pathways from the Hallmarks and Kegg databases, including select T cell signatures from ImmuneSigDB. Enrichment analysis was performed for each group of immune-high or low samples for each treatment group separately. Number of significant treatment/immune conditions for which the signature is enriched (0–7) is shown on the annotations on the right of the heatmap. Broad categories for sets of gene signatures of interest are labelled to the right of the significance annotation; grey bars denote signatures which fall into the given categories. b. Normalized enrichment score for response to anti-PD-1, immune-high. Hypoxia and related signatures are labelled; *: FDR< 0.05, **: FDR< 0.01. c. ssGSEA score of hypoxia response to anti-PD-1 in immune-high and immune-low. *: p < 0.05 by Wilcoxon rank-sum. n.s. denotes “not significant.”

Figure 4.. Single-cell RNA-seq reveals hypoxia associates with immune exclusion and immunosuppressive phenotypes.

a. UMAP showing general cell type clustering for 175,497 cells included in analysis . b. Sample-level hypoxia stratification based on a cell-type normalized sample-level hypoxia score (Methods). Samples are designated as low or high based on tertiles. Mean and standard deviation of normalized hypoxia score shown. c. Pearson correlation between hypoxia score across compartments. d. Immune frequency in each sample by hypoxia state. Significance from Wilcoxon-rank sum, ∗p < 0.05. e. UMAP showing immune cell type clustering . f. Heatmap showing the pre-ranked GSEA normalized enrichment score (NES) for each of 254 signatures plus 4 additional hypoxia signatures. Genes were ranked by the estimated hypoxia score effect for each cell type as (see Methods on “modeling sample-level hypoxia”). Select gene sets of interest are labelled. Signatures and cell types are clustered by hierarchical clustering. Cell type subsets not shown (B cells, plasma cells, pDCs, mast cells, neutrophils) did not have enough representation across low, mid, and high hypoxia samples to include in the analysis. g. Estimated hypoxia score effect for all genes expressed in pseudobulk CD8+ T cells with log2(fold change) greater than 0.1 between high and low hypoxia samples (2,117 genes), as estimated by a linear model. Genes are ranked from lowest to highest based on the estimated effect of the sample-level hypoxia score. Select genes of interest are labelled. h. Normalized ssGSEA score of an LCMV effector versus exhausted CD8+ T cell signature in CD8+ T cell pseudobulks correlated with sample-level hypoxia score (Pearson’s R = −0.52, p < 0.05). Note that sample 194008 (hypoxia score of 2.46), is dropped from this panel because no CD8+ T cells were identified from this sample. i. Same as g. for tumor-associated macrophages (TAMs, 2,107 genes). j. Normalized ssGSEA score of the SPP1+ TAM signatures from Wei et al. 2021 [43] in TAMs correlated with sample-level hypoxia score (Pearson’s R = 0.31, p = 0.14).

Figure 5.. Effect of hypoxia on stromal cell characteristics and stromal-tumor communication.

a. AUCell Hallmark Hypoxia score for all cells in each general cell type. b. Log-normalized expression (counts per ten thousand) of labelled gene for each cell type averaged across all cells of that type per sample. a-b. Cell types with significantly higher scores by Wilcoxon rank-sum are shown, with each cell type compared against all others. *: p < 0.05, **: p < 0.01, ***: p < 0.001. c. Frequency of fibroblasts (out of all cells) in low-, mid-, and high-hypoxia samples. d. Estimated hypoxia score effect for all genes expressed in pseudobulk of labelled cell type with log2(fold change) greater than 0.25 between high and low hypoxia samples (8,410 genes in fibroblasts), as estimated by a linear model (see Methods for details). Genes are ranked from lowest to highest based on the estimated effect of the sample-level hypoxia score. Select genes of interest are labelled. e. Frequency of endothelial cells (out of all cells) in low-, mid-, and high-hypoxia samples (5,871 genes in endothelial cells with log2(fold change) greater than 0.25. f. Like d., but with endothelial cells. g-h. The relationship between the number of stromal-tumor interactions estimated by the CellPhoneDB statistical analysis method [45] and the sample-level hypoxia score for fibroblasts and endothelial cells (ECs). Correlations shown are Pearson correlations. i-j. All stromal-tumor interactions which are significant ≥ 75% of the samples (p < 0.05 in at least 6/8) in each category: low hypoxia (left, pink), mid hypoxia (middle, white), high hypoxia (right, grey). All ≥ 75% significant interactions in multiple categories (e.g. ‘KDR_VEGFA’ in EC-tumor interactions for both mid hypoxia and high hypoxia) are shown in dashed overlapping regions. Interacting pairs are listed as “partnerA_partnerB”, where partnerA is expressed on stromal cells and partnerB is expressed on tumor cells.

Figure 6.. Targeting hypoxia in combination with anti-PD-1 slows tumor growth in pre-clinical models of immune-high tumors

a. Graphic depicting experimental design of combination PT2399/anti-PD-1 treatment experiments. b-c. Growth curves showing tumor growth in response to HIF2α inhibitor (PT2399) and anti-PD-1 blockade combination treatment for (b) B16-OVA and (c) B16.F10. Data representative of two independent experiments (n > 5 per group). d-i. Flow cytometry characterization of T cells isolated from B16-OVA tumors treated with isotype control, anti-PD-1, PT2399, or the combination. Data showing two independent experiments combined (n=4–5 per group). d. Percentage of CD8+ T cells of CD45+ cells. e. Percentage of conventional T cells (Tcon; CD4+FOXP3-) of CD45+ cells. f. Percentage of regulatory T cells (Tregs; CD4+FOXP3+) of CD45+ cells. g. Percent of PD-1+ CD8+ T cells out of total CD8+ T cells. h. Percent of TIM-3+ CD8+ T cells out of total CD8+ T cells. i. Percent of Ki-67+ CD8+ T cells out of total CD8+ T cells. j. Percentage of HIF-1α+ CD31+ cells out of total CD31+ cells in the tumor. Significance levels are denoted as *p < 0.05, **p < 0.01, ***p < 0.001.