Genome-wide DNA methylation profiling in tongue squamous cell carcinoma

Oral Diseases. 2022;00:1–13. | 1wileyonlinelibrary.com/journal/odi Received: 23 June 2022 | Revised: 14 October 2022 | Accepted: 5 November 2022 DOI: 10.1111/odi.14444 ORIGINAL ARTICLE Genomewide DNA methylation profiling in tongue squamous cell carcinoma Óscar RapadoGonzález1,2,3,4 | Nicolás CostaFraga2,5,6 | Aida BaoCaamano2,5,6 | José Luis LópezCedrún7 | Roberto ÁlvarezRodríguez8 | Ana Belén Crujeiras9,10 | Laura MuineloRomay2,3,4 | Rafael LópezLópez2,4,11 | Ángel DíazLagares2,4,5 | María Mercedes SuárezCunqueiro1,2,4,11 1Department of Surgery and MedicalSurgical Specialties, Medicine and Dentistry School, Universidade de Santiago de Compostela (USC), Santiago de Compostela, Spain 2Galician Precision Oncology Research Group (ONCOGAL), Medicine and Dentistry School, Universidade de Santiago de Compostela (USC), Santiago de Compostela, Spain 3Liquid Biopsy Analysis Unit, Translational Medical Oncology Group (ONCOMET), Health Research Institute of Santiago (IDIS), Santiago de Compostela, Spain 4Centro de Investigación Biomédica en Red en Cáncer (CIBERONC), Instituto de Salud Carlos III (ISCIII), Madrid, Spain 5Epigenomics Unit, Cancer Epigenomics, Translational Medical Oncology Group (ONCOMET), Health Research Institute of Santiago (IDIS), University Clinical Hospital of Santiago (CHUS, SERGAS), Santiago de Compostela, Spain 6Universidade de Santiago de Compostela (USC), Santiago de Compostela, Spain 7Department of Oral and Maxillofacial Surgery, Complexo Hospitalario Universitario de A Coruña (CHUAC, SERGAS), A Coruña, Spain 8Department of Pathology, Complexo Hospitalario Universitario de A Coruña (CHUAC, SERGAS), A Coruña, Spain 9Epigenomics in Endocrinology and Nutrition Group, Epigenomics Unit, Health Research Institute of Santiago de Compostela (IDIS), University Clinical Hospital of Santiago (CHUS, SERGAS), Santiago de Compostela, Spain 10Centro de Investigación Biomédica en Red Fisiopatología de la Obesidad y Nutrición (CIBERobn), Instituto de Salud Carlos III (ISCIII), Madrid, Spain 11Translational Medical Oncology Group (ONCOMET), Health Research Institute of Santiago (IDIS), Complexo Hospitalario Universitario de Santiago de Compostela (CHUS, SERGAS), Santiago de Compostela, Spain This is an open access article under the terms of the Creative Commons Attribution-NonCommercial License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited and is not used for commercial purposes. © 2022 The Authors. Oral Diseases published by Wiley Periodicals LLC. Correspondence Ángel DíazLagares, Epigenomics Unit, Cancer Epigenomics, Translational Medical Oncology Group (ONCOMET), Research Institute of Santiago de Compostela (IDIS), Complexo Hospitalario Universitario de Santiago de Compostela (CHUS), C/ Choupana S/N 15706, Santiago de Compostela, Spain. Email: [email protected] María Mercedes SuárezCunqueiro, Translational Medical Oncology Group (ONCOMET), Research Institute of Santiago de Compostela (IDIS), Complexo Hospitalario Universitario de Santiago de Compostela (CHUS), Santiago de Compostela, Spain. Email: [email protected] Funding information Centro de Investigación Biomédica en Red de Cáncer, Grant/Award Number: Abstract Objectives: To provide a comprehensive characterization of DNA methylome of oral tongue squamous cell carcinoma (OTSCC) and identify novel tumorspecific DNA methylation markers for early detection using saliva. Material and Methods: Genomewide DNA methylation analysis including six OTSCC matched adjacent nontumoral tissue and saliva was performed using Infinium MethylationEPIC array. Differentially methylated levels of selected genes in our OTSCC cohort were further validated using OTSCC methylation data from The Cancer Genome Atlas database (TCGA). The methylation levels of a set of tumorspecific hypermethylated genes associated with a downregulated expression were evaluated in saliva. Receiver operating characteristic (ROC) curves were performed to assess the diagnostic value of DNA methylation markers. Results: A total of 25,890 CpGs (20,505 hypomethylated and 5385 hypermethylated) were differentially methylated (DMCpGs) between OTSCC and adjacent nontumoral 2 | RAPADO-GONZÁLEZ et al. 1 | INTRODUCTION Oral tongue squamous cell carcinoma (OTSCC) is the most common subtype of oral cavity cancer, and it is considered one of the most aggressive tumors due to its high potential for local invasiveness and great propensity to metastasize cervical lymph nodes (Ocharoenrat et al., 2003). Consequently, OTSCC involves a poor prognosis even when it is diagnosed at early stages, indicating a different molecular behavior compared to other head and neck cancers (Rusthoven et al., 2008). DNA methylation is a crucial epigenetic modification of the genome which is involved in regulating many cellular processes including embryonic development, transcription, chromatin structure, X chromosome inactivation, genomic imprinting, and chromosome stability (Robertson, 2005). Abnormal distribution of DNA methylation is one of the beststudied epigenetic changes in cancer, where global hypomethylation leading to activation of oncogenes and transposons is often accompanied by focal hypermethylation of the CpG islands (CpGIs) in the promoter region of tumor suppressor genes leading to transcriptional silencing (Kulis & Esteller, 2010). Several genomewide DNA methylation studies have been developed in head and neck squamous cell carcinomas (Basu et al., 2017; Degli Esposti et al., 2017; Furusawa et al., 2011; Jithesh et al., 2013; Langevin et al., 2015); however, there are a few global DNA methylation analyses in OTSCC (Lim et al., 2016; Zhang et al., 2013). Since DNA methylation is an early event in cancer, comprehensive genomewide DNA methylation represents an opportunity to better understand the underlying molecular mechanisms and to identificate of potential DNA methylation markers for early cancer detection and better management of the disease. Liquid biopsies have emerged as a novel strategy for DNA methylation profiling in cancer (Constâncio et al., 2020). In this line, a recent metaanalysis has highlighted the potential of salivary DNA methylation for detecting head and neck cancer (RapadoGonzález, MartínezReglero, et al., 2021). Saliva is the biological fluid of the oral cavity, and oral exfoliated tumor cells released into saliva can reflect the methylation landscape of the tumor, representing an attractive tool for noninvasive detection, monitoring treatment, and recurrence disease detection. In this research, we conducted the first methylomewide study in OTSCC patients using the MehtylationEPIC array which comprises over 850,000 CpG sites. After differential methylation analysis, a set of differentially methylated CpGs (DMCpGs) was identified. Then, we tested the potential of tumorspecific methylated genes as salivary markers for early OTSCC detection. 2 | MATERIAL AND METHODS 2.1 | Study participants and sample collection This observational crosssectional study included six consecutive patients with earlystage OTSCC diagnosed between September 2020 and January 2021 at the Department of Oral and Maxillofacial Surgery from Complexo Hospitalario Universitario of A Coruña (CHUAC, SERGAS) in Galicia, Spain. Inclusion criteria included histologically confirmed OTSCC at earlystage pT1pT2N0M0 according to the 8th edition of the TumorNodeMetastasis (TNM)- staging system promulgated by the American Joint Commission on Cancer (AJCC) and the International Union Against Cancer (UICC). Regarding family history, no relevant data related to cancer were identified. With respect to past medical history, no patient received surgery, chemotherapy, or radiotherapy before sample collection. This study was approved by the Galician Ethics Committee of Clinical Research (Ref. No. 2018/003) and carried out in accordance with principles outlined in the Declaration of Helsinki (World Medical Association, 2013). Written informed consent was obtained from all participants before enrollment in the study and sample collection. Tumor and epithelium without dysplasia adjacent to tumoral tissue were obtained from tissue biopsy (Figure S1). The histologic sections were subsequently confirmed by an expert pathologist (R.A.R). Saliva was also collected from OTSCC patients at the time of cancer diagnosis. Additionally, two salivary samples from healthy individuals were obtained. Unstimulated saliva samples were collected using Danasaliva sample collection kit (Danagene) according to the manufacturer's instructions. The participants were asked to refrain from eating, drinking, smoking, and using oral hygiene products for at least 1 h prior to sampling collection. After collection, saliva was centrifuged at 2600 g for 15 min at room temperature, separating the cellular pellet from the cellfree salivary supernatant. Then, the CB16/12/00328; Universidade de Santiago de Compostela, Grant/Award Number: 2021PU007 tissue. Hypermethylation of 11 tumorspecific genes was validated in OTSCC TCGA cohort. Of these 11 genes, A2BP1, ANK1, ALDH1A2, GFRA1, TTYH1, and PDE4B were also hypermethylated in saliva. These six salivary methylated genes showed high diagnostic accuracy (≥0.800) for discriminating patients from controls. Conclusions: This is the first largest genomewide DNA methylation study on OTSCC that identifies a group of novel tumorspecific DNA methylation markers with diagnostic potential in saliva. KEYWORDS DNA methylation, genomewide methylation, Infinium MethylationEPIC array, Methylome, Oral tongue squamous cell carcinoma, saliva 16010825, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/odi.14444 by Universidade de Santiago de Compostela, Wiley Online Library on [15/06/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License | 3 RAPADO-GONZÁLEZ et al. cellular pellet was resuspended with phosphatebuffered saline (PBS) and preserved at −20°C prior to assay. A description of the clinicopathological characteristics of patients is shown in Table S1. 2.2 | DNA extraction and bisulphite modification DNA extraction from formalinfixed paraffinembedded (FFPE) tissue samples was performed using the AllPrep® DNA/RNA FFPE kit (Qiagen) according to the manufacturer's recommendations. Genomic DNA from tumor samples was eluted in 100 μl elution buffer AE and stored at −20 °C for further use. The QIAmp DNA Blood Mini Kit (Qiagen) was used to extract DNA from cell pellets of saliva samples according to the manufacturer's instructions. Genomic DNA from cell pellets was eluted in 100 μl elution buffer AE and stored at −20 °C for further use. DNA concentration and quality were measured by Nanodrop Spectrophotometer (Thermo Fisher Scientific) and by Qubit 4.0 Fluorometer (Thermo Fisher Scientific USA). Agilent's TapeStation 4200 (Agilent Technologies) was used to assess the integrity of the extracted DNA. Before bisulphite conversion, the DNA from FFPE samples was undergone to Illumina FFPE QC kit for quality control (QC) as per the manufacturer's protocol. Then, 400 ng of extracted DNA from each sample was used for bisulfite conversion using the EZ DNA Methylation Kit (Zymo Research) according to the manufacturer's instructions. After bisulphite conversion, only FFPEisolated DNA was further subjected to restoring with the Illumina Infinium HD FFPE Restoration kit as per the protocol. Genomewide DNA methylation profiling was carried out using the Illumina Infinium MethylationEPIC BeadChip (EPIC array; Illumina), which overspreads the DNA methylation profile across approximately 850,000 CpGs. A total of 4 μl of the bisulphite converted DNA was used as template for target preparation to hybridize on the BeadChip and was processed as per the manufacturer's instructions. To avoid the batch effect, all samples were processed together. The Bead Chip was imaged on HiScan System and intensity values (IDAT files) were extracted. 2.3 | Genomewide DNA methylation analysis Methylation array data were processed in the R 4.1.1 statistical environment (https://www.rproje ct.org/). The R package (R Core Team, 2021) RnBeads was used for quality control and preprocessing (Assenov et al., 2014). Firstly, a greedycut algorithm was used to filter out probes and/or samples. Next, probes overlapping with single nucleotide polymorphisms, probes in sex chromosomes, and probes whose sequence maps to multiple genomic locations (crossreactive) were removed. Raw intensities obtained in the array were normalized using the BMIQ method, which is a modelbased normalization approach to correct β values of type II probes according to the beta distribution of β values of type I probes. For each CpG site, a specific β value was obtained which is defined as the ratio of fluorescent signal between methylated (M) probe relative to the sum of the M and unmethylated (UM) probes (β = M/(M + UM)). The β values range from 0 (no methylation) to 1 (100% methylation of both alleles). Finally, hierarchical linear models were used with limma package to obtain the differences between groups (Ritchie et al., 2015). pValues were corrected for multiple testing (false discovery rate, FDR) using the Benjamini– Hochberg method, and a threshold of p < 0.05 was selected for significance. A similar bioinformatic pipeline was applied to all the comparisons performed. 2.4 | TCGA methylation analysis Methylation data were obtained from The Cancer Genome Atlas database (TCGA, 2013). Infinium Human Methylation 450 K BeadChip (HM450K) IDAT files from 36 earlystage OTSCC based on pathological diagnosis and four healthy controls were downloaded from TCGA database using the TCGAbiolinks R package. The IDAT files were subjected to preprocessing, normalization, and QC steps similar to our EPIC array files. 2.5 | TCGA gene expression analysis To identify the functional relevance of the differential methylated genes in OTSCC, we used normalized RNAsequencing (RNAseq) expression data from 36 earlystage OTSCC and four healthy controls from TCGA database (TCGA, 2013) using the TCGAbiolinks R package. Differentially expressed genes were obtained using generalized linear models with a function included in the package. A FDR < 0.05 and |log2 fold change (FC)| ≥1 was set as the criteria for screening differentially expressed genes. The genes differentially expressed between OTSCC and normal tissue TCGA samples were linked with the differentially methylated genes in our methylation analysis by Venn diagram using the VennDiagram R package. 2.6 | Saliva DNA methylation analysis To identify salivary DNA methylation changes, we obtained salivary DNA methylation data of four healthy controls (HM450K) from the NCBI Gene Expression Omnibus (GEO) under accession number GSE55734. In addition to the IDAT files from these four healthy controls, we also retrieved for salivary methylation analysis the EPIC array data from two saliva samples previously analyzed in our research group. A FDR < 0.05 and pvalue < 0.05 were set as the criteria for screening DMCpG sites. 2.7 | Gene ontology and functional pathway analysis Gene Ontology (GO) enrichment analysis was performed using GeneCodis 4 (TabasMadrid et al., 2012) to classify the differentially methylated genes into categories of cellular component, biological 16010825, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/odi.14444 by Universidade de Santiago de Compostela, Wiley Online Library on [15/06/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License 4 | RAPADO-GONZÁLEZ et al. process, and molecular function (Ashburner et al., 2000). In addition, the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis was performed to detect the potential pathways of the differentially methylated genes (Kanehisa & Goto, 2000). 2.8 | Receiver operating characteristic analysis To assess the diagnostic value of candidate DNA methylation markers, receiver operating characteristic (ROC) analysis was performed using the pROC package in R (Robin et al., 2011). The area under the curve (AUC) and accuracy were calculated with the aforementioned package to assess the performance of each DNA methylation marker. Principal component analysis (PCA) and ROC curve analysis were assessed to obtain the best models for OTSCC diagnosis in saliva. 3 | RESULTS 3.1 | Identification of tumorspecific differentially methylated CpG sites After preprocessing the EPIC array data, a total of 785,613 methylation CpG sites provided a highquality methylation signal. The comparison of the mean methylation values of all CpGs sites analyzed between tumor and adjacent nontumoral tissue samples showed a high correlation across all CpGs (Spearman, r2 = 0.99, p < 2.2 × 10−16), indicating a similar global DNA methylation pattern in both groups (Figure S2). Of them, a total of 25,890 tumorspecific DMCpGs spanning 7210 genes with FDR < 0.05 were identified between OTSCC and adjacent nontumoral tissue samples. The circos representation of DNA methylation locations in the genome showed that these 25,890 DMCpGs were widely distributed throughout all the chromosomes (Figure 1a). In addition, among these DMCpGs, we found a higher number of CpGs hypomethylated (20,505 CpGs; 79% of all DMCpGs) than hypermethylated (5385 CpGs; 21% of all DMCpGs) in OTSCC respect to the adjacent nontumoral tissue, as indicated in Figure 1b. In terms of the genomic distribution, DMCpGs were mostly distributed in the gene body (17,152 hypomethylated and 3875 hypermethylated) followed by the 5´UTR (2875 hypomethylated and 1881 hypermethylated), TSS1500 (2569 hypomethylated and 1314 hypermethylated), TSS200 (533 hypomethylated and 1218 hypermethylated), 1st exon (374 hypomethylated and 987 hypermethylated), 3´UTR (644 hypomethylated and 184 hypermethylated), and ExonBnd (344 hypomethylated and 49 hypermethylated). In terms of the CpG context, DMCpGs were mostly distributed in the opensea (17,592 hypomethylated and 1701 hypermethylated) followed by the island (250 hypomethylated and 2565 hypermethylated), shore (1371 hypomethylated FIGURE 1 Tumorspecific DMCpGs. (a) Circos plot of the differentially methylated CpGs (DMCpGs) in OTSCC respect to adjacent nontumoral tissue analyzed by the Infinium HumanMethylation450 BeadChip platform. The hyperand hypomethylated CpGs in OTSCC across the autosomes are indicated in red and green, respectively. The perimeter of the circular figure represents the human chromosomes. X and Y chromosomes were excluded from the analysis. (b) Graph showing percentages of hypermethylated and hypomethylated DMCpGs between OTSCC and the adjacent nontumoral tissue analyzed on the Infinium HumanMethylation450 BeadChip platform. (a) (b) 16010825, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/odi.14444 by Universidade de Santiago de Compostela, Wiley Online Library on [15/06/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License | 5 RAPADO-GONZÁLEZ et al. and 940 hypermethylated), and shelf (1292 hypomethylated and 179 hypermethylated) (Figure 2). Unsupervised hierarchical clustering and a heatmap using the methylation values of random 1000 DMCpGs, of which 44 were hypermethylated and 956 hypomethylated, were performed. The heatmap showed two DNA methylation clusters indicating that DNA methylation profiling of FIGURE 2 Genomic features of DMCpGs. (a, c, e) graphs showing percentages of hypermethylated and/or hypomethylated DMCpGs according to their functional genomic distribution. (b, d, f) graphs showing percentages of hypermethylated and/or hypomethylated DMCpGs according to their CpG content/neighborhood context. DMCpGs, differentially methylated CpGs; UTR, untranslated region; TSS, transcription start site (a) (c) (e) (f) (d) (b) 16010825, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/odi.14444 by Universidade de Santiago de Compostela, Wiley Online Library on [15/06/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License 6 | RAPADO-GONZÁLEZ et al. OTSCC was significantly different from adjacent nontumoral tissue (Figure S3). A genomewide differential methylation analysis in specific promoter regions including CpGIs and CpG shore was performed between tumor and adjacent nontumoral tissue. A total of 2700 tumorspecific DMCpGs spanning 1597 genes were identified of which 2216 were hypermethylated (82%) and 484 hypomethylated (18%) (Figure 3a). Unsupervised hierarchical clustering analysis and a heatmap using the methylation values of random 1000 DMCpGs within/near to CpGIs, of which 714 were hypermethylated and 286 hypomethylated, was produced. The heatmap showed two DNA methylation clusters indicating that DNA methylation profiling of OTSCC was significantly different from adjacent nontumoral tissue at specific promoter regions (Figure 3b). Enrichment analysis of the 807 genes containing DMCpGs revealed several important functional pathways and ontologies (Figures S4 and S5). 3.2 | Evaluation of the impact of differential methylation on associated gene expression The expression of 19,947 genes was accessible from TCGA data, of which 2494 genes were differentially expressed between OTSCC and normal tissue samples (FDR < 0.05 and |logFC| ≥1). Specifically, the expression of 73 genes from the 807 methylated genes identified in our methylation assay was also available in the TCGA data. A total of 40 genes hypermethylated in our assay showed a downregulated expression in tissue samples from TCGA cohort. By the contrary, 33 genes hypomethylated showed an upregulated expression (Figure 4). 3.3 | Validation of tumorspecific methylated CpGs using TCGA OTSCC data Since promoter methylation can modulate tumorspecific gene expression leading to its transcriptional silencing, the methylation levels of 11 genes, whose hypermethylation in their promoter region was associated with a downregulated gene expression, were further validated using TCGA OTSCC methylation data as an independent sample cohort. We evaluated the overall overlap of the DMCpGs located at the promoter region of 40 hypermethylated and 33 hypomethylated genes selected in our discovery assay with the tumorspecific DMCpGs identified from 36 OTSCC and four normal tissue samples which were available in the TCGA database. Although HM450K array used in the TCGA data was limited in genome coverage with almost half the number of CpGs respect FIGURE 3 Tumorspecific DMCpGs in promoter regions including CpGI and CpG shore. (a) Graph showing percentages of hypermethylated and hypomethylated DMCpGs; (b) Unsupervised hierarchical heatmap clustering of the random 714 hypermethylated DMCpGs and 286 hypomethylated DMCpGs between tumor and adjacent nontumoral tissue samples analyzed on the Infinium HumanMethylation450 BeadChip platform. DNA methylation values were represented as colors, with red representing hypermethylated DMCpGs and green representing hypomethylated DMCpGs. The row represents individual CpGs, and the column represents individual samples. DMCpGs, differentially methylated CpGs; T, tumor; N, nontumoral tissue (epithelium without dysplasia adjacent to tumoral tissue) (a) (b) 16010825, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/odi.14444 by Universidade de Santiago de Compostela, Wiley Online Library on [15/06/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License | 7 RAPADO-GONZÁLEZ et al. to HM850K array, we identified overlapping in 68/88 (77.27%) hypermethylated and 11/34 (32.35%) hypomethylated CpGs. Focused on the hypermethylation, 45 CpGs spinning 22 genes were differentially methylated (FDR < 0.05 and pvalue < 0.05) between tumor and normal tissue of which 11 genes (A2BP1, ACSS3, ALDH1A2, ANK1, CA3, GABRB3, GFRA1, HS3ST1, NDRG2, TTYH1, and PDE4B) were methylated in at least two CpGs. As shown in Figure S6 (Table S2), the differential methylation profiles of each CpG methylation loci of these 11 genes found in our OTSCC cohort were also observed in the TCGA OTSCC cohort, thus validating the findings of our discovery assay. Overall, these results suggest that these tumorspecific hypermethylated CpGs are frequent alterations in OTSCC and could be an early event with a crucial role in carcinogenesis. 3.4 | Accuracy of tumorspecific DNA methylation markers In addition, we performed ROC analysis, estimated AUC, and diagnostic accuracy values of the 11 tumorspecific promoter hypermethylated genes (Figure S7). Overall, these DNA methylation markers showed high sensitivity and specificity values (Table S3) indicating the potential of methylation profiling of these gene promoters for earlyOTSCC detection. 3.5 | Identification of methylated biomarkers in saliva from OTSCC patients in saliva To test the potential of tumorspecific CpG methylation markers in the A2BP1, ACSS3, ALDH1A2, ANK1, CA3, GABRB3, GFRA1, HS3ST1, NDRG2, TTYH1, and PDE4B genes as salivary biomarkers for early detection of OTSCC, we analyzed salivary cell pellets of four OTSCC patients and six healthy individuals. Differential methylation levels between saliva from OTSCC patients and healthy individuals were observed for specificsite methylated CpGs in A2BP1 (p = 0.005), ALDH1A2 (p = 0.033), ANK1 (p = 0.033), GFRA1 (p = 0.033), PDE4B (p = 0.012), and TTYH1 (p = 0.009) genes (Figure 5). The sensitivity observed for the individual markers ranged from 75% to 100%, and the specificity ranged from 83% to 100% (Figure 6 and Table S4). In addition, using PCA and ROC curve analysis, we were able to obtain the best combination of salivary biomarkers for the detection of OTSCC patients (Figure S8 and Table S5). 4 | DISCUSSION To the best of our knowledge, this is the first genomewide DNA methylation study in OTSCC that investigates aberrant DNA methylation at 850,000 CpG sites using the EPIC array. Our study showed FIGURE 4 Venn diagrams showing overlap of differentially hypermethylated CpGs in OTSCC array analysis with downregulated expression in TCGA database (a) and differentially hypomethylated CpGs in OTSCC array with upregulated expression in TCGA database (b) (a) (b) 16010825, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/odi.14444 by Universidade de Santiago de Compostela, Wiley Online Library on [15/06/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License 8 | RAPADO-GONZÁLEZ et al. an OTSCCspecific methylation pattern consisting of 25,890 CpGs that were differentially methylated compared to adjacent nontumoral tissue. We found that OTSCC had more CpGs hypomethylated (79%) than hypermethylated (21%), which is consistent with previously published studies (Basu et al., 2017; Krishnan et al., 2016). Most of the hypermethylated CpG sites (48%) were located in the CpGIs, while hypomethylated CpG sites (86%) were primarily enriched at the opensea region. In agreement with our results, a similar distribution of DMCpGs has been observed in oral squamous cell carcinoma using the HM450K array (Basu et al., 2017). Aberrant DNA promoter hypermethylation has been associated with the transcriptional silencing of different tumor suppressor genes in oral cancer (RapadoGonzález, LópezCedrún, et al., 2021). In this study, we focused on the CpGs that overlapped with the promoter region identifying a methylation pattern of 2700 tumorspecific CpGs mapped 807 unique genes differentially methylated in OTSCC compared to adjacent nontumoral tissue. For 40 of the 807 genes that presented hypermethylated CpGs in their promoter region, we observed a downregulated expression in the available TCGA OTSCC database, suggesting that they may play an important role in OTSCC carcinogenesis. To explore the consistency of our results, we further validated the methylation pattern of the DMCpGs of these 40 genes using genomewide DNA methylation data from the publicly available TCGA head and neck squamous cell carcinoma database. Due to the high heterogeneity of head and neck tumors, we have considered only data from 36 earlystage OTSCC patients in the TCGA dataset. For this approach, we compared data from two different microarray platforms (TCGAHM450K vs. EPIC) finding 45 CpGs sites located at 22 genes available in both databases. After differential methylation analysis, the methylation levels observed at the different CpG sites of 11 genes (A2BP1, ACSS3, ALDH1A2, ANK1, CA3, GABRB3, GFRA1, HS3ST1, NDRG2, TTYH1, and PDE4B) in our OTSCC cohort were also detected in TCGA OTSCC array data, thus confirming the differential methylation levels of these genes in an independent validation cohort. These genes could be functionally important for OTSCC cells and disruption of their control by aberrant DNA methylation could act as a driver event in the OTSCC development. Furthermore, these hypermethylated 11 genes showed a high and comparable diagnostic accuracy in both, our discovery and TCGA OTSCC, cohorts indicating the potential clinical utility of these methylated genes for early OTSCC detection. This study has identified novel convincing sitespecific DNA methylation changes in genes that have not previously been described in OTSCC. However, methylation of various of these genes has been reported in different cancers. Methylationassociated silencing of the ALDH1A2 gene has been observed in head and neck (Seidensaal et al., 2015), prostate (Kim et al., 2005), ovarian (Choi et al., 2019), and cervical (Zhou, Chen, et al., 2021) cancers. In oral cancer, the loss of ANK1 and miR4863p expression has been significantly correlated with ANK1 promoter hypermethylation. Interestingly, the FIGURE 5 Saliva DNA methylation levels (median) of significant hypermethylated DMCpGs in six genes in four OTSCC patients (T) and six healthy individuals (N). *p < 0.05 16010825, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/odi.14444 by Universidade de Santiago de Compostela, Wiley Online Library on [15/06/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License | 9 RAPADO-GONZÁLEZ et al. reexpression of ANK1 and miR4863p yielded to the reduction of the expression of DDR1 oncoprotein, decreasing proliferation and enhancing apoptosis in oral carcinoma cells (Chou et al., 2019). Also, loss of expression of CA3 due to aberrant methylation has been observed both in head and neck squamous cell lines and oral cancer samples (Pereira et al., 2018). Moreover, hypermethylation of CA3 has been reported in bladder cancer (Reinert et al., 2011). NDRG2 promoter hypermethylation has been associated with the transcriptional NDRG2 downregulation in numerous malignancies including meningiomas (Lusis et al., 2005), glioma (Skiriutė et al., 2014; Tepel et al., 2008; Zhou et al., 2014), liver cancer (Gödeke et al., 2016; Lee et al., 2008), colorectal cancer (Feng et al., 2011; Piepoli et al., 2009), gastrointestinal cancer (Chang et al., 2013; Yamamura et al., 2017), and oral squamous cell carcinoma (Furuta et al., 2010). Interestingly, the downregulation of NDRG2 could be modulating the activation of PI3K/Akt signaling in oral carcinogenesis (Furuta et al., 2010). The GFRA1 gene displays an important role in tumor progression and metastasis (Cavel et al., 2012; Dong et al., 2020; Esseghir et al., 2007), and its methylation has been described in rectal (Wei et al., 2016), lung (Sato et al., 2013), and gastric (Liu et al., 2014) cancer. Also, aberrant methylation of GABRB3 and ACSS3 genes has been previously reported in head and neck (GuerreroPreston et al., 2014) and prostate (Zhou, Song, et al., 2021) cancer. Overall, these findings evidence the important role of these epigenetic alterations in the carcinogenesis process; however, future research should be addressed to elucidate their exact biological function in OTSCC carcinogenesis. In the present study, we have also explored the possibility of detecting aberrant promoter hypermethylation in saliva from OTSCC patients. Several groups have examined the potential utility of DNA methylation in saliva for oral cancer diagnosis (GonzálezPérez et al., 2020; Liyanage et al., 2019; Nagata et al., 2012; Puttipanyalears et al., 2018; Srisuttee et al., 2020); however, most of these studies have analyzed the promoter hypermethylation of common tumor suppressor genes such as MGMT, P16, DAPK1, or RASSF1A. By the contrary, only a few studies have interrogated the DNA methylation profile in saliva from head and neck tumors (Langevin et al., 2015; Viet & Schmidt, 2008). Viet et al. identified in saliva a classifier based on 41 gene loci from 34 genes with potential value for oral cancer diagnosis through the Illumina GoldenGate Methylation Array (Viet & Schmidt, 2008). Similarly, Langevin et al. analyzed by the HM450K array the DNA methylation in oral rinses from oral and pharyngeal carcinoma patients, and they reported a classifier based on 22 FIGURE 6 ROC curves of the six salivary DNA methylation markers for discriminating OTSCC patients from healthy controls. ROC curve and AUC values were estimated for A2BP1, ANK1, ALDH1A2, GFRA1, TTYH1, and PDE4B using R software. 16010825, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/odi.14444 by Universidade de Santiago de Compostela, Wiley Online Library on [15/06/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License