Repeated observation of breast tumor subtypes in independent gene expression data sets

Abstract

Gene expression profiling has come into use as a way of defining, at the molecular level, the phenotypes of many kinds of tumors. In the past, we have studied genome-wide expression patterns in several cancers including lymphoma, breast, lung, liver, ovarian cancer and soft tissue sarcomas (1–5). One of the common features of these studies has been the emergence, through hierarchical clustering analysis, of tumor subtypes with distinct gene expression patterns for each of these cancers. The differences in gene expression patterns among these subtypes are likely to reflect basic differences in the cell biology of the tumors. On this basis, one might consider these molecular subtypes as separable diseases.

The molecular differences between the tumor subtypes are often accompanied by differences in clinical features, such as statistically robust differences in relapse-free and overall survival (1, 3, 6, 7). When an alternative approach (i.e., analysis supervised by outcome data) was used, many studies found handfuls of individual genes whose expression is associated with prognosis (8–13). These genes define potential prognostic molecular markers without respect to the biological diversity represented by the subtypes.

Previously we reported a molecular taxonomy of breast cancer based on variation in global gene expression patterns measured by cDNA microarrays (5, 6). In this report, we attempt to reconcile the definition of those subtypes, and the accompanying differences in disease outcome.

Materials and Methods

Tumor samples from two independent studies of response to chemotherapy of locally advanced breast cancer (14, 15) were analyzed on cDNA microarrays as described (6). Details of these and the additional tissue samples analyzed in this study can be found in Table 2, which is published as supporting information on the PNAS web site, www.pnas.org. Altogether, 122 tissue samples were included in the analysis, of which 77 carcinomas and 7 nonmalignant tissues were previously published. Data for all experiments are stored in the Stanford Microarray Database (http://genome-www.stanford.edu/MicroArray/) and can be accessed at (http://genome-www.stanford.edu/breast_cancer/). In addition to these samples, we reanalyzed published data from two independent studies: van't Veer et al. (13) and West et al. (16).

An “intrinsic” gene list was selected (534 genes represented by 552 clones; 500 of these correspond to a single unique UniGene cluster in source (http://source.stanford.edu) (17), consisting of those genes whose expression varied the least in successive samples from the same patient's tumor but which showed the most variation among tumors of different patients. For each data set, as many of these genes as possible were used for clustering (534 for the Norway/Stanford cohorts, 461 for the van't Veer et al. cohort, and 242 for the West et al. cohort).

Centroids (i.e., profiles consisting of the average expression for each of the 500 genes) were computed for each of the five classes found in the Norway/Stanford data. To be conservative, only those tumors that showed the highest correlation with each other within a subtype were used for this calculation (Table 3, which is published as supporting information on the PNAS web site). We then computed the correlation of each sample from the two additional published data sets to each of these five centroids.

Class prediction was performed by using prediction analysis of microarrays (PAM), which is a variant of nearest-centroid classification with an automated gene selection step integrated into the algorithm (18). During cross-validation, a parameter Δ was iteratively increased, and a value that balanced prediction accuracy with a minimal set of genes was chosen for the final model. This value of Δ was used for training on the entire Norway/Stanford set and predicting the class of each sample in the van't Veer et al. and West et al. data sets (Table 4, which is published as supporting information on the PNAS web site). Details of all analyses are discussed in Supporting Materials and Methods, which is published as supporting information on the PNAS web site.

To evaluate the performance of the 231-gene metastasis predictor published by van't Veer et al. (13) on patients from the Norway cohort, we used 77 of the 231 genes (77 genes overlapped with the Norway/Stanford data set) and performed a 10-fold cross validation leave-one-out analysis on the data presented in van't Veer et al. This “training” resulted in 81% accuracy on their own data set, which is similar to the accuracy reported by the authors in a separate validation (13). We then applied this predictor to the Norway cohort (consisting of locally advanced breast tumors from two patient series) by using PAM as was described for the class-predictions above. Univariate Kaplan–Meier analysis was performed by using winstat for excel (R. Fitch Software, Staufen, Germany).

Results

Identification and Validation of Tumor Subtypes in Independent Data Sets. To test the generality of previously proposed subtypes of breast cancer described in Sørlie et al. (6), we analyzed three independent data sets: an extended Norway/Stanford cohort and two published data sets, those of van't Veer et al. (13) and West et al. (16). The genes used for clustering were obtained from an “intrinsic” gene list comprising 534 genes derived from 45 repeated samples in the extended Norway/Stanford cohort. For a comparable analysis of the two other data sets, we used as many of the same genes as possible. Expression centroids were calculated for the core members of each of the five subclasses in the Norway/Stanford data, and correlations with the centroid were calculated for each sample from the two other data sets. Core samples for each class are color-coded similarly in all dendrograms. Table 1 summarizes the distribution of the tumor samples from the different data sets across the subtypes. In an additional, supervised analysis, we trained a predictor on the Norway/Stanford data and used PAM to predict the tumor subtype of a given individual sample in the two other data sets.

Subclassification of Breast Tumors in the Norway/Stanford Data. Our previously published work contained data from 85 tissue samples, 84 of which were reanalyzed in this study. We added 38 more tumor samples from patients with locally advanced breast cancer, which totals 122 breast tissue samples analyzed in this study. The samples fell into five major subgroups, characterized by distinct variation in gene expression pattern, as was previously seen (5, 6) (for the complete cluster diagram, see Fig. 1 and Fig. 6, which is published as supporting information on the PNAS web site). Although the list of genes that was previously used and the gene list that forms the basis for this classification only have ≈200 genes in common, the subtypes are essentially the same and anchored by the same major gene clusters (Fig. 1 C–G). Similarly, with few exceptions, the samples that clustered together in the past still clustered together here (86% within the same groups).

Open in a separate window

Fig. 1.

Hierarchical clustering of 115 tumor tissues and 7 nonmalignant tissues using the “intrinsic” gene set. (A) A scaled-down representation of the entire cluster of 534 genes and 122 tissue samples based on similarities in gene expression. (B) Experimental dendrogram showing the clustering of the tumors into five subgroups. Branches corresponding to tumors with low correlation to any subtype are shown in gray. (C) Gene cluster showing the ERBB2 oncogene and other coexpressed genes. (D) Gene cluster associated with luminal subtype B. (E) Gene cluster associated with the basal subtype. (F) A gene cluster relevant for the normal breast-like group. (G) Cluster of genes including the estrogen receptor (ESR1) highly expressed in luminal subtype A tumors. Scale bar represents fold change for any given gene relative to the median level of expression across all samples. (See also Fig. 6.)

The major distinction seen was between the tumors showing high expression of luminal epithelial specific genes including the ESR1 (Fig. 1G) and all other tumors showing low or no expression of these genes. The basal subtype (red) was the most homogeneous cluster of tumors, as reflected by the relatively short branches linking the tumors in this cluster (node correlation >0.4) and the deep branch separating it from the other subtypes (Fig. 1C).

Breast Cancer Data from van't Veer et al. Gene expression data (log₁₀ ratios) were available for 24,480 genes in a set of 117 tumors from young breast cancer patients (13). Hierarchical clustering was used, exactly as described above for the Norway/Stanford data, to display the expression patterns of 461 “intrinsic” genes in the 97 tumor samples that were obtained from patients diagnosed with sporadic cancer (Fig. 2). Individual dendrogram branches are colored according to the strongest correlation of the corresponding tumor with the subtype centroid as defined for the Norway/Stanford samples.

Open in a separate window

Fig. 2.

Hierarchical clustering of gene expression data from van't Veer et al. (A) The full cluster of 461 genes across 97 sporadic tumors. (B) Experimental dendrogram displaying similarities between the tumors. Branches are color-coded according to the subtype to which the corresponding tumor sample shows the highest correlation. Tumors with low correlation (<0.1) to a specific subtype are indicated by gray branches. (C) Gene cluster associated with the luminal subtype A. (D) Gene cluster containing the ERBB2 oncogene and coclustered genes. (E) Group of genes that tend to be highly expressed in luminal subtype B tumors. (F) Gene cluster characteristic of basal tumors. Scale bar represents fold change for any given gene relative to the median level of expression across all samples. (See also Fig. 7, which is published as supporting information on the PNAS web site.)

As in the Norway/Stanford data, the clearest discrimination was between tumors that expressed genes in the luminal A/ESR1 cluster at high levels (Fig. 2C) and the tumors that were negative for these genes and exhibited expression profiles characteristic of either the basal, the ERBB2+ or the luminal B subtypes (Fig. 2 D–F). It is noteworthy that all samples that showed the strongest correlation with the basal subtype (red branches) are all contained within the left branch of the dendrogram in a tight cluster. The luminal A/luminal B distinction, though less clear than the basal/luminal distinction, is also seen, with many of the luminal B tumors clustering near each other on the right branch of the dendrogram.

BRCA1 Mutations Are Associated with Basal-Type Tumors. The van't Veer et al. study included tumors from 18 carriers of BRCA1 mutations and two carriers of BRCA2 mutations. We did not include this familial subset in the analysis above in order not to risk sample bias in estimating the frequency of different tumor subtypes (Table 1). When we included these 20 tumors along with the 97 samples, we saw little difference in the overall pattern, except for the striking result that all of the tumors from patients carrying BRCA1 mutations fell within the basal subgroup (Fig. 3). This indicates that a mutation in the BRCA1 gene predisposes for the basal tumor subtype, which is associated with lack of expression of the estrogen receptor and poor prognosis. As also reported previously, none of the BRCA1 tumors showed evidence of ERBB2 amplification (20). A distinct expression profile in BRCA1 tumors was also noted by van't Veer et al. as well as by others (21). Two BRCA1 tumors were reported by the authors to be estrogen receptor positive (sample 86 and sample 95); these two nevertheless clustered on the left branch in the dendrogram that contains the basal tumors, but were present on less correlated dendrogram branches. van't Veer et al. also analyzed tumor samples from two BRCA2 carriers, both of which showed luminal A expression patterns.

Open in a separate window

Fig. 3.

BRCA1 tumors associated with a basal tumor profile. (A) Dendrogram showing all tumors from van't Veer et al., including 18 tumors from BRCA1 mutation carriers (black branches) and two tumors from BRCA2 mutation carriers (yellow branches). BRCA1 tumors are indicated with longer arrows; BRCA2 tumors are indicated with shorter arrows. (B) Cluster of genes characteristic of basal tumors and highly expressed in tumors from BRCA1-carriers. (See also Fig. 8, which is published as supporting information on the PNAS web site.)

Breast Cancer Data Set from West et al. Expression levels of a total of 7,129 genes were measured in 49 breast tumor samples by West et al. (16) using Affymetrix oligonucleotide arrays. Data were transformed to a compatible format by normalizing to the median experiment (see supporting information for details). Expression values for 242 “intrinsic” genes were used in a hierarchical cluster analysis exactly as was done for the Norway/Stanford and the van't Veer et al. data sets. Correlation coefficients to the five subtype centroids were calculated for each of the 49 tumors, and the branches of the dendrogram were colored according to the nearest centroid (Fig. 4).

Open in a separate window

Fig. 4.

Hierarchical clustering of gene expression data from West et al. (A) Scaled-down representation of the full cluster of 242 intrinsic genes across 49 breast tumors. (B) Dendrogram displaying the relative organization of the tumor samples. Branches are colored according to which subtype the corresponding tumor showed the strongest correlation with. Gray branches indicate tumors with low correlation (<0.1) to any specific subtype. (C) Luminal epithelial/estrogen receptor gene cluster. (D) Basal gene cluster. (E) ERBB2+ gene cluster. (See also Fig. 9, which is published as supporting information on the PNAS web site.)

Again, the main discrimination seen was between the tumors that highly expressed genes in the luminal A/ESR1 cluster and those that were clearly negative for these genes (Fig. 4 C and D). Data for only about half of the genes from the intrinsic gene list were found in this study, which may explain why a luminal A/luminal B distinction was not seen in this data set. Two tumors showed strong correlation to the ERBB2+ class (Fig. 4E), but they did not cluster together in this analysis.

Prediction of Tumor Subtypes By Using Norway/Stanford Data As a Training Set. Hierarchical clustering analysis is a powerful technique for class discovery; however, we wished also to apply a more supervised analysis that could make a prediction as to what subtype a single sample would belong to when considered only by itself. To accomplish this goal, we used the 115 Norway/Stanford tumor samples and the overlapping “intrinsic” genes for both data sets, respectively, as a training set to develop a breast tumor class predictor using PAM (see Supporting Materials and Methods for details). When we compared these calculations with the results of the hierarchical clustering described above, there was strong agreement, ranging from 79% to 89% for both the van't Veer et al. and West et al. data sets (Table 4). We note that prediction accuracies reported above are somewhat optimistic, as some of the genes used as predictors were used to define the test set groups in the first place.

Tumor Subtypes Are Associated with Significant Difference in Clinical Outcome. In our previous work, the expression-based tumor subtypes were associated with a significant difference in overall survival as well as disease-free survival for the patients suffering from locally advanced breast cancer and belonging to the same treatment protocol (6). To investigate whether these subtypes were also associated with a significant difference in outcome in other patient cohorts, we performed a univariate Kaplan–Meier analysis with time to development of distant metastasis as a variable in the data set comprising the 97 sporadic tumors taken from van't Veer et al. As shown in Fig. 5, the probability of remaining disease-free was significantly different between the subtypes; patients with luminal A type tumors lived considerably longer before they developed metastatic disease, whereas the basal and ERBB2+ groups showed much shorter disease-free time intervals. Although the methodological differences prevent a definitive interpretation, it is notable that the order of severity of clinical outcome associated with the several subtypes is similar in the two dissimilar cohorts. We could not carry out a similar analysis in the West et al. data because the necessary follow-up data were not provided.

Open in a separate window

Fig. 5.

Kaplan–Meier analysis of disease outcome in two patient cohorts. (A) Time to development of distant metastasis in the 97 sporadic cases from van't Veer et al. Patients were stratified according to the subtypes as shown in Fig. 2B. (B) Overall survival for 72 patients with locally advanced breast cancer in the Norway cohort. The normal-like tumor subgroups were omitted from both data sets in this analysis.

Discussion

Breast Tumor Subtypes Represent Distinct Biological Entities. Gene expression studies have made it clear that there is considerable diversity among breast tumors, both biologically and clinically (5, 6, 13, 16, 19, 21, 22). This is not a new idea, as epidemiological studies previously had inferred the existence of two or more subpopulations of breast cancer (23). A straightforward interpretation of the recurrent appearance of several different patterns of gene expression among tumors of similar anatomical origin is to regard each as representing a different biological entity. One possible basis for the consistent differences in these patterns between tumor subtypes might be that they originate from different cell types. Our findings provide some support for this interpretation, as we found breast tumor subtypes with patterns of gene expression similar to those of luminal epithelial cells (the cells that line the duct and give rise to the majority of breast cancers) and patterns of at least one other subtype (termed basal) that resembles the pattern found in basal epithelial cells of the normal mammary gland (characterized by expression of cytokeratins 5/6 and 17).

If indeed luminal and basal tumor subtypes are distinct biological entities, then the cognate expression patterns should be detectable in other genome-scale studies of breast cancer. As shown above, we found strong evidence for the universality of a distinction between basal-like and luminal-like subtypes in two additional, independent data sets comprising different patient populations whose gene expression profiles had been determined by using different microarray technology platforms. We found considerable evidence, in one of the studies, for the distinction between the luminal A and B subtypes. The fact that we could make these distinctions for the basal and luminal subtypes (less so for the luminal B subtype vis a vis luminal A) means that the substantial differences in the characteristics of the patients (e.g., age and tumor stage) are less important determinants of tumor expression phenotypes than intrinsic biology.

The statistical nature of the definition, the differences in the expression technologies, and, more importantly, the limited number of intrinsic genes held in common in particular for the West et al. data set, probably suffices to account for the failure to find coherent clusters for every subtype in each of the cohorts examined.

Another expectation from the concept that the tumor subtypes represent different biological entities is that genetic predispositions to breast cancer might give rise preferentially to certain subtypes. This expectation is amply fulfilled by our finding in the data of van't Veer et al., which shows that the women carrying BRCA1-mutated alleles all had tumors with the basal-like gene expression pattern.

Tumor Subtypes and Clinical Outcome. Consistent with the results previously found in our data (6), we also found differences in clinical outcome associated with the different tumor subtypes in the data set produced by van't Veer et al. The outcomes, as measured here in time to development of distant metastasis, were strikingly similar to what we found previously: worst for basal (and ERBB2+), best for luminal A, and intermediate for luminal B subtypes. Recently, two reports corroborating the poor outcome of the basal subtype solely based on immunohistochemistry with antibodies against keratins 5 and 17 and Skp2, strongly supports our results (24, 25). The finding that our gene cluster profile was of similar prognostic importance in the van't Veer et al. cohort as among our patients is remarkable, taking into account differences regarding disease stage (locally advanced versus stage I primaries) and patient age, but in particular, the fact that the Norwegian patients had presurgical chemotherapy and all patients expressing ESR1 received adjuvant endocrine treatment, whereas the patients from van't Veer et al. in general did not receive any systemic adjuvant treatment.

The observation that BRCA1 mutations are strongly associated with a basal tumor phenotype indicates a particularly poor prognosis for these patients. BRCA1-associated breast cancers are usually highly proliferative and TP53-mutated, and usually lack expression of ESR1 and ERBB2 (20, 26). Status of BRCA1 in familial cancers has failed to be an independent prognostic factor in several studies (reviewed in ref. 27), and is complicated by confounding factors such as frequent screening and early diagnosis.

Molecular Marker Identification. In a mixture of biologically distinct subtypes, it may well be that individual markers derived by supervised analysis will under-perform what is possible if tumor subtypes were separated before searching, in a supervised fashion, for prognostic indicators. Indeed, when we tested the prognostic impact of the 231 markers published by van't Veer et al. on the Norwegian cohort, we found that they performed less well (47%) in predicting recurrences within 5 years (see Materials and Methods). This may in part be due to differences in the patient cohorts and treatments as discussed above.

Both van't Veer et al. and West et al. showed the ability of gene expression profiles to classify tumors into clinically relevant groups and to predict outcome by using supervised statistical analyses. Both reports, however, showed only how the gene expression signatures discriminated tumors based on previously known molecular and clinical parameters, such as ESR1 status, lymph node status, and time to the development of distant metastasis. We have taken a less supervised approach and showed that there are several subtypes of tumors, which may be the result of alterations in different and independent regulatory pathways. The basal subtype was repeatedly recognized as a distinct group in each of three independent data sets, and should be considered as a separate disease with respect to treatment and follow up. The other subtypes are less clear, and require refinement of their molecular definition before they can be reliably defined and diagnosed.

To conclude, classification of breast cancer based on gene expression profiling captures the molecular complexity of tumors. It is for this reason that we believe that the patterns that distinguish subtypes appear to provide a more refined stratification of the patients compared with individual tumor markers. These results imply that the status of the transcriptional programs in the tumor cells and the underlying genetic alterations are major determinants of the tumorigenic potential and ultimately the disease outcome for the patient.

Supplementary Material

Acknowledgments

Notes

Abbreviation: PAM, prediction analysis of microarrays.

References