Spatial transcriptomic analysis of immune checkpoint blockade response in triple negative breast cancers with tertiary lymphoid structures

Abstract

Tertiary lymphoid structures (TLSs) are associated with improved cancer immunotherapy responses. However, TLSs vary in their ability to elicit anticancer immune activity, so it is important to develop databases that allow study of variables that regulate their function. We applied single RNA molecule resolution imaging to longitudinal biopsies taken from women with TLS-enriched triple negative breast cancers prior to therapy, after pembrolizumab and after pembrolizumab plus radiation therapy. We developed a computational framework to align and analyze spatial trajectories between TLSs and tumor beds. Tumors with higher T cell infiltration rates were eradicated after pembrolizumab. In contrast, those with lower malignant cell and T cell interaction rates at baseline showed CXCL9⁺ macrophage infiltration after pembrolizumab, and infiltration of T cells expressing CXCL9-associated programs prior to cancer cell removal after radiation therapy. This manuscript describes single RNA molecule resolution profiling of breast tumors bearing tertiary lymphoid structures throughout an immunotherapy response.

Keywords: Cancer; Cancer systems biology; Immunology.

Graphical abstract

Figure 1

Clustering of cells into spatial niches based on their local environment (A) Tissue region showing segmented cell types (left). myCAF: myofibroblast cancer-associated fibroblast, pDC: plasmacytoid dendritic cell, iCAF: inflammatory cancer-associated fibroblast, cDC: classical dendritic cell, malig. epi: malignant epithelial cell, norm. epi.: normal epithelial cell, endo/PVL: endothelial/perivascular-like cell, macro.: macrophage, corresponding niches (middle) and subniches (right). (B) Bar plot of the percentage of cells assigned to each broad niche (top); Stacked bar plot of the percentage of cells in each niche belonging to each cell type class (bottom). Abbreviations follow those described in (A). (C) Bar plot of the percentage of cells assigned to each subniche (top); Stacked bar plot of the percentage of cells in each subniche belonging to each cell type class (bottom). Abbreviations follow those described in (A). (D) Spatial cross-correlation coefficient for each sub-niche pair at a scale of 500 coordinate units along with subtype proportions of each sub-niche. Spatial sub-niche clusters are annotated at the top. Co-clustering rate at the bottom was determined by computing these spatial cross-correlations on 10,000 samples bootstrapped at the scale of individual FOVs and averaging the rate at which each sub-niche was in the same spatial cluster as each other sub-niche (L); Heatmap of cell subtype proportion in each subniche wherever a proportion exceeded 0.05 (R). myCAF-detox: detoxification-involved myofibroblast cancer-associated fibroblast, norm.-basal epi.: normal basal epithelial cell, norm.-lupr epi.: normal luminal HR + epithelial cell, b5: BLIMP1⁻IGHG1⁺ plasma cell, b6: BLIMP1⁺IGHG1⁺ plasma cell, b2: memory B cell, t2: naive CD4⁺/CD8⁺ T cell, t5: follicular helper T cell, t6: effector/exhausted CD8⁺ T cell, m11: antigen-associated macrophage, b0: naive B cell, b8: BLIMP1⁻IGHG1⁺ plasma cell, t0: naive CD4⁺/CD8⁺ T cell, t1: natural killer cell, BEC-arterial: blood endothelial cell, norm.-luma epi.: normal luminal epithelial cell, m2: Type I Interferon high macrophage, t3: regulatory T cell, m9: classical dendritic cell 2, m4: acute inflammation macrophage, iCAF-IFNG: interferon gamma inflammatory cancer-associated fibroblast, m3: IL-8 macrophage in addition to abbreviations described in (A).

Figure 2

Trajectory analysis in identified TLS and tumor bed regions reveals pathway and signaling gradients moving inward, outward, and between regions (A) Images of B cell zones, T cell zones, and malignant cell zones. Background images show the density of these cells in the same region of tissue, and the white outlines show the location of each connected cluster of cells. epithelial malig.: malignant epithelial cell. (B) Example of the determination of the intermediate zone bridging a B cell zone and malignant cell zone. First, we compute the gradient fields of the convolved images of each zone. Then, we calculate the angle between the gradients at each pixel, setting a minimum cutoff of 0.8π radians for the intermediate zone. BCZ: B cell zone, MCZ: malignant cell zone. (C) Proportions of broad cell type (top), sub-niche (middle), and cell subtype (bottom) moving through each zone from the core of the B cell zone on the left, through the intermediate zone in the center, and ending at the core of the malignant cell zone at the right. We include only those trajectories which account for at least 10% of any one spatial bin. Asterisks indicate statistical significance of a Mann-Whitney U test between proportions in the bracketed regions; ∗ ≤ 0.1, ∗∗ ≤ 0.05, ∗∗∗ ≤ 0.01, ∗∗∗∗ ≤ 10⁻³, ∗∗∗∗∗ ≤ 10⁻⁴. Inter. zone: Intermediate zone, Malig. cell zone: malignant cell zone, macro.: macrophage, malig. epi.: malignant epithelial cell, cDC: classical dendritic cell, T_FH: follicular helper T cell, T1 IFN mac.: type 1 interferon macrophage, migrat. DC: migratory DC, anti. pres. mac: antigen-presenting macrophage, myCAF: myofibrolast cancer-associated fibroblast, PVL: perivascular-like cell. (D) Proportions of broad cell type (top) and subtype (bottom) moving through T cell zones (TCZ) connected to B cell zones (BCZ). Asterisks indicate statistical significance of a Mann-Whitney U test between B cell proportions in the bracketed region. Abbreviations and statistics described in (C). (E) Pathways associated with gene expression programs in each cell type in each zone found through differential expression between each cell type in each zone and the rest. GSEA FDR <0.05. IZ: intermediate zone, TCZ: T cell zone, norm. epi.: normal epithelial cell in addition to abbreviations described in (B–D). (F) Radar plots describing ligand-receptor analysis results for this differential expression. Each plot corresponds to a cell type/zone target population. Color of bars represents cell type sender, radius represents maximal NicheNet score across all ligands sent by cell type, dots represent median NicheNet score across all ligands sent by cell type, and width of bar represents percentage of all interactions coming from cell type. Interactions included here are those with a NicheNet score >0.1. Abbreviations described in (A–C). (G) Example of gradient analysis used throughout this work. We first compute the gradient of the convolved image of one instance of a particular zone. Then, we compute the gradient of expression of each gene in each cell type in this region and divide it by the local density of that cell type. Finally, we take the dot product between these two vector fields to find the angle between them, identifying genes which increase or decrease in expression per cell as you move toward the center of an instance. (H) Pathways associated with gene expression programs found to move together through the previously described gradient analysis. Abbreviations described in (A–E); statistical test described in (E).

Figure 3

Patients separate into two distinct populations based on tumor bed immune infiltration at baseline (A) Box plots of percentage of total cells in each biopsy consisting of niche 4 cells, split by time point. Horizontal lines in the boxes represent the first, second, and third quartiles. Lines extend out to 1.5 times the interquartile range beyond the third quartile and below the first quartile above and below the box, respectively. Points are colored by patient. MCZ: malignant cell zone, pt: patient. (B) Heatmaps showing broad cell type compositions of each instance of malignant cell zones (L) and B cell zones (R). Columns are color-coded by patient and treatment. Cell types highlighted in red indicate where we find a statistically significant difference in infiltrated versus non-infiltrated instances (Mann-Whitney U test p ≤ 0.05). BCZ: B cell zone, macro: macrophage, malig. e.: malignant epithelial cell, endo.: endothelial cell, PVL: perivascular-like cell, iCAF: inflammatory cancer-associated fibroblast, myCAF: myofibroblast cancer-associated fibroblast, cDC: classical dendritic cell, norm. e.: normal epithelial cell, pDC: plasmacytoid dendritic cell. (C) Heatmaps of cell type compositions of each instance of the intermediate zone and T cell zone. IZ: intermediate zone, TCZ: T cell zone in addition to abbreviations and statistics described in (B). (D) Cell type proportions traveling through each zone for infiltrated (solid) and non-infiltrated (dashed) patients at baseline. Trajectories shown are those with a statistically significant difference between infiltrated and non-infiltrated patients in at least two consecutive spatial bins (Mann-Whitney U test; p ≤ 0.05). Asterisks indicate statistical significance of a Mann-Whitney U test between proportions in the bracketed regions (∗ ≤ 0.1, ∗∗ ≤ 0.05, ∗∗∗ ≤ 0.01, ∗∗∗∗ ≤ 10⁻³, ∗∗∗∗∗ ≤ 10⁻⁴). inf: infiltrated, non-inf: non-infiltrated in addition to abbreviations described in (B and C). (E) Sub-niche proportions traveling through the malignant cell zone (MCZ) at baseline for patients of each infiltration status. Abbreviations described in (B and D); statistical tests described in (D). (F) Cell subtype proportions through the malignant cell zone (MCZ). Subtypes shown pass the same test as the previous plot. norm. bas. epi.: normal basal epithelial cell in addition to abbreviations described in (B and D), and statistical tests described in (D). (G) T cell expansion rates from the single cell data of each patient in this study at baseline. NK: natural killer, Treg: regulatory T cell, TFH: follicular helper T cell in addition to abbreviations in (A).

Figure 4

Pathway and signaling differences between infiltrated and non-infiltrated patients at baseline (A) Gene expression programs found from differential expression between immune infiltrated and non-infiltrated biopsies at baseline in each zone. GSEA FDR <0.05. BCZ: B cell zone, IZ: intermediate zone, MCZ: malignant cell zone, TCZ: T cell zone, cDC: classical dendritic cell, endo.: endothelial cell, norm. epi.: normal epithelial cell, malig. epi.: malignant epithelial cell, iCAF: inflammatory cancer-associated fibroblast, myCAF: myofibroblast cancer-associated fibroblast, macro.: macrophage, PVL: perivascular-like cell. (B) Pathway analysis of genes identified as having a difference between infiltrated versus non-infiltrated tumor in their gradient dot products with specific TLS zone principal axes. pDC: plasmacytoid dendritic cell in addition to abbreviations and statistics described in (A). (C) Counts per cell of MHC-related genes in infiltrated (solid) and non-infiltrated (dashed) tumors, going through each zone. inf: infiltrated, non-inf: non-infiltrated in addition to abbreviations described in (A). (D) Heatmap showing primary (outer) and secondary (inner) expressors of each gene in each spatial bin via the isolation forest method. Abbreviations described in (A–C). (E and F) Ligand-receptor analysis for differential expression and (F) gene gradient signatures. All results show interactions corresponding with a NicheNet score >0.1. Abbreviations described in (A–C).

Figure 5

T cell interaction and expansion differences between infiltrated and non-infiltrated patients (A) PDL1/2 counts per cell through each zone at baseline for infiltrated (solid) and non-infiltrated (dashed) patients. Isolation forest results for primary (outer) and secondary (inner) broad cell type expression in each spatial bin shown below. inf: infiltrated, non-inf: non-infiltrated, iCAF: inflammatory cancer-associated fibroblast, myCAF: myofibroblast cancer-associated fibroblast, macro.: macrophage, cDC: classical dendritic cell, epi. m: malignant epithelial cell. (B) Proportion of PD1 interacting T cells moving through each zone. Moving through the intermediate zone (IZ) into the malignant cell zone (MCZ), PD1 interactions drop significantly in non-infiltrated patients. Asterisks indicate statistical significance of a Mann-Whitney U test between proportions in the indicated regions (∗ ≤ 0.1, ∗∗ ≤ 0.05, ∗∗∗ ≤ 0.01, ∗∗∗∗ ≤ 10⁻³, ∗∗∗∗∗ ≤ 10⁻⁴). BCZ: B cell zone, IZ: intermediate zone, MCZ: malignant cell zone in addition to abbreviations used in (A). (C) Box plots of densities of PD1 interacting T cells (T cells expressing PD1 and neighboring a cell expressing PDL1/2) and non-interacting T cells in malignant cell zone and B cell zone instances at baseline, per patient. Horizontal lines represent the first, second and third quartile. Lines extend out to 1.5 times the interquartile range from the first and third quartiles. Asterisks indicate statistically significant differences (Mann-Whitney U test) between the B cell zone and malignant cell zone. pt: patient in addition to significance levels, int. T cell den. cont.: interacting T cell density contrast, in addition to abbreviations described in (B). (D) T cell proportions in the b cell zone in non-infiltrated tumors for base, PD1, and RTPD1. Horizontal lines and whiskers are described in (C). Statistics described in (C) and significance levels in (B). (E) Venn diagram of leading-edge genes associated with cytokine and immune pathways for T cells in the malignant cell zone (GSEA; FDR <0.05). (F) Box plots of overall T cell proportions moving through individual instances of the malignant cell zone at each timepoint in non-infiltrated tumors. Asterisks indicate statistical significance of a Mann-Whitney U test between proportions in the bracketed region between the indicated time points. Horizontal lines and whiskers are described in (C). Abbreviations and significance levels described in (B). (G) T cell subtype proportions moving through the malignant cell zone at each timepoint. Subtypes shown must have a statistically significant difference between at least two timepoints in at least two consecutive spatial bins. (Mann-Whitney U test; p ≤ 0.05). Asterisks indicate statistical significance of a Mann-Whitney U test between proportions in the bracketed region between the indicated time. T_FH: follicular helper T cell in addition to abbreviations and significance levels described in (B). (H) Pathway summaries for T cells in the malignant cell zone at PD1 in non-infiltrated patients. Heatmaps show the number of pathway hits for each summary (GSEA; FDR <0.05), number of unique genes in these pathways, GSEA normalized expression score, and FDR. (I) Upregulated genes in malignant cell zone T cells at PD1 compared to baseline (LFC >0.2, FDR <0.05, and percent expressed in either group >0.25). Abbreviations in (B). (J) T cell subtype expansion rates for each patient at each timepoint from the single cell data. NK: natural killer, Treg: regulatory T cell in addition to abbreviations described in (C and G).

Figure 6

Pembrolizumab induces CXCL9 signaling between myeloid and T cell populations in non-infiltrated tumor beds (A) Sub-niche and cell type proportions through non-infiltrated malignant cell zones that show a statistically significant difference between treatments. Horizontal lines represent the first, second, and third quartiles. Whiskers represent 1.5 times the interquartile range from the first and third quartiles. Asterisks indicate statistical significance of a Mann-Whitney U test between proportions in the bracketed regions between the indicated time points or across the entire zone if no brackets are indicated (∗ ≤ 0.1, ∗∗ ≤ 0.05, ∗∗∗ ≤ 0.01, ∗∗∗∗ ≤ 10⁻³, ∗∗∗∗∗ ≤ 10⁻⁴). malig. epi.: malignant epithelial cell, norm. epi.: normal epithelial cell, endo.: endothelial cell, cDC: classical dendritic cell, MCZ: malignant cell zone. (B) Summaries of pathways found through differential expression of non-infiltrated malignant cell zone at baseline and PD1. GSEA; FDR <0.05. iCAF: inflammatory cancer-associated fibroblast, myCAF: myofibroblast cancer-associated fibroblast, macro.: macrophage in addition to abbreviations described in (A). (C) Bar plot of number of differentially expressed genes in malignant cell zones (LFC >0.2, FDR <0.05, and percent expressed in either group >0.25) between baseline and PD1 assigned to each cell type by the isolation forest analysis. Note the abundance of pathways in cDCs at PD1. Since cDCs are largely absent in non-infiltrated malignant cell zones at baseline, these pathways were missed in the previous differential expression analysis. PVL: perivascular-like cell, pDC: plasmacytoid dendritic cell in addition to abbreviations found in (A and B). (D) Summaries of pathways describing genes (Enrichr; FDR ≤ 0.05) brought into malignant cell zones at PD1 by cDCs (left); Total number of pathways (middle) and genes (right) defining these summaries are indicated in the heatmap. (E) Scatterplot of the log₁₀p-value versus log₂ fold-change of genes brought into malignant cell zones by cDCs after PD1. Secreted ligands are annotated in red. (F) Ligands whose downstream gene programs overlap with gene expression changes seen in malignant cell zone normal epithelial and T cells after pembrolizumab, and the most likely cell type to be overexpressing those ligands after therapy. Abbreviations described in (A and B). (G) Overlap of differentially expressed genes in T cells at PD1 compared to baseline in malignant cell zones, and the genes in the highest 10% of weights in the CXCL9 NicheNet table. P-value results from the hypergeometric test for overlap of these two gene sets. Pathway summaries listed are those that describe both gene sets via Enrichr (FDR ≤ 0.05). (H) Dot plots of log-transformed, normalized CXCL9 expression in all cell types in biopsies containing infiltrated and non-infiltrated malignant cell zones in the single cell data. Abbreviations in (A–C). (I) Dot plots of log-transformed, normalized CXCL9 expression for each zone and infiltration status in cDCs and macrophages. BCZ: B cell zone, IZ: intermediate zone in addition to abbreviations in (A and B). (J) Dot plot of cDC and myeloid CXCL9 expression in all subniches. (K) Upregulated chemokine and cytokine ligands in PD1 malignant cell zones (LFC >0.2, FDR <0.05, and % expression in either group >0.15) where the corresponding receptor is expressed in cDCs (nonzero expression in >30% of cells). Colored boxes show the broad cell type isolation forest results for their expressers. Abbreviations described in (A–C). (L) NicheNet potentials of CXCL9 for each PD1-enriched cytokine and chemokine. (M) Scatter plot of the -log₁₀p-values versus log₂ fold-changes of genes comparing cDCs in the non-infiltrated biopsies in the malignant cell zone from PD1 to baseline. The genes included here are those in the top 5% of weights attributed to IFNG in the NicheNet model and considered differentially expressed based on an FDR <0.05 and log₂ fold-change >0.2.