Understanding and reversing mammary tumor-driven reprogramming of myelopoiesis to reduce metastatic spread

Abstract

Tumor-induced systemic accumulation and polarization of neutrophils to an immunosuppressive phenotype is a potent driver of metastasis formation. Yet, how mammary tumors reprogram granulopoiesis at the molecular level and when tumor imprinting occurs during neutrophil development remains underexplored. Here, we combined single-cell, chromatin and functional analyses to unravel the tumor-driven reprogramming of granulopoiesis in the bone marrow, along with intervention studies aimed at reversing this process. We observe that mammary tumors accelerate commitment to the neutrophil lineage at the expense of lymphopoiesis and erythropoiesis without stimulating the development of a novel myeloid lineage. Moreover, tumor-directed immunosuppressive imprinting of neutrophils starts early in hematopoiesis. Treatment with anti-IL-1β normalizes tumor-induced granulopoiesis, reducing neutrophil immunosuppressive phenotype and mitigating metastatic spread. Together, these data provide molecular insights into the aberrant, tumor-driven neutrophil differentiation pathway leading to metastasis-promoting chronic inflammation and how it can be reversed to reduce metastatic spread.

Keywords: IL-1β; breast cancer; cancer-induced hematopoiesis; immunosuppression; inflammation; metastasis; myelopoiesis; neutrophil programming; neutrophils.

Graphical abstract

Figure 1

Mammary tumors promote myelopoiesis at the expense of erythropoiesis and lymphopoiesis (A) Flow cytometry analysis of peripheral blood neutrophils from healthy controls (gray, n = 63) and patients with metastatic triple negative breast cancer (pink, n = 69) showing neutrophil frequency (left) and neutrophil nitric oxide production (right). (B) Overview of experimental setup for scRNA-seq: cKIT⁺CD34⁺ BM cells were isolated from n = 3 WT and n = 3 KEP mice. Sorted cells were subjected to single-cell RNA sequencing (scRNA-seq), resulting in total n = 11,462 cells after pre-processing. (C) Force-directed graph embedding of local cell neighborhoods from Milo. Neighborhoods were assigned to the cell type which most of the cells within a neighborhood belong to. Differentiation trajectories of various hematopoietic lineages were superimposed on the graph. (D) Differential abundance testing for each cell neighborhood in KEP compared to WT mice from Milo. Color denotes Log₂ ratios of the number of cells in KEP compared to WT, superimposed on the embedding from (C) (left) and summarized for each cell type (right). Neighborhoods with FDR < 0.1 were deemed significant. (E) Flow cytometry analysis of BM population frequencies from WT (n = 9) and end stage KEP (n = 10, tumor size ∼15 × 15 mm²) mice. Data collected across 9 independent experiments. Statistics: Mann-Whitney U test (A), Mann-Whitney with Benjamini, Krieger, and Yekutieli’s multiple comparison correction (E), p and q values <0.05 considered significant. Box and whisker plots: (A and E) boxes represent median and interquartile range; whiskers represent full range. Abbreviations: HSC, hematopoietic stem cell; MPP, multipotent progenitor; GMP, granulocyte monocyte progenitor; GP, granulocyte progenitor; ProNeu, pro-neutrophils; cMoP, common monocyte progenitor; ProMono, pro-monocyte; MHC2-Mo, MHC-II⁺ monocytes; MDP, monocyte dendritic cell progenitor; CLP, common lymphocyte progenitor; MEP, megakaryocyte erythrocyte progenitor; MegaP, megakaryocyte progenitor; EryP, erythrocyte progenitor; Baso/Mast, basophil mast cell progenitor. See also Figures S1–S3.

Figure 2

Mammary tumors accelerate myelopoiesis (A) Validation of inferred cell ordering along the neutrophil trajectory (neutrophil pseudotime). Shown are GAM-smoothed log expression values for Meis1 (HSC-specific gene), Elane (granulocyte-specific gene) and S100a8 (proneutrophil-specific gene) along the inferred trajectory. (B) Kernel density estimates for WT (blue) and KEP (orange) cells across the neutrophil pseudotime. Gray rectangle and triangle denote the extreme values of neutrophil pseudotime where WT cells are absent. (C) Heatmap showing the scaled expression of genes upregulated in cells from KEP mice along neutrophil pseudotime. Genes were grouped into three clusters based on the expression patterns. (D) Enrichment of gene sets for groups of genes in (C), using databases available through Enrichr: GO: Biological Process, Reactome and Tabula muris. Color scale denotes -Log₁₀ transformed adjusted p values. (E) Example expression dynamics of selected genes along neutrophil pseudotime. (F) Incorporation of EdU within BM populations 24 h after EdU injection in WT (blue, n = 6) and KEP (orange, n = 6) mice. Mann-Whitney with Benjamini, Krieger, and Yekutieli’s multiple comparison correction, q value <0.05 considered significant. Data collected across 6 independent experiments. Boxes represent median and interquartile range; whiskers represent full range. (G) Gene set enrichment analysis using gene sets consisting of genes that show expression specific to different hematopoietic lineages, defined from our scRNA-seq data. GSEA was performed on genes ranked based on their log fold changes in KEP vs. WT comparison for merged HSC and MPP cell types, while blocking on the cell type variable. Circle size denotes -Log₁₀ transformed adjusted p values. (H) Proposed model of the aberrant tumor-driven bone marrow hematopoiesis. Tumor signaling causes an inflammatory response at the earliest points of hematopoiesis leading to skewed differentiation toward the granulocyte lineage at the point of transition from MPP to GMP and is accompanied with increased proliferation along the lineage. CD34⁺/cKIT⁺ cells from KEP mice display an accelerated neutrophil signature, upregulating mature neutrophil genes and reaching further along the inferred differentiation trajectory than WT counterparts. HSC, hematopoietic stem cell; CLP, common lymphocyte progenitor; MegaP, megakaryocyte progenitor; EryP, erythrocyte progenitor; ProMono, pro-monocyte; ProNeu, pro-neutrophils. See also Figure S4 and Tables S1, S2.

Figure 3

Dynamics of open chromatin in response to mammary tumors during myelopoiesis (A) Example genomic tracks showing prematurely gained accessibility in KEP progenitor populations at the Prtn3 gene locus. (B) Number of differentially accessible regions (FDR <0.05) in KEP compared to WT progenitors. (C) PCA was performed on 10000 most variable accessible regions for WT (circles) and KEP (triangles) cell populations. Shown are the first two principal components which recapitulate neutrophil differentiation trajectory and suggest accelerated granulopoiesis in KEP. (D) Heatmap showing scaled accessibility dynamics of all differentially accessible regions across. The regions were grouped into five clusters based on their accessibility patterns across cell populations. (E) Average normalized log accessibility for each of the three major clusters in D (clusters 2, 4, 5). (F) Enrichment of differentially accessible regions (clusters 2, 4, 5) in the regulatory domains of cell type marker genes defined using our scRNA-seq data. To define gene domains, we implemented the “basal plus extension” approach as used in GREAT. Significance of the overlap between gene domains and accessible region clusters was tested using Fisher’s exact test with FDR correction. (G) Transcription factor motif enrichment in differentially accessible regions (clusters 2, 4, 5). Highlighted are C/EBP and GATA family transcription factors. Enrichment significance was assessed using Fisher’s exact test with FDR correction. (H) Total red blood cell counts (RBC, left) and red blood cell distribution width standard deviation (RDW-SD) as measured by hemocytometer from WT control (n = 11, blue), KEP control (n = 11, orange). Mann-Whitney U test. Boxes represent median and interquartile range; whiskers represent full range. (I) Normalized bulk-RNA expression of Cebpa and Cebpb genes in WT and KEP cell populations. Statistics represent FDR (q) values from differential expression analysis. See also Figure S5 and Table S3.

Figure 4

Neutrophil immunosuppressive gene signature can be traced back to early progenitors (A) Neutrophil T cell co-culture suppression assay measuring the proliferation of CD8⁺ (above) and CD4⁺ (below) T cells after incubation for 72 h with CD3/CD28 beads either alone (dark gray) or with neutrophils from WT or KEP BM (blue), blood (red), spleen (light gray) or lung (yellow). Data from 4 independent experiments, n = 3–4 mice per genotype. Unpaired t test with Welch-correction. Error bars represent SD. (B) PCA analysis performed on bulk RNA-seq of BM, blood and lung neutrophils from WT (circles) and KEP (triangles) mice, showing separation by tissue along PC1 and by the presence of a tumor along PC3. (C) Volcano plot depicting common up and downregulated differentially expressed genes (DEGs) in BM, blood and lung neutrophils. Genes with FDR < 0.05 were considered significant. (D) Heatmap showing the scaled expression of common differentially expressed genes in (C) across BM, blood and lung neutrophil populations. (E) Neutrophil pseudotime expression dynamics of commonly upregulated genes from (C) which also show upregulation in HSPC populations. (F) Volcano plots showing gene expression changes in KEP mice in each cell population. Genes in the “MDSC signature” obtained from Alshetaiwi H. et al. are highlighted. GSEA was performed using the “MDSC signature” on genes ranked based on their log fold changes. Normalized enrichment scores (NES) and adjusted p values are shown. See also Figure S6.

Figure 5

Treatment with anti-IL-1β normalizes tumor-driven transcriptional re-wiring of progenitor cells (A) Schematic of HSC liquid culture (left) and fold change of gene expression of HSCs cultured with IL-1β over HSCs cultured without IL-1β (n = 4 biological replicates). (B) Schematic overview of experimental set up treating WT and KEP mice with either anti-IL-1β or isotype control antibody twice per week for two weeks (4 injections). (C) Chromatin accessibility log fold changes in KEP mice treated with control antibody (orange) or a-IL-1β (green) in differentially accessible regions within each progenitor cell population. (D) Bulk RNA-seq mRNA expression log fold changes in cells along neutrophil development trajectory from KEP mice treated with control antibody (orange) or a-IL-1β (green) in differentially expressed genes within each cell population. (E) Zoomed-in representation of mRNA expression data from bulk RNA-seq of immunosuppressive genes in neutrophils from KEP mice treated with control antibody (orange) or a-IL-1β (green). (F) Flow cytometry analysis of BM population frequencies from WT control (n = 10, blue), KEP control (n = 10, orange), WT anti-IL-1β (n = 10, dark blue), and KEP anti-IL-1β (n = 10, green) mice. Data collected across 10 independent experiments. (G) Tail vein blood sampling analysis of neutrophil frequency during treatment with anti-IL-1β or isotype control by flow cytometry (n = 10 WT control, n = 10 WT anti-IL-1β, n = 11 KEP control, n = 11 KEP anti-IL-1β). (H) Proportions of blood circulating CD101⁺ neutrophils and cKIT⁺ neutrophils in WT control (n = 9, blue), KEP control (n = 8, orange), WT anti-IL-1β (n = 9, dark blue), and KEP anti-IL-1β (n = 8, green) mice. (I) Neutrophil reactive oxygen species production as measured by ex vivo CellROX green flow cytometry assay (left) and neutrophil nitric oxide production as measured by DAF-FM diacetate fluorescence (right) from KEP control treated mice (orange, n = 7) and KEP anti-IL-1β treated mice (green, n = 7). (J) Total red blood cell counts (RBC, left) and red blood cell distribution width standard deviation (RDW-SD) as measured by hemocytometer from WT (n = 11, blue, repeated from Figure 3H), KEP (n = 11, orange, repeated from Figure 3H), WT anti-IL-1β (n = 9, dark blue) and KEP anti-IL-1β (n = 10, green) mice. (K) Volcano plot depicting differentially accessible chromatin regions in progenitor populations of WT mice treated with a-IL-1β compared to WT mice treated with control antibody. (L) Scaled accessibility dynamics across progenitor populations in all treatment conditions for sites that gain accessibility (upper panel) or lose accessibility (lower panel) in a-IL-1β treated mice as identified in (K). (M) Volcano plot depicting TF motif enrichment in the regions that gain compared to the regions that lose accessibility from (K). Highlighted are SPI and GATA family motifs. Statistics: Mann-Whitney U test (A, G, H, I, and J), Mann-Whitney with Benjamini, Krieger, and Yekutieli’s multiple comparison correction (F) or Fisher's exact test with FDR correction (M). p and q value <0.05 considered significant. Error bars represent SEM (A and G). Box-whisker plots: boxes represent median and interquartile range; whiskers represent full range (C, D, F, H, I, and J). See also Figure S7.

Figure 6

Treatment with anti-IL-1β reduces metastatic spread (A) Schematic of metastasis experimental design including treatment regimen. (B) Neutrophil frequencies in peripheral blood as determined from weekly tail vein bleeding during course of metastasis experiment (n = 3–15 mice per group per time point). (C) Neutrophil reactive oxygen species production as measured by ex vivo CellROX green flow cytometry assay (left) and neutrophil nitric oxide production as measured by DAF-FM diacetate fluorescence (right) from isotype control (gray, n = 7) and anti-IL-1β (red, n = 5) treated mice. (D) Total white blood cell (WBC, left) and red blood cell (RBC, right) counts as measured by hemocytometer (n = 8–10) from isotype control and anti-IL-1β treated mice. (E) Representative lung images of cytokeratin 8-stained lung sections (left) and quantification of lung metastasis and incidence of lymph node metastasis. Statistics are Mann-Whitney U test (B, C, D, and E) and Fisher’s exact test (E) p values <0.05 considered significant. Error bars represent SEM (B, C, and E). Box-whisker plots: boxes represent median and interquartile range; whiskers represent full range (C and D).