Abstract
Aromatase inhibitors are the major first-line treatment of estrogen receptor–positive breast cancer, but resistance to treatment is common. To date, no biomarkers have been validated clinically to guide subsequent therapy in these patients. In this study, we mapped the genome-wide chromatin-binding profiles of estrogen receptor α (ERα), along with the epigenetic modifications H3K4me3 and H3K27me3, that are responsible for determining gene transcription (n = 12). Differential binding patterns of ERα, H3K4me3, and H3K27me3 were enriched between patients with good or poor outcomes after aromatase inhibition. ERα and H3K27me3 patterns were validated in an additional independent set of breast cancer cases (n = 10). We coupled these patterns to array-based proximal gene expression and progression-free survival data derived from a further independent cohort of 72 aromatase inhibitor–treated patients. Through this approach, we determined that the ERα and H3K27me3 profiles predicted the treatment outcomes for first-line aromatase inhibitors. In contrast, the H3K4me3 pattern identified was not similarly informative. The classification potential of these genes was only partially preserved in a cohort of 101 patients who received first-line tamoxifen treatment, suggesting some treatment selectivity in patient classification. Cancer Res; 73(22); 6632–41. ©2013 AACR.
Introduction
Breast cancer is the most frequently diagnosed malignancy among women worldwide, with annually around 1.4 million new cases and half a million patients who die from the disease each year (1). Seventy-five percent of all breast tumors are of the luminal subtype and tumor cell proliferation is thought to depend on activity of the estrogen receptor α (ERα). Inhibition of ERα by endocrine therapy is therefore a major treatment modality for these tumors, either by tamoxifen or the state-of-the-art aromatase inhibitors. Although many studies have focused on defining predictive markers for tamoxifen resistance, relatively little is known about the molecular determinants of aromatase inhibitor response. Such knowledge is essential as intrinsic and acquired resistance to treatment is common (2, 3). Whole-genome sequencing analyses on breast tumor samples revealed a set of 18 genes to be mutated in patients with breast cancer that correlated with differential survival upon aromatase inhibitor treatment (4), including PI3K, TP53, MAP3K1, and GATA3. These data were supported by other reports, indicating differential gene expression pathways to be enriched in poor versus good outcome patients upon aromatase inhibitor treatment (5, 6).
ERα–binding profile assessment in breast cancer cell lines has greatly increased our knowledge of hormonal receptor action. ERα rarely binds promoters and most estrogen receptor/chromatin interactions occur at distal enhancers (7), which are involved in chromatin loop structures to regulate gene expression (8). ERα/chromatin interactions require the functional involvement of other proteins, including FOXA1. FOXA1 is one of the key luminal-defining transcription factors (9, 10) and acts as pioneer factor for ERα by untangling chromatin structures, which enables ERα to bind its targets as shown in cell lines (11) and breast tumor samples (12). Other pioneer factors have been described for ERα function, including AP-2γ (13) and PBX1 (14), which function in synergy with FOXA1 in enabling ERα/chromatin interactions and subsequent activity. These reports tightly link ERα action with chromatin structure.
The epigenetic regulation of chromatin structure is a highly complex interplay between multiple histone modifiers and their downstream histone modifications, each with their intrinsically distinct patterns and functions (15). Two of the best-studied histone marks are the repressive trimethylation on lysine 27 of histone 3 (H3K27me3; ref. 16) and the activating trimethylation on lysine 4 on histone 3 (H3K4me3; ref. 17). H3K4me3 is found enriched at promoter regions (18, 19), whereas H3K27me3 can be found over large genomic regions spanning one or more epigenetically silenced genes (20, 21), which is directly regulated by the EZH2-driven polycomb repressive complex PRC2 (22).
To test for a possible interplay between differential ERα binding and epigenetic regulation of gene activity, we determined the genome-wide chromatin-binding patterns of ERα, H3K4me3, and H3K27me3 in primary breast tumor samples. Following relapse, patients were treated with an aromatase inhibitor for metastatic disease and time to progression (TTP) was correlated with altered ERα/H3K4me3/H3K27me3 binding profiles to identify distinct patterns that could hallmark aromatase inhibitor resistance.
Materials and Methods
Patients, tumor samples, and processing
The Erasmus University Medical Center (EMC; Rotterdam, the Netherlands), the Netherlands Cancer Institute (Amsterdam, the Netherlands), and the Translational Cancer Research Unit (Saint Augustinus Hospital, Antwerpen, Belgium) participated in this study. A detailed description of the tumor samples has been described previously (23) and can be found in the Supplementary Methods.
Chromatin immunoprecipitations
Chromatin immunoprecipitations (ChIP) were performed as described before (25) with minor adjustments. Firstly, as input material, tumor samples were cryosectioned (30× 30 μm sections) before further processing for ChIP-seq. For each ChIP, 10 μg of antibody was used, and 100 μL of Protein A magnetic beads (Invitrogen). ERα (SC-543; Santa Cruz), H3K4me3 (ab8580; Abcam), and H3K27me3 (07–449; Millipore) were used as antibodies. Primer sequences for quantitative PCR (qPCR) analyses are in Supplementary Table S5.
Solexa sequencing and enrichment analysis
ChIP DNA was amplified as described (25). Sequences were generated by the Illumina Hiseq 2000 genome analyzer (using 50 bp reads), and aligned to the Human Reference Genome (assembly hg19, February 2009). Enriched regions of the genome were identified by comparing the ChIP samples to mixed input using the MACS peak caller (26) version 1.3.7.1. Details on the number of reads obtained, the percentage of reads aligned, and the number of peaks called can be found in Supplementary Table S3. Bioinformatic analyses are described in the Supplementary Methods.
RNA isolation and mRNA expression analysis
Total tumor RNA was isolated as described previously (27). mRNA quality was assessed by quantitative real-time PCR (qRT-PCR) and bioanalyzer. Amplification, labeling, and hybridization of samples to 44k mRNA oligonucleotide-arrays (Agilent Technologies) were performed as described (27). Samples were hybridized against a breast cancer reference pool consisting of RNA from 105 primary breast tumors. A detailed description of mRNA expression analyses and progression-free survival analyses can be found in the Supplementary Methods.
Data access
All genomic data are deposited at the NCBI GEO, with accession numbers GSE40867 (ChIP-seq data) and GSE41994 (expression data).
Results
Genome-wide binding patterns of ERα, H3K4me3, and H3K27me3 in primary breast tumor specimens
Fresh frozen primary tumor specimens from a cohort of 84 patients with breast cancer were tested. These patients received aromatase inhibitor treatment of metastatic disease, i.e., anastrozole, letrozole, or exemestane (Supplementary Table S2). TTP was taken as an endpoint (Fig. 1A). Poor outcome patients were defined as patients with a TTP < 12 months, whereas good outcome patients were defined as patients with a TTP > 24 months (Supplementary Tables S1 and S2). To determine the differential patterns in the ERα/chromatin–binding landscape and the epigenetic modifications H3K4me3 and H3K27me3 as well as their possible correlations with outcome after aromatase inhibitor treatment, five good outcome tumor samples and seven poor outcome tumors were randomly selected as “discovery set.” Remaining samples from the entire cohort were used as “validation sets,” as will be discussed later. The tumor cell percentage was consistently high (>70%), and all tumors were ERα and progesterone receptor (PR) positive, whereas negative for HER2. Other clinicopathologic parameters are shown in Supplementary Table S2.
For all 12 tumors from the discovery set, fresh frozen samples were cryosectioned, and chromatin was isolated for ChIP. For each tumor sample, the ChIP for ERα, H3K4me3, and H3K27me3 was performed and isolated DNA fragments were analyzed by high-throughput sequencing (ChIP-seq; ref. 25). The number of aligned reads, unique reads, and number of peaks for each ChIP-seq sample are shown in Supplementary Table S3. Clear and distinct peaks were observed for each condition, as exemplified for one tumor sample (good outcome tumor #5; Fig. 1B). In this example specimen, 13,575 binding events for ERα, 19,012 binding events for H3K4me3, and 33,661 binding events for H3K27me3 were found. Although the vast majority of binding events found from this tissue sample were unique among the markers, overlap was found between ERα and H3K4me3 as well as between H3K4me3 and H3K27me3 (Fig. 1C; for all tumor samples, see Supplementary Fig. S1). Even though the total numbers of binding events greatly varied among tumors, these relative distributions of the three ChIP conditions were consistently found for all tumor samples tested (Supplementary Fig. S1). For the exemplified tumor sample, motif analysis was performed for each ChIP condition, Fig. 1D. As expected, ERα–binding events were enriched for ESR1 motifs, but also for its designated pioneer factors FOXA1 (7) and TFAP2 (13). For H3K4me3 ChIP on this tumor sample, enriched motifs were found for the promoter-selective transcription factors GTF2I (28), ZIC1 (29), and E2F1:TFDP2 (30), whereas only motifs for the mitochondrial transcription factor A (TFAM) were observed for H3K27me3-bound regions.
For all tested tumor samples, example binding events for ERα, H3K4me3, and H3K27me3 are shown (Fig. 1E), illustrating clear and high-quality data for all samples. For one tumor sample, no data could be generated for the H3K4me3 ChIP (poor outcome sample #2). Most ERα–binding events were found at distal enhancers and introns (Fig. 1F), consistent with previous reports in cell lines (7) and breast tumors (12). H3K4me3 was more markedly enriched at promoters, 3′-UTRs and exons, this in contrast to H3K27me3. These distributions were consistent in all analyzed tumors studied. The binding site distributions of all peaks related to the most proximal gene were found not to differ between good and poor outcome tumors.
Distinct genome-wide binding patterns of ERα, H3K4me3, and H3K27me3 correlate with patient survival after aromatase inhibitor treatment
Next, we aimed to determine whether the chromatin-binding patterns of ERα-, H3K4me3-, and H3K27me3–binding events would deviate between patients with a good versus a poor outcome upon aromatase inhibitor treatment. This was achieved through differential binding analysis (DBA; refs. 12, 31), directly comparing good outcome and poor outcome tumor ChIP-seq data. DBA normalizes the sequencing data for each run over the effective library size (reads in peaks) after input background subtraction, and determines relative enrichment of raw sequence reads over two distinct subgroups of tumor samples, in this case good versus poor outcome. Called peaks that were found in at least two tumor samples were considered for DBA to minimize noise, resulting in 14,232 peaks for ERα, 22,587 peaks for H3K4me3, and 35,602 peaks for H3K27me3 (Supplementary Fig. S2). On the basis of these peaks, the raw read counts for all tumor samples were checked for differential binding intensities between the two patient subgroups, resulting in lists of 222 (ERα), 66 (H3K4me3), and 351 (H3K27me3) peaks, with a false discovery rate of < 0.1 (Fig. 2A). Tumors from patients with the same outcome clustered together for ERα and H3K27me3 and, although less pronounced, for the H3K4me3 signals (Fig. 2B). This class separability was also observed in a principal component analysis (Supplementary Fig. S3). The read count and peak caller score did not bias patient stratification (Supplementary Fig. S4).
Peak regions that were significantly differentially enriched between the poor and good outcome patients are shown in a heatmap for all three conditions (Fig. 2C), with a clear distinction for ERα and H3K27me3, whereas case mixing was observed for H3K4me3 signals. The number of poor outcome peaks greatly outweighs the number of good outcome events for both histone marks, but not for ERα (Fig. 2D). A list of all differentially bound regions between good and poor outcome for ERα, H3K4me3, and H3K27me3 is shown in Supplementary Table S4.
The binding events that were differentially enriched between the two patient subgroups were analyzed for their localization relative to the most proximal genes (Fig. 2E). No clear differences in genomic locations of the ERα peaks were observed between good and poor outcome patients. H3K4me3 poor outcome peaks were found enriched at 3′-UTR regions and downstream thereof, whereas the 3′-UTR signals for H3K27me3-binding events were selectively lost in poor outcome tumor samples.
The ERα-, H3K4me3-, and H3K27me3–binding events were analyzed for DNA motifs differentially enriched in the good and poor outcome patient subgroups (Fig. 2F). Although only ESR1 motifs were found in good outcome patients, TFAP2A and TCF4 motifs were enriched next to ESR1 for the poor outcome patients. For H3K27me3, good outcome–enriched binding events were enriched for TERF1 and PPARG:RXRA motifs, whereas for poor outcome–enriched binding events, motifs were found for ERG and SIX1/SIX3. No distinct enriched motifs were found for the good outcome–enriched binding event for H3K4me3, but the poor outcome sites were statistically enriched for LHX2/GBX2, ARID5B, PAX7, and HOXB8 motifs.
Mutual exclusivity of differentially enriched chromatin-binding events between the patient subgroups
A large variation was found in the number (Supplementary Table S2) and overlap of chromatin-binding events between samples as exemplified for three specimens (Fig. 3A; for all samples Supplementary Fig. S2). This is consistent with a previous report that studied ERα ChIP-seq on breast tumor specimens (12). Because ERα/chromatin interactions are dictated by and have a facilitating effect on histone accessibility (11), we determined if altered ERα–binding patterns are accompanied by changes in the epigenetic profile of H3K4me3 and H3K27me3. When ERα–binding events observed in the poor outcome tumor sample were absent in the good outcome tumor, these alterations are not accompanied by a loss or gain of proximal histone marks, as exemplified in Fig. 3B. In line with these data, differences of either histone mark between good and poor outcome patients were virtually mutually exclusive and not shared with changed ERα–binding patterns (Fig. 3C). Consequently, a limited overlap was observed for the genes that were proximal to the altered ERα-, H3K4me3-, and H3K27me3–binding events (Fig. 3D). The altered binding events were mapped over all chromosomes and no clear bias towards distinct regions or chromosomes was observed (Fig. 3E).
Validation, integration with gene expression, and clinical outcome
Next, differential enriched binding patterns between good and poor outcome patients were validated in an independent group of 10 patients (4 good outcome and 6 poor outcome) using qPCR, Fig. 4A (for clinicopathologic parameters, see Supplementary Table S2). For each of the ChIP conditions, four to six primer pairs were designed, detecting randomly picked regions enriched in the good or poor outcome patients from the discovery set. ChIP efficiency can vary among samples due to tumor cell percentage, expression levels, and experimental variations, making inter-sample normalization an essential step in the DBA for the ChIP-seq pipeline. To implement inter-sample normalization in qPCR, the ratios of average good outcome over poor outcome–enriched binding site intensities were calculated. For ERα and H3K27me3, qPCR ratio separated poor outcome patients from good outcome patients, whereas differential enrichment could not be confirmed for H3K4me3.
After qPCR validation, all ChIP-seq–identified differential binding patterns were coupled to genes, based on proximity (gene body plus 20 kb upstream from the transcription start-site), enriching for ERα–binding sites involved in gene regulation (8). For the good outcome binding sites, this resulted in 84 (ERα), 19 (H3K4me3), and 29 (H3K27me3) genes. For poor outcome sites, 99 (ERα), 22 (H3K4me3), and 158 (H3K27me3) genes were found (Supplementary Table S6).
Gene Ontology (GO) analyses were performed on the proximal gene sets to identify plausible functional regulatory networks (Supplementary Table S7). Genes between the different classifiers only marginally overlapped (Fig. 3D), whereas pathways between ERα, H3K4me3, and H3K27me3-based classifiers showed no overlap. The ERα–based classifier was enriched for metabolic processes, whereas nucleoside transport, chemotaxis, and angiogenesis were enriched GO terms for H3K4me3. For the H3K27me3 classifier, developmental processes were strongly enriched. Comparing our data with pathways identified in mutational analysis of aromatase inhibitor–treated patients with breast cancer (4) showed shared processes involved in cell adhesion, cell cycle, chemotaxis, developmental processes, immune responses, metabolism, signal transduction, and transcriptional regulation (Supplementary Table S8).
Next, array-based mRNA expression levels of all proximal genes for ERα, H3K4me3, and H3K27me3 sites were tested for a correlation with time to progression (TTP) in a second independent set of 72 tumors (validation set), of aromatase inhibitor–treated patients with metastatic disease (Fig. 4B, clinicopathologic parameters Supplementary Table S2). Expression of proximal genes for differential ERα (183 genes), H3K4me3 (41 genes), and H3K27me3 (187 genes) binding was used to classify tumors as “poor” and “good outcome.” Classifications for altered ERα (P = 0.0016; HR, 2.42; 95% CI, 1.40–4.19) and H3K27me3 (P = 0.0001; HR,3.19; 95% CI, 1.83–5.57) binding events correlated with a differential TTP, in contrast to H3K4me3 (P = 0.386; HR, 1.38; 95% CI, 0.76–2.49).
Previously reported classifiers, PAM50 (P = 0.0093; HR, 2.11; 95% CI, 1.20–3.71; ref. 32), Oncotype DX (P = 0.0256; HR, 1.87; 95% CI, 1.08–3.23; ref. 33), and the TAM78 Rotterdam classifier (P = 0.0151; HR, 1.98; 95% CI, 1.14–3.42; ref. 34), classified patients in this cohort, performing equally well as the ERα–proximal gene classifier (Fig. 4E). The H3K27me3-proximal gene classifier showed a higher HR, even though the 95% CI did overlap (Fig. 4E). ERα and H3K27me3 ChIP-seq–based classifiers remained significant after multivariate correction analyses (Supplementary Table S9). The classifier genes are listed in Supplementary Table S6. As a second expression-based validation, a cohort of neoadjuvant letrozole treatment patients (N = 54) was used, analyzing samples before and after treatment, where patients were stratified in responsive and non-responsive groups (Fig. 4C; refs. 5, 6). PAM50, Oncotype DX, and TAM78 as well as our ChIP-seq classifiers successfully identified patients with differential outcome.
ERα can affect gene regulation by long-range genomic chromatin-loop interactions, as shown by ERα ChIA-PET (8). Therefore, these published long-range interactions were also considered in our study. In addition, binding sites were analyzed for any transcription factor enrichment, DNAse hypersensitivity, and H3K27Ac, representing active enhancers (35), using ENCODE datasets (Supplementary Table S10). Genes were selected that either had a binding event at a promoter region, or chromatin looping toward the promoter was found, generating three separate lists of 92 (ERα), 38 (H3K4me3), and 150 (H3K27me3) genes. Overlap between the “proximity-based” and “ChIA-PET–based” gene sets was high for H3K4me3 (76%) and H3K27me3 (78%), but considerably lower for ERα (49%; Supplementary Fig. S5). Analogous to our initial classifiers, ChIA-PET–based classifiers for ERα (P = 0.0041; HR, 2.45; 95% CI, 1.41–4.26) and H3K27me3 (P = 0.0002; HR, 3.02; 95% CI, 1.73–5.27) could classify patients, in contrast to H3K4me3 (P = 0.606; HR, 1.30; 95% CI, 0.71–2.37; Supplementary Fig. S6A and Supplementary Table S9). Also for this “ChIA-PET–based” classifier, samples from the neoadjuvant treated letrozole cohort were analyzed (6), identifying patients with differential response (Supplementary Fig. S6B). Conserved genes between both gene classifiers performed equally as the separate classifications (ERα: P = 0.05; HR, 1.96, 95% CI, 1.10–3.49; H3K27me3: P = 0.0001; HR, 3.68, CI 95%, 2.09–6.50; H3K4me3: P = 0.3944; HR, 1.30, 95% CI, 0.711–2.37; Supplementary Fig. S6C and S6D).
A third cohort of patients was analyzed, receiving first-line tamoxifen for metastatic disease instead of aromatase inhibitors (Fig. 4D; ref. 27). The metastatic setup of this cohort enabled a direct comparison between the patient sets. No significant correlation with TTP after tamoxifen was found for the H3K4me3 and H3K27me3 ChIP-seq–based classifiers (H3K27me3: P = 0.0838, HR, 1.53, 95% CI, 0.95–2.47; and H3K4me3: P = 0.0602, HR, 1.556, 95% CI, 0.98–2.47). However, significance was reached for the ERα–based classifier (P = 0.0426, HR, 1.63, 95% CI, 1.02–2.63). These data indicate that while the ERα ChIP-seq–based classifier may identify breast cancer patients with a poor outcome, irrespective of the type of endocrine treatment, the H3K27me3 classification suggests aromatase inhibitor–treatment selectivity in patients with metastatic breast cancer (Fig. 4E).
Discussion
Estrogen receptor genomics has greatly increased our knowledge of hormone receptor functioning in breast cancer cell lines and tumors. The number of ERα/chromatin–binding events greatly exceeds the amount of estradiol-responsive genes (11, 36), suggesting a high level of complexity in ERα–mediated gene regulation (8). ERα/chromatin–binding events have been found to correlate with survival of patients with breast cancer (12), but survival in this heterogeneous study also correlated with traditional pathologic parameters, including PR and HER2, potentially hampering further clinical interpretation of the data. Therefore, we selected ERα+, PR+, and HER2− tumors from a homogeneous cohort of breast cancer patients whose metastasis were all treated with aromatase inhibitors.
Altered ERα–binding patterns between patient groups were not accompanied with altered epigenetic profiles of H3K4me3 and H3K27me3. Our data suggest that any dynamic behavior of ERα uses the accessible regions in the genome that are imprinted and readily accessible in a static epigenetic landscape. Still, H3K27me3-binding events enable the identification of patients with a poor outcome after aromatase inhibitor treatment, in contrast to H3K4me3. The presented data suggest that aromatase inhibitor resistance is accompanied by a specific gain of polycomb-mediated gene repression at distinct sites. Moreover, it suggests an ERα–independent mechanism of therapy resistance in ERα–positive tumors. Our data could be validated in two other cohorts of breast cancer patients, using two different technologic approaches, namely ChIP-qPCR and gene expression analysis.
Because genomic patterns of ERα and H3K27me3 were indicative for patient survival, the enriched motifs may provide clues for transcription factors involved in treatment outcome. For ERα, TFAP2 and TCF4 motifs were selectively enriched for the poor outcome patients. AP-2 can directly guide ERα/chromatin interactions (13), promotes breast cancer cell proliferation (37), and correlates with poor outcome (38). TCF4 enhances breast cancer cell invasion (39) and binds ERα (40), providing a level of cross-control between estradiol and wnt pathways (40).
For H3K27me3, “good outcome” sites were enriched for TERF1 and PPARG:RXRA motifs. TERF1 is a component of the telomere nucleoprotein complex. SNPs in TERF1 have been tested for breast cancer susceptibility and prognosis, but no correlations were found (41). PPARG:RXRA ligands can trigger breast cancer cell apoptosis (42) and PPARγ activation blocks breast cancer cell invasion (43) and induces terminal cell differentiation (44). For poor outcome–enriched H3K27me3 regions, motifs were found for ERG and SIX1/SIX3. ERG and ERα mutually repress each other's activities (45). SIX1 is expressed in breast cancer, stimulating tumor cell proliferation (46), inducing genomic instability and malignant transformation (47), correlating with poor prognosis (48). Collectively, these data highlight possible transcriptional mechanisms that may form the basis for aromatase inhibitor response.
The ERα ChIP-seq classifier was applicable for aromatase inhibitor- and tamoxifen–treated patients with breast cancer with metastatic disease. Because aromatase inhibitors and tamoxifen both affect the functionality of ERα, the genomic downstream signatures for treatment outcome could overlap as well. Because ERα is targeted by these endocrine agents, the ChIP-seq approach could aid in removing noise from expression analyses to exclusively monitor genes that are directly affected by this hormone receptor.
Testing the aromatase inhibitor ChIP-seq classifiers in tamoxifen-treated tumors showed no significant difference in TTP for H3K27me3. Aromatase inhibitor–treated patients cohorts are relatively rare, and our cohort is on metastatic disease. To compare between endocrine treatments, data from a metastatic cohort of tamoxifen-treated patients were used. Performing these ERα ChIP-seq experiments and epigenetic assessments in the adjuvant setting would enable a direct comparison with any of the large (adjuvant-treated) breast cancer patient cohorts (34, 49, 50) for more extensive in silico validations.
Tumor-intrinsic plasticity of ERα and H3K27me3 can be a hallmark of endocrine therapy resistance in breast cancer and may ultimately be applicable to guide endocrine treatment selection for patients with breast cancer.
Disclosure of Potential Conflicts of Interest
I. Simon is employed as Senior Director, R&D in Agendia and has commercial research grant from the same. No potential conflicts of interest were disclosed by the other authors.
Authors' Contributions
Conception and design: M.P.H.M. Jansen, L. Dirix, E.M.J.J. Berns, W. Zwart
Development of methodology: M.P.H.M. Jansen, W. Zwart
Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): M.P.H.M. Jansen, E.A. Reijm, I. Simon, R. Kerkhoven, A. Velds, S. van Laere, L. Dirix, X. Alexi, J.A. Foekens, S.C. Linn, W. Zwart
Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): M.P.H.M. Jansen, T. Knijnenburg, M. Droog, A. Velds, L. Dirix, X. Alexi, L. Wessels, S.C. Linn, E.M.J.J. Berns, W. Zwart
Writing, review, and/or revision of the manuscript: M.P.H.M. Jansen, T. Knijnenburg, E.A. Reijm, S. van Laere, L. Dirix, X. Alexi, L. Wessels, S.C. Linn, E.M.J.J. Berns, W. Zwart
Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): M.P.H.M. Jansen, I. Simon, M. Droog, W. Zwart
Study supervision: M.P.H.M. Jansen, L. Dirix, E.M.J.J. Berns, W. Zwart
Acknowledgments
The authors thank Shan Baban and Marja Nieuwland for sample processing, Kirsten Ruigrok-Ritstier, Paul Roepman, Anieta Sieuwerts, Vanja de Weerd, Anne van Galen, Maxime Look, Marion Meijer-van Gelder, Marleen Kok, Karin Beelen, Sander Canisius, and Hilde Wuyts for their contribution and technical support, Michael Hauptmann for help with statistics, and Jason Carroll (Cambridge, UK) for critically reading the article and valuable suggestions. The authors also thank the surgeons, pathologists, and medical oncologists of the St. Clara Hospital, Ikazia Hospital, St. Franciscus Gasthuis (Rotterdam), and Ruwaard van Putten Hospital (Spijkenisse) for the supply of tumor tissues and/or clinical follow-up data collection.
Grant Support
This work was supported by grants from the KWF Dutch Cancer Society, A Sisters Hope and Top Institute Pharma, grants T3-108 and T3-502.
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.