Clonality, Antigen Recognition, and Suppression of CD8⁺ T Cells Differentially Affect Prognosis of Breast Cancer Subtypes

Numerous immunotherapy approaches are currently being exploited for a variety of human malignancies, including hematologic as well as solid tumors. These approaches generally include oncolytic viral therapy, cancer vaccines, adoptive T-cell therapy, and application of checkpoint inhibitors (CI). Particularly, the latter two have demonstrated impressive objective response rates (OR) of up to 80%, including several complete responses in advanced disease stages (1, 2).

Breast cancer has initially been considered poorly immunogenic due to its low average mutational burden when compared with other tumor types (3). More recently, it has been acknowledged that some breast tumors are extensively infiltrated by immune cells (4), and it became evident that density of tumor-infiltrating lymphocytes (TIL), in particular CD8 T cells, has prognostic value and predicts response to neoadjuvant chemotherapy as well as immune-modulating therapies (5–9). Building on the revisited immunogenicity, several studies are currently exploiting cancer vaccines, adoptive T-cell therapies, or CI for the treatment of breast cancer (10). Unexpectedly, ORs remain variable, and generally do not exceed 20% for CI monotherapy (11).

Breast cancer is a heterogeneous disease comprising several histologic and molecular subtypes. The most well-recognized subtypes include luminal A and normal-like (largely resembling the histologic phenotype: ER⁺, PR⁺ and HER2⁻), luminal B (ER⁺ PR⁺ KI67hi and/or HER2⁺), her2 (ER⁻, PR⁻, HER2⁺), and basal-like subtype [largely resembling the triple-negative (TN) phenotype: ER⁻, PR⁻, HER2⁻; refs. 12, 13]. This subtype classification has clinical relevance with respect to prognosis and choice of targeted therapies (13). Notably, it has been observed that TN tumors respond better to CI treatment when compared with ER⁺ breast cancer (14). Nevertheless, responses to immune therapies are not restricted to TNBC patients, as it has been reported that a metastatic luminal A breast cancer patient showed complete regression following adoptive T-cell therapy (15).

Immune parameters that can be decisive toward an effective antitumor response, such as those reflecting immunogenicity as well as occurrence of T-cell evasion, are poorly characterized in case of breast cancer, and thus critical factors determining tumor immunogenicity are poorly understood. Tumor immunogenicity, which is the extent to which adaptive immune responses are triggered, depends on multiple factors including the expression, processing, and presentation of tumor antigens, and the presence, type, and antigen specificity of TILs. Immunogenic tumor antigens can include oncoviral antigens, cancer germline antigens (CGA), and neoantigens, which all have been reported to elicit T-cell responses in cancer patients (16, 17). Besides immunogenicity, immune evasive mechanisms can affect numbers and antitumor activity of T cells. Several evasive mechanisms have been described that limit influx and migration (e.g., lack or downregulated expression of chemoattractants and/or adhesion molecules), antigen recognition (e.g., lack or downregulated expression of molecules involved in antigen processing and/or presentation), and/or function of CD8 T cells (e.g., presence of immune-suppressor cells, altered expression of immune or metabolic checkpoints, and/or activation of oncogenic pathways; refs. 18 and 19).

Current reports on the prognostic and predictive value of TILs in breast cancer subtypes to a certain extent refer to some of the above mechanisms yet remain inconclusive and sometimes even contradictory (5, 20). In the current study, we have comprehensively assessed immunogenicity as well as T-cell–evasive mechanisms in breast cancer subtypes with respect to characteristics of TILs and the tumor microenvironment (TME). To this end, we applied a series of omics analysis tools on large publicly available data sets and demonstrated that breast cancer molecular subtypes significantly differ with respect to the extent of immunogenicity and occurrence of dominant T-cell–evasive mechanisms, qualities that go beyond the mere presence of TILs. These novel findings may aid future immune therapy trial design and rationalize implementation of differential combination immune therapies for breast cancer subtypes.

Gene-expression data (HG-U133-A array) of lymph node–negative (LNN) breast cancer patients (who did not receive any adjuvant systemic treatment) were retrieved through GEO Series accession numbers GSE2034 (n = 286; ref. 21), GSE5327 (n = 58; ref. 22), GSE11121 (n = 200; ref. 23), GSE2990 (n = 125; ref. 24), and GSE7390 (n = 198; ref. 25). To assess the prognostic value of TILs in different breast cancer subtypes, GSE2034 and GSE5327 constituted a discovery cohort (n = 344), whereas GSE11121, GSE2990, and GSE7390 (n = 523) constituted a validation cohort. Raw CEL files of all GSE entries were downloaded, and data were normalized using fRMA (26) and corrected for batch effects using ComBat (27). Subsequently, to retain statistical power, the combined cohort of 867 LNN primary breast cancer was used for in-depth analyses of TILs and the TME within breast cancer subtypes as described below.

Second, a unique data set was retrieved from primary breast cancer cases with both whole-genome sequencing (WGS; n = 560) and in-depth RNA-seq data (n = 347) and which is accessible through EGAS00001001178 (28, 29). From cases for which both WGS and RNA-seq were available (n = 266), the data were used for predicting neoantigens and from the cases which had just RNA-seq data (n = 347), the data were used for analysis of T-cell clonality.

This study has been approved by the medical ethical committee at Erasmus MC (MEC.02.953), and was performed according the Declaration of Helsinki and the “Code for Proper Secondary Use of Human Tissue in The Netherlands” (FMWV, version 2002, update 2011) of the Federation of Medical Scientific Societies in the Netherlands (http://www.federa.org/), the latter aligning with authorized use of coded spare tissue for research.

To assess TIL abundance, we used a previously published (30), but slightly modified TIL signature. In brief, the signature was built based on the assessment of the proportion of TIL nuclei of total nuclei (including TIL, tumor, and resident stromal cells) for multiple representative areas of hematoxylin and eosin (H&E)–stained slides from 96 familial breast cancer samples by an experienced pathologist (C.H.M. van Deurzen). Gene-expression data (Affymetrix HG-U133_plus_2.0 array) from GEO54219 (31) of corresponding samples were split in two groups based on TIL abundance (high and low TIL count, median split) and tested for differential gene expression, which resulted in a 152-probe signature that highly correlates with TIL percentages in the specimen (r = 0.74, P < 0.001; ref. 30). From the 152 HG-U133_plus_2.0 array immune infiltrate signature probe sets, 120 (resembling 109 genes, see additional file 1) were found to overlap with the HG-U133_A array (cohort A) and were used to classify samples based on TIL abundance according to average linkage hierarchical clustering (correlation as distance metric) in this study. TIL-high samples contained a median stromal TIL count of 35%, whereas TIL-low samples had a median stromal TIL count of 5%. The limitation to 120 probes in this study did not affect TIL scores calculated from either RNA-seq or microarray data (correlations between original and modified signatures, r = 1, P < 0.001).

The prognostic value of this signature according to distant metastasis-free survival (MFS) was tested with cohort A. In the discovery set, the mean age at the time of surgery was 53 years [standard deviation (SD), 12]; 221 patients (64%) were ER-positive; 120 patients were assigned to the TIL-high cluster (35%); 198 patients were premenopausal (58%). T1 tumors (≤2 cm) were present in 168 patients (49%), T2 tumors (>2–5 cm) in 163 patients (47%), T3/4 tumors (>5 cm) in 12 patients (3%), and unknown tumor stage in 1 patient. With respect to disease spread, 226 patients (66%) did not develop metastasis at a distant organ during follow-up (median follow-up time of patients still alive was 101 months; range, 61–171 months). In the validation set, 404 patients were ER-positive (77%); 166 patients were assigned to the TIL-high cluster (32%); and 387 patients (74%) did not have metastasis at a distant organ during follow-up (median follow-up time of 124 months). An overview of clinical characteristics per molecular subtype is given in Supplementary Table S1. To reduce interexperimental variation and minimize biases, both data sets were combined by batch mean centering (32), and subsequently used for hierarchical clustering to divide samples in high and low TIL abundance groups using the 120 probe sets of the TIL signature.

Samples were assigned to subtypes according to breast cancer intrinsic molecular subtypes (AIMS) as described by Paquet and colleagues (12) using the Bioconductor R package genefu (33). AIMS subtyping is considered the most stable subtyping based on expression data, as it is not affected by normalization or subtype frequencies in the cohort. It is noteworthy that AIMS subtyping has good concordance with PAM50 molecular subtypes as well as histologic subtypes and preserves prognosis of subtypes (12). For a breakdown of hormone receptor status (ER, PR, and Her2) per molecular subtype, see Supplementary Fig. S1 and for a comparison of several immune gene sets in AIMS and histologic subtypes see Supplementary Fig. S7.

T-cell receptor-β chain (TRB) reads were extracted from RNA-seq (cohort B, n = 347) using the miXCR (34) algorithm available at https://github.com/milaboratory/mixcr. In brief, the software aligns the sequencing reads to reference V, D, J, and C genes of the T-cell receptors (TCR), assembles the clonotypes, and exports the clones per sample. Because total TCR-Vβ reads significantly differed among breast cancer subtypes, we used only TIL-high samples to compare TCR clonality.

DNA-seq and RNA-seq data (cohort B, n = 266) were used to determine load of neoantigens and expression of CGAs. Neoantigens and epitopes were predicted as described in ref. 29. Briefly, 17-mers of amino acids containing an amino acid arising through a nonsynonymous mutation at the center position were run through the online prediction server Net-MHC to predict EC₅₀ values of all possible 9-mer peptides for HLA class I molecules, and a peptide with a predicted EC₅₀<50 nmol/L was considered a possible neoepitope. A list of CGAs (n = 239) was extracted from the CT database (http://www.cta.lncc.br/), and their expression levels were determined relative to expression of all CGAs per patient.

Microarray samples (cohort A, n = 867) were used for deconvolution of 22 immune cell populations (LM22) using CIBERSORT ((ref. 35; https://cibersort.stanford.edu/). Prognostic values of immune cell populations toward MFS were assessed by Cox regression analysis following classification into “high” and “low” groups, split by the median frequency.

Genes related to T-cell evasion (n = 850, see additional file 2) were selected from reports from others (36–38) as well as the Laboratory of Tumor immunology (PI: R. Debets), Department of Medical Oncology, Erasmus MC (reviewed in refs. 18 and 19). These genes represent different T-cell–evasive mechanisms, and can be divided into three main categories, namely, genes related to (i) influx and migration of T cells, (ii) antigen processing and presentation, and (iii) function of T cells, with each category having various subcategories (see additional file 2). Notably, 4 of 24 genes of the category “immune checkpoints” were also present in the TIL signature. Differential gene-expression analysis was performed for cohort A in R using Limma (39). Genes that were differentially expressed among breast cancer subtypes (Benjamini–Hochberg adjusted P < 0.05) were grouped according to above-mentioned categories of T-cell evasion. Gene sets with a clear direction (>90% upregulated or downregulated) were displayed using violin plots and heat maps. Heat maps were made using average expression of gene sets per breast cancer subtype. For comparison of gene sets expressed by T cells, only TIL-high samples were considered.

Differential gene expression was tested using Kruskal–Wallis tests, and comparisons of gene-expression levels versus basal-like breast cancer were performed with Wilcoxon tests. Distant MFS was used as endpoint for prognosis, and log-rank tests were used to test survivor functions. Cox regression analyses were performed to assess the prognostic value of gene sets (and single genes; see additional file 3) in a univariable as well as multivariable (stratified for tumor size and age) model for the entire cohort A and per individual subtype. An interaction test was performed for all covariates with molecular subtypes. Kaplan–Meier curves were used to plot survival probabilities of subtypes/patients selected based on gene expression. Stata v13 (StataCorp) was used to calculate differences, and two-sided P values of < 0.05 were considered statistically significant.

We have built a 109-gene (120-probe) signature (additional file 1) that captures stromal TIL abundance, according to the guideline method from Salgado and colleagues (ref. 40; Fig. 1A–E). Because TILs can be highly heterogeneous with respect to their composition of T cells and other immune cells, we performed immune fluorescence staining followed by digital image analysis and observed strong correlations between the TIL score and numbers of CD8 T cells (r = 0.82, P < 0.0001) as well as CD4 T cells (r = 0.74, P < 0.0001), but not macrophages (r = −0.031, P = 0.9; Supplementary Fig. S2). We further assessed this signature using the omics tools miXCR (34) and CIBERSORT (35). The first tool, which enumerates the sequence reads of T-cell receptor alpha and beta chains (TCRα and β) from RNA-seq data (cohort B), revealed that the number of TCR-Vβ reads in a specimen strongly correlated with the TIL score (computed as average of all 109 genes per sample) of the same specimen (r = 0.91, P < 0.0001; Fig. 1D). The second tool, a deconvolution algorithm that extracts relative proportions of 22 major immune cell populations from microarray data (cohort A), revealed that proportions of CD8 T cells, activated memory CD4 T cells, naïve B cells, gamma/delta T cells, and M1 macrophages positively correlated with TIL scores, whereas proportions of resting DC, neutrophils, monocytes, M0 macrophages, and activated mast cells negatively correlated with these scores (Supplementary Fig. S3). Collectively, these findings point out that TIL-high samples (median of 35% sTIL) are enriched for activated lymphocytes and reflect frequencies of αβ T cells, and that TIL-low samples (median of 5% sTIL) are enriched for immune-suppressor cells.

Next, we assessed the prognostic value of our TIL signature and observed an association with distant MFS in the discovery set of 344 LNN primary breast cancer (P = 0.009; Supplementary Fig. S4A). To validate our findings, we performed similar analyses in a data set of 523 primary operable, LNN breast cancer and verified a significant association between TIL abundance and MFS (P = 0.016; Supplementary Fig. S4B). Following this validation, and to increase statistical power for subgroup analysis, we combined the discovery and validation sets of samples (cohort A, n = 867) and explored the association between TIL scores and MFS in the different molecular breast cancer subtypes (for hormone receptor status per molecular subtypes, see Supplementary Fig. S1). Basal-like and her2 had equally high TIL scores followed by normal-like breast cancer, whereas luminal A and B subtypes showed the lowest TIL scores (P < 0.0001; Fig. 2A). TIL scores had prognostic value in univariate Cox regression only in basal-like (log-rank P = 0.002), her2 (log-rank P = 0.04), and luminal B (log-rank P = 0.02) subtypes, but not in luminal A or normal-like subtypes (Fig. 3B–F). Notably, TIL scores remained significant for basal-like (P = 0.025) and luminal B (P = 0.02), but not her2 (P = 0.2) in a multivariable Cox regression including tumor size and age. To investigate the biological basis of differential prognostic values of TILs, we subsequently evaluated TIL and TME characteristics among breast cancer subtypes.

Because T-cell clonality can be indicative for tumor reactivity of TILs (41), we assessed TCR-Vβ reads in TIL-high samples of different breast cancer subtypes. Basal-like breast cancer showed the highest number of different T-cell clones (i.e., highest TCR repertoire diversity; average number = 109), which was significantly higher than her2 (61), normal-like (50), luminal A (39), and luminal B breast cancer (35; Fig. 3A). In agreement, basal-like breast cancer showed the highest read counts of the 10 most abundant clones per sample (Fig. 3B). Highly expanded clones (HEC), which were defined as clones with read counts >10% of total clone reads, were present in basal-like, her2, and luminal B tumors. Interestingly, luminal B but not A tumors harbored HECs (Fig. 3C), whereas both tumor types had equally low TCR repertoire diversity. Notably, in all subtypes, individual cases with expanded clones (EC, 5%–10% of total clone reads) were present.

To evaluate antigenicity, we assessed expression of two categories of recognized targets for CD8 T cells, namely, neoantigens and CGAs. We used nonsynonymous mutations (cohort B, see Materials and Methods) and evaluated CGA expression (cohort B, average expression of 239 CGAs per sample). Neoantigen expression was significantly higher in basal-like breast cancer (average number = 21.5) compared with luminal A (7) and B (17), and normal-like breast cancer (15), but not compared with her2 breast cancer (17.5) (P < 0.001, Fig. 3D). CGA expression was again significantly higher in basal-like breast cancer compared with luminal A and B and normal-like breast cancer, but not her2 breast cancer (P < 0.002; Fig. 3F). When correlating antigen expression to TIL scores, we observed a slight positive correlation with neoantigen expression (r = 0.2, P = 0.0017; Fig. 3E). Interestingly, CGA expression was inversely correlated to TIL score (r = −0.21, P = 0.00082; Fig. 3G).

Besides parameters of immunogenicity, we evaluated differential expression of genes related to three main categories of T-cell evasion (influx and migration of T cells; antigen processing and presentation; and function of T cells; Fig. 4). With respect to influx and migration of T cells, basal-like breast cancer showed the highest expression of chemoattractants, which was significantly lower in other subtypes. Normal-like tumors showed the highest expression of T-cell adhesion genes (followed by basal-like breast cancer) and those related to cancer-associated fibroblasts and extracellular matrix products (Fig. 4A). With respect to antigen recognition by T cells, basal-like breast cancer showed the highest expression of antigen processing and presentation (APP) genes, which was significantly lower in all other subtypes (Fig. 4B). Type I IFN gene products, which are recognized for their effects toward antigen priming (42), were equally high in basal-like, her2, and luminal B, but expressed to significantly lower levels in luminal A and normal-like subtypes. With respect to T-cell function, we evaluated the expression of costimulatory ligands and receptors (for costimulatory and inhibitory receptors that are expressed by T cells, we compared TIL-high samples only), immune mediators of the TME, components of oncogenic pathways (see additional file 2 for details), and frequencies of immune (-suppressor) cells (based on CIBERSORT deconvolution). Again, basal-like tumors showed the highest expression of coinhibitory receptors and ligands as well as costimulatory receptors and ligands. Expression levels of coinhibitory/stimulatory receptors were significantly lower in all other subtypes, whereas expression of corresponding ligands was equally high in her2, but not in other subtypes. Expression levels of mediators such as IL10, IDO-1, VEGF, components of adenosine and glycolysis pathways, which can all limit T-cell function, were also highest in basal-like breast cancer, and significantly decreased in all other subtypes. Interestingly, although M0 macrophages were significantly enriched in basal-like breast cancer, frequencies of regulatory T cells, monocytes, and M2 macrophages were highest in luminal B tumors (Fig. 4C).

When interrelating TIL score with above-mentioned gene sets in a subtype-independent manner, we observed that this score correlated significantly with coinhibitory/stimulatory receptors and NFKB pathway (r > 0.7), which may be expected because all of these gene sets can be expressed by T cells. Weaker correlations (r < 0.6) were observed between the TIL score and antigen processing, and IFN signatures. For other gene sets, we observed weak (adenosine, immune mediators), no (glycolysis), or even inverse correlations (M0, WNT pathway) with TIL scores (see Supplementary Fig. S5 for details; interactions between gene sets and TIL score and subtypes are indicated in Fig. 4).

Finally, we performed Cox regression analyses with MFS to determine the prognostic value of the above gene sets in univariable and multivariable settings (Tables 1 and 2; for prognostic value of single genes, see additional file 3). Besides TIL scores, we found significant associations with MFS [hazard ratio (HR) <1, P < 0.05] for the following gene sets: type II IFN in breast cancer (not differentiated per subtype), as well as her2 and luminal B subtypes; coinhibitory/stimulatory receptors and ligands in breast cancer, as well as basal-like, her2 and luminal B subtypes; APP in breast cancer, as well as basal-like and luminal B subtypes; and NFKB pathway in breast cancer as well as basal-like subtype. In addition, activated memory T cells were significantly associated with MFS in breast cancer and basal-like, in the multivariable analysis. Interestingly, we found significant inverse associations with MFS (HR > 1, P < 0.05) for M2 macrophages in breast cancer; M0 macrophages in normal-like subtype; glycolysis and adenosine pathway in breast cancer; and WNT pathway in the breast cancer, basal-like, and luminal A subtype.

CD8-positive TILs in breast cancer are generally associated with a favorable clinical course, yet response to immune therapies is overall low, does not go hand-in-hand with CD8 TIL abundance, and differs significantly across subtypes. Tumor–immune interactions are not well understood and are expected to contribute to lack of response to immune therapies. In this study, we have addressed immunogenicity and T-cell–evasive mechanisms according to multiple parameters in a cohort of 867 LNN, untreated primary breast cancer, and a cohort of 437 primary breast cancer with WGS data. To this end, we have built a 109-gene TIL signature that preserves differential prognostic value of TILs across molecular breast cancer subtypes. We found that T-cell clonality, antigenicity, frequency of immune-suppressor cells as well as expression of genes related to influx/migration, antigen recognition or suppression of T cells provide significant determinants of individual subtypes and may explain the differential prognostic value of TILs.

The TIL signature used in this study was based on quantitative pathologic assessment of TILs on HE-stained breast cancer (30). TILs generally comprise a variety of different immune cell population including subsets with protumor and antitumor activity (5). Nevertheless, we found that this signature correlated very well with numbers of CD8 T cells and frequencies of activated lymphocytes (cell types associated with good prognosis) and inversely correlated with immune-suppressor cells (cells associated with poor prognosis), suggesting that it may be used to assess T-cell abundance when tissues are not available. Moreover, we observed a near-perfect correlation between TIL score and TCR-Vβ read counts, suggesting that TCR-Vβ read counts also act as a surrogate for TIL abundance in silico [correlations between TIL density and either TIL score (r = 0.82, P < 0.0001) or TCR-Vβ read counts (r = 0.76, P < 0.0001) are comparable]. TIL scores were high in basal-like, her2, and normal-like breast cancer. T-cell abundance, however, was not a prerequisite for prognosis because TIL scores were only associated with MFS in basal-like, her2 (univariate analysis only), and Luminal B subtypes, but not in luminal A, nor normal-like subtypes.

Assessment of immunogenicity and T-cell–evasive mechanisms per subtype revealed that basal-like, her2, and normal-like tumors harbor the highest TCR repertoire diversity. TCR clonality, however, was highest in basal-like, her2, and luminal B tumors, which suggests that a tumor-specific T-cell response has taken place in the latter three breast cancer subtypes and to a lesser extent in luminal A and normal-like tumors (despite high TIL score and TCR diversity in normal-like tumors). Even though T-cell clonality assessment based on bulk RNA-seq analyses has several limitations (43) and should be interpreted with caution, it is interesting that a recent study by Scheper and colleagues is in line with our findings and revealed that the tumor-reactive capacity of TILs is highly variable in different cancers and that a proportion of TILs represents true bystander T cells (44), which may be the case for the majority of TILs in normal-like tumors. In line with T-cell clonality, basal-like breast cancer showed the highest expression of neoantigens and CGAs, followed by her2 and luminal B tumors. These observations suggest that a certain level of antigen expression, which is considered a prerequisite for an antitumor T-cell response, and the presence of clonally expanded T cells, which is a consequence of T-cell responses, may at least in part explain the prognostic value of TILs in basal-like, her2, and luminal B subtypes. Nevertheless, when correlating neoantigens to TIL scores or T-cell clonality, we only observed a weak correlation, which did not hold true upon subtype stratification (except for luminal B tumors, which may be due to overrepresentation of APOBEC-driven neoantigens that are considered more immunogenic than other neoepitopes; ref. 29), suggesting that, at least in other subtypes, predicted antigens are not truly immunogenic. Undoubtedly, the landscape of antigens goes beyond CGAs and classic neoantigens (derived from nonsynonymous mutations), to which end we have evaluated correlations between TILs and alternative mutations (i.e., frameshifts, indels, drivers, passengers) which did not improve correlations with TILs (data not shown). The lack of strong correlations between antigenicity and TILs may hint to the occurrence of immune editing or other immune evasive mechanisms.

When evaluating genes related to T-cell evasion, we observed that expression of chemoattractants was highest in basal-like breast cancer followed by her2. T-cell adhesion molecules, on the other hand, were highest in normal-like breast cancer, which may explain the relatively high TIL scores as adhesion molecules may enhance T-cell retention (45) irrespective of antigen specificity. Vice versa, low expression levels of adhesion molecules in luminal B tumors may explain the low TIL score measured for this subtype. Looking into antigen recognition by T cells, we demonstrated that APP was lowest in luminal A, -B, and normal-like tumors. Nevertheless, we observed a significant association of the APP gene set with survival in luminal B tumors, suggesting that APP is functional in at least a subset of these patients. Loss-of-function mutations in APP genes have been reported to result in immune resistance toward CI therapies in melanoma (37, 46). Along this line, it is noteworthy that we found one or more mutations in APP genes in 10% of basal-like patients, a fraction of patients that is 2- to 5-fold higher compared with the other breast cancer subtypes (see Supplementary Fig. S6). Type I interferons, such as IFNα and IFNβ, produced by either tumor cells or dendritic cells, are critical for priming of CD8 T cells and significantly affect natural as well as therapy-induced immune control of tumors (42). Moreover, it has been shown that mice lacking type I IFNs spontaneously develop breast tumors (47). Expression of type I IFN genes was equally high in basal-like, her2, and luminal B tumors and showed significant interaction with the TIL score, which may further support the prognostic value of TILs in these subtypes. Interestingly, we observed that the abundance of memory CD4 T cells was associated with better survival in breast cancer and basal-like breast cancer, which has also been shown recently by others (48). We found that this T-cell subset was enriched in TIL-high samples of basal, her-2, and luminal B tumors, which suggests that not only quantity but also quality of the T-cell infiltrate matters.

Lastly, when analyzing modes of T-cell suppression, we observed that gene sets for both coinhibitory/stimulatory receptors showed significant interaction with TIL scores and were associated with survival in basal, her2, and luminal B tumors (as were TIL scores), which is in agreement with the concept that immune checkpoints become expressed following a tumor-specific T-cell response, because these molecules are often expressed in response to IFNγ (49). M0 and M2 macrophages, cell types that have been associated with various modes of immune suppression and tumor progression in breast cancer (50), were inversely correlated with MFS in breast cancer, which is in agreement with findings by Ali and colleagues (51). Immune inhibitory mediators, which can be secreted by myeloid-derived suppressor cells (MDSC) or tumor cells, were also expressed at the highest level in basal-like breast cancer, and as a gene set demonstrated a high HR in this subtype. Next to suppressor cells, oncogenic and metabolic pathways have also been linked to T-cell evasion. In example, increased glycolysis has been linked to decreased trafficking and cytotoxicity of T cells in other malignancies (52). Interestingly, it has been observed that glycolysis signaling can induce MDSC development, via induction of G-CSF and GM-CFS in breast cancer (53). In line, we observed coexpression of MDSC and glycolysis signatures as well as inverse correlations with TIL scores. Notably, these observations could be caused by a single common process or due to the coexistence of multiple biological processes, indicating that mechanistic conclusions should be drawn with caution. Components of oncogenic pathways, such as WNT and adenosine, were expressed at the highest level in normal-like and basal-like breast cancer, respectively, and remarkably WNT genes were associated with poor prognosis in breast cancer, basal-like breast cancer, and luminal A breast cancer. WNT has been linked to decreased recruitment of T cells in melanoma and it has been shown recently, that WNT signaling can modulate PD-L1 expression in cancer stem cells of breast cancer (54), suggesting that WNT plays a role in multiple modes of T-cell evasion.

In Fig. 5, outcomes of parameters of TILs and TME are summarized per breast cancer subtype, and pointing to the differential qualities of TILs and occurrence of T-cell–evasive mechanisms, which may represent important biomarkers when monitoring patient responses to immune therapies. Based on our multiparameter analyses in relation to subtypes, TIL scores and prognosis, we argue that in particular T-cell clonality, expression of coinhibitory molecules, type I IFN, APP, and frequencies of activated memory T cells better reflect the qualities of TILs. In agreement, recent results of the TONIC trial showed that T-cell clonality correlated with response to pembrolizumab, when combined with chemotherapy, in TNBC (55). A recent trial in which trastuzumab-resistant Her2⁺ patients were treated with pembrolizumab showed clinical benefit for patients whose tumors were positive for PD-L1 (56), a marker that is frequently upregulated following antitumor T-cell responses (49). Moreover, the Impassion130 trial, in which 902 TNBC patients were treated with atezolizumab with/without nab-paclitaxel, revealed that only patients with PD-L1–positive immune infiltrate showed clinical benefit (57). In line with our study, these findings argue that the mere presence of CD8 TILs is not sufficient to accurately chart immune responses, and that markers that reflect the quality of TILs should be included to monitor future trials. In extension to this concept, our findings on the quality of TILs and the immune microenvironment, in particular with regard to subtype-specific differences in metabolic and oncogenic pathways as well as recruitment of different types of suppressor cells, may aid the design and translational testing of subtype-specific (combination) therapies (see examples below). Of note, despite good concordance between subtyping methods (12) and highly similar results for immune gene analysis in histologic and corresponding molecular subtypes (Supplementary Fig. S7), the subtypes presented here do not fully resemble the subtypes used in daily clinical practice. Moreover, we have used gene-expression data at bulk level to interrogate breast cancer subtypes for their immunogenicity and immune evasive mechanisms, which subsequently need to be validated in future studies.

Based on this study, in basal-like breast cancer, there is a strong rationale for testing combinations of checkpoints with inhibitors of oncogenic pathways, MDSC, or immune mediators. Possible combinations include drugs targeting VEGF or adenosine receptors, which are FDA approved for other implications and showed synergistic effects with CI in preclinical models and early clinical studies, including TNBC (58). Furthermore, several WNT-PW inhibitors, which are currently in clinical trials (including studies in breast cancer) or drugs depleting MDSC (59), which are mainly in preclinical development, may increase efficacies of ICI in basal-like breast cancer.

In Her2 breast cancer, even though generally being immunogenic, there may be insufficient numbers of tumor-reactive T-cell clones when compared with basal-like breast cancer (despite equally high TIL scores). Given the high antigenicity and relatively low frequencies of immune-suppressor cells, adoptive T-cell therapy with TCR-engineered T cells may be considered to treat this tumor type. In fact, chimeric antigen receptor (CAR) T cells targeting her2 effectively kill breast cancer in mouse models; however, in patient studies, severe on-target toxicities have been observed (60). Nevertheless, based on our data, CGAs and neoantigens are frequently expressed in her2 breast cancer and may represent safe targets for AT.

In luminal B breast cancer, at least a subset of tumors may very well be immunogenic and may benefit from immune therapies. In fact, TIL scores and expression of APP genes, although being low, were prognostic, and TIL scores correlated to the presence of neoantigens in this tumor type. The generally low frequencies of T cells, which may be due to low expression of chemoattractants and adhesion molecules, argue in favor of treatments that enhance accumulation of T cells (chemotherapy or epigenetic drugs), antigen presentation (CDK4/6 inhibitors; ref. 61), and/or adoptive T-cell therapy. In addition, luminal B tumors were enriched for M2 macrophages and regulatory T cells, which also represent targets in this tumor type for combinatorial treatments (62).

In luminal A breast cancer and normal-like breast cancer, we observed the least signs of immunogenicity and immune evasion implying that these tumor types are unlikely to respond to immune therapies. A recent case study, on the other hand, reported complete regression of a luminal A patient following treatment with tumor-reactive TILs (15). Importantly, however, this patient had an exceptionally high number of mutations, indicating that these results may not be generally translatable to other luminal A tumors.

Finally, normal-like breast cancer demonstrated enrichment for CAF and their products as well as T-cell adhesion molecules, which is suggestive for enhanced T-cell retention and may explain relatively high TIL scores. This, together with the low expression of antigens and markers for ongoing immune responses, suggests that normal-like tumors are also unlikely to respond to immune therapies, but may benefit from therapies that target CAFs and/or their products.

Our data suggest that not frequencies of TILs per se, but rather qualities of TILs, such as T-cell clonality, T-cell subset distribution, APP, expression of type I IFNs and immune checkpoints, and immune microenvironments, in particular oncogenic- and metabolic pathways as well as types and frequencies of suppressor cells, discriminate breast cancer subtypes. Furthermore, our data suggest that evaluation of multiple immune parameters using NGS data enables charting of immunogenicity and immune evasive mechanisms, and provide a guide to select combinatorial approaches that may enhance the efficacy of future immune therapy trials.

R. Debets reports receiving commercial research grants from MSD and speakers bureau honoraria from Cambridge Healthtech Institute, and is an advisory board member/unpaid consultant for Bluebird, Catenion, Genticel, and Intellia. No potential conflicts of interest were disclosed by the other authors.

Conception and design: D. Hammerl, M.P.G. Massink, H.E.J. Meijers-Heijboer, R. Debets, J.W.M. Martens

Development of methodology: J.W.M. Martens

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): C.H.M. van Deurzen, R. Debets, J.W.M. Martens

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): D. Hammerl, M.P.G. Massink, M. Smid, R. Debets, J.W.M. Martens

Writing, review, and/or revision of the manuscript: D. Hammerl, M.P.G. Massink, M. Smid, C.H.M. van Deurzen, H.E.J. Meijers-Heijboer, Q. Waisfisz, R. Debets, J.W.M. Martens

Study supervision: R. Debets, J.W.M. Martens

This study was supported by the Dutch Cancer Society (KWF project nr: 2014-7087).