Integrative and Comparative Genomic Analysis of HPV-Positive and HPV-Negative Head and Neck Squamous Cell Carcinomas

Head and neck squamous cell carcinoma (HNSCC) is the fifth most common non-skin cancer worldwide with an annual incidence of 600,000 cases and a mortality rate of 40% to 50% despite aggressive treatment (1, 2). The major known risk factors are environmental exposure to tobacco products, alcohol, and infection with high-risk human papilloma viruses (HPV). The incidence of HPV-positive tumors is rising rapidly in Western countries and HPV status is the strongest clinically applicable prognostic marker, portending a favorable prognosis (3, 4).

Although HNSCC is widely viewed as comprised of two distinct clinical entities, HPV-positive and HPV-negative tumors, a comprehensive list of differential molecular abnormalities, in particular therapeutically relevant genetic aberrations, has not been reported. A particular problem is the lack of study of HPV-positive HNSCC: Currently no large series of HPV-positive tumors exist and the upcoming TCGA (The Cancer Genome Atlas) cohort is comprised of 85% HPV-negative tumors (5). This bias is likely related to selection of surgically resected, earlier stage oral cavity and laryngeal tumors. This may not be representative for clinically more complex stage IV tumors requiring multimodality or palliative treatments (6–8).

Unlike lung or breast adenocarcinomas, there are currently no defined targetable genetic aberrations for HNSCC, and no approved therapies are tied to genetic alterations as predictive biomarkers. All patients with HNSCC are treated with a largely uniform approach based on stage and anatomic location, typically using surgery, radio-, and chemotherapy alone or in combination (9). Cetuximab, an anti-EGFR antibody, is the only approved targeted therapy for HNSCC with a single-agent response rate of 10% to 13%. Despite the modest response rate. there are no validated predictive biomarkers for benefit from cetuximab (10, 11).

Previous studies have demonstrated frequent mutations of several genes in cohorts of largely HPV-negative HNSCC, most notably TP53, PIK3CA, CDKN2A, the TERT promoter, and NOTCH pathway gene alterations (12–16). However, the genetic makeup of HPV-positive HNSCC remains unclear (15).

In this study, we investigated a fully annotated patient cohort of 120 locoregionally advanced HNSCC (including 42.5% HPV-positive tumors), treated uniformly with organ-preserving chemoradiotherapy using massively parallel sequencing, copy-number profiling, and validation.

We discover distinct mutational and copy-number profiles in HPV-positive and HPV-negative tumors and identify for the first time potentially targetable mutations and copy-number aberrations that are of high translational relevance.

Pretreatment tumor tissues (n = 120) and matched normal DNA for patients with locoregionally advanced HNSCC treated at the University of Chicago (Chicago, IL) were obtained from the HNSCC tissue bank (UCCCC#8980).

An overview of the tissue processing is provided in Supplementary Fig. S1 and described in detail in the Supplementary Methods.

HPV16/18 status was determined by E6/E7–specific qRT-PCR. Results were corroborated by additional tests to increase accuracy, including an E6/E7 DNA–based multiplex PCR for five high-risk HPV types (17) as well as p16/CDKN2A expression, and TP53 mutations (18).

DNA sequencing libraries were prepared following published protocols (19) and enriched using custom capture reagents [Agilent; Nimblegen (validation)]. Of note, 2 × 100 bp paired-end sequencing occurred using Illumina HiSeq 2000/2500 sequencers. Of note, 617 cancer-associated genes (Supplementary Table S1A) were targeted and sequenced at high-depth in tumor/normal pairs. Median coverage at targeted bases was 231× in the tumors and 254× in the matched normal. Mutation and indel calling was done using bioinformatic pipelines as described previously using MuTect (20–22). Furthermore, we used VarWalker to determine potentially relevant genetic aberrations (mutations and copy-number aberrations) using prioritization via a protein–protein interaction network focusing on frequently mutated and cancer gene census genes (23).

We used an established machine learning–based approach—Cancer-Specific High-Throughput Annotation of Somatic Mutations (CHASM)—to predict and prioritize missense mutations leading to functional changes and, thus, likely driving tumorigenesis (24). This approach was previously validated to show high specificity, and to a lesser degree sensitivity to identify “drivers” of oncogenesis (24).

Copy-number analysis was performed using sequencing data and the CONTRA algorithm (25). Results were validated on the Nanostring nCounter using predefined/custom cancer gene panels (Nanostring).

Detailed methods are available in the Supplementary Methods.

The median age of patients was 56 years with a median follow-up of 48 months (Table 1). Fifty-one patients (42.5%) were HPV-positive. Of note, 45% of patients had no or light smoking histories whereas 55% had significant smoking histories, including 35% of HPV-positive tumors (N = 18). The cohort consisted of locoregionally advanced tumors with 67 (55.8%) oropharynx and 115 (95.8%) stage IV tumors (Table 1).

Analysis of targeted hybrid capture sequencing data identified a total of 5,476 point mutations and 4,562 insertion/deletion (indel) events across 120 tumor/normal pairs in targeted 617 genes (Supplementary Table S1A). Of note, 1,536 of the point mutations and 200 of the indel events were in coding regions. Of the 1,536 point mutations in coding regions, 1,209 resulted in an amino acid changes in the corresponding proteins. Individuals displayed a mean point mutation rate of 2.33 mutations per Megabase (Mb) in the targeted regions and 0.22/Mb in the coding territory alone with a range of 0 to 96 mutations in the coding regions of the targeted genes. Among nonsynonymous substitutions, transitions and transversions at CpG sites, which are typically associated with tobacco exposure (26), were the most commonly observed mutational context with a rate of 21/Mb. Mutations at Tp*Cp(A/C/T) sites, a substitution type commonly associated with virally induced cancers (26), were the second most commonly observed context at 13/Mb and enriched in the HPV-positive patients. We did not detect a significant difference in overall mutation rates by HPV status at all targeted bases (HPV-negative 2.22/Mb, HPV-positive 2.46/Mb, P = 0.64) or in the coding territory (HPV-negative 0.23/Mb, HPV-positive 0.21/Mb, P = 0.79). Coding mutations were validated at a rate of 87% (338/390) in a second independent hybrid capture experiment when evaluating statistically enriched mutations (Supplementary Table S1B and S1C).

We applied the MutSig method to identify cancer-relevant genes demonstrating evidence of statistical selection for mutations in our cohort (20–22), (27). MutSig identified 11 genes displaying significant enrichment for mutations as defined by a q value of <0.1 in the overall cohort (Supplementary Fig. S2): TP53, CDKN2A, PIK3CA, MLL2, TPRX1, CUL3, FLG, NSD1, DDX3X, RPIK4, and HRAS. FGFR3 and FBXW7 demonstrated q values between 0.1 and 1. A second MutSig analysis, in which only genes annotated in the COSMIC database were considered, rediscovered TP53, CDKN2A, PIK3CA, HRAS, FGFR3, and FBXW7 and additionally nominated KRAS, MLL3, FGFR2, ZNF217, and RIMS2 as statistically significant with q < 0.1. We performed additional significance analyses on the separate HPV-negative and HPV-positive cohorts using the same MutSig algorithm. Analysis of mutated genes displaying enrichment in the HPV-negative cohorts demonstrated statistical enrichment for mutations of TP53, CDKN2A, MLL2/3, NOTCH1, PIK3CA, NSD1, FBXW7, DDR2, and CUL3, in the HPV-negative samples (Fig. 1A). Genes displaying statistical enrichment in HPV-positive tumors were PIK3CA, MLL3, DDX3X, FGFR2/3, NOTCH1, NF1, KRAS, and FBXW7 (Fig. 1B).

Somatic copy-number alterations were inferred from sequencing data by read-depth analysis using the CONTRA methodology (Supplementary Methods). We applied a modification of GISTIC 2.0 analysis (28) compatible with the focal sequencing data to identify recurrent peaks of amplification and deletion and identified 192 genes with false discovery rate of less than 0.25 (Supplementary Table S1D–S1G). Copy-number alterations were validated at a rate of 83% (45/54) using the Nanostring nCounter, focusing on the genes found to be significantly altered in this dataset (Supplementary Table S1H). Copy-number analysis demonstrated many previously demonstrated regions of amplification and deletion, including focal gains of EGFR, REL, BCL6, PIK3CA, TP63, CCDN1, and MDM2 and losses of ATM, CDKN2A, RB1, NOTCH1, and NF1 (Fig. 1C and D; refs. 12–15, 29). Significantly amplified regions that occurred primarily in HPV-negative tumors were amplifications of 11q13 likely targeting CCND1, 7p11 (EGFR; Fig. 2A). Amplification of 3q26–28, a region containing SOX2/TP63/PIK3CA (Fig. 2A and B) occurred in both HPV-positive and HPV-negative tumors. 3p loss and CDKN2A deletions (Fig. 2C) occurred primarily in HPV-negative tumors whereas ATM deletions occurred primarily in HPV-positive tumors (Fig. 2D).

To identify additional biologic differences between HPV-positive and HPV-negative tumors, we used an unbiased objective protein–protein interaction network-based analysis of genetic aberrations (VarWalker; Fig. 3). The altered networks for HPV-negative tumors contained 84 proteins connected by 76 interactions (Fig. 3A). For HPV-positive tumors, the altered networks contained 88 proteins connected by 83 interactions (Fig. 3B). The two altered networks were substantially different: Alteration of p53 signaling and cell-cycle pathway genes occurred almost exclusively in HPV-negative tumors. The alteration of oxidative stress pathway genes (CUL3, NFE2L2, and KEAP1) occurred more often in HPV-negative tumors (Fig. 3A). The alteration of the DNA damage pathway (BRCA1, BRCA2, FANCG, FANCA, FANCD2, and ATM), FGF signaling (FGFR2, FGFR3, and FGFR4), JAK/STAT signaling (STAT1, JAK1, and JAK2), and immunology related genes (HLA-A and HLA-B) favored HPV-positive tumors (Fig. 3B). HPV-positive and HPV-negative tumors shared network alterations such as PI3K signaling, Notch aberrations, and SMAD signaling (Fig. 3A and B).

Somatic mutations in potentially targetable kinase genes occurred in FGFR2 and FGFR3 uniquely in HPV-positive tumors with an incidence of 17.6% among HPV-positive tumors. FGFR2 mutations included N569D and N569K (both in HPV-positive), and FGFR3 mutations included six S249C mutations, which were validated by Sanger Sequencing (all in HPV-positive), as well as a K413N somatic mutation (Fig. 4A and B). We used a previously validated machine learning–based approach named CHASM to predict and prioritize missense mutations with a high likelihood of being drivers of oncogenesis (24). Both FGFR2 N569K and FGFR3 S249C have been described in several cancer types (30, 31) and showed low CHASM scores suggestive of oncogenic driver character (Fig. 4A and B). Importantly, the FGFR3 S249C mutation was identified recurrently in six HPV-positive tumors. S249C has been reported in one HPV-positive tumor in the TCGA HNSCC cohort and in the TCGA lung squamous cell and bladder carcinoma cohorts. A fraction of DDR2 and EPHA2 mutations also showed low CHASM scores, but were not recurrent, and did not cluster in any domain (Fig. 4A and B).

PIK3CA was the most commonly altered oncogene, including a number of established canonical mutations (E542K, E545K, and H1047R; Fig. 4A and B). Other PI3K pathway mutations included somatic events in PIK3R1, PTEN, TSC1, and TSC2, all of which have been reported to harbor mutations in other tumor types (Fig. 4A). PIK3CA mutations were more commonly observed in HPV-positive individuals, though the difference was not statistically significant in our cohort (P = 0.20).

Mutations in the MAPK pathway genes included seven canonical HRAS (G13V, Q61L, and K117N), KRAS (G12D, G13D, and L19F), and NRAS (Q61R) events (incidence 5.8%), as well as multiple mutations of unclear significance in the tumor suppressor NF1 (Fig. 4A and 4B).

We first examined the mutation difference in smokers versus nonsmokers (Supplementary Table S2). In the overall cohort tumors from smokers displayed a higher overall mutation burden as compared with nonsmokers 11.5 mutations per patient (≥10 pack years) versus 5.6 mutations per patient (nonsmokers; P = 0.0007, t test). We noted a significant association among smoking and mutations in TP53 (P = 0.0009), CSMD3 (P = 0.0062), RB1CC1 (P = 0.0161), THSD7A (P = 0.0319). Mutations in ZFHX4 (P = 0.0167, Fisher exact test) and TRRAP (P = 0.0392, Fisher exact test) were significantly enriched in smokers compared with nonsmokers. Additional subgroup analyses, including drinkers versus nondrinkers, large primary versus small primary, and node negative versus node positive are provided in Supplementary Table S2, including subgroup analyses for oropharyngeal primary tumors.

The association among HPV status and prognosis was evaluated. As expected HPV-positive patients had significantly better prognosis than HPV-negative patients for both overall survival (OS; P = 0.015) and progression-free survival (PFS; P = 0.013), and this holds true for the overall cohort as well as oropharyngeal cancers (OS, P = 4e−05; PFS, P = 0.0001; Supplementary Fig. S3A).

In an exploratory, hypothesis-forming analysis, we correlated gene aberrations, including mutations and copy-number alterations with OS and PFS. In the smaller sub-cohorts of HPV-positive and HPV-negative tumors most associations were not significant, though we did identify in HPV-negative tumors a possible correlation of PIK3CA mutations (OS, P = 0.017; PFS, P = 0.004; Supplementary Fig. S3B), as well as TP53 wild-type status with poor prognosis (OS, P = 0.262; PFS, P = 0.017; Supplementary Fig. S3C). These findings require validation in future studies.

HPV-positive and HPV-negative head and neck cancers are distinct clinical entities. We report for the first time on the underlying differential mutational profiles, generated from a unique patient cohort of poor-prognosis, locoregionally advanced tumors. Unlike prior exome-focused studies, we used a focused 617 cancer gene-candidate approach, allowing deep coverage of cancer-relevant genes. This approach is similar, albeit significantly broader than current clinical grade next-generation sequencing assays (32).

This is the first cohort to include a large number of HPV-positive tumors, including patients with a significant tobacco history (35% of HPV-positive tumors). Although HPV-positive status confers a favorable prognosis, patients with >10 pack-year smoking history have a poorer prognosis (3), and we report a higher mutational burden, including presence of KRAS mutations in such HPV-positive tumors. We identify targetable genetic aberrations, in both HPV-positive and HPV-negative tumors, including frequent FGFR2/3 aberrations in the latter group.

Genetic aberrations correlate closely with HPV status, and less with anatomic site (Supplementary Table S2), suggesting that HPV status (not anatomic site) is the most important factor determining tumor biology. Furthermore, it was recently reported that anatomic site cannot always be accurately determined and misclassification occurs (18).

The mutational spectrum in HPV-negative HNSCC is very similar to lung and esophageal squamous cell carcinomas (SCC) with enrichment for mutations in TP53, CDKN2A, MLL2, CUL3, NSD1, PIK3CA, and NOTCH genes, and copy-number increases in EGFR, CCND1, and FGFR1. In the future, these genetic similarities may potentially allow biologically driven cross-tumor drug development and biomarker discovery.

In contrast, HPV-positive tumors show a distinct genetic profile with unique mutations in DDX3X, CYLD, and FGFR and enrichment for PI3K pathway alterations and rarer KRAS mutations. Somatic aberrations in DNA-repair genes (BRCA1/2, Fanconi anemia genes, and ATM) may contribute to chemo- and/or radiosensitivity of HPV-positive tumors, and interestingly occurred in HPV-positive tumors in non-/light smokers. Both tobacco smoke and defects in DNA repair are known to induce a large number of genetic aberrations, and may be distinct ways to accumulate genetic aberrations required for the emergence of cancer. Of note, a recent RAD51B was reported as an integration target for HPV16 leading to loss of its DNA-repair function (5).

FGFR and PI3K pathway aberrations are potential therapeutic targets in this population. FGFR2/3 mutations are of particular interest as they occurred in 17.6% of HPV-positive tumors, most commonly the S249C mutation, which has been shown to be an oncogenic driven in bladder cancer.

One HPV-negative tumor harbored amplification of FGFR1, which is being explored clinically for lung SCCs (33, 34). Furthermore, oncogenic FGFR3–TACC3 fusions were recently reported in HPV-positive HNC, supporting a prominent role for oncogenic FGFR signaling for HPV-positive HNC whereas other kinases are genetically unaltered (33–37).

In contrast with prior smaller reports, the overall mutational burden in HPV-negative and HPV-positive HNSCC was comparable (15.2 vs. 14.4 somatic exonic mutations in the targeted cancer-associated genes), a finding corroborated when analyzing publically available data from the TCGA HNSCC cohort (5).

HPV-positive tumors in patients with significant tobacco use had a significantly higher mutational burden. Although the number of KRAS (N = 3) and TP53 (N = 2) mutations in HPV-positive oropharyngeal tumors was small, they occurred in patients with significant tobacco history, and showed 4- to 6-fold higher average mutational burden (not significant). Interestingly, the two TP53 mutations (H179R and G361fs) cooccurred with KRAS (G12D and L19F) and one PIK3CA mutation (F872L), suggesting multiple oncogenic drivers in highly mutational altered HPV-positive tumors. The KRAS L19F mutation is thought to be modestly transforming (38).

Low-frequency canonical HRAS, KRAS, and NRAS mutations were identified in both HPV-positive and HPV-negative tumors. RAS mutations have been associated with poor outcome in other cancer types (e.g., lung adenocarcinomas), as well as resistance to cetuximab in colon cancer (39). Further study in HNC will be important, especially about the implications of KRAS mutations in otherwise good-prognosis HPV-positive tumors, which were present in 5.8% of HPV-positive tumors. One additional case of a canonical KRAS mutation in an HPV-positive tumor was previously reported (29).

The specific genetic events in HPV-positive or HPV-negative HNSCC likely contribute to the distinct biologic behavior and may affect therapy in the future. The lack of targeted therapies for squamous cell tumors remains challenging: In addition to FGFRs, we identify multiple potentially targetable genomic alterations, including mutations in the DDR2, EPHA2 kinases, PIK3CA, and PI3K pathway genes (e.g., TSC1/2 and NF2) that are predicted to be oncogenic drivers. Therapeutic implications for these particular changes are being investigated clinically and preclinically (33, 40, 41).

The lack of EGFR aberrations in HPV-positive tumors is consistent with prior reports (5, 42), and clinical data (10, 11). Other potential therapeutic targets include alterations of genes, which regulate the cell cycle [e.g., CCND1 amplification as a target for cyclin-dependent kinase 4/6 (CDK4/6) inhibitors].

We identify a number of additional novel mutated pathways, providing further insight into the molecular pathogenesis of HNSCC: (i) Similar to lung SCC, MLL2 and MLL3 mutations were observed in both HPV-positive and HPV-negative tumors. MLL genes encode histone lysine methyltransferases that are involved in chromatin remodeling and recurrently mutated in several cancer types (43–46). Other studies have demonstrated an association with poor clinical outcomes (47); (ii) CUL3, NFE2L2, and KEAP1 mediate cellular responses to oxidative stress and are a frequent target for mutations in lung SCCs as well as the HPV-negative HNSCC cohort. KEAP1 binds to NFE2L2 and targets NFE2L2 for ubiquitination and subsequent degradation by the CUL3/RBX1 complex (48). Mutations are typically associated with heavy tobacco exposure, and distributed across the gene, which is suggestive of loss of function mutations; (iii) The DDX3X gene was exclusively mutated in HPV-positive tumors. DDX3X encodes an ATP-dependent RNA helicase and is involved in RNA processing (Supplementary Fig. S4). It has been described to be significantly mutated in medulloblastomas, and can potentiate β-catenin activity by promoter transactivation (43, 44); (iv) In addition, HPV-positive tumors harbor mutations in CYLD, a gene that has not been implicated as a hotspot for mutations, but is well known for its role in HPV-positive cervical cancer: CYLD is involved in differentiation and NFκB signaling and in response to hypoxia activates prooncogenic processes (Supplementary Fig. S4; ref. 49); (v) Finally, UBR5 mutations are involved in DNA damage and apoptosis signaling, and in addition to our discovery in HPV-positive and HPV-negative tumors, UBR5 is recurrently mutated in mantle cell lymphomas (Supplementary Fig. S4). Investigations into its role in radio- or chemosensitivity are promising (50).

In conclusion, the mutational landscape of HPV-positive HNSCC differs from HPV-negative HNSCC and may help explain the distinct clinical behavior and prognosis. Further study of clinical implications of the observed mutations will be vital with respect to targeted therapies, radio-, and chemotherapy. For the first time, we identify targetable genetic aberrations in both HPV-positive and HPV-negative HNSCC. In particular, frequent aberrations in the FGFR family are promising as well as PI3K pathway and cell-cycle abnormalities. The significant number of potential targets, including rare aberrations (DDR2, TSC1/2, and EPHA2), or candidate biomarkers such as HRAS/KRAS, indicates a need for tumor profiling of HNSCC in the future—including on clinical trials, as well as functional validation in preclinical studies.

P.S. Hammerman is a consultant/advisory board member for ARIAD, Janssen, and ImClone. No potential conflicts of interest were disclosed.

Conception and design: T.Y. Seiwert, M.K. Keck, E. Blair, E.E.W. Cohen, E.E. Vokes, K.P. White

Development of methodology: T.Y. Seiwert, Z. Zuo, M.K. Keck, M.S. Lawrence

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): T.Y. Seiwert, M.K. Keck, A. Khattri, T. Stricker, A. Langerman, E. Blair, K. Stenson, M.W. Lingen, E.E.W. Cohen, E.E. Vokes, K.P. White, P.S. Hammerman

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): T.Y. Seiwert, Z. Zuo, M.K. Keck, C.S. Pedamallu, T. Stricker, C. Brown, T.J. Pugh, P. Stojanov, J. Cho, M.S. Lawrence, G. Getz, R.R. Weichselbaum, E.E.W. Cohen, E.E. Vokes, K.P. White, P.S. Hammerman

Writing, review, and/or revision of the manuscript: T.Y. Seiwert, Z. Zuo, M.K. Keck, A. Khattri, C.S. Pedamallu, T.J. Pugh, J. Brägelmann, R. DeBoer, R.R. Weichselbaum, A. Langerman, L. Portugal, E. Blair, M.W. Lingen, E.E.W. Cohen, E.E. Vokes, K.P. White, P.S. Hammerman

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): T.Y. Seiwert, Z. Zuo, R. DeBoer, E.E. Vokes, P.S. Hammerman

Study supervision: T.Y. Seiwert

Other (funding): T.Y. Seiwert

This work was supported by an ASCO Translational Professor ship award (to E.E. Vokes), the Grant Achatz and Nick Kokonas/Alinea Head and Neck Cancer Research Fund (to E.E. Vokes and T.Y. Seiwert), a Flight Attending Medical Research Institute Young Clinical Scientist Award (FAMRI, YCSA; to T.Y. Seiwert), and a NCI K08 Award 1K08CA163677 (to P.S. Hammerman).

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.