Genome‐wide DNA methylation analysis identifies potent CpG signature for temzolomide response in non‐G‐CIMP glioblastomas with unmethylated MGMT promoter: MGMT ‐dependent roles of GPR81

Abstract Purposes To identify potent DNA methylation candidates that could predict response to temozolomide (TMZ) in glioblastomas (GBMs) that do not have glioma‐CpGs island methylator phenotype (G‐CIMP) but have an unmethylated promoter of O‐6‐methylguanine‐DNA methyltransferase (unMGMT). Methods The discovery‐validation approach was planned incorporating a series of G‐CIMP−/unMGMT GBM cohorts with DNA methylation microarray data and clinical information, to construct multi‐CpG prediction models. Different bioinformatic and experimental analyses were performed for biological exploration. Results By analyzing discovery sets with radiotherapy (RT) plus TMZ versus RT alone, we identified a panel of 64 TMZ efficacy‐related CpGs, from which a 10‐CpG risk signature was further constructed. Both the 64‐CpG panel and the 10‐CpG risk signature were validated showing significant correlations with overall survival of G‐CIMP−/unMGMT GBMs when treated with RT/TMZ, rather than RT alone. The 10‐CpG risk signature was further observed for aiding TMZ choice by distinguishing differential outcomes to RT/TMZ versus RT within each risk subgroup. Functional studies on GPR81, the gene harboring one of the 10 CpGs, indicated its distinct impacts on TMZ resistance in GBM cells, which may be dependent on the status of MGMT expression. Conclusions The 64 TMZ efficacy‐related CpGs and in particular the 10‐CpG risk signature may serve as promising predictive biomarker candidates for guiding optimal usage of TMZ in G‐CIMP−/unMGMT GBMs.


| INTRODUC TI ON
Glioblastomas (GBMs) are the most frequent and devastating brain malignancy with a reported median survival of only 12 ~ 15 months. 1 Temozolomide (TMZ) stands for the most effective adjuvant chemotherapy for GBMs. 2 It can alkylate genomic DNA of tumor cells at multiple sites and induce cell cycle arrest and apoptosis. 2[5] TMZ was generally beneficial to GBMs that harbor a methylated (me)MGMT promoter and consequently have a low expression of MGMT.However, TMZ was not that sensitive to GBMs harboring an unmethylated (un)MGMT promoter and generally a high expression of MGMT and tended to yield variable outcomes in those tumors.
7][8] Therefore, discovering powerful biomarkers that are predictive of TMZ response for subpopulation with TMZ-resistant unMGMT tumors can be clinically useful. 9,102][13][14][15][16][17][18] Because of the small proportion (10% of all GBMs) and the very distinct molecular and clinical features, the G-CIMP subtype was excluded, and we mainly focused on the majority of GBMs that do not have G-CIMP.In this study, by integrating genome-wide DNA methylation microarray data and clinical information, we identified a panel of 64 CpG candidates with potential linkage to TMZ efficacy in G-CIMP-GBMs with unMGMT promoter, from which a 10-CpG risk signature was constructed and validated to robustly and stably predict TMZ response.Bioinformatic and in vitro experimental analyses further provided biological insights into the 10-CpG predictive signature.

| Determination of G-CIMP and MGMT promoter methylation status
The G-CIMP status and promoter methylation status of MGMT were defined using DNA methylation data from Illumina HumanMethyla-tion27k/450 k microarrays (Illumina Inc.), as previously reported. 19,20thin the Illumina platform, methylation signal of each CpG was summarized as β value which provides a continuous and quantitative measurement of DNA methylation ranging from 0 (completely unmethylated) to 1 (completely methylated). 21The β value was commonly used for intuitive biological interpretation. 21When statistical analysis was employed, β value was transformed to M value, which has a Logit transformation relationship with β value. 21The G-CIMP status was determined by K-means (k = 3) clustering using the 1,503 featured probes as previously described. 19The MGMT promoter methylation was determined by a logistical regression model using two Illumina probes (cg12434587 and cg12981137). 20The logistical regression formula was calculated as logit (y) = 4.3215 + 0.5271 × Mvalue of cg12434587 + 0.9265 × M-value of cg12981137, with a cutoff of 0.358 for stratifying methylated and unmethylated tumors. 20

| CpG selection and model construction
The discovery-validation approach was utilized to identify CpG candidates with potential linkage to TMZ efficacy.The workflow for study design is shown in Figure 1B.In this study, RT/TMZ-treated unMGMT GBMs from TCGA and GEO (GSE22891, GSE50923, and GSE60274) and RT-treated unMGMT tumors from GSE60274 were preset as discovery sets.RT/TMZ-treated unMGMT GBMs from RAUH and RT-treated unMGMT tumors from TCGA were preset as validation sets (Figure 1B).Initial CpG selection was done as reported in our previous study. 29Batch effects across datasets were adjusted using M-value transformation 21 and the Empirical Bayes method (ComBat module, GenePattern). 30CpGs with higher variability in methylation levels (standard deviation [SD] of β value in discovery sets ≥0.15) were kept, and their M-values were used to correlate with overall survival (OS) in each discovery set treated with either RT/TMZ or RT alone.After removing inconsistent results from univariate Cox regression analyses, a panel of 64 CpGs was identified as TMZ efficacy-related candidates.To further reduce data dimensionality, the 64-CpG methylation data were incorporated into a multivariate Cox regression model using Likelihood Ratio (LR) and Forward selection approach, where age, sample source, and G-CIMP status were together adjusted.Finally, we identified a total of ten CpGs to construct a RISK-score formula, which is the sum of the M-value of each CpG weighted by its corresponding multivariate Cox coefficient.The median RISK score value in the RT/TMZ-treated discovery sets was pre-defined as the cutoff for low-risk and high-risk groups.

| Gene set enrichment analysis
Gene set enrichment analysis (GESA) was run through the Gene Ontology Gene Set collection from Molecular Signatures Database (MSigDB), to evaluate the biological profiles of the risk subgroups, with both nominal p-value ≤0.05 and false discovery rate (FDR) ≤ 0.25 as statistical significance. 31The proportion of

MGMT methylation G-CIMP Treatment Methylome/ Transcriptome
Transcriptome (A) (B) estimated using single-sample GSEA (ssGSEA) approach and the 782-gene signature reported by Charoentong et al. 32 The abundance of each TIIC type was summarized as normalized enrichment scores (NES).

| Cell culture and treatment
Two GBM cell lines (A172 and T98G) were obtained from American

| Plasmids and cell transfection
For in vitro gene silencing, RNA interference using plasmids expressing human short hairpin RNAs (shRNAs) was performed.
Transfection efficiency was verified by quantitative real-time PCR (qRT-PCR).

| TMZ cytotoxicity assay
After plasmid transfection, T98G and A172 cells were seeded on 96 well plates (5000 cells per well) and were treated with TMZ at the final concentrations of 7.5, 15, 30, 60, 120, 240, and 480 μM for 48 h.CCK-8 reagent was added to wells (10 μL/well) and incubated at 37°C for 1 h.The absorbance at 450 nm was measured for calculating the half maximal inhibitory concentration (IC 50 ).

| Flow cytometry assay
For cell apoptosis analysis, transfected T98G and A172 cells were treated with TMZ at different concentrations for 48 h.The cells were then washed with cold phosphate-buffered saline (PBS, 4°C).Annexin V-fluorescein isothiocyanate (FITC) and propidium iodide (PI) double staining were used to sort cells in early or late apoptotic phase.S4.In discovery sets, HCA on the 64-CpG methylation data divided all G-CIMP−/unMGMT GBMs into four clusters (Figure 2A).Survival analyses showed that OS of RT/TMZ-treated G-CIMP-GBMs significantly differed while that of RT-treated tumors were undistinguished across the clusters (Figure 2B,C).Similarly, HCA also yielded four clusters of G-CIMP−/unMGMT GBMs in validation sets and confirmed the distinct prognostic correlations with different treatments (RT/TMZ vs. RT alone; Figure 2D-F).Together the results indicated that the 64-CpG panel may serve as an epigenetic biomarker pool that may provide predictive information on TMZ response in G-CIMP−/unMGMT GBMs.

| Identification of a ten-CpG signature with potential linkage to TMZ efficacy in G-CIMP−/unMGMT GBMs
To simplify the HCA-based classification and develop a clinically applicable prediction model, we performed multivariate Cox regression analyses and reduced the 64-CpG panel into a panel of 10 CpGs, each of which was supposed to provide significant, complementary, and independent prediction information for TMZ sensitivity in G-CIMP−/unMGMT GBMs (Table 1).Accordingly, a RISK-score formula was established as follows: Using the median RISK score (0.3028) in the RT/TMZ-treated discovery sets as the cutoff, we divided all G-CIMP−/unMGMT GBMs into low-risk and high-risk subgroups.In RT/TMZ-treated TCGA cohort, low-risk G-CIMP−/unMGMT patients had longer OS than high-risk patients (Figure 3A).Similarly, significant results were also observed in RT/TMZ-treated GEO cohorts (GSE22891, GSE50923, and GSE60274 collectively; Figure 3A).The risk classification in each GEO cohort was also shown in Figure 3C.By contrast, in RT-treated GSE60274 cohort, OS was not significantly different between the risk subgroups (Figure 3B).In addition, survival data of patients with G-CIMP or meMGMT GBMs were shown for comparison (Figure 3).
To further evaluate the performance of the 10-CpG signature, we applied the RISK-score signature to validation sets with different treatments.Expectedly it accurately predicted OS in RT/TMZ-treated RAUH cohort of G-CIMP−/unMGMT GBMs but was not associated with OS in RT-treated TCGA cohort (Figure 3D,E).Patient-level and cohort-level meta-analyses both yielded distinct prognostic correlations of the 10-CpG signature among different treatment subpopulations (Figure 3F-H).Similar results were observed in terms of PFS outcome (Figure 4).Cox regression analyses confirmed the 10-CpG signature as a potent and independent survival predictor for G-CIMP−/unMGMT GBMs treated with RT/TMZ, rather than RT alone (Table S1).The results indicated that, instead of a general prognostic factor regardless of treatment, the 10-CpG signature may have a specific linkage to TMZ efficacy in G-CIMP−/unMGMT GBMs.

| The 10-CpG signature may be a potent predictive factor aiding in TMZ decision-making
To account for potential bias of assigned treatment regimens, only patients receiving standard RT (SRT) with or without (concurrent or adjuvant) TA B L E 1 Genomic and clinical information of the 10 CpGs with potential linkage to TMZ efficacy.TMZ were included for the interaction analyses between different risk subgroups (low-risk vs. high-risk) and different treatments (SRT/TMZ vs.

Relation to CpG Island
SRT).In line with previous report, 33 G-CIMP−/unMGMT GBMs appeared not to benefit much from SRT/TMZ as compared to SRT (Figure 5A,D).
The interaction analyses showed that SRT/TMZ appeared to confer an firmed SRT/TMZ as a better option for low-risk patients, but not for high-risk ones (Table S2).Together, the results indicated that the 10-CpG signature may represent a promising predictive model for TMZ response and be helpful for selecting patients who are likely to benefit from the addition of TMZ. from each dataset were combined by meta-analysis, where the inverse-variance approach was applied using either fixed-or random effect models based on the heterogeneity test, with I 2 ≥ 50% or p value ≤0.05 considered to be statistically significant.Hazard ratio (HR) and 95% confidence interval (CI) for survival curves were presented in Table S4.In TCGA G-CIMP−/unMGMT GBMs, the 10-CpG signature was found to be significantly associated with different gene expression subtypes; low-risk tumors were enriched with proneural subtype, and high-risk tumors were enriched with mesenchymal subtype (Figure 6A).In addition, the treatment cycles of TMZ were significantly higher in low-risk tumors (median cycle: 4) than in high-risk tumors (median cycle: 3; Figure 6A).GSEA showed that low-risk tumors were enriched with gene sets related to normal brain development and function (Figure 6B; Table S3), while high-risk tumors were enriched with a variety of cancer-promoting signatures and especially those relevant to DNA damage response, glucose metabolism, fatty acid metabolism, NF-kB activation, extracellular matrix (ECM) remodeling, immune response, and immune cell function (Figure 6C; Table S3).ssGSEA showed that high-risk tumors were associated with a higher abundance of some TIIC types and particularly those immunosuppressive cells (e.g., regulatory T cells [Treg], myeloidderived suppressor cells [MDSCs]; Figure 6D).The results suggested that the enhanced TMZ resistance in high-risk tumors may be attributable to a complex network of multiple tumorigenic or chemoresistant mechanisms, rather than a single molecular player.

| GPR81 may exhibit distinct impacts on TMZ resistance that depend on MGMT status
The 10-CpG methylation and the expression of their corresponding genes between NTBs and GBMs with each MGMT methylation status were presented in Figure S1.Pearson correlation analyses showed that only one CpG-gene pair (cg13702536 and GPR81) showed stably and significantly negative correlations between CpG methylation and gene expression (Figure 7A).This CpG-gene pair was thus selected for further analyses.GPR81 mRNA and protein levels were not significantly different between NTBs and GBMs (Figure 7B,C and Figure S1).The expression of GPR81 was also not significantly correlated with that of MGMT (Figure 7C,D).Interestingly, meta-analyses showed that both GPR81 methylation and expression appeared to be significantly associated with OS in unMGMT tumors, but not in meMGMT tumors when treated with RT/TMZ, indicting potential MGMT-dependent impacts of GPR81 on TMZ efficacy in G-CIMP-GBMs.
To test this hypothesis in in vitro experiments, we selected two GBM cell lines (A172 and T98G) which were reported to have similarly high GPR81 expression but distinct MGMT methylation and expression status in GSE68379 28 ; specifically A172 was characteristic of high methylation and low expression of MGMT, while T98G was characteristic of low methylation and high expression of MGMT (Figure 7I and Figure S2).The expression of GPR81 and MGMT was validated in our A172 and T98G cells (Figure S3), fol-

| DISCUSS ION
TMZ has long been recognized as the first-choice chemotherapy for treating primary GBMs. 1,2Unfortunately, it cannot benefit every GBM patients and many are resistant to this genotoxic drug.S4.

F I G U R E 5
The predictive performance of the 10-CpG signature in G-CIMP−/unMGMT GBMs; Survival difference between different treatments (SRT/TMZ vs. SRT) in G-CIMP−/unMGMT GBMs from (A) TCGA and (D) GSE60274 cohort; Interaction analysis between treatments (SRT/TMZ vs. SRT alone), and risk subgroups (low-risk vs. high-risk) in (B,C) TCGA and (E,F) GSE60274 cohort.Pooled survival comparisons using patient-level data between different treatments (SRT/TMZ vs. SRT) in (G) low-risk and (H) high-risk subgroups; Metaanalysis using cohort-level data between different treatments (SRT/TMZ vs. SRT) in each risk subgroup incorporating only younger patients (≤70 years).Survival difference of each treatment subgroup was tested by the log-rank test with p value ≤0.05 as statistical significance.Hazard ratios [HR] from each dataset were combined by meta-analysis, where the inverse-variance approach was applied using either fixedor random effect models based on the heterogeneity test, with I 2 ≥ 50% or p value ≤0.05 considered to be statistically significant.   of a straightforward relationship between its detection and TMZ choice in GBMs. 3 Although TMZ yielded much reduced benefits to unMGMT tumors as compared to meMGMT tumors, it is unlikely to withdraw from standard treatment, since there is lack of effective alternative therapies, and TMZ still benefits for some unMGMT cases. 3However, it should be noted that TMZ is not a cost-effective anti-GBM therapy, and its overuse can result in overconsumption of health resources, raise medical cost to caregivers, and increase risk of drug toxicity. 34Therefore, identifying potent predictive biomarkers, other than MGMT methylation, that can be useful for selecting subgroups of unMGMT patients with good sensitivity to TMZ, may represent a promising approach for optimizing decision-making on TMZ.
DNA methylation represented ideal biomarker candidates for precision oncology. 356][37] The last and most appealing advantage is the availability of epigenetic drugs that could reverse aberrant DNA methylation modifications, making it not only an indicator of certain features of a given cancer but also a druggable target to cure the disease. 38 this study, by integrating epigenome data, survival outcome, and treatment information of multi-sourced GBM cohorts, we identified a panel of 64 CpGs that may be specifically linked to TMZ efficacy in G-CIMP−/unMGMT GBMs.To construct a clinically applicable prediction model, we employed a multi-step selection workflow to screen out an optimal combination of a few number of CpGs, each of which not only conferred potent and independent prediction ability but also coordinated with and complemented each other.
Finally, a 10-CpG panel was identified and combined using a RISKscore model.Testing the 10-CpG signature in different cohorts of G-CIMP−/unMGMT GBMs showed that the defined low-risk tumors were stably associated with better OS than high-risk tumors when treated with RT/TMZ but not RT alone.So, it is inferred that the risk signature may be informative of distinct TMZ efficacy in G-CIMP−/ unMGMT GBMs, instead of a treatment-independent prognostic biomarker. 6Furthermore, the interaction analyses revealed that, as compared to RT alone, RT/TMZ was more beneficial to low-risk patients but yielded similar OS outcomes in high-risk patients.0][41] Like our study, Chai et al. 39 reported a 31-CpGs risk signature that predicted survival of TMZ-treated unMGMT GBMs.Ye et al. 40 reported a prognostic 13-gene risk signature that was validated in four RT/TMZtreated cohorts of IDH wild-type (wt) and unMGMT GBMs.Li et al. 41 proposed a 6-lncRNA immune-relevant risk signature that predicted survival in IDHwt/unMGMT GBMs.Table 2 compares the published signatures with our risk signature.The 10-CpG signature appeared to have a good predictive ability than the published models, with the highest AUC values at 1 year, 2 years, and 3 years in TCGA samples (Table 2).Also the present study may have advantages in the following aspects.First, abundant sample sources with relatively large sample size were used for discovery and validation of the risk model.
Second, the treatment information was incorporated into the CpG selection, which is a key variable to distinguish a predictive factor from a prognostic one. 6Third, the predictive value of the risk signature in our study was observed with a prospective objective on building a predictive model for TMZ response, instead of a spurious finding from a post-hoc subgroup analysis.Finally, interaction analyses were performed to compare the survival benefits of different treatment regimens in each risk subgroup, which could provide a direct guide on TMZ usage in specific subpopulations.
The biological implications of the 10-CpG signature may provide molecular clues behind its predictive ability for TMZ response.As highlighted by previous studies 42 that a complex and intertwined network of multiple molecular mechanisms may together determine the therapeutic resistance of GBMs, our bioinformatic analyses showed that the enhanced TMZ resistance observed in high-risk tumors may be partially attributable to the high enrichment of various cancer-promoting or therapy-resistant signatures involving in DNA damage response, energy metabolism, NF-kB activation, ECM remodeling, and tumor immunity, as well as an increased abundance of F I G U R E 6 Molecular and biological correlations of the 10-CpG signature using TCGA multi-omics data; (A) heatmaps of the methylation levels (M-values) of the 10 CpGs; each row represents a CpG and each column represents a sample which is ranked by its risk score.Clinical and molecular features are indicated for each sample, and multivariable Cox coefficients are indicated for each CpG; Representative GSEA enrichment plots of the highly enriched gene sets in (B) low-risk tumors and in (C) high-risk tumors; (D) the abundance of adaptive and innate immune infiltrating cells between low-risk and high-risk tumors.Categorical data (e.g., gene expression, age subgroup, therapies, and the use of bevacizumab) were tested by Chi-square test.Data of TMZ cycles did not pass the normality test and were compared using Mann-Whitney U test.Data of Normalized enrichment scores (NESs) passed the normality test, and were compared using Student t-test.Statistical significance was indicated at the level of ns >0.05, * < 0.05, ** <0.01, *** < 0.001 and **** < 0.0001.ns, non-significant.DNA methylation represents one critical layer of control in gene expression. 35,37,38It is reasonable to assume that the multi-CpG signature may contribute to TMZ resistance via regulating the expression of relevant genes.In our signature, only one CpG-gene pair (cg13702536 and GPR81) was found to show stable and significant correlation between DNA methylation and gene expression across different datasets, suggesting that GPR81 may be epigenetically controlled by DNA methylation.GPR81 has been reported to have multifunctional roles in promoting malignant behaviors of tumor cells by regulating energy metabolism, 43 angiogenesis, 44 therapeutic resistance, 45 and tumor immunity. 46,47 Functional reports on the other CpG-relevant genes may also provide biological clues for the predictive ability of our risk signature.GJB6 (harboring cg03473518) encodes a tumor-suppressive gap junction protein and may prevent GBM growth via rewiring glucose metabolism and inhibiting stemness. 48,49KCNQ1 (harboring cg12578166) encodes a voltage-dependent K+ channel and acts both as a target gene and regulator of the Wnt/β-catenin pathway. 50ss of KCNQ1 has been reported to exert anti-tumor functions via promoting epithelial-to-mesenchymal transition (EMT) and disrupting adheren junctions in epithelial cancers. 50 cg14329157), also called dynein assembly factor with WD repeats 1 (DAW1), belongs to the WD-repeat domain (WDR) family and plays vital roles in cilia motility. 51WDR69 hypermethylation was found to be associated with unfavorable prognosis in hepatocellular carcinoma. 52TNFRSF10D (harboring cg22783363) encodes a plasma membrane-located TNF-related apoptosis-inducing ligand (TRAIL) decoy receptor, and negatively regulates TRAIL-induced apoptosis. 53Hypermethylation and silencing of TNFRSF10D have been reported to occur in multiple cancer types and be associated with poor survival and resistance to DNA-damaging drugs. 53,54FABP6 (harboring cg22861316), encoding a bile acid-binding protein, is physiologically involved in fatty acids metabolism. 556][57] In GBM cells, FABP6 inhibition reversed the malignant phenotypes of tumor cells and increased TMZ sensitivity. 55BARX2 (harboring cg26221631) encodes a member of the homeobox transcription factor family that controls cell adhesion and cytoskeleton remodeling. 58Downregulation of BARX2, partially by CGI hypermethylation, has been reported to correlate with enhanced aerobic glycolysis and aggressive behaviors of tumor cells and be indicative of poor prognosis. 59POMC (harboring cg16302441) encodes a pro-hormone that gives rise to various active peptides such as adrenocorticotropic hormone (ACTH) and melanocyte stimulating hormones (MSHs). 60Some POMC-derived peptides have been reported to have vital roles in neuroendocrine tumors such as guiding optimal choice of chemotherapy. 60By far, little is known about the relevance of C4orf17 (harboring cg26647453) and UNKL (harboring cg26728422) in cancers.The multi-CpG signature may unlikely impact TMZ resistance via direct transcriptional control of the above genes as no significant correlations were observed for the nine CpG-gene pairs.However, in addition to classical epigenetic regulation mechanism, DNA methylation abnormalities may have broader biological effects by affecting heterochromatin structures, leading to loss of epigenetic regulation and resulting in hypervariability of gene expression. 36,37,61Future studies are needed to explore the molecular machinery on TMZ resistance behind the CpG members that may not directly control local gene expression.
Limitations exist in the present study.Our finding should be carefully interpreted due to the following shortcomings, such as lack of validation in a randomized setting or in a prospective manner, potential patient selection bias in retrospectively collected cohorts, very few samples for RT-treated patients, heterogeneous treatment regimens, and incomplete clinical data.Moreover, the risk signature was built on high-throughput DNA methylation detection platform, which is not clinically available and not economical for routine testing.Therefore, the risk signature in its current form is not ready for daily use and should be modified and validated by a more common detection system, such as pyrosequencing. 47Finally, the histochemical characterization of glioma in our study is technically conventional and has limited scope.New histo-methodology, such as tissue clearing and quantitative ultramicroscopy, may provide more comprehensive molecule-histology information. 62 conclusion, we firstly reported a panel of 64 TMZ efficacyrelated CpGs and then built a 10-CpG-based RISK-score signature that may robustly and stably predict response to TMZ in G-CIMP−/ unMGMT GBMs, a subtype characteristic of high TMZ resistance.
Experimental data revealed potential MGMT-dependent roles of GPR81 in TMZ resistance, highlighting the complexity of the chemoresistant mechanisms in GBMs.The 10-CpG signature may be helpful for guiding TMZ choice in such subpopulation.Future studies are needed to explore the molecular mechanisms underlying the risk signature and to translate it into routine practice.
28 tumor-infiltrating immune cell (TIIC) types in tumor bulks was F I G U R E 1 Patient dataset information and study workflow; (A) molecular subtype, available molecular/clinical data, and sample size of included patient datasets; (B) schematic diagram for searching and validating TMZ efficacy-related CpG candidates and a 10-CpG signature.

|
DNA methylation variations of a 64-CpG panel were specifically associated with TMZ efficacy in G-CIMP−/unMGMT GBMsFollowing the selection workflow (Figure1B), we identified a panel of 64 TMZ efficacy-related CpGs from discovery sets, each of which was significantly correlated with OS of G-CIMP−/ unMGMT GBMs treated with RT plus TMZ, rather than RT alone.
benefit to low-risk G-CIMP−/unMGMT GBMs (Figure 5B,E) but was associated with similar OS in high-risk patients (Figure 5C,F).No significant difference was observed in patient baseline information (e.g., surgery, gender, or KPS) between subgroups with different treatment and different risk subgroups except for a high proportion of elderly patients (>70 years) in RT-treated subgroups (data not shown).Patient-level and cohort-level meta-analyses only incorporating younger patients (≤ 70 years) were thus performed, which yielded a statistically significant OS difference between SRT/TMZ-treated and SRT-treated subgroups (Figure 5G-I).Cox regression analyses con-

F I G U R E 3
Prognostic performance of the 10-CpG signature in G-CIMP−/unMGMT GBMs in term of OS outcome; survival comparison among the subgroups characteristic of G-CIMP+, meMGMT, G-CIMP−/unMGMT low-risk or G-CIMP−/unMGMT high-risk in (A) RT/TMZtreated and (B) RT-treated discovery sets, as well as in (C) each GEO cohort; survival comparison among the above indicated subgroups in (D) RT/TMZ-treated and (E) RT-treated validation sets; pooled survival comparison using patient-level data among the above indicated subgroups with (F) RT/TMZ or with (G) RT alone; (H) meta-analysis using cohort-level data between the risk subgroups with RT/TMZ or RT alone.Survival difference of each subgroup was tested by the log-rank test with p value ≤0.05 as statistical significance.Hazard ratios[HR]

3 . 4 |
Clinical and molecular correlations of the 10-CpG signature FigureS4).The results together indicated that the contributions of GPR81 to TMZ resistance may depend on MGMT expression status in GBMs; briefly, GPR81 may enhance TMZ sensitivity when GBM cells highly expressed MGMT, while GPR81 may enhance TMZ resistance when MGMT expression was absent or largely repressed in GBM cells.

2 F I G U R E 4
Molecular biomarkers have been increasingly reported for aiding TMZ Prognostic performance of the 10-CpG signature in G-CIMP−/unMGMT GBMs in term of PFS outcome; survival comparison among the subgroups characteristic of G-CIMP+, meMGMT, G-CIMP−/unMGMT low-risk or G-CIMP−/unMGMT high-risk in (A) RT/TMZtreated TCGA cohort, (B) RT/TMZ-treated RAUH cohort, and (C) RT-treated TCGA cohort.Survival difference of each subgroup was tested by log-rank test with p value ≤0.05 as statistical significance.Hazard ratio (HR) and 95% confidence interval (CI) for survival curves were presented in Table non-G-CIMP, RT/TMZ or RT alone) which MGMT methylation status stands for the most widely validated predictive biomarker. 1,2However, MGMT methylation had limited use in guiding TMZ in clinical practice due to lack

F I G U R E 7
The impacts of GPR81 on TMZ resistance of GBM cells that may depend on MGMT status; (A) Pearson correlation coefficients for each CpG-gene pair from included cohorts; (B).The IHC scores of GPR81 between NTB and GBM samples from Neurosurgery Department, Xijing Hospital; (C) Representative IHC images of GPR81 and MGMT in NTB or GBM samples, with corresponding IHC scores; (D) Pearson correlation between MGMT and GPR81 at protein and mRNA levels in local samples; (E,F) Survival difference between low versus high methylation of GPR81 among RT/TMZ-treated (E) G-CIMP−/unMGMT and (F) G-CIMP−/meMGMT GBMs; the median methylation value (M-value: 2.0708) from RT/TMZ-treated G-CIMP−/unMGMT GBMs was used for stratifying low versus high methylation; (G,H) Meta-analyses for (G) GPR81 methylation-based groups and (H) GRP81 expression-based groups in RT/TMZ-treated non-G-CIMP GBMs with each MGMT methylation status; the median expression value (Z-score: −0.1992) from RT/TMZ-treated G-CIMP−/unMGMT GBMs was used for stratifying low versus high expression.(I) Methylation and expression status of MGMT in A172 and T98G cells from GSE68379; (J) Validation of GPR81 knockdown and MGMT overexpression in A172 cells by qRT-PCR; (K) Validation of GPR81 knockdown and MGMT knockdown in T98G cells by qRT-PCR; (L) GPR81 knockdown increased TMZ sensitivity and cell apoptosis to TMZ treatment in A172 cells originally with no detectable MGMT expression; (M) GPR81 knockdown decreased TMZ sensitivity and cell apoptosis to TMZ treatment in MGMT-overexpressed A172 cells; (N) GPR81 knockdown decreased TMZ sensitivity and cell apoptosis to TMZ treatment in T98G cells originally expressing MGMT; (O) GPR81 knockdown increased TMZ sensitivity and cell apoptosis to TMZ treatment in MGMT-silenced T98G cells; (P) Western bolt results in A172 and T98G cells treated with TMZ; Hazard ratios [HR] from each dataset were combined by metaanalysis, where the inverse-variance approach was applied using either fixed-or random effect models based on the heterogeneity test, with I 2 ≥ 50% or p value ≤0.05 considered to be statistically significant.All continuous data passed normality test except for IHC scores.Statistical significance was indicated at the level of ns >0.05, * < 0.05, ** <0.01, *** < 0.001 and **** < 0.0001.ns, non-significant.
Surprisingly, in in vitro GBM cell experiments, we have revealed the potential MGMT-dependent impacts of GPR81 on TMZ resistance; specifically in GBM cells with high methylation (or low expression) of MGMT, e.g., A172-and MGMTsilenced T98G cells, GPR81 may enhance TMZ resistance while in GBM cells with low methylation (or high expression) of MGMT, e.g., MGMT-overexpressed A172 and T98G cells, GPR81 may increase TMZ sensitivity.In line with experimental data, survival analyses also supported the distinct predictive abilities of GPR81 expression (or methylation) in RT/TMZ-treated G-CIMP-GBMs with different MGMT methylation statues.However, the absence of apparent predictive ability in meMGMT tumors indicated that the tumor intrinsic GPR81 expression may not act as a major contributor to TMZ resistance among the complex molecular mechanisms conferred by the entire tumor microenvironments in meMGMT samples.By contrast, tumor intrinsic GRP81 expression may be a dominant player for TMZ efficacy among the unMGMT GBM microenvironments as supported by the significant predictive ability in clinical samples and the significant impacts on TMZ resistance in GBM cell lines.In summary, the clinical and experimental data of GPR81 may provide additional layer of evidence supporting the predictive ability of the 10-CpG signature in G-CIMP−/unMGMT GBMs.Moreover, the MGMT-dependent roles of GPR81 highlighted the complexity and sophistication of the underlying molecular mechanisms that eventually define the resistant nature of GBMs.However, by far, the data are too preliminary to draw conclusion.Future studies are needed to address how MGMT affects the functions of GPR81 in TMZ resistance of GBMs, and what MGMT status triggers the function transition of GPR81.