Comprehensive analysis of the significance of ferroptosis-related genes in the prognosis and immunotherapy of oral squamous cell carcinoma

doi:10.21203/rs.3.rs-115737/v1

Download PDF

Research Article

Comprehensive analysis of the significance of ferroptosis-related genes in the prognosis and immunotherapy of oral squamous cell carcinoma

https://doi.org/10.21203/rs.3.rs-115737/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

Background: Oral squamous cell carcinoma (OSCC) is a life-threatening disease that emerged as a major international health concern, associated with poor prognosis and the absence of specific biomarkers. Studies have shown that the ferroptosis-related genes (FRGs) can be used as tumor prognostic markers. However, FRGs’ prognostic value in OSCC needs further exploration. Our aim was to construct a novel FRG signature for overall survival (OS) prediction in OSCC patients and explore its role in immunotherapy.

Methods: In our study, gene expression profile and clinical data of OSCC patients were collected from a public domain. FRGs were available from the FerrDb database. We performed univariate and multivariate Cox regression analyses to construct a multigene signature. The Kaplan-Meier (K-M) and receiver operating characteristic (ROC) methods were utilized to test the effectiveness of the FRG signature. A differential gene expression analysis was performed by the limma R package, followed by functional enrichment analyses. CIBERSORT was applied to analyze the tumor microenvironment (TME). Finally, the expression of human leukocyte antigen (HLA) and immune checkpoint molecules were analyzed to confirm the sensitivity of immunotherapy.

Results: A total of 103 FRGs, expressed in OSCC (FRGs-OSCC), were identified from the two datasets by the Venn analysis. The Cox regression analysis identified 5 FRGs-OSCC that were associated with overall survival (all P < 0.01). The FRGs-OSCC risk model was established to classify patients into high risk and low risk groups. Compared with the low risk group, the survival time of the high-risk group was significantly reduced (P < 0.001). According to the multivariate Cox regression analyses, the risk score acted as an independent predictor for OS (HR > 1, P < 0.001). The accuracy of the FRGs-OSCC risk predictive model was confirmed by ROC curve analysis. The results of the Kyoto Encyclopedia of Genes and Genomes (KEGG) showed significant enrichment of immune-related pathways, and a difference in tumor microenvironment between the two groups. The low risk group had the characteristics of higher expression of HLA and immune checkpoints (IDO1, LAG3, PDCD1 and TIGHT), a lower tumor purity and a higher infiltration of immune cells, indicating a more sensitive response to immunotherapy.

Conclusions: The novel FRGs-OSCC risk score system can be used to predict OSCC prognosis. Ferroptosis targeting may be a therapeutic option for OSCC.

Epigenetics & Genomics

Oral squamous cell carcinoma

prognosis

ferroptosis-related genes

tumor immune microenvironment

immunotherapy

OSCC is the most common type of oral cancer, with more than 300,000 newly diagnosed cancer cases worldwide in 2018, according to the GLOBOCAN database (1). OSCC predilection sites include the tongue, the alveolar, the mouth floor, the lips, and the buccal mucosa. OSCC is characterized by high recurrence, metastatic, and mortality rates, especially, in patients with a late diagnosis (2). In recent years, multidisciplinary collaborative diagnosis and treatment has been proposed for OSCC treatment, and which included chemotherapy, biological therapy, and radiotherapy. Although multimodality therapies can improve the prognosis, the 5-year OS rate of OSCC patients remains stable, at approximatively 56%, and the post-treatment local recurrence and distant metastasis rate are 25%-50% (3)(4). Therefore, there is a need to find effective prognostic biomarkers that could guide these management decisions.

Ferroptosis is a new form of programmed cell death, characterized by iron overload and lipid peroxidation that cause lipid reactive oxygen species (ROS) accumulation (5)(6). Numerous studies demonstrated that ferroptosis is involved in cancer initiation, progression, and suppression (7). For instance, the tumor suppressor gene p53 may modulate the susceptibility of cancer cells to ferroptosis in a cell type specific manner (8)(9). Artesunate, a clinically approved drug, can selectively kill OSCC cells by inducing ferroptosis (10). Recently, Kotaro et al. (11) found that non-thermal plasma exposure kills OSCC cells through a specific mechanism that depends on ample catalytic Fe(II). There are two mixed forms of programmed cell death that are caused by this treatment method, including apoptosis and ferroptosis, which suggest that ferroptosis might be closely related to OSCC occurrence. Moreover, ferroptosis is also associated with the efficacy of immunotherapy (12).

It is well known that CD8 + T cells, generally induce tumor cell death, through the pore-forming protein-granzyme and the Fas/FasL pathways (13)(14). However, a new study (15) showed that CD8 + T cells that are activated by immunotherapy, augment ferroptosis in tumor cells, which contributes, therefore, to the antitumor efficacy of immunotherapy. This discovery provided important evidence for the correlation between ferroptosis and antitumor immunity. Besides, Lang et al. (16) also found that interferon-γ that is produced by CD8 + T cells, has a synergistic effect with radiotherapy-activated ataxia-telangiectasia mutated, on promoting lipid oxidation and ferroptosis in tumor cells. These evidences indicate that the induction of ferroptosis is expected to enhance the antitumor efficacy of immunotherapy; however, it is not currently clear whether the immunotherapeutic targeting of ferroptosis, is effective for OSCC patients.

The purpose of this study is to explore potential diagnostic and prognostic markers of OSCC and investigate their biological functions through bioinformatics analysis. We successfully constructed a novel prognostic model that can effectively predict the FRGs-OSCC that is mainly involved in the biological processes of immunity and glycolysis. Finally, we also discussed the prediction of a prognostic model on the sensitivity to immunotherapy.

2.1 Data collection

The mRNA expression profiles and the corresponding clinical characteristics of 273 OSCC patients were obtained from The Cancer Genome Atlas (TCGA) database (https://gdc-portal.nci.nih.gov/) (17). A total of 273 oral samples (254 tumors and 19 controls) were retained, including the tongue, alveolus, buccal mucosa, soft and hard palate, oral cavity, and lips. 108 driver genes and 111 marker genes were downloaded from FerrDb database (http://www.zhounan.org/ferrdb/) (18). The genes tested only in mice (12 genes) and multi-annotated genes in both groups (15 genes) were filtered out. As a result, a total of 192 FRGs was obtained. The accession numbers of all samples were included in Additional file 1. The names of public domains and the direct web links were listed in Additional file 2.

2.2 Cox risk regression establishment

A univariate Cox regression analysis was utilized to filter the prognostic associated factors that were closely related to the OS of OSCC patients. Then, we performed a multivariate Cox regression analysis with a stepwise regression analysis to construct a risk model. Finally, 4 FRGs-OSCC were enrolled in a risk cox regression. The risk score was defined as:

According to the median of the risk score, we divided the patients into two groups: a high risk and a low risk group. The K-M and ROC analyses were performed, based on the risk score, using the survival R package. Moreover, univariate, and multivariate Cox regression analyses were used to analyze whether the risk score was an independent prognostic factor. A nomogram was also established, based on the risk score, the pathologic stage, and the M stage, to obtain survival rates of patients at 1, 3, and 5 years.

2.3 Identification of differentially expressed mRNAs

The limma R package (19) was used to screen the differentially expressed genes (DEGs) between high risk and low risk samples, according to the thresholds of |log2 (fold change)| > 2.0 and p <0.05 (19).

2.4 Functional enrichment analysis

We performed a Gene Ontology (GO) enrichment analysis and KEGG pathway analysis using the clusterProfiler and ggplot2 (20). The GO results were composed of 3 parts: the biological processes (BP), the cellular component (CC), and the molecular function (MF). We considered Functional categories were considered when p-value < 0.05.

2.5 Estimation of immune and stromal scores

ESTIMATE (Estimation of Stromal and Immune cells in Malignant Tumor tissues using Expression data) is a novel algorithm, based on ssGSEA (single sample Gene Set Enrichment Analysis) for predicting the level of tumor tissues’ infiltrating immune and stromal cells, based on gene expression profiles (21) (22). Herein, this method was applied to estimate the immune and stromal scores, for each OSCC patient from the two risk groups using the estimate R package. The tumor purity was inferred according to the formula derived from Prof. Yoshihara’s research (21).

2.6 Estimation of immune cell type fractions

CIBERSORT is a deconvolution algorithm that is used to characterize the cellular constitution of complex tissues (23) (24). The LM22 gene signature contains 547 gene expression signatures that can distinguish 22 human hematopoietic cell phenotypes, including natural killer (NK) cells, T cells, myeloid subsets, B cells, and plasma cells (25). The CIBERSORT R package and the txt files of LM22 are available on the CIBERSORT website (http://cibersort.stanford.edu/).We used the CIBERSORT method and LM22 to compare the proportions of 22 infiltrating immune cell types between the two risk groups. For each sample, the sum of all estimates of immune cell type fractions was equal to 1 (26).

2.7 Differential analysis of immunotherapy between the high risk and low risk groups

Immunotherapy targeting immune checkpoints cytotoxic T lymphocyte antigen 4 (CTLA4) or programmed death 1/ programmed cell death-Ligand 1 (PD/PDL1) has been successfully applied as a first-line treatment of OSCC (27)(28)(29). However, it is disappointing that the response rates of this therapy are low, which suggest that only a subset of OSCC patients respond to immunotherapy, and that the efficacy could be enhanced by determining the type of immune checkpoint inhibitors (30)(31). Besides, the major histocompatibility class 1 (MHC1) complex is essential for presenting endogenous cellular antigens to circulating T cells and in initiating specific anti-cancer immune responses (37). In that case, we checked the expression of HLA-genes and immune checkpoint molecules between the high risk and low risk groups. The HLA signatures were downloaded from Nomenclature (http://hla.alleles.org/genes/index.html). The immune checkpoint molecules included CD27, CD274, CTLA4, HAVCR2, ICOS, IDO1, LAG3, PDCD1, PDCD1LG2, and TIGIT. Tumor Immune Dysfunction and Exclusion (TIDE) is a computational framework for predicting the response of the immune checkpoint blockade response (32). The TIDE scores of OSCC patients were downloaded from the TIDE website (http://tide.dfci.harvard.edu), based on the uploaded transcriptome profiles.

2.8 Statistical analysis

The statistical analyses were undertaken using R v3.6.1 and Bioconductor (https://www.bioconductor.org/). All statistical tests were bilateral, and p < 0.05 was considered significant. Besides, the differences were considered significant if *p < 0.05; **p < 0.01; ***p < 0.001, or ****p < 0.0001.

3.1 A prognosis prediction model based on FRGs in OSCC

Using the TCGA database, 10000 genes that were expressed in the tumor samples, were identified by analyzing the transcriptome data of oral cancer of the tongue, alveolus, buccal mucosa, soft and hard palate, oral cavity, and lips. Subsequently, 192 FRGs were sorted out from the FerrDb database (Additional file 3). The results of the two programs were integrated, and only the intersection was selected to increase the specificity. As a result, 103 FRGs were expressed in the TCGA-OSCC (Fig. 1A, Additional file 4).

To evaluate the prognostic value of the OSCC 103 FRGs, a univariate Cox regression analysis was performed on these genes, and 5 variables (FTH1, FLT3, CDKN2A, DDIT3, and BNIP3; all p < 0.01) were identified (Table 1, Additional file 5). Furthermore, we applied to the 5 FRGs, a multivariate Cox regression with a stepwise regression. As a result, 4 FRGs (FTH1, FLT3, CDKN2A, and DDIT3) were harvested in the Cox regression (p < 0.05; Fig. 1B), and the ferroptosis-related risk model was constructed depending on these 4 FRGs.

Table 1

Univariate Cox Regression of prognostic related genes for OS
id	HR	HR.95L	HR.95H	p-value
FTH1	1.001832	1.000681	1.002985	0.00181
FLT3	0.189899	0.061876	0.582804	0.003689
CDKN2A	0.981546	0.96921	0.994039	0.003896
DDIT3	1.013436	1.0042	1.022756	0.004273
BNIP3	1.021741	1.005669	1.038069	0.007842
p-value < 0.01 was considered statistically significant.HR, hazard ratio; OS, overall survival.

3.2 The prognostic value of the ferroptosis-related signature in OSCC

We calculated the risk score of each patient by the previously described formula and classified the patients into a high risk or a low risk group, based on the median of the risk score. The risk score ranking, the survival status of OSCC patients, and the heatmaps of FRGs’ expressions are shown in Figs. 2A, B and C. The abscissas of each graph represented risk scores of OSCC patients. According to the results of the heatmaps, patients in the high-risk group had much higher expression of DDIT3/FTH1 genes and lower expression of CDKN2A/FLT3 genes. Based on the K-M curves of the TCGA dataset, there was an obvious difference between the two risk groups (P = 1.005e-04; Fig. 2D). Furthermore, the ROC curve confirmed the predictive capacity of the 4-FRG risk signature. The area under the curve (AUC) for the risk signature was 0.749, indicating a considerable predicting power (Fig. 2E).

3.3 The 4-FRG risk model was independent of conventional clinical characteristics

The correlation between the ferroptosis-related risk signature and the clinical traits was analyzed. Although the ferroptosis-related risk score was significantly related to histologic grade and T stages, there was no correlation between genders, pathologic stages, M stages, N stages and risk score (Figs. 3A and B, Additional file 6, p < 0.05). Meanwhile, a stratified analysis of the clinicopathological characteristics was further carried out and the results showed that the risk score level was closely related to prognosis (Additional file 7 and 8).

To analyze the relationship between the OS, the clinico-pathological factors, and the ferroptosis related risk signature, a univariate analysis was performed. Subsequently, a multivariate analysis was utilized to explore the independent prediction of the FRGs signature. The results of univariate and multivariate analyses with the Cox proportional hazard model showed that pathologic stage, risk score, and M stage, were independent prognostic indicators for overall survival (Figs. 3C and D). The ferroptosis-related risk model also had a high accuracy in predicting the patients’ 5-year survival rate (Figs. 4A to C).

3.4 Functional enrichment analysis of differentially expressed genes between the two risk groups

To explore the DEGs between the two groups, |log2FC| > 2 and p < 0.05 were set as the screening criteria using the limma R package. As a result, 631 DEGs were identified, including 244 upregulated and 387 downregulated genes (Additional file 9).

The biological functions of the 631 DEGs were speculated by GO enrichment analyses and KEGG pathway analyses with the help of the clusterProfiler R package (Additional file 10 and 11). The top 10 enriched GO terms and all pathways were listed in Fig. 5. It was revealed that the 631 DEGs were significantly enriched in metabolic processes, including glucose catabolic process to pyruvate (GO:0061718), canonical glycolysis (GO:0061621), NADH regeneration (GO:0006735), and glycolytic process through fructose-6-phosphate (GO:0061615) in the biological process category. As for the cellular component category, the top 4 markedly enriched GO terms, were chaperone complex (GO:0101031), chaperonin-containing T-complex (GO:0005832), collagen-containing extracellular matrix (GO:0062023), and eukaryotic 48S preinitiation complex (GO:0033290). The most significantly enriched molecular function terms included alcohol dehydrogenase (NADP+) activity (GO:0008106), and alditol: NADP + 1-oxidoreductase activity (GO:0004032) (Fig. 5A). Also, 11 enriched pathway terms, including Glycolysis/Gluconeogenesis (hsa00010), Glutathione metabolism (hsa00480) and Hypoxia-inducible factor 1 signaling pathway (hsa04066), were explored by KEGG pathway analysis (Fig. 5B).

Notably, a metabolic switch towards aerobic glycolysis, contributes to the activation of dendritic cells by pro-inflammatory stimuli (33) (34). Besides, the decomposing substance (NAD+) in NADH can increase the activity of macrophages to activate other cells of the immune system, and finally, stimulate the entire immune system (35) (36). This evidence indicated that differentially expressed genes between the two groups are closely related to immunity. Therefore, we turned our attention to the difference in the tumor microenvironment (TME) between the two risk groups.

3.5 Differential analysis of the tumor microenvironment between the two risk groups

We calculated the immune and stromal scores, for each OSCC patient from the two risk groups through the estimate R package, finding that patients from the low-risk group tended to have higher immune scores. However, the stromal score between the two groups showed no significant difference (Fig. 6A, B and C). It was also shown that the patients from the low-risk group had a higher immune-cells content in the tumor microenvironment.

After that, we estimated the relative infiltration of the 22 immune cells in the two risk groups by the CIBERSORT R package and LM22. The relative proportion of 22 immune cells in OSCC samples from two groups is significantly different. The presence of resting memory CD4 + T cells, M0 Macrophages, and M2 Macrophages was positively correlated with the level of risk score. By contrast, the presence of naive B cell, activated memory CD4 + T cells, regulatory T cells, T follicular helper cells, CD8 + T cells, M1 Macrophages, gamma delta T cells, and resting mast cells, negatively correlated with the level of risk score (Fig. 6D). The unique differences in immune infiltration between high and low risk groups suggested that the 4-FRG model can be used to predict the prognosis of immunotherapy.

3.6 Differential analysis of the immunotherapy response between the two risk groups

We further investigated the HLA expression and immune checkpoints between the two different groups. The patients from the low risk group possessed significantly higher expression of most HLA compared with the patients from the high-risk group (Figs. 7A), proving that a higher immune status closely correlates with the OSCC prognosis. The expression of the immune checkpoints IDO1, LAG3, PDCD1, and TIGIT significantly increased in the low-risk patients’ group (Figs. 7B). Based on these results, we speculated that OSCC patients with lower risk score might be promising candidates for immune checkpoint inhibitors. Furthermore, The TIDE score in the high-risk group was much higher (Fig. 7C), suggesting that compared with the high-risk patients, the low-risk patients may be more sensitive to immunotherapy. This also confirmed the previous conclusion that the patients with a lower 4-FRG signature risk score, were more suitable for immunotherapy.

To identify potential diagnostic and prognostic markers of OSCC, we constructed a prognostic model based on FRGs. By applying univariate Cox regression and multivariate Cox regression analyses to 103 FRGs that were downloaded from databases, we found 4 FRGs (FTH1, FLT3, CDKN2A, and DDIT3) that were related to the OS of OSCC patients. The FTH1 protein expression was significantly upregulated in breast cancer cells (37) and the epigenetic silencing of FTH1 and TFRC that is induced by estrogen, reduced liver cancer cell growth and survival (38). FLT3 is a receptor tyrosine kinase that plays a crucial role in the development of hematopoietic progenitor cells (39). Furthermore, FLT3 genetic alterations occurred in up to 30% of cases with acute myelogenous leukemia, and patients with FLT3 mutations, have poor outcomes (40). CDKN2A is a tumor suppressor gene that was reported to be frequently altered in OSCC progression (41). CDKN2A low gene expression is associated with the recurrence of disease in oral cancer patients and could be used as a prognostic marker for OSCC (42). DDIT3 is an endoplasmic reticulum stress-responsive transcription factor which plays an important role in apoptotic execution pathways that are induced by the endoplasmic reticulum stress. The speckle-type POZ protein contributes to prostate cancer by targeting DDIT3 (43). Besides, DDIT3 acts as a transcription factor that enhances the expression of TNFRSF10A and TNFRSF10B, resulting in the initiation of ER stress-mediated apoptosis in human lung cancer cells (44). In summary, all 4 FRGs have been reported to be closely correlated with various cancers.

The prognostic model was constructed based on these 4 FRGs. The patients were separated into high risk and low risk groups, according to the threshold of the median risk score. As mentioned above, the patients with higher expression of CDKN2A/FLT3 genes and lower expression of DDIT3/FTH1 genes, were more likely to be in the low risk group, suggesting that FRGs CDKN2A/FLT3 genes may work as tumor suppressor genes. This result is consistent with previous research results (40) (41) (45). The prognostic value of the 4-FRG risk signature was evaluated using the K-M and Log-rank methods. There were significant differences in the survival curves of patients in the two groups of patients. The prediction capability of the specificity and sensitivity of the FRG risk model was assessed by calculating the AUC of the risk score. Moreover, we found that the ferroptosis-related risk score was an independent prognostic indicator for overall survival when considering conventional clinical characteristics. The results indicated that this risk score model is a firm prognostic tool that can be used to classify patients and guide future targeted therapies.

For further understanding of the biological functions of DEGs, between different risk groups, we performed functional enrichment analyses. The results showed significant enrichment in process, including NADH regeneration, glucose catabolic process to pyruvate, canonical glycolysis, and glycolytic process through fructose-6-phosphate, in the biological process category. NADH regeneration is a metabolic process that consumes NAD + to generate a pool of NADH, which is important to the immune system. It is indicated that NAD + promotes the differentiation of CD4 + T cell without antigens. Furthermore, without relying on antigen-presenting cells, NAD + can also regulate the fate of CD4 + T cell (46). Pro-inflammatory stimuli induce the NAD activation of macrophages and dendritic cells, resulting in a metabolic switch towards glycolysis (33) (34), while inflammatory macrophages depend on NAD + salvage, resulting from ROS-mediated DNA damages (35). For KEGG, 11 pathways, including Glycolysis/Gluconeogenesis, Glutathione metabolism, and the HIF-1 signaling pathway were identified. Several types of cancer, including OSCC highly depend on glycolysis for ATP generation. Zheng et al. (47) found that zeste homolog2 can regulate STAT3 and FoxO1 signaling in human OSCC cells and promote invasion and tumor glycolysis. Another study demonstrated that circMDM2 could promote the proliferation and glycolysis of human OSCC cells by acting as ceRNAs to sponge miR-532-3p (48). The HIF-1 signaling pathway, as a cancer therapy glycolytic target is involved in the regulation of glycolysis at pre-clinical and clinical stages (49) (50). Our results indicated that the DEGs, between the two different groups, may affect OSCC progression by altering these immunity-related biological processes or metabolic pathways.

We also focused on investigating the difference in the TME between the two different risk groups. To this end, we explored the correlation between OSCC and tumor microenvironment. The immune score and stromal score calculated based on the ESTIMATE algorithm can help quantify the immune and stromal components in tumors. The immune scores and ESTIMATE scores of the low risk group were significantly higher than those in the high-risk group; however, no significant differences were found in the stromal scores. A high fraction of T cells gamma delta, Macrophages M1, B cell naive, T cells CD4 memory activated, T cells CD8, T cells regulatory, T cells follicular helper, and mast cells resting, mainly infiltrated the tumors of the low risk OSCC patients. A recent study suggested that the cell density of the high parenchymal CD8+, at the invading tumor edge, was associated with improved overall survival, and therefore could be used as an independent favorable prognostic marker for OSCC (51). Moreover, OSSC patients with high levels of CD4 + CD25highCD127low regulatory T cells (Tregs), were found to have a better survival probability compared to patients with lower Tregs. This result indicated that immune cells might have an important effect on the OSCC TME. What’s more, we analyzed the expression of HLA-related genes, important to the immune system, in the two different risk groups, and found that the expression of most HLAs were significantly higher in the low-risk group, demonstrating that higher immune status was related to the prognosis of OSCC. The HLA molecules on the surface of tumor cells can help T cells recognize new antigens to create opportunities for anti-cancer immune responses (52). The expression of the immune checkpoints IDO1, LAG3, PDCD1, and TIGIT significantly increased in the low risk group. Foy et al. (53) found that OSCC tissues are characterized by a higher level of intratumor T-cells, overexpression of PD-L1 and IDO1, and a higher score of response signature to pembrolizumab, suggesting the inhibition of IDO1 and PD-L1 may have good clinical significance for OSCC. Another research indicated that several immune checkpoint receptors (TIM3, LAG3, IDO, PDL1 and CTLA4) could be considered as biomarkers that reflect the immune status in OSCC patients’ TME during nimotuzumab therapy (54). T cells from peripheral blood mononuclear cells, that were collected from OSCC donors, possessed a high expression level of TIGIT. Moreover, TIGIT blockade can promote the in vitro proliferative ability and effective cytokine secretion capacity of CD4 + T cells and CD8 + T cells isolated from OSCC patients (55). These results provided support for the hypothesis that OSCC patients with lower risk score (patients with higher expression of CDKN2A/FLT3 genes and lower expression of DDIT3/FTH1 genes) might respond better to the IDO1, LAG3, PDCD1 and TIGHT inhibitors.

For the first time, we assessed the effects of ferroptosis-related genes on the prognosis of OSCC and constructed an effective prognostic model to reveal the involved biological processes. We also proved that this model can be used as a criterion for determining whether a patient is suitable for immunotherapy. It is reasonable to believe that ferroptosis will potentially serve as a novel strategy for cancer treatment. However, more data from in vivo/in vitro experiments and clinical trials, are needed to elucidate the mechanisms between FRGs and tumor immunity in OSCC.

AUC: area under the curve; BP: biological process; CC: cellular component; CTLA4: cytotoxic T lymphocyte antigen 4; DEG: differentially expressed gene; ESTIMATE: Estimation of Stromal and Immune cells in Malignant Tumor tissues using Expression data; FC: fold change; FRGs: ferroptosis-related genes; FRGs-OSCC: FRGs expressed in OSCC; GO: Gene Ontology; HIF-1: hypoxia-inducible factor 1; HLA: human leukocyte antigen; KEGG: Kyoto Encyclopedia of; K-M: Kaplan-Meier; lncRNA: long non-coding RNA; MF: molecular function; MHC1: major histocompatibility class 1; miRNA: microRNA; mRNA: messenger RNAs; NK: natural killer; OS: overall survival; OSCC: oral squamous cell carcinoma; PD: programmed death 1; PDL1: programmed cell death-ligand 1; ROC: receiver operating characteristic; ssGSEA: single sample Gene Set Enrichment Analysis; TCGA: The Cancer Genome Atlas; TIDE: Tumor Immune Dysfunction and Exclusion; TME: tumor microenvironment; Tregs: regulatory T cells.

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Availability of data and materials

The case IDs and direct web links of all OSCC samples obtained from the TCGA data portal (http://portal.gdc.cancer.gov/) by Bioconductor package TCGAbiolinks were shown in Additional file 1. Data banks/repositories corresponding to all datasets analyzed in this study were listed in Additional file 2.

Competing interests

The authors declare no potential conflicts of interest.

Funding

This work was supported by the Shanghai Committee of Science and Technology, China [No. 19140904800]. The funder had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Authors' contributions

JY: Drafted the manuscript, XZ& ZZ: Bioinformatics analysis and data mining, ZA& JJ: Statistical analysis, MY& YW: Literature review, SL& JL: Edited the manuscript, corresponding author. All the authors read and approved the final manuscript.

Acknowledgments

The authors would like to express their gratitude to EditSprings (https://www.editsprings.com/) for the expert linguistic services provided.

Authors’ information

1 Department of Oral and Maxillofacial Surgery, School and Hospital of Stomatology, Tongji University, Shanghai 200072, China

2 Shanghai Engineering Research Center of Tooth Restoration and Regeneration, Shanghai 200072, China

Bray F, Ferlay J, Soerjomataram I, Siegel RL, Torre LA, Jemal A. Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin. 2018; 68 (6): 394–424.
Wang B, Zhang S, Yue K, Wang XD. The recurrence and survival of oral squamous cell carcinoma: A report of 275 cases. Chin J Cancer. 2013; 32 (11): 614–618.
Gröbe A, Blessmann M, Hanken H, Friedrich RE, Schön G, Wikner J, Effenberger KE, Kluwe L, Heiland M, Pantel K, Riethdorf S. Prognostic relevance of circulating tumor cells in blood and disseminated tumor cells in bone marrow of patients with squamous cell carcinoma of the oral cavity. Clin Cancer Res. 2014; 20 (2): 425–433.
Cooper JS, Pajak TF, Forastiere AA, Jacobs J, Campbell BH, Saxman SB, Kish JA, Kim HE, Cmelak AJ, Rotman M, Machtay M, Ensley JF, Chao KSC, Schultz CJ, Lee N, Fu KK. Postoperative Concurrent Radiotherapy and Chemotherapy for High-Risk Squamous-Cell Carcinoma of the Head and Neck. N Engl J Med. 2004; 350 (19): 1937–1944.
Dixon SJ, Lemberg KM, Lamprecht MR, Skouta R, Zaitsev EM, Gleason CE, Patel DN, Bauer AJ, Cantley AM, Yang WS, Morrison B, Stockwell BR. Ferroptosis: An Iron-Dependent Form of Nonapoptotic Cell Death. Cell. 2012; 149 (5): 1060–1072.
Galluzzi L, Vitale I, Aaronson SA. Molecular mechanisms of cell death: recommendations of the Nomenclature Committee on Cell Death 2018. Cell Death Differ. 2018; 25 (3): 486–541.
Wang Y, Wei ZH, Pan K, Li J, Chen QM. The function and mechanism of ferroptosis in cancer. Apoptosis. 2020. DOI: https://doi.org/10.1007/s10495-020-01638-w
Jiang L, Kon N, Li TY, Wang SJ, Su T, Hibshoosh H, Baer R, Gu W. Ferroptosis as a p53-mediated activity during tumour suppression. Nature. 2015; 520 (7545): 57–62.
Tarangelo A, Magtanong L, Bieging-Rolett KT, Li Y, Ye JB, Attardi LD, Dixon SJ. p53 Suppresses Metabolic Stress-Induced Ferroptosis in Cancer Cells. Cell Rep. 2018; 22 (3): 569–575.
Roh JL, Kim EH, Jang H, Shin D. Nrf2 inhibition reverses the resistance of cisplatin-resistant head and neck cancer cells to artesunate-induced ferroptosis. Redox Biol. 2017; 11: 254–262.
Huang KJ, Wei YH, Chiu YC, Wu SR, Shieh DB. Assessment of zero-valent iron-based nanotherapeutics for ferroptosis induction and resensitization strategy in cancer cells. Biomater Sci. 2019; 7 (4): 1311–1322.
Tang R, Xu J, Zhang B, Liu J, Liang C, Hua J, Meng QC, Yu XJ, Shi S. Ferroptosis, necroptosis, and pyroptosis in anticancer immunity. J Hematol Oncol. 2020; 13 (1): 110.
Sun C, Mezzadra R, Schumacher TN. Regulation and Function of the PD-L1 Checkpoint. Immunity. 2018; 48 (3): 434–452.
Nguyen LT, Ohashi PS. Clinical blockade of PD1 and LAG3 — potential mechanisms of action. Nat Rev Immunol. 2015; 15 (1): 45–56.
Wang WM, Green MD, Choi JE, Gijón M, Kennedy PD, Johnson JK, Liao P, Lang XT, Kryczek I, Sell A, Xia HJ, Zhou JJ, Li GP, Li J, Li W, Wei S, Vatan L, Zhang HJ, Szeliga W, Gu W, Liu R, Lawrence T, Lamb C, Tanno Y, Cieslik M, Stone E, Georgiou G, Chan TA, Chinnaiyan A, Zou WP. CD8+ T cells regulate tumor ferroptosis during cancer immunotherapy. Nature. 2019; 569 (7755): 270–274.
Lang XT, Green MD, Wang WM, Yu JL, Choi JE, Jiang L, Liao P, Zhou JJ, Zhang Q, Dow A, Saripalli AL, Kryczek I, Wei S, Szeliga W, Vatan L, Stone EM, Georgiou G, Cieslik M, Wahl DR, Morgan MA, Chinnaiyan AM, Lawrence TS, Zou WP. Radiotherapy and immunotherapy promote tumoral lipid oxidation and ferroptosis via synergistic repression of SLC7A11. Cancer Discov. 2019; 9 (12): 1673–1685.
Blum A, Wang P, Zenklusen JC. SnapShot: TCGA-Analyzed Tumors. Cell. 2018; 173 (2): 530.
Zhou N, Bao JK. FerrDb: a manually curated resource for regulators and markers of ferroptosis and ferroptosis-disease associations. Database. 2020; 2020: 1–8.
Ritchie ME, Phipson B, Wu D, Hu YF, Law CW, Shi W, Smyth GK. Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015; 43 (7): e47.
Yu GC, Wang LG, Han YY, He QY. ClusterProfiler: an R Package for Comparing Biological Themes Among Gene Clusters. Omi A J Integr Biol. 2012; 16 (5): 284–287.
Yoshihara K, Shahmoradgoli M, Martínez E, Vegesna R, Kim H, Torres-Garcia W, Trevin˜o V, Shen H, Laird PW, Levine DA, Carter SL, Getz G, Stemke-Hale, Mills GB, Verhaak RGW. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. 2013; 4 (1): 2612.
Hamilton MB, Tartakovsky M, Battocletti A. Speed-ne: Software to simulate and estimate genetic effective population size (Ne) from linkage disequilibrium observed in single samples. Mol Ecol Resour. 2018; 18 (3): 714–728.
Newman AM, Liu CL, Green MR, Gentles AJ, Feng WG, Xu Y, Hoang CD, Diehn M, Alizadeh AA. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods. 2015; 12 (5): 453–457.
Chen B, Khodadoust MS, Liu CL, Newman AM, Alizadeh AA. Profiling Tumor Infiltrating Immune Cells with CIBERSORT. Methods in Molecular Biology. 2018; 1711: 243–259.
Vallania F, Tam A, Lofgren S, Schaffert S, Azad TD, Bongen E, Haynes W, Alsup M, Alonso M, Davis M, Engleman E, Khatri P. Leveraging heterogeneity across multiple datasets increases cell-mixture deconvolution accuracy and reduces biological and technical biases. Nat Commun. 2018; 9 (1): 4735.
Bense RD, Sotiriou C, Piccart-Gebhart MJ, Haanen JBAG, Vugt MATM, Vries EGE, Schröder CP, Fehrmann RSN. Relevance of tumor-infiltrating immune cell composition and functionality for disease outcome in breast cancer. 2017; 109: 1–9.
Champiat S, Ferrara R, Massard C, Besse B, Marabelle A, Soria J-C, Ferté C. Hyperprogressive disease: recognizing a novel pattern to improve patient management. Nat Rev Clin Oncol. 2018; 15 (12): 748–762.
Emens LA. Breast Cancer Immunotherapy: Facts and Hopes. Clin Cancer Res. 2018; 24 (3): 511–520.
Galluzzi L, Chan TA, Kroemer G, Wolchok JD, López-Soto A. The hallmarks of successful anticancer immunotherapy. Sci Transl Med. 2018; 10 (459): eaat7807.
Cramer JD, Burtness B, Ferris RL. Immunotherapy for head and neck cancer: Recent advances and future directions. Oral Oncol. 2019; August: 99.
Groisberg R, Hong DS, Behrang A, Hess K, Janku F, Piha-Paul S, Naing A, Fu SQ, Benjamin R, Patel S, Somaiah N, Conley A, Meric-Bernstam F, Subbiah V. Characteristics and outcomes of patients with advanced sarcoma enrolled in early phase immunotherapy trials. J Immunother Cancer. 2017; 5 (1): 1–8.
Jiang P, Gu SQ, Pan D, Fu JX, Sahu A, Hu XH, Li ZY, Traugh N, Bu X, Li B, Liu J, Freeman GJ, Brown MA, Wucherpfennig KW, Liu XS. Signatures of T cell dysfunction and exclusion predict cancer immunotherapy response. Nat Med. 2018; 24 (10): 1550–1558.
Chang CH, Curtis JD, Maggi LB, Faubert B, Villarino AV, Sullivan DO, Huang SCC, Windt GJW, Blagih J, Qiu J, Weber JD, Pearce EJ, Jones RG, Pearce EL. Posttranscriptional Control of T Cell Effector Function by Aerobic Glycolysis. Cell. 2013; 153 (6): 1239–1251.
Kornberg MD, Bhargava P, Kim PM, Putluri V, Snowman AM, Putluri N, Calabresi PA, Snyder SH. Dimethyl fumarate targets GAPDH and aerobic glycolysis to modulate immunity. Science. 2018; 360 (6387): 449–453.
Cameron AM, Castoldi A, Sanin DE, Flachsmann LJ, Field CS, Puleston DJ, Kyle RL, Patterson AE, Hässler F, Buescher JM, Kelly B, Pearce EL, Pearce EJ. Inflammatory macrophage dependence on NAD+ salvage is a consequence of reactive oxygen species–mediated DNA damage. Nat Immunol. 2019; 20 (4): 420–432.
Minhas PS, Liu L, Moon PK, Joshi AU, Dove C, Mhatre S, Contrepois K, Wang Q, Lee BA, Coronado M, Bernstein D, Snyder MP, Migaud M, Majeti R, Mochly-Rosen D, Rabinowitz JD, Andreasson KI. Macrophage de novo NAD+ synthesis specifies immune function in aging and inflammation. Nat Immunol. 2019; 20 (1): 50–63.
Pinnix ZK, Miller LD, Wang W, Agostino RD, Kute T, Willingham MC, Hatcher H, Tesfay L, Sui GC, Di XM, Torti SV, Torti FM. Ferroportin and iron regulation in breast cancer progression and prognosis. Sci Transl Med. 2010; 2 (43): 43ra56-43ra56.
Muhammad JS, Bajbouj K, Shafarin J, Hamad M. Estrogen-induced epigenetic silencing of FTH1 and TFRC genes reduces liver cancer cell growth and survival. Epigenetics. 2020; 00 (00): 1–17.
Yang T, Hu MS, Qi WY, Yang Z, Tang MH, He J, Chen Y, Bai P, Yuan X, Zhang CF, Liu KJ, Lu YL, Xiang ML, Chen LJ. Discovery of Potent and Orally Effective Dual Janus Kinase 2/FLT3 Inhibitors for the Treatment of Acute Myelogenous Leukemia and Myeloproliferative Neoplasms. J Med Chem. 2019; 62 (22): 10305–10320.
Puente-Moncada N, Costales P, Antolín I. Inhibition of FLT3 and PIM Kinases by EC-70124 Exerts Potent Activity in Preclinical Models of Acute Myeloid Leukemia. Mol Cancer Ther. 2018; 17 (3): 614–624.
Padhi SS, Roy S, Kar M, Saha A, Roy S, Adhya A, Baisakh M, Banerjee B. Role of CDKN2A/p16 expression in the prognostication of oral squamous cell carcinoma. Oral Oncol. 2017; 73: 27–35.
Roshan VGD, Sinto MS, Vargees BT, Kannan S. Loss of CDKN2A and CDKN2B expression is associated with disease recurrence in oral cancer. J Oral Maxillofac Pathol. 2019; 23 (1): 82–89.
Zhang P, Gao K, Tang Y, Jin XF, An J, Yu HX, Wang H, Zhang YY, Wang DJ, Huang HJ, Yu L, Wang CJ. Destruction of DDIT3/CHOP Protein by Wild-Type SPOP but Not Prostate Cancer-Associated Mutants. Hum Mutat. 2014; 35 (9): 1142–1151.
Li TL, Su L, Lei YJ, Liu XF, Zhang YJ, Liu XG. DDIT3 and KAT2A Proteins Regulate TNFRSF10A and TNFRSF10B Expression in Endoplasmic Reticulum Stress-mediated Apoptosis in Human Lung Cancer Cells. J Biol Chem. 2015; 290 (17): 11108–11118.
Smeds J, Berggren P, Ma X, Xu ZJ, Hemminki K, Kumar R. Genetic status of cell cycle regulators in squamous cell carcinoma of the oesophagus: the CDKN2A (p16INK4a and p14ARF) and p53 genes are major targets for inactivation. Carcinogenesis. 2002; 23 (4): 645–655.
Rodriguez Cetina Biefer H, Heinbokel T, Uehara H, Camacho V, Minami K, Nian YQ, Koduru S, Fatimy RE, Ghiran I, Trachtenberg AJ, Fuente MA, Azuma H, Akbari O, Tullius SG, Vasudevan A, Elkhal A. Mast cells regulate CD4+ T-cell differentiation in the absence of antigen presentation. J Allergy Clin Immunol. 2018; 142 (6): 1894-1908.e7.
Zheng M, Cao MX, Luo XJ, Li L, Wang K, Wang SS, Wang HF, Tang YJ, Tang YL, Liang XH. EZH2 promotes invasion and tumour glycolysis by regulating STAT3 and FoxO1 signalling in human OSCC cells. J Cell Mol Med. 2019; 23 (10): 6942–6954.
Zheng Z, Ma XZ, Li HF. Circular RNA circMDM2 accelerates the glycolysis of oral squamous cell carcinoma by targeting miR‐532‐3p/HK2. J Cell Mol Med. 2020; 24 (13): 7531–7537.
Semenza G. HIF-1 Inhibitors for Cancer Therapy: From Gene Expression to Drug Discovery. Curr Pharm Des. 2009; 15 (33): 3839–3843.
Gill KS, Fernandes P, O’Donovan TR, McKenna SL, Doddakula KK, Power DG, Soden DM, Forde PF. Glycolysis inhibition as a cancer treatment and its role in an anti-tumour immune response. Biochim Biophys Acta - Rev Cancer. 2016; 1866 (1): 87–105.
Shimizu S, Hiratsuka H, Koike K, Tsuchihashi K, Sonoda T, Ogi K, Miyakawa A, Kobayashi J, Kaneko T, Igarashi T, Hasegawa T, Miyazaki A. Tumor-infiltrating CD8 + T-cell density is an independent prognostic marker for oral squamous cell carcinoma. Cancer Med. 2019; 8 (1): 80–93.
Roerden M, Nelde A, Walz JS. Neoantigens in Hematological Malignancies—Ultimate Targets for Immunotherapy? Front Immunol. 2019; 10 (December): 1–9.
Foy JP, Bertolus C, Michallet MC, Deneuve S, Incitti R, Bendriss-Vermare N, Albaret MA, Cuaran SO, Thomas E, Colombe A, Gadot CPN, Michot JP, Fayette J, Viari A, Eynde BV, Goudot P, Devouassoux-Shisheboran M, Puisieux A, Caux C, Zrounba P, Lantuejoul S, Saintigny P. The immune microenvironment of HPV-negative oral squamous cell carcinoma from never-smokers and never-drinkers patients suggests higher clinical benefit of IDO1 and PD1/PD-L1 blockade. Ann Oncol. 2017; 28 (8): 1934–1941.
Wang H, Mao L, Zhang T. Altered expression of TIM‐3, LAG‐3, IDO, PD‐L1, and CTLA‐4 during nimotuzumab therapy correlates with responses and prognosis of oral squamous cell carcinoma patients. J Oral Pathol Med. 2019; 48 (8): 669–676.
Liu X, Li Q, Zhou Y. Dysfunctional role of elevated TIGIT expression on T cells in oral squamous cell carcinoma patients. Oral Dis. 2020; odi.13703.

Additionalfile1.xls
Additional file 1. TCGA datasets of all OSCC samples obtained from the TCGA. The data included 254 oral cancer samples and 19 matched normal samples
Additionalfile1.xls
Additional file 1. TCGA datasets of all OSCC samples obtained from the TCGA. The data included 254 oral cancer samples and 19 matched normal samples
Additionalfile2.xls
Additional file 2. Data banks/repositories corresponding to all datasets analyzed in this study.
Additionalfile2.xls
Additional file 2. Data banks/repositories corresponding to all datasets analyzed in this study.
Additionalfile3.xls
Additional file 3. 192 ferroptosis-related genes from the FerrDb database.
Additionalfile3.xls
Additional file 3. 192 ferroptosis-related genes from the FerrDb database.
Additionalfile4.xls
Additional file 4. The transcriptome data of 103 ferroptosis-related genes in TCGA-OSCC.
Additionalfile4.xls
Additional file 4. The transcriptome data of 103 ferroptosis-related genes in TCGA-OSCC.
Additionalfile5.xls
Additional file 5. Univariate Cox Regression of prognostic related genes for OS.
Additionalfile5.xls
Additional file 5. Univariate Cox Regression of prognostic related genes for OS.
Additionalfile6.tif
Additional file 6. Relationships between the risk score and clinicopathological parameters in OSCC. (A) Clinical correlation analysis between risk scores and genders. (B) Clinical correlation analysis between risk scores and stages. (C) Clinical correlation analysis between risk scores and metastasis. (D) Clinical correlation analysis between risk scores and the extent of lymph node involvement. (P < 0.05).
Additionalfile6.tif
Additional file 6. Relationships between the risk score and clinicopathological parameters in OSCC. (A) Clinical correlation analysis between risk scores and genders. (B) Clinical correlation analysis between risk scores and stages. (C) Clinical correlation analysis between risk scores and metastasis. (D) Clinical correlation analysis between risk scores and the extent of lymph node involvement. (P < 0.05).
Additionalfile7.tif
Additional file 7. Kaplan-Meier curves of overall survival (OS) between high-risk and low-risk groups. (A) The K-M analysis of OS based on female patients. (B) The K-M analysis of OS based on male patients. (C) The K-M analysis of OS based on Grade1/2 samples. (D) The K-M analysis of OS based on Grade 3/4 samples.
Additionalfile7.tif
Additional file 7. Kaplan-Meier curves of overall survival (OS) between high-risk and low-risk groups. (A) The K-M analysis of OS based on female patients. (B) The K-M analysis of OS based on male patients. (C) The K-M analysis of OS based on Grade1/2 samples. (D) The K-M analysis of OS based on Grade 3/4 samples.
Additionalfile8.tif
Additional file 8. Kaplan-Meier curves of overall survival (OS) between high-risk and low-risk groups. (A) The K-M analysis of OS based on N + samples. (B) The K-M analysis of OS based on N0 samples. (C) The K-M analysis of OS based on M0 samples.
Additionalfile8.tif
Additional file 8. Kaplan-Meier curves of overall survival (OS) between high-risk and low-risk groups. (A) The K-M analysis of OS based on N + samples. (B) The K-M analysis of OS based on N0 samples. (C) The K-M analysis of OS based on M0 samples.
Additionalfile9.xls
Additional file 9. The list of differentially expressed genes.
Additionalfile9.xls
Additional file 9. The list of differentially expressed genes.
Additionalfile10.xls
Additional file 10. The enriched Gene ontology (GO) terms of differentially expressed genes between high-risk and low-risk groups.
Additionalfile10.xls
Additional file 10. The enriched Gene ontology (GO) terms of differentially expressed genes between high-risk and low-risk groups.
Additionalfile11.xls
Additional file 11. The enriched pathways of differentially expressed genes between high-risk and low-risk groups.
Additionalfile11.xls
Additional file 11. The enriched pathways of differentially expressed genes between high-risk and low-risk groups.

Download PDF

Version 1

posted

You are reading this latest preprint version

Comprehensive analysis of the significance of ferroptosis-related genes in the prognosis and immunotherapy of oral squamous cell carcinoma

Status:

Version 1

Abstract

Figures

1. Introduction

2. Materials And Methods

3. Results

4. Discussion

List Of Abbreviations

Declarations

References

Supplementary Files

Status:

Version 1