Immune gene expression and response to chemotherapy in advanced breast cancer.

BACKGROUND: Transcriptomic profiles have shown promise as predictors of response to neoadjuvant chemotherapy in breast cancer (BC). This study aimed to explore their predictive value in the advanced BC (ABC) setting.
METHODS: In a Phase 3 trial of first-line chemotherapy in ABC, a fine needle aspiration biopsy (FNAB) was obtained at baseline. Intrinsic molecular subtypes and gene modules related to immune response, proliferation, oestrogen receptor (ER) signalling and recurring genetic alterations were analysed for association with objective response to chemotherapy. Gene-set enrichment analysis (GSEA) of responders vs non-responders was performed independently. Lymphocytes were enumerated in FNAB smears and the absolute abundance of immune cell types was calculated using the Microenvironment Cell Populations counter method.
RESULTS: Gene expression data were available for 109 patients. Objective response to chemotherapy was statistically significantly associated with an immune module score (odds ratio (OR)=1.62; 95% confidence interval (CI), 1.03-2.64; P=0.04). Subgroup analysis showed that this association was restricted to patients with ER-positive or luminal tumours (OR=3.54; 95%, 1.43-10.86; P=0.012 and P for interaction=0.04). Gene-set enrichment analysis confirmed that in these subgroups, immune-related gene sets were enriched in responders.
CONCLUSIONS: Immune-related transcriptional signatures may predict response to chemotherapy in ER-positive and luminal ABC.

RCT Entities:   Population Interventions Outcomes

Immunohistochemical staining of tumour tissue for oestrogen receptor (ER) and human epidermal growth factor receptor 2 (HER2) is used to select appropriate candidates for endocrine manipulation and HER2-guided treatment in patients with advanced breast cancer (ABC). However, the development of robust, reproducible markers that predict benefit derived from chemotherapy in this patient population has been more challenging. Clinicopathologic characteristics such as Ki67 staining and chemoresistance assays have been used with variable levels of success (Amadori ; Schrag ). The evolution of the biologic characteristics of the tumour throughout the disease trajectory until the manifestation of clinically overt metastases underscores the value of metastatic lesion biopsies that may accurately capture the temporal heterogeneity of the tumour and thus be a better source of predictive biomarkers (Amir ; Lindstrom ; Kimbung ). Liquid biopsy for enumeration of circulating tumour cells has been shown to be prognostic but not predictive for chemotherapy response in ABC and is not recommended for guiding treatment (Smerage ).

The role of the immune microenvironment as a prognostic and predictive factor in BC has recently been elucidated. An easily accessible marker, tumour-infiltrating lymphocytes (TILs) as assessed by haematoxylin–eosin (H&E) staining, has been extensively studied and shown to harbour significant prognostic power in patients with triple-negative and HER2-positive BC, but not luminal BC, in the neoadjuvant and adjuvant settings (Loi ; Adams ; Denkert ), and in HER2-positive ABC (Luen ). Moreover, TILs predict pathologic complete response after neoadjuvant chemotherapy (Denkert ), but their predictive significance in the metastatic setting is currently unknown (Savas ). In addition, questions remain over the possibility of sampling bias, the inability to distinguish between lymphocyte subpopulations and the interobserver and intraobserver variability, despite guidelines issued by the TILs working group regarding the interpretation of H&E sections (Salgado ).

Gene expression analysis has emerged as a powerful tool that has repeatedly demonstrated the molecular diversity of BC (Perou ; Sotiriou and Pusztai, 2009; Reis-Filho and Pusztai, 2011). Although gene expression signatures have also been extensively used for the prognostic stratification of both early and advanced BC, the prediction of response to chemotherapy at the metastatic setting has been a more complex issue. Factors such as the intratumoural heterogeneity and the ever-increasing number of both genetic and epigenetic aberrations that may affect the chemosensitivity of BC hinder the ability to develop and standardise predictive gene expression signatures (Sotiriou and Pusztai, 2009; Reis-Filho and Pusztai, 2011). Gene expression signatures derived from metastatic lesions have been shown to harbour prognostic information in patients with ABC (Ng ; Tobin ; King ). Here we explore the predictive value for chemotherapy of gene expression signatures from fine needle aspirations obtained in the advanced setting.

Materials and methods

Clinical trial, radiology and metastatic biopsies

TEX was a multicentre, randomised phase III trial (ClinicalTrials.gov identifier NCT01433614) that evaluated the activity of epirubicin and paclitaxel, with or without capecitabine, as first-line treatment for locally advanced inoperable or metastatic BC using tailored doses depending on treatment side effects (see study protocol in the Supplementary Data). A third treatment group with fluorouracil, epirubicin and cyclophosphamide was closed short after the study started. In total, 304 patients were randomised (Hatschek ). Enrolment of patients with HER2-positive disease was initially allowed, but was discontinued after the publication of data regarding the efficacy of trastuzumab in this patient subgroup.

Radiological tumour assessments were performed per protocol after every third 3-week chemotherapy cycle and were evaluated according to the Response Evaluation Criteria In Solid Tumors (RECIST) version 1.0. Best objective response and the size change of target lesions were documented and centrally reviewed for each patient, and were used in this analysis. Progression-free survival (PFS) was defined as the interval from date of randomisation to date of disease progression or death. Time to treatment failure (TTF) was defined as the interval from date of randomisation to date of end of study treatment for any reason (progression, patient’s choice, toxicity or death).

As part of the translational aspect of the TEX trial, 149 patients underwent a biopsy before the start of study treatment from at least one site (36.7% lymph nodes, 22.5% liver, 18.3% skin, 15.8% breast and 6.7% from other sites), by fine needle aspiration biopsy (FNAB) (97.6%) or, in few cases, by core biopsy (2.4%).

The clinical study including the correlative analyses was approved by the ethics committee at Karolinska Institutet, which had jurisdiction for all participating centres and by the Swedish Medical Product Agency. All patients received oral and written information, and consented to participate.

Gene expression profiling and data analysis

RNA was extracted from the obtained biopsies and profiled on Rosetta/Merck Human RSTA Custom Affymetrix 2.0 microarray (Seattle, WA, USA; GEO: GPL10379), as described previously (Tobin ) and are available at the Gene Expression Omnibus (GEO) database under accession number GSE56493.

Assignment of the intrinsic subtype in each tumour according to the PAM50 classification (Parker ) has been described previously (Tobin ). For independent confirmation of the intrinsic subtyping, the Absolute Intrinsic Molecular Subtyping (AIMS) approach was also utilised (Paquet and Hallett, 2015).

Six gene expression signatures were explored for their predictive power in patients with ABC treated with first-line chemotherapy in the TEX trial. These included immune-related (Sota ; Denkert ) and proliferation-related gene signatures (Nielsen ), ER signalling-related gene signatures (Desmedt ) and gene signatures associated with the two most common recurring molecular aberrations in BC (Cancer Genome Atlas Network, 2012), Phosphatidylinositol-4,5-bisphosphate 3-kinase catalytic subunit alpha (PIK3CA) (Loi ) and TP53 mutations (Miller ). Gene module scores were derived as the weighted averages of the expression values of the constituent signature genes, where the weight for each gene is +1 or −1 depending on the direction with the phenotype in the original publication (Supplementary Data). Gene expression data were first collapsed to gene level using a nonspecific filter keeping only the probe sets with highest interquartile range in the case of multiple mappings to the same Entrez Gene ID. Only original probe sets or genes that could be mapped to Entrez Gene IDs were used. The derived gene modules were named as Immune module 1 (derived from Sota ), Immune module 2 (derived from Denkert ), ESR1 module, Proliferation module, PIC3CA module and TP53 module. Additional predictive modelling using all gene-expression data was explored using the R package caret (version 6.0-71; R Foundation for Statistical Computing, Vienna, Austria).

Gene-set enrichment analysis (GSEA) of the Reactome gene sets collection in the Molecular Signatures Database (Broad Institute, version 5.2, Cambridge, MA, USA) in the comparison responders versus non-responders was performed using the GSEA Software (Broad Institute, version 2.2.3) (Subramanian ). Genes were pre-ranked according to moderated t-statistics from a differential expression analysis using the R/Bioconductor package limma (version 3.30.4, R Foundation) (Ritchie ).

Quantitative assessment of immune infiltrate

Lymphocytes and cancer cells were counted using FNAB smears from the same tumour lesions used for gene expression profiling. The smears were stained by standard H&E or Giemsa and were examined by two investigators. To be evaluable, a sample had to contain a minimum of 10 cancer cell clusters of at least 10 cells each. In each cluster, tumour cells and lymphocytes were enumerated and the average percentage of lymphocytes was calculated for each smear.

The absolute abundance of eight immune and two stromal cell populations was estimated using gene expression data with the Microenvironment Cell Populations counter (MCP-counter) method using the R package MCP-counter (version 1.1.0, R Foundation) (Becht ).

Statistical analysis

Objective response rates in two or more groups were compared with Fisher’s exact test. Progression-free survival and TTF outcomes in groups were estimated using Kaplan–Meier curves and compared with the log-rank test. The association between gene expression signatures and objective response was assessed using multivariable logistic regression models with the signature score standardised and as continuous variable. The models included age, recurrence-free interval and treatment group as adjustment variables. Test for interaction was performed by a χ2-test of the reduction in deviance when adding an interaction term to the logistic regression model as compared with a model without. The relation between a gene expression signature score and percentage decrease in target lesions following study treatment was analysed with Spearman’s rank correlation coefficient and exact P-value. Cohen’s κ-coefficient was used to measure agreement between molecular subtype assignments. The performance of the additional exploratory predictive models was assessed by the area under the receiver operating characteristic curve and repeated cross-validation analysis. Multiple testing in the GSEA was controlled by estimating the false discovery rate according to Benjamini and Hochberg. An arbitrary level of 5% statistical significance (two-tailed) was used. All data analysis was done in R/Bioconductor (version 3.3.2, R Foundation), unless otherwise specified.

Results

The characteristics of the patients included in the translational part of the TEX trial have been described in detail previously (Tobin ), as has the radiological assessment according to RECIST (Hatschek ; Suzuki ). For the scope of this analysis, additional quality controls regarding tumour cell purity and cellularity were performed. In total, 109 patients were deemed eligible for this study (Supplementary Figure 1). The median age of the patients at the time of inclusion was 55 years (interquartile range 47–61 years) and the median recurrence-free interval was 2.9 years (interquartile range 0.3–4.6 years). At the time of initial diagnosis, 64 patients had ER-positive and 38 ER-negative disease, whereas in 7 cases the ER status was unknown and was only determined by a biopsy at the time of relapse. The HER2 amplification status was negative in 58, positive in 5 and unknown in 46 cases. Molecular tumour subtypes were determined on the study specific FNAB by PAM50, which classified 11 (10.0%) of the cases as luminal A, 34 (31.2%) as luminal B, 28 (25.7%) as basal-like and 32 (29.3%) as HER2 enriched. The remaining 4 (3.7%) were classified as normal breast-like. Subtyping using AIMS led to similar classification for the majority of tumours (Cohen’s κ=0.644; Supplementary Table 1).

Objective responses to chemotherapy among clinical and molecular subtypes

Overall, objective response to therapy was seen in 56.0% of patients evaluated and did not differ between the treatment groups (47.9% in patients treated with epirubicin and paclitaxel, and 64.6% in those receiving the triplet combination, P=0.149). Thus, the type of chemotherapy was merely considered as an adjustment variable in the analysis. Table 1 shows the best objective response in cases with ER-positive and ER-negative primary tumours, as well as among molecular subtypes of the ABC biopsies defined by PAM50. The response rate did not differ according to ER status; however, PFS was statistically and clinically significantly shorter in patients with ER-negative and non-luminal tumours, median PFS time 14.1 vs 7.6 months for ER-positive and ER-negative tumours, respectively (P=0.017; Supplementary Figure 2A), and 16.5 versus 7.8 months for luminal and non-luminal tumours, respectively (P<0.001; Supplementary Figure 2B).

Table 1

Objective response rates to study chemotherapy in the 109 patients of the translational TEX trial, according to clinical and molecular subtypea

	Best objective response, no. (%)
	All responses (CR+PR)	CR	PR	SD	PD	ND
Oestrogen receptor statusb
Positive	34 (53.1%)	3 (4.7%)	31 (48.4%)	21 (32.8%)	9 (14.1%)	5
Negative	22 (61.1%)	6 (16.7%)	16 (44.4%)	10 (27.8%)	4 (11.1%)	4
Intrinsic molecular subtypec
Luminal A	4 (40.0%)	0 (0.0%)	4 (40.0%)	4 (40.0%)	2 (20.0%)	1
Luminal B	21 (65.6%)	3 (9.4%)	18 (56.2%)	9 (28.1%)	2 (6.2%)	2
HER2-enriched	14 (48.3%)	0 (0.0%)	14 (48.3%)	10 (34.5%)	5 (17.2%)	3
Basal-like	15 (60.0%)	6 (24.0%)	9 (36.0%)	7 (28.0%)	3 (12.0%)	3
Normal breast-like	2 (50.0%)	0 (0.0%)	2 (50.0%)	1 (25.0%)	1 (25.0%)	0

Abbreviations: ABC=advanced breast cancer; CR=complete response; HER2=human epidermal growth factor receptor 2; ND=not determined; PD=progressive disease; PR=partial response; SD=stable disease.

Objective response assessed centrally by Response Evaluation Criteria In Solid Tumors v. 1.0.

Primary tumour, except seven cases where it was only determined in a relapse biopsy.

ABC biopsy, by PAM50.

Gene expression modules in relation to objective response

The six gene modules related to immune response, proliferation, ER signalling, PIK3CA and TP53 mutations were assessed for association with the objective response to chemotherapy, adjusted for age, treatment group and recurrence-free interval (Figure 1). The probability of achieving an objective response to chemotherapy was statistically significantly associated with higher immune module scores for immune module 2 (per s.d. odds ratio (OR)=1.62; 95% CI, 1.03–2.64; P=0.04) but not for immune module 1 (per s.d. OR=1.50; 95% CI, 0.96–2.41; P=0.08). Furthermore, higher immune module scores in both modules were statistically significantly associated with objective response rates in both ER-positive (per s.d. OR=2.05; 95% CI, 1.11–4.13; P=0.02 for Module 1 and per s.d. OR=2.23; 95% CI, 1.21–4.48; P=0.01 for Module 2) and Luminal BC (per s.d. OR=6.91; 95% CI, 2.11–35.05; P=0.006 for Module 1 and per s.d. OR=3.54; 95%, CI 1.43–10.86; P=0.01 for Module 2) but not in ER-negative and non-Luminal BC (Figure 1). When the two immune module scores were assessed as dichotomous variables with the median as the cutoff, the results were essentially unchanged (Supplementary Figure 3).

Figure 1

Odds ratio (responder/non-responder) and 95% confidence interval. Odds ratio (responder/non-responder) and 95% confidence interval per increase in module score of 1 s.d., and adjusted for age and recurrence-free interval, and with stratification by treatment arm, for the whole cohort and biological subgroups. (A) Immune module 1; (B) Immune module 2; (C) Proliferation module; (D) ESR1 module; (E) PI3CA module; and (F) TP53 module. CI=confidence interval; ER=oestrogen receptor.

Exploratory GSEA and predictive modelling

To further explore the biology of chemo-sensitivity in ABC, GSEA was performed to identify sets of genes with differential expression between responders and non-responders. As expected, the highly ranked gene sets were different for ER-positive and ER-negative tumours (Figure 2). In the ER-positive subgroup, immune-related gene sets including the ‘interferon-γ signalling’ and ‘PD-1 signalling’ gene sets were enriched in responders, corroborating the association seen with the published immune modules. Notably, the enrichment of immune-related gene sets was more pronounced in cases where the FNAB was not obtained from a lymph node metastasis (Supplementary Figure 4).

Figure 2

Enrichment of Reactome gene sets in the comparison of responders versus non-responders to chemotherapy. Gene sets (rows) are ordered by the normalised enrichment score (NES) in the oestrogen receptor (ER)-positive (A) or luminal by PAM50 (B) subgroup. Long gene-set names are truncated. In both the ER-positive and luminal subgroups, the immune-related gene sets are enriched in responders. Transparency of the colour of a cell for a gene set indicates below or above a threshold for the false discovery rate (FDR).

Furthermore, we aimed to explore whether a novel gene signature derived from the differentially expressed (nominal P<0.05) genes of the two immune modules would outperform either of them in terms of prediction of response, which was not the case (Supplementary Figure 5A). Advanced predictive modelling starting with all gene-expression data and with stringent evaluation by repeated cross-validation did not either outperform the two immune modules (Supplementary Figure 5B).

Association of immune modules with TTF and with changes in tumour size after study treatment

As the best objective response by RECIST only provides a rough (dichotomous) estimate of chemosensitivity, the correlation of the immune module scores with the continuous variables of TTF and the (%) decrease in target lesions at 4 months (after six 3-weekly cycles of chemotherapy) was also investigated.

For the whole study cohort, TTF did not differ significantly between high and low immune module scores (median TTF 4.3 vs 4.2 months, hazard ratio (HR)=1.12, 95% CI, 0.90–1.39, P=0.32 for immune module 1 and median TTF 4.5 vs 4.1 months, HR=1.04, 95% CI, 0.84–1.30, P=0.71 for immune module 2). However, higher immune module 1 scores were associated with a shorter TTF in ER-negative (median TTF 3.8 vs 5.3 months, HR=1.7, 95% CI, 1.15–2.53, P=0.008) and non-luminal (median TTF 4.0 vs 5.0 months, HR=1.34, 95% CI, 1.01–1.78, P=0.046) disease, whereas there were trends for a positive association in both ER-positive (median TTF 4.7 vs 3.8 months, HR=0.9, 95% CI, 0.68–1.19, P=0.46) and Luminal tumours (median TTF 5.5 vs 3.6 months, HR=0.74, 95% CI, 0.0.49–1.10, P=0.14), further corroborating the relationship between the immune modules and chemosensitivity in these subgroups (Figure 3). Similar results were seen for immune module 2 (data not shown).

Figure 3

Kaplan–Meier plots of time to treatment failure. Kaplan–Meier plots of time to treatment failure for high vs low immune module 1 scores (using median as the cutoff) in patient groups defined by (A) ER status of the primary tumour or (B) intrinsic subtype of the metastasis/relapse. ER=oestrogen receptor.

Furthermore, immune module 1 was correlated with decrease in tumour size in the ER-positive (Spearman’s ρ=0.31; P=0.03) or luminal (ρ=0.53; P=0.001) subgroup of patients (Figure 4A and C, respectively). Similar correlation was seen between immune module 2 and tumour decrease in the luminal subgroup; however, in the ER-positive tumours the association was weak and non-statistically significant (Figure 4B and D).

Figure 4

Correlation of the Immune module 1 and Immune module 2. Correlation of the Immune module 1 (A and C) and Immune module 2 (B and D) at baseline with the percentage decrease in target lesions at 4 months. A local polynomial regression (loess) smoother with 95% confidence band is superimposed in each scatterplot. ER=oestrogen receptor.

Quantification of lymphocytes in FNAB smears and by in silico analysis

Fifty of the 109 tumours had FNAB smears that were evaluable for lymphocyte counts. In general, lymphocytes were scarce with a median of 1.3% of the total cells (interquartile range 0.6–2.3%). The correlation between lymphocyte counts and the immune gene modules was low (Spearman’s ρ=0.26 and 0.18 between lymphocytes and Immune module 1 and 2, respectively, Supplementary Figure 6). There was no statistically significant association between the lymphocyte counts and response, whereas the immune signatures retained their predictive value also in this subset of patients (Supplementary Figure 7).

Using the MCP-counter method, the absolute abundancies of immune and stromal cell populations were calculated in silico using the gene expression profiles of the 100 patients that were evaluable for response. As expected, immune cell lineages were correlated to each other and to the immune gene modules (Figure 5). In addition, lymphocytes but also monocytes were correlated with response to chemotherapy in the ER-positive and luminal subgroup, however without outperforming the predictive value of the two immune modules (Figure 5).

Figure 5

Correlation of estimated absolute abundances. Correlation of estimated absolute abundances of immune and stromal cell populations (MCP-counter), immune modules 1 and 2, single-sample GSEA scores for the Reactome interferon-gamma and PD-1 signalling gene sets, and objective response. Colourgram of the Spearman’s rank correlation matrix for the ER-positive (A), ER-negative (B), Luminal (C) and Non-luminal (D) subgroups. The variables are arranged based on clustering using the one minus absolute correlation coefficient as distance metric and the average-linkage hierarchical clustering algorithm, and ordered according to absolute correlation with the objective response, while maintaining the constrains of the dendrogram from the hierarchical clustering. Objective response is encoded as 0 for non-responders and 1 for responders; hence, the Spearman’s rank correlation is equivalent to the rank-biserial correlation. ER=oestrogen receptor; GSEA=gene-set enrichment analysis.

Discussion

In this correlative analysis of a phase 3 randomised clinical study of chemotherapy in the first-line setting, the predictive power of gene expression profiling of ABC biopsies was comprehensively evaluated. To our knowledge, this is the first study that demonstrates the predictive value of immune-related gene signatures in patients with ER-positive ABC treated with chemotherapy, similar to previous observations made in early disease (Ignatiadis ). Our findings suggest that the immune-based gene modules derived by the previous work of Sota and Denkert statistically significantly and consistently outperformed other gene signatures in predicting sensitivity to chemotherapy in ER-positive and luminal BC. In contrast, they did not predict response to chemotherapy in ER-negative and non-luminal tumours as shown in early BC, possibly due to low number of ER-negative tumours in this study, or presumably due to the underrepresentation of stroma in FNAB, which has been shown to exert a strong immunogenic effect in TNBC and potentially drives the association between chemosensitivity and immune function in this subgroup (Bonsang-Kitzis ).

It should be noted, that the two immune modules consisting of 19 and 12 genes, respectively, only have a minor overlap of 3 genes. However, in an exploratory analysis, a combination of the two as well as data-derived novel signatures did not outperform any of the original modules, indicating that the predictive information is dependent on the presence of an activated immune microenvironment and not on the expression of specific genes. The importance of the immune infiltrate in this context was confirmed by an independent and unbiased GSEA, showing that gene sets related to cross-talk between immune cells and cancer cells, PD-1 signalling and interferon-γ signalling were among the most highly enriched in chemo-sensitive tumours. Furthermore, the findings were accentuated when FNAB from lymph nodes were excluded, in which the presence of non-activated lymphocytes may have masked the signals seen in the remaining cases.

The role of TILs in early BC is well established (Savas ); however, using a single marker is a relatively rudimentary approach that does not fully grasp the complexity of the tumour–host interactions. An effort was made to quantify lymphocytes in relapse and metastatic sites using FNAB smears; however, the number of cells was too low in order to have any clinical usefulness. It is unclear whether this is a result of a lower number of lymphocytes infiltrating metastatic lesions per se, as has been previously shown (Ogiya ; Luen ), or whether FNAB is not an appropriate method for this analysis, as it cannot distinguish between intratumoural and stromal TILs.

The study has some limitations that should be acknowledged. It remains unclear whether the predictive role of the immune infiltrate and its interactions with cancer cells in ER-positive ABC can be generalised to other chemotherapeutic drugs, later treatment lines, the early BC setting or in tissue obtained with core instead of FNA biopsies. In addition, objective response by RECIST and TTF but not PFS was used as efficacy outcome in this analysis. Progression-free survival could not be utilised, as almost half the patients enrolled at the TEX trial either discontinued treatment due to toxicity or switched to other treatments, mainly endocrine, after a period of disease stabilisation. Nevertheless, rate of objective response is a widely accepted measure of the activity of a regimen and it has also been correlated with overall survival in patients with ABC treated with chemotherapy (Bruzzi ), while TTF also exhibited trends for the same associations in ER-positive and Luminal tumours. Furthermore, this was a retrospective analysis of prospectively collected samples from a relatively small number of patients, with missing HER2 status in a considerable percentage of those, which could have influenced our results. Prospective evaluation of these findings in a randomised trial is warranted, as no similar cohort in ABC are currently available for external validation. Ideally, for confirming the predictive value of immune activity in this setting, a comparator group receiving other treatment than chemotherapy should be included. Given these limitations, the results should be seen as hypothesis generating.

In conclusion, overexpression of immune-related genes was found to predict chemosensitivity in patients with luminal ABC enrolled in a prospective randomised trial. With the availability of therapies that can effectively modulate immune activity in the tumour microenvironment, these findings may have implications not only in patient selection but also in the design of combination therapies of immunotherapies and chemotherapy, as has been successfully implemented in other malignancies (Langer ) and is being actively pursued in phase 3 trials of ABC (ClinicalTrials.gov identifiers NCT02819518 and NCT02425891).