Mitophagy-Related Gene CHDH Predicts Prognosis and Immune Response and Inhibits Proliferation and Migration in vitro and in vivo of Oral Squamous Cell Carcinoma

Introduction

OSCC is a frequent malignancy of the head and neck. In some areas of the world, both the incidence and death rates are high.¹ In 2018, it caused over 170,000 fatalities and about 350,000 new cases.² The development of OSCC is typically linked to heavy alcohol use, tobacco use, and chewing betel quid.³ OSCC is generally seen in people over 40 years of age. Men are more likely than women to have OSCC globally, at 5.8 and 2.3 per 100,000, respectively.⁴ Surgical excision followed by adjuvant radiation therapy or chemotherapy with radiation therapy are the standard treatments for OSCC. Conformal radiation therapy (CRT) has been the mainstay of treatment for pharyngeal cancer.⁵ Patients’ needs have not been met by current medicines, and just half of OSCC patients make it through five years.⁶ Therefore, it’s critical to identify potential treatment targets, prognostic indicators, and ways to lower medication resistance in OSCC.

Mitophagy is crucial in tumorigenesis, especially cancer therapy.⁷ To maintain intracellular and mitochondrial homeostasis, damaged mitochondria are selectively degraded via autophagic clearance in a process called mitophagy.⁸ According to reports, in tumor cells, aerobic glycolysis and metabolic reprogramming are tightly regulated by mitophagy.⁹ It has a great deal of potential for treating different types of cancer.¹⁰ Dysfunctional mitophagy permits the survival of impaired mitochondria, initially suppressing but eventually accelerating tumor development through mechanisms that shift from protective to oncogenic during cancer progression.¹¹ It has been pointed out that mitophagy triggers mitochondrial damage-induced apoptosis in oral cancer cells,¹² that inducing mitophagy-dependent cellular senescence inhibits oral carcinogenesis,¹³ and that up-regulation of the expression of the proto-oncogene ROS 1 (located in mitochondria) enhances the invasive nature of OSCC,¹⁴ suggesting that mitophagy plays a vital role in OSCC. Although previous studies have found that a reliable prognostic model based on MRGs can be constructed to accurately predict the individualization of OSCC patients,^15,16 more research is needed to develop new therapeutic targets.

Therefore, this study investigated the prognostic and immunologic value between mitophagy and OSCC by exploring the TCGA and GEO databases. Using ex vivo and in vivo research, we selected independent prognostic genes and investigated their expression and functional influence on OSCC. The function of mitophagy in OSCC was further investigated in this research, which may offer fresh treatment options for OSCC patients.

Methods And Materials

Data Download and Pre-Processing

The TCGA database (https://portal.gdc.cancer.gov/) provided us with transcriptome information and clinical profiles of 358 OSCC patients on April 26, 2024 (Supplementary Table 1). The initial count data was transformed into transcripts per million (TPM) and subsequently converted to log2 for further analysis. We obtained the GSE41613 dataset from the GEO (https://www.ncbi.nlm.nih.gov/gds) database. It includes 97 tumor specimens with available clinical data. This dataset underwent both quantile normalization and log2 transformation before being analyzed.

MRGs Identification and Mutation Frequency Analysis

To select MRGs from the GeneCards database (https://www.genecards.org/), a correlation score of more than 1 was applied (Supplementary Table 2). Genes that displayed differential expression between the tumor and normal samples (|logFC| ≥ 1, p-value < 0.05) were chosen from the dataset using the “limma” package. The transcriptome data and the GSE41613 dataset were merged using the “limma” and “sva” packages in R. When merging, we deleted the normal samples from the transcriptome data and used the “limma” package to obtain expression levels for differential analysis. The clinical data from TCGA and GEO were then merged and organized. Using the “survival” and “survminer” packages, survival evaluation was carried out on the combined clinical data and differential gene expression files. Prognostic network graphs were obtained using the “igraph” R package. Using the copy number proof file, we can get the copy number compilation frequency by R software, and then we can get the graphs by using the “RCircos” package, which is able to display the MRGs copy number variation (CNV) locations across 23 chromosomes.

Consistent Cluster Analysis

The prognostic MRGs were subjected to a consistent cluster analysis with a maximum of K = 9 using the R package “ConsensusClusterPlus” and the clusters were validated 1000 times. Following that, based on the best classification of K = 2–9, patients with OSCC were divided into several genetic subgroups for further research. Graphs were created and survival rates across several clusters were compared using the “survival” and “survminer” programs. To see if it was possible to distinguish between distinct clusters, PCA analysis was conducted using the “ggplot2”, “Rtsne” and “umap” packages. To verify the two clusters’ gene expression, we can also create box plots and heatmaps using the R packages “ggpubr”, “reshape2” and “pheatmap”. Using the “GSEABase”, “GSVA” and “pheatmap” packages, GSVA analysis was carried out. It revealed the pathway enrichment difference between cluster B and cluster A. Using the “limma”, “org.Hs.eg.db”, “clusterProfiler” and “enrichplot” packages, GSEA analysis was carried out. The GSEA pathway analysis shows the pathways’ functional activity in clusters A and B.

Constructing a Prognostic Model

Univariate Cox regression analysis was used to select MRGs associated with overall survival (OS) using the “glmnet” R software package. A risk model was created based on the MRGs’ predictive capacities after the MRGs that independently predicted the prognosis of OSCC were selected using multivariate Cox regression analysis. Patients with OSCC were divided into low-risk and high-risk groups based on their median risk score. The OS rates of OSCC patients in the low-risk and high-risk groups were calculated using the R packages “survival” and “survminer” which were based on the survival curves. The clinicopathologic characteristics of OSCC and risk score diagnostic accuracy were evaluated using the R package “timeROC”. Using the “pheatmap” package, a risk heatmap was created to determine whether the genes are high- or low-risk. The R packages “ggplot2”, “ggalluvial”, “dplyr” and “ggpubr” were used to carry out the risk variance analysis. For risk variance analysis, we can get box plots and Sankey plots.

Tectonic Nomogram

We created a nomogram using the R packages “rms” and “regplot” combining risk scores with additional clinicopathologic data. Each variable receives a score in the nomogram scoring system. The scores of all variables were summed to calculate the overall rating for every sample, estimating the survival rates for OSCC patients at one, three and five years. The association between the actual observations and the 1-, 3- and 5-year survival events that were anticipated was shown using calibration plots. Using the R package “ggDCA” decision curves were created to show the accuracy of patient survival forecasts at 1, 3 and 5 years.

Drug Sensitivity Analysis and the Immunological Microenvironment Landscape

Using the R packages “ggpubr”, “pheatmap”, “vioplot” and “corrplot” we get histograms for high and low risk, heatmaps for correlation and violin plots for variance. The risk file and the immune cell infiltration file were processed using the R packages “limma”, “reshape2”, “tidyverse”, “ggplot2”, “ggpubr” and “ggExtra” to produce immune cell correlation heatmaps. ESTIMATE assessed stromal cell abundance (Stromal Scores) and immune cell abundance (Immune Scores). Using the R software “estimate” we calculated the stromal, immune, and ESTIMATE scores of OSCC. Utilizing the R packages “ggpubr” and “reshape2”, we produced violin plots for difference analysis as well. We obtained the drug sensitivity scoring file using the R software “oncoPredict” for any drug sensitivity difference analysis using R package “ggplot2” and “ggpubr” to get the drug box plot. Sensitivity difference analysis to get the box plot of the drug.

Looking for Separate Genes That Predict Outcomes

To identify independent prognostic genes, both univariate and multivariate Cox regression analyses were performed. The results, including hazard ratios (HRs), 95% confidence intervals (CIs), and corresponding p-values, were visualized using forest plots. The R package “ggplot2” was employed to analyze and compare the expression patterns of these prognostic genes between OSCC tissues and adjacent normal tissues. Additionally, Kaplan-Meier survival analysis combined with the Log rank test and univariate Cox regression was conducted to assess survival differences between groups with high and low expression of the identified genes. Statistical significance was defined as p < 0.05, with HR and 95% CI reported for significant findings.

Tissue Samples

Between January and December 2024, twelve matched pairs of OSCC and adjacent normal tissue samples were collected from the Affiliated Cancer Hospital of Xiangya School of Medicine, Central South University.

Immunohistochemistry

The tissue samples were processed into 6-μm-thick paraffin-embedded sections and fixed in 4% paraformaldehyde. Following deparaffinization and rehydration, the sections were prepared for immunohistochemical staining. To unmask antigens and block endogenous peroxidase activity, the rehydrated sections were treated with Tris-EDTA buffer (1 mm EDTA, 10 mm Tris-HCl, pH 9.0) and subjected to heat-induced epitope retrieval (HIER) by boiling in a pressure cooker for 5 minutes. After three washes, the sections were blocked with BSA for 30 minutes and then incubated overnight at 4°C with primary antibodies against CHDH, Ki67, PCNA, or E-Cadherin. The next day, the sections were incubated with a secondary antibody for 1 hour at room temperature. Finally, the nuclei were counterstained with hematoxylin, and the sections were dehydrated, cleared, and mounted using neutral resin.

Quantitative Real-Time PCR

Total RNA was extracted using the RNeasy Mini Kit (QIAGEN, Beijing, China). Subsequently, cDNA synthesis was performed with the High-Capacity cDNA Reverse Transcription Kit (Thermo Fisher Scientific, Shanghai, China). For quantitative PCR (qPCR), reactions were carried out using TaqMan Gene Expression Master Mix (Bio-Rad, Shanghai, China) along with TaqMan probes specific for human CHDH and GAPDH (Sangon Biotech, Shanghai, China), following the manufacturer’s protocols.

Western Blot

Following 48 hours of transfection with either the control lentiviral plasmid (OC) or the CHDH-overexpressing plasmid (OE), total cellular protein was extracted. The lysates were resolved by SDS-PAGE using dual-color protein standards (Bio-Rad, Shanghai, China, Cat. #1610374) for molecular weight reference. Subsequently, proteins were transferred to PVDF membranes for immunoblotting. Membranes were probed with a primary antibody against CHDH, followed by incubation with a peroxidase-conjugated secondary antibody. Protein-antibody interactions were visualized using the ChemiDoc Go Imaging System (Bio-Rad), and band intensities were quantified using ImageJ software (NIH).

MTT

Prior to transfection, 2×10⁵ cells per well were seeded in 6-well plates and cultured for 24 hours until reaching 30–40% confluency. Lentiviral particles carrying either the control vector (OC) or CHDH-overexpressing construct (OE) were added to the wells in 1 mL of medium supplemented with 40 μL of transfection reagent. After 12 hours of incubation, the viral medium was replaced with fresh complete medium. To select successfully transfected cells, puromycin treatment (concentration not specified) was applied for 72 hours to eliminate untransfected cells. Post-selection, cells were trypsinized and plated in a 96-well plate at a density of 6,000 cells/well for an additional 72 hours. For MTT cell viability assessment, 50 μL of MTT solution (2 mg/mL, final concentration 0.5 mg/mL) was added to each well, followed by 4 hours of incubation at 37°C. The medium was then carefully aspirated, and 150 μL of DMSO was added to solubilize the formazan crystals. Plates were gently shaken for 10 minutes, and absorbance was measured at 490 nm using a microplate reader.

Clonogenic Assay

For lentiviral transfection, 2×10⁵ cells were seeded per well in 6-well plates and cultured for 24 hours to achieve 30–40% confluency. Cells were then transduced with CHDH-overexpressing lentivirus (OE) or control vector (OC) in 1 mL medium containing 40 μL transfection reagent. After 12 hours, the viral supernatant was replaced with fresh complete medium. Following 72 hours of incubation, untransfected cells were eliminated through puromycin selection (note: “piphofomycin” appears to be a typographical error and should be corrected to puromycin). The successfully transfected cells were subsequently plated in 24-well plates at a density of 1,000 cells per well for colony formation analysis. During days 6–8 post-plating, cells were fixed with 10% neutral buffered formalin for 15 minutes, washed twice with PBS, and stained with 0.1% crystal violet solution for 30 minutes at room temperature. After thorough washing with distilled water, the stained colonies were solubilized in 10% acetic acid and absorbance was measured at 550 nm to quantify colony formation capacity.

Scratch Assay

For the wound healing assay, cells were seeded in 6-well plates (2×10⁵ cells/well) and cultured for 24 h to reach 30–40% confluency. Lentiviral transduction was performed using 1 mL medium containing CHDH-overexpressing vectors and 40 μL transfection reagent. After 12 h incubation, the viral supernatant was replaced with fresh complete medium. Following 72 h of transduction, untransfected cells were removed by puromycin selection (2 μg/mL). Selected cells were then replaced in 12-well plates at 3×10⁵ cells/well. After 24 h, a uniform scratch wound was created using a sterile pipette tip, washed twice with PBS, and maintained in serum-free medium. Cell migration was assessed by capturing images at 0 h and 24 h post-scratching using phase-contrast microscopy (40×). Migration distance was quantified by measuring wound width at three predetermined locations per well using ImageJ software.

Transwell Assay

For the migration assay, 2×10⁵ cells were seeded in 6-well plates and cultured for 24 h to reach 30–40% confluency before lentiviral transduction with CHDH-overexpressing vectors (1 mL medium containing 40 μL polybrene). After 12 h, the viral supernatant was replaced with complete medium, followed by 72 h puromycin selection (2 μg/mL) to eliminate untransfected cells. Subsequently, 40,000 viable cells were suspended in serum-free medium and loaded into the upper chamber of transwell inserts (8 μm pore), while the lower chamber contained 10% FBS as chemoattractant. Following 24 h incubation, migrated cells on the membrane’s lower surface were fixed with 4% paraformaldehyde, stained with 0.5% crystal violet, and quantified by counting five random fields per insert under 100× magnification. Non-migratory cells were removed from the upper surface prior to imaging.

In vivo Tumor Models

We acquired nude mice (Strain: BALB/c, Age: five to six weeks old, Gender: female) from Skajingda Biotechnology Company, ensuring they were raised in a meticulously maintained, pathogen-free environment. Temperatures were from 22 to 26 °C, humidity was about 55%, and there was a 12-hour light/dark cycle. To prevent bias, the mice were randomly assigned to several experimental groups, with five mice housed in each cage. They were subcutaneously injected with 2×10^6 SCC-9-OC or SCC-9-OE cell suspensions to create xenograft animals. Once tumors reached 70–100 cm^3, tumor size was determined using the formula 1/2 × (long diameter × short diameter²), and mouse weights were measured every other day. China’s legal laboratory animal usage and care regulations were adhered to by the experiments.

CD8⁺ T Cell Cytotoxic Assay

The procedures outlined in the literature were used to acquire CD8⁺ T cells.¹⁷ After transfection with either control or overexpressed CHDH plasmid packed in a lentivirus, tumor cells were cultured for a full day. They were then plated at a density of 10,000 cells per well in a 96-well plate after this incubation time. CD8⁺ T cells were added to the 96-well plate after 12 hours. The culture media was withdrawn after the tumor cells had been co-incubated for 48 hours, and the T cells were eliminated by giving them two PBS washes. Finally, the MTT test was used to determine the killing power of CD8⁺ T lymphocytes on tumor cells.

Statistical Analysis

All gene data in this study were log-transformed to achieve normalization. PERL annotated and organized transcriptomic, gene and clinical expression data. R software was used for some statistical operations; specifically, versions 4.2.1 and 4.4.0 were used. GraphPad Prism 8.0 software was used to analyze the RT-PCR data of this study. The quantitative data were analyzed using the “mean ± standard deviation” method. The data between the two groups were analyzed using a t-test. One-way analysis of variance was used for comparison between multiple groups. All experiments were repeated three times, and p-value less than 0.05 were considered statistically significant.

Results

Research Process

We downloaded 358 OSCC data from TCGA-HNSCC, and in addition, we added 97 clinical data from GEO-GSE411613. 1686 MRGs were obtained from the GeneCards database. Combining the 2 sets of data for difference analysis, we received 298 differentially expressed MRGs. Univariate Cox analysis of differentially expressed genes allowed us to obtain 108 genes associated with OSCC prognosis. These 108 prognostically related genes were subjected to consensus cluster analysis, then GSVA, GSEA and single-sample gene set enrichment analysis (ssGSEA). To verify MRG features’ effectiveness in creating a risk model, 8 genes were chosen from these 108 genes using multivariate Cox analysis and LASSO regression analysis, followed by K-M analysis, column-line graph analysis and immune relationship analysis. Finally, these 8 genes were screened for independent prognostic genes, and we discovered that CHDH was an independent OSCC predictive gene. As a result, we validated OSCC’s functionality in vitro, such as IHC, PCR, WB, clonogenic assay, MTT, transwell assay, scratch assay, in vivo tumor models and CD8⁺ T cell cytotoxic assay-related processes, as Supplementary Figure 1 shows.

Identification of MRGs and Mutation Frequency Analysis

1686 MRGs were acquired from the Genecards database. We obtained 358 samples in OSCC and nearby normal tissues from the TCGA database, together with 97 samples obtained from the GEO database. With differential expression analysis, we were able to identify 298 genes involved in mitophagy that are differentially expressed. The heatmap and volcano diagram showed the top 100 genes, of which 28 were down-regulated and 72 were up-regulated in OSCC (Figure 1A and B). Considering the contribution of these DEGs to OSCC, we performed GO and KEGG enrichment analysis on these DEGs. In GO enrichment analysis, these genes were mainly involved in ribose phosphate metabolic process, cell-substrate junction, ATP hydrolysis activity and other processes (Supplementary Figure 2A). In KEGG enrichment analysis, these genes were mainly involved in Viral carcinogenesis, Cell cycle, Human immunodeficiency virus 1 infection and other processes, which play an important role in the development of OSCC (Supplementary Figure 2B). We conducted a univariate Cox regression analysis on these 298 differentially expressed genes to develop a more precise characterization of MRGs, and we found 108 genes linked to the prognosis of OSCC (Figure 1C). Nine out of the 108 genes were indicated as low-risk genes by the network diagram (SPATA18, RNF31, CORO1A, CHDH, CDKN2A, FYCO1, IGKC, ACADSB, XRN1), and there was a positive regulatory relationship among the other genes (Supplementary Figure 3). The magnitude of the relationship between the frequency at which copy numbers are acquired and lost can be seen in Figure 1D and Supplementary Figure 4 shows the locations of 108 gene mutations on 23 pairs of chromosomes.

Figure 1 Identification and mutation frequency analysis of MRGs and screening of prognostic-related genes and mutation frequency analysis. (A) Heatmap showing 100 differentially expressed genes in OSCC and adjacent normal tissues. (B) Volcano diagram shows the differentially expressed genes with the threshold set at |logFC| ≥ 1 and p-value < 0.05. (C) The forest map displays the top 108 MRGs with high or low risk. (D) Copy number variation frequency of 108 MRGs in TCGA-OSCC.

Consensus Clustering and Immune Microenvironment Landscape Analysis

Divide OSCC patients into multiple subgroups based on 108 MRGs. According to the results, K = 2 was the best classification for patients with OSCC (Figure 2A–C). Group A included 186 samples, while Group B included 263 samples. After analyzing the two subgroups’ survival rates, group A outperformed group B regarding survival duration (Figure 2D). Therefore, it was known that groups A and B could be distinguished based on the expression of mitophagy. According to differential analysis, these 108 genes had high expression levels in typology B (Figure 2F). Additionally, the heatmap displayed the high expression of these 108 genes in typology B (Figure 2E). The PCA of 108 samples revealed a significant separation between groups A and B (Figure 2G). Immunocyte differential analysis was performed on both sets of samples. The ssGSEA revealed a significant difference of 17 immune cells between the two typologies. T follicular helper cell, Type 17 T helper cell, Activated CD8 T cell, Type 1 T helper cell, Monocyte, MDSC, Mast cell, Activated B cell, Activated CD4 T cell, Eosinophil, Immature B cell, Plasmacytoid dendritic cell, and Immature Dendritic cell in cluster A had high expression levels; In cluster B, there was a high expression of CD56dim natural killer cell, Gamma delta T cell, CD56bright natural killer cell and Natural killer T cell (Figure 2H). GSVA identified disparities in pathway enrichment between clusters B and A. KEGG-BLADDER-CANCER, KEGG-PATHOGENIC-ESCHERICHIA-COLI-INFECTION, KEGG-NOTCH-SIGNALING-PATHWAY, KEGG-GLYCOSAMINOGLYCAN-BIOSYNTHESIS-CHONDROITIN-SULFATE, KEGG-REGULATION-OF-ACTIN-CYTOSKELETON, KEGG-FOCAL-ADHESION, KEGG-SMALL-CELL-LUNG-CANCER and KEGG-RENAL-CELL-CARCINOMA showed a strong enrichment in cluster B. The other pathways, however, were significantly more enriched in cluster A (Figure 2I). GSEA pathway analysis showed that KEGG-FOCAL-ADHESION and KEGG-ECM-RECEPTOR-INTERACTION were functionally active in cluster B. At the same time, KEGG-PRIMARY-IMMUNODEFICIENCY, KEGG-DRUG-METABOLISM-CYTOCHROME-P450, and KEGG-METABOLISM-OF-XENOBIOTICS-BY-CYTOCHROME-P were functionally silent in cluster B (Figure 2J). These findings imply that MRGs could offer fresh perspectives on immune infiltration and response in OSCC.

Figure 2 Analysis of consensus cluster and the immune microenvironmental landscape. (A–C) The consensus matrixes for the legend, k = 2 and k = 3, were obtained by applying consensus clustering. (D) Overall survival of two subtypes (p < 0.001). (E) A heatmap showing the expression of MRGs and the associated clinicopathological characteristics of two subgroups. (F) MRGs expression in two subtype clusters. (G) The PCA of 108 samples showed significant separation between cluster A and cluster B. (H) The results of single sample gene set enrichment analysis (ssGSEA) between two subtype clusters. (I) Gene Set Variation Analysis (GSVA) reveals enrichment differences in pathways between Cluster B and Cluster A. (J) Gene Set Enrichment Analysis (GSEA) pathway analysis shows the functional activity of pathways in clusters A and B. The significance levels were set as *p < 0.05, **p < 0.01, and ***p < 0.001.

Constructing a Prognostic Model

In order to create the risk model using multivariate Cox regression analysis, eight MRGs (SHMT2, GFPT2, CHDH, FADS3, CISD2, SFXN1, XRN1, and PYGL) that could independently forecast OSCC prognosis were chosen. These eight prognostic genes’ regression analysis results (Supplementary Table 3) were applied to develop the prognostic model with a risk score (Figure 3A and B). Patients with OSCC were split into the Train and Test groups at random, and for survival analysis, patient groups were categorized as high-risk or low-risk based on their median risk ratings. The findings demonstrated that patients in the low-risk group had a greater survival rate as compared to those in the high-risk group, with p < 0.05 suggesting that the developed model could correctly distinguish between patients in the two categories (Figure 3C–E). The 1-, 3- and 5-year survival rates for the Train group had AUC values of 0.769, 0.774, and 0.738, respectively, based on the ROC curves. The 1-, 3- and 5-year survival rates for the Test group showed AUC values of 0.599, 0.697, and 0.606 respectively. AUC values for the All group’s 1-, 3- and 5-year survival rates were 0.679, 0.732, and 0.671, respectively. The developed model could correctly predict the patients’ survival duration, as shown by the AUCs of the Train group, the Test group, and the All group being > 0.5 (Figure 3F–H). The heat map of risks indicated that XRN1 and CHDH were low-risk genes, whose expression decreased as the patient’s risk increased, and SFXN1, SHMT2, CISD2, PYGL, GFPT2, and FADS3 were high-risk genes, whose expression increased as the patient’s risk increased (Figure 3I). When the difference between categories A and B was analyzed, it was found that there was a significant difference (p < 0.05) between the risk scores of the two patient groups. The risk score of the patient was higher in typology B (Figure 3J). Figure 3K illustrates the relationship between patient type and survival status. These results suggest that the model based on eight risk genes can accurately predict the prognosis of OSCC patients.

Figure 3 Construction of a predictive model. (A and B) LASSO regression analysis shows the minimum lambda and optimal coefficients of the prognostic MRGs. (C–E) The K-M curves showed different prognoses in the different risk groups. (C) Train group, (D) Test group, (E) All group. (F–H) The time-dependent ROC curves for OS at 1-, 3- and 5-years. (F) Train group, (G) Test group, (H) All group. (I) The risk heatmap divides 8 MRGs into low-risk and high-risk genes. (J) Risk scores for A cluster and B cluster. (K) The correspondence between patient type and survival status.

Constructing a Nomogram

We created a nomogram utilizing risk scores for eight genes together with additional clinicopathologic variables to predict the likelihood of survival in OSCC, given the significance of risk scores in evaluating patient survival. There was a score for each feature; the patient’s risk score was 135 for that reason. The patient survival rates at one, three, and five years were 0.570, 0.194, and 0.0847, as indicated in Figure 4A. Compared to the optimal model for the total cohort, the 1-, 3-, and 5-year OS curves displayed quite good results, according to the prediction plots (Figure 4B). Over time, the risk to patients rose, with the risk being higher for those in the high-risk group compared to those in the low-risk group (Figure 4C). The developed nomogram performed better in predicting patient survival than other clinical parameters, according to decision curves (Figure 4D–F). These findings imply that OSCC patients’ survival can be reliably predicted by a risk score derived from eight risk genes.

Figure 4 Constructing a nomogram. (A) Nomogram plot based on score and clinicopathological factors. (B) Calibration plot for the validation of the nomogram. (C) The cumulative hazard curve represented the probability of survival over time progression. (D–F) DCA curves of the nomogram for 1-, 3- and 5-year OS in OSCC patients.

Analysis of the Tumor Microenvironment and Immune Cell Infiltration

Since the immune response is essential to the fight against cancer, we contrasted the immune cell makeup of the groups at high and low risk. The results showed a substantial variation in immune cell content (Figure 5A). Figure 5B displays the correlation heat map between immune cells. Supplementary Figure 5 provides evidence that Dendritic cells activated, Mast cells activated, B cells naive, Macrophages M0, Dendritic cells resting, Neutrophils, NK cells resting, T cells follicular helper, Mast cells resting, Plasma cells, T cells CD4 memory activated, T cells CD8, T cells CD4 memory resting, T cells gamma delta, and T cells regulatory were notably dissimilar across the two groups. Analysis of the relationship between immune cells and risk scores for patients revealed that Mast cells activated, NK cells resting, Macrophages M0, Neutrophils, and T cells CD4 memory resting had a positive correlation with risk ratings. In contrast, Plasma cells, B cells naive, T cells CD4 memory activate, T cells follicular helper, Mast cells resting, T cells gamma delta, Dendritic cells resting, T cells CD8 and T cells regulatory had a negative correlation with risk ratings (Figure 5C). The high- and low-risk groups’ immune scores varied significantly, based on tumor microenvironment (TME) difference analysis. Moreover, in comparison to the low-risk group, the high-risk group’s TME score was much lower (Figure 5D). The interaction between immune cells and these eight MRGs is depicted in the heat map (Figure 5E). Finally, MRGs might offer fresh perspectives on how the immune system responds and infiltrates the OSCC.

Figure 5 Immune cell infiltration and tumor microenvironment analysis. (A) Immune cell content between high-risk and low-risk groups. (B) Correlation heatmap between immune cells. (C) Correlation analysis between immune cells and patient risk score. (D) Estimate the expression profile scores of high-risk and low-risk groups. (E) The heatmap shows a correlation between 8 MRGs and immune cells. The significance levels were set as *p < 0.05, **p < 0.01, and ***p < 0.001.

Drug Sensitivity Analysis

Medication therapy is a crucial component in treating tumors, so we analyzed 198 drugs’ sensitivity associated with prognostic genes in OSCC patients (Supplementary Table 4). The outcomes demonstrated that, for 81 medications, the high- and low-risk categories differed significantly from one another. The high-risk group had substantially greater IC50 values for 72 medications, suggesting that their drug sensitivity was less than the group at low risk (Figure 6A–E). The high-risk group showed significantly lower IC50 values for nine medications (AZD1332, Dasatinib, Docetaxel, Entospletinib, BPD-00008900, Paclitaxel, Sepantronium bromide, Luminespib, and Staurosporine) in contrast to the group at low risk. This implies that these drugs had a greater effect on the high-risk group (Figure 6F–N). According to these findings, medications with high sensitivity in the high-risk population could be able to offer fresh suggestions for clinical dosage.

Figure 6 Drug sensitivity analysis. IC50 values were calculated for patients in the high- and low-risk groups to assess the sensitivity of chemotherapeutic agents. (A–E) Drugs that are highly sensitive in high-risk groups. (F–N) Drugs that are highly sensitive in low-risk groups.

Finding Independent Prognostic Genes and Validating Them Initially

Then, for eight predictive gene expressions and clinical features, we carried out univariate and multivariate Cox analyses on the risk coefficients, p-values, CI, and HR. Figure 7A and B show how CHDH, Grade, and pN staging may each function as separate predictors of OSCC. In terms of survival duration, the CHDH high-expression group outlasted the low-expression group (Figure 7C). The expression of CHDH varies between nearby normal tissues and OSCC tissues, and Figure 1A shows that CHDH was significantly expressed in neighboring normal tissues but poorly expressed in OSCC tissues. It implies that as a stand-alone prognostic gene in OSCC, CHDH may be important.

Figure 7 Screening for independent prognostic genes. (A and B) Hazard ratio and p-value of the constituents involved in univariate and multivariate Cox regression considering clinical parameters and 8 prognostic MRGs in OSCC. (C) Comparison of OS of CHDH in OSCC high expression group and low expression group.

Validation of CHDH Expression and Prognostic Function in OSCC

In order to investigate CHDH expression in OSCC, we identified CHDH in three adjacent tissues as well as three OSCC tissues. When comparing OSCC to nearby normal tissues, the IHC data revealed a considerable reduction in CHDH (Figure 8A). Next, we used WB and PCR to assess CHDH expressions in two OSCC and one normal cell line. We discovered OSCC had much lower amounts of CHDH protein and mRNA than normal cells, as Figure 8B and C illustrates. According to WB data, OSCC had much less CHDH than surrounding tissues (Figure 8D), which was in line with what IHC revealed. We looked into CHDH’s function in OSCC in more detail. MTT, clone formation, scratch, and transwell assays were used to find the effects of CHDH on OSCC after WB and PCR confirmed the successful overexpression of CHDH (OE) (Figure 9A and B). The findings demonstrated that CHDH overexpression dramatically reduced cell migration, clone formation, and proliferation (Figure 9C–G). Activated T cells were more likely to destroy human OSCC cells when CHDH was overexpressed, according to an interesting finding from co-incubating these cancer cells with CD8⁺ T cells (Figure 9H). Further studies revealed that overexpression of CHDH significantly enhanced the expression of LC3II/ LC3I (Figure 9I). These results highlight the pivotal function of CHDH in OSCC cell proliferation, metastasis, and immune evasion by activating autophagy.

Figure 8 CHDH is low expressed in OSCC tissue and cells. (A) Immunohistochemistry determined the expression of CHDH in OSCC and adjacent tissues. (B) Western Blot determined the expression of CHDH in two OSCC cells and one normal cell line. (C) The expression of CHDH in two OSCC cells and one regular cell line was determined by PCR. (D) The expression of CHDH in OSCC and adjacent tissues was determined by Western Blot.

Figure 9 Validation of CHDH prognostic function in OSCC. (A) Western Blot measured the expression of CHDH after overexpressing CHDH in tumor cells. (B) PCR measured the expression of CHDH after overexpressing CHDH in tumor cells. (C) MTT was used to evaluate cell proliferation after overexpressing CHDH in tumor cells. (D) Clonogenic assay was used to evaluate colony suppression after overexpressing CHDH in tumor cells. (E–G) Scratch and transwell assay evaluated migration suppression after overexpressing CHDH in tumor cells. (H) OC and OE tumor cells’ cell survival was treated with CD8+ T cells. (I) Western Blot measured the expression of LC3 after overexpressing CHDH in tumor cells. The significance levels were set as *p < 0.05, ***p < 0.001.

CHDH Overexpression Inhibits OSCC Cell Growth and Migration in vivo

With a nude mouse subcutaneous xenograft model, we also looked into how overexpression of CHDH affects tumor development in living mice. Tumor growth was considerably suppressed by overexpressing CHDH, as seen by the injection of cells overexpressing CHDH (OE) and non-targeted control (OC) (Figure 10A and B). Furthermore, the OE group tumor weight was significantly lower than the OC group tumor (Figure 10C). Notably, CHDH overexpression did not affect the mice’s average body weight (Figure 10D). The number of proliferation-related proteins Ki-67 and PCNA-positive cancer cells compared to the OC group, was considerably lower in the OE group, according to immunohistochemical results (Figure 10E). Furthermore, the results of the immunohistochemistry analysis verified that the overexpression of CHDH inhibited the in vivo expression of the migration-related protein E-Cadherin (Figure 10E). Our data provides compelling evidence that OSCC cell growth and migration are inhibited in vivo by overexpressing CHDH, which does this via downregulating the expression of PCNA, E-Cadherin, and Ki67.

Figure 10 CHDH overexpression inhibited the growth and migration of OSCC cells in vivo. (A and B) 2×10^6 SCC-9 -OC and SCC-9 -OE cell suspensions were injected into the right flank of mice. When the tumor volume reached 70–100 mm^3 (10 days after cells were injected into the mice flank), each mouse’s tumor volume was measured every two days. 24 days after cells were injected into the flank of mice, tumors were removed and photographed. (C) Statistical analysis of tumor weight at 24 days. (D) Changes in each group of mice’s weight. (E) The expression of CHDH, Ki67, PCNA and E-Cadherin in tumor tissues was detected by immunohistochemistry.

Discussion

Clinically speaking, the most prevalent sites of OSCC are the lips, tongue, and floor of the mouth.¹⁸ Because the disease progresses almost entirely without symptoms, some lesions in the early stages may go unrecognized or ignored, despite the fact that the oral cavity is a highly accessible area for both medical assessment and self-examination. Consequently, approximately half of OSCC patients were found in advanced stages (stage III or IV), suggesting a bleak prognosis and elevated death rate.^19–21 Therefore, the early diagnosis and risk assessment of OSCC patients are necessary to increase their survival rate. Eight MRGs were found to be connected to OS in OSCC in our investigation, and these were then applied to develop a new risk model for OSCC prognosis. Using immunological cell infiltration and drug susceptibility analysis, the prognosis-related genes in OSCC were evaluated further. Basic research further established the expression and operation of the independent prognostic gene CHDH. It suggests that CHDH is an effective therapeutic target for OSCC, which may help solve the complex problem of poorer prognosis and higher mortality in OSCC.

Mitochondria are the primary site of aerobic respiration in eukaryotic cells. Mitophagy is a selective form of autophagy. The PINK1-PRKN/Parkin pathway is the primary regulatory mechanism of mitophagy.²² Although mitophagy has been linked to several malignant cancers, including pancreatic cancer,²³ lung cancer,²⁴ and ovarian cancer,²⁵ discoveries of its role in OSCC are scarce. We created a predictive model of OSCC using MRGs in this study, and we also looked into how OSCC’s mitophagy genes functioned.

We employed the LASSO technique to identify these genes through univariate and multivariate Cox regression analysis in order to reliably predict the prognosis of OSCC patients. We then constructed a risk model and identified 8 genes independently predicting the OSCC prognosis. Previous investigations have identified some connections between these indicators and the development and progression of tumors. An essential enzyme in the metabolism of single-carbon units that catalyzes the conversion of serine to glycine in OSCC is serine hydroxymethyltransferase (SHMT2).²⁶ Hypoxia is associated with mitochondrial metabolic activity in the OSCC microenvironment and is indicative of a poor prognosis, has been demonstrated to encourage immune escape in the tumor microenvironment.²⁷ Up-regulation of SHMT2 may predict a poor prognosis for OSCC patients.²⁸ Glutamate-fructose-6-phosphate transaminase 2 (GFPT2) has been linked to the development of cancer, including stomach, colon, and pancreatic cancers, according to an increasing number of studies.^29–31 GFPT2’s function in OSCC, however, is still unclear. Much research has been documented on CHDH. Transmembrane protein CHDH is found in mitochondria and is responsible for mitophagy following mitochondrial injury. It is a significant catalytic enzyme that catalyzes the transformation of betaine from choline.³² CHDH may be linked to head and neck squamous cell carcinoma (HNSCC).³³ Using total RNA sequencing (RNA-seq) data from TCGA and confirmed in the GEO dataset, the findings revealed that lysophospholipid acyltransferase 1 (LCLAT1) and CHDH in HNSCC were related to multiple prognostically relevant genes, with CHDH having the lowest risk.³³ OSCC is a small branch of HNSCC. Thus, it is possible that CHDH and the emergence of OSCC are connected. In our study, CHDH was confirmed to be a prognostically relevant gene in OSCC, which is additional evidence that the development of OSCC depends on CHDH expression. Fatty acid desaturase 3 (FADS3) has been identified as one of the critical genes influencing the production and buildup of polyunsaturated fatty acids in mammals.³⁴ FADS3 is statistically associated with critical genes related to the tumorigenesis of OSCC.³⁵ It has been determined that CDGSHiron-sulfurdomain2 (CISD2) is an important gene in human aging, a mammalian lifespan-regulated regulator, and an oncoprotein that controls the development of cancer. CID2 mediates the emergence of diseases associated with human aging and the proliferation, differentiation, metastasis, and invasion of different cancer cells through several mechanisms.³⁶ For the first time, we discovered in this research that CDGSH was connected to the OSCC prognosis. In non-viral hepatocellular carcinoma, sideroflexin1 (SFXN1) expression reduced lipid buildup and ROS production, avoided the harmful consequences of fat overload, and forecasted the clinical prognosis of these patients.³⁷ When interferon-stimulated gene (ISG) expression is high in cancer cells, Ribonucleic acid exonuclease XRN1 acts as a negative regulator of the dsRNA sensor protein kinase R (PKR). When PKR is activated due to XRN1 deficiency, cancerous cells become lethal. Thus, for cancer cells that have an active interferon cellular state, XRN1 has been proposed as a possible therapeutic target.³⁸ One gene involved with glycogen degradation that is frequently up-regulated in malignancies is glycogen phosphorylase L (PYGL).³⁹ In several studies, PYGL is a prognostically relevant gene in OSCC.⁴⁰

The fast-development idea of cancer immunotherapy has been dubbed the “fifth pillar” of cancer treatment. However, the four other pillars of surgery- radiation, chemotherapy, and targeted therapy- continue to exist. It stimulates the immune system to control tumor growth and selectively targets tumor cells rather than healthy cells. These days, monoclonal antibodies, chimeric antigen receptor T cells (CART cells), checkpoint inhibitors, cytokines, radioimmunotherapy and cancer vaccines are among the therapeutic approaches that use the immune system.^41–43 Studies have found that high expression of CHDH will inhibit the infiltration of M1 macrophages and reduce the expression of CD86, indicating that CHDH plays an important role in tumor immune response.⁴⁴ In our study, we also found that when CHDH is overexpressed, activated T cells increase the possibility of destroying human OSCC cells. In the battle against cancer, the immune system is essential. It was investigated how immune cells differed in high- and low-risk groups and how immune cells and patient risk scores related to each other. Our research indicates that the MRGs risk score model has a substantial predictive value for immune response in OSCC, as evidenced by the significant correlation our model found between it and particular immune cells.

Our drug sensitivity analysis showed 9 drug sensitivities in the at-risk population. Among them, AZD1332 has been studied in various cancers.^45–47 but few studies have been done in OSCC. Dasatinib-induced LCK inhibition prevented the highly aggressive oral cancer cell line SAS cells from invading the body. LCK is a crucial regulator of oral cancer cell motility.⁴⁸ In conjunction with the drug sensitivity analysis of this study, for patients with oral cancer, dasatinib appears to be a viable anti-metastatic therapeutic choice. Docetaxel treats cancer through antimitotic activity, but its effects on the tumor microenvironment have recently received more attention.⁴⁹ In chronic lymphocytic leukemia, entospletinib is a B-cell receptor signaling kinase that targets B-cell receptor kinases,⁵⁰ and no studies have pointed to its therapeutic efficacy in OSCC. BPD-00008900 has been noted to be effective in hepatocellular carcinoma,⁵¹ but has not been explored in OSCC. Paclitaxel, a commonly used antitumor agent, is currently a first-line therapeutic choice for individuals with OSCC and is essential in the treatment of a variety of malignancies.⁵² Sepantronium bromide has been described as the first “survivin inhibitor.” Survivin proteins play a crucial role in OSCC, suggesting that Sepantronium bromide is of great value in OSCC.^53,54 The drug Luminespib is used in the treatment of ovarian cancer,⁵⁵ pancreatic cancer,⁵⁶ and non-small cell lung cancer,⁵⁷ and this study is the first to suggest that it is effective for OSCC treatment. In dopaminergic neurons, staurosporine triggers autophagy and mitophagy during cell death.⁵⁸ However, its potential as a treatment for OSCC is still largely unknown. In our investigation, we examined the medication’s sensitivity to offer additional viable choices for the pharmacological treatment of OSCC.

Our study still has several limitations, even if the MRG characterization model we built is effective in predicting the prognosis of OSCC patients. Since most of the OSCC samples used in our study were from public databases and were limited in number, some key information might be missing, resulting in inaccurate prediction results. For example, the constructed risk model had a 1-year AUC value of 0.599 in the test set. Therefore, more clinical data from the future are needed to demonstrate the effectiveness of our predictive risk model. Second, while collecting clinical tissue samples to support fundamental research, we only chose CHDH genes with independent predictive values for the initial validation of basic experiments, which may have missed other genes with fundamental values. To further validate our results, more research is needed. These elements could all present fresh chances for OSCC treatment.

Conclusion

The risk model constructed in this study can accurately predict the prognosis of OSCC patients and evaluate their immune status. CHDH is an important independent prognostic gene that is involved in the proliferation, migration, and invasion of OSCC. For physicians treating OSCC patients, our findings may open up new treatment options and new ideas for future research.