Prognostic Correlation of an Autophagy-Related Gene Signature in Patients with Head and Neck Squamous Cell Carcinoma

Abstract

Considerable evidence indicates that autophagy plays a vital role in the biological processes of various cancers. The aim of this study is to evaluate the prognostic value of autophagy-related genes in patients with head and neck squamous cell carcinoma (HNSCC). Transcriptome expression profiles and clinical data acquired from The Cancer Genome Atlas (TCGA) database were analyzed by Cox proportional hazards model and Kaplan-Meier survival analysis to screen autophagy-related prognostic genes that were significantly correlated with HNSCC patients' overall survival. Functional enrichment analyses were performed to explore biological functions of differentially expressed autophagy-related genes (ARGs) identified in HNSCC patients. Six ARGs (EGFR, HSPB8, PRKN, CDKN2A, FADD, and ITGA3) identified with significantly prognostic values for HNSCC were used to construct a risk signature that could stratify patients into the high-risk and low-risk groups. This signature demonstrated great value in predicting prognosis for HNSCC patients and was indicated as an independent prognostic factor in terms of clinicopathological characteristics (sex, age, clinical stage, histological grade, anatomic subdivision, alcohol history, smoking status, HPV status, and mutational status of the samples). The prognostic signature was also validated by data from the Gene Expression Omnibus (GEO) database and International Cancer Genome Consortium (ICGC). In conclusion, this study provides a novel autophagy-related gene signature for predicting prognosis of HNSCC patients and gives molecular insights of autophagy in HNSCC.

Figure 1

A schematic flowchart of this study.

Figure 2

Differentially expressed autophagy-related genes (ARGs). (a) The heat maps of the 38 differently expressed ARGs. The red color indicates high gene expression while the green color indicates low gene expression. N indicates nontumor tissues; T indicates tumor tissues. (b) The volcano plot for the 232 ARGs from TCGA database. Red indicates high expression, and green indicates low expression. Black indicates that those genes show no difference between HNSCC and paired nontumor tissues. (c) The boxplot of the differentially expressed ARGs.

Figure 3

The barplot and heat map of gene ontology (GO) enrichment analysis. (a) The relationship between enriched GO terms and differentially expressed autophagy-related genes. BP indicates biological process; CC indicates cellular component; MF indicates molecular function. (b) The color of each block depends on the logFC value of the specified gene.

Figure 4

The bubble plot and circle plot of Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis. (a) The relationship between KEGG pathways and differentially expressed autophagy-related genes. (b) The outer circle shows a scatter plot for each term of the logFC value of the specified gene. Red circles display upregulation, and blue means downregulation. The higher the Z-score indicates, the higher expression of the enriched pathway.

Figure 5

(a) The Kaplan–Meier plot of the overall survival (OS) for high-risk and low-risk patient cohorts divided by the 6-ARGs signature in the training set (N = 495). The OS differences are determined by the two-sided log-rank test. (b) Receiver operating characteristic analysis for the 6-ARG signature in predicting the patients of 1, 3, and 5 years OS in the training set. (c) The distribution of the risk score. Each dot represents a specific patient. The color red indicates high-risk level, and green indicates low-risk level. (d) The distribution of patients' survival status. Each dot represents a specific patient. The color red indicates dead, and green indicates alive. (e) The heat map displays the relationship between risk score and respective gene expression level of this signature in the training set.

Figure 6

Different expression of the six key genes between the high-risk group and low-risk group.

Figure 7

The correlation between the six genes included in the signature and HNSCC patients' overall survival. Kaplan–Meier plots demonstrate results from the analysis of correlation between each gene expression level and OS, all using the best separation.

Figure 8

(a) The 6-ARG signature in the cohort stratified by clinical stages. (b–e) CDKN2A expression in the cohorts stratified by sexes, ages, grades, and anatomic subdivisions. (f, g) HSPB8 expression in the cohorts stratified by anatomic sites and alcohol histories. (h, i) PRKN expression in the cohorts stratified by sexes and clinical stages. (j, k) FADD expression in the cohorts stratified by clinical stages and smoking statuses. (l–n) ITGA3 expression in the cohorts stratified by sexes, anatomic subdivisions, and smoking statuses.

Figure 9

(a) EGFR expression in the cohorts stratified by HPV status. (b) CDKN2A expression in the cohorts stratified by HPV status. (c) ITGA3 expression in the cohorts stratified by HPV status.

Figure 10

The forest plots of univariate (a) and multivariate (b) Cox analyses display the correlation of different indexes and overall survival of HNSCC patients. (c) Receiver operating characteristic analysis for risk signature and clinicopathological features in predicting HNSCC patients' OS in TCGA cohort.

Figure 11

The forest plots of univariate (a) and multivariate (b) Cox analyses display the correlation of different indexes and overall survival of HNSCC patients with available HPV status information. (c) Receiver operating characteristic analysis for risk signature and clinicopathological features in predicting HNSCC patients' OS.

Figure 12

(a) Top ten genes with the highest mutation rates in HNSCC patients from TCGA database. (b) The results of Kaplan–Meier analysis in HNSCC patients with TP53 mutation. (c) The results of Kaplan–Meier analysis in HNSCC patients without TP53 mutation. The OS differences are determined by the two-sided log-rank test.

Figure 13

The Kaplan–Meier curves of the overall survival (OS) for high-risk and low-risk patient cohorts divided by the 6-ARG signature in the GEO dataset GSE41613 (a) and ICGC dataset ORCA-IN (c). (e) The Kaplan–Meier curves of the progression-free survival (PFS) for high-risk and low-risk patient cohorts divided by the 6-ARG signature in GSE117973. Receiver operating characteristic analysis in the GEO datasets GSE41613 (b) and GSE117973 (f) and the ICGC dataset ORCA-IN (d).