Abstract
Purpose: The complete molecular basis of the organ-specificity of metastasis is elusive. This study aimed to provide an independent characterization of the transcriptional landscape of breast cancer metastases with the specific objective to identify liver metastasis–selective genes of prognostic importance following primary tumor diagnosis.
Experimental Design: A cohort of 304 women with advanced breast cancer was studied. Associations between the site of recurrence and clinicopathologic features were investigated. Fine-needle aspirates of metastases (n = 91) were subjected to whole-genome transcriptional profiling. Liver metastasis–selective genes were identified by significance analysis of microarray (SAM) analyses and independently validated in external datasets. Finally, the prognostic relevance of the liver metastasis–selective genes in primary breast cancer was tested.
Results: Liver relapse was associated with estrogen receptor (ER) expression (P = 0.002), luminal B subtype (P = 0.01), and was prognostic for an inferior postrelapse survival (P = 0.01). The major variation in the transcriptional landscape of metastases was also associated with ER expression and molecular subtype. However, liver metastases displayed unique transcriptional fingerprints, characterized by downregulation of extracellular matrix (i.e., stromal) genes. Importantly, we identified a 17-gene liver metastasis–selective signature, which was significantly and independently prognostic for shorter relapse-free (P < 0.001) and overall (P = 0.001) survival in ER-positive tumors. Remarkably, this signature remained independently prognostic for shorter relapse-free survival (P = 0.001) among luminal A tumors.
Conclusions: Extracellular matrix (stromal) genes can be used to partition breast cancer by site of relapse and may be used to further refine prognostication in ER positive primary breast cancer. Clin Cancer Res; 22(1); 146–57. ©2015 AACR.
This article is featured in Highlights of This Issue, p. 1
Although metastasis is the principal cause of cancer-related deaths, the scarcity of clinical breast cancer metastases has impeded their characterization in large genomic and transcriptomic studies. Although metastases may be genetically similar to their seeding primary tumors, distinct differences that could be exploited to improve disease control may nonetheless exist. We performed global transcriptional profiling of 91 clinical breast cancer metastases, aiming to identify genes associated with liver metastases, given the inferior outcome associated with liver recurrence. We identified a set of 17 liver metastasis–selective genes of prognostic relevance in early breast cancer. Importantly, this signature showed an independent ability of identifying patients at higher risk of recurrence and death within the luminal A molecular subtype. These patients may benefit from closer disease monitoring and may, in addition, be amenable to enrollment into clinical trials investigating novel antineoplastic therapeutics targeting features other than increased proliferation.
Introduction
Metastasis is a significant clinical and socio-economic problem, accounting for over 90% of cancer-related deaths (1). After diagnosing metastatic breast cancer (MBC), the site of recurrence is an important feature for estimating the patient's prognosis. Liver metastasis is associated with the poorest survival relative to locoregional, bone, and lung colonization (2–7). Noteworthy, the diagnosis of liver metastases is on the rise (8, 9), suggesting that available adjuvant therapies may have limited efficacy in preventing liver colonization compared with metastases at other sites. Consequently, the increasing numbers of patients presenting with these adverse events warrants a better understanding of the molecular attributes of site-specific metastases to enable the identification of novel biomarkers to guide surveillance and improve personalization of therapy.
The selection of metastatic sites is not a random process. Once disseminated, circulating tumor cells exhibit tissue-specific tropisms beyond what can be explained by normal circulatory patterns. Tissue selectivity for breast cancer metastatic colonization has been associated with primary tumor pathologic characteristics such as estrogen receptor (ER) expression and tumor molecular subtypes (10, 11). However, a marked redundancy of metastatic site selectivity prevails between these molecularly heterogeneous groups, limiting their accuracy as site-specific predictive markers.
Conventionally, at time of primary breast cancer diagnosis, the prognosis for a favorable outcome and decision for the exemption from chemotherapy is based on a combination of factors, including ER positivity, negative nodal status, small tumor size, and low histologic grade (4). Tumors displaying these favorable prognostic factors are significantly enriched within the luminal A intrinsic subtype. However, intrinsic or acquired resistance to hormonal therapy and disease recurrence to distant sites, including the liver, may eventually occur in a clinically relevant number of patients with luminal A tumors, underlining the heterogeneity even within this favorable subtype. Metastases remain the main cause of breast cancer-related mortality. It is therefore necessary to identify better prognostic biomarkers, and if possible subtype-specific prognostic biomarkers to improve individualization of therapy.
A few studies have shown that primary tumors and their metastases generally share similar copy number aberrations (12, 13) and gene expression profiles (14, 15), but these studies were underpowered by the scarcity of metastatic biopsies, limiting the identification of differences between these matched tumor pairs. By using experimental mouse models and a limited series of clinical metastatic biopsies, genes associated with the propensity of breast cancer relapse to the bone (16), lung (17, 18), and brain (19) have been published. Furthermore, we (7) and others (20, 21) have shown an association between claudin-2 expression and liver recurrence. However, because experimental mouse models incompletely capture the relevant genetic complexity of tumor progression within the human host, studies using patient-derived biopsies from metastases may reveal additional clinically relevant site-specific attributes to complement and/or validate these preliminary reports.
The aim of this study was to provide an independent characterization of the transcriptional landscape of breast cancer metastases with the specific objective to identify genes selective for breast cancer liver metastases with prognostic potential at time of primary tumor diagnosis.
Materials and Methods
Patients and tumors
The study cohort consisted of 304 women diagnosed with locally advanced (inoperable) or MBC, enrolled in a randomized phase III trial (TEX) conducted between 2002 and 2007 in Sweden. As first-line treatment for metastatic disease, patients received a combination of epirubicin and paclitaxel alone (ET) or with the addition of capecitabine (TEX). Patients presenting with brain metastases, approved for first-line HER2-targeted therapy, or diagnosed with other malignancies within 5 years of the trial commencement were exempted. Complete clinical and pathologic data were recorded in a central clinical trial database. The median follow-up for postrecurrence survival was 45 months (range, 9–135 months) for patients alive at last update (July 2013). Detailed information regarding the design and outcome of the trial has been published (22). Fine-needle aspirates (FNA) of at least one metastatic lesion were collected before commencement of treatment whenever possible. In addition, archival formalin-fixed paraffin-embedded primary tumor blocks were collected for tissue microarray (TMA) construction and central reassessment of biomarkers by immunohistochemistry and in situ hybridization techniques where applicable. Table 1 and Supplementary Table S1 show the distribution of clinicopathologic factors in the cohort.
. | All tumors . | ER-positive tumors . | ||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|
. | . | Metastatic category . | . | . | Metastatic category . | . | ||||||
Primary tumor characteristic . | N . | Locoregional . | Bone . | Lung . | Liver . | P . | N . | Locoregional . | Bone . | Lung . | Liver . | P . |
ER status | ||||||||||||
Negative | 68 | 20 | 7 | 19 | 22 | 0.002 | ||||||
Positive | 213 | 29 | 43 | 40 | 101 | |||||||
PR status | ||||||||||||
Negative | 81 | 18 | 10 | 18 | 35 | 0.21 | 44 | 8 | 5 | 7 | 24 | 0.32 |
Positive | 110 | 16 | 25 | 21 | 48 | 107 | 16 | 25 | 20 | 46 | ||
HER2 status | ||||||||||||
Negative | 179 | 32 | 32 | 37 | 78 | 0.76 | 143 | 23 | 29 | 24 | 67 | 0.39 |
Positive | 17 | 2 | 2 | 5 | 8 | 8 | 0 | 1 | 3 | 4 | ||
Number of metastatic sites | ||||||||||||
Oligo (n = 1) | 76 | 25 | 23 | 9 | 19 | <0.001 | 57 | 17 | 19 | 5 | 16 | <0.001 |
Multiple (n > 1) | 226 | 25 | 33 | 54 | 114 | 156 | 12 | 24 | 35 | 85 | ||
Histologic grade | ||||||||||||
Grade 1/2 | 80 | 9 | 17 | 12 | 42 | 0.03 | 72 | 9 | 15 | 10 | 38 | 0.58 |
Grade 3 | 105 | 24 | 15 | 26 | 40 | 67 | 10 | 13 | 15 | 29 | ||
Adjuvant endocrine therapy | ||||||||||||
No | 147 | 27 | 30 | 37 | 53 | 0.05 | 73 | 8 | 19 | 17 | 29 | 0.17 |
Yes | 154 | 22 | 26 | 26 | 80 | 140 | 21 | 24 | 23 | 72 | ||
Adjuvant chemotherapy | ||||||||||||
No | 152 | 24 | 30 | 32 | 66 | 0.93 | 112 | 15 | 22 | 21 | 54 | 0.99 |
Yes | 148 | 25 | 25 | 31 | 67 | 101 | 14 | 21 | 19 | 47 | ||
Age at primary diagnosis | ||||||||||||
<50 y | 152 | 17 | 28 | 37 | 70 | 0.08 | 107 | 12 | 19 | 24 | 52 | 0.38 |
≥50 y | 149 | 32 | 28 | 26 | 63 | 106 | 17 | 24 | 16 | 49 | ||
Metastasis-free interval | ||||||||||||
≤24 mo | 80 | 15 | 15 | 16 | 34 | 0.91 | 55 | 7 | 11 | 11 | 26 | 0.99 |
>24 mo | 221 | 34 | 41 | 47 | 99 | 158 | 22 | 32 | 29 | 75 | ||
Nodal status | ||||||||||||
N0 | 91 | 14 | 17 | 22 | 38 | 0.69 | 65 | 7 | 12 | 16 | 30 | 0.47 |
N+ | 202 | 35 | 37 | 37 | 93 | 145 | 22 | 30 | 23 | 70 | ||
Tumor size | ||||||||||||
≤20 mm | 119 | 14 | 26 | 21 | 58 | 0.21 | 82 | 9 | 20 | 12 | 41 | 0.47 |
>20 mm | 178 | 33 | 30 | 40 | 75 | 128 | 18 | 23 | 27 | 60 | ||
Molecular subtypea | ||||||||||||
Luminal A-like | 65 | 9 | 19 | 11 | 26 | 0.01 | 65 | 9 | 19 | 11 | 26 | 0.05 |
Luminal B-like | 81 | 13 | 9 | 16 | 43 | 81 | 13 | 9 | 16 | 43 | ||
HER2 positive | 9 | 2 | 1 | 2 | 4 | |||||||
Triple negative | 24 | 8 | 2 | 9 | 5 |
. | All tumors . | ER-positive tumors . | ||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|
. | . | Metastatic category . | . | . | Metastatic category . | . | ||||||
Primary tumor characteristic . | N . | Locoregional . | Bone . | Lung . | Liver . | P . | N . | Locoregional . | Bone . | Lung . | Liver . | P . |
ER status | ||||||||||||
Negative | 68 | 20 | 7 | 19 | 22 | 0.002 | ||||||
Positive | 213 | 29 | 43 | 40 | 101 | |||||||
PR status | ||||||||||||
Negative | 81 | 18 | 10 | 18 | 35 | 0.21 | 44 | 8 | 5 | 7 | 24 | 0.32 |
Positive | 110 | 16 | 25 | 21 | 48 | 107 | 16 | 25 | 20 | 46 | ||
HER2 status | ||||||||||||
Negative | 179 | 32 | 32 | 37 | 78 | 0.76 | 143 | 23 | 29 | 24 | 67 | 0.39 |
Positive | 17 | 2 | 2 | 5 | 8 | 8 | 0 | 1 | 3 | 4 | ||
Number of metastatic sites | ||||||||||||
Oligo (n = 1) | 76 | 25 | 23 | 9 | 19 | <0.001 | 57 | 17 | 19 | 5 | 16 | <0.001 |
Multiple (n > 1) | 226 | 25 | 33 | 54 | 114 | 156 | 12 | 24 | 35 | 85 | ||
Histologic grade | ||||||||||||
Grade 1/2 | 80 | 9 | 17 | 12 | 42 | 0.03 | 72 | 9 | 15 | 10 | 38 | 0.58 |
Grade 3 | 105 | 24 | 15 | 26 | 40 | 67 | 10 | 13 | 15 | 29 | ||
Adjuvant endocrine therapy | ||||||||||||
No | 147 | 27 | 30 | 37 | 53 | 0.05 | 73 | 8 | 19 | 17 | 29 | 0.17 |
Yes | 154 | 22 | 26 | 26 | 80 | 140 | 21 | 24 | 23 | 72 | ||
Adjuvant chemotherapy | ||||||||||||
No | 152 | 24 | 30 | 32 | 66 | 0.93 | 112 | 15 | 22 | 21 | 54 | 0.99 |
Yes | 148 | 25 | 25 | 31 | 67 | 101 | 14 | 21 | 19 | 47 | ||
Age at primary diagnosis | ||||||||||||
<50 y | 152 | 17 | 28 | 37 | 70 | 0.08 | 107 | 12 | 19 | 24 | 52 | 0.38 |
≥50 y | 149 | 32 | 28 | 26 | 63 | 106 | 17 | 24 | 16 | 49 | ||
Metastasis-free interval | ||||||||||||
≤24 mo | 80 | 15 | 15 | 16 | 34 | 0.91 | 55 | 7 | 11 | 11 | 26 | 0.99 |
>24 mo | 221 | 34 | 41 | 47 | 99 | 158 | 22 | 32 | 29 | 75 | ||
Nodal status | ||||||||||||
N0 | 91 | 14 | 17 | 22 | 38 | 0.69 | 65 | 7 | 12 | 16 | 30 | 0.47 |
N+ | 202 | 35 | 37 | 37 | 93 | 145 | 22 | 30 | 23 | 70 | ||
Tumor size | ||||||||||||
≤20 mm | 119 | 14 | 26 | 21 | 58 | 0.21 | 82 | 9 | 20 | 12 | 41 | 0.47 |
>20 mm | 178 | 33 | 30 | 40 | 75 | 128 | 18 | 23 | 27 | 60 | ||
Molecular subtypea | ||||||||||||
Luminal A-like | 65 | 9 | 19 | 11 | 26 | 0.01 | 65 | 9 | 19 | 11 | 26 | 0.05 |
Luminal B-like | 81 | 13 | 9 | 16 | 43 | 81 | 13 | 9 | 16 | 43 | ||
HER2 positive | 9 | 2 | 1 | 2 | 4 | |||||||
Triple negative | 24 | 8 | 2 | 9 | 5 |
aMolecular subtyping using immunohistochemical staining for ER, PR, HER2, and Ki67 according to the 2013 St Gallen consensus guidelines. Patients were categorized according to the most advanced metastatic site affected (locoregional, locally advanced, or regional metastases in the lymph nodes or skin; bone, skeletal metastases with or without locoregional metastases; lung, pleural metastases with or without skeletal and locoregional metastases; liver, hepatic metastases with or without pleural, skeletal or locoregional metastases). P values are from Fisher exact tests. Values in bold are significant.
Ethics statement
This substudy was approved by all the regional ethics committees at the participating hospitals [Karolinska Institutet Stockholm (Karolinska, Sweden; KI 02-205 and 02-206); Sahlgrenska University Hospital (Gothenburg, Sweden; M090-02 and M091-02); Linköping University Hospital (Linköping, Sweden; 02-519 and 02-339); Örebro University Hospital (Örebro, Sweden; 308/02 and 308/03); Umeå University Hospital (Umeå, Sweden; Um 02-336 and Um 03-03), and Lund University Hospital (Lund, Sweden; LU 290-02 and LU 291-02)]. All patients provided written informed consent to participate in the clinical trial and translational studies. This study adheres to the REMARK guidelines for reporting prognostic biomarker studies (23).
RNA extraction and gene expression microarrays
Tumor cellularity of FNAs was assessed by a cytologist (L. Skoog) on Giemsa stained, ethanol-fixed, cytospin preparations and total RNA was extracted from samples with high (>50%) tumor cell content using the Qiagen RNA Mini Kit (Qiagen) following the manufacturer's recommendations. RNA quantity and integrity were analyzed on the NanoDrop spectrophotometer (NanoDrop Technologies) and the Agilent 2100 Bioanalyzer (Agilent), respectively, and cDNA was generated and biotin-labeled using the NuGen 50 ng amplification protocol (Covance Genomics Laboratory). Labeled cDNA was hybridized onto custom-made whole-genome Affymetrix HuRSTA-2a520709 gene chips following the GeneChip Hybridization, Wash, and Stain Kit protocol (Affymetrix). Data preprocessing and normalization were performed using the robust multichip average (RMA) algorithm. After normalization, a presence filter was applied to select only features present in ≥90% of assays, and features with low intensities (below the median intensity for Y chromosome gene probes) were filtered out. The data were log2 transformed and only transcripts showing high variance across assays were selected (variance filter SD ± 1), leaving a final dataset with 8,339 features representing 5,232 unique gene variants for further analyses. All processes were performed using packages in R (24) and the TM4 microarray software suite (25). The final dataset included 91 samples from 85 patients [liver (n = 16), bone (n = 5), lung (n = 2), lymph node (n = 39), local [breast (n = 11), skin (n = 17), and ascites (n = 1)]. The distribution of baseline clinicopathologic features in the original study cohort (n = 304) and the subpopulation included in the transcriptional profiling study (n = 85) is presented in Supplementary Table S1. Raw and processed data have been deposited in the Gene Expression Omnibus (GSE46141).
Multivariable data analyses
Unsupervised analyses.
Principal component analysis (PCA) was performed using SIMCA P version 13.0.2 software package (Umetrics AB). The dataset was mean-centered across rows (genes), unit variance scaled and model complexity was estimated by leave-one-out cross-validation. Unsupervised hierarchical clustering (HCL) was performed using the Pearson correlation distance metric and average linkage.
Supervised analyses.
The intrinsic molecular subtypes of the metastases were determined using the research-based PAM50 algorithm as previously described (26). A two-class significance analysis of microarray (SAM; ref. 27) analysis was performed to identify significant differentially expressed genes in liver metastases compared with metastases from other sites. The liver-selectivity of the identified genes was verified in an external dataset of 36 breast cancer metastases (GSE14018; ref. 28). The biologic processes and pathways enriched among the liver metastasis-selective genes were uncovered by gene ontology analysis using the DAVID (29, 30) database. Furthermore, the activity of eight gene expression-based modules representing relevant breast cancer-specific biologic processes (stroma, lipid, immune response, mitotic progression, mitotic checkpoint, basal, early response, and steroid response; ref. 31) was assessed in the metastases.
In a final step, candidate liver metastasis–selective genes that may serve as biomarkers for predicting the liver metastatic potential of a primary tumor were identified using an external dataset of 192 primary breast tumors (GSE12276; ref. 15) and associations between these candidate genes and outcome in early breast cancer were independently tested using Gene expression-based Outcome for Breast cancer Online (GOBO; ref. 32), an online tool for validation of the prognostic value of single genes or sets of genes in primary breast cancer (n = 1,881).
Survival analyses.
Kaplan–Meier plots were generated and the log-rank test was used to check for statistically significant differences between target groups. Cox-proportional hazards models were used to evaluate the independent prognostic significance of biomarkers, adjusting for conventional prognostic factors. P values correspond to two-sided statistical tests and values <0.05 were considered significant.
Results
Associations between primary tumor clinicopathologic factors and the first site(s) of recurrence
Because many patients with MBC present with relapses in more than one anatomic site at time of first metastasis diagnosis, we classified patients into four metastatic categories reflecting the most advanced site affected at first clinical presentation. These categories were: locoregional (locally advanced or regional metastases in the lymph nodes or skin), bone (skeletal metastases with or without locoregional disease), lung (lung parenchymal/plural metastases with or without bone and locoregional disease), and liver (hepatic metastases with or without lung, bone or locoregional metastases). Associations between primary tumor clinicopathologic factors and the first site(s) of recurrence are shown in Table 1. ER positivity was found to be associated with bone and liver recurrences, while negative ER status correlated with locoregional and lung relapses (Fisher exact, P = 0.002). Liver recurrence was also common among patients with HER2 positive tumors (8/17), but this association was not statistically significant (probably due to the limited number of HER2-positive tumors in the study). Locoregional and bone metastases were often detected as oligo-metastases, while liver and lung metastases were often diagnosed in parallel with deposits at other sites (P < 0.001). Furthermore, low histologic grade (grades 1 and 2) was associated with bone and liver recurrences, while high grade (grade 3) correlated with locoregional and lung relapses (P = 0.03). However, no significant association between histologic grade and metastatic site was observed when ER-positive tumors were analyzed separately (P = 0.58). When the surrogate (IHC-based) molecular subtype of the primary tumor was considered, bone and hepatic recurrences were found to be associated with luminal-like (A and B) tumors, while relapses to the lung and locoregional sites were associated with the triple-negative subtype (P = 0.01). Subanalyses within ER-positive tumors revealed a borderline association of bone metastases with the luminal A–like and liver metastases with the luminal B–like subtypes, respectively (P = 0.05). Overall, these results confirm that conventional tumor pathologic biomarkers provide important insights into a primary tumor's metastatic propensity, with liver relapse commonly associated with poor prognostic pathologic features. Nevertheless, in this cohort, a remarkably high prevalence of liver metastases was noted among patients who presented with primary tumors with favorable prognostic features; 40% of luminal A and 53% of histologic grade 1 and 2 tumors progressed to the liver.
Liver-only relapse is associated with a relatively better outcome
We recently reported that liver relapse was associated with inferior survival after recurrence in the present study cohort (7). However, some studies suggest that patients with liver-only metastatic disease may experience longer survival compared with patients harboring liver metastases in parallel with metastases in other organs (5, 9). We found a similar trend in this cohort [Fig. 1; log-rank, P = 0.01, Multivariable Cox model P < 0.001 adjusting for age (>50 years or ≤50 years), metastasis-free interval (≤2 years or >2 years), nodal status, adjuvant endocrine therapy, and adjuvant chemotherapy], emphasizing the significance of tumor burden in addition to metastatic site for postrelapse survival.
Identification of shared and distinct transcriptional portraits of site-specific metastases
The transcriptional landscape of breast cancer metastases has generally been inferred from primary tumors due to scarcity of clinical biopsies from metastases to perform independent studies. PCA analyses revealed that the first three principal components partitioned breast cancer metastases into groups that were strongly associated with primary tumor ER expression (Fig. 2A) and the intrinsic molecular subtypes of the metastases (Fig. 2B). Remarkably, liver metastases were the only class that was tightly clustered in the PCA score plot (Fig. 2C), indicating a transcriptional distinction relative to other metastases. A similar tight clustering pattern for liver metastases was observed by unsupervised HCL of the samples using the top 3,000 most variable probes (Fig. 2D). Of note, all biologic replicates (independent metastatic biopsies from the same patient) clustered together pair-wise and adjacent to each other in the sample dendrogram, confirming that transcriptional profiles of intraindividual tumors are more similar than interindividual profiles.
To provide more insight into the biology of site-specific metastases, the activity of eight gene modules representing key biologic aspects associated with breast cancer (31) was compared between the specific metastatic sites. Four modules were found to be significantly differentially expressed between the metastatic sites (Fig. 3A–D). Liver metastases displayed a significantly lower expression of the “stroma” (relative to the skin and lymph nodes, adjusted P = 0.01 and P = 0.001, respectively); “basal” (relative to skin, P = 0.043) and “early response” (relative to bone, P = 0.005) modules, and a higher expression of the “steroid response” module relative to metastases in the skin (P = 0.003). Considering that the transcriptional profiles of independent tumors from the same individual are highly similar, the activities of the eight gene modules were next compared between samples classified according to the four metastatic categories as previously defined in Table 1. Similarly, differential expression of the same four modules was observed (Fig. 3E–H). Low expression of the “stroma” module was observed in the liver category relative to bone (adjusted P = 0.015). In addition, the “basal” module was elevated in the lung category relative to the liver (P = 0.018), while “steroid response” was higher in the liver and bone categories relative to the lung category (P = 0.013 and P = 0.017, respectively). These results further confirm the association between the metastatic site and ER expression or molecular subtype because the “basal” and “steroid response” modules were shown to be strongly associated with the basal-like (ER−) and luminal subtypes (ER+), respectively (31).
Identification of liver metastasis–selective genes
To further dissect the transcriptional distinctiveness of liver metastases, a two-class SAM analysis was performed comparing liver versus other metastases. This analysis was restricted to ER-positive primary tumors, because the liver metastases were mainly of this phenotype (12/16; with 3 missing ER status; Fig. 2). We found 358 genes to be differentially expressed (309 upregulated and 49 downregulated, FDR = 0.1; Supplementary Table S2); henceforth, referred to as “breast cancer liver metastasis–selective genes.” HCL of an independent set of 36 breast cancer metastases (GSE14018; ref. 28) using only these 358 genes revealed a similar expression pattern in liver metastases (Supplementary Fig. S1), thus confirming their liver selectivity.
Gene set enrichment analysis (29, 30) showed significant upregulation of biologic processes, including endopeptidase inhibitor activity, complement activation, blood coagulation, immune response, and steroid metabolism in liver metastases (Supplementary Table S3). Conversely, processes associated with extracellular matrix, biologic adhesion, skeletal system development, and blood vessel development were enriched among the downregulated genes in liver metastases (Supplementary Table S3). To ascertain that the enriched upregulated biologic processes, which are also common biologic processes occurring in normal liver, was not a reflection of normal tissue contamination, we performed unsupervised HCL of all samples in the test cohort using previously reported normal breast and liver tissue-specific genes (33), as well as breast cancer-selective genes (34), respectively. Reassuringly, even though the liver metastases formed a distinct cluster in the sample dendrogram when clustered using the normal liver genes (Supplementary Fig. S2A), no separation of the samples based on metastatic site was seen upon clustering with normal breast (Supplementary Fig. S2B) or breast cancer–specific genes (Supplementary Fig. S2C). Instead, clustering correlated with other biologic characteristics, such as ER expression and molecular subtype. These results suggest that breast cancer liver metastases maintain a transcriptional profile consistent with the site of origin of the tumor cells (breast), and in addition adopt other transcriptional features associated with the metastatic microenvironment (liver) which may be important for their survival at this foreign site.
Associations between breast cancer liver metastasis–selective genes and primary tumor clinicobiologic factors and clinical outcome
Robust tissue-specific metastasis biomarkers may be detectable in the primary tumors of patients who eventually develop metastases in the corresponding target organ. Gene signatures with the potential to predict breast cancer metastasis to the lung, bone, and brain (16, 17, 19) have been reported. Using an external primary breast cancer dataset including only patients with metastatic disease and for whom the annotation of the site(s) of metastasis was recorded (15, 19), we performed a restricted analysis using only the 358 liver metastasis–selective genes. Only ER-positive tumors (n = 119) were interrogated. Of note, 347 of the 358 genes could be mapped across datasets. We found 17 genes to be significantly (FDR < 0.05) differentially overexpressed in tumors relapsing in the liver. This list was enriched for genes involved in cadherin and integrin signaling pathways, as well as in skeletal system development. Of note, 6 of 17 (CDH11, COL11A1, FBN1, MFAP5, SFRP4, and SPON1) genes overlapped with the previously described “stroma” module (Spearman correlation coefficient, 0.7). Figure 4 shows the expression of these 17 genes in both the test and the validation cohorts. Surprisingly, while all 17 genes were upregulated in primary tumors with liver metastatic potential, 14 of 17 genes were downregulated in liver metastases in both the test and the validation (GSE14018) metastasis datasets.
Finally, in GOBO (32), a database containing 1,881 annotated primary breast tumors, we aimed to identify relevant associations between the 17-gene signature and other primary tumor pathologic features and prognosis. The expression of the signature was heterogeneous between different molecular subtypes and histologic grades (Fig. 5A–D; ANOVA, P < 0.00001), with a significantly lower expression in luminal B and basal-like tumors compared with the other subtypes (adjusted P < 0.0001 for all pairwise comparisons of luminal B or basal tumors vs. other subtypes). In addition, low expression was significantly correlated with high histologic grade (Fig. 5C and D; adjusted, P < 0.0001 for pairwise comparisons between grade 3 tumors vs. grades 1 and 2). Exploratory analyses revealed significant differential expression of the 17-gene signature across the recently described IntClust subgroups (Supplementary Fig. S3A; refs. 35, 36). Decreased expression was observed in IntClust subgroups 10, 1, and 9 relative to IntClust 3 and 4 (adjusted P < 0.0001 for all pairwise comparisons). Furthermore, among ER-positive tumors, a significantly lower expression was also noted in subgroups 7 and 8 relative to subgroups 3 and 4 (P < 0.0001).
Remarkably, low expression of the 17-gene signature was significantly associated with shorter recurrence-free survival (RFS; Fig. 5E; log-rank, P = 3 × 10−5; Supplementary Table S4; multivariable Cox model HR, 1.5; P = 0.001) and overall survival (OS; Fig. 5G, log-rank, P = 0.00927; Supplementary Table S5; multivariable Cox model HR, 1.4; P = 0.026) in patients with ER-positive tumors. More importantly, the 17-gene signature remained significantly and independently prognostic for RFS when the subset of luminal A tumors was analyzed separately (Fig. 5F; log-rank, P = 0.00097; Supplementary Table S5; multivariable Cox model HR, 2.2; P = 0.004). A trend toward poor OS for patients with luminal A tumors with low expression of the 17-gene signature was observed in univariable analysis (Fig. 5H; log-rank, P = 0.083; Supplementary Table S4; multivariable Cox model HR, 1.4; P = 0.29). In subanalyses restricted to tumors in IntClust subgroups 3, 7, and 8 (those highly enriched for luminal A tumors), low expression of the 17-gene signature was associated with an inferior RFS (Supplementary Fig. S3B; log-rank, P = 0.001) and OS (log-rank, P = 0.06). Exploratory analyses confirmed the association of the signature with poor prognosis when all tumors were included in the analysis, irrespective of ER status (Supplementary Fig. S3C and S3D; RFS, P = 0.00012; OS, P = 0.01872).
Discussion
In this study, we identified a 17-gene signature enriched for extracellular matrix or stroma genes, the majority of which were selectively downregulated in breast cancer liver metastases. Furthermore, downregulation in primary tumors, irrespective of site of relapse, was associated with aggressive tumor biologic features and inferior prognosis. Liver metastases are deleterious, leading to the early demise of MBC patients (2–7). We observed significant positive associations between liver recurrence and poor tumor biologic characteristics, including luminal B subtype, high histologic grade, and large tumor burden. However, despite the statistically significant associations, the prevalence of liver relapses was notable in all subgroups, indicating a low specificity and sensitivity of these factors for accurate metastatic site prediction. Because liver relapse is indicative of inferior postrecurrence survival, there is a need for more specific and independent biomarkers to identify patients at risk. Recently, we demonstrated that CLDN2, which is significantly upregulated in liver metastases, is an independent prognostic factor for early liver recurrence in breast cancer (7). Here, we show that downregulation of various genes involved in cell adhesion is characteristic of liver metastases.
Patients with liver-only metastatic disease had a better postrecurrence survival compared with those harboring liver metastases in parallel to metastases in other organs. This finding corroborates results from other studies (5, 9). There is great interest in evaluating local treatment options such as surgery or stereotactic radiotherapy in patients with oligo-metastases in the liver but randomized studies are needed to evaluate the efficacy of these treatment options.
Transcriptional profiling has increased our understanding of the biology of organ-specific metastases and has led to the identification of site-specific metastasis genes and signatures (16, 17, 19, 20). PCA and unsupervised HCL analyses reported herein revealed that the major variation across breast cancer metastases was strongly associated with ER status and molecular subtype, an observation consistent with the conventional understanding of breast cancer biology. This similarity underscores that primary tumor molecular traits are conserved across stages of tumor progression. Interestingly, we observed minor but significant site-specific differences at the transcriptional level, which reflects additional alterations acquired by breast cancer cells to thrive and evolve into overt metastases in the foreign milieu. Interestingly, our data suggest that mimicry of “normal processes” of the new microenvironment may be a necessary adaptation. An enrichment of genes and biologic processes commonly observed in normal liver was noted among upregulated genes. Most of these genes code for signaling peptides commonly found in the extracellular space, further highlighting the importance of the microenvironment in metastatic colonization. The deregulation of genes that mimic target organ functions has previously been observed in other studies investigating the organ-specificity of metastases. Differential expression of genes important for ossification in bone metastases (16, 37), brain metabolism in brain metastases (19), pulmonary function in lungs (17, 18), and liver function in liver metastases (20) have been reported. This phenomenon can be interpreted within the confinements of the “seed and soil theory” of tumor invasion and metastatic colonization (38). Of note, mimicry of target-organ properties was observed even when pure tumor cell line populations displaying distinct site-specific preferences were studied (16, 17, 19, 20), suggesting that part of this expression profile is indeed intrinsic to the tumor cells. Furthermore, we did not observe any segregation of our samples according to metastatic site when subjected to HCL on normal breast (33) or breast cancer (34) selective genes, confirming that all samples were enriched for breast cancer cells and that the transcriptional profiles observed are most likely mainly tumor cell intrinsic. Nonetheless, the possibility of normal tissue contamination cannot be completely ruled out. On the other hand, downregulation of extracellular matrix genes and genes involved in cell adhesion and the development of blood vessels and the skeletal system, which are all processes that have been linked with invasion and metastasis in breast cancer (39) was seen in liver metastases. Of note, the top downregulated gene was the epithelial mesenchymal transition inducer PRRX1, recently reported to play an important role in metastatic colonization through repression of its expression to favor reversion of the mesenchymal phenotype that is necessary for the outgrowth of metastases (40). However, analyses performed in the external dataset that included many more lung metastases and in addition included brain metastases, suggested that downregulation of these genes may also be a trait of lung and brain, but not bone metastases. Further studies are necessary to investigate this phenomenon.
Predicting the future metastatic site(s) of a primary breast cancer is multifaceted and challenging. In their recent study aimed at unraveling how bone-specific metastatic traits arise in the primary tumor, Zhang and colleagues (37) showed that stromal signals resembling those of the distant target organ play important roles at the primary tumor site to prime cells for colonizing of a specific metastatic niche. Also, three independent gene modules enriched for extracellular matrix (i.e., stroma) genes were among the 11 gene modules recently identified to shape the transcriptional landscape of primary breast cancer (41). Interestingly, in this study (41), only expression of the ECM modules showed significant associations with the site of recurrence, although liver metastases were not annotated in this study. Our 17-gene signature was enriched for stroma-related genes and was significantly correlated to the stroma module described by Fredlund and colleagues (31). Consistent with our results, they found that low expression of the stroma module was associated with shorter distant metastasis-free survival among patients with luminal A primary tumors (31). Furthermore, an independent study by Bergamaschi and colleagues (42) identified four extracellular matrix gene modules (ECM1–ECM4) with prognostic significance in ER-positive (luminal) breast cancer, but their survival analyses were not stratified to assess differences between the luminal subtypes. Of note, downregulation of several genes in our signature was characteristic of the ECM1 module (42), which was associated with the poorest outcome. Taken together, these studies highlight the possibility of harnessing the heterogeneity in the expression of extracellular matrix (stroma) genes to improve prognostication in hormone receptor positive disease.
Currently, prediction of the prognosis in ER-positive breast cancer at the transcriptional level is limited to the expression of proliferation-related genes, but high proliferative rate alone is not sufficient to account for all the recurrences observed among patients with ER-positive breast cancer, especially among patients with luminal A tumors which are generally of a low proliferating phenotype. Downregulation of the 17-gene signature was indirectly associated with high proliferation, because features such as high histologic grade and luminal B subtype are common to proliferative tumors. Consequently, low expression was independently prognostic of shorter time to recurrence and shorter OS among patients with ER-positive tumors. Remarkably, the 17-gene signature and tumor size were the only independently prognostic factors for early recurrence among patients with luminal A (low proliferative) tumors in multivariable analyses. Importantly, the luminal A tumors in this cohort were mostly of histologic grades 1 and 2. The significantly lower expression of the 17-gene signature in IntClust subgroups 3, 7, and 8, which are predominantly composed of luminal A tumors, confirms that the IntClust subtypes may also be used to further stratify luminal A tumors into groups with distinct outcome. InClust 3 is mainly characterized by low genomic instability, while IntClust 7 and 8 harbor the characteristic (“luminal”) 16p gain/16q loss and 1q gain/16q loss aberrations, respectively. Interestingly, the 17-gene signature captures the diversity in prognosis even within these well-characterized subgroups. Metastases remain the main cause of death from cancer. The goal of individualizing therapy for breast cancer can only be achieved if all patients at risk can be accurately identified. The prognostic relevance of the 17-gene signature in luminal A breast cancer holds great promise in this context and needs to be independently validated.
The fact that all 17 genes in our signature were found to be overexpressed in the group of primary tumors from patients who subsequently developed liver metastases is surprising because the majority of the genes showed low expression in the liver metastases. The SAM analysis comparing metastases from specific organs, that is, liver versus other sites, disregards the fact that the same patient from whom the liver metastasis was collected may have metastases in other organs. Also, in this study, we confirm that paired tumors from the same individual have highly similar global transcriptional profiles. Taken together, the high concordance in global transcription and the fact that SAM analysis only detects differences in levels of gene expression between groups and not an absolute presence or absence thereof, argues that the genes we identified are more liver-selective and therefore likely not uniquely liver-specific per se. Furthermore, searching for the expression of site-selective genes in primary tumors, which represent a heterogeneous mix of clones with diverse site-specific metastatic propensities is also complex. In the primary tumor cohort used to identify the subset of liver-selective genes differentially expressed at this early time point during tumor progression, the subcategorization of patients was also confounded by intraindividual overlap of several metastatic sites. Nonetheless, the inverse correlation in the direction of expression of many of the genes between primary tumors and metastases is intriguing and requires further functional investigation. However, importantly, low expression as observed in the liver metastases was prognostic of an inferior outcome. Because decreased expression of most of the genes (as observed in the liver metastases) is associated with inferior outcome, we hypothesize that the lower expression may be a stronger marker of overall inferior prognosis rather than only a marker for liver-specific recurrence. This is in line with the understanding that liver metastasis is an indicator of poor prognosis. The scarcity of datasets with annotations for the metastatic site(s) hindered an independent evaluation of the ability of this signature to specifically predict breast cancer liver recurrence.
Diagnosis of liver-only metastases in breast cancer is not common and liver metastases are frequently diagnosed in tandem with other sites as can be seen in our patient cohort where only 19 of 133 (14%) patients presented with liver-only disease at first diagnosis of metastatic disease. This suggests that liver metastases and tumor burden are strongly associated and our signature may to some extent be associated with tumor burden. Of importance however, the liver metastases clustered together and displayed similar transcriptional profiles regardless of whether they were diagnosed as oligo-metastases or in parallel to other known metastatic deposits, supporting the liver selectivity of the identified gene signature. Identification of independent site-specific signatures would therefore require a well-annotated and sufficiently large cohort of patients with oligo-metastatic disease, which is challenging given the scarcity of patients presenting with oligo-metastases as well as the fact that biopsies are seldom taken from patients presenting with oligo-metastatic disease. Biopsies of metastases are now routinely collected whenever possible for reassessment of biomarkers to guide treatment for MBC. Ultimately, the gap of scarcity of these samples will be bridged and larger collections of metastases will become available for research purposes, enabling, for example, validation of the data presented herein. Notwithstanding this limitation, our analysis pipeline enabled us to identify a biologically important gene set, the clinical relevance of which was independently validated in a large cohort of primary breast cancer.
In conclusion, we have identified a 17-gene signature enriched for genes selectively underexpressed in breast cancer liver metastases, with a remarkable ability to independently identify patients with luminal A primary breast cancers who may benefit from closer disease monitoring and may in addition be candidates for enrollment into clinical trials investigating novel targeted therapies. Further studies are warranted to validate our results especially in more recently diagnosed patient series to adjust for modern advances in adjuvant breast cancer management.
Disclosure of Potential Conflicts of Interest
No potential conflicts of interest were disclosed.
Disclaimer
None of the funding agencies were involved in study design, data collection, analysis or interpretation, or in the writing and submission of the report for publication.
Authors' Contributions
Conception and design: L. Carlsson, Z. Einbeigi, P.-O. Malmström, T.M. Walz, M. Fernö, T. Hatschek, I. Hedenfalk
Development of methodology: S. Kimbung, L. Carlsson
Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): S. Kimbung, S.E. Brage, M.F. Stolt, L. Skoog, L. Carlsson, Z. Einbeigi, B. Linderholm, N. Loman, P.-O. Malmström, M. Söderberg, T.M. Walz, M. Fernö, T. Hatschek, I. Hedenfalk
Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): S. Kimbung, I. Johansson, A. Danielsson, S. Veerla, L. Carlsson, N. Loman, T. Hatschek, I. Hedenfalk
Writing, review, and/or revision of the manuscript: S. Kimbung, I. Johansson, A. Danielsson, S.E. Brage, L. Carlsson, B. Linderholm, N. Loman, P.-O. Malmström, M. Söderberg, T.M. Walz, M. Fernö, T. Hatschek, I. Hedenfalk
Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): S. Kimbung, A. Danielsson, S.E. Brage, M.F. Stolt, L. Carlsson, Z. Einbeigi, M. Fernö, T. Hatschek, I. Hedenfalk
Study supervision: L. Carlsson, I. Hedenfalk
Other (discussed and introduced some changes in the article): E. Lidbrink
Other (Principal Investigator of the TEX trial): T. Hatschek
Acknowledgments
Microarray experiments were performed at Merck Inc. The authors are also indebted to the TEX Trialists Group [Coordinating Investigator: Thomas Hatschek; Translational research: Mårten Fernö, Linda Lindström, Ingrid Hedenfalk; QoL: Yvonne Brandberg; Statistics: John Carstensen; Laboratory: Suzanne Egyházy, Marianne Frostvik Stolt, Lambert Skoog; Clinical Trial Office: Mats Hellström, Maarit Maliniemi, Helene Svensson; Radiology: Gunnar Åström; Karolinska University Hospital (Stockholm, Sweden): Jonas Bergh, Judith Bjöhle, Elisabet Lidbrink, Sam Rotstein, Birgitta Wallberg; Sahlgrenska University Hospital: Zakaria Einbeigi, Per Karlsson, Barbro Linderholm; Linköping University Hospital: Thomas Walz; Malmö University Hospital (Malmö, Sweden): Martin Söderberg; Lund University Hospital: Niklas Loman, Per Malmström; Helsingborg General Hospital (Helsingborg, Sweden): Martin Malmberg; Sundsvall General Hospital (Sundsvall, Sweden): Lena Carlsson; Umeå University Hospital: Birgitta Lindh; Kalmar General Hospital (Kalmar, Sweden): Marie Sundqvist; Karlstad General Hospital (Karlstad, Sweden): Lena Malmberg] for providing samples and clinical data.
Grant Support
This work was supported by grants from the Swedish Cancer Society, the Swedish Research Council, the Gunnar Nilsson Cancer Foundation, the Berta Kamprad Foundation, the Gyllenstierna Krapperup's Foundation, the Swedish Cancer and Allergy Foundation, the Research Funds at Radiumhemmet, the Swedish Breast Cancer Association (BRO), ALF/FOU research funds at the Karolinska Institutet and Stockholm County Council, and unrestricted grants from Bristol-Myers Squibb Sweden AB, Pfizer Sweden AB, and Roche Sweden AB.
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.