Immunogenomic analysis of human brain metastases reveals diverse immune landscapes across genetically distinct tumors

Abstract

Brain metastases (BrMs) are the most common form of brain tumors in adults and frequently originate from lung and breast primary cancers. BrMs are associated with high mortality, emphasizing the need for more effective therapies. Genetic profiling of primary tumors is increasingly used as part of the effort to guide targeted therapies against BrMs, and immune-based strategies for the treatment of metastatic cancer are gaining momentum. However, the tumor immune microenvironment (TIME) of BrM is extremely heterogeneous, and whether specific genetic profiles are associated with distinct immune states remains unknown. Here, we perform an extensive characterization of the immunogenomic landscape of human BrMs by combining whole-exome/whole-genome sequencing, RNA sequencing of immune cell populations, flow cytometry, immunofluorescence staining, and tissue imaging analyses. This revealed unique TIME phenotypes in genetically distinct lung- and breast-BrMs, thereby enabling the development of personalized immunotherapies tailored by the genetic makeup of the tumors.

Keywords: T cells; brain metastasis; cancer immunology; genomics; immunogenomics; microglia; monocyte-derived macrophages; neutrophils; transcriptomics; tumor microenvironment.

Graphical abstract

Figure 1

The genomic landscape of lung-to-brain metastases (A) Schematic overview of the study. BrM, brain metastasis; PBMC, peripheral blood mononuclear cells; WES, whole-exome sequencing; LP-WGS, low-pass whole-genome sequencing; RNA-seq, RNA sequencing; IF, immunofluorescence; CD45-neg, CD45-negative cells; MDM, monocyte-derived macrophages. (B) Oncoplot summarizing genomic features of lung-BrMs in the immunogenomics cohort (n = 21, in red), including mutational load (upper panel), mutations (central panel), and mutation prevalence (central left panel) of mutated lung cancer drivers; copy number variation (CNV) status of selected genes (second central panel, only CNV alterations identified by both WES and LP-WGS are reported here); summary of individual nucleotide changes in all mutated genes (third central panel); mutational signatures (fourth central panel) and clinical information, including time from the end of the last chemo/radio/immunotherapy treatment prior to surgical resection of the BrM (i.e., > 6 months indicates that the corresponding treatment was finished more than 6 months before surgery), duration of dexamethasone treatment, sex, and histological subtype of the primary tumor (bottom panels). Radiation treatment modality corresponds to stereotactic radiosurgery (SRS) in all radiotherapy-treated samples. (C) Comparison of tumor mutational burden (TMB) between lung-BrM samples in this study (n = 21) and 33 primary cancer types from The Cancer Genome Atlas (TCGA). Each dot represents the mutational burden of one tumor; total number of tumors per primary cancer is indicated above the plot. (D) Total number of mutations in lung-BrMs (TP53^mut n = 9; TP53^mutKRAS^mut n = 5; KRAS^mut n = 2; NP/NK n = 5 biological replicates; one-way ANOVA test p = 0.0002; Dunnett’s multiple comparisons test ∗∗adjusted p < 0.01, ∗∗∗∗adjusted p < 0.0001). (E) Total number of neoantigens in lung-BrMs (TP53^mut n = 9; TP53^mutKRAS^mut n = 5; KRAS^mut n = 2; NP/NK n = 5 biological replicates; ANOVA/Kruskal-Wallis test p = 0.025; Dunn’s multiple comparisons test ∗adjusted p < 0.05, ∗∗∗adjusted p < 0.001). See also Table S1. Table S1A, Table S1B, and Table S1C, Table S2. Annotated variants across the immunogenomics cohort (n = 30), related to Figures 1, 2, and S1, Table S3. Annotated variants across the validation cohort (n = 18), related to Figures 1, 5, and S6, Table S4. Copy number variation estimation from WES data, related to Figures 1, 2, and S3, Table S5. Copy number variation validation by LP-WGS, related to Figures 1, 2, and S3, Table S10. Gene set enrichment analysis (GSEA) and overrepresentation analysis (ORA) in the lung-BrM (n = 21) and breast-BrM (n = 9) cohorts, related to Figures 3, 4, 6, and S7.

Figure 2

The genomic landscape of breast-to-brain metastases (A) Oncoplot summarizing genomic features of breast-BrMs (n = 9) in the immunogenomics cohort (panels as in Figure 1B). (B) Comparison of TMB between breast-BrMs from this study (n = 9, in red) and 33 primary cancers from the TCGA. (C) Rainfall plots of kataegic breast-BrMs. Y axis depicts the distance between consecutive mutations; red arrows indicate genomic regions where kataegis was detected. (D) Total number of mutations (left panel) and neoantigens (right panel) in breast-BrMs (kataegic, n = 4; non-kataegic, n = 5 biological replicates; unpaired two-tailed Mann-Whitney test ∗p < 0.05; data are represented as mean ± SEM). (E) Estimation of intratumoral heterogeneity in breast-BrMs (see STAR Methods for details; kataegic, n = 4; non-kataegic, n = 5 biological replicates; two-tailed unpaired Mann-Whitney test ∗p < 0.05). See also Table S1. Table S1A, Table S1B, and Table S1C, Table S2. Annotated variants across the immunogenomics cohort (n = 30), related to Figures 1, 2, and S1, Table S3. Annotated variants across the validation cohort (n = 18), related to Figures 1, 5, and S6, Table S4. Copy number variation estimation from WES data, related to Figures 1, 2, and S3, Table S5. Copy number variation validation by LP-WGS, related to Figures 1, 2, and S3, Table S10. Gene set enrichment analysis (GSEA) and overrepresentation analysis (ORA) in the lung-BrM (n = 21) and breast-BrM (n = 9) cohorts, related to Figures 3, 4, 6, and S7.

Figure 3

Characteristic immune landscapes in genetically distinct lung-BrMs (A) Flow cytometry (FCM) quantification of non-immune cells (CD45-neg), lymphocytes (CD45⁺, CD11B⁻), and myeloid cells (CD45⁺, CD11B⁺) in lung-BrMs (NP/NK n = 5; KRAS^mut n = 2; TP53^mut n = 9; TP53^mut;KRAS^mut n = 5 biological replicates; Kruskal-Wallis plus Dunn’s multiple comparisons test ∗adjusted p < 0.05). (B) Mean of 14 different immune cell populations analyzed in lung-BrMs as percentage of total CD45⁺ cells. NK, natural killer; DN, double-negative; Treg, regulatory T cells; mono, monocytes; MDM, monocyte-derived macrophages; DC, dendritic cells; iMC, immature myeloid cells; see STAR Methods for details on the markers used to define individual populations. (C) Summary of differential gene expression analyses by DESeq2 (contrasts: TP53^mut vs. NP/NK; KRAS^mut vs. NP/NK; TP53^mut;KRAS^mut vs. NP/NK). Only genes with an adjusted p value ≤0.05 and absolute fold-change ≥2 were defined as differentially expressed. Black, total number of differentially expressed genes (DEGs); red, number of upregulated genes; blue, number of downregulated genes. (D) Visualization of intersects of DEGs between genetically distinct BrMs in sorted immune (MDM, microglia, neutrophils, CD8⁺ T cells, CD4⁺ T cells) and non-immune (CD45-neg) populations. Left panel indicates total number of DEGs (up, upregulated genes; down, downregulated genes); number of intersecting genes for each comparison, as encoded in the combination matrix, are indicated above individual bars; numbers below the combination matrix show percentage of intersecting genes for each comparison. (E) Circos representation of selected significantly enriched pathways in sorted immune and non-immune populations (gene set enrichment analysis (GSEA) on MSigDB HALLMARK (H) and Gene Ontology Biological Process (BP) gene sets; adjusted p ≤ 0.05). Outer rings show genotype of the BrM samples (left half), and whether enrichment scores are positive (red, indicates enrichment in upregulated genes) or negative (blue, indicates enrichment in downregulated genes) (right half). See also Tables S6A and S7A.

Figure 4

Transcriptional analysis of the immune microenvironment of TP53 mutant lung-BrMs (A) Volcano plot representation of DEG in TP53^mut vs. NP/NK lung-BrMs in sorted immune and non-immune cell populations. Horizontal and vertical red lines indicate adjusted p value (≤0.05) and fold-change thresholds (≥2) respectively for genes to be considered differentially expressed. (B) Normalized enrichment score (NES) of selected gene sets from the MSigDB hallmark (H) and Gene Ontology BP collections in sorted CD45-neg cells, MDMs, neutrophils, CD8⁺ T cells, and microglia (adjusted p value ≤ 0.05). (C and D) Normalized counts (TPM, transcripts per million) of T cell activation-related genes (C) and interferon response-related genes (D) in sorted CD8⁺ T cells. (E) Normalized counts of selected genes in sorted microglia; please note log₁₀ scale (C–E, TP53^mut n = 9; NP/NK n = 5 biological replicates; adjusted p values from differential expression analysis by DESeq2; ∗adjusted p < 0.05; data are represented as mean ± SEM). See also Tables S6A and S7A.

Figure 5

Kataegic breast-BrMs present a distinct immune microenvironment (A) Unsupervised clustering of different immune cell fractions (of total CD45⁺ immune cells) as quantified by FCM. (B) FCM quantification of the abundance of selected immune populations in kataegic and non-kataegic breast-BrMs (kataegic, n = 4; non-kataegic, n = 5 biological replicates; unpaired two-tailed t test, ∗p < 0.05, ∗∗p < 0.01; data are represented as mean ± SEM). (C–E) Representative immunofluorescence staining of (C) CD8⁺ T cells, (D) CD4⁺ T cells, and (E) CD15⁺ neutrophils in breast-BrM tumors. (F) Quantification of CD8⁺ T cells, CD4⁺ T cells, and neutrophils from immunofluorescence stainings of kataegic and non-kataegic breast-BrMs from the immunogenomics and validation (Val) cohorts (CD8⁺ T cells: kataegic, n = 6; non-kataegic, n = 10 biological replicates; CD4⁺ T cells: kataegic, n = 6; non-kataegic, n = 14 biological replicates; neutrophils: kataegic, n = 7; non-kataegic, n = 10 biological replicates; unpaired two-tailed Mann-Whitney test, ∗p < 0.05; data are represented as mean ± SEM). Only samples with a false to true positive ratio <2% were included in these analyses (see STAR Methods for details). The two samples closest to the mean of each group were selected as representative for the images shown in (C)–(E). Scale bars represent 100 μm.

Figure 6

Transcriptional analysis of the immune microenvironment of kataegic and non-kataegic breast-BrM tumors (A) NES of selected gene sets from the MSigDB hallmark (H) and Gene Ontology BP collections in sorted CD45-neg cells, MDMs, microglia, neutrophils, and CD8⁺ T cells (adjusted p value ≤ 0.05). (B) Normalized counts (TPM) of selected genes in sorted neutrophils (kataegic, n = 4; non-kataegic, n = 5 biological replicates; adjusted p values from differential expression analysis by DESeq2; ∗adjusted p < 0.05; data are represented as mean ± SEM). (C) GSEA of interferon alpha response genes (MSigDB H collection; upper left panel) and genes involved in positive regulation of type I interferon production (Gene Ontology Biological Process collection; lower left panel) in sorted CD8⁺ T cells; top 40 interferon alpha response genes sorted by fold change (kataegic vs. non-kataegic; numbers in pink indicate log₂ fold changes) from DESeq2 analysis (right panel). (D) Gene set variation analysis (GSVA) of “Early CD8 activation” (as defined in Andreatta et al.⁶¹), AKT, NFKB, and apoptosis pathway-related gene sets from the Reactome database in sorted CD8⁺ T cells. Each dot corresponds to one BrM sample (kataegic, n = 4; non-kataegic, n = 5 biological replicates; unpaired two-tailed t test, ∗p < 0.05; ∗∗p < 0.01; data are represented as mean ± SEM). (E and F) (E) Representative immunofluorescence staining and (F) quantification of CD103⁺CD8⁺ T cells in breast-BrM tumors from the validation cohort (kataegic, n = 3; non-kataegic, n = 8 biological replicates; unpaired two-tailed t test, ∗∗p < 0.01; data are represented as mean ± SD). Only samples with a false to true positive ratio <2% were included in these analyses (see STAR Methods for details). Scale bars represent 100 μm. (G) Heatmap of the average immune cell Danaher scores (computed as in Danaher et al.⁶²) from bulk RNA-seq of breast-BrM samples from the validation cohort (n = 13). See also Tables S6B and S7B.