Survival-related indicators ALOX12B and SPRR1A are associated with DNA damage repair and tumor microenvironment status in HPV 16-negative head and neck squamous cell carcinoma patients

Abstract

Objectives: To investigate prognostic-related gene signature based on DNA damage repair and tumor microenvironment statue in human papillomavirus 16 negative (HPV16-) head and neck squamous cell carcinoma (HNSCC).

Methods: For the RNA-sequence matrix in HPV16- HNSCC in the Cancer Genome Atlas (TCGA) cohort, the DNA damage response (DDR) and tumor microenvironment (TM) status of each patient sample was estimated by using the ssGSEA algorithm. Through bioinformatics analysis in DDR_high/TM_high (n = 311) and DDR_high/TM_low (n = 53) groups, a survival-related gene signature was selected in the TCGA cohort. Two independent external validation cohorts (GSE65858 (n = 210) and GSE41613 (n = 97)) with HPV16- HNSCC patients validated the gene signature. Correlations among the clinical-related hub differentially expressed genes (DEGs) and infiltrated immunocytes were explored with the TIMER2.0 server. Drug screening based on hub DEGs was performed using the CellMiner and GSCALite databases. The loss-of-function studies were used to evaluate the effect of screened survival-related gene on the motility of HPV- HNSCC cells in vitro.

Results: A high DDR level (P = 0.025) and low TM score (P = 0.012) were independent risk factors for HPV16- HNSCC. Downregulated expression of ALOX12B or SPRR1A was associated with poor survival rate and advanced cancer stages. The pathway enrichment analysis showed the DDR_high/TM_low samples were enriched in glycosphingolipid biosynthesis-lacto and neolacto series, glutathione metabolism, platinum drug resistance, and ferroptosis pathways, while the DDR_high/TM_low samples were enriched in Th17 cell differentiation, Neutrophil extracellular trap formation, PD - L1 expression and PD - 1 checkpoint pathway in cancer. Notably, the expression of ALOX12B and SPRR1A were negatively correlated with cancer-associated fibroblasts (CAFs) infiltration and CAFs downstream effectors. Sensitivity to specific chemotherapy regimens can be derived from gene expressions. In addition, ALOX12B and SPRR1A expression was associated with the mRNA expression of insulin like growth factor 1 receptor (IGF1R), AKT serine/threonine kinase 1 (AKT1), mammalian target of rapamycin (MTOR), and eukaryotic translation initiation factor 4E binding protein 1 (EIF4EBP1) in HPV negative HNSCC. Down-regulation of ALOX12B promoted HPV- HNSCC cells migration and invasion in vitro.

Conclusions: ALOX12B and SPRR1A served as a gene signature for overall survival in HPV16- HNSCC patients, and correlated with the amount of infiltrated CAFs. The specific drug pattern was determined by the gene signature.

Keywords: Cancer-associated fibroblasts; DNA damage response; Head and neck squamous cell carcinoma; Human papilloma virus 16-negative; Tumor microenvironment.

Fig. 1

Patient sample stratification based on DNA damage response (DDR) level in human papillomavirus 16 negative (HPV16-) head and neck squamous cell carcinoma (HNSCC) TCGA discovery cohort (n = 433). A ssGSEA matrix plot of two subtypes identified in HPV16- HNSCC TCGA discovery cohort according to five DDR-associated genesets. B Kaplan-Meier (K-M) plot of overall survival probability for patients in the above two subtypes. One patient belonging to DDR_high subtype lacked survival data as shown in risk table

Fig. 2

The Kaplan–Meier (K-M) overall survival (OS) curve of patients in the TCGA final discovery cohort (n = 364). A Two clusters were obtained from TCGA discovery cohort (n = 433) by dichotomizing tumor microenvironment (TM) score (“ESTIMATE” package in R). B K-M plot of OS probability in high- and low- TM score group. C K-M plot of OS probability in DDR_high/TM_low, DDR_high/TM_high, and DDR_low/TM_high groups. D K-M plot of OS probability for patients in the TCGA final discovery cohort including only DDR_high/TM_high and DDR_high/TM_low groups. One patient in the TCGA final discovery cohort lacked survival data as mentioned in Fig. 1

Fig. 3

Somatic mutant genes in DDR_high/TM_low and DDR_high/TM_high groups of the TCGA final discovery cohort. A Mutant genes (top 20) in DDR_high/TM_low group. B Mutant genes (top 20) in DDR_high/TM_high group. C Differential mutant genes between DDR_high/TM_low and DDR_high/TM_high groups

Fig. 4

Identification of DDR- and TM- related hub genes in the TCGA final discovery cohort. A Heatmap of clustered top 100 differentially expressed genes (DEGs) between DDR_high/TM_high and DDR_high/TM_low groups. B-D WGCNA algorithm screened out seven eigengenes module (blue, green, red, pink, black, brown, and grey) based on DDR and TM status. E-F The top 50 hub genes selected in black and brown modules by protein-to-protein network method; respectively

Fig. 5

Biological characteristics of hub genes in black and brown modules. A Biological process, molecular function, and cellular component terms enriched by hub genes in black module. B KEGG pathway enriched by hub genes in black module. C Biological process, molecular function, and cellular component terms enriched by hub genes in brown module. D KEGG pathways enriched by hub genes in brown module

Fig. 6

KEGG pathways analysis of hub genes in blue, green, and red modules. A KEGG pathways enriched by top 100 hub genes in blue module. B KEGG pathways enriched by top 100 hub genes in green module. C KEGG pathways enriched by top 100 hub genes in red module

Fig. 7

Wilcoxon rank test of prognostic-related hub gene expression between early and advanced TNM stage groups. A-C Correlation between the expression of ALOX12B and N stage in the TCGA final discovery cohort and GSE65858 data set, respectively. Association between the expression of ALOX12B and TNM stage in GSE41613 data set. D-F Correlation between the expression of SPRR1A and N stage in the TCGA final discovery cohort and GSE65858 data set, respectively. Association between the expression of SPRR1A and TNM stage in GSE41613 data set

Fig. 8

Immunogenicity analysis based on the expression of prognostic-related hub genes. Differences of IPS, IPS-PD1/PD-L1/PD-L2 blocker, IPS-CTLA4 blocker, and IPS-PD1/PD-L1/PD-L2 + CTLA4 blocker between the groups with high and low expression of ALOX12B (A-D) and SPRR1A (E-H), respectively. I Drug patterns based on ALOX12B and SPRR1A in CellMiner database, and panels were arranged according to correlation coefficient

Fig. 9

Association between the mRNA expression of ALOX12B (A-D) or SPRR1A (E-H) and the expression of insulin like growth factor 1 receptor (IGF1R), AKT serine/threonine kinase 1 (AKT1), mammalian target of rapamycin (MTOR), and eukaryotic translation initiation factor 4E binding protein 1 (EIF4EBP1) in HPV- HNSCC samples at TIMER2.0 online database

Fig. 10

ALOX12B suppresses the migration and metastasis of HN6 cells in vitro. A Quantification mRNA expression of ALOX12B in HN6 and CAL27 cell lines by qRT-PCR. B HN6 cells were transfected with siNC and siALOX12B, and the expression of ALOX12B was examined using qRT-PCR. C CCK8 assay was conducted to assess cell viability in HN6-siNC and HN6-siALOX12B cells. D-E Transwell assays were performed to examine the invasion and migration of HN6-siNC and HN6-siALOX12B cells. Column, mean; Error bars, S.D.; *p < 0.05; **p < 0.01; ALOX12B, Arachidonate 12-Lipoxygenase, 12R Type; qRT-PCR, quantitative real-time polymerase chain reaction