Spatially distinct tumor immune microenvironments stratify triple-negative breast cancers

Triple-negative breast cancer (TNBC) lacks targeted therapies, has high rates of distant recurrence, and poor overall survival for patients. Large-scale gene expression and sequencing studies have revealed high heterogeneity within the TNBC subtype and few common actionable targets (1–5). Despite the overall poor outcomes, a subset of patients with TNBC respond well to standard-of-care chemotherapy, implying the existence of distinct TNBC phenotypes (1–5). However, identifying those patients most likely to respond to therapy remains an important clinical challenge.

Increasing evidence indicates that interactions between tumor cells, tumor stroma, and the tumor immune microenvironment (TIME) evolve during the course of disease and play a key role in the response to therapies (6). Some tumors evade immune control, enabling tumor progression. Others are subject to immune attack mediated by the establishment of a Th cell type 1 response that can be subsequently modulated through expression of immunomodulatory ligands (7, 8). Such heterogeneity in the TIME and its evolution throughout tumor progression are still poorly understood.

Importantly, the presence and localization of tumor-infiltrating lymphocytes (TILs) correlate with a better prognosis and an improved response to neoadjuvant chemotherapy in TNBC (9, 10). This understanding has led to the development of guidelines for TIL scoring to harmonize TIL evaluation in breast cancer. However, these guidelines do not elaborate on the importance of TIL location. While stromal TILs (sTILs) constitute the most reproducible parameter when assessed on H&E-stained sections (9), studies suggest an important role for the infiltration of lymphocytes and specifically CD8⁺ T cells, which represent the cytotoxic arm of the adaptive immune response, into the epithelial compartment (10, 11). This suggests that sublocalization of TILs might contribute to the prognostic and even predictive stratification of patients with TNBC.

Elevated expression of the coinhibitory immune ligand programmed cell death ligand 1 (PD-L1, also known as CD274) is associated with the presence of infiltrating lymphocytes (12), supporting the therapeutic value of immune modulation in the TNBC setting through immune checkpoint blockade (ICB). However, in the setting of advanced disease, only 8% to 20% of TNBCs preselected for tumor or immune cell expression of PD-L1 respond to therapy targeting PD-L1 or its receptor programmed cell death 1 (PD-1) (13). To improve these response rates, a better understanding of the various factors that influence differential lymphocyte infiltration and/or activation in TNBC is needed. In addition, gene expression studies have revealed the existence of a paradoxical subset of patients with high expression of immune-associated signatures, who nevertheless experience poor outcomes (14), indicating that additional complexity exists beyond the current information provided by bulk tumor immune signatures.

The tumor microenvironment (TME) plays a crucial role in tumor progression (15), including modulation of the local immune environment (12, 16). Previous studies are limited to bulk tumor–derived expression profiles of TNBC (1–3) and have only identified an immunomodulatory subtype of TNBC associated with good outcomes. Studies of bulk tumors lose information on compartment-specific signals within the tumor core and therefore do not reflect the spatial landscape of the TIME. A recent study identified heterogeneity in TIME architecture using a predetermined set of biomarkers (17). Here, we identify distinct TIME subtypes defined by spatial patterns of CD8⁺ T cell localization and gene expression signatures in therapy-naive TNBC tumors. By integrating spatial characterization of the immune response, using laser capture microdissection–derived (LCM-derived) gene expression profiling data from matched stromal and epithelial tumor compartments, we discovered that each TIME was associated with distinct metasignatures of the TME, prognosis, and biomarkers. The biological processes identified here allowed us to stratify and characterize each TIME subtype, which we believe will support the development of TIME-dependent targeted therapeutic approaches to treat TNBC.

Distinct CD8⁺ T cell localization profiles are observed in TNBC. To define patterns of CD8⁺ T cell localization in TNBC, we performed immunohistochemical analysis to assess the spatial distribution of CD8⁺ T cells in a cohort of 38 therapy-naive TNBC patients using whole sections of formalin-fixed, paraffin-embedded (FFPE) samples (Supplemental Table 1 and Supplemental Figure 1A; supplemental material available online with this article; https://doi.org/10.1172/JCI96313DS1). For each tumor, CD8⁺ T cell density was quantified in 4 distinct compartments: the tumor margin (marCD8) and tumor core (corCD8), the latter of which includes the tumor stroma (strCD8) and tumor epithelium (epiCD8) (Figure 1, A and B, and Supplemental Figure 1B). Importantly, CD8⁺ T cell quantification in the stromal (strCD8) and epithelial (epiCD8) compartments was highly correlated with sTILs and intratumoral TILs (iTILs) (i.e., in the tumor epithelial compartment), as evaluated on matched H&E-stained sections and as defined by the International Immuno-Oncology Biomarker Working Group (9), respectively (Figure 1C). This demonstrates the robustness of our CD8⁺ T cell quantification scheme with respect to standard approaches for clinical TIL evaluation (9).

Figure 1

Therapy-naive TNBC tumors are classified into subtypes on the basis of distinct spatial localization of CD8⁺ T cells. (A) Representative images of CD8⁺ T cell staining at the tumor margins (top panels, dotted lines) and in the tumor core (bottom panels) (n = 38). Scale bars: 100 μm. (B) Quantification of CD8⁺ T cell densities at the tumor margins (marCD8) and in the tumor core (corCD8) (n = 38). (C) Comparison of strCD8 with sTILs and epiCD8 with iTILs (n = 38). Data were analyzed using Spearman’s correlation. (D) Working model of TNBC grouping based on CD8⁺ T cell localization. Black, green, blue, and red represent ID, MR, SR, and FI tumors, respectively. marCD8, corCD8, strCD8, and epiCD8 are the CD8⁺ T cell densities in the tumor margin, core, and stromal and epithelial compartments, respectively. Data represent the mean ± SEM.

Using the spectrum of infiltration of CD8⁺ T cells into the different compartments, we defined TNBC subgroups according to the presence and/or differential localization of CD8⁺ T cells, referred to as TIME subtypes (Figure 1A and stratification diagram in Figure 1D). Tumors were first divided into 2 groups on the basis of their corCD8⁺ T cell infiltration: corCD8^hi and corCD8^lo. The majority of tumors in the corCD8^lo group had an accumulation of CD8⁺ T cells at the tumor margins (marCD8^hi) and were designated as “margin-restricted” (MR) tumors (12 of 16), while a few tumors (4 of 16) displayed a low abundance of CD8⁺ T cells at the margins (marCD8^lo) and were defined as “immune desert” (ID) tumors (Figure 1, A, B, and D, and Supplemental Figure 2). Alternatively, corCD8^hi tumors (n = 22) were divided into 2 subgroups consisting of “fully inflamed” (FI) (11 of 22) tumors, which exhibited significant CD8⁺ T cell infiltration into the tumor epithelial compartment (epiCD8^hi) in addition to their presence in the stroma, and “stroma-restricted” (SR) (11 of 22) tumors, which showed CD8⁺ T cell accumulation in the stroma (strCD8^hi) and exclusion from the tumor epithelial compartment (epiCD8^lo) (Figure 1, A, B, and D, and Supplemental Figure 2).

No significant differences in clinical variables, including tumor size, grade, and lymph node status, were observed among these groups (Supplemental Table 1). On the basis of gene expression profiling of matched bulk tumor specimens (n = 37), we found that 31 of the 37 TNBC samples in this data set belonged to the Prediction Analysis of Micrarray 50–defined (PAM50-defined) basal-like subtype (ref. 5 and Figure 2A). Consistent with the TNBC subtypes (TNBC types) defined by Lehmann et al. (1), we found that corCD8^hi tumors were significantly enriched in the immunomodulatory subtype of TNBC. In contrast, corCD8^lo tumors were significantly enriched in the mesenchymal subtype (Figure 2, B–D).

Figure 2

TIME subtypes compared with PAM50 and Lehmann breast cancer subtype stratifications. (A) Comparison of CD8⁺ T cell grouping (TIME subtypes) with PAM50 molecular subtyping of our TNBC cohort (n = 37). (B–D) Comparison of CD8⁺ T cell grouping (TIME subtypes) with Lehmann molecular subtyping of TNBC (n = 37), showing enrichment of the mesenchymal and immunomodulatory subtypes in corCD8^lo (ID + MR) and corCD8^hi (SR + FI), respectively. Data were analyzed by Fisher’s exact t test.

Metasignatures derived from immune microenvironments reflect different biologies. To assess whether specific biological processes are linked to differential CD8⁺ T cell localization, we profiled gene expression in matched samples derived from bulk tumor (n = 37) and tumor stromal and epithelial compartments, isolated by LCM (n = 38) (bulk: Gene Expression Omnibus [GEO] GSE88847; LCM: GSE88715. To understand the contributions of stromal and epithelial compartments to the biological differences between the TIME-based TNBC subgroups identified, we used GSEA-based metasignatures (MSigs) derived from bulk tumor gene expression data to interrogate LCM-derived tumor stromal and epithelial gene expression data sets (Figure 3). Through this analysis, we were able to more precisely decipher the source of the various biological pathways identified (Figure 4).

Figure 3

CD8⁺ T cell localization–derived metasignature methodology. (A) Analysis workflow for metasignatures and associated biological processes discovery in our data set (n = 37) and validation of the external data set (n = 578). (B) Heatmap on the left shows clustering of the pathway scores, determined according to a positive or inverse correlation with CD8⁺ T cell density in the tumor core (corCD8), and identifies 4 corCD8 metasignatures (corCD8 MSig) in whole-tumor gene expression for all TNBC tumors (n = 37). Heatmap on the right shows clustering of the pathway scores, determined according to a positive or inverse correlation with CD8⁺ T cell density in the tumor epithelium (epiCD8), and identifies 6 epiCD8 metasignatures (epiCD8 MSig). epiCD8 MSigs were generated from whole-tumor gene expression in tumors showing CD8⁺ T cell infiltration into the tumor core (corCD8^hi, n = 22, SR and FI tumors).

Figure 4

CD8⁺ T cell localization–derived metasignatures identify distinct biological processes. (A and B) Cellular pathways positively or inversely correlated with (A) CD8⁺ T cell density in the core (corCD8) of all TNBC tumors (n = 37) and (B) with CD8⁺ T cell density in the epithelium (epiCD8) of TNBC corCD8^hi (SR and FI) tumors (n = 22). Data were analyzed using Spearman’s correlation. Cellular pathway score FDR values are represented in whole tumor (black), tumor stroma (gray), and tumor epithelium (light gray). The top 5 significant pathways per metasignature are represented (except for corCD8 Msig1, for which the top 10 pathways are represented to reflect the dominance of the MSig). MODY, maturity-onset diabetes of the young; CTL, cytotoxic T lymphocyte; ESC, embryonic stem cell; GGPP, geranylgeranyl pyrophosphate; GVHD, graft versus host disease. Pathway names have been abbreviated; full names of the pathways can be found in Supplemental Tables 2 and 3.

To identify pathways associated with each specific pattern of CD8⁺ T cell localization, we followed a 2-step classification scheme (Figure 1D). We first identified pathways that differed between corCD8^hi (SR and FI) and corCD8^lo (ID and MR) tumors by clustering all significant pathways that were correlated (positively or inversely) with corCD8 (FDR <5%). Clustering analysis of pathway enrichments identified 4 metasignatures (corCD8 MSigs; Figure 3B, Figure 4A, and Supplemental Table 2). As expected, the predominant metasignature enriched in corCD8^hi bulk tumors reflected elevated immune signaling (corCD8 MSig1; “Immune”) (Figure 3B, Figure 4A, and Supplemental Table 2). These pathways were equally enriched in both the tumor stromal and tumor epithelial compartments (Figure 3B, Figure 4A, and Supplemental Table 2), consistent with the distribution of CD8⁺ T cells within the tissue. By contrast, the predominant metasignature associated with corCD8^lo bulk tumors included fibrosis and matrix remodeling pathways (corCD8 MSig 3; “fibrosis”) (Figure 3B, Figure 4A, and Supplemental Table 2). Consistent with signatures of reactive stroma, we observed that corCD8^lo tumors were enriched in fibrotic foci (Figure 5), defined as scar-like areas associated with reactive tumor stroma and a poor prognosis (18). Moreover, signaling pathways linked to TGF-β, a key regulator of fibrosis, were enriched in the tumor stroma of corCD8^lo tumors (Supplemental Figure 3). Accordingly, corCD8^lo tumors were enriched for the mesenchymal Lehmann TNBC subtype (Figure 2, B and C). This subtype is characterized by extracellular matrix (ECM) receptor interaction, as well as actin remodeling by Rho and TGF-β signaling (1). Collectively, these results indicate the presence of elevated reactive stroma and matrix remodeling in the poorly infiltrated corCD8^lo (ID and MR) TNBC tumors.

Figure 5

Tumors with poor infiltration of CD8⁺ T cells are enriched for fibrotic foci. (A) Representative image of fibrotic focus in H&E-stained images. Original image size: 7 × 10 mm; enlarged inset size: 0.5 × 0.5 mm. Dotted line identifies the fibrotic focus area. (B) Fibrotic focus presence represented as a 2D plot of marCD8 over corCD8 densities (n = 38). (C) Fibrotic focus presence was enriched in corCD8^lo compared with corCD8^hi tumors (n = 38). Data were analyzed using the Fisher’s exact test. (D) corCD8 was higher in tumors that did not contain a fibrotic focus (n = 38). Data were analyzed using the Mann-Whitney U test. Data represent the mean ± SEM.

Within corCD8^hi tumors, we identified 6 metasignatures (epiCD8 MSigs) reflecting biological processes associated with epiCD8 status (Figure 3B, Figure 4B, and Supplemental Table 3). The metasignatures enriched in FI tumors included JAK/STAT signaling (epiCD8 MSig 1), IFN signaling and cytotoxic activity (epiCD8 MSig 2), and stem cell–linked transcription (epiCD8 MSig 3) (Figure 4B and Supplemental Table 3). In contrast, SR tumors showed enrichment for signatures of mTOR signaling (epiCD8 MSig 4), cholesterol biosynthesis (epiCD8 MSig 5), and IL-17 signaling (epiCD8 MSig 6) (Figure 4B and Supplemental Table 3). Together, these data indicate the association of specific patterns of CD8⁺ T cell infiltration into SR and FI tumors with distinct TMEs.

FI TNBCs are characterized by a type 1 IFN proinflammatory environment. To identify pathways implicated in enhanced CD8⁺ T cell recruitment into the tumor epithelium (epiCD8^hi), we analyzed epithelium- and stroma-specific gene expression to find biological pathways distinguishing tumors with the FI TIME subtype from those of the SR subtype. As shown by the epiCD8 MSig2 signature (IFN signaling and cytotoxic activity), multiple pathways associated with various immune response activities are elevated in FI compared with SR tumors (Supplemental Table 3). This includes genes associated with a type I IFN response (e.g., OASL, ISG15), antigen presentation (e.g., TAP, B2M), cytotoxic activity (e.g., GZMB, FASLG), as well as cell death (e.g., CASP and PARP genes) (Figure 6, A and B, and Supplemental Figure 4).

Figure 6

Fully inflamed TNBC are associated with a proinflammatory TIME. (A and B) Heatmap depicting the expression of genes associated with a type 1 IFN response and cytotoxic activity in the tumor stroma (A) and epithelium (B) (n = 22). (C) Representative images showing a higher number of GzmB⁺ CD8⁺ T cells in FI tumor epithelium compared with numbers in SR tumor. Pan-cytokeratin (Pan-CK) staining (pink) identifies tumor cells, and DAPI (blue) identifies nuclei. White squares outline the position of the zoomed area in the stromal region for SR tumor and the epithelial region for FI tumor. Scale bars: 50 μm (merge, CD8, and GzmB) and 5 μm (enlarged insets showing CD8 and GzmB colocalization). n = 22. (D) Quantification of GzmB⁺CD8⁺ T cells in the tumor core, tumor stroma, and tumor epithelium (n = 32). Data were analyzed using the Kruskal-Wallis test. Data represent the mean ± SEM. (E) GzmB⁺CD8⁺ T cell density in tumor epithelium was positively correlated with epiCD8 (n = 20). Data were analyzed using Spearman’s correlation. Green, blue, and red dots represent MR, SR, and FI tumors, respectively.

Effector CD8⁺ T cells with cytotoxic activity were distinguished by elevated expression of GzmB compared with expression in memory CD8⁺ T cells. To validate the presence of cytotoxic CD8⁺ T cells in FI tumors, we quantified GzmB staining in the tumor stroma and epithelium in samples from our 38 patients. We found that the density of GzmB⁺CD8⁺ T cells was highest in the epithelial compartments of FI tumors, while it was decreased and predominantly restricted to the stroma in SR tumors and absent in corCD8^lo tumors (Figure 6, C–E). These results indicate increased cytotoxic activity of CD8⁺ T cells in FI tumors compared with SR tumors, as confirmed by the gene expression metasignatures. Consistent with this, we observed a significant accumulation of proinflammatory CD68⁺CD206^– macrophages in the epithelium of FI tumors when compared with tumor epithelium in the other groups (Figure 7), further indicating the presence of an active immune response in the epithelium of FI tumors. Together, these results identify a distinct antitumorigenic immune microenvironment mediated by CD8⁺ T cell cytotoxicity in FI compared with SR TNBC.

Figure 7

FI TNBC tumors are infiltrated with proinflammatory macrophages. (A) Representative images of CD68⁺CD206^– and CD68⁺CD206⁺ macrophages show proinflammatory CD68⁺CD206^– accumulation in FI tumors (n = 15). White squares indicate the position of the zoomed-in stromal regions for MR and SR tumors and epithelial region for the FI tumor. Pan-CK (pink) identifies tumor cells and DAPI (blue) identifies cell nuclei. Scale bars: 50 μm (merge, CD68, and CD206) and 5 μm (enlarged insets showing CD68 and CD206 colocalization). (B and C) Quantification of CD68⁺CD206^– (B) and CD68⁺CD206⁺ (C) macrophages in each tissue compartment across groups (n = 15). Green, blue, and red dots represent MR, SR, and FI tumors, respectively. Data were analyzed using the Kruskal-Wallis test and represent the mean ± SEM.

Stromal CD8⁺ T cell–restricted tumors have a distinct TIME. To determine pathways enriched in SR versus FI tumors, we examined metasignatures specifically associated with SR tumors. One of these metasignatures (epiCD8 MSig5, Figure 4B) was dominated by the “superpathway (SPP) of cholesterol biosynthesis” (https://targetexplorer.ingenuity.com/pathway/ING/ING:8h0v2#!/api/rest/v1/client/searchPathwayNodes?pathwayId=ING:8h0v2&rows=0&facetLimit=5000&responseType=default). This is the only pathway that is significantly inversely correlated with epiCD8 in the tumor stroma as well as bulk tumor (Supplemental Table 3). Accordingly, we found that key leading-edge genes involved in cholesterol biosynthesis were elevated in the SR samples compared with the FI samples of bulk tumor (Figure 8, A and B, and Supplemental Figure 5, A and B) as well as in tumor stroma (Supplemental Figure 5, C and D) and, to a lesser extent, in tumor epithelium (Supplemental Figure 5, E and F). Hence, compared with FI tumors, SR tumors were enriched for signatures of cholesterol biosynthesis.

Figure 8

SR TNBC tumors are defined by a cholesterol gene expression signature and a distinct TIME. (A) Heatmap depicting the expression of genes of the SPP of cholesterol biosynthesis in bulk tumor from SR and FI tumors (n = 22). (B) Signature score of genes depicted in (A). Data were analyzed with Spearman’s correlation. (C) Signature score of ISGs repressed by SREBP2 showing decreased expression in SR versus FI tumors (n = 22). Data were analyzed with Spearman’s correlation. (D) Representative IHF images showing the presence of IL-17–producing cells in the tumor stroma of SR and FI tumors (n = 22). Blue, DAPI; pink, pan-CK; yellow, IL-17F. White squares represent the zoomed position in the images. Scale bars: 50 μm (merge) and 20 μm (enlarged insets). (E and F) Density of IL-17–producing cells (E) and neutrophils (F) across SR and FI tumors (n = 10; 5 patients with the lowest and 5 patients with the highest epiCD8, respectively, for SR and FI tumors). Data were analyzed by Mann-Whitney U test and represent the mean ± SEM.

Cholesterol metabolism and the type 1 IFN response have been shown to be inversely coregulated (19, 20). Indeed, type 1 IFN signaling, reflective of the SREBP2-regulated IFN-stimulating genes (ISGs) highly expressed in FI compared with SR tumors (Figure 8C and Supplemental Figure 5, E and H), negatively regulates SREBP2, the transcription factor controlling expression of cholesterol biosynthesis genes (19, 20). Our data consistently showed mutual exclusion between the cholesterol biosynthesis signature identified in SR tumors and a type 1 IFN response in FI tumors (Figure 4B and Supplemental Table 3). These data support the idea that the observed cholesterol signature associated with exclusion of CD8⁺ T cells from the tumor epithelium (SR tumors) is inversely related to a type 1 IFN signature identified in FI tumors.

A second metasignature associated with SR samples (epiCD8 MSig6) contains the “IL-17A in psoriasis” pathway, which is associated with autoimmunity (21) and shown to be protumorigenic and immunosuppressive in cancer (22). This includes genes expressed in IL-17–producing cells and in response to IL-17, including psoriasin (S100A7) (ref. 23 and Supplemental Figure 6). Immunofluorescence analyses revealed that SR tumors with the lowest epiCD8 had higher infiltration of IL-17–producing cells than did FI tumors with the highest epiCD8 profile (Figure 8, D and E, and Supplemental Figure 7A). The majority of the infiltrating IL-17–producing cells were negative for CD4 expression by immunofluorescence (Supplemental Figure 7, B–G). Thus, enrichment of IL-17–producing cells in SR tumors is probably due to the infiltration of γδ T cells or other IL-17–producing cells. IL-17 and γδ T cells are associated with an immunosuppressive microenvironment in part through their ability to recruit neutrophils (22, 24). In support of this, elevated levels of neutrophil-associated genes, including AMICA1 (also known as junction adhesion molecule like [JAML]) (Supplemental Figure 6B) as well as elevated levels of neutrophils were observed in SR tumors displaying strong stromal restriction when compared with levels in FI tumors (Figure 8F, Supplemental Figure 6C, and Supplemental Figure 7H). Together, these results link extreme stromal restriction of CD8⁺ T cells within the SR TIME subtype with the presence of IL-17–producing cells and neutrophils.

Since CD8⁺ T cells recognize antigens presented by MHC class I (MHC-I) molecules (heterodimers composed of HLA class I [HLA-I]and 2M subunits) via interaction with their T cell receptor (TCR), the levels of MHC-I expressed by tumor cells may influence CD8⁺ T cell localization (25). To investigate alternative mechanisms that could determine the spatial pattern of CD8⁺ T cell localization, we examined the level of tumoral HLA-I by IHC and found that a subset of ID, MR, and SR TIME subgroups contained tumors with decreased expression or loss of HLA-I (Supplemental Figure 8). In contrast, we detected no loss of HLA-I in FI tumors (Supplemental Figure 8). The consistent positivity for MHC-I expression in all FI tumors supports a capacity for antigen presentation while the decreased levels of HLA-I in subsets of ID, MR, and SR TIME subtypes, which may contribute to diminished infiltration of CD8⁺ T cells in these TNBCs. While SR tumors are equipped with different potential immune evasion mechanisms (HLA loss, neutrophil infiltration, IL-17–producing cell infiltration), all of these tumors are characterized by a high cholesterol biosynthesis signature consistent with their low IFN signature.

Metasignatures are prognostic in an independent validation cohort. Immune signatures can stratify TNBC patients and predict outcomes (1–3), yet some patients with high expression of immune-based signatures still have a poor outcome (14). Therefore, we aimed to assess the prognostic value of our metasignatures (Figure 9). To do so, we first identified metasignatures that best discriminated TIME subtypes in our discovery cohort (Figure 9) and applied Cohen’s κ statistic. The corCD8 MSig 1 “immune” and corCD8 MSig 3 “fibrosis” signatures were most enriched in and accurately predicted corCD8^hi and the corCD8^lo tumors, respectively (κ = 0.55) (Figure 4A and Figure 9B). Moreover, a combination of corCD8 MSig 1 “immune” and corCD8 MSig 3 “fibrosis” (referred to as “immune/fibrosis”) had a better predictive value (κ = 0.71) than did either metasignature alone (Figure 9B). Similarly, the epiMSig 2 “IFN” and epiMSig 5 “cholesterol” metasignatures were enriched in both LCM compartments and accurately predicted FI and SR TIME subtypes (κ = 0.45 and κ = 0.46), respectively, whereas combining these 2 metasignatures (referred to as “IFN/cholesterol”) demonstrated a better predictive value (κ = 0.52) (Figure 9B).

Figure 9

TIME metasignatures show prognostic value in the external TNBC cohort. (A) Analysis pipeline showing (a) Discovery of metasignatures in a 2-step process (corCD8 and then epiCD8 stratification); (b) prediction in our training set; and (c) validation of the external set. (B) Cohen’s κ statistics measuring the prediction accuracy of each metasignature and combinations. (C) Recurrence-free survival curves using the identified combinations of corCD8 MSig (up; n = 337) and epiCD8 MSig (down; n = 196). For the corCD8 MSig immune^hi fibrosis^lo versus immune^lo fibrosis^hi, the HR is 0.63 (0.456, 0.887), P = 0.006. For the epiCD8 MSig IFN^hi cholesterol^lo versus IFN^lo cholesterol^hi, the HR is 0.52 (0.330, 0.843), P = 0.01. The log-rank P value is shown on the plots. (D) Working model of TNBC stratification into immune subgroups on the basis of the metasignatures identified.

To assess the predictive value of this 2-step approach, we first applied the “immune/fibrosis” metasignature to an independent external data set of chemotherapy-naive TNBC (n = 579) (ref. 3 and Figure 9C) for which recurrence-free survival (RFS) data were available (Supplemental Figure 9). This stratified a poor outcome immune^lo fibrosis^hi (MR-like) (log-rank P = 0.04) and a good outcome immune^hi fibrosis^lo tumor subset. As a second step, the “IFN/cholesterol” metasignature, when applied to the good outcome subset, further stratified this subgroup into patient subsets with intermediate outcomes (IFN^lo cholesterol^hi, SR-like) and good outcomes (IFN^hi cholesterol^lo, FI-like) (log-rank P = 0.04) (Figure 9C). Hence, stratification by sequential use of “immune/fibrosis” and “IFN/cholesterol” metasignatures showed a prognostic value in independent gene expression data sets derived from TNBC bulk tumor samples.

Importantly, when applied using the same 2-step process (Figure 9, A and D), both combination signatures (“immune/fibrosis” and “IFN/cholesterol”) had prognostic value. Thus, the corMSig1/3 “immune/fibrosis” combination could be used as a first step to stratify MR-like versus SR- and FI-like tumors. As a second step, the epiMSigs 2 and 5 “IFN/cholesterol” combination could then be applied to distinguish between SR-like versus FI-like subtypes (Figure 9D). Together, these findings demonstrate that CD8⁺ T cell localization pattern–derived metasignatures, when applied in a 2-step approach, capture distinct aspects of TNBC patient prognosis that cannot be detected using existing approaches.

Specific coinhibitory molecules are associated with distinct immune microenvironments. Tumors can escape immune attack using various immunosuppressive mechanisms. These include the recruitment of immunosuppressive cell types such as Tregs as well as expression of negative regulators such as PD-1 and PD-L1 that can lead to a progressive decrease in T cell effector activity and functional hyporesponsiveness (12, 26, 27). By examining gene expression, we found that multiple regulators, known to act as negative feedback loops following immune activation, were elevated in the epithelial compartment of FI tumors (Figure 10, A and B, and Supplemental Figure 10). These regulators include members of the Ig superfamily of B7 coinhibitory receptors (PDCD1 [PD-1], CTLA, TIGIT) and the PD-1 ligand PD-L1 (CD274), as well as other checkpoint receptors including lymphocyte-activating 3 (LAG3), T cell immunoglobulin and mucin domain 3 (TIM3, encoded by the gene HAVCR2), and the immunoregulatory enzymes indoleamine 2,3-dioxygenase 1 and indoleamine 2,3-dioxygenase 2 (IDO1 and IDO2) (27). In contrast, expression levels of the B7 family members B7-H4 (VTCN1) and B7-H3 (CD276), which can mediate inhibition of T cell activity and infiltration (28, 29), were elevated in corCD8^lo (ID and MR) tumors (Supplemental Figures 10 and 11) and were inversely correlated with CD274 (PD-L1) expression (Figure 10C and Supplemental Figure 10B). Immunostaining revealed an inverse correlation between PD-L1 and B7-H4 protein expression (Figure 11) as well as an association with distinct TIME subgroups. Notably, PD-L1 and other markers of inflammation, such as IDO1, were both enriched in the epithelial compartment of FI tumors and in the stromal compartment of SR tumors, but were low or absent in corCD8^lo (ID and MR) tumors (Figure 11, A and B, and Supplemental Figures 12 and 13). Similarly, CD4⁺ T cells that express the transcription factor FOXP3 and can function as Tregs were elevated in SR and FI tumors compared with that seen in MR and ID tumors and followed the spatial distribution patterns of CD8⁺ T cells, with accumulation in stroma in SR tumors and infiltration into the epithelial compartment in FI tumors (Supplemental Figure 14). Hence, expression and localizatin of PD-L1 and IDO1 as well as precence of FOXP3⁺CD4⁺ T cells were positively correlated with presence of CD8+ T cells. Conversely, expression of B7-H4 was predominantly observed in the epithelial compartment of corCD8^lo (ID and MR) tumors (Figure 11, A and B). Together these results demonstrate mutually exclusive expression of known negative regulators of T cells (such as PD-L1 or IDO1) which are positively correlated with presence of CD8⁺ T cells and the PD-L1 family member, B7-H4 which is negatively correlated with presence of CD8⁺ T cells. Our data thus reveals that expression of immune checkpoint targets is distinct between ID, MR, SR, and FI TNBC tumors. These results highlight the importance of delineating the localization of PD-L1 expression as well as of other immune checkpoints to better understand and target mechanisms of response to immune checkpoint inhibitors and improve their efficacy.

Figure 10

TIME TNBC subtypes express distinct markers of immune suppression. (A) Heatmap depicting the expression of classical immunosuppressive genes in tumor stroma and epithelium for each patient (n = 38). (B) Signature scores for the immunosuppressive gene list from A (n = 38). Data were analyzed using the Kruskal-Wallis test and represent the mean ± SEM. (C) VTCN1 (B7-H4) and CD274 (PD-L1) gene expression was inversely correlated (n = 37). Data were analyzed using Spearman’s correlation.

Figure 11

TIME TNBC subtypes display mutual exclusion and distinct localization of PD-L1 and B7-H4. (A) Representative images of PD-L1 (green), B7-H4 (pink), and DAPI (blue) IHF-stained sections. Scale bars: 10 μm and 20 μm (enlarged insets). (B) Quantification of staining in the tumor epithelial compartment for B7-H4 (pink) and PD-L1 (green) for each patient (n = 35). (C) B7-H4 and PD-L1 quantification after IHF show an inverse correlation in the tumor epithelial compartment. This excludes tumors with less than 1% staining in the epithelial compartment for both markers. Data were analyzed using Spearman’s correlation.

The immune context of TNBC has gained acceptance as an important clinical correlate, raising hopes that modulating immune responses via immunotherapies may constitute an effective therapeutic strategy. However, only 8% to 20% of preselected patients with TNBC benefit from anti–PD-L1 or anti–PD-1 immunotherapy (13), highlighting the need for a better understanding of how the TIME architecture influences outcomes in TNBC and responses to current treatment modalities (30). In this study, we provide a deeper understanding of complex TIMEs. By combining immune cell identification and localization in matched clinical samples with gene expression profiling from matched tumor epithelial and stromal compartments, we identified 4 distinct TIMEs (ID, MR, SR, and FI) associated with disease outcomes in TNBC (summarized in Figure 12). These TIMEs are defined by distinct CD8⁺ T cell localization patterns and harbor distinct GzmB positivity, gene expression metasignatures that predict outcome in independent whole tumor data sets, and expression of distinct patterns of immune checkpoint proteins and immunomodulatory cell types. Each TIME class represents a substantial fraction of TNBC cases, a finding that translates well to patient stratification approaches. Metasignatures associated with each TIME subtype have prognostic value, allow for improved stratification, and support TIME-dependent therapeutic strategies for TNBC.

Figure 12

Schematic of TIME TNBC subtype stratification. Poorly infiltrated tumors (ID and MR) are characterized by signatures of fibrosis, enrichment of fibrotic foci, and expression of the immune checkpoint B7-H4. Tumors that display significant infiltration of CD8+ T cells into the tumor core but that showaccumulation specifically in the tumor stroma, display signatures of cholesterol and infiltration of IL-17–producing cells and neutrophils. These tumors also display stromal expression of the immune checkpoint PD-L1. Tumors of the FI TIME subtype are consistently positive for MHC-I, whereas MHC-1 loss is observed in tumors from other TIME subtypes. Tumors with significant CD8+ T cell infiltration in the tumor epithelium are characterized by type 1 IFN signatures as well as activated CD8⁺ T cell (GzmB⁺) and PD-L1 expression in the tumor epithelial compartment. Scale bar: 100 μm.

Our study demonstrates that tumors with low CD8⁺ T cell expression (corCD8^lo, including MR and ID TIME subtypes) have the poorest prognosis. These tumors are negative for expression of the immune checkpoint PD-L1 and the immune modulator IDO1 and display elevated signatures of fibrosis and fibrotic foci. Similar signatures characterized by desmoplasia and matrix remodeling are associated with low immune content and resistance to chemotherapy in ovarian cancer (31–34). These CD8^lo TNBC tumors display elevated TGF-dependent signatures associated with activated stroma and immunosuppressive signals (35) as well as expression of the B7 family coinhibitory molecule B7-H4. B7-H4 can promote an immunosuppressive environment by negatively regulating T cell effector function and infiltration (28, 36, 37). Although originally characterized by its expression on hematopoietic cells, elevated levels of B7-H4 on tumor cells, as was observed here, correlate with poor clinical outcome across multiple types of solid tumors (38). Consistent with our data, B7-H4 expression is upregulated by TGF-1 signaling in colorectal cancer (39). Notably, for TNBC with poor outcome, expression of B7-H4 was inversely correlated with PD-L1 protein and signatures of inflammation indicating that, in these TNBCs, B7-H4 may actively suppress immune infiltration. A similar inverse correlation between B7-H4 and PD-L1 was described in lung carcinoma (40). This indicates the utility of B7-H4 as a biomarker for “immune-cold” tumors and as a potential target for TNBC and other PD-L1–negative solid tumors that display low infiltration of CD8⁺ T cells, and indicates the potential use of inhibitors of TGF-1 signaling to sensitize immune-cold tumors to PD-1 or PD-L1 immunotherapy (41).

In contrast, the FI TIME TNBC subtype had a proinflammatory microenvironment defined by a type I IFN gene signature, the presence of GzmB⁺CD8⁺ T cells and proinflammatory CD68⁺CD206^– macrophages in the tumor epithelium, as well as good outcomes. This observation supports previous studies showing that elevated immune infiltration (10, 42, 43) and immune signatures (1–3) predict good outcomes and responses to chemotherapies.Although the FI TIME TNBC subtype displays good outcome, it has the highest gene expression levels of several immune checkpoints (e.g. LAG-3, TIM-3, TIGIT, CTLA-4, and PD-L1) as well as positivity for PD-L1 and IDO1 and elevated infiltration of FoxP3⁺CD4⁺ T cells in the tumor epithelial compartment. This likely reflects a negative feedback loop consistent with the proinflammatory nature of these tumors and indicates that patients with this TNBC subtype are good candidates for immunotherapy (12, 44). Such negative feedback may be selected for during the progression from breast ductal carcinoma in situ to invasive ductal carcinoma (45).

Using our compartment-specific analysis, we identified a unique TIME characterized by CD8⁺ T cell accumulation in the stroma (SR TIME). Patients with this TIME had worse overall survival in large independent breast cancer cohorts (10), yet factors contributing to this have been poorly described. SR tumors with the highest accumulation of CD8⁺ T cells in tumor stroma had elevated levels of IL-17–producing CD4^– γδT cells as well as neutrophils when compared with levels in FI TIME tumors. IL-17–producing γδT cells are associated with the recruitment of neutrophils with protumorigenic activity and tumor progression (22, 24). In SR tumors, the stromal, but not the epithelial, compartment displays positivity for PD-L1 and IDO1 as well as FOXP3⁺CD4⁺ T cell infiltration, indicating the development of a potentially immunosuppressive stromal microenvironment. In this context, macrophages in peritumoral stroma can foster immune privilege and disease progression through expression of PD-L1 (46).

In many solid tumors, a clinical response to anti–PD-L1/PD-1 therapy occurs most often in patients with tumors identified as inflamed (12). Clinical trials using PD-1/PD-L1 immune checkpoint inhibitors showed responses in up to 8% to 20% of PD-L1⁺ TNBCs in the advanced disease setting (13). Although FI TIME tumors that harbor PD-L1 expression in the tumor epithelial compartment may be expected to derive the greatest benefit from ICB therapies, our study and others indicate that patients with FI TNBC have the best outcomes with standard-of-care chemotherapies (12). Hence, patients with FI TNBC may not be well represented in most of the anti–PD-1/PD-L1 clinical trial cohorts of advanced disease and may constitute a patient population suitable for immune checkpoint therapies in the neoadjuvant or adjuvant setting.

Selection of patients for anti–PD-1/PD-L1 therapies is based on PD-L1 protein positivity (47), regardless of its spatial localization. However, PD-L1 expression has been shown to be enriched in tumor cells or immune cells in distinct TNBCs (17). Our study highlights the importance of stromal PD-L1 positivity in patients with SR TNBC, raising the possibility that these patients may respond differently to anti–PD-L1 therapies. In addition, although PD-L1 positivity, as determined by IHC, corresponds with enrichment of cell populations associated with clinical benefit, PD-L1 testing alone is not sufficient to accurately predict the response to ICB therapy (48). Hence, further analysis of spatial expression patterns of immune checkpoint biomarkers (including PD-L1 but also B7-H4, IDO1, and others) is critical to assess the predictive value, clinical relevance, and optimal combinations of such biomarkers and improve patient stratification for clinical trials of immune checkpoint therapies.

In our study, identification of the SR TIME subtype, an immune signature subtype with a poorer outcome in TNBC (10), demonstrates the importance of including the spatial pattern of CD8⁺ T cell localization when characterizing TNBC subgroups. Clinically, SR tumors with accumulation of CD8⁺ T cells in the tumor stroma would be scored as “immune-positive” using current pathology guidelines. This includes the “immunoscore” for bulk immune infiltration (49) or current guidelines for TIL assessment, in which sTILs are preferentially scored (9). Hence, using the current guidelines, SR tumors would not be fully distinguished from FI tumors, despite their association with a poorer prognosis (10) and the idea that these subtypes are expected to respond differently to therapies.

Currently, bulk gene expression analyses have failed to integrate spatial information of immune cells and do not accurately predict the TIME subtypes identified here (1–3). Using data from bulk tumors, we found that our patients with the SR TIME subtype had elevated levels of “pan-immune” signatures shown by some studies to be predictive of a good outcome (10, 14). In support of this finding, the Lehmann immunomodulatory subtype, derived using bulk tumor, was unable to differentiate between FI and SR TNBCs. In contrast, our metasignatures, derived from gene expression analysis of matched tumor stroma and tumor epithelium with additional stratification based on CD8⁺ T cell localization, allowed for discrimination between SR and FI tumors. Moreover, our data showed that the inverse correlation of cholesterol biosynthesis and type 1 IFN signatures can distinguish TNBCs with FI and SR TIME subtypes and thus provide an approach for identifying the SR TIME subtype TNBC in bulk tumor retrospective data sets. This is consistent with reciprocal negative regulation of the cholesterol and IFN pathways following viral infection, whereby the cholesterol biosynthesis pathway (low in FI TNBC) is transcriptionally downregulated by type 1 IFN signaling (19, 20, 50), which is characteristic of FI tumors. Hence, the decrease of the cholesterol biosynthesis signature in FI versus SR tumors is consistent with enhanced type 1 IFN signatures in FI tumors over SR tumors (Figure 6, A and B). In support of a negative feedback loop between the IFN and cholesterol biosynthesis pathways, the combination of “IFN/cholesterol” metasignatures (epiMSigs 2 and 5) more accurately predicted SR and FI subtypes in our cohort and had a better prognostic value than did using either metasignature alone in independent bulk tumor data sets. In our study, the prognostic value of these signatures in external gene expression data sets derived from bulk tumor TNBCs demonstrates that the TIME expression subtype and the biological processes underlying differential T cell localization govern TNBC progression and responses to standard-of-care therapy. Considering that patients with SR TNBC have an immunosuppressive TIME with reduced type 1 IFN signaling, they may benefit from emerging vaccine-based approaches and/or IFN-stimulating therapeutic strategies such as stimulator of IFN genes (STING) agonists to enhance T cell trafficking to the tumor site (51) or demethylating agents that activate a viral mimicry and type I IFN signaling (52).

To our knowledge, this is the first study to integrate specific spatial distribution of CD8⁺ T cells within whole sections of TNBC tumors with compartment-derived gene expression profiling of tumor stroma and epithelium to identify distinct TIME subtypes. The TIME subtypes identified here (ID, MR, SR, and FI) depict distinct immune landscapes and potential escape strategies involving differential patterns of immune checkpoint proteins (PD-L1 and B7-H4), immune modulators (IDO1), immunomodulatory cell type infiltration (macrophages, Tregs, neutrophils, and IL-17–producing cells), as well as HLA-I loss. These distinct spatial TIMEs could enable an enhanced stratification of TNBC patients for optimal standard-of-care therapy as well as for differential immunotherapy. The approach developed here sheds light on the limitations of current TNBC stratification, immune infiltration assessments, and use of bulk tumor gene expression data sets, all of which fail to fully integrate the heterogeneity of CD8⁺ T cell spatial distribution and the distinct underlying biologies. We believe that a better understanding of the TIME subtypes identified here will contribute to the understanding of the mechanisms by which TNBCs evade immune surveillance and that the integration of these TIME subtypes into associated clinical studies has the potential to aid in the development of new therapeutics and biomarkers.

Sample collection and selection

Detailed protocols and procedures are available in the Supplemental Methods. Samples were collected from patients undergoing breast surgeries at the MUHC between 1999 and 2012. For the purposes of this study, samples were selected according to the following criteria: therapy naive at the time of surgical excision; clinically documented lack of expression or amplification of estrogen receptor (ER), progesterone receptor (PR), and human epidermal growth factor receptor 2 (HER2); a histological subtype assignment of IDC (not otherwise specified [NOS]). LCM and gene expression profiling methodology has been previously published (15), and details are provided in the Supplemental Methods. Agilent Technologies SurePrint G3 Human GE 8×60K Microarrays (catalog G4851A) were used for gene expression profiling.

Gene expression normalization

We complemented our set of LCM samples with sample-matched bulk tumor gene expression data for 37 of 38 samples extracted from a previously published data set (GEO GSE58644) (14). Loess normalization was used to correct within arrays, and quantile normalization was used to correct between arrays using the R package Limma (version 3.22.7). Mean expression values were used to aggregate probes, and the most variable probe was used to summarize transcripts with multiple probes mapping to it. Gene expression from stroma and epithelium was normalized separately. Raw and normalized microarray data have been deposited in the NCBI’s GEO database (GEO GSE88715, for stromal and epithelial gene expression, and GSE88847, for bulk tumor gene expression).

IHC and immunohistofluorescence

IHC directed against HLA-ABC and neutrophil elastase was performed on a Ventana Benchmark XT (Roche Diagnostic) automated system. Briefly, for Neutrophil elastase, deparaffinization was followed by cell conditioning CC1 40 minutes, primary antibody incubation for 60 minutes, followed with the Optiview DAB detection kit (Ventana). HLA was processed in a similar manner with a cell conditioning CC1 16 minutes. Other IHC procedures were performed manually. Briefly, sections were deparaffinized and conditioned, and antigens were retrieved using proprietary buffers (pH6 or pH9). After blocking, primary antibodies were applied at optimized concentrations overnight at 4°C, followed by 30 minutes of incubation with a secondary HRP antibody and then DAB revelation and counterstaining. For immunohistofluorescence (IHF), samples were processed the same manner until incubation of the primary antibody. Detection was performed with tyramide signal amplification. Slides were counterstained with DAPI. Details on staining protocols, antibodies and quantification procedures are provided in the Supplemental Methods.

Pathological assessment

Scoring of TILs and evaluation of fibrotic focus on H&E-stained sections were performed by 2 trained pathologists, who were blinded to the clinical and experimental data, following the proposed guidelines for TILs in breast cancer (9) and fibrotic foci, respectively (18). Details are provided in the Supplemental Methods.

Grouping of patients

All tumors with a CD8⁺ T cell density in the tumor core of fewer than 100 cells/mm² were assigned to the corCD8^lo group as opposed to the corCD8^hi group. The majority of corCD8^lo tumors showed accumulation of CD8⁺ T cells at the tumor margins (marCD8) (corCD8^hi marCD8^hi) and were named MR tumors. The corCD8^lo tumor group also includes a small group of tumors with no accumulation of CD8⁺ T cells at the tumor margins (marCD8 <200 cells/mm²) and were named ID tumors. Tumors having a significant infiltration of CD8⁺ T cells in the tumor core (corCD8^hi) were then divided depending on the infiltration of CD8⁺ T cells into the epithelium (epiCD8). Tumors with epiCD8 infiltration below the median (of 204.5 cells/mm²) were classified as corCD8^hi epiCD8^lo and named SR tumors. Tumors with epiCD8 infiltration above the median were classified as corCD8^hi epiCD8^hi and named FI tumors.

PAM50 and Lehmann subtyping analyses

For PAM50 subtypes, classification of samples by PAM50 subtypes was performed using the Genefu R package (version 1.16.0). For Lehmann subtypes, TNBC patients were labeled according to “TNBC type” via the web-based tool (http://cbc.mc.vanderbilt.edu/tnbc/) (1).

Pathway analyses, signature score, and development of metasignatures

Pathway analyses. To identify genes and pathways associated with CD8⁺ T cell tumor core density, we first ranked all the genes on the microarray platform on the basis of their correlation with CD8⁺ T cell tumor core density (corCD8) for bulk tumor (n = 37 for the bulk tumor data set) and tumor stromal and epithelial compartments (n = 38). To identify genes and pathways associated with CD8⁺ T cell epithelium density (epiCD8) in tumors with core infiltration (corCD8^hi), we computed the correlation between the gene expression profiles and epiCD8 in the corCD8^hi tumor set (n = 22), and the resulting rankings were subjected to gene set enrichment analysis (GSEA). Spearman’s correlations were performed using R (version 3.1.3), and GSEA was performed using the Piano R package (version 1.12.0). Nominal P values obtained for each pathway were corrected for multiple testing using the FDR approach.

Development of metasignatures. To develop whole tumor pathway–based signatures associated with the variable of interest (namely, corCD8 and epiCD8) in each patient, we performed computed GSEA analysis on our bulk tumor cohort (n = 37) with pathways defined by QIAGEN’s Ingenuity Pathway Analysis (IPA) (www.qiagen.com/ingenuity). The signature score for each pathway was calculated as a signed average of the leading-edge genes (the sign being defined as the sign of the correlation between the expression of the gene and the phenotype of interest) using the sig.score function of the Genefu R package (version 1.16.0). Hierarchical clustering was performed using the correlation measure as the distance, and complete linkage was used to cluster the pathway scores for all patients. To group highly correlated pathway–based signatures, we cut the dendrogram at a height of 0.70 to define metasignatures (MSig) for both corCD8 and epiCD8. These metasignatures are associated with different biological processes, as indicated by their underlying pathways.

Analysis for the enrichment of publicly available TGF signatures. Analysis of the enrichment of publicly available TGF signature pathways correlated or inversely correlated with CD8⁺ T cell core density was carried out using TGF signatures available from the following subsets: chemical and genetic perturbation (C2 CGP) and oncogenic signatures (C6) from the Molecular Signatures Database (MSigDB), version 4.0. (see the Supplemental Methods for additional details). To calculate the TGF signature score, all genes from each signature represented in the corresponding heatmaps were combined. All duplicates were removed. Gene expression values were first normalized by Z score, and signature scores were derived using R, version 3.1.3.

Statistics

Spearman’s correlation analyses were performed for assessment of correlations. Nonparametric Mann–Whitney U tests were applied for comparisons between 2 different groups of patients. Ordinary 1-way ANOVA was performed for comparisons of all 3 groups. The statistics described above were computed using GraphPad Prism, version 6.03. Correlations in the generation of pathway analyses and signature scores were derived using R, version 3.1.3.

Prediction by metasignatures

For stratification by corCD8 MSig and then epiCD8 MSig, we first evaluated the predictive value of each metasignature in our training set (n = 37). To quantify the strength of association between the metasignature and the immune groups in the discovery cohort, we used the κ coefficient statistic to measure agreement. The range of the κ statistic is from –1 to +1, where the value 0 represents the amount of agreement that can be expected from random chance, and +1 represents perfect agreement. A κ value of less than 0 implies that the agreement is worse than expected by chance and indicates disagreement. In other words, a κ value of less than 0 indicates no agreement; 0 to 0.2, slight agreement; 0.21 to 0.40, fair agreement; 0.41 to 0.60, moderate agreement; and 0.61 to 0.80, substantial agreement. We performed a meta-analysis on the most biologically relevant metasignatures that represent the corCD8 and epiCD8 phenotypes.

Survival analysis

The prognostic value of gene signatures was assessed using a log-rank test for Kaplan-Meier survival curves, as implemented in the survcomp R package. To assess the prognostic value of the gene signatures, we used the data set from Rody et al. (GEO GSE31519; n = 578) (53), which is a compendium of normalized TNBC samples. Briefly, this compendium collected data from a single platform (Affymetrix U133A and U133 Plus 2.0 chips) and included only samples that were defined as triple-negative on the basis of the mRNA expression of ER, PgR, and HER2. Overall relapse-free survival was used as the endpoint. We first used the combination of metasignatures (immune and fibrosis) from the corCD8 phenotype to mimic the split observed in the discovery cohort, i.e., we divided the TNBC patients into 2 groups: 60% with immune^hi fibrosis^lo status and 40% with immune^lo fibrosis^hi status. We then assessed the association of the combination metasignatures (IFN and cholesterol) from the epiCD8 phenotype with relapse-free survival within the immune^hi fibrosis^lo tumor subcohort.

Data and code availability

Gene expression and clinical data. Raw and normalized microarray data have been deposited in the GEO database under accession numbers GSE88715 (for stromal and the epithelial compartment gene expression) and GSE88847 (for bulk tumor gene expression).

Analysis code. The code and data links required to reproduce this analysis are publicly available on GitHub (Branch name: master; commit ID: 6e03b4b4ebea477ac0d4ba0ef318fdc4f0dd5cee; https://github.com/bhklab/EpiStromaImmune). All software dependencies are available in the Comprehensive Repository R Archive Network (CRAN) (https://cran.r-project.org/) or on Bioconductor (https://www.bioconductor.org/). A tutorial describing how to run our analysis pipeline to generate the figures and tables is provided on GitHub (github.com/bhklab/EpiStromaImmune). The procedure to set up the software environment and run our analysis pipeline is also provided. This work complies with the guidelines proposed by Sandve et al. (54) in terms of code availability and replicability of results.

Study approval

The human subject study was approved by the IRB of the MUHC, and all patient data and biological samples from patients at the MUHC were obtained with their written informed consent.