# A five-CpG signature of microRNA methylation in non-G-CIMP glioblastoma

En-Ming Kang, An-An Yin, Ya-Long He, Wei-Jun Chen, Amandine

Etcheverry, Marc Aubry, Jill Barnholtz-Sloan, Jean Mosser, Wei Zhang, Xiang
Zhang

► To cite this version:

En-Ming Kang, An-An Yin, Ya-Long He, Wei-Jun Chen, Amandine Etcheverry, et al.. A five-CpG sig-nature of microRNA methylation in non-G-CIMP glioblastoma. CNS Neuroscience and Therapeutics, 2019, 25 (9), pp.937-950. 10.1111/cns.13133 . hal-02121050

HAL Id: hal-02121050

https://univ-rennes.hal.science/hal-02121050

Submitted on 7 Jul 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
ORIGINAL ARTICLE

CNS Neuroscience & Therapeutics
WILEY

WILEY

# A five-CpG signature of microRNA methylation in non-G-CIMP glioblastoma

En-Ming Kang<fnref n="1" /> ID | An-An Yin<fnref n="1" />,<fnref n="2" /> ID | Ya-Long He<fnref n="1" /> | Wei-Jun Chen<fnref n="3" /> |

Amandine Etcheverry<fnref n="4" />,<fnref n="5" />,<fnref n="6" /> | Marc Aubry<fnref n="5" />,<fnref n="7" /> | Jill Barnholtz-Sloan<fnref n="8" /> | Jean Mosser<fnref n="4" />,<fnref n="5" />,<fnref n="6" />,<fnref n="7" /> |

Wei Zhang<fnref n="1" /> | Xiang Zhang<fnref n="1" />

Correspondence

Xiang Zhang, Department of Neurosurgery,
Xijing Institute of Clinical Neuroscience,
Xijing Hospital, Air Force Medical University,
Xi'an, China.

Emails: xzhang@fmmu.edu.cn, xzhangneurosurgery@aliyun.com

Funding information

National Natural Science Foundations of
China, Grant/Award Number: 81471266,
81671302 and 81672909

Summary

**Aims:** DNA methylation has been found to regulate microRNAs (miRNAs) expression, but the prognostic value of miRNA-related DNA methylation aberration remained largely elusive in cancers including glioblastomas (GBMs). This study aimed to investigate the clinical and biological feature of miRNA methylation in GBMs of non-glioma-CpG island methylator phenotype (non-G-CIMP).

**Methods:** Prognostic miRNA methylation loci were analyzed, with TCGA and Rennes cohort as training sets, and independent datasets of GBMs and low-grade gliomas (LGGs) were obtained as validation sets. Different statistical and bioinformatic analysis and experimental validations were performed to clinically and biologically characterize the signature.

**Results:** We identified and validated a risk score based on methylation status of five miRNA-associated CpGs which could predict survival of GBM patients in a series of training and validation sets. This signature was independent of age and O-6-methyl-guanine-DNA methyltransferase (MGMT) promoter methylation status. The risk subgroup was associated with angiogenesis and accordingly differential responses to bevacizumab-contained therapy. MiRNA target analysis and in vitro experiments further confirmed the accuracy of this signature.

**Conclusion:** The five-CpG signature of miRNA methylation was biologically relevant and was of potential prognostic and predictive value for GBMs. It might be of help for improving individualized treatment.

KEYWORDS

angiogenesis, DNA methylation signature, glioblastoma, miRNA, prognostication

<footnote n="1">¹Department of Neurosurgery, Xijing
Institute of Clinical Neuroscience, Xijing
Hospital, Air Force Medical University, Xi'an,
China</footnote>

<footnote n="2">²Department of Neurosurgery, The 88th
Hospital of the People's Liberation Army,
Taian, China</footnote>

<footnote n="3">³Department of Emergency Medicine,
Jinling Hospital, Medical School of Nanjing
University, Nanjing, China</footnote>

<footnote n="4">⁴CNRS, UMR 6290, Institut de Génétique et
Développement de Rennes (IGdR), Rennes,
France</footnote>

<footnote n="5">⁵UEB, UMS 3480 Biosit, Faculté de
Médecine, Université Rennes1, Rennes,
France</footnote>

<footnote n="6">⁶CHU Rennes, Service de Génétique
Moléculaire et Génomique, Rennes, France</footnote>

<footnote n="7">⁷Plate-forme Génomique Santé
Biosit, Université Rennes1, Rennes, France</footnote>

<footnote n="8">⁸Case Comprehensive Cancer Center, Case
Western Reserve University, Cleveland,
Ohio</footnote>
# 1 | INTRODUCTION

Glioblastomas (GBMs) are the most common and devastating sub-types of primary central nervous system tumors.¹ Unfortunately, despite the multimodal treatment of surgical resection, radiotherapy, and chemotherapy, the reported median survivals of GBM patients were only 16-19 months.¹⁻³

Cancer-specific DNA methylation changes play important roles in cancer development and progression. The best-known epigenetic abnormality in cancers is promoter-specific CpG island (CGI) hypermethylation of tumor suppressor genes which consequently cause transcriptional silencing.⁴ Altered DNA methylation affected the expressions of not only protein-coding genes but also noncoding RNAs (ncRNAs).⁵ Among those ncRNAs, microRNAs (miRNAs), the 20-22 nucleotides small ncRNAs, have been demonstrated to have multiple roles in the pathogenesis of cancers.⁶ It has been reported that miRNAs could be regulated by DNA methylation and abnormal methylation in miRNAs was closely correlated with cancer progression.⁶,⁷ However, the biological and clinical implications of miRNA methylation abnormality were largely unstudied in cancers including GBMs. Glioma-CpG island methylator phenotype (G-CIMP) represents a distinct subgroup of glioma which is featured by genome-wide hypermethylated CGIs and favorable prognosis.⁸ The G-CIMP+ tumors have already been broadly studied, while the relevance features of non-G-CIMP GBMs remain largely unclear.

In this study, we analyzed miRNA methylation data of non-G-CIMP GBMs from The Cancer Genome Atlas (TCGA), Gene Expression Omnibus (GEO), and Rennes cohort⁹ to reveal the relationship between miRNA methylation and GBM survival. Bioinformatic methods and in vitro experiments were used to validate our results.

# 2 | MATERIALS AND METHODS

## 2.1 | GBM datasets

Rennes cohort of 77 newly diagnosed non-G-CIMP GBMs with clinical and genome-wide DNA methylation microarray data by Infinium HumanMethylation450k BeadChip (Illumina Inc, San Diego, CA, USA) was obtained from the ArrayExpress under the accession number "E-MTAB-4969."⁹ All patients received standard adjuvant treatment of radiotherapy (RT) and concurrent temozolomide (TMZ). Public DNA methylation datasets of non-G-CIMP GBM samples were also downloaded from The Cancer Genome Atlas (TCGA) data portal,¹⁰ and Gene Expression Omnibus (GEO) under the accession number "GSE60274."¹¹ (Detailed clinical data of and relative CpG information are listed in the Supporting Information S1) We also obtained clinical and DNA methylation data of LGGs from TCGA¹² and GSE48462.¹³ Among the heterogeneous datasets, only those with age over 18 years old and a molecular diagnosis of non-G-CIMP tumors were included in this study. For survival analysis, patients with a follow-up data >1 month were included, in order to reduce the bias caused by noncancer death.¹⁰ In addition, nontumor brain tissues were obtained from apparently healthy individuals or chronic epilepsy patients with

<figure><img src="image_2.png" /><figcaption>**FIGURE 1** Schematic diagram of the probe selection workflow for the study</figcaption></figure>

pathological evidence of other neurological or psychiatric diseases in each dataset. The G-CIMP status was determined by K-means (k = 3) clustering on the 1503 probes reported by Noushmehr et al¹⁴ MGMT (O-6-methylguanine-DNA methyltransferase) promoter methylation status was determined by a logistic regression model using two CpGs, that is, cg12434587 and cg12981137.¹⁵ Batch effects from different datasets and platforms were adjusted by a nonparametric empirical Bayes approach (ber package).¹⁶ Methylation level of each integrated CpGs was summarized as M-value.¹⁷

## 2.2 | Construction and validation of a miRNA methylation-based risk score model

CpG probes were filtered by removing those targeting the X and Y chromosomes, containing a single nucleotide polymorphism (SNP) within five base pairs of the targeted CpG. We then selected probes
annotated with miRNAs (n = 2448) for this study. The discovery phase was performed within TCGA and Rennes cohort (training sets). Univariate Cox regression analysis with permutation test was performed within each training set. Potential prognostic CpGs with consistent survival correlation (permutation *P* < 0.2) in each training set were subjected to multivariate Cox regression analysis within the combined training set (TCGA and Rennes collectively), and those with a *P* value < 0.05 were finally selected for risk score modeling. The risk score formula was constructed by integrating the *M*-values of all selected CpGs which were weighted by their multivariate Cox regression coefficients after adjusted by patient age and MGMT promoter methylation status.¹⁸,¹⁹ Patients were then classified into high-risk or low-risk groups with the cutoff point as the median risk score from the combined training set. The validation phase was performed in GSE60274 and datasets of LGGs and in particular those with wide-type IDH.

## 2.3 | Gene Set Enrichment Analysis (GSEA)

GSEA was performed to evaluate the functional gene expression profiles between the risk subgroups on reported gene sets from

Molecular Signature Database (MSigDB), with nominal *P* value ≤ 0.05 for significance.²⁰

## 2.4 | MiRNA target gene prediction and pathway analysis

The online databases TargetScan (http://www.targetscan.org/vert_72/). Accessed November 11, 2018), miRanda (http://www.microrna.org/microrna/home.do. Last update: 2010-11-01), and miRDB (http://mirdb.org/. Accessed November 11, 2018) were used to identify the target genes of the relative miRNAs. Genes appeared in all three databases were included for the following analysis.²¹ GO analysis was then performed with PANTHER (version 14.0 Released 2018-12-03) based on the GO database (http://www.geneontology.org/ Accessed January 11, 2019) for biological process (BP), cellular component (CC), molecular function (MF), and pathway enrichment of the predicted target genes.²² Fisher's two-side exact test was used to classify the GO categories, and the Bonferroni correction for multiple test was calculated to correct the *P* values. Bonferroni-corrected for *P* < 0.05 was considered to

<table><caption>**TABLE 1** Patient characteristics of included patient cohorts of non-G-CIMP GBMs</caption><thead><tr><th rowspan="2">Variables</th><th colspan="2">Training set</th><th>Validation set</th></tr><tr><th>Rennes cohort</th><th>TCGA</th><th>GSE60274</th></tr></thead><tbody><tr><td>Sample size</td><td>77</td><td>102</td><td>59</td></tr><tr><td>Clinical factors</td><td></td><td></td><td></td></tr><tr><td>Age</td><td></td><td></td><td></td></tr><tr><td>Median</td><td>60</td><td>63</td><td>52</td></tr><tr><td>Range</td><td>36-75</td><td>23-85</td><td>26-70</td></tr><tr><td>Pre-operative KPS<sup>a</sup></td><td></td><td></td><td></td></tr><tr><td>Median</td><td>80</td><td>80</td><td>NA</td></tr><tr><td>Range</td><td>40-100</td><td>40-100</td><td>NA</td></tr><tr><td>Gender</td><td></td><td></td><td></td></tr><tr><td>Male/Female</td><td>55/22</td><td>58/44</td><td>45/14</td></tr><tr><td>Extent of surgery</td><td></td><td></td><td></td></tr><tr><td>Surgery (total/partial)/ Biopsy</td><td>72 (55/17)/4</td><td>101 (NA/NA)/1</td><td>57 (NA/NA)/2</td></tr><tr><td>Adjuvant Treatments</td><td></td><td></td><td></td></tr><tr><td>RT + TMZ/RT</td><td>77/0</td><td>71/31</td><td>32/27</td></tr><tr><td>BVZ/non-BVZ/UN</td><td>29/32/16</td><td>NA</td><td>NA</td></tr><tr><td>Molecular factors</td><td></td><td></td><td></td></tr><tr><td>MGMT methylation status</td><td></td><td></td><td></td></tr><tr><td>Methylated/ Unmethylated</td><td>26/51</td><td>37/65</td><td>26/33</td></tr><tr><td>Gene expression subtype</td><td></td><td></td><td></td></tr><tr><td>P/N/C/M</td><td>18/6/24/27</td><td>20/13/37/30</td><td>8/4/17/20</td></tr><tr><td>TCGA methylation clusters</td><td></td><td></td><td></td></tr><tr><td>Clusters 2/3</td><td>29/48</td><td>35/67</td><td>23/36</td></tr></tbody></table>

KPS, Karnofsky performance score; NA, not available; RT, radiotherapy; TMZ, temozolomide; UN, unknown.

ªKPS was available for only a small subset of patients from TCGA cohort.
be significant. Enrichment analysis based on Kyoto Encyclopedia of Genes and Genomes (KEGG) was performed and visualized using ClueGO (Version 2.5.3),²³,²⁴ a Cytoscape (version 3.7.1, http://cytoscape.org/) plug-in. The main parameters for constructing the network were as follows: ontologies/pathways: KEGG (321 terms/ pathways with 7454 available unique genes, December 7, 2018), GO term/pathway selection: Min Gene = 3 & Min Percentage = 3.0%, Kappa Score Threshold = 0.5, Statistical Test Used = Enrichment/ Depletion (Two-sided hypergeometric test), Correction Method Used = Bonferroni step down. Only pathways with *P* value < 0.05 were considered to be significant.

## 2.5 | Cell culture and transfection

The human GBM cell lines U87MG, U251, T98G, and SHG44 were obtained from the cell bank of the Air Force Medical University (Xi'an, China) and cultured in Dulbecco's modified Eagle's medium (DMEM; Gibco, USA) supplemented with 10% fetal bovine serum (FBS; Shanghai BioSun Sci&Tech Co., Ltd., Shanghai, China) and incubated at 37°C with 5% CO₂. miR-1284 mimic (UCU AUA CAG ACC CUG GCU UUU C) and mimic negative control (mimic NC; UUC UCC GAA CGU GUC ACG UTT) were synthesized by Sangon Biotech Co., Ltd. (Shanghai, China). Cells for transfection were planted into 60-mm dishes at 4 × 10⁵ cells/well. After 48 hours incubation, miR-1284 mimic, mimic NC, or control (only treated with Polymer) was transfected into cells at 50 pmol/mL using Xfect RNA Transfection Reagent (Takara Bio, USA). The transfection efficiency was verified by real-time quantitative PCR (qPCR).

## 2.6 | RNA extraction and Real-time quantitative PCR

Total RNA was extracted by TRIzol reagent (Invitrogen, USA) and reverse transcribed with PrimeScript RT reagent kit (Takara, Tokyo, Japan). PCR amplification was performed in triplicate with SYBR Premix Ex Taq II (Takara) using CFX96 Real-Time PCR Detection System (Bio-Rad, Hercules, CA, USA). The expression values of miR-1284 were normalized to the levels of small nuclear U6. The primer sequences were listed as follows: (a) miR-1284: Reverse

transcription primer: 5'-CTC AAC TGG TGT CGT GGA GTC GGC AAT TCA GTT GAG GAA AAG-3'; (b) U6 Reverse transcription primer: 5'-CGC TTC ACG AAT TTG CGT GTC AT-3'; miR-1284-F: 5'-CGT CTA TAC AGA CCC TGG CTT TTC-3'; miR-1284-R: 5'-CTC AAC TGG TGT CGT GGA-3'; U6-F: 5'-CTC GCT TCG GCA GCA CAT A-3'; U6-R: 5'-CGC TTC ACG AAT TTG CGT G-3'.

## 2.7 | Pyrosequencing

Pyrosequencing was performed by Pyromark Q96 ID platform and analyzed by PyroMark CpG software (Qiagen, Germany). The following primers were used: miR-1284-F 5'-ATT TTT ATT GGT TAA ATT AAT ATT ATA GG-3', miR-1284-R biotin-5'-AAC TTA TTA CAT TAA ATA CAA ACA ACA AC-3', miR-1284-seq 5'-TTT TTA GTT TTT AAG TAT ATT-3'. The DNA methylation value for each sample was calculated as the average methylation value of the interrogated CpGs.²⁵

## 2.8 | 5-Aza-2'-deoxycytidine (5-Aza-dC) demethylation treatment

U251 and U87MG cells were grown for 4 days in the presence of 10 µmol/L 5-Aza-dC (Sigma-Aldrich, St. Louis, MO, USA). Fresh 5-Aza-dC was added every 24 hours.

## 2.9 | Cell proliferation assay

Cells with different treatments were implanted in 96-well plates at 5 × 10³ per well. At indicated time points, CCK-8 kit (Yeasen, Shanghai, China) was assayed for cell viability measurement.

## 2.10 | Cell cycle and apoptosis analysis

For cell cycle analysis, cells were harvested, fixed in 70% ethanol on ice, and stained with propidium iodide in phosphate-buffered saline containing RNase. The DNA contents were analyzed by flow cytometry. For cell apoptosis analysis, Annexin V-fluorescein isothiocyanate and propidium iodide double staining (Roche Diagnostics, Germany) was used to sort cells in early or late apoptotic phase.

<table><caption>TABLE 2 The five prognostic CpGs associated with miRNA</caption><thead><tr><td>Probes</td><td>Chr.</td><td>miRNA name</td><td>miRNA region</td><td>Relation to CpG island</td><td>Methylation status in GBM</td><td>Average M-value of high-risk GBMs<sup>a</sup></td><td>Cox regression coefficients</td></tr></thead><tbody><tr><td>cg05744073</td><td>17</td><td>miR-132</td><td>Body</td><td>Island</td><td>Hypermethylated</td><td>-4.073</td><td>-0.534</td></tr><tr><td>cg08244382</td><td>14</td><td>miR-127; miR-433</td><td>TSS1500;TSS200</td><td>Island Shore</td><td>Hypermethylated</td><td>3.185</td><td>-0.446</td></tr><tr><td>cg20382675</td><td>3</td><td>miR-1284</td><td>TSS200</td><td>Open sea</td><td>NS</td><td>0.287</td><td>-0.263</td></tr><tr><td>cg24082174</td><td>3</td><td>miR-1248</td><td>TSS1500</td><td>Island Shore</td><td>NS</td><td>0.991</td><td>0.255</td></tr><tr><td>cg13767001</td><td>13</td><td>miR-759</td><td>TSS1500</td><td>Open sea</td><td>Hypomethylated</td><td>-2.223</td><td>0.368</td></tr></tbody></table>

NS, no significance; TSS, transcription start sites.

Methylation level assessed with M-value: low (-Inf, -2), middle [-2, 2], and high (2, Inf).

<footnote>aIncluded all high-risk samples of three datasets.</footnote>
<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

**FIGURE 2** The survival correlation of the five-CpG signature in each dataset. A, The five-CpG signature predicted overall survival (OS) in training sets. B, The signature was validated by yielding apparent OS difference in GSE60274. C, The five-CpG signature was also able to predict PFS in Rennes cohort. D, The signature could not identify patients with different prognoses in IDH wide-type LGG (grade III or II)

## 2.11 | Wound-healing assay

Cell motility was assessed by wound-healing assay as described previously.²⁶ A scratch wound was generated by a 200 µL pipette tip on the confluent cell monolayers in 6-well plates. The spread of the wound closure was observed after 48 hours of the scratch.

## 2.12 | Statistical analysis

The distributions of known molecular and clinical features with respect to the risk subgroups were tested by Fisher's exact or chi-square test. Survival data, for example, overall survival (OS) and progression-free survival (PFS), were estimated by the Kaplan-Meier method and compared by log-rank test. Univariate and multivariate Cox regression models were performed to evaluate the correlation and independence of potential prognostic factors. For in vitro experiments, data were expressed as mean ± SEM from three independent experiments and analyzed by Student's t test. All the calculations were done within SPSS19.0 (IBM Corporation,

New York, NY, USA) and R software (version 3.2.5; https://www.r-project.org/), and a difference was considered significant when $P \le 0.05$.

# 3 | RESULTS

## 3.1 | Identification of prognostic miRNA methylation loci from the training sets

The included cohorts and the workflow of probe selection were schematically presented in Figure 1, and patient characteristics were summarized in Table 1. By employing a multistep selection criterion, we identified a five-CpG panel of miRNA methylation that showed consistent prognostic significance in both training sets (Table 2). Among the panel, two CpGs (eg, cg05744073 and cg08244382) were hypermethylated and one CpG (eg, cg13767001) was hypomethylated in GBMs, while the other two were not differentially methylated in GBMs (Table 2). Upon the correlation with prognosis, three CpGs (eg, cg05744073, cg08244382, and cg20382675) showed
<figure><img src="image_11.png" /><figcaption>**FIGURE 3** Molecular and clinical characteristics of the 5-CpGs miRNA methylation signature. A, the heat maps of K-means (k = 2) clustering on the 5-CpGs methylation signature according to the M-value from all GBM groups; each column represented a sample; for each sample (n = 238), subgroup correlation was indicated; P values for Fisher' exact test and chi-square test were accordingly shown; B, GSEA enrichment plots for representative functional gene sets enriched in high-risk tumors from TCGA. C, High-risk but not low-risk tumors conferred significant OS benefits when treated with bevacizumab in Rennes cohort with available second-line therapies</figcaption></figure>

<table><caption>**TABLE 3** Results of the miRNA methylation signature in Cox regression analysis</caption><thead><tr><th rowspan="2">Variables</th><th colspan="3">Univariate Cox model</th><th colspan="3">Multivariate Cox model</th></tr><tr><th>HR</th><th>95% CI</th><th>P value</th><th>HR</th><th>95% CI</th><th>P value</th></tr></thead><tbody><tr><td>Rennes (n = 61)<sup>a</sup></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Patient age</td><td>1.046</td><td>1.015-1.078</td><td><strong>0.003</strong></td><td>1.040</td><td>1.003-1.078</td><td><strong>0.033</strong></td></tr><tr><td>miRNA methylation signature</td><td>2.926</td><td>1.733-4.942</td><td><strong>&lt;0.001</strong></td><td>3.129</td><td>1.782-5.493</td><td><strong>&lt;0.001</strong></td></tr><tr><td>MGMT methylation status</td><td>0.438</td><td>0.236-0.813</td><td><strong>0.009</strong></td><td>3.047</td><td>0.140-0.569</td><td><strong>&lt;0.001</strong></td></tr><tr><td>DNA methylation clusters</td><td>0.849</td><td>0.492-1.465</td><td>0.557</td><td></td><td></td><td></td></tr><tr><td>Proneural subtype</td><td>0.905</td><td>0.483-1.695</td><td>0.754</td><td></td><td></td><td></td></tr><tr><td>BVZ treatment</td><td>0.607</td><td>0.357-1.031</td><td>0.065</td><td>0.536</td><td>0.273-1.049</td><td>0.069</td></tr><tr><td>Gender</td><td>0.918</td><td>0.522-1.614</td><td>0.767</td><td></td><td></td><td></td></tr><tr><td>Extent of surgery</td><td>0.957</td><td>0.623-1.469</td><td>0.840</td><td></td><td></td><td></td></tr><tr><td>TCGA + GSE60274 + Rennes (n = 238)<sup>b</sup></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Patient age</td><td>1.028</td><td>1.012-1.044</td><td><strong>0.001</strong></td><td>1.034</td><td>1.018-1.051</td><td><strong>&lt;0.001</strong></td></tr><tr><td>Treatments (RT/TMZ vs RT)</td><td>0.479</td><td>0.345-0.666</td><td><strong>&lt;0.001</strong></td><td>0.438</td><td>0.314-0.609</td><td><strong>&lt;0.001</strong></td></tr><tr><td>DNM methylation clusters</td><td>0.995</td><td>0.732-1.351</td><td>0.973</td><td></td><td></td><td></td></tr><tr><td>miRNA methylation signature</td><td>2.207</td><td>1.704-2.859</td><td><strong>&lt;0.001</strong></td><td>2.368</td><td>1.838-3.050</td><td><strong>&lt;0.001</strong></td></tr><tr><td>MGMT methylation status</td><td>0.627</td><td>0.455-0.863</td><td><strong>0.004</strong></td><td>0.589</td><td>0.427-0.812</td><td><strong>0.001</strong></td></tr><tr><td>Gender</td><td>1.009</td><td>0.732-1.392</td><td>0.956</td><td></td><td></td><td></td></tr></tbody></table>

KPS, Karnofsky performance score; NA, not available; RT, radiotherapy; TMZ, temozolomide.

<footnote><sup>a</sup>Rennes cohort excluded 16 patients with insufficient treatment information.</footnote>

<footnote><sup>b</sup>Including all patients from TCGA, Rennes cohort, and GSE60274.</footnote>

<footnote><sup>c</sup>The significance of bold values indicate P value < 0.05.</footnote>
<figure><img src="image_12.png" /><figcaption>**FIGURE 4** The survival correlation of the five-CpG signature within current GBM classification. A, The five-CpG signature predicted overall survival (OS) in both MGMT promoter methylated and unmethylated patients treated with both radiotherapy (RT) and temozolomide (TMZ). B, It was also correlated with different OS in subgroups of ≤60 or >60 y. C, The correlation between five-CpG signature and different prognoses was significant in proneural and neural subtypes and marginally significant in the classical and mesenchymal subtypes</figcaption></figure>

negative correlation with OS, while two CpGs (eg, cg24082174 and cg13767001) with positive correlation (Table 2).

Accordingly, the risk score model was constructed as follows: risk score = (-0.534 × M-value of cg05744073) + (-0.446 × M-value of cg08244382) + (-0.263 × M-value of cg20382675) + (0.254 × M-value of cg24082174) + (0.368 × M-value of cg13767001). With the cutoff as the median risk score from the combined training sets (-0.382), patients were divided into low-risk groups (with lower risk scores) and high-risk groups (with higher risk scores). In the combined training sets, the assigned low-risk patients (n = 89) were significantly associated with longer overall survival (OS) than those high-risk ones (n = 90; log-rank *P* < 0.0001; Figure 2A). The 5-CpG signature also showed consistent prognostic value in each training set (Figure 2A).

## 3.2 | Validation of the five-GpG miRNA methylation signature for prognostication

To validate the prognostic performance of the 5-CpG miRNA methylation signature, we applied it to the independent validation set of GSE60274. With the prespecified cutoff, patients were classified into a low-risk group (n = 28) and a high-risk group (n = 31). Consistent with the training sets, low-risk patients were associated with longer OS than high-risk ones (log-rank *P* = 0.013; Figure 2B). We also observed a significant correlation between progression-free survival (PFS) and the miRNA methylation-based risk groups in Rennes cohort (Figure 2C).

In addition, we applied the GBM-derived signature to independent validation cohorts of IDH wild-type LGGs. The miRNA methylation signature failed to yield significant OS differences between the risk subgroups within both grade III and II gliomas using their median risk scores as cutoffs, respectively, which supported the signature as a GBM-specific prognostic model (Figure 2D).

## 3.3 | Molecular and clinical correlation of the 5-CpG miRNA methylation signature

Correlation with current established molecular features showed that the 5-CpG signature appeared not to be significantly correlated with TCGA gene expression subtypes, DNA methylation clusters, and MGMT promoter methylation status (Figure 3A). Also, the signature seemed not to be correlated with treatments, gender, and age (Figure 3A). GSEA on gene expression data of TCGA samples showed that the high-risk tumors were enriched with pro-oncogenic gene sets such as ErbB signaling pathway (*P* < 0.0001), MAPK signaling pathway (*P* = 0.029), pro-angiogenic gene sets such as hypoxia (*P* = 0.035), and VEGF pathway (*P* = 0.029), which might biologically explain the inferior survival of those high-risk tumors (Figure 3B).

## 3.4 | High-risk patients appeared to be beneficial for bevacizumab therapy

As reported by GSEA, high-risk tumors seemed to be featured by upregulation of various pro-angiogenic gene sets (Figure 3B). Accordingly, we tested the potential survival benefits conferred by the anti-angiogenic agent bevacizumab as combined therapy within each risk subgroup. In Rennes cohort with available second-line therapies, we found that the addition of bevacizumab did confer significant OS benefits in high-risk tumors, but was associated with similar OS to bevacizumab-free therapy (Figure 3C).

## 3.5 | The 5-CpG signature was an independent prognostic factor in non-G-CIMP GBMs

Within Rennes cohort with RT/TMZ, univariate Cox regression analysis showed that patient age, MGMT promoter methylation
<figure><img src="image_13.png" /><figcaption>**FIGURE 5** Target prediction results of signature associated miRNAs</figcaption></figure>

status, and our miRNA methylation signature were significantly associated with OS (Table 3). Multivariate Cox regression analysis further demonstrated the prognostic independence of the above-mentioned factors (Table 3). Multivariate Cox regression model within the combined cohorts of TCGA, GSE60274, and Rennes cohort not only confirmed the prognostic independence of our miRNA methylation signature but also supported its treatment in-dependence (Table 3).

## 3.6 | The prognostic value of the miRNA methylation signature with respect to current GBM classification

We also tested the prognostic interrelationship of the 5-CpG sig-nature with known prognostic factors within available patients from the combined training and validation sets. We found that it could consistently predict OS within the subtypes of unmethylated
<figure><img src="image_14.png" /></figure>

<table><caption>TABLE 4 PANTHER analysis for predicted target genes</caption><thead><tr><th>Terms</th><th>Target gene</th><th>Expected gene Nr</th><th>Fold enrichment</th><th>P value</th></tr></thead><tbody><tr><td colspan="5">PANTHER GO-slim molecular function</td></tr><tr><td>MAP kinase activity</td><td>17</td><td>5.04</td><td>3.38</td><td>0.028</td></tr><tr><td>→Protein serine/threonine kinase activity</td><td>44</td><td>20.45</td><td>2.15</td><td>0.005</td></tr><tr><td>→Protein kinase activity</td><td>80</td><td>38.61</td><td>2.07</td><td>&lt;0.001</td></tr><tr><td>→Catalytic activity, acting on a protein</td><td>170</td><td>99.8</td><td>1.7</td><td>&lt;0.001</td></tr><tr><td>→Catalytic activity</td><td>396</td><td>324.05</td><td>1.22</td><td>0.012</td></tr><tr><td>Ubiquitin-like protein transferase activity</td><td>45</td><td>22.2</td><td>2.03</td><td>0.017</td></tr><tr><td>RNA polymerase II transcription factor activity, sequence-specific DNA binding</td><td>56</td><td>27.85</td><td>2.01</td><td>0.002</td></tr><tr><td>→DNA-binding transcription factor activity</td><td>125</td><td>75</td><td>1.67</td><td>&lt;0.001</td></tr><tr><td>→Transcription regulator activity</td><td>140</td><td>83.93</td><td>1.67</td><td>&lt;0.001</td></tr><tr><td>DNA binding</td><td>92</td><td>56.39</td><td>1.63</td><td>0.008</td></tr><tr><td>→Nucleic acid binding</td><td>193</td><td>125.21</td><td>1.54</td><td>&lt;0.001</td></tr><tr><td>→Binding</td><td>509</td><td>404.85</td><td>1.26</td><td>&lt;0.001</td></tr><tr><td>→Organic cyclic compound binding</td><td>196</td><td>129.71</td><td>1.51</td><td>&lt;0.001</td></tr><tr><td>Unclassified</td><td>738</td><td>859.06</td><td>0.86</td><td>0</td></tr><tr><td colspan="5">PANTHER GO-slim biological process</td></tr><tr><td>Regulation of transcription by RNA polymerase II</td><td>52</td><td>26.71</td><td>1.95</td><td>0.044</td></tr><tr><td>→Regulation of transcription, DNA-templated</td><td>55</td><td>28.84</td><td>1.91</td><td>0.038</td></tr><tr><td>→regulation of biological process</td><td>390</td><td>287.8</td><td>1.36</td><td>&lt;0.001</td></tr><tr><td>→Biological regulation</td><td>420</td><td>312.6</td><td>1.34</td><td>&lt;0.001</td></tr><tr><td>→Regulation of metabolic process</td><td>182</td><td>106.9</td><td>1.7</td><td>&lt;0.001</td></tr><tr><td>→Regulation of macromolecule metabolic process</td><td>157</td><td>91.18</td><td>1.72</td><td>&lt;0.001</td></tr><tr><td>→Regulation of gene expression</td><td>114</td><td>64.47</td><td>1.77</td><td>&lt;0.001</td></tr><tr><td>Transcription by RNA polymerase II</td><td>128</td><td>68.59</td><td>1.87</td><td>&lt;0.001</td></tr><tr><td>→Transcription, DNA-templated</td><td>162</td><td>93.7</td><td>1.73</td><td>&lt;0.001</td></tr><tr><td>→Cellular macromolecule biosynthetic process</td><td>186</td><td>118.42</td><td>1.57</td><td>&lt;0.001</td></tr><tr><td>→Metabolic process</td><td>460</td><td>364.94</td><td>1.26</td><td>&lt;0.001</td></tr><tr><td>→Biosynthetic process</td><td>192</td><td>123.53</td><td>1.55</td><td>&lt;0.001</td></tr><tr><td>→Macromolecule biosynthetic process</td><td>188</td><td>120.02</td><td>1.57</td><td>&lt;0.001</td></tr><tr><td>→Organic substance biosynthetic process</td><td>192</td><td>123.07</td><td>1.56</td><td>&lt;0.001</td></tr><tr><td>Cellular protein modification process</td><td>84</td><td>50.21</td><td>1.67</td><td>0.032</td></tr><tr><td>→Protein modification process</td><td>84</td><td>50.36</td><td>1.67</td><td>0.033</td></tr><tr><td>Unclassified</td><td>670</td><td>820.69</td><td>0.82</td><td>0</td></tr><tr><td colspan="5">PANTHER GO-slim cellular component</td></tr><tr><td>Nuclear chromatin</td><td>66</td><td>36.62</td><td>1.8</td><td>0.008</td></tr><tr><td>→intracellular part</td><td>398</td><td>311.92</td><td>1.28</td><td>&lt;0.001</td></tr><tr><td>→cell part</td><td>540</td><td>448.19</td><td>1.2</td><td>0.001</td></tr><tr><td>→cell</td><td>544</td><td>450.25</td><td>1.21</td><td>&lt;0.001</td></tr><tr><td>Unclassified</td><td>809</td><td>909.65</td><td>0.89</td><td>0</td></tr><tr><td>Extracellular space</td><td>22</td><td>56.31</td><td>0.39</td><td>&lt;0.001</td></tr><tr><td>→Extracellular region part</td><td>28</td><td>62.41</td><td>0.45</td><td>0.001</td></tr><tr><td>→Extracellular region</td><td>35</td><td>70.2</td><td>0.5</td><td>0.002</td></tr><tr><td colspan="5">PANTHER pathways</td></tr><tr><td>EGF receptor signaling pathway</td><td>29</td><td>10.15</td><td>2.86</td><td>0.001</td></tr><tr><td>Cadherin signaling pathway</td><td>33</td><td>11.98</td><td>2.75</td><td>&lt;0.001</td></tr><tr><td>FGF signaling pathway</td><td>24</td><td>9.16</td><td>2.62</td><td>0.015</td></tr><tr><td>CCKR signaling pathway</td><td>31</td><td>13.2</td><td>2.35</td><td>0.008</td></tr><tr><td>PDGF signaling pathway</td><td>26</td><td>11.37</td><td>2.29</td><td>0.046</td></tr><tr><td>Wnt signaling pathway</td><td>52</td><td>23.73</td><td>2.19</td><td>&lt;0.001</td></tr><tr><td>Unclassified</td><td>1304</td><td>1404.08</td><td>0.93</td><td>0</td></tr></tbody></table>
<figure><img src="image_15.png" /><figcaption>**FIGURE 6** Bioinformatic analysis of predicted target genes. A, PANTHER GO-Slim biological process. B, PANTHER GO-Slim molecular function. C, PANTHER GO-Slim cellular component. D, PANTHER pathway enrichment. E, KEGG pathway enrichment analysis, relative genes were shown as well</figcaption></figure>

or methylated MGMT tumors (Figure 4A), and the subgroups of ≤60 or >60 years (Figure 4B). Regarding the TCGA expression subtypes, the signature was significantly associated with OS in the proneural and neural subtypes, and also yielded marginally significant OS difference in the classical and mesenchymal subtypes (Figure 4C).

## 3.7 | Target gene prediction of the 5-CpG signature-related miRNAs

Targetscan, miRanda, and miRDB databases were used to predict the target genes of miR-132, miR-127, miR-433, miR-1284, miR-1248, and miR-759. To ensure the specificity and sensitivity of our prediction, we kept the identical targets predicted in all three databases without setting additional criterion. Totally 1578 target genes left for further functional analysis (Figure 5). Among them,

139 genes were candidate targets for two miRNAs, seven genes were common targets of three miRNAs, and one was target of four miRNAs.

## 3.8 | Biological characteristics of predicted target genes

Predicted target genes were further analyzed using PANTHER GO-slim tools on MF, BP, and CC (Table 4). The most enriched MF terms were MAP kinase activity, ubiquitin-like protein transferase activity, DNA binding, and RNA polymerase II transcription factor activity (Figure 6A). The most enriched BP terms were regulation of transcription by RNA polymerase II, transcription by RNA polymerase II, and cellular protein modification process (Figure 6B). The most enriched CC terms were nuclear chromatin and extracellular space (Figure 6C). PANTHER pathway analysis showed EGF receptor
<table><caption>TABLE 5 KEGG pathway enrichment analysis of predicted target genes</caption><thead><tr><th>GOID</th><th>GOTerm</th><th>Term P value<sup>a</sup></th><th>Group P value<sup>a</sup></th><th>Nr genes</th><th>Associated genes found</th></tr></thead><tbody><tr><td>KEGG:04520</td><td>Adherens junction</td><td>&lt;0.001</td><td>&lt;0.001</td><td>21</td><td>ACP1, ACTN4, BAIAP2, CREBBP, CTNND1, EP300, INSR, MAPK1, MAPK3, MET, NLK, SMAD4, SNAI1, SNAI2, SORBS1, SRC, SSX2IP, TGFBR1, VCL, WASL, YES1</td></tr><tr><td>KEGG:04110</td><td>Cell cycle</td><td>&lt;0.001</td><td>&lt;0.001</td><td>29</td><td>CCND2, CDC14A, CDC16, CDC27, CDC6, CDK2, CDKN2B, CDKN2D, CHEK2, CREBBP, CUL1, E2F3, E2F5, EP300, GADD45A, MAD2L1, MCM5, RAD21, RB1, SKP1, SMAD4, SMC1B, SMC3, STAG1, TTK, WEE1, YWHAB, YWHAG, YWHAQ</td></tr><tr><td>KEGG:04350</td><td>TGF-beta signaling pathway</td><td>0.001</td><td>&lt;0.001</td><td>21</td><td>ACVR1, ACVR1C, BMPR1B, CDKN2B, CREBBP, CUL1, E2F5, EP300, GDF5, ID2, ID4, MAPK1, MAPK3, NODAL, PPP2R1B, RPS6KB1, SKP1, SMAD4, SMAD5, SMAD7, TGFBR1</td></tr><tr><td>KEGG:04012</td><td>ErbB signaling pathway</td><td>0.004</td><td>&lt;0.001</td><td>20</td><td>CAMK2G, ERBB4, EREG, GRB2, HBEGF, KRAS, MAP2K1, MAP2K4, MAPK1, MAPK10, MAPK3, NCK2, PAK2, PAK4, PAK6, PIK3R1, PRKCB, RPS6KB1, SHC3, SRC</td></tr><tr><td>KEGG:04360</td><td>Axon guidance</td><td>0.005</td><td>&lt;0.001</td><td>32</td><td>ABLIM1, ABLIM3, ARHGEF12, BMPR1B, CAMK2G, EPHA5, EPHA7, GNAI1, KRAS, MAPK1, MAPK3, MET, NCK2, NFATC2, NRP1, NTNG1, PAK2, PAK4, PAK6, PIK3R1, PTCH1, ROCK2, RYK, SEMA3C, SEMA4A, SEMA4G, SLIT3, SRC, SGAP1, SGAP2, SSH2, UNC5B</td></tr><tr><td>KEGG:05211</td><td>Renal cell carcinoma</td><td>0.009</td><td>&lt;0.001</td><td>17</td><td>CREBBP, EP300, EPAS1, ETS1, GRB2, HIF1A, KRAS, MAP2K1, MAPK1, MAPK3, MET, PAK2, PAK4, PAK6, PIK3R1, RAP1A, SLC2A1</td></tr><tr><td>KEGG:04114</td><td>Oocyte meiosis</td><td>0.011</td><td>&lt;0.001</td><td>25</td><td>ADCY7, AURKA, CAMK2G, CDC16, CDC27, CDK2, CPEB2, CPEB4, CUL1, FBXW11, MAD2L1, MAP2K1, MAPK1, MAPK3, PPP1CB, PPP2R1B, PPP2R5C, RPS6KA3, SKP1, SMC1B, SMC3, SPDYE1, YWHAB, YWHAG, YWHAQ</td></tr><tr><td>KEGG:04218</td><td>Cellular senescence</td><td>0.048</td><td>&lt;0.001</td><td>28</td><td>CCND2, CDK2, CDKN2B, CHEK2, E2F3, E2F5, ETS1, FBXW11, FOXO3, GADD45A, HIPK1, HIPK3, KRAS, LIN52, LIN9, MAP2K1, MAP2K6, MAPK1, MAPK3, NBN, NFATC2, PIK3R1, PPP1CB, RASSF5, RB1, RRAS2, TGFBR1, TSC1</td></tr></tbody></table>

signaling pathway, cadherin signaling pathway, FGF signaling pathway, CCKR signaling pathway, PDGF signaling pathway, and Wnt signaling pathway were the most enriched pathways (Figure 6D). Then, we utilized ClueGO to make a KEGG pathway enrichment analysis (Figure 6E, Table 5). The most enrichment terms were adherens junction, cell cycle, TGF-beta signaling pathway, ErbB signaling pathway, axon guidance, renal cell carcinoma, oocyte meiosis, and cellular senescence.

## 3.9 | MiR-1284 suppressed glioma cell proliferation and migration

To further validate the functional relevance of this miRNA methylation signature, we selected miR-1284 for in vitro experiments. Pyrosequencing of cg20382675 showed that the miR-1284-associated CpG was consistently associated with high DNA methylation status in glioma cells, that is, U87MG, U251, T98MG, and SHG44 (Figure 7A). Accordingly, the expression levels of miR-1284 were

comparable in those glioma cells (Figure 7B). However, after treated with 5-Aza-dC, we found that the expressions of miR-1284 were significantly decreased in U251 and U87MG, indicating a positive impact of DNA methylation on miRNA expression (Figure 7C). By transferring miR-1284 mimics into U251, we established a miR-1284-overexpressed U251 model, which was validated by qPCR (Figure 7D). CCK-8 analysis showed that over-expression of miR-1284 reduced cell viability of U251 (Figure 7E). Flow cytometry analysis showed that over-expression of miR-1284 was also associated with lower frequency of tumor cells in S and G2 phase (Figure 7F), and higher frequency of apoptotic cells (Figure 7G). Finally, wound-healing assay showed that migration was inhibited by over-expression of miR-1284 (Figure 7H).

## 4 | DISCUSSIONS

This study investigated genome-wide DNA methylation microarray data of miRNA-associated CpGs to explore the clinical value of

<footnote><sup>a</sup>Corrected with Bonferroni step down.</footnote>
<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /></figure>

<figure><img src="image_18.png" /></figure>

<figure><img src="image_19.png" /></figure>

<figure><img src="image_20.png" /></figure>

<figure><img src="image_21.png" /></figure>

<figure><img src="image_22.png" /></figure>

<figure><img src="image_23.png" /></figure>

<figure><img src="image_24.png" /><figcaption>**FIGURE 7** Characteristics of miR-1284 in GBM cell lines. A, Methylation level of miR-1284 in glioma cell lines (U251, U87MG, T98G, and SHG44). B, Relative expression of miR-1284 compared with U6 in glioma cell lines. No difference was found between each cell line. C, Expression of miR-1284 by qRT-PCR in U251 and U87MG cells treated with 5-Aza-2-deoxycytidine (AZA). D, Expression of miR-1284 transfected with mimic and mimic NC for 48 h (P < 0.001). E, CCK-8 assay testing cell viability from 1 to 5 d. F, Flow cytometry detecting cell cycle of U251 and PI values in different groups (G). Flow cytometry testing cell apoptosis after transfection. H, Representative results of wound-healing assay and the percentage of healing area determined using the ImageJ. \*P < 0.05, \*\*P < 0.01, \*\*\*P < 0.001</figcaption></figure>

miRNA methylation in non-G-CIMP GBMs. We identified a 5-CpG signature of miRNA methylation which could predict survival of non-G-CIMP GBM patients. This signature showed consistent and robust prognostic values within each subgroup of different ages, molecular subtypes, and treatments and was validated in independent patient cohort. Notably, different risk groups distinguished by this signature showed different bevacizumab therapy outcomes. These findings suggest that methylation status of this 5-CpG relative miRNAs is closely correlated with GBM malignancies, especially with tumor angiogenesis, and this 5-CpG methylation signature has good potential to be an indicator for bevacizumab therapy in non-G-CIMP GBMs.

As GBMs are characterized by high heterogeneity and massive vessels, anti-VEGF therapy was expected to improve the outcome of GBM patients.²⁷ Bevacizumab, a humanized monoclonal antibody against VEGF, has been the most promising anti-angiogenic agents for treating GBMs and was approved for recurrent GBM treatment.²⁸ However, in newly diagnosed GBM patients, recent randomized trials failed to yield clear survival benefits when applied bevacizumab plus Stupp regimen,²⁹ implying that proper indicators are needed for bevacizumab treatment. In this study, the miRNA methylation-based risk subgroups were associated with differential enrichments of pro-angiogenic gene sets (eg, hypoxia or VEGF pathways), suggesting the possibility of differential responses to bevacizumab within the risk
subgroups. Accordingly, distinct survival benefits were observed in Rennes cohort with the use of bevacizumab at progression: High-risk patients seemed to benefit more from the bevacizumab-contained therapy. These results suggest that the 5-CpG signature is of potential use to optimize bevacizumab therapy by identifying appropriate patient candidates.

For biological features of this 5-CpG signature, target genes of relative miRNAs were analyzed with bioinformatic methods. Results showed these miRNAs regulated a great many genes and might cooperate with each other to regulate specific genes. The profound characteristics of the target genes were analyzed based on GO database and KEGG database. From the results of GO analysis of BF of these target genes, we can infer that these miRNAs might greatly participate in biological regulation, especially regulation of transcription by RNA polymerase II, also in metabolic processes such as macromolecule biosynthetic process and cellular protein modification process. For MF, these miRNAs might affect ubiquitin-like protein transferase activity, RNA polymerase II transferase activity, DNA-binding process, and catalytic activity especially in MAP kinase activity. For CC results, these miRNAs might regulate synthesis of nuclear chromatin and extracellular space components. Pathway enrichment analysis implied these miRNAs were correlated with differentiation, cell motility, immunology, cell proliferation, and migration. Of note, the pathway enrichment analysis showed target genes were enriched in pathways of renal cell carcinoma (KEGG:05211), TGF-beta signaling pathway (KEGG:04350), and ErbB signaling pathway (KEGG:04012). The renal cell carcinoma pathway includes HIF-α pathway and strongly correlates with VEGF and PDGF production. This is consistent with the above GSEA results on gene expression data of TCGA and reminds the effects of these signature relative miRNAs.

To further explore the biological relevance of the 5-CpG signature, we selected one miRNA (miR-1284) for functional experiments. MiR-1284 has been reported to have tumor-inhibiting roles in lung, gastric and ovarian cancers,³⁰⁻³² but its roles in GBMs were still unknown. The in vitro experiments confirmed the anti-tumor role of miR-1284 as inhibiting glioma cell proliferation and migration and inducing glioma cell apoptosis. Interestingly, when treated with demethylation agent, the expression of miR-1284 was significantly decreased, which indicated a positive correlation between its DNA methylation and expression. These results supported the biological implications of the 5-CpG signature that higher methylation status of miR-1284 was positively correlated with patient survival. In general, hypermethylation of promoter inhibits transcription processes and decreases miRNA expression. However, there have been emerging evidences reporting that DNA methylation could facilitate the expression in some situations despite not knowing the exact mechanisms.³³,³⁴

Among the other panel-associated miRNAs, some have been reported to be implicated in glioma biology. MiR-132 was upregulated in GBMs and was a potential indicator of poor prognosis.³⁵,³⁶ MiR-127 and miR-433 are both derived from an overlapping gene locus and colocalized within a cancer-associated genomic region.³⁷ MiR-127 was reported to promote GBM cell migration and invasion

by targeting tumor-suppresser gene SEPT7.³⁸ MiR-433 was reported to be commonly dysregulated in GBMs and suppressed glioma cell proliferation, migration, invasion, and enhanced sensitivity to TMZ therapy.³⁹,⁴⁰ Regarding miR-759 and miR-1248, no biological or clinical evidences have been reported in cancers so far.

Our study has several limitations. First, as a retrospective study, the identification and validation of the signature were based on open source databases which had already been uploaded before. The follow-up information of these researches could not be considered in our study. Also, clinical information such as drug data and recurrent therapy of some cases was not detailed enough, which made it hard to make more subtle analysis. Second, bias caused by the differences among these selected trial platforms should be considered even with compensatory statistical measure. More proof should be collected before conducting further trials. Third, we only performed in vitro study on one miRNA, more in vitro and in vivo studies are needed, especially those on GBM angiogenesis.

In conclusion, by analyzing genome-wide DNA methylation microarray data of miRNAs-associated CpGs, we presented the initial report on the prognostic relevance of aberrant DNA methylation in miRNA regions in GBMs. The identification of the biologically and clinically relevant miRNA methylation signature may represent a promising approach for optimizing prognostication of GBMs, and be of potential value for improving individualized treatment and anti-angiogenic therapy in particular.

## ACKNOWLEDGMENTS

This work was supported by grants from the National Natural Science Foundation of China (No: 81471266; 81671302; and 81672909).

## CONFLICT OF INTEREST

No potential conflicts of interest were disclosed.

## ORCID

En-Ming Kang ID https://orcid.org/0000-0003-1281-3458

An-An Yin ID https://orcid.org/0000-0002-1888-5551

## REFERENCES

1. Omuro A, DeAngelis LM. Glioblastoma and other malignant gliomas: a clinical review. JAMA. 2013;310(17):1842-1850.

2. Gilbert MR, Wang M, Aldape KD, et al. Dose-dense temozolomide for newly diagnosed glioblastoma: a randomized phase III clinical trial. J Clin Oncol. 2013;31(32):4085-4091.

3. Thomas AA, Brennan CW, DeAngelis LM, Omuro AM. Emerging therapies for glioblastoma. JAMA Neurol. 2014;71(11):1437-1444.

4. Dawson MA, Kouzarides T. Cancer epigenetics: from mechanism to therapy. Cell. 2012;150(1):12-27.

5. Baylin SB, Jones PA. A decade of exploring the cancer epigenome - biological and translational implications. Nat Rev Cancer. 2011;11(10):726-734.
6. Baer C, Claus R, Plass C. Genome-wide epigenetic regulation of miRNAs in cancer. *Cancer Res.* 2013;73(2):473-477.

7. Lujambio A, Calin Ga, Villanueva A, et al. A microRNA DNA methylation signature for human cancer metastasis. *Proc Natl Acad Sci USA.* 2008;105(36):13556-13561.

8. Malta TM, de Souza CF, Sabedot TS, et al. Glioma CpG Island Methylator Phenotype (G-CIMP): Biological and Clinical Implications. *Neuro Oncol.* 2017; 20 (5): 608-620.

9. Yin AA, Lu N, Etcheverry A, et al. A novel prognostic six-CpG signature in glioblastomas. *CNS Neurosci Ther.* 2018;24:167-177.

10. Brennan C, Verhaak R, McKenna A, et al. The somatic genomic landscape of glioblastoma. *Cell.* 2013;155(2):462-477.

11. Kurscheid S, Bady P, Sciuscio D, et al. Chromosome 7 gain and DNA hypermethylation at the HOXA10 locus are associated with expression of a stem cell related HOX-signature in glioblastoma. *Genome Biol.* 2015;16:16.

12. Ceccarelli M, F, Malta T, et al. Molecular profiling reveals biologically discrete subsets and pathwayBarthelys of progression in diffuse glioma. *Cell.* 2016;164(3):550-563.

13. van den Bent MJ, Erdem-Eraslan L, Idbaih A, et al. MGMT-STP27 methylation status as predictive marker for response to PCV in anaplastic Oligodendrogliomas and Oligoastrocytomas. A report from EORTC study 26951. *Clin Cancer Res.* 2013;19(19):5513-5522.

14. Noushmehr H, Weisenberger DJ, Diefes K, et al. Identification of a CpG island methylator phenotype that defines a distinct subgroup of glioma. *Cancer Cell.* 2010;17(5):510-522.

15. Bady P, Sciuscio D, Diserens A-C, et al. MGMT methylation analysis of glioblastoma on the Infinin methylation BeadChip identifies two distinct CpG regions associated with gene silencing and outcome, yielding a prediction model for comparisons across datasets, tumor grades, and CIMP-status. *Acta Neuropathol.* 2012;124(4):547-560.

16. Johnson WE, Li C, Rabinovic A. Adjusting batch effects in microarray expression data using empirical Bayes methods. *Biostatistics.* 2007;8(1):118-127.

17. Du P, Zhang X, Huang C-C, et al. Comparison of Beta-value and M-value methods for quantifying methylation levels by microarray analysis. *BMC Bioinformatics.* 2010;11:587.

18. Alizadeh AA, Gentles AJ, Alencar AJ, et al. Prediction of survival in diffuse large B-cell lymphoma based on the expression of 2 genes reflecting tumor and microenvironment. *Blood.* 2011;118(5):1350-1358.

19. Lossos IS, Czerwinski DK, Alizadeh AA, et al. Prediction of survival in diffuse large-B-cell lymphoma based on the expression of six genes. *N Engl J Med.* 2004;350(18):1828-1837.

20. Subramanian A, Tamayo P, Mootha VK, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. *Proc Natl Acad Sci USA.* 2005;102(43): 15545-15550.

21. Peterson SM, Thompson JA, Ufkin ML, Sathyanarayana P, Liaw L, Congdon CB. Common features of microRNA target prediction tools. *Frontiers in genetics.* 2014;5:23.

22. Mi H, Muruganujan A, Ebert D, Huang X, Thomas PD. PANTHER version 14: more genomes, a new PANTHER GO-slim and improvements in enrichment analysis tools. *Nucleic Acids Res.* 2019;47(D1): D419-D426.

23. Bindea G, Galon J, Mlecnik B. CluePedia Cytoscape plugin: pathway insights using integrated experimental and in silico data. *Bioinformatics.* 2013;29(5):661-663.

24. Bindea G, Mlecnik B, Hackl H, et al. ClueGO: A Cytoscape plug-in to decipher functionally grouped gene ontology and pathway annotation networks. *Bioinformatics.* 2009;25(8):1091-1093.

25. Joo JK, Kim SH, Kim HG, et al. CpG methylation of transcription factor 4 in gastric carcinoma. *Ann Surg Oncol.* 2010;17(12):3344-3353.

26. Rodriguez LG, Wu X, Guan JL. Wound-healing assay. *Methods Mol Biol.* 2005;294:23-29.

27. Weathers SP, de Groot J. Resistance to antiangiogenic therapy. *Curr Neurol Neurosci Rep.* 2014;14(5):443.

28. Field KM, Jordan JT, Wen PY, Rosenthal MA, Reardon DA. Bevacizumab and glioblastoma: scientific review, newly reported updates, and ongoing controversies. *Cancer.* 2015;121(7):997-1007.

29. Diaz RJ, Ali S, Qadir MG, De La Fuente MI, Ivan ME, Komotar RJ. The role of bevacizumab in the treatment of glioblastoma. *J Neurooncol.* 2017;133(3):455-467.

30. Li J, Jin H, Yu H, Wang B, Tang J. miRNA1284 inhibits cell growth and induces apoptosis of lung cancer cells. *Mol Med Rep.* 2017;16(3):3049-3054.

31. Cao W, Wei W, Zhan Z, Xie Y, Xiao Q. MiR-1284 modulates multidrug resistance of gastric cancer cells by targeting EIF4A1. *Oncol Rep.* 2016;35(5):2583-2591.

32. Pan C, Wang D, Zhang Y, Yu W. MicroRNA-1284 inhibits cell viability and induces apoptosis of ovarian cancer cell line OVCAR3. *Oncol Res.* 2016;24(6):429-435.

33. Tilley SK, Kim WY, Fry RC. Analysis of bladder cancer tumor CpG methylation and gene expression within The Cancer Genome Atlas identifies GRIA1 as a prognostic biomarker for basal-like bladder cancer. *Am J Cancer Res.* 2017;7(9):1850-1862.

34. Yang X, Han H, De Carvalho DD, Lay FD, Jones PA, Liang G. Gene body methylation can alter gene expression and is a therapeutic target in cancer. *Cancer Cell.* 2014;26(4):577-590.

35. Liu Q, Liao F, Wu H, et al. Upregulation of miR-132 expression in glioma and its clinical significance. *Tumour Biol.* 2014;35(12):12299-12304.

36. Parker NR, Correia N, Crossley B, Buckland ME, Howell VM, Wheeler HR. Correlation of MicroRNA 132 up-regulation with an unfavorable clinical outcome in patients with primary glioblastoma multiforme treated with radiotherapy plus concomitant and adjuvant temozolomide chemotherapy. *Transl Oncol.* 2013;6(6):742-748.

37. Song G, Wang L. A conserved gene structure and expression regulation of miR-433 and miR-127 in mammals. *PLoS One.* 2009;4(11):e7829.

38. Jiang H, Hua D, Zhang J, et al. MicroRNA-127-3p promotes glioblastoma cell migration and invasion by targeting the tumor-suppressor gene SEPT7. *Oncol Rep.* 2014;31(5):2261-2269.

39. Hua D, Mo F, Ding D, et al. A catalogue of glioblastoma and brain MicroRNAs identified by deep sequencing. *OMICS.* 2012;16(12):690-699.

40. Sun S, Wang X, Xu X, et al. MiR-433-3p suppresses cell growth and enhances chemosensitivity by targeting CREB in human glioma. *Oncotarget.* 2017;8(3):5057-5068.

## SUPPORTING INFORMATION

Additional supporting information may be found online in the Supporting Information section at the end of the article.

**How to cite this article:** Kang E-M, Yin A-A, He Y-L, et al. A five-CpG signature of microRNA methylation in non-G-CIMP glioblastoma. *CNS Neurosci Ther.* 2019;25:937-950. https://doi.org/10.1111/cns.13133