Abstract
The abundance and functional orientation of tumor-infiltrating lymphocytes in breast cancer is associated with distant metastasis-free survival, yet how this association is influenced by tumor phenotypic heterogeneity is poorly understood. Here, a bioinformatics approach defined tumor biologic attributes that influence this association and delineated tumor subtypes that may differ in their ability to sustain durable antitumor immune responses. A large database of breast tumor expression profiles and associated clinical data was compiled, from which the ability of phenotypic markers to significantly influence the prognostic performance of a classification model that incorporates immune cell–specific gene signatures was ascertained. Markers of cell proliferation and intrinsic molecular subtype reproducibly distinguished two breast cancer subtypes that we refer to as immune benefit-enabled (IBE) and immune benefit-disabled (IBD). The IBE tumors, comprised mostly of highly proliferative tumors of the basal-like, HER2-enriched, and luminal B subtypes, could be stratified by the immune classifier into significantly different prognostic groups, while IBD tumors could not, indicating the potential for productive engagement of metastasis-protective immunity in IBE tumors, but not in IBD tumors. The prognostic stratification in IBE was independent of conventional variables. Gene network analysis predicted the activation of TNFα/IFNγ signaling pathways in IBE tumors and the activation of the transforming growth factor-β pathway in IBD tumors. This prediction supports a model in which breast tumors can be distinguished on the basis of their potential for metastasis-protective immune responsiveness. Whether IBE and IBD represent clinically relevant contexts for evaluating sensitivity to immunotherapeutic agents warrants further investigation. Cancer Immunol Res; 4(7); 600–10. ©2016 AACR.
Introduction
As cancer cells emerge and propagate, they must acquire evasive tactics to escape immune destruction. These tactics must endow malignant cells with the ability to avoid elimination by the immune system early in their formation, establish equilibrium with the immune system as the tumor grows, and eventually escape from immune control, enabling tumor progression (1). However, results from both mechanistic and correlative studies indicate that some established tumors, with or without therapeutic intervention, can eventually succumb to immune control by establishment of a helper T-cell type 1 (Th1)–type immune response that mediates tumor rejection and gives rise to immunologic memory that guards against future metastases (2–4). Consistent with these observations, numerous histologic studies have demonstrated that the abundance and functional orientation of tumor-infiltrating effector cells within primary tumors can significantly predict recurrence-free survival of patients with cancer (5). This is particularly true of cancers of the breast, colon, ovary, lung, liver, head/neck, skin, and brain, whereby abundant tumor-infiltrating lymphocytes (TIL), namely cytotoxic T cells (CTL), memory T cells, natural killer (NK) cells, and NKT cells, have been associated with the suppression of metastatic recurrence (5–7). Together, these findings suggest that (i) abundant tumor infiltration by effector cells is indicative of an antitumor immune response, and (ii) the induction of tumor-reactive immunity capable of guarding against metastatic progression is a fairly common event in multiple cancer types.
Microarray tumor profiling studies are now providing compelling evidence that the relative abundance of tumor-infiltrating immune cells can be quantified by the surrogate measure of intratumoral transcript levels of immune cell–specific genes that are coordinately expressed (8–13). Through hierarchical clustering analysis, we and others have demonstrated that these genes self-organize into multiple gene clusters that embody signaling networks fundamental to distinct and identifiable immune cell subpopulations (12–18). These immune gene signatures or “metagenes,” consistent with the effector cell types they reflect, have been shown by multiple groups to be associated with patient outcomes such as disease-free and overall survival (2, 8, 12, 13, 16, 19–21), response to neoadjuvant chemotherapy (22, 23), and adjuvant molecular therapy (anti-Her2/neu; ref. 24), paralleling historical immunohistochemical observations.
We characterized the underlying biology and prognostic ramifications of several distinct immune metagenes in breast cancer (12). The term “metagene” is defined here as a cluster of coordinately expressed genes that together comprise a single cognate expression vector. By averaging the expression levels of these genes, a composite expression score can be computed for each tumor (25). We showed that the immune metagene scores quantify the relative abundance of distinct effector cell populations, namely cytotoxic T and/or NK cells (the T/NK metagene), antibody-producing plasma B cells (the B/P metagene), and antigen-presenting myeloid/dendritic cells (the M/D metagene). The expression scores of all three metagenes were significantly associated with prolonged distant metastasis-free survival (DMFS) of patients with breast cancer and displayed additive prognostic information when considered in multivariate Cox regression models. We have also observed that these same prognostic immune metagenes are predictive of taxane and anthracycline efficacy in the neoadjuvant setting (22).
Although the density of immune effector cell infiltrates represents a significant prognostic marker, it does not comprehensively characterize the immune system's functional orientation toward cancer. For even in the presence of abundant TIL, some tumors still escape immune destruction and progress to metastatic disease. Why effector immune infiltrates equate with protective antitumor immunity in some cancers but not others could depend on one or more immunomodulatory mechanisms acting alone or in concert. A better understanding of which tumors will and will not give rise to productive immunity could serve as a valuable predictor of immunotherapy responsiveness as well as creating opportunities to discover and model causal mechanisms of pro- and antitumor immunity in specific cancer types.
In breast cancer, we previously observed that the immune metagenes exhibited clear prognostic power in some tumors, but not others. We found that a gene signature highly correlated with tumor cell mitotic index and Ki67 staining (termed the “proliferation,” or “P” metagene) enabled the stratification of breast tumors into subgroups that either permitted or prohibited prognostication by the different immune metagenes (12). Specifically, in tumors with a high proliferative capacity (defined by the upper tertile of the P metagene scores, denoted as PH), high immune metagene scores were equated with significantly prolonged DMFS compared with the poor DMFS associated with low immune metagene scores. By contrast, tumors with intermediate and low proliferative capacities (denoted as the PI and PL tertiles, respectively) could not be stratified by the immune metagenes into significantly different survival groups. Intriguingly, these observations suggest that the successful recognition and elimination of breast cancer by the immune system may be governed, in part, by certain definable tumor characteristics associated with intrinsic tumor immunogenicity.
In this study, we address the question of how to optimally stratify breast tumors into subsets that differ in their potential for protective immune responsiveness. Although a number of immunohistochemical and microarray-based studies have analyzed the protective effects of relative levels of immune infiltration in breast cancer, none have addressed the question of how to distinguish tumors in which elevated immune involvement does or does not confer protective benefit. We hypothesized that the statistical association between the immune metagenes and DMFS of patients could be used as a means to extrapolate the underlying breast tumor phenotypes that differ in their ability to potentiate long-term, immune-mediated tumor rejection. Using an immune metagene model that combines the prognostic attributes of the B/P, T/NK, and M/D metagenes, we define and validate proliferative and intrinsic subtype attributes of breast cancer that indicate the existence of immune benefit-enabled (IBE) and immune benefit-disabled (IBD) tumor subtypes, and shed light on the underlying molecular pathways that may contribute to their differing immunogenic dispositions.
Materials and Methods
Microarray datasets and data processing
All microarray data were obtained from MIAME-compliant (26) data repositories and were associated with IRB-approved specimen collection protocols compliant with REMARK criteria (27). The meta-cohort #1 (MC1) dataset comprised of 1,954 tumor expression profiles of primary invasive breast cancer generated on the Affymetrix U133 series GeneChip microarrays. Queried transcripts included the 22,268 probe sets common to the U133A, U133A2, and U133 PLUS 2.0 array platforms. Assembly of the dataset is detailed in Nagalla and colleagues (12). Briefly, the data represent diverse patient populations spanning 16 medical centers in the United States, Europe, and Asia. Raw data (CEL files) were extracted from the NCBI's Gene Expression Omnibus (GEO; ref. 28; accessions: GSE1456, GSE2034, GSE5327, GSE12093, GSE7390, GSE6532, GSE9195, GSE2603, GSE7378, GSE8193, GSE4922, GSE11121, and GSE45255), the NCI's caArray database (accession: mille-00271), and the EMBL-EBI's ArrayExpress database (accession: E-TABM-158). The meta-cohort #2 (MC2) dataset consisted of 616 breast tumor expression profiles derived from patients in the United States, Asia, Europe, and New Zealand. Samples were analyzed on the same Affymetrix platforms as MC1, and CEL files were accessed via GEO accessions GSE19615, GSE20685, GSE31519, and GSE36771. For MC1 and MC2, all CEL files were preprocessed and normalized using the R software package (29) and library files were provided by the Bioconductor project. Array data were MAS5.0 normalized using the justMAS function in the simpleaffy library from Bioconductor (30) using a trimmed mean target intensity of 600 without background correction. Nonbiologic batch effects were corrected using the COMBAT empirical Bayes method (31) with estrogen receptor status used as a categorical covariate. The Cancer Genome Atlas (TCGA) meta-cohort originally included 532 breast tumors derived from 13 TCGA Tissue Source Sites and profiled on the Agilent 244K G4502A-07 custom microarray platform. TCGA Level-3 lowess-normalized data (17,813 genes) were downloaded from the TCGA Data Portal (https://tcga-data.nci.nih.gov/tcga/). Three metastatic tumors (nonprimary) and 5 male breast tumors were omitted for a total of 524 cases/profiles. Intrinsic molecular subtype calls for the TCGA data were utilized as previously generated (32). For MC1 and MC2 meta-cohorts, molecular subtype calls were assigned according to ref. 12.
Patient and tumor clinical annotations
De-identified clinical data corresponding to tumor expression profiles and including pathology, demographic, and clinical follow-up data were compiled from several sources, including published supplemental tables, microarray data repositories, and through direct inquiry with the original authors. Patient-tumor clinical characteristics are presented in Supplementary Table S1. See Supplementary Methods for further description.
Metagenes, tertiles, and the immune metagene model
The proliferation and immune metagenes were constructed as described in ref. 12, based on the probe set composition as listed in Additional File 6 of that publication. The immune metagenes, consistent with their transcriptional origins in tumor-infiltrating leukocytes (i) significantly correlate with histologic TIL abundance in primary breast tumors (12); (ii) display immune cell–specific expression patterns in FACS-sorted leukocyte populations; and (iii) exhibit negligible and uncorrelated expression in cultured breast cancer cell lines that lack an immune component, as illustrated in Supplementary Fig. S1. Tumors were grouped into metagene tertiles by ranking tumors by metagene scores (gene averages) and identifying the 33rd and 66th percentile thresholds. The rationale behind the selection of the tertile combinations used to define the favorable immunogenic disposition (FID), weak immunogenic disposition (WID), and poor immunogenic disposition (PID) subclasses of the Immune Metagene Model is described by Nagalla and colleagues (12). In the MC2 and TCGA cohorts, metagene scores were computed as described for MC1. See Supplementary Methods and Supplementary Fig. S2 for analysis of comparisons of metagene tertile thresholds across the three tumor cohorts.
Survival analyses
Statistical associations between gene expression levels (metagene scores, probe set signal intensities) and the DMFS (or overall survival) of patients were evaluated by Cox proportional hazards regression or the Kaplan–Meier method using SigmaPlot 11.0 software or R software (survival package). Clinical endpoints and censoring criteria related to DMFS were applied as described in ref. 12. In the multivariable models, covariates were entered as semicontinuous variables. A stepwise covariate selection procedure was used with an entry significance level of 0.05. Likelihood ratio test P values were reported. For analysis of cell type–specific marker genes in IBE and IBD populations, Affymetrix probe sets were processed as described in Supplementary Methods.
Differential gene expression analysis
Microarray probe sets (n = 22,268) were analyzed for differential gene expression between the IBE-FID (n = 179) and IBD-FID (n = 135) tumor subtypes. The two-tailed Student t test (unequal variance) with false discovery correction (Benjamini–Hochberg) was used to assess the significance of gene differential expression. Probe sets having false discovery rates <1% (q < 0.01) and absolute value of log2 fold change >0.50 were selected for downstream analyses. Pathway activation scores were generated as described previously (33), with further details provided in Supplementary Methods.
Results
Defining phenotypic subclasses of immune responsiveness
We performed a pooled analysis of 1,954 breast tumor expression profiles and corresponding recurrence data for MCI (Supplementary Table S1), to derive an optimal classification strategy for distinguishing tumors that are and are not amenable to prognostic stratification by an integrated immune metagene model (IMM). The IMM combines the prognostic attributes of the three immune metagenes (B/P, T/NK, and M/D) into a simplified classification system using unsupervised tertile-based metrics. Specifically, the IMM classifies tumors into one of three outcome-related subclasses according to different combinations of the immune metagene expression-level tertiles (i.e., low, intermediate, and high gene expression; Fig. 1A and B). Based on metastatic recurrence rates in highly proliferative tumors (12), we refer to these subclasses as FID, WID, and PID, as defined in Materials and Methods (favorable, weak, and poor, respectively). The FID subclass is associated with excellent long-term DMFS in highly proliferative tumors and comprises tumors ranked into the upper tertiles of all three immune metagenes, reflective of a histopathologic prominent immune infiltrate (Supplementary Fig. S3A and S3B). In contrast, the PID subclass comprises tumors ranked into the lower tertile of any one or more of the immune metagenes (consistent with an absent or reduced immune infiltrate; Supplementary Fig. S3C and S3D) and is significantly associated with a short interval to distant metastasis. The WID subclass captures tumors of the intermediate/high immune metagene tertiles and corresponds to patient DMFS rates intermediate to those from the FID and PID subclasses.
In addition to tumor proliferative capacity, we previously found that the intrinsic molecular subtypes of breast cancer defined by the PAM50 classifier (34) further influenced, to variable degrees, the prognostic potential of the individual immune metagenes (12). Therefore, we used Cox proportional hazards regression to determine the prognostic relevance of the IMM within each molecular subtype stratified according to proliferative capacity via proliferation metagene tertiles (Supplementary Table S2). The IMM achieved significant associations (P < 0.05, FDR < 10%) with DMFS in four subgroups: intermediate and high proliferating basal-like tumors [designated as Basal (PI, PH)], highly proliferative HER2-enriched tumors [HER2-E (PH)], and highly proliferative luminal B tumors [LumB (PH)]. As illustrated in Fig. 1C, we refer to this group as IBE in reference to the potential of the FID subclass to exhibit a significantly better DMFS rate as compared with the WID and PID subclasses. In contrast, the IMM was not significantly associated with DMFS in luminal A (LumA) and claudin-low (CL) subtypes, regardless of proliferation tertile, nor was an association observed in low and intermediate proliferating HER2-E and LumB tumors, or low proliferative basal tumors. We refer to these tumors as IBD in reference to the inability of the FID subclass to exhibit a significantly better DMFS rate than the WID and PID subclasses (Fig. 1D). Because 56% of IBD tumors belong to the LumA subtype, we considered the possibility that LumA tumors may mask an otherwise significant association of the IMM with DMFS in non-LumA IBD tumors. Thus, we analyzed the FID, WID, and PID survival curves in the non-LumA and LumA tumors separately. By Kaplan–Meier analysis, neither the non-LumA IBD (comprising basal-like, CL, HER2-E, and LumB subtypes) nor the LumA IBD tumors showed a statistically significant stratification of the IMM (Supplementary Fig. S4A and S4B). This indicates that the LumA subtype is not driving the negligible immune-prognostic characteristic of IBD. Together, these observations suggest that prominent immune infiltrate associated with FID tumors can guard against metastatic progression in some tumor subtypes (IBE), but not others (IBD).
To discern the prognostic value of the IMM in the presence of conventional variables, we fit a stepwise Cox proportional hazards model within IBE and IBD to elucidate independent associations with DMFS (Table 1). In IBE tumors, the IMM and lymph node status were each significant univariably, with only the IMM and tumor size remaining significant after adjusting for other variables. In contrast, in IBD tumors, histologic grade, tumor size, lymph node status, age, and ER status were each individually significant, with histologic grade, tumor size, lymph node status and adjuvant treatment remaining significant after adjustment. These findings suggest that in IBE tumors the IMM adds substantial prognostic information to conventional prognostic markers, and that the statistically significant prognostic stratification conferred by the IMM in IBE tumors is not a spurious event driven by other known associations.
. | Univariable . | Multivariable . | ||
---|---|---|---|---|
Variablesa . | HR (95% CI)b . | Pc . | HR (95% CI) . | P . |
IBE | ||||
IMM (PID, WID, FID) | 0.46 (0.38–0.55) | <0.0001 | 0.49 (0.39–0.61) | <0.0001 |
Histologic grade (I, II, III) | 1.00 (0.77–1.30) | 0.99 | — | — |
Tumor size (T1, T2, T3) | 1.20 (0.93–1.54) | 0.16 | 1.37 (1.02–1.84) | 0.04 |
LN status (−, +) | 1.41 (1.05–1.88) | 0.03 | — | — |
Age (≤40 y, >40 y) | 0.95 (0.64–1.43) | 0.82 | — | — |
ER status (+, −) | 1.15 (0.86–1.53) | 0.34 | — | — |
Adjuvant treatment (no, yes) | 0.88 (0.67–1.13) | 0.31 | — | — |
IBD | ||||
IMM (PID, WID, FID) | 1.07 (0.89–1.28) | 0.47 | — | — |
Histologic grade (I, II, III) | 2.23 (1.76–2.84) | <0.0001 | 2.14 (1.65–2.76) | <0.0001 |
Tumor size (T1, T2, T3) | 1.57 (1.24–1.99) | 0.0003 | 1.61 (1.18–2.20) | 0.003 |
LN status (−, +) | 1.62 (1.21–2.16) | 0.002 | 2.50 (1.44–4.31) | 0.001 |
Age (≤40 y, >40 y) | 0.49 (0.33–0.74) | 0.002 | — | — |
ER status (+, −) | 0.49 (0.35–0.68) | <0.0001 | — | — |
Adjuvant treatment (no, yes) | 0.89 (0.68–1.16) | 0.39 | 0.56 (0.33–0.96) | 0.03 |
. | Univariable . | Multivariable . | ||
---|---|---|---|---|
Variablesa . | HR (95% CI)b . | Pc . | HR (95% CI) . | P . |
IBE | ||||
IMM (PID, WID, FID) | 0.46 (0.38–0.55) | <0.0001 | 0.49 (0.39–0.61) | <0.0001 |
Histologic grade (I, II, III) | 1.00 (0.77–1.30) | 0.99 | — | — |
Tumor size (T1, T2, T3) | 1.20 (0.93–1.54) | 0.16 | 1.37 (1.02–1.84) | 0.04 |
LN status (−, +) | 1.41 (1.05–1.88) | 0.03 | — | — |
Age (≤40 y, >40 y) | 0.95 (0.64–1.43) | 0.82 | — | — |
ER status (+, −) | 1.15 (0.86–1.53) | 0.34 | — | — |
Adjuvant treatment (no, yes) | 0.88 (0.67–1.13) | 0.31 | — | — |
IBD | ||||
IMM (PID, WID, FID) | 1.07 (0.89–1.28) | 0.47 | — | — |
Histologic grade (I, II, III) | 2.23 (1.76–2.84) | <0.0001 | 2.14 (1.65–2.76) | <0.0001 |
Tumor size (T1, T2, T3) | 1.57 (1.24–1.99) | 0.0003 | 1.61 (1.18–2.20) | 0.003 |
LN status (−, +) | 1.62 (1.21–2.16) | 0.002 | 2.50 (1.44–4.31) | 0.001 |
Age (≤40 y, >40 y) | 0.49 (0.33–0.74) | 0.002 | — | — |
ER status (+, −) | 0.49 (0.35–0.68) | <0.0001 | — | — |
Adjuvant treatment (no, yes) | 0.89 (0.68–1.16) | 0.39 | 0.56 (0.33–0.96) | 0.03 |
Abbreviations: ER, estrogen receptor; LN, lymph node.
aEntered as semi-continuous variables: IMM (PID, WID, and FID) as 1, 2, or 3; histologic grade (I, II, III) as 1, 2, or 3; tumor size (T1, T2, and T3) as 1, 2, or 3; LN status (−, +) as 0 or 1; age (≤40 years, >40 years) as 0 or 1; ER status (+, −) as 0 or 1; adjuvant treatment (no, yes) as 0 or 1.
b95% confidence interval (CI).
cLikelihood ratio test P value.
IBE and IBD are reproducible disease states
To validate our parameters for defining IBE and IBD, we sought to confirm their existence using microarray datasets from two independent patient populations: MC2 and TCGA. MC2 (Meta-Cohort #2) comprises 616 breast tumors profiled on the Affymetrix U133 series GeneChips, whereas the TCGA cohort comprises 524 breast tumors profiled on an Agilent custom microarray as part of the Cancer Genome Atlas project (32). Using the same procedures used in MC1, we assigned tumors of MC2 and TCGA to metagene tertiles, IMM subclasses and PAM50 subtypes. We then categorized tumors as IBE or IBD according to MC1-defined thresholds and examined the prognostic performance of the IMM subclasses by Kaplan–Meier analysis. In the MC2 cohort, the prognostic performance of the IMM was highly significant in the IBE-designated tumors, stratifying patients into markedly different DMFS risk groups (Fig. 2A). By contrast, the IMM failed to risk stratify patients in the IBD-designated tumors (Fig. 2B). In the TCGA cohort, the IMM again exhibited statistical significance in the tumors classified as IBE (Fig. 2C), but not IBD (Fig. 2D) despite potential limitations including patient overall survival as the primary clinical endpoint (rather than DMFS), >75% of patients lost to follow-up at 5 years, and the different microarray technology platform used in the TCGA study (Agilent versus Affymetrix). Together, these findings confirm the reproducibility of the parameters used to distinguish the IBE and IBD subtypes and demonstrate the robust association of the IMM with both distant recurrence and overall survival in IBE tumors, but not IBD tumors.
The IMM classifies tumors according to a combination of immune-derived gene expression patterns. To better understand how common immunologic mechanisms relate to the IBE and IBD designations, we investigated how individual genes with more clearly defined roles in antitumor immunity associate with clinical outcome in IBE and IBD. CTLs and NK cells are the primary effectors of cell-mediated tumor rejection. Expression of the CD8 transmembrane receptor contributes to synaptic stability during tumor-responsive antigenic activation and is a defining feature of CTLs. Expression of the CD16 (FCGR3A/B) transmembrane receptor is a defining feature of NK cells (though it may also be expressed by neutrophils and macrophages) and contributes to antibody-dependent cell-mediated cytotoxicity (ADCC) through binding of the Fc portion of IgG-type antibodies bound to tumor-associated antigens expressed on cancer cells. Activated CTLs and NK cells alike secrete cytolytic proteins that induce apoptosis in their cancer cell targets. Central proteins of this class include GZMB, GNLY, and PRF1. Finally, professional antigen-presenting cells—including dendritic cells, macrophages, and B cells—present tumor antigens to T cell effectors (e.g., helper T cells) through MHC (HLA) class II molecules (e.g., HLA-DR) and provide costimulatory signals necessary for T-cell activation through expression of B7-1/B7-2 type I membrane proteins (e.g., CD86). DMFS survival curves based on expression tertiles of the genes coding for CD8, CD16, GZMB, HLA-DR, and CD86 are shown in Fig. 3. In all instances, significant positive associations were observed in IBE tumors, but not in corresponding IBD tumors. In contrast, some genes such as GZMB, CD16, and CD86, which showed positive associations with DMFS in IBE tumors consistent with antitumorigenic functions, displayed small but statistically significant negative associations with DMFS in IBD tumors – consistent with possible protumorigenic tendencies in IBD.
Immunomodulatory pathways distinguish IBE and IBD
The FID subclass is defined as the subset of tumors that sort into the upper tertiles of all three immune metagenes, simultaneously. Yet, IBE-FID and IBD-FID tumors reflect different DMFS associations relative to the WID and PID subclasses. This observation is consistent with the biologic hypothesis of a differential metastasis-protective potential between IBE-FID and IBD-FID, with IBE-FID tumors potentiating a metastasis-protective immune phenotype whereas IBD-FID tumors do not. Moreover, patients with IBD-FID tumors exhibit a worse prognosis than those with IBE-FID (Supplementary Fig. S4C and S4D), despite the lower proliferative capacity of IBD tumors in general. Thus, we investigated the transcriptional differences that distinguish IBE-FID from IBD-FID. First, we considered the expression of a predefined panel of genes with known roles in tumor immunology (Table 2). Central to this panel were the genes comprising the immunologic constant of rejection (ICR; ref. 2), a group of genes repeatedly observed to participate in different forms of immune-mediated tissue destruction (2, 35–37). Consistent with the prognostic observations, a number of these genes were found to be statistically significantly overexpressed in IBE-FID tumors, including genes involved in Th1 polarization (IFNG, IRF1, STAT1), cell-mediated cytotoxicity (GZMB, GNLY, PRF1, GZMH), antitumor chemoattraction (CXCL9, CXCL10, CX3CL1, CCL5) and leukocyte adhesion (ICAM1). Positive regulators of CTL activation—including CTLA4, CD80, and CD86—were also overexpressed in IBE-FID. Conversely, we observed that TGFB1, encoding the immunosuppressive cytokine TGFβ, was overexpressed in IBD-FID tumors. Using an unsupervised approach, we applied Gene Ontology (GO) enrichment analysis to the genes significantly differentially expressed between IBE-FID and IBD-FID using the DAVID Bioinformatics Resource (38, 39). Controlling for false discoveries, we identified 460 and 569 probe sets significantly overexpressed in IBE-FID and IBD-FID tumors, respectively. As anticipated, we found that many genes overexpressed in IBE-FID tumors were significantly associated with the GO terms cell-cycle transit, cell division and proliferation, consistent with the IBE subtype comprising mostly of highly proliferative tumors. Specifically, we found that 38% of the IBE-FID–overexpressed probe sets were correlated with the proliferation metagene at R>0.5. After adjusting for these cell-cycle-associated genes and analyzing for GO enrichment, one annotation cluster of significantly enriched terms was related to immunological processes, namely chemokine activity (P = 7.3 × 10−7) and defense response (P = 2.9 × 10−5). In addition to the previously-described ICR genes, these terms comprised a number of additional chemokine ligands (CCL7, CCL8, CCL13, CCL18, CXCL11, CXCL13, XCL1) as well as genes involved in acute phase responses (ORM1, DEFB1, LBP). In contrast, genes overexpressed in the IBD-FID tumors were associated with highly significant enrichment of GO terms such as secreted (P < 0.0001), extracellular matrix (P < 0.0001), skeletal system development (P < 0.0001), and cell adhesion (P < 0.0001). Within these categories were genes representing more distinct biological processes such as angiogenesis (BMP4, CTGF, SLIT2, SLIT3, LEPR, EMCN, FGF18, PLXDC1, EGFL7, CXCL12), proteolysis (ADAM12, ADAM23, MMP28, PLAT, CPE, CPZ, C4A, C7, CFD, HTRA1, AEBP1, PCSK5, TPSAB1), TGF-beta signaling (TGFBR3, BMP4, BMPR1B, THBS4, DCN, CHRD, COMP) and response to wounding (C7, C4A, CFD, CTGF, EPHA3, IGF1, IGF2, IGFB4, LAMB2, LYVE1, MMRN1, THBD, TIMP3, VWF).
. | IBE-FID vs. IBD-FID . | |
---|---|---|
Marker name/functional type . | Mean log differencea . | Pb . |
ICR: Th1 polarization | ||
IFNG | 210354_at | 0.52 | 7.6E−07 |
IRF1 | 202531_at | 0.24 | 7.6E−05 |
STAT1 | 209969_s_at | 0.56 | 4.3E−07 |
IL12B | 207901_at | 0.03 | 8.5E−01 |
TBX21 | 220684_at | 0.13 | 9.6E−02 |
IL12A | 207160_at | 0.09 | 5.2E−01 |
ICR: Cytotoxic mechanisms | ||
GZMB | 210164_at | 0.96 | 1.3E−12 |
GNLY | 205495_s_at | 0.65 | 3.4E−06 |
PRF1 | 214617_at | 0.36 | 4.8E−04 |
GZMH | 210321_at | 0.52 | 4.2E−03 |
GZMA | 205488_at | 0.22 | 8.9E−02 |
ICR: Tissue rejection chemokines | ||
CXCL10 | 204533_at | 0.90 | 7.6E−11 |
CXCL9 | 203915_at | 0.83 | 2.7E−08 |
CX3CL1 | 203687_at | 0.46 | 3.5E−03 |
CCL5 | 1405_i_at | 0.43 | 6.9E−04 |
CCL2 | 216598_s_at | 0.20 | 9.5E−02 |
ICR: Adhesion molecules | ||
ICAM1 | 202637_s_at, 215485_s_at | 0.18 | 1.2E−02 |
VCAM1 | 203868_s_at | 0.14 | 1.4E−01 |
MADCAM1 | 208037_s_at | −0.12 | 2.4E−01 |
CTL activation and regulatory markers | ||
CD8A | 205758_at | 0.15 | 2.8E−01 |
CD8B | 207979_s_at, 215332_s_at | 0.11 | 4.1E−01 |
CD28 | 206545_at | 0.02 | 9.0E−01 |
CD86 | 205686_s_at | 0.25 | 8.3E−04 |
CD80 | 207176_s_at | 0.15 | 2.1E−02 |
CTLA4 | 221331_x_at | 0.29 | 1.7E−02 |
PDCD1 | 207634_at | 0.03 | 8.1E−01 |
FOXP3 | 221333_at | −0.15 | 2.5E−01 |
NK cell activation and regulatory markers | ||
KLRK1 (NKG2D) | 205821_at | 0.02 | 9.0E−01 |
MICA | 205904_at | −0.08 | 5.0E−01 |
MICB | 206247_at | −0.08 | 5.2E−01 |
ULBP1 | 221323_at | 0.05 | 7.9E−01 |
ULBP2 | 221291_at | −0.02 | 9.0E−01 |
Protumorigenic cytokines | ||
TGFB1 |203085_s_at | −0.25 | 1.4E−03 |
IL8| 202859_x_at | 0.36 | 1.0E−01 |
IL10 | 207433_at | −0.13 | 4.2E−01 |
TGFB2 | 220407_s_at | 0.12 | 4.1E−01 |
. | IBE-FID vs. IBD-FID . | |
---|---|---|
Marker name/functional type . | Mean log differencea . | Pb . |
ICR: Th1 polarization | ||
IFNG | 210354_at | 0.52 | 7.6E−07 |
IRF1 | 202531_at | 0.24 | 7.6E−05 |
STAT1 | 209969_s_at | 0.56 | 4.3E−07 |
IL12B | 207901_at | 0.03 | 8.5E−01 |
TBX21 | 220684_at | 0.13 | 9.6E−02 |
IL12A | 207160_at | 0.09 | 5.2E−01 |
ICR: Cytotoxic mechanisms | ||
GZMB | 210164_at | 0.96 | 1.3E−12 |
GNLY | 205495_s_at | 0.65 | 3.4E−06 |
PRF1 | 214617_at | 0.36 | 4.8E−04 |
GZMH | 210321_at | 0.52 | 4.2E−03 |
GZMA | 205488_at | 0.22 | 8.9E−02 |
ICR: Tissue rejection chemokines | ||
CXCL10 | 204533_at | 0.90 | 7.6E−11 |
CXCL9 | 203915_at | 0.83 | 2.7E−08 |
CX3CL1 | 203687_at | 0.46 | 3.5E−03 |
CCL5 | 1405_i_at | 0.43 | 6.9E−04 |
CCL2 | 216598_s_at | 0.20 | 9.5E−02 |
ICR: Adhesion molecules | ||
ICAM1 | 202637_s_at, 215485_s_at | 0.18 | 1.2E−02 |
VCAM1 | 203868_s_at | 0.14 | 1.4E−01 |
MADCAM1 | 208037_s_at | −0.12 | 2.4E−01 |
CTL activation and regulatory markers | ||
CD8A | 205758_at | 0.15 | 2.8E−01 |
CD8B | 207979_s_at, 215332_s_at | 0.11 | 4.1E−01 |
CD28 | 206545_at | 0.02 | 9.0E−01 |
CD86 | 205686_s_at | 0.25 | 8.3E−04 |
CD80 | 207176_s_at | 0.15 | 2.1E−02 |
CTLA4 | 221331_x_at | 0.29 | 1.7E−02 |
PDCD1 | 207634_at | 0.03 | 8.1E−01 |
FOXP3 | 221333_at | −0.15 | 2.5E−01 |
NK cell activation and regulatory markers | ||
KLRK1 (NKG2D) | 205821_at | 0.02 | 9.0E−01 |
MICA | 205904_at | −0.08 | 5.0E−01 |
MICB | 206247_at | −0.08 | 5.2E−01 |
ULBP1 | 221323_at | 0.05 | 7.9E−01 |
ULBP2 | 221291_at | −0.02 | 9.0E−01 |
Protumorigenic cytokines | ||
TGFB1 |203085_s_at | −0.25 | 1.4E−03 |
IL8| 202859_x_at | 0.36 | 1.0E−01 |
IL10 | 207433_at | −0.13 | 4.2E−01 |
TGFB2 | 220407_s_at | 0.12 | 4.1E−01 |
aPositive values reflect log fold increase in IBE–FID relative to IBD–FID.
bAdjusted for false discovery.
Next, we utilized the Ingenuity Pathway Analysis (IPA) Upstream Regulator Analysis tool to identify potential transcriptional regulators that may explain the differential expression of genes observed between IBE-FID and IBD-FID. This tool leverages the Ingenuity Knowledge Base of regulatory relationships among genes, proteins, microRNAs, protein complexes and chemicals to statistically rank potential transcriptional regulators causally involved in activation and inhibition of differentially expressed genes (40). Table 3 shows the top two most significant protein regulator candidates identified as activated in IBE-FID or IBD-FID. In IBE-FID tumors, tumor necrosis factor–α (TNF; TNFα) and interferon-γ (IFNG; IFNγ) were identified as the top two most significant gene regulators with largest activation scores. Biologically, TNFα and IFNγ are pro-inflammatory cytokines known to play central roles in the activation of Th1-type tissue rejection and antitumor immune responses. The many differentially expressed genes that support the activation z-scores of TNFα and IFNγ are shown in the network schematic of Supplementary Fig. S5A. By contrast, among the top most significant regulator to emerge in IBD-FID with greatest activation z-score was transforming growth factor β1. TGFβ is a potent immunosuppressive cytokine that antagonizes Th1-type immune responses. Shown in Supplementary Fig. S5B are the genes contributing to the significant TGF β activation z-score. The identification of estrogen receptor α -1 (ESR1) as a top regulator may be explained by the significant over-representation of ER-positive (luminal) tumors comprising the IBD subtype. To further corroborate the IPA findings of TNFα/IFNγ and TGFβ activity in IBE-FID and IBD-FID, respectively, we investigated the corresponding pathway activation signatures previously defined by Gatza and colleagues (33). Briefly, these signatures were empirically deduced from expression profiles of cell culture experiments and utilize a binary regression algorithm that measures the probability of pathway activation for each sample in a microarray dataset. Consistent with the IPA Upstream Regulator analysis results, we found that the IBE-FID tumors were characterized by significantly higher activation scores for TNFα (P = 9.2 × 10−05) and IFNγ (3.0 × 10−07), whereas IBD-FID tumors displayed significantly higher scores for TGFβ pathway activation (P = 3.5 × 10−06). Together, these findings suggest that major regulators of pro- and anti-Th1 immune responses distinguish the IBE and IBD breast tumor subtypes, which could explain their differential metastasis-protective associations.
Upstream regulator . | Molecule type . | Activation z-score . | P . |
---|---|---|---|
Activated in IBE-FID | |||
TNF | Cytokine | 4.68 | 7.4E−16 |
IFNG | Cytokine | 4.15 | 3.1E−11 |
Activated in IBD-FID | |||
TGFB1 | Growth factor | 3.40 | 1.6E−14 |
ESR1 | Ligand-dependent nuclear receptor | 2.62 | 1.5E−15 |
Upstream regulator . | Molecule type . | Activation z-score . | P . |
---|---|---|---|
Activated in IBE-FID | |||
TNF | Cytokine | 4.68 | 7.4E−16 |
IFNG | Cytokine | 4.15 | 3.1E−11 |
Activated in IBD-FID | |||
TGFB1 | Growth factor | 3.40 | 1.6E−14 |
ESR1 | Ligand-dependent nuclear receptor | 2.62 | 1.5E−15 |
Discussion
The contemporary view of breast cancer as poorly immunogenic is predicated on the absence of an observed increase in breast cancer incidence in immunocompromised patients and lack of responsiveness to otherwise clinically effective immunotherapy (41). Nonetheless, mounting evidence from studies of breast cancer-specific immune responses suggest variable immunogenic activity (42, 43). Recently, the immune contexture of breast cancer has gained broad acceptance as an important clinical correlate in both patient prognosis and therapy prediction. An understanding of how the prognostic value of immune infiltrates varies across the spectrum of breast tumor subtypes is relevant to disease management and can provide much needed biologic insight into the underlying variation in tumor immunogenicity. To this point, we hypothesized that breast tumor immunogenicity, like other biologic properties of cancer, is subject to phenotypic heterogeneity, varying from tumor to tumor according to molecular cues intrinsic to cancer cells or to other components of the tumor microenvironment. Our findings suggest that the immunogenic heterogeneity of human breast cancer may be resolved, to a substantial degree, by the molecular classification of breast tumors into two discernible subtypes that vary on both clinical and biologic scales. We show that these subtypes, IBE and IBD, can be distinguished by phenotype and differ widely with respect to the ability of an immune-dependent model to stratify patients into prognostic subclasses with significantly different DMFS. Our data suggest the IMM subclasses (FID, WID, and PID) reflect a differential ability of the immune system to mount durable antitumor immune responses that are protective of distant recurrence in IBE tumors, but not IBD tumors. Furthermore, we show that the distinguishing biologic features of IBE and IBD are differential measures of proliferative capacity and composition of intrinsic molecular subtype, and as such provide a reproducible classification system for delineating IBE and IBD (Fig. 2).
Our data indicate that IBD tumors are least likely to be controlled by immune infiltrates, perhaps because their functional orientation is less clearly polarized toward a Th1/effector phenotype as supported by evidence in Tables 2 and 3. Expression of CD8 did not significantly differ between IBE and IBD (Table 2), suggesting that the functional molecular orientation is the key determinant of these opposite phenotypes. The results of a gene network analysis suggest that some of the phenotypic characteristics may be related to the production by IBD tumors of potent immune-suppressant factors, such as TGFβ, that may alter the functional status of the immune response in the tumor microenvironment, though this hypothesis remains to be confirmed. Recent evidence suggests that in some breast tumors, TGFβ signaling (in association with Smad3) contributes to the maintenance of an antiproliferative effect manifesting as a reduced tumor proliferative capacity (44). Thus, the potential duality of TGFβ function, as an intrinsic tumor suppressor and potent immunosuppressive cytokine, could explain, in part, the low proliferative phenotype of IBD tumors in general, and an impaired immune response associated with IBD-FID specifically.
In contrast, our findings suggest that IBE tumors are endowed with the potential to achieve a metastasis-protective immune response (in the FID subclass) that may be facilitated by a proclivity toward Th1 polarization (Table 2) via activation of TNFα and/or IFNγ signaling (Table 3). However, as indicated by the reduced DMFS associated with WID and PID tumors, this immunogenic potential may be negatively regulated by the selection of oncogenic mechanisms that induce immune tolerance or block effector cell trafficking into the tumor mass. Indeed, this hypothesis is consistent with recent reports that breast and other solid tumors evade immune destruction and the control of metastatic spread through expression of secreted factors that fortify the tumor endothelial barrier, thereby limiting tumor infiltration by effector cells and impeding immunotherapeutic efficacy in vivo (45, 46).
In reality, the causality for these distinct phenotypes of breast cancer remains elusive as it is for other types of cancer such as colon cancer and melanoma. Although the presence of a Th1 polarization appears to be relevant for improved survival in IBE, it remains unclear why some tumors are endowed with this immune phenotype independent of their ontogeny (47). Genetic determinants of intratumoral immune reactions are only now beginning to be elucidated (48).
Our findings suggest that biologic characteristics of tumor cells, and in particular their proliferative connotation, may be related to the immune phenotype. This possibility may be of particular interest because it suggests that intrinsic characteristics of tumor cells may be driving the recruitment, and shaping the function, of leukocyte populations in the tumor microenvironment. This proliferation component may complement the notion of genetic instability of cancer cells as a potential cause for the production of mutated protein products recognized as foreign by the immune system, leading, therefore, to higher immunogenicity (49). This hypothesis, on the other hand, would be in contrast with the alternative explanation that the genetic background of the host determines tumor immune responsiveness, which in turn would preferentially explain why such immune phenotypes are seen in tumors of different histology, as recently reviewed (5). To date, there remains no valid explanation for this phenomenon, and this study may provide a clue toward dissecting the genetic and epigenetic mechanisms leading to immune responsiveness.
Independent of the causality, however, the current subclassification strategy may provide the practical benefit of allowing a better prognostic stratification of patients if the association with better prognosis will prove to also reflect a phenotype predictive of responsiveness to immunotherapy, a hypothesis that warrants clinical testing. This would not be completely surprising because the expression of ICR genes that was observed in IBE-FID tumors has been shown to span a continuum of immune surveillance that includes immune responsiveness (prediction of response) and prognostic significance (improved survival independent of treatment) and has mechanistic implications (activation of genes during active rejection; ref. 2). These observations lead to the hypothesis that all such phenomena related to immune surveillance are part of the same process and only quantitative variations determine the ultimate outcome (2).
In conclusion, this work suggests a new immunogenic perspective for breast cancer, in which tumor proliferative capacity and intrinsic subtype are organizing principles that distinguish a tumor's capacity for clinically beneficial immunogenicity. Future studies aimed at optimizing the IBE/IBD tumor classification, as well as the prognostication of patients using more comprehensive measures of immune responsiveness, should increasingly improve our contextual understanding of breast tumor immunogenicity and complement ongoing efforts to standardize histopathologic assessment of TIL in breast cancer (50). In addition, such studies should help to bring about a better understanding of the factors that regulate antitumor immunity, with the potential to spark novel therapeutic strategies that target specific immunogenic subtypes.
Disclosure of Potential Conflicts of Interest
No potential conflicts of interest were disclosed.
Authors' Contributions
Conception and design: L.D. Miller, A. Alistar, F.M. Marincola
Development of methodology: L.D. Miller, J.A. Chou
Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): T. Putti
Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): L.D. Miller, J.A. Chou, M.A. Black, C. Print, J. Chifman, X. Zhou, D. Bedognetti, W. Hendrickx, A. Pullikuth, J. Rennhack, E.R. Andrechek, F.M. Marincola
Writing, review, and/or revision of the manuscript: L.D. Miller, J.A. Chou, M.A. Black, A. Alistar, T. Putti, X. Zhou, D. Bedognetti, W. Hendrickx, J. Rennhack, E.R. Andrechek, S. Demaria, E. Wang, F.M. Marincola
Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): J A. Chou, T. Putti
Study supervision: F.M. Marincola
Grant Support
This work was supported in part by the American Cancer Society (RSG-12-198-01-TBG; to L.D. Miller) and the Mary Kirkpatrick Professorship for Breast Cancer Research (to L.D. Miller). The Cancer Genomics Shared Resource and the Biostatistics and Bioinformatics Shared Resource of the Comprehensive Cancer Center of Wake Forest University (supported by NCI CCSG P30CA012197) provided resources for this study. Breast Cancer Cure NZ provided funds for tissue collection and array analysis of the New Zealand cohort.
The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.