Development and validation of a molecular classifier of meningiomas

Abstract

Background

Meningioma, the most common primary intracranial tumor,¹ is associated with notoriously heterogeneous clinical behavior. While the majority exhibit benign growth patterns and can be cured with surgical resection in the setting of growth or symptoms, a critical subset of patients suffer from aggressive disease resulting in significant morbidity and mortality, both from the disease itself and the sequelae of surgery. Importantly, World Health Organization (WHO) grade, which is mainly based on histopathology, has been shown to insufficiently capture the biological and clinical heterogeneity of meningiomas,² prompting the recognition of a need for molecular subgroups which are able to better resolve disease heterogeneity and predict clinical outcomes.^3–7

We previously described 4 consensus molecular groups of meningioma by integrating DNA methylation, transcriptomic, and copy number alterations (CNAs) derived from whole-exome sequencing data³: immunogenic (MG1), benign NF2-wildtype (MG2), hypermetabolic (MG3), and proliferative (MG4). Others have described similar groups in parallel with varying nomenclature and using various molecular platforms.^4–9 Importantly, these molecular groups more accurately predict clinical outcomes, provide an understanding of the underlying biological processes that drive clinical behavior, and yield group-specific therapeutic vulnerabilities. Despite these important advances, the operationalization of these groups into clinical practice is lacking, as their discovery was based on high-dimensional clustering approaches which are inherently stochastic and do not provide a clear avenue through which to classify prospective cases. Additionally, the need for multiple genomics platforms in group discovery is prohibitively expensive for routine clinical use.

To address this critical gap and translate the biologic insights of these novel molecular groups into clinical practice, we focused on adapting a multiplatform description of molecular classification into a more convenient single platform approach by building a classifier that predicts meningioma molecular group classification using only raw DNA methylation files as input. We used DNA methylation to anchor this classifier since CNAs can be inferred from these data, and because methylation-based classifiers have gained considerable traction in the clinical care of brain tumor patients.^2,10,11 We show that independent, prospective cases can be reliably classified into 1 of 4 molecular groups, with biology and clinical outcomes that closely resemble those described in their original discovery. Overall, this model, which we make publicly available, could easily be integrated into existing DNA methylation-based workflows that are already in place and will add considerable clinically-actionable granularity, thereby allowing for the possibility of robust molecularly-informed clinical trials for meningioma.

Materials and Methods

Study Cohort

The cohort for this study was generated by using a combination of previously published cases analyzed by our group, datasets made publicly available, and additional cases we obtained for this new analysis. Our total study cohort was divided into a discovery cohort and a validation cohort. The discovery cohort included a previously published cohort used to describe the original 4 molecular groups.³ The molecular group membership labels of the discovery cohort, generated using DNA methylation, RNA sequencing, and exome sequencing as previously published, were considered the ground-truth. The validation cohort consisted of samples from multiple institutions including Indiana University (n = 305), Case Western/University Hospitals (n = 49), Northwestern University (n = 37), Fred Hutchison Cancer Center (n = 173), Vanderbilt University (n = 47), Vancouver General Hospital (n = 14), the University Health Network, Toronto (n = 210), Baylor (n = 110, previously published⁷), University of California San Francisco (n = 540, previously published⁵), and prospectively-collected samples from the RTOG-0539 clinical trial (n = 92). This study was approved by the University Health Network Institutional Review Board (#18-5820).

Sample Processing

Sample acquisition and processing have been described in detail previously.³ Our combined cohort consisted of a combination of fresh frozen tissue, which was immediately frozen at the time of surgery and stored at −80 °C, and formalin-fixed embedded tissue. DNA was extracted from all samples using the DNeasy Blood and Tissue Kit or QIAamp DNA FFPE Tissue Kit (Qiagen). A representative subset of samples with available tissue had RNA extracted using the RNEasy Kit (Qiagen).

DNA Methylation

Between 250 and 500 ng of extracted DNA was bisulfite converted using the EZ DNA Methylation Kit followed by methylation profiling using the Illumina 850k EPIC array. Raw output (.idat files) were imported, preprocessed, and ssNoob normalized using the minfi¹² package. Filtering was done to remove probes which failed to hybridize (detection P-values > .05), probes overlapping with known single-nucleotide polymorphisms, cross-reactive probes, and probes located on the X or Y chromosome. Output beta-values were used for all downstream analysis.

Inferred Copy Number Alterations

We inferred chromosome arm-level gains and losses using the conumee R package.¹³ Normalized DNA methylation data (beta-values) were aligned to the human genome, binned, and segmented by chromosome arm. Median segment values were computed, and a threshold of |0.2| was used to identify chromosome arm gains and losses.

RNA Sequencing

Extracted RNA was paired-end sequenced on the Illumina HiSeq 2500 to target 70 million reads per sample. Output fastq files were aligned to the GRCh38 human reference genome using the STAR aligner.¹⁴ Aligned reads were sorted and deduplicated using SamTools¹⁵ and subsequently converted to raw counts using the Rsubread¹⁶ package with default parameters. Finally, the edgeR¹⁷ package was used to convert raw counts to log-transformed counts per million (CPM) and apply trimmed mean of M (TMM) normalization. The resultant processed data was used for all downstream analysis.

Mutational Data

All samples in the discovery cohort, as well as prospective samples from the prospective RTOG-0539 clinical trial, underwent whole-exome sequencing. Mutational data from a publicly available cohort⁷ within our combined validation cohort was also used with mutation calls as per their prior publication. As previously described,³ exome libraries were prepared using 100 ng of DNA from the tumor or matched normal DNA (plasma) prior to pair-ended sequencing on a HiSeq 2500 platform. Fastq files were aligned to the hg19 genome using BWA-MEM¹⁸ (v.0.7.12). Quality control was performed using Picard (v1.72) and GATK¹⁹ (v3.6.0). Mutect²⁰ (v1.1.7), and Strelka²¹ (v1.0.13) were used to detect somatic mutations, and known canonical group-specific mutations (NF2, TRAF7, KLF4, AKT1, and POLR2A) were retained for statistical analysis.

In addition, amplification of the TERT promoter region containing the hotspots C228T and C250T was performed on 712 meningioma tumor samples in the validation cohort. The Platinum SuperFi II PCR Master Mix (Thermo Fisher cat no. 12368010) with primers 5’-AGTGGATTCGCGGGCACAGA-3’ and 5’- CAGCGCTGCCTGAAACTC-3’²² was used to produce a 235 base pair amplicon. PCR products were separated by gel electrophoresis to verify a product at the correct size. After purification with the ZR-96 DNA clean-up kit (Zymo Research, cat no. D4018), PCR products were sent for Sanger sequencing at The Centre for Applied Genomics. Chromatograms were analyzed in Geneious Prime, and all mutations were confirmed by sequencing a second PCR reaction.

Model Training

To predict molecular groups, we used elastic net classifiers which were optimized using 10-fold cross validation repeated 3 times on the discovery cohort, where the true label was defined as the original molecular group assignment. For each sample in this cohort, a predicted molecular group (MG) was generated by training models on all remaining samples (n–1), such that no data from the sample being predicted was used in the feature selection or model training process. For each model, input features were centered and scaled and alpha and lambda were allowed to vary from 0–1 to 0–0.1, respectively. Once a prediction was generated for each sample, model configurations (the number of input features) were compared using the area under the receiver operating characteristic curves.

For feature selection, differential methylation analysis was performed in a pairwise fashion between each input group label using the limma package,²³ and the most differentially methylated regions ranked by moderated t-statistic. The number of hyper- and hypo-methylated probes used in the model was determined based on the leave-one-out analysis above and only probes with adjusted P < .05 were considered (10, 50, 100, 250, 500, 1000, and 2500 hyper- and hypo-methylated probes for each pairwise comparison were tested for all models, and the configuration with the highest area under the curve (AUC) chosen for the final model). Probes that were common to the Illumina EPIC 850k and the EPICv2 array were considered as possible features, such that our final model is fully compatible with current versions of the EPICArray.

Model Validation

Once the model configuration was determined, an elastic net classifier (a form of regularized logistic regression) was trained using the full discovery cohort and used to predict the molecular group on our combined validation cohort. The training parameters were identical to those used in the model optimization. The biological profiles of these predicted molecular groups were compared to the molecular descriptions of the original cohort using methylation, clinical, and orthogonal genomic data (RNA sequencing and somatic point mutation data from whole-exome sequencing). Inferred copy number changes were also compared between cohorts, with a specific focus on canonical MG-specific gains/losses such as loss of chromosome 1p, gain of chromosome 1q, and loss of chromosome 22q.

RNA Sequencing Data

A subset of the validation cohort with available RNA sequencing data was used as further validation of the biological recapitulation of our original molecular groups. Single sample gene set variation analysis (ssGSVA),²⁴ a method for computing the expression of a gene set within individual samples, was used to compute the relative expression of each MG-specific gene signature derived from our original cohort for each sample. Expression of each signature was compared between predicted group assignments in the validation cohort using an analysis of variance (ANOVA) with significance set at P = .05.

RNA Pathway Analysis

Within the validation cohort, differential gene expression analysis was performed between each predicted molecular group and the remaining samples (immunogenic versus others, benign NF2-wildtype versus others, hypermetabolic versus others, and proliferative versus others) using the limma package.²³ Ranked gene lists were used to generate pathway enrichment scores as previously described.³ Network maps were generated with nodes representing pathways with associated P-value < .001, and separate pathways were connected when gene overlaps were associated with Jaccard coefficient of >0.25. Nodes were finally grouped based on shared keywords to allow high-level pathway enrichment annotation and plotted using Cytoscape.

Cell Deconvolution Analysis

As further validation of our predicted molecular groups, we compared the inferred cell proportions in our discovery and validation cohorts using CIBERSORTx.²⁵ We applied the same single-cell RNA reference matrix defining neoplastic cells, macrophages, fibroblasts, T-cells, and endothelial cells as was used for the cell deconvolution in the original discovery of the groups.³ All samples with available bulk RNA data were then deconvolved to infer a relative proportion of each cell type. Notably, inputs for the deconvolution were raw counts, since this was the input used to generate the single-cell RNA signature matrix in the discovery cohort.

Statistical Analysis

For continuous variables, a student’s t-test (2 variables) or ANOVA (>2 variables) was used for statistical comparisons. Chi-squared or Fisher’s exact tests were used for comparisons of categorical variables. Outcomes were compared using the log-rank test. For all comparisons, alpha is set to 0.05 unless otherwise specified. All analysis was performed using R²⁶ version 4.3.1.

Results

Study Cohort

Full details of our discovery (n = 121) and validation (n = 1577) cohorts are presented in Table 1. The median age of our validation cohort was 58 years (interquartile range [IQR] 48–68 years) and 1023/1527 (67%) of patients with annotated sex were female. There were 974 WHO grade 1 tumors, 487 WHO grade 2 tumors, and 113 WHO grade 3 tumors of those with annotations available. Our validation cohort was enriched in WHO grade 3 cases compared to the discovery cohort (P < .001) and in cases that underwent adjuvant radiotherapy after surgery (P = .03, Table 1). There were otherwise no statistical differences between the demographic parameters of our cohorts. Of the 1514 patients with comprehensive follow-up data, 428 (28.3%) had recurrence during follow-up after a median of 2.3 (IQR 0.8–3.9) years from surgery, reflective of a cohort enriched in aggressive tumors (which are particularly in need of molecularly-informed care). Median follow-up among cases that did not recur during follow-up was 2.9 years (IQR 0.2–5.7 years), and median follow-up overall, censored at the time of recurrence, was 2.6 years (IQR 0.2–5.0 years). Notably, our validation cohort includes newly generated data from the only prospective clinical trial samples in meningioma (RTOG-0539).

Model Construction

Given the previous finding that hypermetabolic and proliferative tumors are more strongly distinguished by RNA sequencing and mutational burden than DNA methylation, reinforced by the fact that previous work using DNA methylation alone to identify subgroups has yielded three groups,⁵ one of which can subsequently be further clustered into 2 groups which parallel the hypermetabolic and proliferative groups,²⁷ we constructed a 2-layer model to predict molecular group (Figure 1A). In the first layer (model 1), hypermetabolic and proliferative tumors were considered as one group, yielding three possible outputs with associated probabilities. The second layer (model 2) was trained to specifically distinguish hypermetabolic and proliferative tumors. For model 1, the optimal configuration used the top 250 hyper- and hypo-methylated probes for each pairwise comparison as features (a total of 1383 unique probes) and yielded AUCs of 0.99 (immunogenic versus others), 0.96 (benign NF2-wildtype versus others), and 0.97 (hypermetabolic/proliferative versus others) in a leave-one-out analysis (discovery cohort, probability threshold 0.5; Supplementary Figure 1). Notably, cases that were discordant between the predicted and actual molecular group (n = 10) had similar group-specific outcomes, suggesting these cases may not be fully represented by a single molecular group (Supplementary Figure 2). In model 2, the top 100 hyper- and hypo-methylated probes were included (total of 200 unique probes), with an associated AUC of 0.91 in a leave-one-out analysis (discovery cohort, threshold 0.8; Supplementary Figure 3). The discordant cases (n = 10) also had similar group-specific outcomes (log-rank P = .2, Supplementary Figure 4).

Figure 1.

Predicted molecular groups recapitulate the methylation signatures, inferred copy number profiles, and outcomes of the original molecular groups. (A) Model diagram, depicting information flow from model 1 to model 2 with additional CNA-based refinement. (B) Distribution of methylation model probability scores in both models based on final predictions (excluding cases that were reassigned based on CNV profile). (C) Barplots depicting the proportion of samples in each molecular group with canonical CNAs among the validation cohort (predicted, top) and the discovery cohort (original, bottom). (D) Group-specific Kaplan Meyer survival curves are plotted for the validation cohort (left), the prospective RTOG-0539 cohort (a subset of the validation cohort, middle), and the discovery cohort (right). Vertical lines represent median progression-free survival.

Open in new tabDownload slide

To apply these trained models to our validation cohort, processed beta-values were input into model 1, and preliminary molecular groups were assigned based on the highest output probability. Cases assigned as “hypermetabolic or proliferative” were then input into model 2 to separate them into 1 of these 2 groups. Therefore, if the output of model 1 was either immunogenic or benign NF2-wildtype, the prediction probability was simply the corresponding output probability from model 1; if the output of model 1 was “hypermetabolic or proliferative,” the prediction probability was the product of the probability in model 1 and the highest probability in model 2 (corresponding to its associated group prediction). Cases with a final output probability <0.5 were not assigned to a single molecular group to ensure confidence in classification. Cases that were reassigned based on copy number variations (CNVs) profiles are not assigned a final probability for clinical purposes, but otherwise, the output probability is provided for assessment of the degree of confidence of a particular prediction.

Model Refinement With Copy Number Alterations

Given that stereotypical CNAs were specific to some molecular groups in their original description, we refined the methylation model output using CNA data inferred from DNA methylation (effectively creating a 3-layer model using only DNA methylation as input data). Based on findings from our discovery cohort, cases with the presence of both 1p loss and either 1q gain or 10 (p or q) loss, while rare, were assigned as proliferative. Tumors labeled as hypermetabolic in the methylation model with the presence of benign NF2-wildtype-specific polysomies identified on original discovery (gain of chromosome 5, 12, 13, 17, or 20) and the absence of chromosome 1p loss were assigned benign NF2-wildtype, and those with either a completely neutral copy number profile or isolated loss of chromosome 22q were not assigned to an individual molecular group (n = 134 cases).

The Methylation-Only Model Faithfully Reproduces the MG-Specific Methylation and CNA Patterns of Original MGs

We applied the complete 3-layer classifier to all samples in the validation cohort, classifying 1378/1577 (87.4%) of cases. Of these 1378 cases, 68 (4.9%) cases had final assignments that were refined by their inferred CNV profiles. Two hundred and seventy-three cases (19.8%) were classified as immunogenic, 550 (39.9%) as benign NF2-intact, 272 (19.7%) as hypermetabolic, and 283 (20.5%) as proliferative. Notably, among the prospectively-collected samples from the RTOG-0539 clinical trial which classified as 1 of the 4 molecular groups with sufficient confidence, 14 (15.7%) were classified as immunogenic, 29 (32.6%) as benign NF2-wildtype, 23 (25.8%) as hypermetabolic, and 23 (25.8%) as proliferative. Model output scores were highly specific to each molecular group in both layers of the methylation model (Figure 1B).

Given the observation that specific chromosome arms are stereotypically gained or lost in an MG-specific manner, we also compared inferred copy number changes of our final MG predictions in the validation cohort to those described in the original discovery cohort (Figure 1C). Loss of chromosome 1p was found in 1.1%, 7.1%, 50%, and 88% of immunogenic, benign NF2-wildtype, hypermetabolic, and proliferative tumors in the validation cohort (P < .001), respectively. This is concordant with the proportions inferred from the discovery cohort that also showed 1p loss was enriched in hypermetabolic and proliferative tumors (0%, 18.8%, 62.8%, and 82.7%, respectively, P < .001). Gain of chromosome 1q was highly specific for proliferative meningiomas and found in 14.1% of predicted proliferative tumors, no immunogenic tumors, and 0.3% of both benign NF2-wildtype and hypermetabolic tumors in the validation cohort (P < .001). Similarly, gain of 1q was present in no immunogenic tumors, 3.1% of benign NF2-wildtype tumors, 2.3% of hypermetabolic tumors, and 31.0% of proliferative tumors in the discovery cohort (P < .001). Loss of chromosome 22q was present in 65.2% versus 64.7% of immunogenic tumors, 5.6% versus 12.5% of benign NF2-wildtype tumors, 77.6% versus 69.8% of hypermetabolic tumors, and 89.0% versus 75.9% of proliferative tumors in the validation and discovery cohorts, respectively. Finally, CDKN2A hetero/homozygous deletions were highly enriched in the proliferative group (5.7/12.3% of proliferative cases versus 0%/0% of immunogenic cases, 0.2/0.4% of benign NF2-wildtype cases, and 1.5/0% of hypermetabolic cases in the validation cohort [P < .001]). In the discovery cohort, they were present in no immunogenic cases, 6.2/3.1% of benign NF2-wildtype cases, 2.3/2.3% of hypermetabolic cases, and 6.9%/10.3% of proliferative cases. Overall, we demonstrate that canonical MG-specific copy number changes are upheld in our predicted molecular groups.

Most Undifferentiated Cases are Hypermetabolic-Like With Bland Copy Number Profiles

Of the 199 “undifferentiated” cases which were not confidently assigned to a single group, 65 (32.7%) had final output model probabilities <0.5 and 134 (67.3%) were assigned as hypermetabolic in the 2-layer model but had neutral copy number profiles (or isolated loss of 22q) and were therefore not assigned to a particular group. Undifferentiated cases were slightly more likely to be male (P = .01), but there were no other clear differences in demographics (Supplementary Table 1). Taking the highest output probability to represent the best possible group assignment, 9 undifferentiated cases (4.5%) could be considered immunogenic-like, 16 (8.0%) as benign NF2-like, 162 (81.4%) as hypermetabolic-like, and 12 (6.0%) as proliferative-like. Interestingly, outcomes based on these group labels did not follow the same trends as the cases that were confidently assigned, suggesting that these tumors are not fully representable by a single group. In particular, the hypermetabolic-like undifferentiated tumors, making up the majority of undifferentiated cases, had considerably longer PFS than expected (median 10.28, 95% CI 8.46–Inf years), suggesting the presence of a benign hypermetabolic-like, copy number neutral subgroup which may not have been captured on smaller discovery cohorts to date (Supplementary Figure 5A). Additionally, unsupervised clustering of the top 10 000 variable probes in the validation cohort demonstrates that the undifferentiated cases tend to cluster between 2 or more groups, suggesting they are not fully encapsulated by a single group using methylation alone (Supplementary Figure 5B). Nevertheless, they do partially retain the MG-specific gene signatures of their “nearest” MG, suggesting that these cases could be classifiable if additional data modalities, such as RNA sequencing and exome sequencing, were used in addition to DNA methylation (Supplementary Figure 5C).

Predicted Molecular Groups Have Outcomes That are Concordant With the Original MGs

We next sought to compare molecular group-specific outcomes in the validation cohort to the discovery cohort (Figure 1D). Importantly, the discovery of the original molecular groups was agnostic to outcome, but each group was associated with distinct survival patterns. For immunogenic tumors, median PFS was not reached in either cohort, and outcomes were overall most favorable. For benign NF2-wildtype tumors, median PFS was also favorable and not reached in the validation cohort; in the discovery cohort, it was 18.5 (95% CI 16.0–Inf) years. Median PFS for hypermetabolic tumors was 7.38 (95% CI 4.89–Inf) years in the validation cohort and 6.47 (95% CI 3.32–Inf) years in the discovery cohort. Finally, proliferative tumors had the shortest median PFS in both cohorts: 2.53 (95% CI 2.31–3.52) years in the validation cohort and 2.23 (95% CI 1.35–6.05) years in the discovery cohort. The distribution of outcomes was similar in the prospective RTOG-0539 samples within the validation cohort (n = 89, 5-year PFS 92.9% in immunogenic, 89.7% in benign NF2-wildtype, 73.9% in hypermetabolic, and 52.2% in proliferative). Importantly, we show that the molecular group remains strongly predictive of outcome even when stratifying by WHO grade, extent of resection, and primary versus recurrent status (Supplementary Figure 6). Additionally, MG-specific outcomes were similar whether or not the final prediction probabilities were >0.8, suggesting confidence in outputs with even for probabilities of 0.5–0.8 (Supplementary Figure 7). While the classification of molecular groups was agnostic to the outcome in the original discovery and the classifier we describe in this paper, we found remarkable concordance between the outcome patterns of the original molecular groups and the newly predicted groups in our validation cohort.

Molecular Group-Specific Mutations are Preserved in Our Predicted Groups

We next sought to further validate the biological fidelity of the predicted molecular groups using orthogonal genomic platforms, which were not used in our classifier’s training or predictions. We first investigated the somatic point mutation profiles of a subset of the validation cohort with available whole-exome sequencing (n = 165, Figure 2A). This includes newly generated whole-exome sequencing from the prospective RTOG-0539 clinical trial samples. Notably, a total of 176 cases had available exome sequencing data including the undifferentiated cases, of which 109 were obtained from a publicly available cohort.⁷ Importantly, we found that none of the cases predicted to be benign NF2-wildtype had mutations in NF2, reinforcing the appropriate naming of this group of tumors. Conversely, known canonical non-NF2 meningioma driver mutations (KLF4, TRAF7, AKT1, and POLR2A), were demonstrated to be specific to this group as originally described. Specifically, we found at least one of these mutations in 40% of benign NF2-wildtype tumors in the discovery cohort (compared to 0%, 2.4%, and 3.6% of immunogenic, hypermetabolic, and proliferative tumors, respectively) and in 51.4% of predicted benign NF2-wildtype tumors in the validation cohort (compared with 2.7%, 0%, and 8.0% of immunogenic, hypermetabolic, and proliferative tumors, respectively). Finally, TERT promoter mutations were not present in any immunogenic cases, 1/296 (0.3%) of benign NF2-wildtype cases, 4/149 (2.7%) of hypermetabolic cases, and 8/147 (5.4%) of proliferative cases in the validation cohort. This is compared to 1/17 (5.8%), 2/32 (6.2%), 5/43 (11.6%), and 3/29 (10.3%) in the discovery cohort, respectively, though the numbers are small in this cohort. Overall, we show that the mutation patterns of our predicted molecular groups follow the same trends as those identified in their original discovery despite the current classifier being agnostic to somatic point mutation data and whole-exome sequencing data in its entirety. This further validates the biologic robustness of our methylation-only model.

Figure 2.

Validation of predicted molecular groups using orthogonal genomic platforms. (A) Whole-exome sequencing reveals similar molecular group-specific patterns of mutations between cohorts, with benign NF2-intact (MG2) cases harboring the majority of canonical non-NF2 mutations TRAF7, KLF4, AKT1, and POLR2A, but comparatively few NF2 mutations in both cohorts. (B) Pathway analysis reveals expected transcriptomic signatures in each molecular group, with enrichment of immune pathways in immunogenic (MG1, and, to a lesser extent, hypermetabolic, MG3, cases), metabolic signaling pathways in hypermetabolic (MG3) cases, and cell cycle pathways in proliferative (MG4) cases. In this network, nodes represent pathways and edges represent shared genes between pathways. Nodes colored in red represent upregulated pathways and those colored in blue represent downregulated pathways. Boxplots below depict the expression of each group-specific gene signature identified on the discovery cohort (ANOVA P-values plotted for each). (C) Bulk deconvolution analysis demonstrates consistent MG-specific cell populations between discovery and validation cohorts. Notably, immunogenic tumors are associated with greater macrophage infiltration while hypermetabolic and proliferative tumors are associated with a higher proportion of neoplastic cells.

Open in new tabDownload slide

Predicted Molecular Groups Recapitulate the Transcriptomic Signatures of Original Molecular Groups

Given that the transcriptional profiles of the originally described consensus molecular groups were highly distinct, we used RNA sequencing to further characterize the biology of predicted molecular groups in the validation cohort (n = 557 classified cases including newly generated data from prospective clinical trial samples, Figure 2B). Notably, a total of 631 cases had available RNA sequencing data including the undifferentiated cases, of which 255 were obtained from publicly available cohorts.^5,7 Using pathway analysis, we showed that predicted molecular groups were associated with upregulation of their namesake processes, such as immune-related pathways in immunogenic tumors, metabolic/mitochondrial pathways in hypermetabolic tumors, and cell cycle pathways in proliferative tumors. Notably, the hypermetabolic group shares a subset of immune pathways with the immunogenic group, which aligns with our previous finding that these tumors are moderately immune-enriched when compared to the benign NF2-wildtype and proliferative groups.

As further confirmation of the transcriptome-level validity of our predicted molecular groups, we computed the expression of group-specific gene signatures generated during the initial discovery of these groups³ for each sample in the validation cohort. The expression of each signature was computed using GSVA and compared between predicted molecular groups using ANOVA. Notably, we found that each predicted group was associated with the highest median expression of its corresponding gene signature, and this difference in expression was highly significant in all groups (P < .0001, Figure 2B). Therefore, we demonstrate that the RNA profiles of our predicted molecular groups recapitulate the signatures from the original groups despite these models being agnostic to RNA data.

Finally, we inferred the relative proportions of neoplastic cells, macrophages, T-cells, fibroblasts, and endothelial cells by deconvoluting the transcriptome of samples from the validation cohort. To do this, we applied a cell-type-specific single-cell signature matrix that was previously generated using the original discovery cohort to samples from the validation cohort using CIBERSORTx (Figure 2C). We found that the proportion of neoplastic cells increased with the molecular group in both cohorts (median 45.4%/47.0% [discovery/validation] in immunogenic tumors, 73.2%/60.0% in benign NF2-wildtype tumors, 80.3%/72.8% in hypermetabolic tumors, and 90.5%/86.6% in proliferative tumors. Fibroblasts were comparatively enriched in immunogenic [P < .001 discovery and validation] and benign NF2-wildtype tumors [P < .001 discovery and validation]. Meningiomas predicted to be immunogenic had the highest proportion of macrophages in both the discovery [median 28.4%, P < .001]), and validation (median 22.1%, P < .001) cohorts, again reinforcing their defining immune-enriched phenotype.

Discussion

The inclusion of TERT promoter mutations and CDKN2A deletions to classify meningioma represents a significant update to the longstanding histopathology-based WHO criteria.²⁸ The molecular taxonomies of meningioma have been expanded with multiple classifications described in recent years, each consistently resolving the biological and clinical outcome heterogeneity of the disease. We previously demonstrated that the combination of DNA methylation, RNA sequencing, and CNAs derived from whole-exome sequencing data can robustly classify meningiomas into 4 consensus molecular groups with defining biological processes, outcomes, and therapeutic vulnerabilities. Notably, while attempts have been made to define immunohistochemical surrogates to these groups to allow for easy implementation into the laboratory setting,^3,29 a 1:1 correspondence with molecular groups has yet to be established and molecular approaches continue to be required. In this follow-up, we have built a publicly available, easy-to-use point-and-click classifier which can predict molecular group with high fidelity using only DNA methylation data as input. We showed that a multi-layered DNA methylation model with additional refinement based on inferred CNAs (Supplementary Figure 8) achieves suitable performance with confident classification and low group-specific error rates (Supplementary Figure 9). Remarkable concordance was achieved when comparing the methylome profiles, CNAs, canonical somatic mutations, RNA signatures, and clinical outcomes between predicted molecular groups in the validation cohort and the ground-truth molecular group assignments identified through the integration of multiple molecular platforms in the discovery cohort, reinforcing the validity of our model.

While methylation profiling is not currently part of the standard of care in meningioma management and is not available at many laboratories, it is becoming increasingly accessible as the cost of genomic profiling decreases, further mitigated by using microarray-based technologies³⁰ such as the EPICArray used here. Further, centralized labs such as our institution are increasingly offering methylation-based services to external centers on a referral basis, allowing widespread access to this valuable technology. Additionally, we are aware that billing codes are being developed for methylation profiling in the United States, which will further reduce the barrier to implementation.

Our institution routinely applies a previously published DNA methylation-based nomogram to stratify meningiomas into high and low risk groups,² thereby helping inform adjuvant treatment decisions in a multidisciplinary fashion. For example, it is used to select patients for adjuvant RT among completely resected WHO grade 2 cases, for whom there remains ongoing clinical equipoise. We (and many other institutions) also routinely apply the DKFZ brain tumor classifier¹¹ to classify central nervous system tumors with unclear histopathology using DNA methylation. Based on this precedent, our novel molecular group predication model could be seamlessly integrated into already functional workflows, allowing for increasingly meaningful output for each patient diagnosed with meningioma. Importantly, since we have recently shown that molecular group is strongly predictive of response to radiotherapy,³¹ the ability to classify prospective cases using DNA methylation alone will have immediate practice clinical impact to patient care. The ready integration into existing workflows will allow for immediate data accrual on a prospective basis across multiple centers, allowing for increasingly rigorous model validation and generalizability over time. Finally, and perhaps most importantly, this will allow for molecularly-informed clinical trials, whereby recruited patients are stratified into homogeneous molecular subgroups in order to address the heterogeneity that may be a source of the lack of positive meningioma trials to date. This could allow for the identification of increasingly personalized, MG-specific treatments in the future.

The results of our study should be considered in the context of some limitations. First, the labels in our training are based on assignments from prior clustering analyses, which are associated with inherent stochasticity. Basing our classifier on a different datatype such as RNA sequencing, furthermore, may have led to differing results; notably, DNA methylation was chosen given its increasing role in tumor classification and its ability to robustly infer copy number profiles. Tumors which are not fully encapsulated by a single molecular group in either cohort, for example, could unpredictably skew results. To mitigate this in our classifier, we have built in confidence scores that are returned to the user and return the result of “undifferentiated” for cases where prediction probabilities are low, suggesting a case that is not fully encapsulated by a single molecular group. Of course, this comes with the caveat that a small minority of patients will not have a robust molecular group assignment, and other approaches such as outcome prediction models and copy number profiling will be needed to help prognosticate those cases. Secondly, most samples collected and included in this study were done so in a retrospective manner and so the limitations of retrospective data collection remain. However, the inclusion of valuable prospectively-collected, uniformly treated samples from the RTOG-0539 trial, and the large size and multi-institutional nature of our validation cohort addresses some of these limitations.

Conclusions

Overall, we present a novel meningioma molecular classifier using only DNA methylation data as input and rigorously demonstrate its validity using a large, multi-centered cohort of 1698 patients including samples from the prospective RTOG-0539 clinical trial. This represents an important step toward increasingly personalized care for patients with meningioma and will help usher in a new era of molecularly-stratified clinical trials.

Funding

This study was funded by the Canadian Institutes of Health Research (CIHR) and Brain Tumour Charity (BTC). This study was supported in part with funding from the National Cancer Institute, National Institutes of Health, under the Cancer Moonshot Initiative (contract number HHSN261201500003I; task order number HHSN26100039). The content of this publication does not necessarily reflect the views or policies of the US Department of Health and Human Services nor does mention of trade names, commercial products or organizations imply endorsement by the US Government. Additional funding for this study was provided by the Canadian Institutes of Health Research (CIHR) Project Fund (RN482811-481519 to G.Z., F.N., and J.Z.W.), Brain Tumor Charity United Kingdom (GN-00043 and GN-000693 to G.Z. and F.N.), the UHN Foundation, Mary Hunter Meningioma Research Fund, the V Foundation for Cancer Research (G.Z.), the CIHR Vanier Scholarship (F.N., J.Z.W., and C.G.), the American Association of Neurological Surgeons Neurosurgery Education & Research Foundation Research Fellowship (F.N. and J.Z.W.), the Congress of Neurological Surgeons Tumor Section (F.N. and J.Z.W.), and the Princess Margaret Hospital Foundation Hold ‘em For Life Oncology Fellowship (F.N., J.Z.W., A.P.L., and C.G.).

Conflict of interest statement

The authors of this manuscript have no conflicts of interest to disclose.

Acknowledgments

We thank the Princess Margaret Genomics Centre (https://pmgenomics.ca), as well as The Centre for Applied Genomics (https://tcag.ca) for their expertise in generating some of the sequencing data used in this study. For the NRG Oncology RTOG-0539 cases, total RNA sequencing was performed by the High Throughput Sequencing Facility at the University of North Carolina School of Medicine and whole-exome sequencing was generated at the Broad Institute. Additionally, we would like to thank key members of NRG Oncology including Stephanie Pugh and Minhee Won, the Cancer Therapy Evaluation Program (CTEP), and Leidos Biomedical Research for their collaborative efforts as part of the Molecular Profiling to Predict Response to Treatment (MP2PRT) program.

Authorship Statement

A.P.L. and J.Z.W. drafted the manuscript; A.P.L., J.L., V.P., and J.S. performed all analysis; A.P.L., J.Z.W., C.G., Z.P., A.A., Y.E., C.W., Q.W., O.S., and S.M. performed all experimental procedures; A.A.C.-G., M.A.Z., G.T., M.T., F.B., J.S.B.-S., A.E.S., S.C., L.B.C., A.D.R., S.M., S.Y., A.M., A.G., and D.S.T. contributed valuable samples and data from their respective institutions; K.A., A.G., F.A., and G.Z. oversaw all aspects of the project; A.P.L., J.Z.W., F.N., and G.Z. conceived the study. All authors critically revised the manuscript and approved its submission.

Data Availability

Our classifier is freely available for public use (https://www.meningiomaconsortium.com/models/). Source code and data are also available for download on Zenodo (https://zenodo.org/records/13236998).

References

Author notes