Primary Colorectal Cancers and Their Subsequent Hepatic Metastases Are Genetically Different: Implications for Selection of Patients for Targeted Treatment

This article is featured in Highlights of This Issue, p. 593

The main challenge for the future of cancer treatment is to provide every individual patient with the most effective drug tailored to their specific cancer. Personalized cancer treatment is still in the very early stages, but several recent studies have shown the potential of this approach. For example, trastuzumab (Genentech) is approved for patients with breast and gastric cancer who have HER2-amplified tumors (1). For patients with metastatic colorectal cancer (mCRC), the anti-epidermal growth factor receptor (EGFR)-antibodies panitumumab (Amgen) and cetuximab (Merck) are approved for patients with wild-type KRAS tumors (2; 3); patients with mCRC with a mutated KRAS gene do not benefit from treatment with these antibodies (4–8). Unfortunately, treating patients with anti-EGFR antibodies, while they express (wt) KRAS remains beneficial to only a subset of patients. In addition, it is becoming increasingly clear that patients with loss-of-function mutations in PTEN or BRAF are also unresponsive to cetuximab (2; 3). These examples illustrate that solitary gene analyses in personalized anti-cancer treatment have limitations (9; 10) and emphasize the need for more complete genetic profiling of relevant pathways within the tumor to optimize treatment strategies (11). In addition, it is important to determine which tumor tissue will best predict treatment outcome. Current clinical practice is to use archived material of the primary tumor to determine the constitution of the molecular target to select patients for treatment (2–6). However, there are several biologic reasons why this may not be optimal. First, genomic instability is a hallmark of cancer and caused by constant selection pressure tumors rapidly change their genetic makeup over time. Second, specific populations of tumor cells may be more prone to metastasis than others, which is likely to result in an enrichment of these cells and consequently their genetic aberrations in the metastases. Third, systemic treatment may induce selection pressure toward a specific genetic phenotype or induce additional genetic changes. Taken together, tumors are genetically dynamic, which suggests that selecting patients for targeted treatments based on the characteristics of the primary tumor and not their metastases may not be optimal. It has been suggested that cancer-initiating mutations will be present in every tumor localization but that cancer-driving mutations may be enriched or depleted in the metastasis (12). Sequencing 189 candidate cancer genes (CAN) in breast and colorectal cancer showed substantial differences within tumors, indicating that each tumor type is heterogeneous and that tumorigenesis is likely tumor-specific (13). Importantly, these analyses were conducted on primary tumors and did not address the question of genetic differences between primary tumors and their metastases (14).

Next-generation sequencing (NGS) is an extremely powerful technology for genetic analysis of complete signaling pathways in large patient cohorts (11). We therefore embarked on a study to comprehensively compare the genetic constitution of the 1,264 genes comprising the most relevant cancer-associated pathways and processes in 21 primary colorectal adenocarcinomas (CRC) and their subsequent hepatic metastases by using targeted deep-sequencing.

We selected 21 patients with CRCs, from whom both formalin-fixed, paraffin-embedded (FFPE) samples of the primary tumor and their sequential hepatic metastasis were available (Table 1). We included only patients with a minimal time of 6 months between primary tumor resection and metastasectomy. Patients who received no treatment (n = 10, chemo-naïve group) and patients who received chemotherapy (n = 11, chemo-treated group) between removal of the primary tumor and metastasis were included. Eight patients in the chemo-treated group received 5-fluorouracil, oxaliplatin, and leucovorin as chemotherapy and 3 patients received capecitabine, oxaliplatin, and bevacizumab between the surgical interventions. The primary CRC tumor of the majority of patients were localized in the colon (n = 16), where 5 patients were diagnosed with rectal adenocarcinoma. At time of primary CRC resection, patients were classified according to the TNM-staging system (Table 1), showing that all patients were diagnosed with stage II, III, or IV disease. In the chemo-treated group, 2 patients were diagnosed with synchronous hepatic metastasis and the other 8 patients with metachronous hepatic metastasis. All procedures were approved by the University Medical Center Utrecht ethics committee.

A pathologist confirmed the CRC diagnosis of each sample and demarcated tumor areas from normal tissue. To obtain samples consisting of 80% or more tumor cells, we microdissected tumor tissue using a laser-capture microdissection microscope (PALM, Carl Zeiss). Dissected tumor cells were incubated with sodium thiocyanate to dissolve DNA and protein cross-linking. Subsequently, samples were treated overnight at 55°C with 20 mg/mL Proteinase-K, and DNA was isolated using the QIAamp DNA Micro Kit (Qiagen). At least 2 μg DNA was isolated per sample, of which more than 1 μg was stored at −20°C for validation purposes and 1 μg was used for library preparation.

Sequencing libraries were prepared as described previously (15), with minor modifications. Briefly, 1 μg of tumor DNA was sheared using a Covaris S2 sonicator (Covaris; duty cycle, 20%; intensity, 5,200 cycles/burst for 10 minutes) and end-repaired using End-It DNA end-repair kit (Epicenter Technologies). After ligation of short sequencing adaptors, the library was amplified using a truncated version of sequencing primers, and size was selected on 2% agarose gel for a 125- to 175-bp fraction.

To interrogate differences in genetic makeup between primary tumor and metastasis with possible therapeutic consequences, we composed a Cancer Mini-Genome consisting of all the exons (ensembl v56) of 1,264 genes (totaling ∼7 Mbp), including known oncogenes, tumor suppressor genes, identified colorectal CAN genes (13), all 518 kinases (16), and important pathways related to tumorigenesis and anti-cancer treatment (e.g., angiogenesis; apoptosis and EGFR, PIK3CA, TGFβ, mTOR, and VEGF pathways were all included; ref. 17). All genes of the Cancer Mini-Genome are listed in Supplementary Table S1. A custom-made PERL script was used to design 60-mer capture probes for the target regions with a 15-bp moving window on both genomic strands. The best probes in each window was selected on the basis of T_m, GC%, and monomer stretches as described previously (18). Probes with more than one additional location in the genome with a similarity more than 60% were discarded and the remaining probes were ordered on a custom 1M CGH array (Agilent). A list of sequences of the capture probes used for the Cancer Mini-Genome is mentioned in Supplementary Table S2.

Enrichment was conducted as described previously (18). Briefly, size-selected libraries were PCR amplified using 10 additional cycles to produce amounts necessary for hybridization and subsequently purified and hybridized to CGH-capture arrays to enrich for exonic regions of our designed Cancer Mini-Genome. After 72 hours of hybridization at 42°C, arrays were extensively washed and specifically hybridized library molecules were eluted and amplified by 13 cycles of PCR. After purification, barcodes (SOLiD barcode primers; Applied Biosystems) unique to each sample were incorporated into each library by carrying out 4 additional PCR cycles. Ten pre-barcoded individual libraries were pooled and sequenced on the SOLiD 3+ system (Life Technologies) according to the manufacturer's instructions.

Sequence data were mapped against the human reference genome (GRCh37/hg19) using BWA (19). Single nucleotide variations (SNV) were called with a custom PERL script as previously described (15) with the settings mentioned below and the results were further processed by custom-made PERL scripts, with subsequent prediction of the consequence of genetic variations at the amino acid level. Bioinformatic analyses were conducted to assess frequency of variants, affected pathways, and functional relevance of identified somatic variations. Normal nonneoplastic tissue from 5 patients–-of whom also the primary CRC tumor and hepatic metastasis—was sequenced to annotate germ line variants and to exclude common existing SNVs in accordance with previous studies (20, 21). Identified variations in the Cancer Mini-Genome were annotated according to existing databases (ensembl59). To limit the influence of SNV caller artifacts, we applied 3 different SNV callers [custom PERL script (refs. 15, 18), VarScan (ref. 22), and the GATK toolkit (ref. 23)] and considered the overlapping variants as high-confidence variants. The following settings were applied:

1. Custom PERL filter: coverage between 20 and 2,000×; 3 reads supporting each allele on each strand, minimum of 3 independent reads per variant, no strand imbalance greater than 1:10, call quality ≥ 10, reads should map uniquely, clonality filtering: no more than 5 identical reads counted per allele, alleles counted if present ≥ 3×.

2. VarScan: Variants with a P ≤ 10⁻⁸, minimal coverage of 20×, and 2 non-reference reads and a variant frequency ≥ 0.20.

3. GATK toolkit: Minimal 20× coverage, variant frequency ≥ 0.20.

In general, the overlap was 33% ± 10% of all the variable positions called by the SNV callers, indicating significant bias depending on the used SNV caller. In cases where a difference between the tumors was discovered, we examined the raw data of the primary sample for any indication of the variant that might have escaped the strict SNV caller filters. If the allele was detected only minimally (presence of 3 non-reference reads), this variant was not called as de novo but as equal between tumors. To predict the possible impact of the identified variation on amino acid structure and thus the protein function, PolyPhen-2 [Polymorphism Phenotyping version2 (ref. 24)] was used.

To validate NGS-identified variants with conventional Sanger sequencing, genomic DNA of all samples was whole genome–amplified (WGA) using the Repli-g FFPE Kit (Qiagen) with 100 ng of input DNA according to the company's protocol. Primers were designed to amplify an about 200 bp amplicon containing a random number of putative variants and amplicons were resequenced using conventional Sanger sequencing according to standard protocols with WGA-amplified DNA as input. In addition, the samples were analyzed through routine diagnostics on KRAS (HRM and Sanger sequencing in duplo).

Basic patient characteristics are shown in Table 1 and a schematic representation of followed sequencing procedure is presented in Fig. 1. The number of sequenced nucleotides of all included samples was greater than 100 GB. The mean percentage of mapped reads was 62%, of which 65% on average was on target for our Cancer Mini-Genome. Per sample, an average of about 789 Mbp (range, 250—2,150) of the target genes were sequenced, resulting in an average mean and median coverage of 163 and 87, respectively (Table 2). In the exons of our Cancer Mini-Genome of all 42 tumor samples combined, 8,405 high-confidence genetic variations were detected when compared with the reference genome (hg37), including 8,001 SNVs (96%) and 404 (4%) small insertions or deletions (Fig. 2). Of these identified SNVs, 84% (N = 6,696) were known variants when compared with ensembl59 (1,174 genes) and 1,305 were novel (667 genes). The complete catalog of all scored variations is shown in Supplementary Table S3. To obtain an estimate of the inherited background of the genetic variations, we sequenced “normal” nonneoplastic colorectal tissue of 5 patients of whom also the paired primary CRC tumor and hepatic metastasis was analyzed. A total of 192 variants (2%) were found to be common in all these samples and were therefore excluded from the comparison. To further validate the NGS results, we assigned a random subset of candidate genetic variations for conventional sequencing techniques. In total, 108 variations in 32 genes were validated with a true-positive rate of 84%. Supplementary Table S4 lists all re-validated variants using Sanger sequencing and the results of the diagnostic sequencing of KRAS.

Using stringent quality control, we filtered nonrelevant genetic variations and subsequently compared the primary tumor and its metastasis to determine differences in somatic mutational constitution. On an individual basis, on average 83 (SD = 69) variants were gained in the metastasis and 70 (SD = 48) variants were lost as presented in Fig. 3. Next, we investigated the mutational status of a selected set of genes including KRAS, HRAS, NRAS, EGFR, PI3KCA, FLT1, KDR, PTEN, and BRAF, which are well studied in relation to therapeutic responsiveness (7, 8, 25–30). Table 3 shows a brief overview of the current clinically relevant codons of BRAF, EGFR, HRAS, KRAS, NRAS, and PIK3CA. In almost all 21 patients with CRC, we found aberrations in KRAS and EGFR (both ≥90%) either in the primary tumor and/or in liver metastasis as shown in Table 4, which includes an overview of variants per chromosome position for each gene. We identified genetic variants in KRAS codon 12 and 13 and in codon 61 of NRAS, which are currently biomarkers of resistance to anti-EGFR therapy. No variants were detected in BRAF codon 600, EGFR codon 790 and 858, HRAS codon 12, 13, and 61, KRAS codons 61 and 146, NRAS codon 12 and 13, and PIK3CA codons 542, 545, and 1,047. In the majority of patients (67% and 52%, respectively), PI3KCA and FLT1 genes were affected; in a smaller but substantial part of our patient population, HRAS, NRAS, KDR, PTEN, and BRAF were mutated (10%, 24%, 19%, 24%, and 38%, respectively). Subsequently, mutational status differences per individual gene between primary tumor and liver metastasis were investigated. Taking all variations into account, dissimilarities of the KRAS and EGFR mutational status between both tumor entities were detected in 52% and 86% of the 21 patients with CRC, respectively. Modest variability was observed for HRAS (24%), PIK3CA (19%), FLT1 (10%), NRAS (10%), and BRAF (14%). KDR and PTEN showed a more or less stable pattern of mutational status between both tumor identities, with only 5% of patients presenting deviations. Supplementary Table S5 presents on overview of all variant information for these genes per patient. To investigate the impact of chemotherapy on the mutational status of the metastasis, patients with CRC were categorized into 2 equally matched groups of chemotherapy-naive (N = 10) and chemotherapy-treated (N = 11) patients. In this relatively small data set, we could not show a significantly increased number of variants as a result of chemotherapy. Accordingly, a significant relation with genetic burden and time between primary tumor and metastasis resection was not observed.

Prediction of possible impact on protein function of all nonsynonymous variations (1,455) identified 41 nonsense, splice site variants, and indels (small insertions/deletions). Of all other 1,414 nonsynonymous variants, 281 were predicted as probably damaging (P > 0.85, using PolyPhen-2). Moreover, 170 were classified as possibly damaging (0.20 >> P << 0.85), 740 as benign (P < 0.2), and 261 not classified at all. The entire catalog of all genes potentially affected at the protein level is listed in Supplementary Table S3. Novel potentially protein-changing mutations were identified throughout our Cancer Mini-Genome which contains genes that are relevant for a variety of cellular functions and responses (i.e., mTOR, TGFβ pathways, tyrosine kinases, apoptosis, etc.; Table 5).

Given the current significance of KRAS status for treatment decision making, we specifically addressed the EGFR pathway (25). First, we analyzed the variants in codons 12 and 13 by conventional sequencing and found mutations in either codon 12 or 13 in 6 of 21 patients (Supplementary Table S4). One patient converted to wild-type KRAS in the metastasis. These data were identical to the data generated by our NGS platform. Next, we analyzed individual components of the EGFR pathway. None of the 21 patients with CRC had a predicted protein-changing EGFR mutation in the primary tumor or liver metastasis. Ten of the 21 patients with CRC had a relevant genetic variation in BRAF or KRAS relevant for possible therapeutic responsiveness. Interestingly, every patient showed one or more protein-altering genetic aberration in genes downstream of EGFR. These genes included MAP3K1 (MEKK1), MAPK10 (JNK3), PLCG1 (PLCγ), PIK3CG (PI3K), PRKCB (PKC), and RAF1 (Table 5). Furthermore, we also found KSR1, MOS, NF1, and RAPGEF-3/4—all potentially relevant variants in other genes strongly associated to the EGFR pathway.

The purpose of the study was to determine which tissue best reflects the presence of a target for therapy: the primary CRC or their metastasis. The most important finding was that substantial genetic differences existed between the primary tumor and its metastasis. Significant numbers of either loss or gain of functionally relevant genetic variations were found, suggesting that metastases may provide a better predictive window for targeted therapy than the primary cancer.

In the era of targeted therapy and personalized cancer treatment, individual genetic mutations in the primary tumor drives patient selection and therapeutic strategies (7, 8, 31). If we focus on the codon 12 and 13 KRAS mutations commonly used to determine whether patients with metastatic CRC will benefit from treatment with antibodies against EGFR, our data are in line with recent reports where a concordance of more than 95% for the 12G mutation was found between the primary tumor and metastasis (27, 29, 32). However, we analyzed all exons of the KRAS gene and found a substantial number of additional genetic changes outside codons 12 and 13. Their clinical relevance has yet to be determined but our findings underline the importance of having a complete overview of the genetic alterations. However, the power of the NGS analysis used in this study is that it provides a complete overview of all known components of the targeted pathway. Using this technology, we showed that in addition to the well-known diagnostic KRAS codon 12 and 13 mutations, several variants that could potentially alter protein function were present in other components of the pathway. This may explain why the clinical benefit of anti-EGFR antibodies remains limited, even in preselected patients on the basis of codon 12/13 KRAS analysis. We also clearly showed a substantial loss or gain of variants from the primary tumor to its metastasis in our selected set of 1,264 genes. The finding that chemotherapy did not change the genetic profile significantly was rather surprising but can have several explanations. This study was not primarily designed to investigate chemotherapy-induced changes and therefore lacked sample power to draw any definitive conclusions.

Metastases are generally heterogeneous, and clonality may be a confounder in the analysis (33–36). In CRCs, the cancer-initiating genes of the primary tumor were investigated in the metastasis; this showed heterogeneity of identified variants in the primary tumor and its paired metastasis, which supports our findings (34; 37). Similarly, in lobular breast cancer, the mutational evolution between the metastasis and its primary origin was studied in a single patient, illuminating 19 novel nonsynonymous coding mutations of the 32 variations discovered in the metastasis (38). In a patient with basal-like breast cancer, the metastasis displayed two de novo mutations and an additional large deletion (39). These numbers are lower than observed in our analysis, which may be explained by the enrichment for specific genes and higher sequencing coverage in our study and the stringency of the analysis that allows the detection of lower frequency mutations.

Recently, the genomes of the primary origin and metastasis have been sequenced in patients with pancreatic cancer in two independent investigations (40, 41). Both studies showed that clonal expansions found in the metastasis primarily originated from the primary tumor, showing that clones evolve over time. Mutations in cancer-driving genes may occur later and help the tumors to become more malignant under a certain selection pressure. Equally, Goranova and colleagues showed that multipoint microsampling elucidated discrepancies in genetic variations in APC, KRAS, and TP53 between primary CRC tumors and its subsequent hepatic metastasis and also showed subpopulations for these genes in most individual tumors (42).

These studies show the potentially very important contribution of genetic heterogeneity. Future prospective studies including the analysis of multiple specimen of the same tumor may provide a better perspective of the relevance of genetic heterogeneity to clinical decision making. This holds true for both current diagnostic testing of tumors and the use of NGS. In our study, chemotherapy did not further increase the number of genetic variations; however, the number of patients analyzed is too small to draw any firm conclusions. Furthermore, the selective pressure does not necessarily act through quantity of mutations but rather their quality. In the present study, we went to great lengths to make sure that the quality of the sequence data generated was high. Samples were sequenced until we reached a predetermined level of 40× mean coverage, which is generally accepted in the field to allow reliable heterozygote variant detection. Using this strategy, we were able to cover our intended target regions for more than 83% with sufficient coverage to allow for confident variant calling.

Our study does reveal some future challenges. First, differentiation between low- and high-frequency variations was not our primary aim but could potentially identify “initiating” and “driving” variants (12, 43). Complicating factors as tumor heterogeneity, genome or gene copy number status, and quantitative sequencing need also to be addressed in further studies. Second, to acquire a basic idea of the fraction of germ line variants, we sequenced normal nonneoplastic tissue of 5 patients to determine an approximate level of germ line contribution. However, for clinical interventions of metastatic disease, only the differences in genetic makeup between primary CRC tumor and hepatic metastasis are therapeutically relevant. Third, for target-enriched NGS, we used genomic DNA from selected tumor-rich areas (>80%) from FFPE material. It is debated whether the FFPE material procedure introduces variation by the depurination or fragmentation. However, solid sequencing is not limited by the shorter fragment sizes in an FFPE DNA sample. Schweiger and colleagues compared FFPE tissue with snap-frozen tissue and showed identical results for copy number variants and SNVs (44). We compared the mutational spectra with non-FFPE material and did not observe significant difference as well (data not shown). However, as DNA yields are generally low and highly variable, systematic exploitation of large collections of archived material for retrospective studies will still be challenging. Fourth, functional confirmation of the novel identified variants is necessary to discern real pathogenic from benign variants and to validate bioinformatic prediction algorithms. Systematic collection of these large-scale data in open-source databases enables the possibility to address these issues. Over the next few years, enormous amounts of genetic data will likely be generated. Improved tools to analyze cellular pathways and network analyses will become readily available (45–51). Linking clinical outcome to these genetic data and enforcing open-source data sharing will help to accelerate the development of clinical decision-making algorithms (52, 53). Taken together, our study indicates that obtaining a biopsy from the metastasis and reevaluating the genetic makeup at the time of treatment may be preferred over an archived primary tumor. This requires a different mindset for oncologists because biopsies from metastases are not commonly taken at the start of treatment. Obviously, the clinical benefits need to outweigh the risk that patients take by undergoing a biopsy. However, in the metastatic setting, the complication rate of ultrasound-guided biopsies is 0% to 5%, depending on tumor location (54–56). Biopsies should be guided by stringent quality control, and items such as the percentage of tumor cells in the biopsy, tumor inflammation, and necrosis should be noted. If patients are allocated to clinical studies with targeted agents, an NGS readout is currently not sufficient and the target should be validated by other certified methods. Similarly, drug approval by the regulatory authorities will also require stringent quality control measures. Regarding therapeutic responsiveness, it is also crucial to investigate differences in genetic profiles between primary tumors and synchronous metastases. In conclusion, this study shows that the differences in potentially relevant genetic variations between the primary CRC and hepatic metastases in important pathways are of such magnitude that an impact on treatment outcome is realistic. This indicates that genetic pathway analysis of the metastasis may have more predictive power when patients are selected for specific treatment modalities, thus allowing for further refinement of treatment algorithms.

S.J. Scherer has stock ownership of Roche. E.E. Cuppen received unrestricted research funding from Applied Biosystems. E.E. Voest received unrestricted research funding from Roche. No potential conflicts of interest were disclosed by the other authors.

Literature search was done by J.S. Vermaat, I.J. Nijman, E.E. Cuppen, and E.E. Voest. Statistical analysis was conducted by J.S. Vermaat, I.J. Nijman, F.L. Gerritse, E.E. Cuppen, and E.E. Voest. J.S. Vermaat, I.J. Nijman, P.J. van Diest, R. van Hillegersberg, E.E. Cuppen, and E.E. Voest gave the study concept and design. Acquisition of data collection and interpretation and drafting and critical revision of the manuscript were done by all authors. Technical/material support was provided by M.J. Koudijs, F.L. Gerritse, W.M. Roessingh, M. Mokry, N. Lansu, and E. de Bruijn. Data analysis was done by J.S. Vermaat, I.J. Nijman, M.L. Koudijs, F.L. Gerritse, E.E. Cuppen, and E.E. Voest. E.E. Cuppen and E.E. Voest obtained unrestricted funding. Study was supervised by J.S. Vermaat, P.J. van Diest, R. van Hillegersberg, E.E. Cuppen, and E.E. Voest.

This study was supported by unrestricted grants from Roche and the Dutch “Barcode for Life” Foundation.

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.