Cancer cell genetics shaping of the tumor microenvironment reveals myeloid cell-centric exploitable vulnerabilities in hepatocellular carcinoma

Abstract

Myeloid cells are abundant and plastic immune cell subsets in the liver, to which pro-tumorigenic, inflammatory and immunosuppressive roles have been assigned in the course of tumorigenesis. Yet several aspects underlying their dynamic alterations in hepatocellular carcinoma (HCC) progression remain elusive, including the impact of distinct genetic mutations in shaping a cancer-permissive tumor microenvironment (TME). Here, in newly generated, clinically-relevant somatic female HCC mouse models, we identify cancer genetics' specific and stage-dependent alterations of the liver TME associated with distinct histopathological and malignant HCC features. Mitogen-activated protein kinase (MAPK)-activated, Nras^G12D-driven tumors exhibit a mixed phenotype of prominent inflammation and immunosuppression in a T cell-excluded TME. Mechanistically, we report a Nras^G12D cancer cell-driven, MEK-ERK1/2-SP1-dependent GM-CSF secretion enabling the accumulation of immunosuppressive and proinflammatory monocyte-derived Ly6C^low cells. GM-CSF blockade curbs the accumulation of these cells, reduces inflammation, induces cancer cell death and prolongs animal survival. Furthermore, GM-CSF neutralization synergizes with a vascular endothelial growth factor (VEGF) inhibitor to restrain HCC outgrowth. These findings underscore the profound alterations of the myeloid TME consequential to MAPK pathway activation intensity and the potential of GM-CSF inhibition as a myeloid-centric therapy tailored to subsets of HCC patients.

Fig. 1. Genetically-distinct murine HCCs display distinct transcriptomic and histopathological features.

a Schematic diagram of experimental design: Hydrodynamic tail vein injection (HDTVi) method used to deliver the Sleeping beauty (SB) transposon and CrispR-Cas9 constructs that enforce expression of the depicted oncogenic drivers in hepatocytes in order to generate genetically-distinct HCCs. b Kaplan–Meier survival curves of HCC-bearing mice HDTV injected with the indicated combinations of oncogenic drivers (Myc^OE/Trp53^KO n = 19, median survival of 23 days; Myc^OE/Pten^KO n = 14, median survival of 29.5 days; Nras^G12D/Pten^KO n = 15, median survival of 62 days; Nras^G12V/Pten^KO n = 14, median survival of 118 days). Nras^G12V/Pten^KO vs Nras^G12D/Pten^KO p < 0.001; Nras^G12V/Pten^KO vs Myc^OE/Trp53^KO p < 0.001; Nras^G12V/Pten^KO vs Myc^OE/Pten^KO p < 0.001; Nras^G12D/Pten^KO vs Myc^OE/Trp53^KO p < 0.001; Nras^G12D/Pten^KO vs Myc^OE/Pten^KO p < 0.001; Myc^OE/Trp53^KO vs Myc^OE/Pten^KO n.s. (non-significant). ****p < 0.0001 c Longitudinal HCC volumes of individual tumor-bearing mice (Myc^OE/Trp53^KO n = 9, Myc^OE/Pten^KO n = 6, Nras^G12D/Pten^KO n = 15, Nras^G12V/Pten^KO n = 9) determined from weekly/bi-weekly magnetic resonance imaging (MRI), and represented relative to each animal initial tumor volume. d Representative images of H&E, Masson Trichrome, and Oil Red staining performed on sectioned livers collected from end-stage HCC-bearing mice. (Scale bars = 100 µm; representative of n = 4 mice). e Principal Component Analysis (PCA) plot depicting the transcriptome differences between genetically-distinct HCC bulk tumors (Myc^OE/Trp53^KO n = 3, Myc^OE/Pten^KO n = 3, Nras^G12D/Pten^KO n = 4, Nras^G12V/Pten^KO n = 3) and control livers injected with empty vectors, hereafter referred to as control (n = 3), following RNA-seq analyses (Supplementary Data 1). f Heatmap of unsupervised hierarchical clustering depicting the geneset enrichment of the MAPK, PI3K, and MYC signaling pathways in genetically-distinct HCCs (Myc^OE/Trp53^KO n = 3, Myc^OE/Pten^KO n = 3, Nras^G12D/Pten^KO n = 4, Nras^G12V/Pten^KO n = 3) compared to control (n = 3). The color scale represents the significance of the enrichment in –log₁₀(FDR). FDR False Discovery Rate. Statistical significance was determined by log-rank (mantel-cox) test (b), two-sided hypergeometric test with multiple testing correction using Benjamini-Hochberg (f). Source data are provided as a Source Data file.

Fig. 2. Transcriptional gene signatures derived from genetically-distinct HCC models recapitulate human HCC heterogeneity and predict patient prognosis.

a Heatmap of unsupervised hierarchical clustering depicting the global gene expression level as TPMs of genetically-distinct HCCs (Myc^OE/Trp53^KO n = 3, Myc^OE/Pten^KO n = 3, Nras^G12D/Pten^KO n = 4, Nras^G12V/Pten^KO n = 3). Top annotations represent the classification of HCC mouse models’ transcriptomic profiles based on the human molecular HCC subtypes^–. *FDR < 0.01. b Heatmap of unsupervised hierarchical clustering depicting the link between the MAPK, PI3K, and MYC signaling pathway enrichment (Supplementary Fig. 2c) and patients correlating to genetically-distinct HCC-derived transcriptional signatures across TCGA: Liver Hepatocellular Cancer (LIHC) patients (n = 423) (Supplementary Data 2). c Kaplan–Meier survival curves displaying TCGA:LIHC patients segregated according to their high/low correlation with the transcriptional signatures of each genetically-distinct HCC relative to control. Lines at survival probability = 0.5 depict median survival. Risk tables show the number of patients at the indicated time points (in months). d Representative IHC image for p-ERK1/2 and MYC performed on human HCC TMA sections from the Wu et al. dataset. Pie charts depict the percentage of patients positive for p-ERK1/2 and MYC. (Scale bars = 100 µm; representative of n = 99 p-ERK1/2 positive patients and n = 78 MYC positive patients). e Kaplan–Meier curves displaying the overall survival (OS) and recurrence-free survival (RFS) of HCC patients (from the Wu et al. dataset) segregated according to p-ERK1/2 negative (n = 368) or positive (n = 99) and c-MYC negative (n = 389) or positive (n = 78) staining in cancer cells. Risk tables show the number of patients at the indicated time points (in months) (see Supplementary Data 4 for median OS and RFS time). Statistical significance was determined by one-sided Fisher’s test using Bonferroni multiple testing correction (a), log-rank test (c, e). The shading represent 95% confidence interval (c, e). TPMs Transcripts per Million. Source data are provided as a Source Data file.

Fig. 3. Distinct oncogenic drivers’ activation impact HCC immune landscape.

a Myeloid and lymphoid cell contents relative to total CD45⁺ leukocytes in end-stage HCC and age-matched control livers Control-4-weeks n = 5, Myc^OE/Trp53^KO n = 6, Myc^OE/Pten^KO n = 5, Nras^G12D/Pten^KO n = 4, Control-15-weeks n = 5, and Nras^G12V/Pten^KO n = 4). b Quantification of circulating myeloid cells relative to total CD45⁺ leukocytes at the indicated timepoints for each of the HCC models (Control n = 10, Myc^OE/Trp53^KO n = 20, Myc^OE/Pten^KO n = 18, Nras^G12D/Pten^KO n = 30, Nras^G12V/Pten^KO n = 16). Arrows indicate timepoints when tumors were detectable. Statistical significance: final analyzed blood samples of tumor-bearing mice versus aged-matched controls. c Percentage of the immune cell populations relative to total CD45⁺ leukocytes in control and HCCs at intermediate and end-stage (Control-4-weeks n = 5, Myc^OE/Trp53^KO intermediate n = 4 and end-stage n = 5, Myc^OE/Pten^KO intermediate n = 5 and end-stage n = 5, Control-15-weeks n = 5, Nras^G12D/Pten^KO intermediate n = 7 and end-stage n = 7 (for CD8T and CD4T n = 6); and Nras^G12V/Pten^KO intermediate n = 4 and end-stage n = 5). Statistical analyses of each model are performed comparing each tumor stage to its relative control group (Control 4-weeks for Myc^OE/Trp53^KO and Myc^OE/Pten^KO; and Control-15-weeks for Nras^G12D/Pten^KO and Nras^G12V/Pten^KO). d Quantification of CD11B, CD15, CD204 and S100A9 positive cells by immunohistochemistry in paraffin-embedded HCC patient samples from the Wu et al. dataset segregated according to p-ERK1/2 (positive n = 99, negative n = 369 for CD11B, CD15 and CD204; positive n = 98, negative n = 367 for S100A9) and MYC (positive n = 78, negative n = 389 for CD11B, CD15 and CD204; positive n = 77, negative n = 386 for S100A9) expression in cancer cells. e Unsupervised hierarchical clustering of the transcriptome of CD45⁺ cells isolated from genetically-distinct HCCs (n = 3 for all genotypes) and their classification according to the human HCC immune subtypes clustering. f Unsupervised hierarchical clustering of the ssGSEA enrichment scores per control liver and HCC models (end-stage) using immune-related pathways presented in the Biocarta database. The color scale represents the z-score normalized enrichment per pathway (row) between HCC models (Control n = 5, Myc^OE/Trp53^KO n = 3, Myc^OE/Pten^KO n = 3, Nras^G12D/Pten^KO n = 4, Nras^G12V/Pten^KO n = 5). Graphs show mean ± SEM (a, c, d), +SEM (b). Statistical significance was determined by unpaired two-sided Student’s t-test (a–c), two-sided Mann-Whitney U test (d), one-sided Fisher’s test using Bonferroni multiple testing correction (e). See Supplementary Fig. 9j, k for gating strategy (a–c). Source data are provided as a Source Data file.

Fig. 4. The Nras^G12D/Pten^KO TME displays transcriptomic heterogeneity with mixed pro-tumorigenic and inflammatory signatures.

a Uniform Manifold Approximation and Projection (UMAP) representation of CD45⁺ immune cells from control liver (n = 1 mouse; 7447 cells) and Nras^G12D/Pten^KO HCC (n = 1 mouse; 4814 cells) (left) with annotated populations identified by scRNA-seq (right, Supplementary Data 7). b K-nearest neighbor (KNN) graph (left) and dotplot (right) depicting the differential abundance of cell types between Nras^G12D/Pten^KO HCC relative to control. Each dot represents a group of cells clustered in ‘neighborhoods’. Colors represent significant logFC (FDR ≤ 0.05), whereas white is non-significant difference in abundance. Edge thickness represents the number of overlapping cells between neighborhoods. c Violin plots depicting the normalized expression levels of inflammatory and immunosuppressive genes in the indicated immune cell subsets in control liver (n = 1) and Nras^G12D/Pten^KO HCC (n = 1) (Supplementary Data 8). d UMAP representation of the ‘Monocytic cell’ subset (2616 cells) grouped with cDC1 (115 cells) from (a). e K-nearest neighbor (KNN) graph (left) and dotplot (right) depicting the differential abundance of myeloid cell subsets between Nras^G12D/Pten^KO HCC and control liver from cells grouped as ‘neighborhoods’ in (d). The colored dots represent significant changes in abundance using a threshold of FDR ≤ 0.05, whereas white is non-significant difference in abundance. Edge thickness represents the number of overlapping cells between neighborhoods. f Violin plots depicting the normalized expression levels of inflammatory and immunosuppressive genes in the indicated myeloid cell subsets in control livers and Nras^G12D/Pten^KO HCC from the ‘Monocytic cell’ population (Supplementary Data 8). Nhood = neighborhood. Statistical significance was determined by two-sided, Wilcoxon rank sum test with Bonferroni multiple testing correction (c, f).*** FDR ≤ 0.001. Source data are provided as a Source Data file.

Fig. 5. Nras^G12D/Pten^KO cancer cells activate GM-CSF signaling in the TME.

a–c Differential expression analyses from bulk tumor (Fig. 1e) and cancer cell RNA-seq datasets comparing Nras^G12D/Pten^KO (n = 3) and Nras^G12V/Pten^KO HCC (n = 3). a DEGs between bulk tumors (y-axis) and cancer cells (x-axis). Colored dots represent genes significantly deregulated in both datasets. b Enriched signaling pathways for genes in (a). c Overlap of up-regulated genes and proteins identified by RNA-sequencing and cytokine arrays, respectively. d GM-CSF levels in conditioned media of genetically-distinct HCC cells (Myc^OE/Trp53^KO n = 4, Myc^OE/Pten^KO n = 4, Nras^G12D/Pten^KO n = 6, Nras^G12V/Pten^KO n = 4). e GM-CSF signature enrichment (Csf2rb^KO-WT) in DEG from transcriptional signatures of Nras^G12D/Pten^KO (n = 4) and Nras^G12V/Pten^KO (n = 3) bulk tumors (top), and Nras^G12D/Pten^KO (n = 3) and Nras^G12V/Pten^KO (n = 3) CD45⁺ cells (bottom). Colors represent DEG overlap from Nras-driven models and GM-CSF signature. f DEGs from CD45⁺ cell RNA-seq of Nras^G12D/Pten^KO (n = 7) relative to Nras^G12V/Pten^KO (n = 7) overlapping with GM-CSF signature. Colors represent direction of GM-CSF regulation, relevant genes from GM-CSF signature are depicted (Supplementary Data 9). g Csf2ra, Csf2rb and Ccl6 expression from scRNA-seq (Fig. 4a). h GM-CSF signature enrichment for the indicated immune cell subsets identified by scRNA-seq. i. Survival of TCGA:LIHC patients segregated into high correlation with Nras^G12D/Pten^KO or Nras^G12V/Pten^KO transcriptional signatures (Fig. 2c) (Nras^G12D/Pten^KO n = 65, Nras^G12V/Pten^KO n = 66). high- and low-correlating patients are unique to each signature. Tables show patient numbers at the indicated time points. j Enrichment of GM-CSF positively-regulated genes (n = 628) from the GM-CSF signature in TCGA:LIHC patients segregated into high correlation with Nras^G12D/Pten^KO or Nras^G12V/Pten^KO transcriptional signatures (Fig. 2c) (Nras^G12D/Pten^KO n = 65, Nras^G12V/Pten^KO n = 66). k Survival of TCGA:LIHC patients segregated into high/low enrichment of GM-CSF positively-regulated genes (n = 628) from the GM-CSF signature (High n = 103, Low n = 82). Tables show patient numbers at the indicated time points. Graph show mean ± SEM (d). Statistical significance: two-sided hypergeometric test with Benjamini-Hochberg multiple testing correction (b, e, h), one-way ANOVA with Tukey’s multiple comparison test (d), log-rank test (i, k), unpaired two-sided Student’s t-test (j). Vertical lines at −log₁₀(FDR) = 2.5 indicate significance threshold (b, e, h). The shading represents 95% confidence interval (i, k). Source data are provided as a Source Data file.

Fig. 6. GM-CSF expression is regulated through activation of the ERK1/2 pathway and SP1 transcription factor.

a Barplot depicting the expression levels (signal intensity) of RAS/MAPK-associated phospho-proteins in Nras^G12D/Pten^KO and Nras^G12V/Pten^KO HCC cell lysates (n = 1 per cell line). b Barplot depicting GM-CSF protein levels quantified in Nras^G12D/Pten^KO cancer cell conditioned media after 24 h of treatment with Trametinib (MEK1/2 inhibitor; 0 nM n = 5, 5 nM n = 5, 10 nM n = 4), Temuterkib (ERK1/2 inhibitor; 0 μM n = 3, 5 μM n = 3, 10 μM n = 3) and Vx-11e (ERK2 inhibitor; 0 μM n = 6, 5 μM n = 5, 10 μM n = 4) at the indicated drug concentrations. c Barplots depicting the GM-CSF protein level quantified in the supernatant of scramble (shScr, n = 9), shErk2 (n = 4), and shErk1 (n = 4) Nras^G12D/Pten^KO HCC cell lines. d Venn diagram depicting the overlap of DEG between Nras^G12D/Pten^KO and Nras^G12V/Pten^KO cancer cells treated or not with Vx-11e (n = 3 per cell line per condition) (Supplementary Data 9). e Graphical representation of the gene expression pattern (y-axis) across conditions (x-axis) that follow Csf2 regulation (d; 1270 overlapping DEG) shown (d). “+” indicates up-regulated genes and “-“ indicates down-regulated genes. f Motif enrichment analysis in promoters from the 373 genes, identifying SP1 as a top candidate (Supplementary Data 9). g Barplots depicting the GM-CSF protein levels quantified in the supernatant of shScr (n = 9, shown in (c)) and shSp1 (n = 3) Nras^G12D/Pten^KO HCC cell lines. h Scatterplot depicting the correlation between the enrichment of ERK signaling pathway (y-axis) and GM-CSF signature (x-axis) for each TCGA cancer type. Correlation analyses were performed on the median pathway enrichment score of all patients per cancer type. Graph shows mean ± SEM (b, c, g). Statistical significance was determined by two-sided Student’s t-test (b, c, g), two-sided Fisher’s exact test, followed by Benjamini-Hochberg multiple test correction (f) and two-sided test for association between paired samples, using Pearson’s product moment correlation coefficient (h). The shading represents 95% confidence interval (h). Source data are provided as a Source Data file.

Fig. 7. GM-CSF blockade curbs monocyte-derived inflammatory Ly6C^low cell accumulation in the Nras^G12D/Pten^KO TME.

a, b Percentage of Ly6C^lowF4/80^low cells relative to total myeloid cells (a) obtained from bone marrow (BM) cells differentiated in either recombinant M-CSF (n = 6) or GM-CSF (n = 6), or in conditioned media (CM) prepared from distinct HCC cell lines, with or without GM-CSF neutralizing antibody (a-GM-CSF) (Nras^G12D/Pten^KO n = 7, Nras^G12D/Pten^KO + a-GM-CSF n = 6, Nras^G12V/Pten^KO n = 6, Nras^G12V/Pten^KO + a-GM-CSF n = 5) and (b) analyzed for the indicated phenotypic markers relative to the total Ly6C^low F4/80^low population. c IL-1R signaling activity in HEK reporter cells exposed to CM from BM cells differentiated with either recombinant M-CSF (n = 7), GM-CSF (n = 7), or to CM from Nras^G12D/Pten^KO or Nras^G12V/Pten^KO cell lines, in presence or not of a-GM-CSF (Nras^G12D/Pten^KO n = 7, Nras^G12D/Pten^KO + a-GM-CSF n = 7, Nras^G12V/Pten^KO n = 7, Nras^G12V/Pten^KO + a-GM-CSF n = 3). d Experimental design: mice were HDTV-injected to induce Nras^G12D/Pten^KO HCC, monitored by weekly MRI starting from 3 weeks post-injection, and enrolled into treatments with a-IgG2a (12.5 mg/kg three times per week) or a-GM-CSF (12.5 mg/kg three times per week) for 2 weeks (time point T1) or until end-stage (time point T2) for flow cytometry analyses. e Percentage of intratumoral Ly6C^high monocytes, Ly6G^high neutrophils, Ly6C^lowF4-80^high/int macrophages and Ly6C^lowF4-80^low subsets relative to total CD45⁺ leukocytes in HDTVi-induced Nras^G12D/Pten^KO HCC-bearing mice 2 weeks post treatment (T1) with a-IgG2a (n = 5) or a-GM-CSF (n = 5). f Percentage of Ki67⁺ intratumoral Ly6C^lowF4/80^low cells in HDTVi-induced Nras^G12D/Pten^KO HCC upon 2 weeks of treatment (T1) with a-IgG2a (n = 5) or a-GM-CSF (n = 5). g Mean fluorescence intensity of the depicted markers in four different Ly6C^lowF4-80^low subsets (CD206^lowCCR2^-, CD206^lowCCR2⁺, CD206⁺CCR2^high, and CD206^highCCR2^high) identified by FlowSOM analysis performed on HDTVi-induced Nras^G12D/Pten^KO end-stage (T2) HCCs (n = 3). h Proportions of Ly6C^lowF4-80^low subsets (CD206^lowCCR2^-, CD206^lowCCR2⁺, CD206^intCCR2^high, and CD206^highCCR2^high) identified by FlowSOM analysis performed on HDTVi-induced Nras^G12D/Pten^KO end-stage (T2) HCCs (n = 3). Graphs show mean ± SEM (a, c, e, f) and median (b). Statistical significance was determined by unpaired two-sided Student’s t-test in (a–c, e, f). See gating strategy in Supplementary Fig. 7b (a) and 9j (e-f). Source data are provided as a Source Data file.

Fig. 8. GM-CSF blockade hampers tumor growth and synergizes with VEGF inhibition to extend Nras^G12D/Pten^KO HCC-bearing mice survival.

a Kaplan–Meier survival curves of LOI-induced Nras^G12D/Pten^KO HCC-bearing mice treated with a-IgG2a (n = 8; median survival of 14 days) or a-GM-CSF (n = 10; median survival of 23 days). b Kaplan–Meier survival curves of LOI-induced Nras^G12D/Pten^KO HCC-bearing mice treated with a-IgG2a (n = 9; median survival of 14 days), a-GM-CSF (n = 14; median survival of 25 days), a-VEGF + a-PD-L1 (n = 8; median survival of 21 days), a-GM-CSF + a-PD-L1 (n = 4; median survival of 24 days), a-GM-CSF + a-VEGF (n = 6; median survival of 39.5 days). Graph includes a-IgG2a (n = 8) and a-GM-CSF (n = 10) treated mice from Fig. 8a. c Representative H&E and cleaved caspase 3 (CC3⁺) IHC staining performed on liver sections from end-stage LOI-induced Nras^G12D/Pten^KO HCC-bearing mice treated with a-IgG2a, a-GM-CSF or a-GM-CSF + a-VEGF. (Scale bars = 2 mm; representative CC3⁺ for data shown in d and representative H&E for data shown in (e)). d Barplot depicting the percentage of cleaved caspase 3 (CC3⁺) positive cells per area in tumor nodules (end-stage) isolated from LOI-induced Nras^G12D/Pten^KO HCC-bearing mice treated with IgG2a (n = 9), a-GM-CSF (n = 5) or a-GM-CSF + a-VEGF (n = 6). e Barplot depicting the percentage of necrotic area in tumor section (end-stage) from LOI-induced Nras^G12D/Pten^KO HCC-bearing mice treated with IgG2a (n = 7 mice), a-GM-CSF (n = 11 mice) or a-GM-CSF + a-VEGF (n = 6 mice). Graph shows mean ± SEM (d, e). Statistical significance was determined by log-rank test (a, b) and unpaired two-sided Student’s t-test (d, e). Source data are provided as a Source Data file.