Gene Expression Signatures and Immunohistochemical Subtypes Add Prognostic Value to Each Other in Breast Cancer Cohorts

The past 15 years have seen an exponential growth in the number of cancer-specific gene expression signatures aiming to describe tumor biology, estimate patient survival, or predict how likely a malignancy is to respond to specific treatment modalities (1). Unfortunately, few of these signatures have demonstrated validity, utility, and prognostic or predictive capacity in large randomized clinical trials. Breast cancer, however, represents a clear exception to this generalization with a number of signatures warranting discussion and recommendation in treatment guidelines (2–4). Moreover, results from the MINDACT and TAILORx clinical trials have provided level 1A clinical utility evidence supporting the use of the Mammaprint and RS signatures, respectively, to determine which breast cancer patients may be safely spared from chemotherapy (5, 6).

Despite the emergence of gene signatures, targeted breast cancer treatment continues to be administered on the basis of immunohistochemical (IHC) analysis of the estrogen and progesterone receptors (ER and PR), human epidermal growth factor receptor 2 (HER2), and the proliferation marker Ki67 (7). Given the clear importance of both transcriptional and IHC classifiers, studies comparing the prognostic and predictive capacity of both are vital. We have previously demonstrated that the prognostic capacity of Ki67 alone was similar to that of binary gene expression signatures (8) and that Ki67 and gene signatures differ in prognostic classification of tumors into good/poor prognosis subgroups in 18% to 33% of patients (8). This led us to hypothesize that a combination of Ki67 on its own or combined with ER, PR, and HER2 (to form the IHC subtypes) and gene signatures could provide more prognostic information than either classifier alone when considering long-term breast cancer–specific survival (BCSS).

To comprehensively address this supposition, we determined the change in likelihood ratio (LR) when adding gene expression signatures to Ki67 or IHC subtypes in all, ER-positive (ER⁺) lymph node–positive (LN⁺), ER⁺ LN-negative (LN⁻) and ER-negative (ER⁻) subgroups of two breast cancer cohorts with 379 and 209 patients, respectively. Furthermore and uniquely, we performed the reverse analysis, that is, adding Ki67/IHC subtypes to gene signatures, thereby addressing the question of whether gene signatures have the potential to compete well with pathology-based classifications of breast tumors.

Cohort 1 is described in further detail in the Supplementary Methods. Briefly, this cohort consists of 621 individual patients diagnosed with primary breast cancer between January 1, 1997, and December 31, 2005, with available gene expression profiles. Patients were originally selected for inclusion to a nested case–control study design (totaling 768 study subjects) on the basis of development of metastatic disease; however, this design is not being employed here, and as such, it is important to note an enrichment of metastatic events in this cohort relative to a typical breast cancer population. Of these 621 patients, 379 were included in our analysis; reasons for exclusion were bilateral tumor (n = 2), missing ER (n = 13), PR (n = 147), HER2 (n = 96), or Ki67 (n = 55) information or tumors that could not be classified into an IHC subgroup (n = 13); see “IHC subtypes” section below; a CONSORT diagram for this cohort is shown in Fig. 1.

Cohort 2 has been previously extensively described (8, 9). Briefly, this cohort is commonly referred to as the Uppsala cohort and consists of 484 breast cancer patients who received primary therapy in the Uppsala region of Sweden between 1987 and 1989. A CONSORT diagram for this cohort is shown in Fig. 1.

Follow-up information for both cohort 1 (complete to January 10, 2015) and cohort 2 (complete to December 31, 2008) was retrieved from the Stockholm-Gotland Breast Cancer Registry and Swedish National Board of Health and Welfare (Socialstyrelsen), respectively. Median BCSS, defined as patients who have not died from breast cancer in the study period from date of surgery to end of follow-up, ICD 10 code C509) was 12.4 and 12.5 years for cohort 1 and cohort 2, respectively.

ER and PR (assessed as a continuous variable by IHC, >10% cutoff for positivity) were collected from pathology reports, while HER2 and Ki67 were assessed using chromogenic in situ hybridization (CISH) and MIB-1 antibody (1:100 dilution, DAKO) on tissue microarrays (TMA). TMAs were assessed by a pathologist (P.K. Wright), and in the case of Ki67, tumors were split into seven bins based on positive nuclear expression of the protein.

ER and PR were determined using ligand-binding assay where positivity was defined as ≥0.05 and 0.1 fmol/ng DNA, respectively. This assay provides similar prognostic and treatment-predictive information to IHC analysis (10, 11) and the chosen cutoffs are analogous to a >10% IHC cutoff (12). HER2 and Ki67 expression were assessed by a pathologist (J.W. Carlson) on whole tumor sections using HER2/neu (CB11, 1:300, NovoCastra Laboratories Ltd.) and MIB-1 (1:100 dilution, DAKO) antibodies, respectively. Nuclear Ki67 expression was scored at the invasive tumor edge as a continuous variable (13). BRISQ criteria for cohort 1 are shown in Supplementary Methods and have been previously published for cohort 2 (8). This article was performed and is reported in accordance with REMARK guidelines (14). Both gene expression studies were approved by the ethics committee at Karolinska Institutet (Stockholm, Sweden).

IHC subtypes were constructed as follows; Luminal A-like: ER⁺ and PR^+/− (hormone receptor⁺, abbreviated to HR⁺), HER2⁻, Ki67 low, Luminal B-like: HR⁺, HER2⁻, Ki67 high; HER2⁺: any HR, HER2⁺, any Ki67; and triple negative: HR⁻, HER2⁻, any Ki67. Ki67 cutoff was chosen as the median value in cohort 1 (low/high, <16/≥16) and the same cutoff was applied to cohort 2. We saw no survival differences between HR⁺/HER2⁺ and HR⁻/HER2⁺ tumors using Kaplan–Meier analysis and as such grouped HER2⁺ tumors together regardless of HR status (data not shown).

Research versions of the Genomic Grade Index (GGI), 70-gene (commercially Mammaprint), cell-cycle score (CCS), recurrence score (RS, commercially OncotypeDx), and PAM50 signatures were applied to both cohorts as described in the original publications, and we have previously published our R code for these classification calls (8). Array data for both cohorts can be retrieved using NCBI GEO accession numbers GSE48091 and GSE3494. Signatures were chosen on the basis of their relevancy in both a routine clinical setting and in an ongoing Swedish clinical trial (15). The CCS was derived by adding the expression of cell-cycle genes identified from three different databases (KEGG, HGNC, Cyclebase, 463 in total, see Supplementary Table S1; refs. 16–18) and splitting the resulting continuous variable into tertiles of low, intermediate, and high cell-cycle activity. A full description of the methodology for this signature is provided in the Supplementary Materials and Methods.

All statistical analyses were performed using R statistical software version 3.3.3 (19). The LR was used to determine the additional prognostic information gene signatures added to Ki67/IHC and vice versa. The LR can be interpreted as the goodness-of-fit of a model and allows us to compare biomarkers with two (e.g., Ki67) to five (e.g., PAM50) levels. LR and concordance index (c-Index) values were taken from the coxph function of the survival package, whereby each classifier (Ki67/IHC subtypes or gene signatures) was added alone or in combination to a Cox proportional hazards model with BCSS as the clinical endpoint. To test for differences in clinicopathologic variables between ER^+/− patients, χ² and Mann–Whitney tests were used as indicated.

On the basis of previous studies (8) we hypothesized that a combination of Ki67/IHC subtypes with gene signatures could provide more prognostic information than either classifier alone. To test this, we used two breast cancer datasets with matching IHC and gene expression profiles totaling 379 (cohort 1) and 209 (cohort 2) patients, respectively. Exclusion criteria for both cohorts are shown in the CONSORT diagram in Fig. 1, and patient/tumor characteristics split by ER status for both cohorts are displayed in Table 1.

As cohort 1 is derived from a nested case–control study (20) where cases developed distant breast cancer metastasis and controls did not during the same time period, we anticipated an overrepresentation of aggressive tumors. This is most clearly demonstrated by the high number of patients with lymph node metastases (>50%, Table 1, cohort 1) and grade 3 disease (45%, Table 1, cohort 1). Moreover and as expected, ER⁺ tumors were more likely to have received endocrine treatment (Table 1, compare ER⁺ vs. ER⁻ columns, cohort 1). Using the median Ki67 value (16%) for this dataset as a cutoff, we found the rate of prognostic discordance between Ki67 and gene expression signatures to be between 14% and 22% (Supplementary Table S2), in line with our previously published discordance rates (8).

In keeping with our aim to perform a comprehensive analysis in clinically relevant subgroups, we next assessed the additional prognostic information provided by gene signatures when added to Ki67/IHC subtypes in all (n = 379), ER⁺ LN⁻ (n = 104), ER⁺ LN⁺ (n = 167), and ER⁻ patients (n = 103), using the LR test (LR-χ²) as a measure of prognostic power. These groupings represented 100%, 27%, 44%, and 27% of the cohort 1 dataset, respectively, and patient/tumor characteristics for each group are shown in Supplementary Table S3. Kaplan–Meier curves for Ki67, IHC subtypes, and gene signatures are shown for all patients in Supplementary Fig. S1.

In all patients, addition of RS and PAM50 signatures provided statistically significant information beyond Ki67 (Table 2, all patients column, see “G.Sigs addition to Ki67/IHC”, Δ LR-χ² RS = 12.8 and PAM50 = 20.7; P = 0.001 and P < 0.001, respectively) and IHC subtypes (Table 2, all patients, G.Sigs addition to Ki67/IHC, Δ LR-χ² RS = 12.9 and PAM50 = 11.7; P = 0.001 and P = 0.020, respectively). A similar trend for RS and PAM50 was found in ER⁺/LN⁻ patients but did not reach statistical significance.

All signatures added significant prognostic information beyond Ki67 and IHC subtypes in ER⁺/LN⁺ patients, suggesting a prognostic capacity of gene signatures not captured by ER, PR, HER2, or Ki67 in patients with nodal positive breast cancer (Table 2, ER⁺/LN⁺ column, G.Sigs addition to Ki67/IHC). Of note, removal of HER2⁺ tumors from ER⁺/LN⁻ and ER⁺/LN⁺ subgroups did not change these conclusions (data not shown). The prognostic capacity of all gene classifiers was minimal in ER⁻ patients, with no signature reaching statistical significance when added to Ki67 or IHC subtypes (Table 2, ER⁻ patients, see “G.Sigs addition to Ki67/IHC”). This was anticipated given the highly proliferative, poor prognostic nature of ER⁻ tumors.

To determine whether gene signatures can compete well with Ki67/IHC subtypes from a prognostic perspective, we repeated our analysis but this time added Ki67/IHC subtypes on top of gene signatures in the same patient subgroups. In general, Ki67 alone did not add significant prognostic information to any gene signature in any patient subgroup, with the exception of in ER⁺/LN⁻ patients when Ki67 was added to PAM50 (Table 2, ER⁺/LN⁻ patients, see “Ki67/IHC addition to G.Sigs,” LR- Δχ² PAM50/Ki67 = 4.2, P = 0.040). On the other hand, IHC subtypes proved to be more robust and in general added statistically significant prognostic information to gene signatures in all and ER⁻ subgroups, but not in ER⁺/LN⁺ or ER⁺/LN⁻ subgroups (Table 2). One notable exception to this trend in all patients was the reduced prognostic capacity of IHC subtypes when added to PAM50 (Table 2, all patients, Ki67/IHC addition to G.Sigs, LR- Δχ² PAM50/IHC subtypes = 7.1, P = 0.068), highlighting the strength of this gene signature and its ability to potentially provide analogous/superior prognostic information to that found in the pathology-based IHC subtypes. A flow diagram with a simplified summary of all results from cohort 1 is shown in Supplementary Fig. S2, and a more visual representation of all results in Table 2 is shown in Supplementary Figs. S3 and S4.

As a second measure of the prognostic capacity provided by gene signatures, Ki67 and IHC subtypes, we examined the c-index of all classifiers alone or in combination in all subgroups (Supplementary Table S4, cohort 1) and found that results were concordant with LR analysis.

We conducted an identical set of analyses in a second cohort of 209 patients. First, we visually inspected whole-section Ki67 staining in tumors from four representative patients where Ki67 and gene signature classifications were discordant. These samples show instances where (i) Ki67 staining is low (good prognosis) but signature classifications in general point to an aggressive, poor prognosis tumor (Fig. 2A and B) and conversely (ii) where Ki67 is high but signatures designate the tumor as good prognosis (Fig. 2C and D). This indicates that the discordance we see is truly a difference in classification rather than a technical artifact of Ki67 scoring or microarray analysis.

Patient numbers along with clinicopathologic variables and gene signature distributions of the cohort 2 dataset for all (n = 209), ER⁺ LN⁻ (n = 115), ER⁺ LN⁺ (n = 65), and ER⁻ patient (n = 24) subgroups are shown in Supplementary Table S5. These subgroups account for 100%, 55%, 31%, and 11% of the entire cohort, respectively, and Kaplan–Meier curves for Ki67 and gene signatures in this cohort have been previously published (8).

Using the same cutoff for Ki67 (<16/≥16), we found predominantly comparable results for cohort 1; however, reduced patient numbers particularly in the ER⁺/LN⁻ subgroup led to reduced or no statistical significance. The main exception to this was in all patients when gene signatures were added to Ki67 and IHC subtypes and when IHC subtypes were added to signatures (Table 3, all patients), which could be attributed to both a large reduction in patient numbers (379 vs. 209, cohort 1 and cohort 2, respectively) and a weaker performance of IHC subtypes in this dataset (IHC subtypes LR-χ² = 22.9 and 6.2, cohort 1 and cohort 2, respectively). These results are also presented in bar plots in Supplementary Figs. S5 and S6, and c-indices for cohort 2 dataset are shown in Supplementary Table S4, cohort 2. Finally, as the number of patients who received endocrine therapy and chemotherapy differ greatly between cohorts 1 and 2, we reanalyzed all data shown in Tables 2 and 3 adjusting for therapy received in multivariate analysis. No significant changes in results were observed (data not shown).

In this study, we aimed to evaluate whether a range of gene expression signatures add prognostic information beyond that provided by Ki67 or IHC subtypes and then, reversing the analysis, determine whether the same IHC biomarkers added prognostic information to gene signatures. We found that in general, neither binary signatures (GGI and 70-gene) nor the CCS provided additional prognostic capacity in all, ER⁺/LN⁻ or ER⁻ patient subgroups where Ki67 or IHC subtypes have been assessed. This is in line with previous publications from us and others showing the partial reliance of gene expression signatures on proliferation-related genes (21–23). In all patients, both the RS and PAM50 signatures added statistically significant prognostic information to Ki67 and IHC subtypes, suggesting the ability of these signatures to more accurately separate tumors into prognostically relevant groups. All signatures demonstrated prognostic strength beyond Ki67 or IHC subtypes in the subgroup of ER⁺/LN⁺ patients, implying that signatures may contain additional relevant biological information beyond that of proliferation and ER/HER2 signaling. Indeed, the genes of these 70-gene signature encapsulate various aspects of the hallmarks of cancer (24), and the RS genes include those representative of cellular invasion (MMP11, CSL2) and a monocyte/macrophage marker (CD68); however, precisely which genes contribute to the signatures strong performance in this patient subgroup is beyond the scope of the current study.

In the second and highly novel part of our analysis, we examined whether Ki67 or the IHC subtypes added prognostic information to gene signatures. This could be interpreted as an exploratory analysis of whether gene signatures have the capacity to replace IHC profiling with ER, PR, HER2, and Ki67. When considering all patients, IHC subtypes but not Ki67 alone provided significant prognostic information beyond all tested signatures except PAM50, a result analogous to that of Nielsen and colleagues when comparing PAM50 and IHC subtyping in ER⁺ patients (25). Studies comparing both PAM50 and IHC subtyping methods have found clear discordances in classifications (26–28) and notably, all of the IHC subtypes (Luminal A/B-like, HER2⁺, triple-negative) have been shown to be present within each of the PAM50 subtypes (29). This has important clinical implications as among tumors classified as HER2⁺ by IHC, those that are also HER2 enriched by PAM50 benefit the most from anti-HER2 therapy (30, 31). Although our results imply that PAM50 has the potential to provide superior prognostic information to IHC subtyping, it is important to consider that IHC does demonstrate a trend toward significance when added to PAM50 (Table 2, all patients, PAM50 + IHC, LR-χ² for IHC = 7.1, P = 0.068). This implies that with a greater number of samples, IHC may add significant prognostic information beyond PAM50. It is also worth emphasizing that treatment-predictive rather than prognostic information is more relevant in a clinical setting. Additional trials comparing treatment responses for patients with matching versus discordant IHC and PAM50 subtypes will be required to clarify which one, if any, is superior in a treatment-predictive capacity.

Others have also sought to determine signature prognostic capacity in the presence of Ki67 and IHC markers. Niikurra and colleagues and Reyal and colleagues found strong correlations between Ki67 as a continuous variable and the GGI signature (32, 33), and between Ki67 as a binary variable (with a median 20% cutoff) and the GGI G1-like and G3-like subgroups (33). The latter result in particular is important, as it demonstrates similar results to ours while using a different Ki67 cutoff. Conversely, with disease-free survival as clinical endpoint, Bertucci and colleagues found that GGI outperforms Ki67 as a continuous variable in multivariate analysis (34). In ER⁺ HER2⁻ patients, a poor correlation has been demonstrated between Ki67 as a continuous variable and the Endopredict gene signature (a clinically relevant signature not tested here). Similarly, the Endopredict score has been shown to retain statistical significance in multivariate analysis adjusting for ER, PR, and Ki67 (cutoff of 11%; ref. 35). Cuzick and Dowsett and colleagues have examined the additional prognostic capacity of commercial versions of the RS and PAM50 risk of recurrence score beyond IHC4, in the presence of clinicopathologic information (clinical treatment score), using the TransATAC material (36, 37). Direct comparisons between the TransATAC studies and ours are complicated owing to differences in chosen study population and treatment received; however, despite these differences, the prognostic superiority of PAM50 and the ability of IHC4/IHC subtypes to compete well against gene classifiers are common to both analyses. The novelty in our study, however, arises from first, the clear delineation of the prognostic capacity of Ki67 from IHC subtypes through examination of both independently in all analyses, and second the determination of whether IHC can add any prognostic information to gene signatures to understand whether signatures have the potential to compete well (from a prognostic perspective) with routine clinical IHC subtype classifications.

The limitations of our study are as follows: First, our analysis is retrospective in nature; second, we are using the research version of gene expression signatures rather than commercial versions; third, Ki67 was scored at discrete intervals on TMAs in cohort 1, and as a continuous variable on whole tumor sections in cohort 2; fourth, our choice of Ki67 cutoff of 16% in cohort 1 was based on a median value rather than the suggested 20% to 29% cutoff (2), however, in this regard, we have previously shown a 20% cutoff to be a poor indicator of prognosis in our hands (8). Fifth, although results were similar between cohort 1 and cohort 2, we did not see complete agreement in all analyses. We believe that this is most likely owing to lower patient numbers in cohort 2 (particularly in the ER⁺ LN⁺ subgroup) and that an increase in sample size would make results more similar between cohorts. However, to our knowledge, no other dataset of greater size with matching gene expression, ER, PR, HER2, and Ki67 (as a continuous variable) plus full long-term clinical follow-up is currently available in public databases.

In summary, we show that in general, the RS and PAM50 signatures can provide additional prognostic information beyond Ki67 alone or the IHC molecular subgroups. Conversely, the IHC subtypes, but not Ki67, add prognostic information to all gene signatures except PAM50, highlighting the potential of this signature to convey more prognostic information than ER, PR, HER2, and Ki67 combined. Further independent studies are required to confirm these results.

C.M. Perou holds ownership interest (including patents) in and is a consultant/advisory board member for Bioclassifier LLC. No potential conflicts of interest were disclosed by the other authors.

Conception and design: L.S. Lindström, T. Foukakis, N.P. Tobin

Development of methodology: A. Lundberg, C.M. Perou, N.P. Tobin

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): J.W. Carlson, P.K. Wright, T. Foukakis, J. Bergh

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): A. Lundberg, L.S. Lindström, C. Falato, P.K. Wright, C.M. Perou, K. Czene, J. Bergh, N.P. Tobin

Writing, review, and/or revision of the manuscript: A. Lundberg, L.S. Lindström, J.C. Harrell, C. Falato, T. Foukakis, C.M. Perou, K. Czene, J. Bergh, N.P. Tobin

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): N.P. Tobin

Study supervision: N.P. Tobin

This work was supported by BRECT, the Swedish Cancer Society, the Cancer Society in Stockholm, the King Gustaf V Jubilee Foundation, the Swedish Breast Cancer Association (BRO) and the Swedish Research Council (to J. Bergh), and by the Swedish Research Council (grant no: 521-2014-2057; to L.S. Lindström). C.M. Perou and J.C. Harrell were supported by funds from the NCI Breast SPORE program (P50-CA58223-09A1), by R01-CA195754-01, and the Breast Cancer Research Foundation (to C.M. Perou).

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.