Profiling the Tumor Immune Microenvironment of HPV-Associated Base of Tongue Squamous Cell Carcinoma

Introduction

Head and neck squamous cell carcinomas (HNSCCs) rank as the sixth most common cancer globally, contributing to significant morbidity and mortality, with over 900,000 new cases and 450,000 deaths reported in the past two years.^1,2 Despite notable advancements in diagnostic and therapeutic strategies, the survival rates for HNSCC remain stagnant, with a five-year survival rate of only 40% to 50%.^3,4 Base of tongue squamous cell carcinoma (BOTSCC), a subset of HNSCC, is particularly challenging due to its propensity for occult metastasis, high invasiveness, local recurrence, and poor prognosis. Human papillomavirus (HPV) has been identified as the primary etiological agent in BOTSCC, with HPV type 16 emerging as a major risk factor for oropharyngeal squamous cell carcinoma (OPSCC) according to findings by the International Agency for Research on Cancer.⁵

Among head and neck squamous cell carcinomas (HNSCCs), the tonsils and base of the tongue are the primary sites commonly affected by this malignancy. The identification of human papillomavirus (HPV) type 16 as a significant risk factor has provided critical insights into the etiological link between HPV infection and the development of oropharyngeal squamous cell carcinoma (OPSCC). This discovery underscores the importance of early detection and preventive measures in mitigating disease progression. HPV vaccination has become a cornerstone of cancer prevention strategies, particularly for individuals at high risk of developing base of tongue squamous cell carcinoma (BOTSCC). Immunization targeting nine HPV genotypes has demonstrated significant efficacy in reducing the incidence of OPSCC and its associated tumorigenesis.⁶

The global incidence of base of tongue squamous cell carcinoma (BOTSCC) has risen markedly in recent years.⁷ This increase is primarily driven by the growing prevalence of HPV-positive tongue squamous cell carcinoma (TOSCC) and BOTSCC metastasis.^8,9 Notably, approximately 70% of cases occur in men.⁵ Historically, HPV vaccination programs have focused predominantly on females, but this concerning trend underscores the critical need to expand vaccination efforts to include males. Broader immunization strategies targeting both genders are essential to mitigate the rising burden of HPV-associated malignancies and address this escalating public health challenge.^5,10

The tumor microenvironment (TME) of oropharyngeal squamous cell carcinoma (OPSCC) is a heterogeneous and complex milieu composed of various cell types, including cancer cells, diverse immune cells, stromal cells, and non-cellular extracellular matrix components.^11–14 The dynamic and intricate process of crosstalk and interaction between tumor cells and the surrounding microenvironment plays a pivotal role in the pathogenesis, prognosis, and metastatic potential of head and neck squamous cell carcinoma (HNSCC).^15–17 Moreover, the TME has a profound impact on immune cell activation and infiltration, which can significantly influence the efficacy of cancer immunotherapy.¹⁸

We have previously revealed a significant heterogeneity in gene expression within the TME in both HPV-associated and HPV-non-associated tonsillar squamous cell carcinoma (TSCC), highlighting the roles of cellular mediators and factors involved in various processes, including inflammation, immune crosstalk, transcription regulation, immune signaling pathways, signal transduction, oncogenesis, tumor suppression, angiogenesis, and apoptosis.⁸ Building on this foundation, the current study aims to explore the mRNA expression landscape of the HPV-positive base of tongue squamous cell carcinoma (BOTSCC). We conducted a differential expression analysis of genes specific to inflammation and immunity, oncogenes and performed cancer hallmark enrichment analysis on all upregulated genes.

Materials and Methods

FFPE Tissue Collection Patients Characteristics

We gathered eleven formalin-fixed paraffin-embedded (FFPE) tissue samples from BOTSCC for our multicenter investigation. In blind confirmation, two competent pathologists verified BOTSCC. The majority of the samples were invasive squamous cell carcinoma. The tumor size of all samples ranged from 0.2 to 0.5 cm in greatest dimension (microscopic measurements). Participants in this study were required to be at least eighteen years old, male, or female, have enough material for HPV detection and genotyping, and have access to BOTSCC FFPE tissue samples that were of acceptable quality. Individuals under the age of eighteen, immunocompromised, using chemotherapy, immunotherapy, or radiation, or with incomplete data were not included. From each block of paraffin, we cut a series of tissue slices that were 10 μm thick. The tissue was sectioned and processed using standard procedures to prevent cross-contamination and guarantee reliable results.

Nucleic Acid Extraction and cDNA Synthesis

The AllPrep DNA/RNA FFPE Kit from Qiagen (Hilden, Germany) was used to extract viral DNA total RNA in accordance with the instructions given. The purity and concentration of extracted DNA and RNA were assessed using Thermo Fisher Scientific’s NanoDrop2000 (Waltham, MA, USA). Following the manufacturer’s instructions, an RT2 first-strand synthesis kit (Qiagen, Germantown, MD, USA) was used for cDNA synthesis and genomic DNA removal. The synthesis procedure used 1.5 µL of RNA in total. Using the internal housekeeping gene GAPDH, the isolated nucleic acid’s purity and integrity were evaluated. When handling RNA, stringent precautions were taken to guarantee precise readings and avoid contamination, including operating in RNase-free conditions and within a biological safety cabinet.

Detection and Genotyping of HPV

SYBR Green/ROX Master Mix was used in a quantitative polymerase chain reaction (qPCR) approach for HPV detection and genotyping, following the manufacturer’s instructions. 8.5 µL of RNase-free water, 1 µL of forward primer, 1 µL of reverse primer, 2 µL of cDNA, and 12.5 µL of 2× RT2 SYBR Green Master Mix made up the qPCR reaction mixture. Supplementary Table 1 provides specific primers used in the qPCR tests. Utilizing an ABI 7500 Fast device (Applied Biosystems, Waltham, MA, USA), amplification was carried out. An initial denaturation at 95 °C for 10 minutes was followed by 40 cycles of denaturation at 95 °C for 15 seconds, annealing at 57 °C, and extension at 72 °C for 45 seconds for the cycling conditions.

Cellular RNA Processing and Cancer Inflammation, Immunity Crosstalk, Oncogenes and Tumor Suppressor Genes Transcriptomes

The following methods were used to examine mRNA expression. First, the RT2 First Strand Synthesis Kit (Qiagen, Germantown, MD, USA) was used to accomplish cDNA synthesis and genomic DNA removal in compliance with the manufacturer’s instructions. SABiosciences (Frederick, MD, USA) provided two different kits for RT2-PCR array analysis of the generated cDNA: the Oncogenes and Tumor Suppressor Genes Array (PAHS-502ZC-12) and the Cancer Inflammation and Immunity Cross-talk Array (PAHS-181Z). These arrays made it possible to evaluate the expression levels of 168 genes linked to angiogenesis, apoptosis, transcription factors, tumor suppressor functions, oncogenesis, inflammatory mediators, and tumor cell communication. The following methods were used to examine mRNA expression. First, the RT2 First Strand Synthesis Kit (Qiagen, Germantown, MD, USA) was used to accomplish cDNA synthesis and genomic DNA removal in compliance with the manufacturer’s instructions. SABiosciences (Frederick, MD, USA) provided two different kits for RT2-PCR array analysis of the generated cDNA: the Oncogenes and Tumor Suppressor Genes Array (PAHS-502ZC-12) and the Cancer Inflammation and Immunity Crosstalk Array (PAHS-181Z). These arrays made it possible to evaluate the expression levels of 168 genes linked to angiogenesis, apoptosis, transcription factors, tumor suppressor functions, oncogenesis, inflammatory mediators, and tumor cell communication. The master mix, which had a total volume of 2700 µL for one 96-well plate, was made up of 1248 µL of RNase-free water, 102 µL of cDNA, and 1350 µL of 2× RT2 SYBR Green Master Mix, as previusly described.¹⁹ After that, each PCR array well received 25 µL of the qPCR component mix. The ABI 7500 Fast device (Applied Biosystems, Waltham, MA, USA) was then used to amplify the 96-well array. To enable extension and fluorescence data collection, the cycling conditions comprised a 10-minute initial denaturation at 95 °C, 40 cycles of denaturation at 95 °C for 15 seconds, and annealing at 60 °C for 1 minute. Qiagen’s online RT2-PCR Profiler data analysis program (www.qiagen.com) was used to examine the relative gene expression data. The ΔΔCt method was utilized to quantify the levels of gene expression, and the results were standardized against five housekeeping genes.

Cancer Hallmark Enrichment Analysis

All elevated genes found in HPV-associated BTOSCC patients were subjected to enrichment analysis using the Cancer Hallmarks database (CancerHallmarks.com). This database provides a proven framework that describes basic biological concepts that are common to different types of cancer. An existing gene set and the cancer hallmark pathways listed in the database were integrated for the enrichment analysis. Genes and molecular changes associated with each hallmark were systematically annotated in order to perform a survival analysis based on the expression status.

Functional Analysis and Biological Network Construction

A biological network was constructed by doing a functional analysis of genes that were found to be considerably over-expressed. This study made use of the Retrieval of Interacting Genes/Proteins (STRING) database, version 11.0, which makes it easier to investigate known and anticipated interactions between proteins. With the use of these instruments, we sought to clarify the intricate connections between the over-expressed genes, offering information on their possible functions in biological functions and disease mechanisms. Through the STRING database, we were able to create a comprehensive interaction network by combining several kinds of data, such as computational predictions, database annotations, and experimental outcomes.

Survival Prognosis Analysis

The Gene Expression Profiling Interactive Analysis (GEPIA) database (http://gepia.cancer-pku.cn/index.html) was used to analyze the relationship between the survival prognosis and the expression levels of the ten most highly upregulated genes in HPV-associated BTOSCC patients. The GEPIA provided extensive data from the Cancer Genome Atlas (TCGA), allowing the analysis of gene expression patterns across various cancer types and normal tissues.

Analysis of Real-Time RT² Profiler PCR Array

The online RT2 Profiler Data study program was utilized for the study of real-time PCR array data and gene expression. Using the relative quantification 2-ΔΔCt technique, fold regulation values for every gene were ascertained. Five housekeeping gene mean values were used to standardize ΔCt results.

Data Analysis

Data analysis was conducted using GraphPad 10 software (GraphPad Software, San Diego, CA, USA). A t-test was used for group comparison. Pearson’s correlation test was used to evaluate the relationships between each variable. Statistical significance was defined as a p-value of less than 0.05.

Results

Basic BOTSCC Patients’ Characteristics

This study includes all FFPE from patients with head and neck cell cancer as well as healthy persons. Table 1 shows the demographic details of the HPV-associated BTOSCC patients, whose ages ranged from 46 to 77 years old, with a median of 57.

Table 1 Demographics and Prevalence of HPV Genotypes

Prevalence of HPV in the Base of Tongue Squamous Cell Carcinoma (BOTSCC)

In BOTSCC patients, HPV DNA was found in the PEPP tissue-extracted DNA. All of the BOTSCC patients had three hr-HPV genotypes identified. In terms of genotype detection, HPV-35 8 (88%), HPV-16 7 (77%), and HPV-45 1 (11%) were the most common. There was variation in the HPV infection across BOTSCC patients. Table 1 shows that among BOTSCC patients, double HPV infections (HPV16+HPV35), 6 (85%), and 1 (14%) were found, whereas single HPV-16 and HPV-35 infections were found in 1 patient (50%) and 1 patient (50%) respectively. The average age of the BOTSCC patients was around 43 years old, with a range of 21 to 77 years. HPV-16, HPV-35, and HPV-45 were the three kinds of HPV found. Each type of HPV has a different viral burden. The viral load of HPV-16 ranged from 21 to 36, with an average of 32.14 ± 4.79. The viral load of HPV-35 ranged from 18 to 35, with an average of 31.75 ± 5.56. A viral load of 36 was found in HPV-45. It is significant to remember that only three HPV genotypes had data on viral loads. Our findings demonstrated that BOTSCC patients were more evenly distributed, with half being Saudi and the other half not. There may be differences in the viral loads of different HPV varieties, which could affect how the disease develops and spreads. Our findings indicated that the viral load was higher for some HPV genotypes. Increasing viral loads are frequently linked to enhanced infectivity and a higher chance of developing related HNSCC malignancies. In contrast, certain HPV genotypes were either undetectable or had a decreased viral load. This is by no means an indication that the illness is less serious or less dangerous. Diseases may emerge even in people with low viral loads, especially those with weakened immune systems or other risk variables.

Transcriptome Analysis of Inflammation and Immunity in HPV-Associated BOTSCC

The Human Cancer Inflammation and Immunity Crosstalk RT2 Profiler PCR Array was used to analyze the transcriptome profile of the interaction between inflammation and immunity. According to the analysis, 42 of the 84 immune genes in this qPCR array that are involved in the connection between tumor cells and the biological mediators of inflammation and immunity crosstalk were differently dysregulated mRNA in individuals with HPV-positive BOTSCC (Figure 1). Patients with HPV-positive BOTSCC had elevated expression of 20 of their genes (Figure 1). In individuals with HPV-associated BOTSCC, the mRNA of four genes (IL6, IRF1, TLR2, and TNF) was substantially elevated (Figure 1). Patients with HPV-associated BOTSCC showed downregulation of 22 genes. We discovered that the levels of VEGFA and CCL18 mRNA were low. Four genes (BCL2, CXCL9, IGF1, and MYC) were found to be elevated in all HPV-associated BOTSCC patients, whereas six genes (VEGFA, CSF3, CCR2, CSF2, CXCL12, and MIF) were frequently found to be downregulated. In patients with HPV-associated BOTSCC, there was significant variation in TME gene expression. Of the patients, one (p1) tested positive for two genes (4%), three (p2, p3, and p6) tested positive for seven genes (16%), one (p4) tested positive for thirteen genes (30%), one (p5) tested positive for five genes (11%), one (p7) tested positive for four genes (9%), and one (p8) tested positive for six genes (14%), as shown in Figure 1. According to inter-age compression, patients between the ages of 42 and 55 showed more gene expressions than those between the ages of 55 and 97 (Figure 2). A number of variables, including age, can cause physiological alterations in the immune system that may affect both innate and adaptive immune responses. These results imply that the gene expression profile may be influenced by age. The gene expression of cancer inflammation and immunity interaction in the HPV-associated BOTSCC microenvironment showed distinct heterogeneity.

Figure 1 Heatmap of the 42 dysregulated genes’ expression profiles in patients with BOSTCC linked to HPV. Dysregulated tumor inflammation and immune crosstalk gene expression are displayed in the ring heatmap cluster gram. HPV-associated BOSTCC patients are represented by rings P1 through P8. While the color gradient depicts expression from low (blue) to high (red), the heatmap color scale shows the range of expression values. Co-expression is indicated by genes with comparable colors across samples. The hierarchical clustering of genes according to their patterns of expression is displayed by the dendrogram, or tree structure.

Figure 2 Mean fold-change expression of genes involved in tumor cell-to-tumor cell contact and cellular mediators of inflammation and immune crosstalk in two age groups of patients with HPV-associated BOTSCC. The data emphasize the effect of age on gene expression patterns in the context of HPV-associated BOTSCC by showing the difference in expression levels between patients under the age of 55 and those over 55.

Analysis of mRNA Expression of Oncogenes and Tumor Suppressor Genes of TME in HPV-Associated BOTSCC

A key factor in the development, spread, and clinical results of malignancies is the dysregulation of the molecular landscape. In the microenvironment of HPV-associated BOTSCC, we used the RT² Profiler™ PCR Array for human oncogenes and tumor suppressor genes to examine the expression levels of 84 oncogenes and tumor suppressor genes. According to our study, patients with HPV-positive BOTSCC had dysregulated 64 genes (Figure 3). Of these, 33 genes were found to be increased, with nine genes exhibiting significantly higher levels of mRNA expression: APC, CDKN2B, E2F1, JUN, KRAS, MYB, RAF1, RUNX1, and XRCC1. In contrast, three genes—CDKN1A, MCL1, and MDM2—that showed noticeably low expression were among the 31 genes that were downregulated. Additionally, we found that six generally downregulated genes (BAX, CDKN1A, MCL1, MDM2, RAR, and STK11) and seven commonly upregulated genes (BCR, CASP8, NF1, PIK3C2A, S100A4, VHL, and WWOX) were shared by the majority of HPV-positive BOTSCC patients (Figure 3). Two patients (p1 and p15) received positive results for 15 genes (23%), one patient (p2) received positive results for 18 genes (28%), one patient (p4) received positive results for 7 genes (10%), one patient (p5) received positive results for 13 genes (23%), one patient (p6) received positive results for 9 genes (14%), one patient (p7) received positive results for 4 genes (6%), and one patient (p8) received positive results for 5 genes (7%). Furthermore, compared to patients aged 55–97 years, those aged 42–55 years showed more gene expressions, according to interage comparisons (Figure 4). A number of variables, including age, can cause physiological alterations in the immune system that could affect both innate and adaptive immune responses.

Figure 3 The 64 dysregulated gene expression profiles in HPV-associated BOSTCC patients are shown in a heatmap. The cluster-gram of the ring heatmap of the expression of tumor suppressors and dysregulated oncogenes. Patients with BOSTCC linked to HPV are represented by rings (P1 to P8). The range of expression values is shown by the color scale of the heatmap. Expression is represented by the color gradient, which ranges from low (blue) to high (red). Co-expression is indicated by genes with comparable colors across samples. The hierarchical clustering of genes according to their expression patterns was displayed by the dendrogram, or tree structure.

Figure 4 Genes linked to tumor suppression and oncogenesis expressed at mean fold-change in two age groups of patients with HPV-associated BOTSCC. Age has an impact on gene expression profiles in the setting of HPV-associated BOTSCC, as seen by the data showing the varied expression levels between patients under and over the age of 55.

Clustering Analysis of Inflammation, Immunity, and Oncogene/Tumor Suppressor Gene Expression Profiles

Six different clusters of gene expression profiles were identified in the examination of the crosstalk between inflammation and immunity (Figure 5). The interleukins and chemokines (IL6, IL15, IL17A, CCL4, CXCL9, and CXCL10) that are part of Cluster 2 showed the greatest levels of overall expression. ACKR3, CCR1, CCR9, CXCR2, and TLR2 are part of Cluster 3 (chemokine receptors and Toll-like receptor signaling), which came after Cluster 4 (growth factors and receptors), which includes IGF1 and TGFB1. Significant expression of BCL2 and MYC was also observed in Cluster 6 (anti-apoptotic factors). TNF, GZMB, and MICB are among the enzymatic modulators of inflammation, immunostimulatory factors, and antigen presentation that make up Cluster 1, which, in contrast, showed comparatively low gene expression. The expression profile of FOXP3 and IRF1 in Cluster 5 (transcription factors) was in the middle. Likewise, tumor suppressor genes and oncogenes were divided into five groups (Figure 5). CASP8, CDK4, ETS1, JAK2, JUND, KITLG, MOS, MYB, PML, PRKCA, RAF1, RUNX1, and STAT3 are all members of Cluster 1 (oncogenes), which had the greatest expression levels. Following this, significant expression levels were shown by Cluster 3 (both oncogenic and tumor suppressor genes), which includes BCR, ESR1, JUN, KRAS, MYC, PIK3C2A, RB1, and RET, and Cluster 5 (epithelial-to-mesenchymal transition), which includes S100A4 and VHL. The expression profile of Cluster 2 (tumor suppressor genes), which comprises NF1, XRCC1, and IGF2R, was moderate. On the other hand, when compared to the other clusters, Cluster 4 (cell cycle), which included CDKN2A, ATM, BRCA1, BRCA2, CDKN2B, and E2F1, showed the lowest expression of tumor suppressor genes and oncogenes.

Figure 5 Analysis of Gene Expression Profile Clustering. Tumor suppressor genes, oncogenes, inflammation, and immunity-related gene expression profiles are grouped into different clusters by the study. The relationship between inflammation and immunity Groupings: Cluster 1: Antigen presentation, inflammatory enzymatic modulators, and immunostimulatory factors. Interleukins and chemokines make up Cluster 2. Cluster 3: Signaling from Toll-like receptors and chemokine receptors. Growth factors and receptors make up Cluster 4. Transcription factors are in cluster five. Cluster 6: Factors that prevent death. Oncogenes and Suppressors of Cancer Groupings: Group 1: Oncogenes. Cluster 2: Genes that inhibit tumor growth. Cluster 3 includes both tumor suppressor and carcinogenic genes. Cell cycle-related genes make up Cluster 4. Transition from epithelium to mesenchymal tissue is Cluster 5.

Dysregulated Gene Expression Was Negatively Correlated With Age

One of the variables influencing a patient’s immune level is their age. Thus, we looked into a possible relationship between age and over-expressed genes in BOTSCC patients. The following genes have a positive association with age, according to correlation analysis, however, it was not statistically significant: MYC (r=0.669; p = 0.331) and CAPS8 (r=0.982; p = 0.122). Four genes—BCL2 (r=−0.70; p=0.125), CXCL9 (r=−0.87; p=0.317), PIK3C2A (r=−0.970; p= 0.157), VHL (r=−0.458; p=0.543), and WWOX (r=−0.432; p=0.285)—showed a negative connection with age, but this link was again not deemed significant. While more research is required to ascertain the significance of these associations, these results imply that age may influence the expression levels of these genes.

Functional Analysis and Biological Network Construction

The online STRING site (cn.string-db.org) was used to build the PPI network. Humans were chosen when the most elevated genes of HPV-associated BOTSCC (IL6, IRF1, TLR2, TNF, APC, CDKN2B, E2F1, JUN, KRAS, MYB, RAF1, RUNX1, and XRCC1) were added to the gene list. 0.400 was selected as the median degree of confidence. By using the k-means clustering option, the PPI network was divided into three clusters (Figure 6). This PPI included direct (physical) or indirect (functional) interactions, as well as several gene linkages, including gene neighborhood, gene fusions, co-occurrence, co-expression, text mining, protein homology, and biochemical/genetic data (experiments).

Figure 6 Protein-protein interactions with STRING PPI. Different colored lines indicate string PPI interactions. The protein-protein association is represented by the edges. Interactions with known co-occurrence evidence and experimental evidence are indicated by blue and purple borders. Text mining and co-expression are indicated by black and yellow margins. Various hues correspond to different kinds of interactions.

Cancer Hallmark Enrichment Analysis

The integrated gene set of the cancer Hallmarks database, which includes 6763 genes associated with different hallmark capabilities, was used to map all elevated genes to the database. By methodically annotating the genes and molecular changes connected to every characteristic (Figure 7 and Supplementary Table 2), we found 51 genes (CXCL10, MIF, TLR2, ESR1, CDKN2B, AKT1, MYC, IGF2R, ACKR3, MCL1, CDK4, NF1, CCR2, CCR1, ATM, APC, BAX, IL17A, ETS1, CXCL9, KRAS, CCR10, TNF, CSF2, BCR, JAK2, CXCR2, IL15, BRCA1, CASP8, CSF3, STAT3, IGF1, IL13, RAF1, RET, JUND, JUN, IL6, MYB, CXCL12, TSC1, CCL5, CCR9, IL12B, KITLG, PRKCA, BCL2, CCL4, TGFB1 and MOS) which were linked to sustaining proliferative signaling. 30 genes (CXCL10, TLR2, ESR1, AKT1, CCR2, CCR1, ATM, IL17A, CXCL9, KRAS, CCR10, TNF, IRF1, CSF2, JAK2, CXCR2, IL15, CASP8, STAT3, RAF1, RET, JUN, IL6, CXCL12, CCL5, CCR9, IL12B, PRKCA, BCL2 and CCL4) were involved in tumor-promoting inflammation. Genome instability was correlated with 16 genes (TNF, BRCA2, IGF1, XRCC1, CDK4, APC, ATM, CDKN2A, ESR1, BAX, RB1, AKT1, MYC, BRCA1, GZMB, CASP8). 37 genes (CXCL10, TLR2, ESR1, AKT1, MYC, ACKR3, NF1, CCR2, CCR1, IL17A, ETS1, CXCL9, KRAS, CCR10, TNF, CSF2, CXCR2, IL15, BRCA1, CASP8, CSF3, STAT3, IGF1, IL13, RAF1, VHL, IL6, CDKN2A, CXCL12, CCL5, CCR9, IL12B, KITLG, PRKCA, RUNX1, CCL4 and TGFB1) were implicated in sustained angiogenesis. 32 genes (CXCL10, TLR2, ESR1, RB1, AKT1, MYC, CDK4, CCR2, CCR1, IL17A, CXCL9, KRAS, CCR10, TNF, IRF1, CSF2, JAK2, CXCR2, GZMB, CASP8, MICB, STAT3, RAF1, RET, JUN, IL6, CXCL12, CCL5, CCR9, IL12B, PRKCA and CCL4) helped to evade immune destruction. Replicative immortality involved 30 genes (MIF, ESR1, CDKN2B, RB1, MYC, AKT1, E2F1, BRCA2, CDK4, APC, ATM, ETS1, CXCL9, KRAS, TNF, JAK2, PML, CSF3, STAT3, IGF1, RAF1, JUND, JUN, IL6, CDKN2A, CXCL12, RUNX1, PRKCA, BCL2, TGFB1). 54 genes (CXCL10, MIF, TLR2, ESR1, CDKN2B, RB1, MYC, AKT1, IGF2R, E2F1, ACKR3, BRCA2, MCL1, CDK4, NF1, CCR2 CCR1, ATM, APC, BAX, IL17A, CXCL9, KRAS, CCR10, TNF, CSF2, JAK2, CXCR2, IL15, PML, BRCA1, GZMB, CASP8, CSF3, STAT3, IGF1, IL13, RAF1, VHL, RET, JUN, IL6, MYB, CDKN2A, CXCL12, TSC1, CCL5, CCR9, IL12B, KITLG, PRKCA, BCL2, CCL4 and TGFB1) were associated with resisting cell death. Evading growth suppressors were linked to 52 genes (CXCL10, TLR2, RUNX3, ESR1, CDKN2B, RB1, MYC, AKT1, IGF2R, E2F1, ACKR3, MCL1, APC, CDK4, NF1, CCR2, CCR1, ATM, BAX, ETS1, CXCL9, KRAS, CCR10, TNF, FOXP3, IRF1, BCR, JAK2, CXCR2, IL15, PML, BRCA1, CSF3, STAT3, IGF1, IL13, RAF1, VHL, RET, JUN, IL6, CDKN2A, CXCL12, TSC1, CCL5, CCR9, KITLG, PRKCA, BCL2, CCL4, TGFB1 and MOS). 19 genes (CDKN2B, RB1, MYC, AKT1, E2F1, ATM, KRAS, JAK2, BRCA1, STAT3, VHL, IGF1, RAF1, RET, JUND, IL6, TSC1, BCL2 and TGFB1) were implicated in the reprogramming of energy metabolism. Finally, 43 genes (CXCL10, TLR2, ESR1, CDKN2B, RB1, MYC, AKT1, E2F1, ACKR3, CDK4, NF1, ATM, CCR1, APC, BAX, ETS1, S100A4, CXCL9, KRAS, TNF, BCR, JAK2, CXCR2, BRCA1, GZMB, CASP8, STAT3, VHL, IGF1, RAF1, RET, JUND, JUN, PIK3C2A, CDKN2A, CXCL12, IL12B, KITLG, PRKCA, BCL2, CCL4, TGFB1 and MOS) had a role in metastasis and tissue invasion. Multiple cancer signature pathways were mapped to the individual genes in HPV-associated BOTSCC, illustrating the intricate interdependence of the molecular mechanisms underlying tumor growth.

Figure 7 Using CancerHallmarks.com, the elevated gene set in HPV-associated BOTSCC is compared to established hallmark gene sets in a cancer hallmark enrichment plot. The red dotted line indicates an adjusted p-value threshold of less than 0.05, while colored slices indicate significantly enhanced characteristics. The distribution of elevated genes across several cancer hallmark pathways in HPV-associated BOTSCC is depicted by the dot plots. Markers with more mapped genes are indicated by black dots, and those with less mapped genes are indicated by grey dots.

The Association of Upregulated Genes in HPV-Associated BTOSCC Patients With a Poor Prognosis

As determined by the GEPIA2 tool, the Ten most strongly elevated genes (CCL4, BCR, FOXP3, IGF1, IL17A, S100A4, MYB, MYC, NF1, and PIKC2A) were evaluated for prognostic significance. Several of these genes were shown to be strongly linked to worse patient outcomes and lower overall survival, according to the data (Figure 8). These results imply that in HPV-related BTOSCC patients, these genes are associated with a poor prognosis.

Figure 8 The survival prognosis of cervical cancer is correlated with the expression levels of the 10 genes that are most substantially elevated in CSCC. Overall survival rate: The inadequate survival rate is shown by the red line; the opposite is shown by the blue line.

Discussion

Recent years have seen a significant increase in studies into the TME because of its well-established involvement during cancer and its importance in assessing therapy response. Gaining a thorough grasp of the distinct elements of the TME, such as immunogens, oncogenes, and tumor suppressor genes, is still very difficult. To overcome this difficulty, this study profiled the expression of 168 genes in the TME of HPV-associated BOTSCC that are involved in communication between tumor cells and the cellular mediators of inflammation and immunity, oncogenesis promotion, tumor suppressor genes, transcription factors, angiogenesis, and apoptosis.

The tumor microenvironment (TME) of oropharyngeal squamous cell carcinoma (OPSCC) is a complex and diversified environment made up of stromal cells, cancer cells, immune cells, and non-cellular extracellular matrix components.^11–13 One of the most important factors in tumor growth, progression, and metastasis in OPSCC is the dynamic and intricate process of crosstalk and interaction between tumor cells and the surrounding microenvironment.¹⁷ Modulating immune cell infiltration and/or activation is another way that the TME can impact the results of cancer immunotherapy.¹⁸ By altering immune cell activation and invasion, the TME can also affect the results of cancer immunotherapy.⁹ The TME is heterogeneous, as demonstrated by its characterization in OPSCC; a variety of causes can cause variations in the patterns of gene expression within the TME.¹² Important gene signatures that are used as prognostic and immunotherapeutic targets have been found in earlier research,¹⁴ and OPSCC patients who have dysregulated genes have a lower overall survival rate.¹²

The most common kind of cancer in the oral cavity is called BOTSCC, and it is classified as an oropharyngeal cancer (OPC). In addition to having aggressive clinical behavior,²⁰ BOTSCC is becoming more common.²¹ BOTSCC exhibits a high proportion of hidden metastases, a high rate of local recurrence, and poor prognostic outcomes.²¹ The development of cancer may be influenced by genomic aberrations, which are modifications in the DNA sequence. The malignant phenotype may be caused by a variety of genetic defects in addition to well-known risk factors like alcohol, smoking, and virus-induced carcinogenesis.

In OPSCC, for example, oncogenes—genes that might cause cancer—have been shown to be activated. To improve tumor clinical outcome prediction and treatment planning, it is crucial to characterize the BOTSCC microenvironment and find informative prognostic biomarkers. Tumor start and progression are caused by genetic stability. We sought to investigate 84 genes associated with oncogenesis promotion, tumor suppressor genes, transcription factors, angiogenesis, and apoptosis to highlight the significance of the BOTSCC microenvironment.

The TME of HPV-associated BOTSCC was identified as unique, exhibiting dysregulation of several genes linked to angiogenesis, apoptosis, tumor suppression, and oncogenesis. Ten genes that were shown to have dysregulated expression profiles were BAX, CASP8, MCL1, S100A4, STK11, BCR, NF1, PIK3C2A, CDKN1A, MCL1, MDM2 VHL, RARA, and WWOX. In the development of head and neck squamous cell carcinoma, genes linked to angiogenesis, transcription factors, tumor suppressor genes, oncogenesis promotion, and the apoptosis pathway are crucial.²² There is evidence that the casp8 gene is important in tumor biology, as its overexpression is linked to the onset and spread of several cancer types.^23–25

Cell proliferation, apoptosis, and metastasis are critical processes in the development of cancer, and their dysregulation has been connected to the increased expression of casp8.²⁴ In addition to its role in the primary tumor cells, CASP8 also demonstrates its importance in the TME.²³ Its ability to maintain a healthy and functional immune system is one of its essential roles in controlling the immunological response. Additionally, Caspase-8 is essential for angiogenesis and is resistant to treatment in glioblastoma.²⁶ The study discovered that HPV-associated BOTSCC had higher levels of CASP8 expression. Potential therapeutic targets for the treatment of cancer may become clearer if the molecular mechanisms underlying CASP8 overexpression in malignancies are understood.

Diffuse large B-cell lymphoma (DLBCL) and chronic lymphocytic leukemia (CLL) are two examples of B-cell malignancies that are associated with the B-cell receptor gene (BCR).^27,28 Uncontrolled cell growth and proliferation in cancer can result from the expression of the BCR gene, according to mounting evidence.²⁹ A predictive biomarker associated with immune infiltrates in kidney renal clear cell carcinoma, phosphoinositide 3-kinases (PIK3C2A) expression shows a strong association with the infiltrating levels of primary immune cells.³⁰ According to an analysis of the study’s findings, HPV-associated BOTSCC had increased expression levels of PIK3C2A and, more significantly, MCL-1. There is a correlation between cancer and MCL-1 overexpression.

The development of drug resistance and a poor prognosis have been linked to elevated MCL-1 levels, which further complicates the treatment of cancer patients.^31,32 MCL-1 is a desirable target for cancer treatment due to its crucial function in controlling the mitochondrial apoptotic pathway.³² MCL-1 is a desirable possible target for cancer treatment because of its crucial function in controlling the mitochondrial apoptotic process. Analysis revealed that in around 75% of cases of HPV-associated BOTSCC, there is a notable elevation in the gene expression levels of calcium-binding protein A4 (S100A4). A significant contributor to the development and metastasis of several cancer types is S100A4, a member of the S100 protein family.^33,34 Interestingly, tumor pathophysiology is not the only area where S100A4 is upregulated. In a variety of disease disorders unrelated to tumors, it has also been noted.³³ This implies that S100A4 may have more general effects than only cancer biology.

In addition to its role in malignancies, S100A4 plays a complex role in controlling immunological responses and inflammation.³⁵ S100A4 may be a suitable therapeutic target for the treatment of HPV-associated BOTSCC, according to these findings. The fact that HPV-associated BOTSCC has upregulated S100A4 highlights the protein’s significance in cancer biology, its potential as a target for cancer therapy, and its ability to regulate immunological and inflammatory responses.

Moreover, retinoic acid (RARA) mRNA levels were shown to be downregulated. RARA may stop human tumor cell lines from growing, according to earlier research.^36,37 In soft agar tests, RARA treatment was found to decrease DNA synthesis, cause morphological alterations, lengthen the time it takes for cells to double, and decrease colony formation and saturation density. Importantly, tumor cells are arrested during the G1 stage of the cell cycle rather than dying as a result of retinoid reduction of growth.³⁷ Bax is upregulated in various cancer types, according to numerous studies.³⁸ All patients with HPV-associated BOTSCC, with the exception of one, had downregulated BAX. The main negative regulatory element for p53 is MDM2.

Several solid tumors including 40–60% of human sarcomas have over-expressed and/or up-regulated MDM2 proteins.³⁹ MDM2 overexpression was discovered in the HPV-associated BOTSCC. OPSCC linked to HPV has a unique TME. Immunological cell infiltration and the overexpression of many immunogens and immunological signaling pathways are characteristics of HPV-associated OPSCC.^40,41 This study demonstrated that the HPV interaction may be crucial in modifying the TSCC microenvironment and resulting in suboptimal clinical outcomes. One new approach to treating BOTSCC may be to use antagonized drugs to target those genes. To classify individuals who might benefit from immunotherapy, more differentiation of the immune-related genes of BOTSCC is required.

The precise function of those genes in HPV-associated BOTSCC may be fully and deeply understood molecularly, which may help find a possible target for efficient cancer therapy strategies. In patients with BOTSCC, the expression patterns of 84 different genes were analyzed. It is commonly known that these genes are essential for many biological functions, such as signal transduction pathways, transcriptional control, immunological crosstalk, and inflammation. Forty-two genes that are essential for interactions between tumor cells and cellular mediators of inflammation and immunity were found to be dysregulated.

The intricate interactions between the immune system and the tumor microenvironment in these patients were clarified by this discovery. Twenty of these forty-two genes were found to be elevated. According to this, these genes might be implicated in both immune response evasion and tumor growth promotion. On the other hand, a downregulation of 22 genes indicated a decrease in their expression. These downregulated genes may contribute to immune response suppression or tumor growth inhibition. CXCL9, IGF1, and anti-apoptotic (BCL2 and MYC) genes were discovered to be elevated. The complex relationship between BCL2 gene expression and cancer can change based on the setting and kind of disease.

According to earlier research, patients with BCL2 expression who have breast cancer typically had higher overall and disease-free survival rates.⁴² BCL2 expression’s predictive value in colorectal cancer is unknown. According to certain research, there is no predictive value for BCL2 expression in colorectal cancer.⁴³ Both favorable and unfavorable prognoses for many malignancies have been linked to BCL2 expression.⁴⁴ The results of this investigation are consistent with a number of studies that have systematically documented elevated levels of CXCL9 expression in HNC squamous cell carcinoma, which have been strongly linked to adverse clinical outcomes.¹³

A possible function of CXCL9 as a biomarker for determining the prognosis of HNC patients is thus suggested. These findings also demonstrated the increase in IGF1 mRNA levels. By inhibiting cell death and promoting cell proliferation, IGF-1 fuels cancer. Research has shown a connection between elevated IGF-1 levels and malignancies such as prostate, colorectal, and breast cancer.^45,46 Demonstrating even more how it might be connected to the pathophysiology of HNCC. VEGFA, CSF2, and CSF3 are growth factor-related genes that were down-expressed in most HPV-associated BOTSCC patients.

The present study’s results contradict those of other investigations that demonstrated increased expression of angiogenic genes, such as vascular endothelial growth factor (VEGF), which has been linked to poor prognosis and resistance to treatment in head and neck cancer.^47,48

The chemokine CXCL12 is essential for several physiological and pathological processes, such as tumor metastasis, immune cell migration, and stem cell mobilization.⁴⁹ Cancer cell migration, invasion, and metastasis are linked to CXCL12 downregulation, which has been seen in a variety of cancer types, including breast, lung, and colon cancer.^49,50

Those with HPV-associated BOTSCC had down-regulated levels of CXCL12 and CSF3 mRNA, according to the current study. Targeting CSF3 may be a viable cancer therapeutic approach, as some research has revealed biochemical pathways and immunological responses in the colo-rectal TME.⁵¹ A noteworthy finding of this study was that patients with HPV-associated BOTSCC also had significantly lower levels of CCR2 and Macrophage Migration Inhibitory Factor (MIF). It has been shown that CCR2 downregulation contributes significantly to the onset and spread of cancer. Research has demonstrated that a decrease in tumor development and metastasis is linked to a decrease in CCR2 expression.^52,53

One promising treatment approach for cancer is to target the down-regulation of CCR2.⁵⁴ Increased survival rates, better patient outcomes, and the creation of innovative therapeutic strategies could result from this. Several cancers, including breast, lung, and colon cancers, have been linked to MIF downregulation.⁵⁵ MIF is a cytokine involved in inflammation and immunological control, and its downregulation may promote the migration and proliferation of cancer cells.⁵⁶ Research has indicated that the downregulation of MIF may represent a new mechanism of resistance to anti-angiogenic treatment.⁵⁷ To develop targeted medicines and improve treatment results for individuals with HPV-associated BOTSCC, it is imperative to comprehend the differential expression of these genes. The collected information highlights the TME’s varied and intricate role in OPSCC. Effective therapeutic approaches for cancer require an awareness of the TME’s heterogeneity and its critical role in the disease’s course and response to treatment. By focusing on the channels of communication between tumor cells and the immune system, immunotherapies may work better or new therapeutic approaches may be developed.

The inter-age compression study revealed that patients between the ages of 42 and 55 had a considerably greater number of gene expressions than patients between the ages of 55 and 97. This finding draws attention to the possible impact of aging on the gene expression profile. Age is known to be one of the main variables that might cause physiological alterations in the immune system and numerous other tissues.^58,59 Genes that are essential for preserving cellular function and response to different stimuli, such as tumor suppressor genes, oncogenes, and immune-related genes, may be impacted by these alterations.

The results of this analysis showed that patients of various ages exhibited a unique pattern of gene expressions. These results provide insight into how aging may affect the immune system and its gene expression profile. The mechanisms behind age-related variations in gene expression and their effects on immune function require more investigation in this field. The results of the study, however, offer important new information about the molecular processes that underlie the interaction of immune cells and tumor cells.

The most common hr-HPV genotypes found in HNSCC worldwide are HPV-16, which is followed by HPV-31, HPV-18, HPV-56, HPV-52, HPV-33, and HPV-35, among other hr-HPV genotypes.⁶⁰ The most frequently discovered genotypes in this investigation were HPV-35, HPV-16, and HPV-45. In individuals with BOTSCC, HPV infection levels ranged from one to two infections. More HPV infections are associated with larger viral loads, which raises the incidence of cervical lesions in comparison to single HPV infections, according to several studies.^61–63 However, according to some studies, compared to several infections, a single HPV infection considerably raises the risk of developing HNSCC. The increasing incidence of HPV+ tonsillar (TSCC) and BOTSCC has prompted study into prognostic markers and driver genes that could guide targeted therapy. These cancers generally exhibit a more favorable prognosis than their HPV-negative counterparts. Identifying predictive markers and understanding the genetic landscape of HPV+ tumors is crucial for refining treatment strategies. Recent studies focus on specific driver genes that could inform tailored therapies, aiming to improve patient outcomes while minimizing the burden of aggressive treatments.⁶⁴

Achieving successful cancer treatment depends on the TME of HNSCC, which is diverse and complex. Comprehending the interplay between cancer inflammation and immunology, oncogenes, and tumor suppressor genes is crucial for creating efficacious cancer treatments. The development of a three-dimensional collagen-based scaffold model offers a novel approach to studying the microenvironment and drug-resistance mechanisms in oropharyngeal squamous cell carcinomas (OSCC).⁶⁵ In this study, the TME of HPV-associated BOTSCC and HPV-associated TSCC showed variability in gene expression. The mRNA levels of multiple genes, including STK11, RARA, and MDM2, were found to be upregulated in HPV-associated TSCC. However, it was discovered that HPV-associated BOTSCC had downregulated versions of these same genes. For these two kinds of oropharyngeal cancer subsites, a reverse pattern of gene expression was observed. While the mRNA levels of genes such as WWOX, BCR, BCL2, and CXCL9 were elevated in HPV-associated BOTSCC, they were downregulated in HPV-associated TSCC.⁸ These results point to the different molecular profiles found in the TME of HPV-associated BOTSCC and HPV-associated TSCC, which may indicate variations in the underlying biological processes and potential treatment targets. All things considered, the work emphasizes how varied and intricate the gene expression profiles are in HPV-associated BOTSCC, illuminating possible areas for further investigation and treatment.

In conclusion, our study highlights the complex interplay between HPV genotypes and gene expression profiles in patients with HPV-positive BOTSCC. The presence of multiple high-risk HPV types, particularly HPV-16 and HPV-35, indicates a significant viral burden, which may influence disease progression and clinical outcomes. The observed dysregulation of immune and oncogenic gene expression, particularly among younger patients, suggests that age-related physiological changes could modulate immune responses and tumor behavior. Despite variations in HPV genotypes and gene expression among patients, the consistent elevation of certain inflammatory and oncogenic markers underscores the potential of HPV in shaping the tumor microenvironment. These findings emphasize the need for further investigation into the molecular mechanisms underlying HPV-related malignancies, as well as the implications for targeted therapeutic strategies in BOTSCC patients. Understanding these dynamics could ultimately enhance patient management and treatment outcomes in this vulnerable population.

This study has limitations. First off, the sample size was small, which could limit the findings’ generalizability and have an effect on their statistical significance and generalizability. Furthermore, the study did not look into the functional implications of the identified over-expressed genes in HPV-associated BOTSCC, which could have shed more light on the underlying mechanisms. To obtain a more thorough grasp of this field, it is advised that future research build on these constraints and investigate the functional implications of the found over-expressed genes.