Molecular Characterization of Locally Relapsed Head and Neck Cancer after Concomitant Chemoradiotherapy

Head and neck squamous cell carcinomas (HNSCC) are among the most common incident cancers worldwide, with more than 600,000 newly diagnosed cases each year (1, 2). HNSCC arises in the mucosal lining of the upper aerodigestive tract encompassing the oral cavity, oropharynx, hypopharynx, and larynx. Most tumors arise in the oral cavity (30%–40%), followed by the larynx (30%–35%; ref. 3). Classic risk factors for HNSCC are smoking and excessive alcohol consumption. Also, patients with inherited genetic predispositions such as Fanconi anemia have an increased risk for HNSCC (4). For tumors arising in the oropharynx (OP) infection with high-risk human papillomavirus (HPV), in particular HPV16, is a more recently emerged risk factor (5, 6).

Unfortunately, the majority (∼60%) of patients with HNSCC present with advanced stage disease, meaning that the malignancies have invaded neighboring anatomic structures and/or have spread to the lymph nodes in the neck. Clinical management of these patients consists of surgery followed by postoperative (chemo)radiotherapy, locoregional radiotherapy, or concurrent chemoradiotherapy with surgical salvage for residual or relapsed cancers. Upfront chemoradiotherapy is applied when surgical resection is considered too invasive, the tumor is deemed very sensitive to nonsurgical modalities or when severe problems with swallowing and speech are expected after surgery. This multimodality approach has led to improved outcomes in terms of quality of life in patients with advanced stage disease, but prognosis still leaves much to be desired (2, 7). The most critical event in the course of the disease is the development of local relapses (LR; refs. 7–10). After definitive treatment with concurrent chemoradiotherapy, the local failure rate at 3 years of follow-up is between 15% and 50% (7, 11, 12).

There are several pathobiological origins for LRs (13, 14). A local relapse is clinically defined as a lesion that develops within 3 years after and within 2 cm from the treated primary tumor, relapses developing spatially or temporally more distinct are considered to be second primary tumors (SPT; ref. 15). Our current knowledge on the pathobiological mechanisms of the development of relapses is mainly based on molecular studies of resected oral tumors and their surgical margins (13, 14, 16). Two mechanisms seem to play a role. First, even when the surgical margins are histologically diagnosed as tumor-free, tumor cells may have remained behind, also called minimal residual cancer (MRC), and these could give rise to an LR. A genetic relationship between tumor and LR is then expected. Second, fields of preneoplastic cells that preceded and surrounded the primary tumor might stay behind unnoticed after tumor excision, and develop in LRs. These LRs have also been depicted as “second field tumors” (SFT; refs. 13, 14, 16). When studied for genetic alterations, these second field tumors share the early changes with the index tumor, such as TP53 mutations, CDKN2A mutations and 9p losses, but not the late changes. Finally, patients with HNSCC are also at risk for second primary tumors, and as these often originate independently from the index tumor, no genetic relationship is expected.

To date, the genetic origin of local relapse after chemoradiotherapy has not been studied very well, and we wondered which mechanisms might play a role. As surgical margins are not available from these patients, the studies have to rely on genetic characterization of the index tumors and local relapses to investigate genetic relationships. Caveats in these analyses are both intratumor genetic heterogeneity, which can be determined by multiple biopsies or single-cell analyses, and treatment-induced genetic changes.

Recent research has demonstrated a role for intratumoral genetic heterogeneity in HNSCC development and progression (17–19), which may indeed contribute to LR (20). Therefore, genetic assays to study genetic relationships that are validated by analyzing multiple biopsies of primary tumors, would be of great value to study genetic relationships between index tumors and LRs. Both chromosomal changes and point mutations in bona fide cancer driver genes can be used as markers of genetic relationship. Chromosomal changes are easiest determined, but their discriminating power to assess genetic relationships is more limited. Algorithms that consider chromosome breakpoints and other genetic changes in HNSCC often relate to losses of complete chromosome arms, and are consequently comparable in many different tumors, losing discriminative power. Somatic mutations in cancer driver genes are therefore considered to be most accurate markers to assess genetic relationships as they are very unique (21, 22). However, low-quality DNA from formalin-fixed paraffin-embedded (FFPE) specimens and PCR errors may cause false positive findings, sometimes with a high variant allele frequency (VAF), which could hamper the interpretation of the data. Therefore, a combination of both approaches might be the preferred strategy to assess genetic relationships.

The aims of this study are to evaluate copy number alterations (CNAs) as biomarker of genetic relationship by studying intratumor genetic heterogeneity using multiple spatially distinct biopsies obtained from individual tumors. Second, genetic relationships between LRs and index tumors in patients treated with chemoradiotherapy were investigated with this approach, and findings extended with mutation analysis of head and neck cancer driver genes.

To study the suitability of copy number alterations to assess genetic relationships we had to take into account intratumor genetic heterogeneity within distinct samples from one tumor (intratumor) and comparable alterations (e.g., frequent 3p losses) between different tumors (intertumor). We therefore collected and analyzed 44 biopsies from 13 resected primary oral tumors (3–5 biopsies/tumor). The samples were obtained from the resection specimen available at the department of Pathology at Amsterdam UMC (Amsterdam, the Netherlands), location VUmc, or at the University Hospital Parma (Parma, Italy), and were snap-frozen. All primary tumors were large enough to allow for widely spaced biopsies; care was taken to avoid areas of necrosis and ulceration; all sampled areas were documented histologically. None of the patients had received treatment prior to definitive surgery.

We included all patients who developed an LR after cisplatin-based chemoradiotherapy with curative intent for an advanced stage HPV-negative oropharyngeal, hypopharyngeal, or laryngeal squamous cell carcinoma at Amsterdam UMC, location VUmc, between 2009 and 2014. Criteria for LR were residual or recurrent tumor within 2 cm from and within 3 years after the index tumor. All patients were evaluated by medical history, physical examination, examination under general anesthesia/panendoscopy, and imaging (e.g., CT and/or MRI), staging was according to the 7th edition of the American Joint Committee on Cancer staging manual. Only HPV-negative OPSCC tumors were included, determined by p16 immunohistochemical (IHC) staining as described by Smeets and colleagues (23). Patients received intensity modulated radiotherapy (IMRT) concurrent with cisplatin (100 mg/m²every 3 weeks), with a curative intent.

In total, 113 patients with an HPV-negative oropharyngeal, hypopharyngeal, or laryngeal squamous cell carcinoma were treated within the given timespan with curative-intent chemoradiation at our clinic. In total, 20 (17.7%) of these 113 patients developed a local relapse (12) according to the clinical criteria described above or a regional relapse in the neck (8). From the 10 of 12 patients with an LR, we could retrieve the paired FFPE biopsies from the pathology archive. Clinical characteristics and follow-up period of the patients of whom the tumor and LRs were analyzed did not significantly differ from patients who developed a regional relapse or of whom tumor–LR pairs could not be analyzed (Supplementary Table S1). FFPE material with sufficient tumor purity, assessed by inspection of the hematoxylin/eosin–stained slide, was used for DNA isolation. When required, tumor percentage was enriched by macrodissection and was between 45% and 90%. Definitive tumor purity was estimated by ACE (see below).

Studies were carried out in accordance to the Declaration of Helsinki. According to the decisions of the Institutional Review Board and when patients consented, these studies were performed following the guidelines of the Code of Conduct for Human Tissue and Medical Research (https://www.federa.org/codes-conduct) and the EU General Data Protection Regulation.

Genomic DNA was obtained from the snap-frozen samples using the Life PureLink Kit (Thermo Fisher Scientific) and from the FFPE samples with either the Life PureLink Kit or the QIAamp DNA Micro Kit (Qiagen). DNA yield was analyzed on a Qubit 2.0 (Thermo Fisher Scientific). In total, a minimum of 200 or 250 ng DNA from snap-frozen or FFPE material, respectively, was used as input for library preparation. First, DNA was sheared on a Covaris S2 (Covaris). Library preparation for DNA isolated from FFPE was conducted using the TruSeq Nano Kit (Illumina), following the manufacturer's instructions. In short, sheared DNA was end-repaired, the 3′ ends were adenylated, indexed adapters were ligated, and the library was amplified with 10 PCR cycles. Libraries of DNA from frozen specimen were prepared with SOLiD reagents (Applied Biosystems) and amplification with 12 PCR cycles. The quality of the libraries was verified on a Bioanalyzer DNA 7500 chip (Agilent Technologies). Whole-genome sequencing libraries of 20 (FFPE) or 24 (fresh frozen) samples were pooled equimolarly and sequenced on a single lane of a HiSeq 2500 (Illumina) in a single-read 50-cycle run mode (SR50). Raw sequencing reads were mapped to the human reference genome (build GRCh37/hg19) with BWA (24). Data analysis was performed using the “QDNAseq” R package (version 1.16.0; ref. 25). Data were filtered against a blacklist of regions known to be germline copy number variants (26), and log₂ ratios were median-normalized and segmented with the Circular Binary Segmentation algorithm (27). Estimates of tumor purity and absolute copy numbers were obtained through the “ACE” R package (version 0.99.6). Calling of the segments was done using CGHcall (version 2.42.2; ref. 28).

The segmented data of CGHcall were used to compare copy number profiles. To compare copy number profiles and take the specific breakpoints into consideration, we calculated the log-likelihood ratio (LLR) using the “Clonality” R package (version 1.28.0), which quantifies the likeliness of two tumors being clonally related (29). In addition, we did a correlation analysis of segments as described by Sie and colleagues (30). Significance of both the LLR and correlation score for discriminating pairs and nonpairs was assessed by a permutation test. Under the null-hypothesis, the scores (LLR or correlation) of pairs and nonpairs are exchangeable. Therefore, we created a permutation null-distribution for the difference in mean score between pairs and nonpairs, where the pairs were created by permutation of the sample profiles. The P value for the observed difference in mean scores was then computed using this null-distribution.

Calculation of both assays were plotted in 2D and the optimal cutoff determined to define genetic relationships. We used the multi-biopsied OSCC samples to develop a classification model for genetic relatedness. For evaluating the performance of the predictive model, leave-2-out cross-validation (CV) was used as other analyses such as bootstrapping demand larger sample sizes. For evaluation by CV, it is essential that the left-out samples are completely independent of the training samples. Because biopsies from the same patient are related, we did not perform CV on the pairs, but in fact on the patients. In short, all samples of 2 patients were left out for testing. The profiles of all remaining patients were used for training the logistic regression model with the pairing indicator as response variable and LLR and correlation score as covariates. The model was then applied to all paired and nonpaired samples of the two left-out patients, and this was repeated for all possible sample couples. The predicted probability scores and true pairing indicators were then used to produce the ROC –curve (31), and the corresponding AUC. The optimal cutoff was determined by the highest sensitivity and specificity combination using the Youden index (32). This threshold was applied to the testing cohort to test the performance of the model.

In addition, mutation analysis of the index tumors and relapses of chemoradiotherapy-treated patients was performed. We selected a panel of 12 cancer driver genes which are frequently altered in HNSCC according to TCGA data (Supplementary Table S2; ref. 33). The sequencing libraries as described above were used for target-enrichment using the SeqCap EZ (Roche) protocol. We aimed for equimolar pooling of 50 ng of all twenty samples to a combined DNA amount of 1 μg. Because of low yields in some samples, input of individual samples ranged from 6.6 to 83.1 ng, which obviously impacted coverage. The library was hybridized to the SeqCap Oligo Pool. Streptavidin beads were used to capture the complex of oligos and DNA fragments, and unbound fragments were washed away. The captured fragment pool was amplified by 14 PCR cycles, and sequenced on a single lane of a HiSeq 2500 (Illumina) in a paired-end 150-cycle run mode (PE150). Raw sequencing data were mapped to the human genome (build GRCh37/hg19) with BWA and duplicate reads were marked for removal using Picard (http://broadinstitute.github.io/picard). Variants were called using multiple tools, including GATK MuTect2 (34) and Samtools mpileup (35) with VarScan 2 (36). To assess mutations, we required a minimum coverage of 70 unique reads. We filtered all variants with a frequency higher than 0.01% in the European population of the 1000 Genome Project cohort (37). In addition, we manually filtered variants that were apparent in both the primary and the relapse with a VAF of around 50%, independent of the tumor purity (germline SNPs). Estimated tumor purity, absolute copy numbers (computed with ACE), and VAF obtained with MuTect2 was used for calculation of the 95% confidence interval (CI) of the VAF and the absolute number of variant (mutant) alleles (aVA) in the tumor genome. To determine the genetic relation of primary and relapse, we looked at the cooccurrence of variants found in the primary, which were present in at least 33.3% of the tumor-derived alleles as determined by ACE, implying that the variant occurs in the larger part of the tumor cell population. Furthermore, we reviewed each variant manually using the IGV viewer to ensure that no false positives were called by the algorithms. Finally, variants that appeared to be specific for either the primary tumor or the relapse sample were cross-checked for low coverage reads to ensure that minor subclones in either the index tumor or relapse were not missed.

As a first step, we developed and tested an algorithm for assessing genetic relationships on basis of copy number alterations, taking intratumoral genetic heterogeneity into account. We collected multiple biopsies of surgical specimen and isolated DNA. We performed whole-genome low-coverage sequencing to generate copy number profiles for 44 samples of 13 patients. Similarity of profiles was assessed by (i) the calculation of a LLR, which is based on the concordance of gains and losses in two samples and (ii) a correlation score of segments of all intra- and interpatient combinations of tumor pairs. Calculations of the two methods were plotted in a 2D graph. We made a total number of 946 comparisons, of which 56 were intratumor comparisons and 890 intertumor comparisons. The patient-matched intratumor comparisons had a higher mean correlation coefficient (0.83 ± 0.18 vs. 0.25 ± 0.18) and a higher mean log-likelihood ratio (17.3 ± 21.7 vs. −5.5 ± 3.5) than the comparisons between different tumor pairs. Both means differed significantly (P < 0.001) analyzed by permutation test (10,000 iterations). To classify sample pairs as genetically related or genetically unrelated, we used a logistic regression model analysis on the basis of the LLR and the correlation score together. The probability of genetic relation by CNA profiles could be predicted by the equation: −13.11 + 19.87 (correlation score) + 0.50 (LLR) in our training cohort. Across patients, the AUCs ranged from 93% to 100%, as assessed by comparing biopsy couples within a given patient with couples when biopsies were from different patients. An optimal cutoff at P(Y) > 0.094 was determined by the highest combination of sensitivity and specificity (Youden index) of the classification. Both the correlation score and LLR were significant predictors of genetic relationship (P < 0.001). Our model could identify sample pairs in our testing cohort as being biopsies from the same tumor with a sensitivity of 95%, a specificity of 96%, and an overall accuracy of 96%, (Fig. 1A and B)

Because the developed algorithm using CNAs appeared to be sufficiently accurate to investigate genetic relationships within tumors and between independent tumors, we studied paired primary tumors and local relapses of chemoradiation-treated HNSCC patients likewise. From a consecutive cohort of 113 patients with HPV-negative tumors treated with chemoradiotherapy between 2009 and 2014, 20 cases developed a local or regional relapse within a median follow-up time of 35.7 months (1.05–85.3; Supplementary Table S3). In 10 of 12 cases with a local relapse, the paired biopsies could be retrieved from the pathology archive. The studied cohort of these 10 cases with local relapse consisted of predominantly male patients (80%) all with locally advanced disease (T3–T4), 60% of whom had a primary tumor in the oropharynx (all HPV-negative). On the basis of imaging during follow-up, two patients were diagnosed with residual disease after treatment, which was confirmed by histopathologic assessment of a biopsy or the specimen after salvage surgery (Table 1). In the remaining 8 patients, locoregional control seemed achieved but relapses occurred nonetheless. The median time to diagnosis of residual or recurrent disease, both considered as local relapses, was 8.2 months (5.8–35.1). As expected for patients with relapsed HNSCC, they had a poor overall survival.

From both primary tumors and local relapses, copy number profiles were established by low-coverage whole-genome sequencing. The copy number profiles showed the characteristic alterations for HNSCC with losses and gains at 3p, 8p, 9p, and 11q. Furthermore, a high frequency of amplifications and gains was observed in regions containing PIK3CA (85%), EGFR (65%), and CCND1 (60%; Fig. 2)

We applied the genetic relationship algorithm described above to the corresponding pairs of primary tumors and relapses. Primary tumor and relapse of two patients were designated by the algorithm as genetically related, VUMC0905 and VUMC1406, while the eight other pairs were identified as genetically unrelated, although two were borderline and within the 0.95 CI of the cut-off value (Fig. 1C).

The result that six or maybe even eight of the local relapses in chemoradiotherapy-treated patients seemed not genetically related to the index tumors by copy number profile analysis was unexpected. This could imply that despite the results reported above, CNA profiling is less suited for assessment of genetic relationships of index tumors and LRs after chemoradiotherapy, also because the biological context is very different. Therefore, we decided to complement the data by mutation analysis of genes commonly mutated in HNSCC and that are assumed cancer driver genes (Supplementary Table S2). All 10 paired tumor and relapse samples were sequenced at a median sequencing coverage of 527×, but with a large range (19–3,790) depending on the sample quality and nucleotide position. In one sample, DNA quality was too poor and coverage too low to pass our quality checks, and we lost that pair for mutation analysis, but it was convincingly genetically related by copy number analysis. In the nine remaining tumor pairs, sequence variants were found in 10 of 12 tested genes, and the median number of variants were 4 (0–8) and 3 (2–8), in respective primary tumors and relapses.

As examples, we displayed two cases in Table 2, one that we consider as genetically related and one apparently not. We applied the ACE algorithm on the data that corrects for stromal contamination and ploidy, and provides an estimate of the absolute number of alleles in the tumor cells, and the number of variant alleles. When we analyzed the variants with high number of absolute variant alleles (aVA) and that thus were present in the very large majority of the tumor/relapse cell population, we found that these were shared in case VUMC1406 between primary and relapse, in accordance with the result of the CNA profiling. Case VUMC1406 had a TP53 (c.716A>G) variant that occurred in all of its tumor alleles of the primary tumor (absolute tumor alleles/variant alleles; 0.98/1.02) and in an equal genetic composition in the relapse. This suggests that in this particular case, the original population of cells in the index tumor had recurred. However, we also identified somatic variants in apparently a smaller fraction of tumor cells (range: 0.05%–27.9%) in either the primary tumor or the relapse, which likely reflects slight variations in tumor cell populations. The same was observed in all other cases (VUMC1122, −1419, −1102, and −1322) that showed a genetic relationship within tumor and relapse pairs based on one or more cancer driver mutations with a high aVA in the index tumor (Fig. 2; Supplementary Table S4).

In five of the nine cases, there seemed to be no relation at all between the tumor and relapse. As an example, in case VUMC1417 the primary tumor had two mutations in TP53, one with an aVA of 83.4% (c.1022_1023delTC) and one with 7.0% (c.1013A>G). The mutation with a frequency of 83.4% should certainly be considered as a cancer driver mutation but it lacks completely in the relapse with 2,779 reads on that position. Likewise, the relapse of VUMC1417 had double mutations in FAT1 and TP53, but these variants all lacked in the primary tumor, although it should be mentioned that the coverage in the tumor was somewhat lower on the respective positions but still 96, 68, 181, and 151 reads. All four remaining cases followed this pattern of mutually exclusive mutations: VUMC1004, −0923, −1417, and −1123. Hence, these relapses originated from very minor subclones in the tumor, or they had an independent genetic origin. A striking example is case VUMC1123 with a high VAF mutation in TP53 (c.833C>T) in the tumor, and no mutant reads in 1,495 reads on that position in the relapse, while the relapse had another high VAF TP53 (c.574C>T) mutation with no mutant reads in 984 reads on that position in the index tumor.

Assessment of the genetic landscape of tumors, relapses, and metastases, is technically challenging when frozen specimens are not available. Often only archived FFPE tissue material is available, which gives poor DNA yields and quality. This negatively impacts genetic analyses as the prepared sequencing libraries are generated on only few DNA strands, hampering coverage and introducing more PCR errors. For CNA analysis this is not too problematic; only the noise will increase, which can be reduced by increasing the bin size. For mutations, it is more problematic as a basic coverage of at least 30 deep is required to identify a mutation in 10% of the tumor cells in 3 reads. Moreover, formalin fixation may cause nucleotide deamination, causing high sequencing error rates. These limitations can in part be overcome if mutation analysis is restricted to established cancer driver genes that usually show high VAFs and are biologically relevant. However, clonal evolution during tumor progression may impact also bona fide cancer drivers.

Although many HPV-negative HNSCC tumors share specific genomic alterations with frequent patterns of losses in chromosomal regions 3p and 8p, and gains in 3q, 5p, and 8q (33), most tumors exhibit a unique combination of CNAs, particularly when breakpoints are considered. The combination of the log-likelihood ratio (29) and the correlation of segments (30) on the 44 biopsies of 13 oral cavity samples could make an accurate classification between different biopsies from one tumor and independent tumors using copy number profiles. The high accuracy (96%) of our classification model reflects the much higher degree of inter- than intratumoral heterogeneity, when considering copy number alterations. This result is in line with a recent study by Tabatabaeifar and colleagues (38). Unfortunately, analysis of genetic relationships by CNA profiling seemed less conclusive for tumors and local relapses of chemoradiotherapy-treated patients. This might relate to the assay itself, but more likely it reflects the complex pathobiological origin of the relapses. Therefore, we extended the biomarker panel with somatic mutation analysis of the known and established HNSCC driver genes. When we consider mutations as gold standard, accuracy of CNA profiling is indeed less accurate to assess genetic relationships as it was in 6 of 9 cases correct (67%). As indicated above, this likely relates to the nature of the relapse after chemoradiotherapy treatment, and comparable studies to lymph node metastases may shed more light on this issue. Moreover, the analysis of 3-5 biopsies might not capture the genetic heterogeneity of a tumor completely, and we may have underestimated intratumor heterogeneity. Nonetheless, CNA profiles are extremely helpful to correct tumor percentages and ploidy using algorithms such as ACE, and provides quantitative estimates of the mutant allele numbers (39).

We showed by mutation analysis that 5 of 10 relapses of chemoradiotherapy-treated patients analyzed seemed genetically related to the tumor based on one or more cancer driver mutations with a high frequency of mutant alleles. This is not unexpected as tumor cell populations may vary between tumor and relapse, which may relate to the extra cell divisions in relapses but which is also in line with the view that some populations are more resistant to particular treatments than others (40).

Most remarkable are the five cases of which the relapses seemed not to be genetically related to the index tumor. Whether this truly refers an independent origin remains a challenging question. The availability of only FFPE specimen hampered DNA sequencing, and in some pairs the coverage was somewhat low. Moreover, chemoradiotherapy treatment might have selected resistant clones that were present in the index tumor but with frequencies of less than 1% or even 0.1%, and that are not picked up with the current sequencing coverage of a single randomly taken biopsy. Moreover, treatment may induce genetic changes, although we would expect additional mutations and copy number alterations, and not a disappearance. Hence, this cannot explain why cancer driver mutations in the primary tumor are not present in the relapse, and it does not really support that treatment induced changes underlie the observed phenomenon.

LRs remain a major obstacle in the success rate of HNSCC treatment. Most of our knowledge on relapses is based on studies performed on surgically treated OSCC (13, 14, 16, 17). Despite the high frequency of local failure, studies to investigate the local relapses after treatment with cisplatin-based chemoradiotherapy are scarce, due to the lack of available biological material in the absence of surgery. Incidence of locoregional failure after concurrent chemotherapy varies between 15% and 50% (7, 11, 12). In our cohort of 113 patients, only 20 (17.7%) patients developed a local and/or regional relapse, of whom 12 (10.6%) had a local relapse. These numbers are low but still comparable to other studies on cisplatin-based concomitant chemoradiotherapy in the IMRT era (12). Nonetheless, the relatively effective treatment in this cohort might relate to the remarkably high number of genetically unrelated relapses.

We evaluated whether the applied treatment protocols played a role. Development of local recurrences from residual cells can be explained by insufficient treatment or an inadequate response to treatment. Patients receiving less than a cumulative dose of 200 mg/m²cisplatin have a significant higher risk on developing a locoregional relapse (41). Also, a lower total irradiation dose or greater tumor volume can potentially lead to the outgrowth of residual cells (42). All our patients received the planned dose of 70 Gy on their primary tumor. However, two patients (VUMC1123 and VUMC1322) did not reach a cumulative dose of >200 mg/m²cisplatin due to inacceptable nephrotoxicity, although patient VUMC1322 switched to a weekly carboplatin regiment. Interestingly, this somewhat less optimal primary tumor management did apparently not lead to outgrowth of the index tumor as the origin of the relapses of both VUMC1123 and VUMC1322 seemed genetically unrelated.

Several molecular mechanisms of action and resistance to chemoradiotherapy are proposed to explain inadequate response to treatment, and preclinical studies show many different mechanisms for resistance (43). In HNSCC the fraction of DNA-bound platinum and genes in the DNA repair pathway (especially Fanconi anemia/BRCA pathway) are associated with the response to treatment (44, 45). In clinical studies specific mutations (e.g. TP53, NOTCH1; refs. 42, 46) or expression of specific genes (e.g., hypoxia-related genes or SDF-1 and CXCR4; refs. 42, 47) have been associated with local treatment failure. These studies generally assume that the development of local relapses relates primarily to the outgrowth of residual tumor cells. The idea that so many relapses seem to be genetically unrelated to the bulk of the index tumor may impact the interpretation of these studies and could explain why we failed to identify clinically meaningful predictive biomarkers. A previous study by Hedberg and colleagues also described a case where there was no genetic relationship between the index tumor and the recurrence, and the authors suggested that this recurrence might have originated from a second independent field or developed as a true second primary tumor (17). Irrespective of the precise origin of these genetically unrelated relapses, this observation warrants further studies. The impact is tremendous for future studies, as it will imply that gene profiling or radiomics studies on the bulk index tumor will not always be able to indicate the accurate risk for relapse. Associations may improve when the genetically unrelated relapses are filtered out. Also, the evaluation of novel experimental agents in clinical trials might benefit from the precise analysis of the relapse type that occurred. In addition, relapses may be treated as new tumors. Given these potential impacts, similar studies should be performed but in larger cohorts, and preferably making use of multiple frozen biopsies to allow ultradeep sequencing and better estimates of intratumoral heterogeneity.

Our study provides insight in the complex biology of relapsed HNSCC after chemoradiotherapy, and might have large consequences for prognostic modeling, the use of predictive and prognostic biomarkers and therapeutic innovations.

C.R. Leemans is a consultant/advisory board member for MSD Global Advisory Board HNSCC. No potential conflicts of interest were disclosed by the other authors.

Conception and design: R.H. de Roest, S.W. Mes, R.H. Brakenhoff

Development of methodology: R.H. de Roest, J.B. Poell, R.H. Brakenhoff

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): R.H. de Roest, S.W. Mes, A. Brink, E. Bloemena, E. Thai, T. Poli, R.H. Brakenhoff

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): R.H. de Roest, S.W. Mes, J.B. Poell, A. Brink, M.A. van de Wiel, C.R. Leemans, R.H. Brakenhoff

Writing, review, and/or revision of the manuscript: R.H. de Roest, S.W. Mes, J.B. Poell, A. Brink, M.A. van de Wiel, E. Bloemena, T. Poli, C.R. Leemans, R.H. Brakenhoff

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): A. Brink, R.H. Brakenhoff

Study supervision: R.H. Brakenhoff

Other (sampling of surgical piece from the Surgery Unit of Parma; diagnosis and staging of the tumor; collection and dispatch of biopsy samples): E. Thai

We would like to thank D. Sie for providing his “clonality” tool to assess the correlation and log-likelihood analysis of CNAs and Prof Enrico Silini for providing patient material. This work was supported by the Design Project (KWF-A6C7072).

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.