Mutation-associated transcripts reconstruct the prognostic features of tongue squamous cell carcinoma

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

Download PDF

Article

Mutation-associated transcripts reconstruct the prognostic features of tongue squamous cell carcinoma

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

This work is licensed under a CC BY 4.0 License

Journal Publication

published 03 Jan, 2023

Read the published version in International Journal of Oral Science →

You are reading this latest preprint version

Tongue squamous cell carcinoma is highly malignant with a poor prognosis. This study aimed to combine whole-genome sequencing, whole-genome methylation, and whole transcriptome analyses to better understand the molecular mechanisms of this cancer. Cancerous and paraneoplastic tissues from five patients with tongue squamous cell carcinoma were included as five paired samples. After multi-omics sequencing, differentially methylated intervals, methyl loop sites, methylated promoters, and transcripts were screened for variation in all paired samples. Correlations between them were analyzed to determine biological processes in tongue squamous cell carcinoma. We found five mutated methylation promoters that were significantly associated with the expression levels of mRNAs and lncRNAs. Functional annotation of these transcripts revealed their involvement in triggering the mitogen-activated protein kinase cascade, which is associated with cancer progression and the development of drug resistance during treatment. The prognostic signature models constructed based on the WDR81 and HNRNPH1 genes and combined clinical phenotype-gene prognostic signature models have shown high predictive efficacy and can be applied to predict patient prognostic risk in clinical settings. We identified biological processes in tongue squamous cell carcinoma that are initiated by mutations in the methylation promoter and are associated with the expression levels of specific mRNAs and lncRNAs. Ultimately, changes in the transcript levels affect the prognosis of tongue squamous cell carcinoma patients.

Tongue squamous cell carcinoma

Multi-omics

Methylation

Oral tongue squamous cell carcinoma (OTSCC) is the most common cancer in the oral cavity, characterized by insidious metastasis and high lymphatic metastasis. As a result, OTSCC has a higher risk and worse prognosis than other oral cancers (1). Initially, it was thought that OTSCC incidence was higher in the elderly population due to the accumulation of genetic mutations as well as risk factors such as long-term smoking, alcohol, and betel nut use (2–4). However, the lack of precancerous staging and effective early diagnostic markers for OTSCC has prevented establishing an efficient and accurate early warning system. The early warning system, which should be non-invasive or minimally invasive, can be used in high-risk groups so that a diagnosis can be made before the lesions are fully formed or after surgery. It is of great clinical significance to detect the cancer risk before the lesions are fully formed or before the metastases become established after surgery.

Recent reports have suggested an increase in OTSCC incidence among young people; no studies have explained this phenomenon (5). Currently, no effective diagnostic technology can meet the needs of early clinical diagnosis. The biggest problem is the lack of molecular diagnostic markers specific to OSCC. Therefore, it is necessary to identify key molecular markers with high sensitivity and specificity that can be monitored, screened, and diagnosed in a non-invasive or minimally invasive manner to accurately assess the disease status, improve the prognosis, and provide a better understanding of OSCC. Developing new targets is crucial for early diagnosis, precise drug use, accurate prognosis, and understanding OSCC pathogenesis. Therefore, searching for effective OSCC-specific molecular diagnostic markers and establishing rapid, sensitive, simple and non-invasive diagnostic tests have become the focus of research in OSCC prevention and treatment.

High-throughput sequencing technology has enabled researchers to obtain complete expression profiles, genome-wide data, and genome-wide methylation profiles. Whole transcriptome sequencing provides access to numerous differentially expressed genes and regulatory metabolic pathways, but genes do not always represent the entire molecular mechanism, and key signaling pathways are challenging to identify with too many differential genes. Therefore, transcriptome analysis often falls short of the intended research purpose (6).

DNA methylation is a common alteration at the molecular level and can be readily detected in various states of cell differentiation, especially among cancer cells (7). Analyses of DNA methylation have the potential to predict differences in survival and can help detect susceptibility to therapeutic approaches (8). In recent years, mRNA markers in serum or tissues as diagnostic or therapeutic targets for OSCC have been gaining attention because of their effectiveness, utility, and ability to identify mutations with high throughput screening. Circulating mRNA markers in serum and plasma have been extensively studied as tumor markers (9). Long-stranded non-coding RNAs (lncRNAs) play an important role in the development and prognosis of cancer, but their pathology has been poorly studied (10). Transcriptome-wide analysis has shown that 90% of human genetic DNA is transcribed into non-coding RNAs (ncRNAs) lacking protein-coding potential, while lncRNAs are ncRNAs ranging from 200 b to 100 kb in length. In malignant tumors, the abnormal expression of numerous lncRNAs is associated with cancer development, including lung, breast, and prostate cancer (11, 12). It has been suggested that abnormal lncRNA expression has been found in patients with OSCC and OSCC metastasis and that lncRNA detected in saliva can be used as a potential marker for OSCC. lncRNA in saliva has the potential to be a non-invasive rapid diagnostic marker for oral cancer (12). DNA methylation plays a vital role in normal mammalian development, but aberrant methylation has been associated with various differentiation-related diseases, including several human cancers. Early epigenetic alterations may contribute to the abnormalities in cellular genes that cause tumorigenesis. Thus, identifying methylation genes may provide a means of preventing and treating OSCC (13).

Therefore, in addition to analyzing transcriptomic data from tongue cancer tissues and paraneoplastic tissues, this study combined whole genome sequencing (WGS) and whole genome bisulfite sequencing (WGBS) to investigate variations in expression profiles. This joint analysis allows for greater precision in targeting key regulatory genes associated with tongue cancer development and progression. We also explored the upstream and downstream regulatory relationships of key genes.

2.1 Sample Information

Five patients with HE-stained tissue sections including sample information are shown in supplementary material figure S1. All the 5 samples were squamous cell carcinoma of the floor of the mouth or tongue, and 3 of them had lymph node metastasis. The details of patient’s gender, age, weight, smoking history, alcohol consumption history, and tumor grade are shown in Table S1.

2.2 Results of the differentially methylated site (DMS) and differentially methylated promoter (DMP) screening

The analysis of different cytosine methylation sites (CG, CHG, CHH, and C) identified 291 mC-, 2262 mCG-, 1 mCHH-, and 0 mCHG-types among the five groups of samples. The CHH-type DMS was not located in the coding or promoter region, while some of the C- and CG-type DMSs were. Statistical analysis identified 82 shared C-type DMSs in the coding region and nine in the promoter region, as well as 1,138 shared CG-type DMSs in the coding region and 160 in the promoter region (Table 1). Significant differential methylation sites shared by the 5 groups of cancer and paracancer samples are shown in Table S2. Significantly different methylation sites appearing in the coding regions of the genes common to the 5 groups of cancer and paraneoplastic samples are shown in Table S3. Significantly differentially methylated sites in the promoter region common to the 5 cancer groups and paraneoplastic samples are shown in Table S4.

The DMP analysis for different types of cytosine methylation modification sites showed a shared significant difference in methylation sites in all five cancer and paracancerous groups, yielding a total of 5837 mC-, 1804 mCG-, 633 mCHG-, and 5872 mCHH-types (Table 1). Detailed information regarding all DMPs is available in the Supplemental File (Table S5).

2.3 Results of the differential transcriptional analysis

The results of differential expression analysis of tumor and paracancerous transcripts showed that significant differences in tumor tissue expression profiles compared with paracancerous tissue (Fig. 2a). 1213 mRNAs were significantly upregulated, and 1768 mRNAs were significantly downregulated (Fig. 2b). 93 lncRNAs were significantly upregulated and 259 lncRNAs were significantly downregulated (Fig. 2c); 128 miRNAs were significantly upregulated, and 117 miRNAs were significantly downregulated; (Fig. 2d).

2.4 Methylation promoter and transcriptome correlation

Among the different types of cytosine methylation sites with common DMPs in the five groups of samples, corresponding mRNA transcripts were also significantly different in all the groups; there were five CG-type, seven CHG-type, 40 CHH-type, and 39 C-type DMPs. The results showed that SDR9C7 and MAPK8IP2 expression was significantly correlated with the corresponding C-type DMPs; HAND2 and SEPP1 expression levels were significantly correlated with the corresponding CG-type DMPs; and GALNT2 expression was significantly correlated with the corresponding CHH-type DMPs (Fig. 3a). No genes had expression significantly correlated with the corresponding CHG-type DMPs. Based on the human transcription factor target gene data included in the TRRUST database, HAND2 was identified as a transcription factor with four well-defined target genes (DBH, GATA4, NPPA, and PHOX2A) and was positively correlated with all of them except GATA4.

Among the different types of cytosine methylation sites with common DMPs in the five sample groups, lncRNA transcripts were also significantly different in all groups; there were 0 CG-type, three CHG-type, eight CHH-type, and six C-type DMPs. The result showed that LINC00885 expression was significantly correlated with the corresponding C-type DMPs (Fig. 3a). LINC00885 had 45 target genes and was significantly correlated in expression with all of them (Fig. 3b). Table 2 provides correlations of these six C-types DMPs with gene transcripts and their corresponding promoters.

2.5 Functional enrichment analysis of shared DMP-related transcripts and differential mRNA

GO analysis of mRNA transcripts significantly associated with shared DMPs (Fig. 4a) revealed six pathways containing more than two genes: positive regulation of stress-activated mitogen-activated kinase (MAPK) cascade, positive regulation of stress-activated protein kinase signaling cascade, regulation of stress-activated MAPK cascade, regulation of stress-activated protein kinase signaling cascade, stress-activated MAPK cascade, and stress-activated protein kinase signaling cascade. KEGG analysis resulted in only two signaling pathways; target genes of lncRNA transcripts correlated with shared DMPs were significantly enriched in arginine and proline metabolism (Fig. 4b).

The results of mRNA transcriptome analysis identified 2981 differentially expressed genes between the five paraneoplastic tissue samples and the five cancer tissue samples. According to GO analysis, the main functions of these genes were the formation of extracellular matrix tissues and structural tissues, promotion of myoblast development and formation, as well as the mediation of myofiber movement (Fig. 4c). KEGG analysis revealed the interaction between tumor cells and the extracellular matrix, which constitute the tumor metastasis channel, the formation of adhesion spots, and the formation of proteoglycans in the extracellular matrix of tumor cells (Fig. 4d).

2.6 TCGA cohort validation

Among the LINC00885 target genes, the expression levels of HNRNPH1, SEMA6D, and NTRK2 were significantly associated with prognosis (Fig. 5a–c). For TCGA-HNSCC cohort data, we selected eight genes with the largest significant differences in the distribution of differentially expressed genes between paracancerous and cancerous tissues, as well as in Kaplan–Meier survival curve results. Among them, TMPRSS11B and GAS2 were significantly overexpressed in paracancerous tissues compared with cancerous tissues. Accordingly, patients with high TMPRSS11B and GAS2 expression had a better prognosis than those with low expression. MMP11, TMBIM6, NOMO2, LAMC2, HMGA2, and CSF2 were significantly downregulated in paracancerous tissues compared with cancerous tissues. Accordingly, patients with low expression of these genes had a better prognosis than patients with high expression (Fig. 5d, e).

In the TCGA-HNSCC cohort, the difference in LINC00885 expression between cancer and paraneoplastic tissues was statistically significant (Fig. 6a). Twelve of LINC00885's target genes were upregulated, and 16 were downregulated in cancer tissues (Fig. 6b). The upregulated genes were PLK4, RHOC, CHEK1, P4HA2, LEMD1, TMEM244, MICAL2, SLC15A3, FOXI3, PNCK, and NELL2, and the down-regulated genes were FAM149A, BOLL, RBPMS, ALDH3A2, GALNT1, B9D1, ZNF85 SMAD3, NEDD4L, WDR81, PGAP1, PRPF40A, PALM2, KIAA1429, CLPB, PCYT1A, and HNRNPH1 (Fig. 6c). The prognostic signature was developed based on 45 target genes, using the least absolute shrinkage and selection operator Cox (LASSO-Cox) analysis (Fig. 6d, e). Furthermore, 2 survival-associated target genes, WDR81 and HNRNPH1, were selected in the final prognostic signature, and the coefficients were obtained from the LASSO algorithm. The signature calculated the risk score for each patient using the function predict, and each patient was grouped into a high or low-risk group according to the median risk score. The Kaplan-Meier survival analysis showed that the HNSCC patients in the high-risk group had a significantly shorter OS than HNSCC patients in the low-risk group (Fig. 6f). The area under the ROC curve (AUC) of the prognostic signature model equal to 0.75 indicated an acceptable prediction efficiency (Fig. 6g). Single-factor random forest plots showed that patient age, pathological grade, pathological TNM stage, and RS reduced survival time and promoted adverse prognostic events (Fig. 6h). Multi-factor random forest plots showed that the model's RS significantly reduced survival time and promoted adverse prognostic events after association with the patient's clinical phenotype (Fig. 6i). The combined clinical phenotype-genetic prognostic risk model divided patients into high-risk and low-risk groups based on the median RS. Table 3 shows the risk coefficient corresponding to each clinical phenotype and gene. Kaplan-Meier survival curves showed a significant difference in survival time between the high-risk and low-risk groups, with a significantly higher overall 5-year survival rate for patients in the low-risk group than in the high-risk group (Fig. 6j). The AUC of the prognostic signature model equal to 0.817 indicated a higher prediction efficiency than the prognostic risk model built on genes alone (Fig. 6k). The final nomogram was constructed based on the factors included in the gene-clinical phenotype-based prognostic characteristics model (Fig. 7a). The closer the red line matches the black diagonal line, the closer the predicted result is to the actual situation (Fig. 7b–d).

2.7 Summary of CCLE lineage analysis

WDR81 expression in head and neck squamous cell carcinoma tissues was lower than in many other diseases. In addition, WDR81 expression was high in skin cancer and myeloma and, conversely, low in diseases such as cervical cancer and teratoma (Fig. 8a). HNRNPH1 expression in head and neck squamous cell carcinoma tissues was significantly lower than in various other diseases but higher in gallbladder cancer and teratoma. Furthermore, HNRNPH1 was significantly highly expressed in leukemia and embryonal cancer tissues (Fig. 8b). Our sequencing results for WDR81 and HNRNPH1 are consistent with the mRNA expression interval of the CCLE pan-cancer spectrum results.

The expression levels of WDR81 and HNRNPH1 differed significantly in various cell lines of head and neck squamous cell carcinoma. We selected 20 cell lines with the highest and lowest expression for each gene separately for demonstration (Fig. 8c, d).

Studies on tongue cancer are severely lacking; with no clear marker genes, prognosis and progression are difficult to predict, and low drug sensitivity during treatment is a challenge. Prior research addressing the role of single transcripts in the treatment of OTSCC has failed to elucidate how their effects at the genetic level. Reviews examining the relationship between clinical phenotype and tumor have also been unsuccessful at answering this question (14–17). Although the incidence of this cancer is increasing each year among young people (1–4), the reason for the increase remains unclear. No credible association between smoking and OTSCC has been found among young or old patients (18, 19).

Here, we demonstrated the MAPK cascade involvement of all target genes corresponding to mutant methylation promoters simultaneously present in the tumors of five patients with OTSCC. Thus far, no other studies have found abnormalities in the expression of MAPK cascade pathways in OTSCC tissues. The MAPK pathway mediates cell proliferation, differentiation, and chemotaxis. Its negative feedback regulation of the MAPK cascade in cancer cells reduces the sensitivity and efficacy of cancer therapeutic agents (20, 21). Moreover, the pathway influences important physiological processes (e.g., neuronal function, immune response, and embryonic development) through regulating gene expression, cytoskeletal protein dynamics, and cell proliferation or apoptosis pathways (22, 23). Based on our findings and previous reports, we put forward a preliminary hypothesis that the MAPK cascade is deeply involved in the biological variation of oral cancer development, migration, and drug resistance. The results of this study suggest that mutation of methylated promoters triggers aberrant expression of mRNA transcripts, ultimately activating the MAPK cascade. Therefore, if we can target and block specific methylated promoter mutations and the resulting MAPK cascade, we may be able to reduce the likelihood of adverse events in OTSCC.

The occurrence of extracellular matrix heterogeneity is inextricably linked to tumors. The precipitation and mechanical sclerosis of the extracellular matrix are considered key factors leading to tumor infiltration and metastasis (24–27). Here, we confirmed that in tumor samples of patients with OTSCC, an abnormal extracellular matrix is associated with cancer progression. Recent research suggests that activation of the PI3K-Akt signaling pathway promotes epithelial-mesenchymal transition (EMT), ultimately resulting in tumor invasion, metastasis, and drug resistance (28). Studies have also confirmed that some specific inhibitors of the PI3K-Akt signaling pathway (e.g., marine drugs) can reverse EMT and thus reduce the degree of drug resistance in tumor tissue during treatment (29, 30). Similar to these previous results, we also noted activation of the PI3K-Akt signaling pathway in OTSCC. Thus, further investigation is needed to determine whether drugs that specifically act on this pathway can inhibit EMT and improve patient prognosis.

Human papillomavirus (HPV)-negative tumors are believed to be associated with the development of oral squamous cell carcinoma, including of the tongue, and are also predictive of poor prognosis and treatment resistance (31, 32). In our study, HPV infection was the KEGG pathway with the second-highest number of aberrantly expressed transcripts enriched in OTSCC tissue. Another pathway of importance was the calcium signaling pathway, which is involved in crosstalk with reactive oxygen species (ROS) signaling pathways leading to tumorigenesis (33, 34). Inhibition of calcium signaling can inhibit cancer cell proliferation and metastasis in some cancers (35). Consistent with previous findings, our study confirmed the activation of the calcium signaling pathway in OTSCC tissues.

Data from the TCGA-HNSCC cohort were used to verify whether our results aligned with the results of previous large sample analyses. Notably, the most differentially expressed transcripts in OTSCC were clearly correlated with prognosis in the TCGA-HNSCC cohort. Nevertheless, the potential heterogeneity between HNSCC and OTSCC, ethnic differences in the patients studied, and differences in sequencing methods can explain why some of our results are inconsistent with the results of the HNSCC cohort analysis. With a larger sample size in future research, the credibility of our findings will also increase.

It has been demonstrated that LINC00885 promotes tumor cell proliferation and invasion (36). Current research has focused on breast and cervical cancers, but its expression in OTSCC remains unclear (37, 38). Our findings revealed that LINC00885 was also upregulated in OTSCC tissues. The prognostic signature model based on the target gene of LINC00885 has high predictive efficacy in predicting the patients' prognosis after surgery or treatment. Notably, since both HNRNPH1 and WDR81 have risk coefficients less than 0 in the risk profile model, the patients’ risk scores are the absolute value of the actual risk score. In this study, the HNRNPH1 and WDR81 genes were highly expressed in OTSCC tumor samples and serve as protective genes for predicting prognosis. Compared to gene-based models, the combined clinical phenotype-gene model has a more reliable predictive efficacy, but more complex information needs to be collected.

Our study found no correlation between mutations occurring in the genome and differences in the transcriptome. However, this does not mean that genes are not mutated in OTSCC tissues or that gene mutations do not affect the transcriptome, thereby leading to functional changes. On the contrary, genomic sequencing of OTSCC tissue samples from five patients showed a large number of mutations in the genome of the tumor tissue compared to the paracancerous tissue samples.

Our study shows that two mC-, two mCG-, and one mCHH-type methylation mutations cause aberrant expression of the transcriptome in OTSCC. Mapping such molecular changes to cellular functions reveals differences in MAPK cascade pathways. Further, lncRNAs as well as their target genes in the variants were used to predict the prognostic risk of patients. Ultimately, such changes lead to cancer development, increased drug resistance, and suboptimal prognosis in patients with OTSCC. Analysis of transcripts showed that five patients with OTSCC had differential genes mainly clustered in pathways with multiple functions. These pathways include deposition and mechanical sclerosis of extracellular matrix tissue, PI3K-Akt signaling pathway leading to EMT, HPV infection, and interaction of calcium signaling with ROS signaling. In conclusion, this study can provide a theoretical basis for follow-up research on experimental etiology or interventions. Targeted blockade of specific methylated promoter mutations and the resulting MAPK cascade may be a new direction for reducing adverse events in OTSCC. The prognostic signature models constructed based on the WDR81 and HNRNPH1 genes and the combined clinical phenotype-gene prognostic signature models have shown high predictive efficacy and can be constructed to predict patient prognostic risk in the clinical setting.

4.1 Sample and data collection

Samples were surgically excised from tongue cancer tissues and paracancerous tissues. During tissue processing, we washed off bloodstains with saline while removing non-essential tissues. Samples were dried with gauze, then cut into tissue blocks less than 0.5 cm thick and placed into labeled RNase-free cryotubes or EP tubes. The tubes were snap-frozen in liquid nitrogen and stored at − 80°C. The study included 10 samples from five patients with tongue cancer, divided into tumor and control groups based on histological characteristics. Paired samples (tumor and control) were obtained for each patient. Next, we sequenced the whole transcriptome via RNA sequencing (RNA-seq), WGS, and WGBS for each of the 10 samples.

4.2 Differential methylation sites and promoter methylation screening

Based on the WGS results, we screened for mutations co-existing in five sets of paired samples. We next screened and genetically annotated the differentially methylated regions (DMRs) present in all paired samples. We then screened for differentially methylated sites (DMSs) and differentially methylated promoters (DMPs) common to the five sets of paired samples based on the WGBS results. We defined cytosine loci with P < 0. 05 as statistically significant. Loci with statistically significant differences between tumor samples and comparison samples were considered DMSs. Similarly, promoters with P < 0.05 were selected as DMPs.

4.3 Differential expression of the transcriptome

We further analyzed and collated whole transcriptome sequencing data, including long non-coding RNA (lncRNA), mRNA, and micro-RNA (miRNA) sequences. Differential expression analysis between groups was performed with the tumor and paracancerous tissues separately using DESeq2, a software used to detect differentially expressed genes with duplicate samples (39). The screening was conditioned on differential ploidy of ≥ 2 and P < 0. 05.

4.4 Correlation between DMPs and the transcriptome and gene function prediction

Based on the DMPs shared by the five sets of paired samples, we further analyzed the effect of their modifying effects on the transcripts. Transcript data were obtained through differential genes analysis. During the analysis, we extracted mRNA and lncRNA transcripts with significantly different expression in five groups of samples. Their expression was subjected to Pearson correlation analysis with their corresponding shared DMPs. The relationships between DMPs and transcripts were considered significant at P < 0.05, with the absolute value of correlation coefficients greater than 0.9. Screened transcripts were significantly correlated with their corresponding shared DMPs, indicating that DMPs regulated the expression of these transcripts in tumor tissues.

We used the R packages “org.Hs.eg.db” (40) and “clusterProfiler” (41) to perform Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses of transcripts associated with shared DMPs to determine their functional pathways. Similarly, we performed GO and KEGG enrichment analyses on differentially expressed transcripts co-occurring in the five sets of paired samples to understand how biological functions vary across cancerous and paraneoplastic tissues.

4.5 The Cancer Genome Atlas (TCGA) cohort validation

We downloaded the head and neck squamous cell carcinoma (HNSCC) cohort mRNA transcriptome data, as well as follow-up data, using the TCGA database. We used Perl software to initially organize the data. We organized the data using the R package “survival” (42, 43) and performed survival analysis with “survminer” (44). We verified whether the expression of lncRNA target genes associated with shared DMPs has an impact on the prognosis of patients with HNSCC. The Least Absolute Shrinkage and Selection Operator (LASSO) can retain the most representative variables, which were valuable for the model's accuracy. Consequently, it is considered by researchers as an effective high-dimensional predictive regression method to avoid over-fitting of the model variables. [16]. Multivariate Cox analysis was performed to determine whether the selected gene is a prognostic factor in patients with HNSCC. We used the Lasso-Cox regression model to show the ideal risk coefficient of each prognostic feature for the genes in the HNSCC prognostic signature model. In addition, the risk score (RS) was calculated by the function PREDICT. Patients were defined as high risk if their RS was above the median and low risk if their RS was below the median. Risk coefficients were used to distinguish whether factors were prognostic protective or risk factors, and the ability of each factor to affect prognosis. A factor greater than 0 was considered a risk factor, and a factor less than 0 was a protective factor. A larger absolute value indicated that the factor had a more significant impact on prognosis. The next step was to observe whether the prognosis between the two groups was different over time. The Kaplan-Meier survival analysis with a two-sided log-rank test was performed to assess the difference in the prognosis between the two groups. Next, the Receiver Operating Characteristic (ROC) curve was used to identify the accuracy of the model’s prediction. Univariate and multivariate Cox regression analyses were performed to identify the independent prognostic factors for the HNSCC cohort. The survival difference between the high-risk and low-risk groups was stratified based on age, gender, histologic grade, tumor stage, and pathological T/N/M stage. Finally, the prognostic risk of a combined gene-clinical phenotype was modeled in a similar manner.

The prognostic impact of transcripts was verified using box plots showing differential gene expression in cancer and paracancer tissues and also using Kaplan–Meier survival curves. Ultimately, the nomogram was constructed based on the factors included in the model. The calibration curve was used to determine how well the predicted results match the formal situation.

4.6 Cancer Cell Line Encyclopedia (CCLE) shows pan-cancer lineage and cell lineage

We download the expression of genes in the prognostic risk model across the cancer spectrum from the CCLE database. We visualized the results so that the expression levels of these genes in each cancer or in different cells of head and neck squamous cell carcinoma are clearly presented.

Additionally, while some studies have investigated the molecular mechanism of this cancer, few have combined WGS, WGBS, and transcriptome sequencing. The lack of such a joint analysis prompted us to conduct this study. Given the complexity of our analyses and the results obtained, we have created a flow chart for easy interpretation (Fig. 1).

4.7 Statistical analysis and data processing

For all statistical analyses in this study, a P value less than 0.05 was defined as statistically significant. All statistical analyses were performed by the software Perl and R. The graphs were plotted based on the packages "ggpubr" and "ggplot2" (45, 46). Analysis of variance was performed with the program package "edgeR" (47). The package "survival" was used for the integration of survival times and ending events (48), the package "glmnet" was used for the Lasso regression analysis (49, 50), and the package "survminer" was used to plot Kaplan-Meier survival curves (51). The final patient RS was calculated by the function predict, while the function coef algorithm calculated the coefficients of each factor in the prognostic risk profile model. The plotting of ROC curves and the calculation of AUC values were implemented by the package "survivalROC" (52). The plotting of the random forest plot was performed with the package "forestplot," and the function forestplot was used to plot (53). Finally, the program package "rms" was used to construct the nomogram. Since quantitative numerical data are required in the random forest analysis, we transformed the TNM staging according to the staging values and kept only numerical data. For gender, we defined female as 0 and male as 1.

Acknowledgements

This work was supported by funding from the National Nature Science Foundations of China (grant numbers 81772275 and 32071462).

Ethics approval statement

The study received approval from the Ethics Committee of West China Hospital of Stomatology, Sichuan University.

Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.

Bello IO, Soini Y, Salo T. Prognostic evaluation of oral tongue cancer: means, markers and perspectives (I). Oral Oncol 2010;46:630–635.
Khowal S, Naqvi SH, Monga S, Jain SK, Wajid S. Assessment of cellular and serum proteome from tongue squamous cell carcinoma patient lacking addictive proclivities for tobacco, betel nut, and alcohol: Case study. J Cell Biochem 2018;119:5186–5221.
Sankaranarayanan R, Masuyer E, Swaminathan R, Ferlay J, Whelan S. Head and neck cancer: a global perspective on epidemiology and prognosis. Anticancer Res 1998;18:4779–4786.
VanMonsjou HSAW, V.B.%A VanDenBrekel,M.W.M.%A Balm,A.J.M. Head and neck squamous cell carcinoma in young patients. Oral oncology%J 2013;49:1097–1102.
Paderno A, Morello R, Piazza C. Tongue carcinoma in young adults: a review of the literature. Acta Otorhinolaryngol Ital 2018;38:175–180.
Reuter JAAS, Damek V.%A Snyder, Michael P. High-Throughput Sequencing Technologies. Molecular cell%J 2015;58:586–597.
PA SSAKTAJ. Epigenetics in cancer. Carcinogenesis%J 2010;31:27–36.
Jones PA, Baylin SB. The fundamental role of epigenetic events in cancer. Nat Rev Genet 2002;3:415–428.
Li Y, Elashoff D, Oh M, Sinha U, St John MA, Zhou X, et al. Serum circulating human mRNA profiling and its utility for oral cancer detection. J Clin Oncol 2006;24:1754–1760.
Zhang H, Zhao L, Wang YX, Xi M, Liu SL, Luo LL. Long non-coding RNA HOTTIP is correlated with progression and prognosis in tongue squamous cell carcinoma. Tumor Biology 2015;36:8805–8809.
Kong XP, Yao J, Luo W, Feng FK, Ma JT, Ren YP, et al. The expression and functional role of a FOXC1 related mRNA-lncRNA pair in oral squamous cell carcinoma. Mol Cell Biochem 2014;394:177–186.
Tang HK, Wu ZY, Zhang JP, Su B. Salivary lncRNA as a potential marker for oral squamous cell carcinoma diagnosis. Molecular Medicine Reports 2013;7:761–766.
Li YF, Hsiao YH, Lai YH, Chen YC, Chen YJ, Chou JL, et al. DNA methylation profiles and biomarkers of oral squamous cell carcinoma. Epigenetics-Us 2015;10:229–236.
Jia LF, Wei SB, Gan YH, Guo Y, Gong K, Mitchelson K, et al. expression, regulation and roles of miR-26a and MEG3 in tongue squamous cell carcinoma. International Journal of Cancer 2014;135.
Jiang L, Yang D, Liu X, Wang C, Wang A, Chen Z, et al. Identification and experimental validation of G protein alpha inhibiting activity polypeptide 2 (GNAI2) as a microRNA-138 target in tongue squamous cell carcinoma. Human Genetics 2011;129:189–197.
Zhang K, Wang Q, Zhong B, Gong Z. LUCAT1 as an oncogene in tongue squamous cell carcinoma by targeting miR-375 expression. Journal of Cellular and Molecular Medicine 2021.
Poonam R. Clinicopathological correlation of tumor-stroma ratio and inflammatory cell infiltrate with tumor grade and lymph node metastasis in squamous cell carcinoma of buccal mucosa and tongue in 41 cases with review of literature. Journal of cancer research and therapeutics;16:445–451.
Khadijah M. Meta-analysis on risk factors of squamous cell carcinoma of the tongue in young adults. Journal of oral and maxillofacial pathology: JOMFP;23:450–457.
Maura L G. Distinct risk factor profiles for human papillomavirus type 16-positive and human papillomavirus type 16-negative head and neck cancers. Journal of the National Cancer Institute 2008;100:407–420.
M M M. Integrating signals from RTKs to ERK/MAPK. Oncogene 2007;26:3113–3121.
David L. Negative feedback regulation of the ERK1/2 MAPK pathway. Cellular and molecular life sciences: CMLS 2016;73:4397–4413.
R J D. Signal transduction by the JNK group of MAP kinases. Cell 2000;103:239–252.
András Z. JNK Signaling: Regulation and Functions Based on Complex Protein-Protein Partnerships. Microbiology and molecular biology reviews: MMBR 2016;80:793–835.
Michael W P. The extracellular matrix modulates the hallmarks of cancer. EMBO reports 2014;15:1243–1253.
FuiBoon K. The Extracellular Matrix Modulates the Metastatic Journey. Developmental cell 2019;49:332–346.
Paolo P P. Collagen density promotes mammary tumor initiation and progression. BMC medicine 2008;6:11.
Caroline B. Remodelling the extracellular matrix in development and disease. Nature reviews Molecular cell biology 2014;15:786–801.
Wenting X. A new role for the PI3K/Akt signaling pathway in the epithelial-mesenchymal transition. Cell adhesion & migration 2015;9:317–324.
Junli D. Inhibition of PI3K/Akt/mTOR signaling pathway alleviates ovarian cancer chemoresistance through reversing epithelial-mesenchymal transition and decreasing cancer stem cell marker expression. BMC cancer 2019;19:618.
Jiaen W. Marine compounds targeting the PI3K/Akt signaling pathway in cancer therapy. Biomedicine & pharmacotherapy = Biomedecine & pharmacotherapie 2020;129:110484.
K Kian A. Human papillomavirus and survival of patients with oropharyngeal cancer. The New England journal of medicine 2010;363:24–35.
Ashley H. Targeted Therapy in Oropharyngeal Squamous Cell Carcinoma: The Implications of HPV for Therapy. Rare cancers and therapy;3:89–117.
Nadine H. Crosstalk between calcium and reactive oxygen species signaling in cancer. Cell calcium 2017;63:70–96.
Sarah J R-T. The Calcium-Signaling Toolkit in Cancer: Remodeling and Targeting. Cold Spring Harbor perspectives in biology 2019;11.
Peng C. Erianin, a novel dibenzyl compound in Dendrobium extract, inhibits lung cancer cell growth and migration via calcium/calmodulin-dependent ferroptosis. Signal transduction and targeted therapy 2020;5:51.
Liu YW, Tu HY, Zhang LL, Xiong JP, Li L. FOXP3-induced LINC00885 promotes the proliferation and invasion of cervical cancer cells. Molecular Medicine Reports 2021;23.
Liu YL, Chen JR, Zhou LZ, Yin CH. LINC00885 promotes cervical cancer progression through sponging miR-3150b-3p and upregulating BAZ2A. Biology Direct 2022;17.
Abba MC, Canzoneri R, Gurruchaga A, Lee J, Tatineni P, Kil H, et al. LINC00885 a Novel Oncogenic Long Non-Coding RNA Associated with Early Stage Breast Cancer Progression. International Journal of Molecular Sciences 2020;21.
Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol 2014;15:550.
M. C. org.Hs.eg.db: Genome wide annotation for Human. 2015.
Yu G, Wang LG, Han Y, He QY. clusterProfiler: an R package for comparing biological themes among gene clusters. Omics 2012;16:284–287.
Lin H, Zelterman D. Modeling Survival Data: Extending the Cox Model. Technometrics 2000;44:85–86.
Therneau T. A package for survival analysis. R package 2.37-2. 2012.
Kassambara A. Drawing Survival Curves using 'ggplot2' [R package survminer version 0.2.0]. 2017.
Ginestet C. ggplot2: Elegant Graphics for Data Analysis. J R Stat Soc a Stat 2011;174:245–245.
ggpubr: 'ggplot2' Based Publication Ready Plots. 2017.
Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 2010;26:139–140.
Li JCA. Modeling survival data: Extending the Cox model. Sociol Method Res 2003;32:117–120.
Friedman JH, Hastie T, Tibshirani R. Regularization Paths for Generalized Linear Models via Coordinate Descent. Journal of Statistical Software 2010;33.
Simon N, Friedman J, Hastie T, Tibshirani R. Regularization Paths for Cox's Proportional Hazards Model via Coordinate Descent. Journal of Statistical Software 2011;39:1–13.
Kassambara A, Kosinski M, Biecek P. survminer: Drawing Survival Curves using 'ggplot2'. 2016.
Heagerty PJ. survivalROC: Time-dependent ROC curve estimation from censored survival data. 2013.
Gordon M, Lumley T. forestplot: Advanced Forest Plot Using 'grid' Graphics. 2016.

Tables 1 to 3 are not available with this version.

(Not answered)

Download PDF

Journal Publication

published 03 Jan, 2023

Read the published version in International Journal of Oral Science →

Editorial decision: revise
09 Sep, 2022
Review #1 received at journal
07 Sep, 2022
Review #2 received at journal
21 Aug, 2022
Reviewer #2 agreed at journal
16 Aug, 2022
Reviewer #1 agreed at journal
16 Aug, 2022
Reviewers invited by journal
16 Aug, 2022
Submission checks completed at journal
16 Aug, 2022
Editor assigned by journal
14 Aug, 2022
First submitted to journal
14 Aug, 2022

You are reading this latest preprint version

Mutation-associated transcripts reconstruct the prognostic features of tongue squamous cell carcinoma

Status:

Journal Publication

Version 1

Abstract

Figures

1. Introduction

2. Results

2.1 Sample Information

2.2 Results of the differentially methylated site (DMS) and differentially methylated promoter (DMP) screening

2.3 Results of the differential transcriptional analysis

2.4 Methylation promoter and transcriptome correlation

2.5 Functional enrichment analysis of shared DMP-related transcripts and differential mRNA

2.6 TCGA cohort validation

2.7 Summary of CCLE lineage analysis

3. Discussion

4. Materials And Methods

4.1 Sample and data collection

4.2 Differential methylation sites and promoter methylation screening

4.3 Differential expression of the transcriptome

4.4 Correlation between DMPs and the transcriptome and gene function prediction

4.5 The Cancer Genome Atlas (TCGA) cohort validation

4.6 Cancer Cell Line Encyclopedia (CCLE) shows pan-cancer lineage and cell lineage

4.7 Statistical analysis and data processing

Declarations

References

tables

Additional Declarations

Supplementary Files

Status:

Journal Publication

Version 1