Introduction

In 2020, there were approximately 840,000 newly diagnosed head and neck squamous cell (HNSCC) cases worldwide, and HNSCC ranked eighth in incidence among malignant tumors.1 HNSCC can be caused by genetic changes induced by carcinogens (such as alcohol, tobacco, and betel nut) or high-risk human papillomavirus (HPV) infection. The etiology varies significantly between regions. For example, in South and Southeast Asia, areca and tobacco irritation are the leading causes, whereas in the United States, approximately 50% of HNSCC cases are attributed to high-risk human papillomavirus (HPV) infections.

HNSCC is cancer with the highest leukocyte fractions compared to other cancers and is associated with HPV infection status.2 Over 65% of individuals diagnosed with HNSCC are likely develop to recurrent/metastatic (R/M) HNSCC 3,4. The primary therapeutic approaches for HNSCC include surgical intervention, radiotherapy, and chemotherapy.5 With the development of targeted therapy and immunotherapy, the survival time of (R/M) HNSCC patients had been significantly prolonged. Unfortunately, it is estimated that 70–80% of HNSCC patients do not respond to immunotherapy or only respond initially. 6–12 Therefore, it is very important to identify and evaluate tumor markers for predicting the benefit of immunotherapy.

Tumor cells, inflammatory cells, and chemical factors together orchestrate a sophisticated signaling network to form an inflammatory tumor microenvironment (TME). 11,13,14. Many studies have shown that tumor-related inflammation can promote angiogenesis and metastasis, and the persistent inflammatory microenvironment can also induce tumor production by triggering certain gene mutations.15 Additionally, several conventional drugs known for their anti-inflammatory effects have demonstrated protective roles against cancer. A meta-analysis indicated that the use of nonsteroidal anti-inflammatory drugs(NSAIDs) reduced the specific mortality rate by 30% and the recurrence rate by 40% in HNSCC patients.16

The present study analyzed genomic information from 913 HNSCC samples intending to identify patterns of inflammatory response in HNSCC and comprehensively assess the immune microenvironment of different patterns of inflammatory response-related genes. Then, we identified ten genes to establish the prognostic signature in the TCGA cohort. In addition, we evaluated the relationship between prognostic genes, immune microenvironment, and chemotherapy sensitivity. And lastly, the mRNA expression levels of 10 genes in adjacent tissues and HNSCC tissues were verified by qPCR and single-cell RNA sequencing.

Materials and MethodsData Collection (TCGA-HNSC Cohort and GEO Cohort)

The RNA sequencing (RNA-seq) data, and corresponding clinical information of 546 cases of HNSCC were downloaded from https://portal.gdc.cancer.gov/. In addition, the RNA-seq data and clinical information data of 367 tumor samples were obtained from the Gene Expression Omnibus (GEO) (GSE41613 and GSE65858). The expression profile (FPKM value) of TCGA-HNSCC dataset was converted to transcripts per million (TPMs), which was the same as the microarray results. ComBat algorithm was used to eliminate batch effects by “sva” package, and three datasets were combined into one dataset.17

Consensus Clustering for Inflammatory Response-Related Genes

200 inflammatory response-related genes (IRGs) were identified in the Molecular Signatures database. The k-means algorithm was used to perform unsupervised clustering of 200 IRGs, which performed using the “ConsensusClusterPlus” package.18 Principal component analysis (PCA) was employed to investigate the distributions of various groups.19

Gene Set Variation Analysis (GSVA) and Functional Annotation

GSVA is a nonparametric, unsupervised algorithm that was employed to access the activity and pathway changes of biological processes in various clusters.20 This process was performed using “GSVA” and “clusterProfiler” packages.21

Tumor Microenvironment Analysis

By using single-sample gene set enrichment analysis (ssGSEA), we were able to analyze the relative abundance of immune infiltration within the TME, which encompassed 16 distinct immune cell subtypes and 13 immune cell functions.22,23 In this process, the packages “GSVA” and “limma” packages were used.

Establishment and Evaluation of an IRGs Signature

Within the TCGA-HNSCC dataset, the “limma” package was applied to identify differentially expressed genes (DEGs) between malignant and non-malignant tissues. The screening criteria were an absolute value of the fold change greater than 1 and a false discovery rate (FDR) was ≤0.05. To identify inflammatory response genes with potential prognostic implications, univariate Cox regression analysis was conducted. Subsequently, a prognostic model was formulated using LASSO-penalized Cox regression analysis.24 The risk score of each patient was calculated as following: risk score , where “C” is the coefficient of gene and “Exp” is the gene expression level. Based on the median risk score, the HNSCC patients were stratified into low- and high-risk groups. Kaplan-Meier (K-M) survival curves were generated to contrast the overall survival (OS). Additionally, time-dependent receiver operating characteristic (ROC) curve analysis was performed, and the relationship between the risk model and clinicopathological characteristics was evaluated. This process was performed using the R survival, pHeatmap, glmnet, timeROC, and ggplot2 packages.

Evaluation of the Predictive Function of the Risk Model for Common Chemotherapeutic Drugs

In order to access the clinical application of the risk model, we calculated the half maximal inhibitory concentration (IC50) values of common chemotherapeutic drugs for the TCGA-HNSCC dataset. The differences in IC50 among two group were calculated using the Wilcoxon signed rank test. This process was undertaken using the “pRRophetic” and “ggplot2” packages.25

qRT–PCR

To validate the mRNA expression levels of prognostic genes in HNSCC tissues, we collected five sets of matched malignant and non-malignant tissues from the Stomatological Hospital of Shandong University in China (Table S1). The study protocol was approved by the Human Research Ethics Committee. qRT–PCR was employed to detect the expression levels of 10 inflammatory response-related mRNAs in the samples, malignant and non-malignant tissues.

scRNA‑seq Data Processing and Analysis

scRNA-seq samples were downloaded from GSE181919, including the tissue types of normal tissue (NL, n=9) and primary cancer (CA, n=20). For normalizing and analyzing the single-cell data, we used the Seurat (version 3.1.5) package (https://satijalab.org/seurat/). Dimension reduction and cell cluster visualization were performed using t-DistributedStochastic Neighbor Embedding (t-SNE). Using Seurat’s FindAllMarkers function, unique cluster-specific markers genes were identified with the cutoff values of |log2 fold change (FC)| > 1 and adjusted p-value < 0.01. The crosstalk between IRGs was explored using CellChat and the patterns of incoming and outgoing signals were identified.

Protein-Protein Interaction and Molecular Docking

We conducted a protein-protein interaction (PPI) analysis utilizing the STRING database, which can be accessed at http://string-db.org. Download the structures of paraplatin, docetaxel, fluorouracil and cetuximab small molecules from the Pubchem database. Search the Uniprot database for CCR7 proteins, species human. Protein preparation was performed using the Protein Preparation Wizard module of Schrödinger software. The protein was prepared by assigning bond levels, hydrogenation, and completing missing side chains. Molecular docking was performed using the Glide module in Schrödinger.

ResultsIdentification of HNSCC Patterns Based on IRGs

The gene-expression profiles of 200 IRGs from TCGA-HNSCC and GEO (GSE41613, GSE65858) were used for a consensus clustering analysis of HNSCC (Figure 1A). The 913 HNSCC patients were divided into three patterns: type A (n=218), type B (n=379), and type C (n=269). We then compared the prognosis of each subgroup. The results showed that subgroup C had a significant survival advantage (p=0.002, Figure 1B). PCA indicated that patients categorized into distinct subgroups exhibited three separate clustering patterns. (Figure 1C). Figure 1D shows that the three patterns of HNSCC patients had different clinicopathological characteristics, in which stage and survival showed a remarkable difference (p < 0.01). In brief, type C was associated with a favorable prognosis.

Figure 1 Identification of HNSCC subtypes based on IRGs. (A) The k-means algorithm clustering using 200 IRGs. (B) K-M curves for the three inflammatory-related molecular patterns of HNSCC patients in the TCGA-HNSCC and GEO cohorts. (C) Confirming the discrimination of k = 3 using PCA analysis, cluster A (red), cluster B (blue), and cluster C (Orange). (D) The heat map depicts the correlation between the pattern and different clinicopathological characteristics.

Characteristics of Cell Infiltration in the TME in the Different Inflammatory Response Patterns

GSVA was used to explore the differences between the biological behaviors of different inflammatory response patterns. As shown in Figure 2A and B, the type A inflammatory response was characterized by enriched tumor metabolism; the type B inflammatory response was characterized by enrichment of the oncogenic activation of the Janus kinase/signal transducers; and the type C inflammatory response was characterized by enrichment of immune response-related pathways (Figure 2A–C).

Figure 2 Characteristics of cell infiltration in the TME in the different inflammatory response patterns. (A–C) GSVA enrichment analysis showing the activation states of KEGG in distinct inflammatory response subtypes. In the heat map, red represents activated pathways, and blue represents inhibited pathways. A cluster A vs cluster B; B cluster B vs cluster C; C cluster A vs cluster C (D, E) Box plots showed the scores of immune infiltrations (D) and immune functions (E) among the inflammatory response patterns. ***p<0.001. (F–H) The comparison of the expression of PD-1, PD-L1 and CTLA4 in the different inflammatory response patterns. ***p<0.001.

Subsequently, ssGSEA data showed that there were significant differences among the three inflammatory response patterns, which was consistent with the results of the GSVA enrichment analysis; that is, Type C had the most infiltrated immune cells and the most active immune functions (Figure 2D and E). The expression levels of CTLA4 and PD1/PD-L1 in immune response type C were the highest, followed by the immune response types B and A, which showed sequential decreases in expression levels (Figure 2F–H). In summary, the consistency of the prognosis of different immune response patterns with their TME suggests that a classification method based on inflammatory response is reasonable.

Establishment and Assessment of an Inflammatory Response-Related Signature

From 200 IRGs, 18 prognostic-related genes were identified by univariate Cox analysis; 143 IRGs were differentially expressed in tumor tissues and paracancerous tissues (p<0.01, logFC>1). Figure 3A). Finally, 14 IRGs were selected as candidate prognostic indicators. Figures 3B delineated the expression level of the 14 genes in tumors and paracancerous tissues. Figures 3C showed the prognosis characteristics associated with 14 genes. For instance, the risk ratio of ADGRE1, ITGA5, OLR1, PVR, SERPINE1 and SRI were different from the others, indicating that elevated expression of these genes was associated with better prognosis. In addition, Figure 3D shows the correlations of the 14 genes.

Figure 3 Identification of IRGs in the TCGA-HNSCCC cohort. (A) Venn diagram shows that 14 genes were identify the prognosis-related DEGs. (B) The heat map depicts 14 overlapping genes expression levels in HNSCC tissues and adjacent normal tissues. (C) Forest plots shows the hazard ratio and 95% confidence intervals (95% CI) of 14 prognosis-related DEGs by univariate Cox regression analysis. (D) The correlation network of 14 overlapping genes.

The prognostic characteristics of 10 genes were determined by LASSO-Cox regression. Subsequently, the TCGA-HNSCC cohort was stratified into low- and high-risk groups. Figure 4A–C shows the heat map of risk value, survival status, and gene expression for each patient in the TCGA-HNSCC cohort. The results suggested that the clinical pathological staging of HNSCC patients, such as T-stage, Stage, and survival status, was significantly correlated with the risk score (Figure 4C). The K-M curve indicated a much more favorable prognosis for the patients in the low-risk group (Figure 4D). PCA suggested that the patients in different subgroups were distributed in different directions (Figure 4E). The 1, 2, and 3 years of AUC rates were 0.703, 0.724, and 0.742, respectively (Figure 4F). The 3-year AUC curve was also compared with other clinicopathological characteristics (Figure 4G). The two GEO cohort (GSE41613 and GSE65858) were used to verify the predictive effect of the risk score (Figure 4H–K). Both low-risk groups in the two datasets exhibited significantly better prognosis than the high-risk groups. And the AUC analysis at one, two, and three years showed that the risk score had good predictive efficiencies, which indicated that IRGs risk model had a strong performance for predicting overall survival.

Figure 4 Construction and assessment of ten inflammatory response-related signature. (A)The risk score distribution of HNSCC patients. (B) Survival time and survival status of HNSCC patients. (C) The heat map depicts the expression of inflammatory response-related and the correlation between the risk group and different clinicopathological characteristics. (D) K-M curves were used to analyze OS in the low and high-risk groups in the TCGA cohort. (E) Principal component analysis (PCA) of the inflammatory response-related signature in TCGA cohort. (F) The 1-year, 2-year, and 3-year ROC curves with AUC values of 0.671, 0.684, and 0.710. (G) A comparison of 3-year AUC values with different clinicopathological characteristics: age, sex, grade, T, N, and M. (H–K) K-M curves were used to analyze OS in the low and high-risk groups in 2 GEO cohort(H), GSE41613; (J), GSE65858. AUC time-dependent ROC curves for OS (I), GSE41613; (K), GSE65858.

Next, the patients were categorized into distinct subgroups based on their diverse clinicopathological characteristics, including age, sex, T stage, HPV status, smoking history, anatomical location (floor of mouth, tonsil, tongue, larynx, hypopharynx, and base of tongue), and N stage. The results of the K-M curve suggested that the low-risk group was significantly prolonged OS in the presence of different clinicopathological characteristics, which further suggested that our risk model was dependable and accurate (Figure 5 and Figure S1).

Figure 5 K-M analysis of overall survival in different subgroups. Patients were classified according to age (age <= 65 and age >65), sex (female and male), T (T1-T2, T3-T4), N (N0, N1-N3).

Different Immune Infiltration Characteristics and Enrichment Pathways in the High- and Low-Risk Groups

GSEA using the TCGA data for the hallmark genome suggest that the enrichment of cancer-related pathways such as glycolysis, MYC, and epithelial-mesenchymal transition (EMT) in the high-risk group may be associated with a poor prognosis (Figure 6A–B). The infiltrating immune-cells and immune-cells functions were evaluated by ssGSEA. Our analysis revealed that the abundance of the majority of immune cell types was considerably greater in the low-risk cohort compared to the high-risk group, and the low-risk group’s immune checkpoint, T-cell co-stimulation and inflammation promoting were also relatively high (Figure 6C–D). In summary, the inflammatory response score exhibited a strong correlation with the TME of HNSCC patients.

Figure 6 Different immune infiltration characteristics and enrichment pathways in the high- and low-risk groups. (A and B) GSVA enrichment analysis showing the activation states in the high and low-risk group. In the heat map, red represents activated pathways, and blue represents inhibited pathways. (A) HALLMARK; (B) KEGG. (C and D) Box plots showed the scores of immune infiltrations (C) and immune functions (D) among the high and low-risk groups. Ns: not significant; *p<0.05; **p<0.01; ***p<0.001.

The Role of the Inflammatory Response in Predicting Responses to Common Chemotherapeutic Drugs and Immune Responses

We attempted to determine the relationship between the efficacy of common chemotherapeutic drugs and risk scores. We found that low-risk scores were associated with a low IC50 for methotrexate and a high IC50 for docetaxel. These results indicate that the risk score derived from IRGs could serve as a potential predictor of sensitivity to common chemotherapeutic drugs (Figure 7A–D). Additionally, we assessed the potential correlation between the risk model and the expression level of immune checkpoint inhibitor (ICI) markers. The results suggested that a low risk score was positively correlated with high expression of CTLA4 and PDCD1 (p<0.001). However, there were no correlation between the expression level of CD274 and the risk score (p=0.19) (Figure 7E–J). As depicted in Figure 7K, the IPS-CTLA4-Positive/PD1-Negative scores and IPS-PD1-Negative/CTLA4-Negative scores, were elevated in the low-risk group. This suggests that individuals classified as low-risk may exhibit a more favorable response to immunotherapeutic treatments.

Figure 7 The role of the inflammatory response in predicting responses to common chemotherapeutic drugs and immune responses. (A–D) The model acted as a potential predictor for common chemotherapeutic drugs. The low-risk scores were associated with a low IC50 for methotrexate and a high IC50 for docetaxel. (E–G) The comparison of the expression of PD-1, PD-L1, and CTLA4 in the different risk scores. (H–J) Scatter plot showing the correlation between risk score and CTLA4, PD-1, and PD-L1 expression. (K) IPS score reflecting the response to ICIs.

The mRNA Expression of IRGs Between HNSCC Tissues and Adjacent Normal Tissues

To confirm the variation in expression levels of 10 IRGs genes (CCL22, CCR7, ADGRE1, F3, HBEGF, ITGA5, KCNA3, LCK, OLR1, and SRI) between HNSCC tissue and adjacent non-tumor tissues, qRT-PCR was conducted to measure mRNA levels. The results showed that compared with neighboring non-tumor tissues, 10 IRGs genes were all highly expressed in HNSCC tissues. The qRT-PCR results are consistent with the expression of 10 IRGs genes in the TCGA-HNSCC (Figure 8).

Figure 8 The mRNA expression of IRGs between HNSCC tissues and adjacent normal tissues. The mRNA expression analysis by qRT-RCR.

The Docking Conformation and Interaction Force Analysis of 10 IRGs Genes

In our study, we employed molecular docking techniques to investigate the function of SCC7 in chemotherapy. Initially, we conducted a PPI analysis among 10 inflammation-related genes. From this analysis, SCC7 emerged as a central hub gene due to its highest connectivity. Next, we molecularly docked SCC7 with each of the common chemotherapeutic agents (paraplatin, docetaxel, fluorouracil, and cetuximab). Figure 9 demonstrates the docking conformation and interaction force analysis between SCC7 and four common chemotherapeutic agents. The results showed that the binding energies formed by the four common chemotherapeutic agents (paraplatin, docetaxel, fluorouracil, and cetuximab) with SCC7 amino acid residues were −8.880 kcal/mol, −8.241 kcal/mol, −8.153 kcal/mol, −9.79 kcal/mol, respectively (Figure 9).

Figure 9 The docking conformation and interaction force analysis between SCC7 and paraplatin, docetaxel, fluorouracil and cetuximab. Color symbols: yellow sticks for drug molecules, blue sticks for amino acid residues, blue lines for hydrogen bonding, and gray dashed lines for hydrophobic interaction.

scRNA-Seq Analysis of Inflammatory Response-Related Gene in HNSCC

We performed the “ScaleData” function to scale all genes extracted from the scRNA-seq dataset GSE181919. Finally, 19 clusters were found in NL and CA (Figure 10A). We screened the cell markers by the “FindAllMarkers” function (logFc = 0.5, Minpct −0.35). Marker genes were used to label the identified clusters as different cell types (Figure 10A). Following this, the features of the nine markers in different clusters were displayed using violin and UMAP plots in the sc-RNA sequencing profile (Figure S2 and S3). For example, CCR7 was specifically expressed in T cells, Plasma cells, Macrophages, and Dendritic cells and was significantly higher in CA tissues than in NL tissues.

Figure 10 scRNA-seq analysis of IRGs in HNSCC. (A) UMAP clustering colored by groups in normal tissue (NL) and primary cancer (CA). (B and C) Comparisons of overall changes in cell-cell communication between NL and CA, including the number or weight of interactions (left) in NL, number or weight of interactions (medium) in CA and differential interaction strength (right) between CA and NL. (D) Differences in the overall signaling pathway between NL and CA. (E and F) The heatmap demonstrated the incoming signaling patterns and outgoing f signaling patterns in CA. (G) Communication probabilities of important ligand-receptor pairs in CA.

CellChat was employed to conduct a comprehensive assessment of the interactions between NL and CA, with a particular focus on the number and weight of cell communications (Figure 10B–C). A summary of communication probabilities in the information network was conducted to facilitate a comparison of the overall information flow between NL and CA (Figure 10D). The results revealed that inflammatory response-related signalling pathways, WNT, IL-6, IL-10, IL16, and TGFb, were more abundant in CA (blue). The heatmap demonstrated an increase in the incoming signaling patterns from T cells and macrophage in CA. The most significant increase in outgoing signaling patterns was observed in the endothelial cells and fibroblast, which LCK signaling was most significant outgoing signaling patterns in T cell (Figure 10E–G).

Discussion

The inflammatory response affects the entire process of tumor development and treatment.26 In fact, Rudolf Virchow proposed the hypothesis that tumors originate from chronic inflammation as early as the 19th century. In this study, we examined inflammatory response-related genes in HNSCC samples and stratified the risk levels of HNSCC patients. To our knowledge, this is the first time that bioinformatics methods have been used to comprehensively evaluate the inflammatory response mode and biological functions in HNSCC, with special attention to its immune characteristics.

According to the profiles of IRGs genes expression, HNSCC was divided into three patterns with different prognoses, immune microenvironments, and enrichment pathways. Specifically, the OS of type C was significantly longer than that of the other types. The GSEA results showed that type C had enriched immune response-related pathways, including primary immunodeficiency, T cell receptor signaling pathway, and intestinal immune network. Therefore, we hypothesized that the inflammatory response-based patterns had different immune microenvironments. The ssGSVA results suggested that type C had significantly more infiltrated immune cells, including B cells and CD8+ T, than the other patterns. Additionally, immune function was more active in type C, and this subtype showed higher expression of ICIs genes, such as PD-1, PD-L1 and CTLA4, than the other two patterns.

Although the immunotherapy strategies represented by ICIs have achieved good clinical effects, their response rates are low, and primary or secondary drug resistance can develop during treatment.27,28 Many studies have shown that the current “dilemma” of immunotherapy is inseparable from tumor-related inflammation.29 Currently, the anti-inflammatory drugs and inflammation-targeted drugs used in clinical practice mainly include anti-infective drugs, nonsteroidal anti-inflammatory drugs, metformin, statin, natural supplements, and cytokine-specific drugs. In animal experiments, the preoperative administration of the nonsteroidal anti-inflammatory drugs ketorolac and/or resolvin inhibited micrometastasis in a variety of tumor resection models and could prolong survival time.30 Jayaprakash et al pointed out that aspirin usage is correlated with a lower incidence of HNSCC: aspirin users had a 25% reduction in the risk of HNSCC.31 In addition, Becker et al showed that patients who took more than 7 tablets of aspirin per week had a significantly reduced risk of HNSCC (approximately 32%).32 Regarding the effect of anti-inflammatory drugs on the prognosis of HNSCC patients, Macfarlane et al reported that the use of nonsteroidal anti-inflammatory drugs was associated with increased OS in HNSCC patients in a dose-dependent manner.33

We used 10 IRGs genes to construct prognostic characteristics. Some of these genes have been reported in HNSCC. One study showed that high expression of CCR7 was associated with prolonged disease-free survival.34 Our study showed that the four commonly used chemotherapeutic drugs (cisplatin, docetaxel, fluorouracil, and cetuximab) formed strong binding affinities with the amino acid residues of SCC7. HBEGF may activate tumor growth, invasion, and metastasis through EGFR.35

Currently, ideal and reliable ICI-related biomarkers for predicting efficacy remain to be explored. In patients with HNSCC, higher CD8+ T cell infiltration was observed in responders to PD1/PD-L1 inhibitor treatment.27,36 A recent meta-analysis indicated that CD8+ cells could serve as a standalone predictor of the prognosis of HNSCC and a therapeutic target.16 In our study, the immune-related pathways in inflammatory response subtype C and the low-risk group were activated and had the advantage of immune infiltration, and the OS was better, which is consistent with the previously described phenomenon.

Our findings indicated that patients classified as high-risk exhibited a greater sensitivity to docetaxel, whereas those in the low-risk showed a higher sensitivity to methotrexate. This indicates that the risk score based on the IRGs genes can be serve as a potential predictor of sensitivity to common chemotherapeutic drugs.

Summary

In summary, our study has delineated three distinct inflammatory response patterns in HNSCC, each with unique prognostic and immunological features. Moreover, we have formulated a risk score that serves as a predictive tool for gauging the sensitivity of chemotherapy and their potential response to immunotherapy. This research underscores the pivotal role of genes associated with the inflammatory response within the tumor microenvironment (TME) of HNSCC, offering a valuable framework for the development of comprehensive and personalized treatment approaches for this disease.

Data Sharing Statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Ethics Statement

The studies involving human participants were reviewed and approved by The Ethics Committee of Stomatological Hospital of Shandong University (NO.20210227). The patients/participants provided their written informed consent to participate in this study. Our study complies with the Declaration of Helsinki. No animal studies are presented in this manuscript. No potentially identifiable human images or data is presented in this study.

Acknowledgment

Yong Zhu and Yongzhe Zhang share the first authorship.

Funding

This work was supported by the National Natural Science Foundation of China (82201092), and the Natural Science Foundation of Shandong Province (ZR2021MH086).

Disclosure

The authors report no conflicts of interest in this work.

References

1. Muzaffar J, Bari S, Kirtane K, Chung C. recent advances and future directions in clinical management of head and neck squamous cell carcinoma. Cancers. 2021;13(2):338. doi:10.3390/cancers13020338

2. Thorsson V, Gibbs D, Brown S, et al. The immune landscape of cancer. Immu. 2018;48:812–830. doi:10.1016/j.immuni.2018.03.023

3. Kitamura N, Sento S, Yoshizawa Y, Sasabe E, Kudo Y, Yamamoto T. Current trends and future prospects of molecular targeted therapy in head and neck squamous cell carcinoma. Int J Mol Sci. 2020;22(1):240. doi:10.3390/ijms22010240

4. Fulcher C, Haigentz M, Ow T. AHNS series: Do you know your guidelines? Principles of treatment for locally advanced or unresectable head and neck squamous cell carcinoma. Head Neck. 2018;40:676–686. doi:10.1002/hed.25025

5. Dlamini Z, Alaouna M, Mbatha S, et al. Genetic drivers of head and neck squamous cell carcinoma: Aberrant splicing events, mutational burden, hpv infection and future targets. Genes. 2021;12(3):422. doi:10.3390/genes12030422

6. Jung K, Narwal M, Min S, Keam B, Kang H. Squamous cell carcinoma of head and neck: What internists should know. Korean j Intern Med. 2020;35(5):1031–1044. doi:10.3904/kjim.2020.078

7. Ebnoether E, Muller L. Diagnostic and therapeutic applications of exosomes in cancer with a special focus on head and neck squamous cell carcinoma (HNSCC). Int J Mol Sci. 2020;21(12):4344. doi:10.3390/ijms21124344

8. Patel J, Levy D, Nguyen S, Knochelmann H, Day T. Impact of PD-L1 expression and human papillomavirus status in anti-PD1/PDL1 immunotherapy for head and neck squamous cell carcinoma-Systematic review and meta-analysis. Head Neck. 2020;42:774–786. doi:10.1002/hed.26036

9. Cohen EEW, Bell RB, Bifulco CB, et al. The society for immunotherapy of cancer consensus statement on immunotherapy for the treatment of squamous cell carcinoma of the head and neck (HNSCC). J ImmunoTher Cancer. 2019;7:184. doi:10.1186/s40425-019-0662-5

10. Hellmann M, Nathanson T, Rizvi H, et al. Genomic features of response to combination immunotherapy in patients with advanced non-small-cell lung cancer. Cancer Cell. 2018;33:843–852. doi:10.1016/j.ccell.2018.03.018

11. Mlecnik B, Bindea G, Kirilovsky A, et al. The tumor microenvironment and immunoscore are critical determinants of dissemination to distant metastasis. Sci trans med. 2016;8:327ra26. doi:10.1126/scitranslmed.aad6352

12. Ferris R, Blumenschein G, Fayette J, et al. Nivolumab for recurrent squamous-cell carcinoma of the head and neck. New Engl J Med. 2016;375:1856–1867. doi:10.1056/NEJMoa1602252

13. Lin Z, Xu Q, Miao D, Yu F. An inflammatory response-related gene signature can impact the immune status and predict the prognosis of hepatocellular carcinoma. Front Oncol. 2021;11:644416. doi:10.3389/fonc.2021.644416

14. Lan M, Lu W, Zou T, et al. Role of inflammatory microenvironment: potential implications for improved breast cancer nano-targeted therapy. Cel mol life sci. 2021;78:2105–2129. doi:10.1007/s00018-020-03696-4

15. Bonomi M, Patsias A, Posner M, Sikora A. The role of inflammation in head and neck cancer. Adv Exp Med Biol. 2014;816:107–127.

16. Huang Z, Xie N, Liu H, et al. The prognostic role of tumour-infiltrating lymphocytes in oral squamous cell carcinoma: A meta-analysis. J Oral Pathol Med. 2019;48:788–798. doi:10.1111/jop.12927

17. Teschendorff A, Zhuang J, Widschwendter M. Independent surrogate variable analysis to deconvolve confounding factors in large-scale microarray profiling studies. Bioinfo. 2011;27:1496–1505. doi:10.1093/bioinformatics/btr171

18. Wilkerson M, Hayes D. ConsensusClusterPlus: A class discovery tool with confidence assessments and item tracking. Bioinfo. 2010;26:1572–1573. doi:10.1093/bioinformatics/btq170

19. Raychaudhuri S, Stuart J, Altman R. Principal components analysis to summarize microarray experiments: application to sporulation time series. Pac Symp Biocomp. 2000;455–466. doi:10.1142/9789814447331_0043

20. Hänzelmann S, Castelo R, Guinney J. GSVA: Gene set variation analysis for microarray and RNA-seq data. BMC Bioinf. 2013;14:7. doi:10.1186/1471-2105-14-7

21. Wu T, Hu E, Xu S, et al. clusterProfiler 4.0: A universal enrichment tool for interpreting omics data. Innov. 2021;2:100141.

22. Barbie D, Tamayo P, Boehm J, et al. Systematic RNA interference reveals that oncogenic KRAS-driven cancers require TBK1. Nature. 2009;462:108–112. doi:10.1038/nature08460

23. Bindea G, Mlecnik B, Tosolini M, et al. Spatiotemporal dynamics of intratumoral immune cells reveal the immune landscape in human cancer. Immu. 2013;39:782–795. doi:10.1016/j.immuni.2013.10.003

24. Wang H, Lengerich B, Aragam B, Xing E. Precision lasso: Accounting for correlations and linear dependencies in high-dimensional genomic data. Bioinfo. 2019;35:1181–1187. doi:10.1093/bioinformatics/bty750

25. Geeleher P, Cox N, Huang R. pRRophetic: An R package for prediction of clinical chemotherapeutic response from tumor gene expression levels. PLoS One. 2014;9:e107468. doi:10.1371/journal.pone.0107468

26. Dolan R, Laird B, Horgan P, McMillan D. The prognostic value of the systemic inflammatory response in randomised clinical trials in cancer: A systematic review. Crit Rev Oncol/Hematol. 2018;132:130–137. doi:10.1016/j.critrevonc.2018.09.016

27. Tang L, Liu T, Chen J, Dang J, Li G. Immune-checkpoint inhibitors versus other systemic therapies in advanced head and neck cancer: A network meta-analysis. Immunother. 2021;13:541–555. doi:10.2217/imt-2020-0070

28. Okada T, Okamoto I, Sato H, Ito T, Miyake K, Tsukahara K. Efficacy and safety of paclitaxel combined with cetuximab for head and neck squamous cell carcinoma. In vivo. 2021;35:1253–1259. doi:10.21873/invivo.12376

29. Hou J, Karin M, Sun B. Targeting cancer-promoting inflammation - have anti-inflammatory therapies come of age? Nature reviews. Clin Oncol. 2021;18:261–279. doi:10.1038/s41571-020-00459-9

30. Panigrahy D, Gartung A, Yang J, et al. Preoperative stimulation of resolution and inflammation blockade eradicates micrometastases. J Clin Invest. 2019;129:2964–2979. doi:10.1172/JCI127282

31. Jayaprakash V, Rigual N, Moysich K, et al. Chemoprevention of head and neck cancer with aspirin: A case-control study. Arch Otolaryngol. 2006;132:1231–1236. doi:10.1001/archotol.132.11.1231

32. Becker C, Wilson J, Jick S, Meier C. Non-steroidal anti-inflammatory drugs and the risk of head and neck cancer: A case-control analysis. Int j Cancer. 2015;137:2424–2431. doi:10.1002/ijc.29601

33. Macfarlane T, Murchie P, Watson M. Aspirin and other non-steroidal anti-inflammatory drug prescriptions and survival after the diagnosis of head and neck and oesophageal cancer. Cancer Epidemiol. 2015;39:1015–1022. doi:10.1016/j.canep.2015.10.030

34. González-Arriagada W, Lozano-Burgos C, Zúñiga-Moreta R, González-Díaz P, Coletta R. Clinicopathological significance of chemokine receptor (CCR1, CCR3, CCR4, CCR5, CCR7 and CXCR4) expression in head and neck squamous cell carcinomas. J Oral Pathol Med. 2018;47:755–763. doi:10.1111/jop.12736

35. Gao J, Ulekleiv C, Halstensen T. Epidermal growth factor (EGF) receptor-ligand based molecular staging predicts prognosis in head and neck squamous cell carcinoma partly due to deregulated EGF- induced amphiregulin expression. J Exp Clinic Cancer Res. 2016;35:151. doi:10.1186/s13046-016-0422-z

36. Davis A, Patel V. The role of PD-L1 expression as a predictive biomarker: An analysis of all US food and drug administration (FDA) approvals of immune checkpoint inhibitors. J ImmunoTher Cancer. 2019;7(1):278. doi:10.1186/s40425-019-0768-9