Signatures of somatic mutations and gene expression from p16INK4A positive head and neck squamous cell carcinomas (HNSCC)

Abstract

Human papilloma virus (HPV) causes a subset of head and neck squamous cell carcinomas (HNSCC) of the oropharynx. We combined targeted DNA- and genome-wide RNA-sequencing to identify genetic variants and gene expression signatures respectively from patients with HNSCC including oropharyngeal squamous cell carcinomas (OPSCC). DNA and RNA were purified from 35- formalin fixed and paraffin embedded (FFPE) HNSCC tumor samples. Immuno-histochemical evaluation of tumors was performed to determine the expression levels of p16INK4A and classified tumor samples either p16+ or p16-. Using ClearSeq Comprehensive Cancer panel, we examined the distribution of somatic mutations. Somatic single-nucleotide variants (SNV) were called using GATK-Mutect2 ("tumor-only" mode) approach. Using RNA-seq, we identified a catalog of 1,044 and 8 genes as significantly expressed between p16+ and p16-, respectively at FDR 0.05 (5%) and 0.1 (10%). The clinicopathological characteristics of the patients including anatomical site, smoking and survival were analyzed when comparing p16+ and p16- tumors. The majority of tumors (65%) were p16+. Population sequence variant databases, including gnomAD, ExAC, COSMIC and dbSNP, were used to identify the mutational landscape of somatic sequence variants within sequenced genes. Hierarchical clustering of The Cancer Genome Atlas (TCGA) samples based on HPV-status was observed using differentially expressed genes. Using RNA-seq in parallel with targeted DNA-seq, we identified mutational and gene expression signatures characteristic of p16+ and p16- HNSCC. Our gene signatures are consistent with previously published data including TCGA and support the need to further explore the biologic relevance of these alterations in HNSCC.

Fig 1. Summary of somatic mutations identified.

A) Barplots showing numbers of mutations detected in the HNSCC cohort using GATK tumor-only method after each filtering step (see Methods section). MAFtools was used to generate the list of top mutated genes from 27 subjects, including p16+, p16- and p16-unknown status. B) Top 15 mutated genes. C) PIK3CA mutations identified in p16+ samples. Schematic of protein domains, displaying sites of mutations identified in the most frequently mutated gene PIK3CA in p16+ tumor samples.

Fig 2. Summary of SNV identified.

A) Individual number of SNVs identified from all 3 classes of HNSCC samples. For example, 72, 3 and 1 C>T substitutions were identified from 22 p16+, 4 p16- and 1 p16-unknown samples, respectively B) Percent contribution of individual SNVs from all 3 classes of HNSCC samples.

Fig 3. Enrichment of known oncogenic pathways based on 26 HNSCC samples including p16+ (N = 22) and p16- (N = 4) groups.

Oncogenic signaling pathways are derived from TCGA cohorts [17]. A-B) Mutation affected pathways C) 25 mutated genes of the RTK-RAS and PI3K signaling pathways (22 p16+ and 4 p16-). See S1 Fig for their known/reported drug-gene interactions and druggable categories. All gene-drug interactions and drug claims are compiled from the Drug Gene Interaction Database [16].

Fig 4. Differentially Expressed Genes (DEGs) from a cohort of HNSCC.

A) Number of DEGs in p16+ vs p16- and p16+ vs p16-unknown comparisons B-C) Heatmaps of 8 and 58 DEGs with FDR<0.01, respectively, between p16+ vs p16- and p16+ vs p16-unknown D) Heatmap generated based on TCGA [14] HNSCC data confirming the gene expression consistency. Yellow indicates reduced expression while red indicates increased expression.

Fig 5. Significantly affected pathways (p < = 0.05) based on Ingenuity Pathway Analysis (IPA).

Fig 6. Kaplan-Meier curves of the molecular subtypes of HNSCC cohort based on p16 status.