Mathematical Modeling and Mutational Analysis Reveal Optimal Therapy to Prevent Malignant Transformation in Grade II IDH-Mutant Gliomas

Diffuse gliomas account for approximately 80% of malignant central nervous system tumors and are classified by the World Health Organization (WHO) into grades II–IV, according to their histopathology and clinical behavior (1). WHO grade II diffuse gliomas are usually referred to as low-grade gliomas (LGG), while WHO grade III and IV gliomas are high-grade gliomas. In general, LGGs are characterized by no significant enhancement in post-contrast T1 MRI, unlike WHO grades III and IV (2). LGGs are generally less aggressive tumors with a longer and indolent clinical course but often undergo malignant transformation (MT) and recur as a higher grade with a poor prognosis (1). Historically, in most LGG cases, the “wait and see” approach was often favored because of the lack of symptoms. However, several reports support the early use of surgery to influence the MT-free survival and overall survival (OS; refs. 3–5). Because of these tumors' infiltrative nature, surgery alone cannot eradicate tumors, and other treatments are necessary to control the disease (6). Essentially, radiotherapy (RTx) or RTx plus chemotherapy (CTx; chemoradiotherapy, hereafter called CRTx) with temozolomide or the regimen of procarbazine, lomustine (CCNU), and vincristine (known as PCV) is widely used for the treatment of patients with LGGs (7–12).

The advent of genome-wide analysis technologies has revealed the characteristic genetic lesions in various tumors, including LGGs, and those results imply clinical outcomes (13–15). Mutations in the genes encoding isocitrate dehydrogenase (IDH) 1 and 2 (IDH1 and IDH2) account for approximately 80% of LGGs and represent a fundamental and earliest molecular event in tumorigenesis (16–19). The IDH-mutant LGGs (IDHmut-LGG) are genetically classified into two subtypes depending on the presence of codeletion of the short arm of chromosome 1 (1p) and the long arm of chromosome 19 (19q) (1p/19q): IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel. We previously reported that primary IDH-mutant grade III gliomas have a greater number of gene mutations and broad copy-number variations (CNV) than IDHmut-LGGs (20). Longitudinal studies of IDHmut-LGGs have identified several genetic and epigenetic alterations more common in recurrent tumors with MT than primary tumors (21, 22). While RTx and/or CTx inhibit tumor growth, they induce genomic (or epigenomic) abnormalities that can promote MT (23, 24). Because MT affects the OS of patients with IDHmut-LGGs, therapeutic efforts should focus on preventing MT to improve OS (25). Although a randomized controlled trial revealed the long-term efficacy of early RTx compared with delayed RTx, there have been few studies that have comprehensively examined the relationship between the timing of treatment and MT (26).The association between treatment and MT and therapeutic strategies to prolong the MT-free survival in IDHmut-LGGs remains elusive.

Mathematical and computational approaches have been applied to interpret experimental and clinical data in various cancers, including LGGs and glioblastomas (27–31). Among several interdisciplinary studies, two mathematical models accounted for the MT of LGGs (32, 33). These studies analyzed a small number of patients who underwent surgery alone and estimated growth parameters, leading to the approximate prediction of tumor progression. However, treatment options to prevent MT using adjuvant therapies in IDHmut-LGGs with a large dataset have not yet been investigated.

Here, we designed a mathematical model based on time-series data, including tumor volume and matched clinical information of 276 IDHmut-LGGs classified by the 1p/19q codeletion to investigate tumor progression dynamics, calculate the risk of MT, and suggest optimal therapeutic strategies. We also performed genome-wide mutational analysis, including targeted deep sequencing (n = 111) and whole-exome sequencing (WES) with serial multisampling [n = 100 (45 patients)] to elucidate the relationship between genetic alterations and the risk of MT. This study is the first mathematical analysis to clarify the relationship between treatment and MT in IDHmut-LGGs using a relatively large dataset, enabling the proposal of the optimal treatment strategy for IDHmut-LGGs.

In total, 276 patients diagnosed with supratentorial IDHmut-LGGs at 18 years of age or over, followed by MRI scans at 10 hospitals in Japan between 1990 and 2018 were considered. Patients who underwent supratotal resection with no apparent residual tumors during follow-up were excluded. Tumors were classified into IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel according to the 2016 revision of the WHO Classification of Tumors of the Central Nervous System (1). All tumors demonstrated the absence of tumor contrast enhancement in MRI at the time of diagnosis. The median follow-up time was assessed for individuals with censored data. Written informed consent was obtained from all patients. This study was approved by the ethics committee or Institutional Review Boards of all participating institutes (approval number: 2012-0067).

Tumor volumes were measured on either fluid-attenuated inversion recovery or a T2-weighted image by two board-certified neurosurgeons. Quantitative assessment was performed using Horos, a free, open-source research platform for medical image analysis (https://horosproject.org/). Tumor segments were traced using a semi-automatic region-growing tool and manually corrected on axial sections. Tumor volumes were retrieved in cubic centimeters (cm³). MT was defined as pathologic confirmation of grade III–IV glioma and/or new contrast enhancement with a growth pattern consistent with MT as determined by multidisciplinary consensus.

Our genome-wide mutational analysis dataset included 100 WES for serial multisampling from 45 patients (IDH^mut/1p19q^codel: 38 samples from 17 patients; IDH^mut/1p19q^noncodel: 62 samples from 28 patients) and 111 targeted sequencing (IDH^mut/1p19q^codel: n = 46; IDH^mut/1p19q^noncodel: n = 65), with overlapping data from 28 patients (IDH^mut/1p19q^codel: n = 8; IDH^mut/1p19q^noncodel: n = 20). Sequencing data, including WES from 21 samples (15 patients) with matching germline sequences and all targeted sequencing, were published previously (13). WES for 79 tumor samples and 30 blood samples was performed using targeted capture of all exon sequences, as previously described with minor modifications (13). See Supplementary Materials and Methods for additional details.

In targeted sequencing, we selected 185 genes that included recurrently mutated genes in LGGs and related disorders as described previously (34). In brief, somatic mutation calling was performed using the empirical Bayesian mutation calling method, in which we adopted variants with variant allele frequencies ≥ 0.05 in tumor samples.

For 100 patients, we analyzed SNP array data to assess broad and focal CNVs based on a hidden Markov model using CNAG, as described previously (13, 35). Using WES data, copy-number identification was performed, according to the GATK Best Practice (https://gatk.broadinstitute.org/hc/en-us). See Supplementary Materials and Methods for additional details.

OS was calculated from the time of diagnosis until death or the last follow-up and MT-free survival was calculated from the time of diagnosis until MT or the last follow-up. OS and MT-free survival was evaluated using the log-rank test. The Pearson correlation coefficient between OS and MT-free survival was evaluated on patients with both MT and death. Kaplan–Meier curve analysis for tumor volume at diagnosis was performed by classifying into three categories (<20, 20–59, and ≥ 60 cm³) according to the recursive partitioning and regression-tree analysis. The Wilcoxon rank-sum test (for nonpaired samples) and Wilcoxon signed-rank test (for paired samples) were used to compare diagnostic values between the two groups. Comparisons of frequencies were made using the Fisher exact test. The associations of the cumulative MT risk (cMTrisk) and tumor mutational burden (TMB) were analyzed using the Pearson correlation analysis. Statistical analyses were performed using R version 3.6.1 (https://www.r-project.org/). Additional details can be found in Supplementary Materials and Methods. The Q values (FDR) were adjusted using the Benjamini and Hochberg method. Differences were considered statistically significant at P and Q values < 0.05.

WES, targeted sequencing data and SNP-array data are available in the European Genome-phenome Archive (EGA) under accession EGAS00001001044 and the Japanese Genotype-phenotype Archive (JGA, https://www.ddbj.nig.ac.jp/jga), under accession number JGAS000268.

We utilized a dataset of 276 patients with IDHmut-LGGs, which were classified into IDH^mut/1p19q^codel (n = 123) or IDH^mut/1p19q^noncodel (n = 153; Table 1; Supplementary Table S1). The median age at diagnosis of IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel was 40 (first and third quartiles: 33 and 51) and 36 (29 and 45), respectively, and the median follow-up periods were 1,932 days (first and third quartiles: 1,235 and 2,896) and 1,964 days (1,278 and 3,175), respectively. During the follow-up, 34 patients with IDH^mut/1p19q^codel (28%) and 73 patients with IDH^mut/1p19q^noncodel (48%) experienced MT. Patients with IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel tumors underwent a median of 9 (first and third quartiles: 6 and 13) and 8 (6 and 12) MRI scans, respectively. The median tumor volume at diagnosis of IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel was 48.4 cm³ (first and third quartiles: 21.7 and 73.6) and 52.9 cm³ (19.8 and 90.1), respectively. Until MT, RTx alone, CTx alone, or both, were administered to 10 (8%), 23 (19%), and 34 (28%) patients with IDH^mut/1p19q^codel tumors, and 31 (20%), 12 (8%), and 45 (29%) patients with IDH^mut/1p19q^noncodel tumors. Approximately half of the adjuvant CTx was temozolomide, and the remaining protocol was procarbazine, vincristine, and/or nimustine (ACNU), which was frequently substituted for CCNU in Japan (the regimen called PAV; ref. 36).

Patients with IDH^mut/1p19q^codel tumors had a significantly longer OS and MT-free survival [OS: median 7,405 days (95% confidence interval, CI, 6,960–not reached); MT-free survival: 4,079 days (95% CI, 3,791–not reached)] than those with IDH^mut/1p19q^noncodel tumors [OS: median 3,843 days (95% CI, 3,182–not reached), MT-free survival: 2,508 days (2,333–3,511)] (both P < 0.0001, log-rank test; Fig. 1A; Supplementary Fig. S1A). Comparing the duration before and after MT, the median survival after MT was 1,818 days (95% CI, 723–not reached) for patients with IDH^mut/1p19q^codel tumors and 614 days (482–791) for patients with IDH^mut/1p19q^noncodel tumors, which was less than half of the MT-free survival (Fig. 1B). However, MT-free survival correlated significantly with OS in both subtypes [r = 0.98 (95% CI, 0.93–1.00) for IDH^mut/1p19q^codel; 0.90 (95% CI, 0.83–0.95) for IDH^mut/1p19q^noncodel, Pearson correlation coefficient] (Fig. 1C and D). Thus, efforts should be made to prolong MT-free survival. In log-rank tests using clinical factors (sex, age, and tumor volume at diagnosis), tumor volume at diagnosis was significantly associated with MT-free survival, as reported previously (Fig. 1E and F; Supplementary Fig. S1B and S1C; ref. 37). In our cohort, MT-free survival did not vary significantly according to the type of the adjuvant therapies [P = 0.53 (IDH^mut/1p19q^codel), 0.40 (IDH^mut/1p19q^noncodel); Supplementary Fig. S2A]. However, the relationship between adjuvant therapy and MT-free survival should be carefully considered because the tumor volume at the time of surgery significantly varies with and without adjuvant therapy in our dataset (Supplementary Fig. S2B).

Using tumor volume datasets from patients who underwent at least three times head MRIs for each treatment (pre-surgery, post-surgery, post-RTx, during-CTx, post-CTx, during-CRTx, and post-CRTx), we first examined whether linear or exponential regression models fit better by comparing the Akaike information criterion (AIC) in IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel tumors. The exponential model had a better fit than the linear model (P = 0.008 for IDH^mut/1p19q^codel; P = 0.042 for IDH^mut/1p19q^noncodel, Wilcoxon signed-rank test; Fig. 2A). Other growth models, such as a logistic model, did not converge for most patients because of data sparsity. Next, we computed the net proliferation rate in each treatment period using the cases that fit significantly better to the exponential model (P < 0.05). In both IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel tumors, there was a significant decrease in the net proliferation rate of post-RTx, during CTx, during CRTx, and post-CRTx, compared with that of surgery alone, which included both pre- and post-surgery cases (P < 0.01, Wilcoxon rank-sum test; Fig. 2B). In contrast, the net proliferation rate of post-CTx did not differ from that of surgery alone. Notably, those of pre-surgery (P = 0.50 for IDH^mut/1p19q^codel; P = 0.31 for IDH^mut/1p19q^noncodel, Wilcoxon rank-sum test), and the post-second and subsequent surgery (P = 0.40; IDH^mut/1p19q^codel; P = 0.12 for IDH^mut/1p19q^noncodel, Wilcoxon rank-sum test) were similar to those of the post-initial surgery. These results suggest that the type of adjuvant therapy affects the suppression of coefficient of net proliferation rates in both IDHmut-LGG subtypes. When we compared the net proliferation rate between two subtypes, only the post-surgery net proliferation rate exhibited a significant difference (median net proliferation = 0.00073 for IDH^mut/1p19q^codel, 0.00125 for IDH^mut/1p19q^noncodel, P = 0.00026, Wilcoxon rank-sum test; Supplementary Table S2). The net proliferation rate of LGGs treated with surgery alone has been estimated by the velocity of diameter expansion (VDE) as an average of about 4 mm/year in previous reports (38, 39). In our cases, the estimates also showed consistency with the reports as the VDE [deduced from volume by the formula (2 × tumor volume)^1/3] of IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel tumors with a volume of 20 cm³ at diagnosis was 3.1 and 5.4 mm/year, respectively.

When we examined the differences among CTx (temozolomide vs. PAV), there were no statistically significant differences in net proliferation rate under each CTx protocol (Supplementary Table S3) or MT-free survival for both genetic subtypes of IDHmut-LGGs (P = 0.17 for IDH^mut/1p19q^codel; P = 0.60 for IDH^mut/1p19q^noncodel, log-rank test). Therefore, for simplicity, we did not distinguish the type of CTx used in this study.

Because we determined the most appropriate growth kinetics as an exponential growth (Fig. 2A), we developed an exponential growth model of glioma cells. The dynamics are given by

Here y(t), α, β, and t were denoted by tumor volume at time t, a net proliferation rate, an initial tumor volume, and time, respectively. Because a net proliferation rate varied among treatments (Fig. 2B), we considered different growth parameters according to treatment conditions. Then, the dynamics between t₁ and t₂ under treatment option i is given by

Here α_i and β₁ were a net proliferation rate between t₁ and t₂ under treatment option i and the tumor volume at t₁, respectively. Treatment options were divided into seven: pre-surgery, post-surgery, post-RTx, during CTx, post-CTx, during CRTx, and post-CRTx (Fig. 2B).

Next, let us consider the evolutionary path of IDHmut-LGGs to MT. We assumed that the occurrence of MT was a consequence of one or more “MT events” accumulated in the tumor with a small probability per unit volume per unit time. Here, the events were supposed to be unknown (epi)genetic alterations within the tumor. The standard alteration rate ignoring treatment effects is denoted by m. We hypothesized that adjuvant therapies changed the rates additionally with coefficient, |$\varepsilon $| , based on the previous reports that treatments induced (epi)genetic alterations (23). The additional coefficients under each therapy—RTx, CTx, and CRTx were denoted by |${\varepsilon _{RTx}}$| , |${\varepsilon _{CTx}}$| , and |$\ {\varepsilon _{CRTx}}$| , respectively, where |${\varepsilon _{CRTx}}$| were defined to be higher than |$\ {\varepsilon _{RTx}}$| and |${\varepsilon _{CTx}}$| . Then the expected number of “MT events” named the cMTrisk was given by

Here, cumulative tumor volumes under surgery alone, RTx, CTx, and CRTx were denoted by S_s, S_RTx, S_CTx, and S_CRTx, respectively (Fig. 3A). The cumulative tumor volume (S) between t₁ and t₂ under treatment option i was calculated as

When we used time-series data of tumor volume, the equation for calculating the cumulative tumor volume between t₁ and t₂, S, was approximately given by

Here, V_j is the tumor volume at time j. All information about V_j and time was included in the datasets. We used Eq. C and Eq. E for calculating the cMTrisk of each patient.

The model considered an exponential expansion of IDHmut-LGG cells starting from a single cell. The number of tumor cells in 1 cm³ tumor bulk was estimated as 1 billion (40). The initial time of tumorigenesis was calculated using tumor volume at diagnosis and the estimated net proliferation rate of the pre-surgery period. We also estimated the tumor volume at the start and end of each treatment using Eq. B.

On the basis of the time-series data, including the tumor volumes and treatment history of IDHmut-LGGs, we estimated the additional coefficients of alteration rates under each treatment, |$\varepsilon $| . The clinical dataset was randomly and evenly divided into training and validation sets. These two subsets did not show significant differences in clinical characteristics (Supplementary Table S4). By using a training set, the cMTrisk of each patient was calculated with a large variety of parameter combinations of |${\varepsilon _{RTx}}$| , |${\varepsilon _{CTx}}$| , and |${\varepsilon _{CRTx}}$| using Eq. C. According to the MT occurrence in each patient, we estimated the coefficients to maximize the AUC of the ROC curve for MT. As a result, the |${\varepsilon _{RTx}}$| , |${\varepsilon _{CTx}}$| , and |${\varepsilon _{CRTx}}$| were found to be 1.8, 2.2, and 2.6 in IDH^mut/1p19q^codel tumors and 2.7, 1.8, and 2.8 in IDH^mut/1p19q^noncodel tumors, respectively (Fig. 3B).

Furthermore, we estimated the number of MT events required for MT. Using Eq. C, the expected number of “MT events” was calculated. Because the events were supposed to happen rarely, a Poisson distribution was adopted to predict the MT occurrence with mean cMTrisk. Then, if the number of MT events required for transformation was more than k, the probability of MT denoted by P was given by

Here, λ denotes cMTrisk. By comparing AIC with different MT events required for MT, the best number and the standard alteration rate m were estimated. The numbers of required MT events that showed the lowest AIC value in the training set, suggesting the best-fit model, were three for IDH^mut/1p19q^codel and two for IDH^mut/1p19q^noncodel (Fig. 3C and D; Supplementary Fig. S3), where the standard alteration rate m was identified to be 8.99·10^–6/cm³/day for IDH^mut/1p19q^codel and 8.77·10^–6/cm³/day for IDH^mut/1p19q^noncodel. Notably, clinical outcomes could be accurately predicted for the validation and training sets using these parameters (Fig. 3C).

Simultaneously, there was a simple idea to regard the MT occurrence as a time-dependent phenomenon rather than developing a mathematical model as shown above. Thus, we also calculated the probability of MT occurrence (P) by using Eq. F with λ proportional to time instead of cMTrisk, assuming the “MT events” were accumulated just over time. As a result, our mathematical model demonstrated a better fit and lowered residual sum of squares than the time-dependent one in both genetic subtypes (Fig. 3E and F; Supplementary Fig. S4A–S4C).

To delineate the association between cMTrisk and genetic alterations, we performed WES of samples obtained by multiple temporal sampling, including 38 samples from 17 patients with IDH^mut/1p19q^codel tumors and 62 samples from 28 patients with IDH^mut/1p19q^noncodel tumors (Fig. 4A). Excluding DNA hypermutators (n = 2), the TMB was significantly associated with cMTrisk in both IDHmut-LGG subtypes (Fig. 4B).

Several markers for shorter OS or MT-free survival were previously proposed: PI3K pathway alteration and mutation of notch receptor 1 (NOTCH1) in oligodendroglial tumors and focal amplification of platelet-derived growth factor receptor alpha (PDGFRA), mutation of PI3K, RB pathway alteration [including homozygous deletion of cyclin-dependent kinase inhibitor 2A and 2B (CDKN2A and CDKN2B; hereafter called CDKN2A/B) and focal amplification of cyclin-dependent kinase 4 (CDK4) in astrocytic tumors (20, 41–46). Mutations of PI3K [including phosphatidylinositol-4,5-bisphosphate 3-kinase catalytic subunit alpha (PIK3CA) and phosphoinositide-3-kinase regulatory subunit 1 (PIK3R1)] and NOTCH1 were identified in 4 (24%) and 3 (18%) cases, respectively, in recurrent IDH^mut/1p19q^codel tumors even before MT, but not in primary tumors (P = 0.0072, Fisher exact test). These genetic mutations did not necessarily coincide with MT; in particular, all PI3K mutations had been identified before (Fig. 4C). Next, we investigated the correlation between clinical characteristics and genetic alterations of IDH^mut/1p19q^codel tumors using targeted sequencing of samples obtained at initial surgery. Net post-surgery proliferation rate in the PI3K-mutant IDH^mut/1p19q^codel tumors was significantly higher than that of the PI3K wild-type IDH^mut/1p19q^codel tumors (P = 0.042, Wilcoxon rank-sum test), unlike the NOTCH1 mutation (Fig. 4D). At the same time, the presence of PI3K mutations was associated with a decrease in MT-free survival (P = 0.0056, log-rank test; Fig. 4E). These results suggest that PI3K mutations are not involved in tumorigenesis but are acquired during tumor expansion and result in a significant increase in the net proliferation rate and decreased MT-free survival.

Focal amplification of CDK4 (n = 2; 7%) and PDGFRA (n = 4; 14%) and homozygous deletion of CDKN2A/B (n = 5; 19%) were identified in recurrent IDH^mut/1p19q^noncodel tumors, but not in primary tumors (P = 0.0044, Fisher exact test; Fig. 4A). All these focal CNVs coincided with MT (P = 0.029, Fisher exact test; Fig. 4C). Neither was observed in the SNP array data of the IDHmut-LGG samples obtained at initial surgery (n = 57). These results suggest that these three focal CNVs are not involved in tumorigenesis but are acquired during tumor expansion and are significantly associated with MT.

We determined the optimal therapeutic strategy that would prolong the MT-free survival using the mathematical model. Tumor volume at diagnosis was a significant prognostic factor for MT-free survival in both IDHmut-LGG subtypes, as described above (Fig. 1C). Tumor volume at surgery and the extent of resection (EOR) strongly influence cMTrisk in patients with IDHmut-LGGs because cMTrisk strongly depends on cumulative tumor volume. Here, we hypothesized that these two factors would substantially impact the optimal therapeutic strategy to prolong MT-free survival. We simulated the scenario with IDHmut-LGGs classified by tumor volume at surgery and the EOR and investigated the optimal therapeutic strategy for each case. Figure 5A shows the estimated MT-free survival in 50% of patients with IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel tumors who underwent surgery alone. We also estimated the MT-free survival of patients with the PI3K-mutant IDH^mut/1p19q^codel tumors. In all genetic subtypes, the smaller the tumor volume at surgery and the higher the EOR, the longer the MT-free survival. To elucidate optimal adjuvant therapy and timing, simulated cases received adjuvant therapies: RTx, CTx, or CRTx, at 30 days or later after surgery. When the tumor was small (≤ 50 cm³ in IDH^mut/1p19q^codel), prompt initiation of CRTx was the optimal treatment regardless of the EOR (Fig. 5A). In contrast, when the IDHmut-LGG tumor was large (> 50 cm³), optimal treatment differed according to genetic alteration and the EOR; adjuvant therapies prolonged MT-free survival in IDH^mut/1p19q^noncodel tumors, but it was often shortened in IDH^mut/1p19q^codel tumors (Fig. 5A). PI3K-mutant IDH^mut/1p19q^codel tumors exhibited increased net proliferation rate and decreased MT-free survival compared with PI3K wild-type tumors (Fig. 4D and E), which led to the essential need for adjuvant therapy to delay MT (Fig. 5A). When it comes to the relationship between the timing of adjuvant therapy and MT-free survival in cases where adjuvant therapy was effective, it was most effective to initiate adjuvant therapy as soon as possible to delay MT regardless of tumor volume, the EOR, and genetic subtype (Supplementary Fig. S5). In our dataset, administration of CRTx within 50 days after surgery was significantly associated with increased MT-free survival among IDH^mut/1p19q^noncodel patients with tumor volume ≤ 100 cm³, [HR (95% CI): 0.29 (0.10–0.83); P = 0.021, tumor volume-stratified Cox regression analysis] (Fig. 5B). In IDH^mut/1p19q^codel patients with tumor volume ≤ 100 cm³ [or ≤ 50 cm³ when EOR was partial resection (PR)], administration of CRTx within 50 days after surgery exhibited longer MT-free survival, but this difference was not significant [HR: 0.30 (0.038–2.36); P = 0.25, tumor volume-stratified Cox regression analysis] (Fig. 5B). IDH^mut/1p19q^codel patients who did not receive CRTx within 50 days after surgery were significantly more likely to undergo reoperation than those who did (P = 0.04, Wilcoxon rank-sum test), which may have influenced the difference in MT-free survival (Supplementary Fig. S6). The differences in MT-free and reoperation-free survival between the patients who received CRTx within 50 days after surgery and those who did not were significant [HR: 0.11 (0.015–0.86); P = 0.036 in IDH^mut/1p19q^codel; HR: 0.21 (0.074–0.59); P = 0.0033 in IDH^mut/1p19q^noncodel, tumor volume-stratified Cox regression analysis] (Fig. 5C).

Our analyses revealed that small-sized tumors in IDHmut-LGGs had prolonged MT-free survival and were more effectively treated with adjuvant therapy than large ones. We also examined MT-free survival with an assumption of earlier diagnosis of the same tumors. Comparing the case where the IDHmut-LGG with a tumor volume of 50 cm³, the maximal tumor volume for which prompt initiation of CRTx was optimal treatment, was surgically removed and immediately treated by CRTx to the case where the same tumor had been diagnosed when the tumor volume was 10 cm³, surgically removed and followed by prompt initiation of CRTx, the latter scenario had substantially delayed MT (Fig. 5D; Supplementary Fig. S7). These results suggest the importance of early diagnosis and intensive treatment, if possible, before the onset of neurologic symptoms.

We developed a mathematical model of tumor progression based on serial tumor volume data, treatment history, and genome-wide mutational analysis of 276 IDHmut-LGGs to understand the tumor dynamics and estimate the risk of MT and suggest optimal therapeutic strategies. Our analyses successfully predicted the MT-free survival with validation, revealed a link between genetic alterations and the increased risk of MT or rapid tumor expansion, and emphasized the importance of early diagnosis and treatment, and therapeutic strategy based on genetic alterations for IDHmut-LGGs.

The required number of MT events to undergo MT was different for IDH^mut/1p19q^codel and IDH^mut/1p19q^noncodel tumors: three and two, respectively. For example, when cMTrisk was 2, 32% of IDH^mut/1p19q^codel tumors were predicted to undergo MT, whereas approximately twice as many cases of IDH^mut/1p19q^noncodel tumors (59%) would become malignant in IDH^mut/1p19q^codel cases. The difference in the required number of MT events and the net proliferation rate of surgery alone is consistent with the decreased MT-free survival of IDH^mut/1p19q^noncodel compared with IDH^mut/1p19q^codel tumors. The mutation of PI3K pathway, which is an intracellular signal transduction pathway that promotes metabolism, proliferation, and cell survival, have been reported in many carcinomas, including WHO grade II and III (II/III) IDH-mutant gliomas with 1p19q codeletion (13). The PI3K-mutant WHO grade II/III IDH-mutant gliomas with 1p19q codeletion successfully propagated patient-derived xenografts, unlike the PI3K wild-type, suggesting an association with a malignant phenotype (41). In this study, although the number of patients harboring PI3K mutations was not so large, these mutations were predominantly found in recurrent tumors compared with paired primary tumors and were associated with increased post-surgery net proliferation rate and decreased MT-free survival in IDH^mut/1p19q^codel cases. However, these PI3K mutations did not coincide with MT. These results are reasonable, given that IDH^mut/1p19q^codel requires three or more MT events for transformation.

It is established that temozolomide inordinately induces somatic mutations in LGGs (23). In the current study, the TMB was significantly associated with cMTrisk in both subtypes (Fig. 4B), suggesting the possibility that the additional coefficients of (epi)genetic alteration rates under therapy, such as CTx ( |${\varepsilon _{CTx}}$| ) and CRTx ( |${\varepsilon _{CRTx}}$| ), could be different between temozolomide and PAV. We estimated the additional coefficients of temozolomide and PAV separately using the whole dataset. The CRTx with temozolomide ( |${\varepsilon _{CRT\_TMZ}}\ = \ 3.4$| ) was approximately 1.5 times higher than the CRTx with PAV ( |${\varepsilon _{CRT\_PAV}}\ = \ 2.2$| ) in IDH^mut/1p19q^codel tumors. Consequently, the prompt initiation of CRTx with PAV was the optimal treatment for most patients with IDH^mut/1p19q^codel tumors and those with IDH^mut/1p19q^noncodel tumors (Supplementary Fig. S8; Supplementary Materials and Methods), which is consistent with previous reports that CRTx with PCV significantly prolongs PFS and OS in patients with IDH^mut/1p19q^codel tumors and those with IDH^mut/1p19q^noncodel tumors (47). Contrastingly, temozolomide, which has not been shown to increase MT-free survival in IDH^mut/1p19q^codel in a randomized controlled study, was not generally effective for patients with IDH^mut/1p19q^codel tumors. Notably, these results were not validated in an independent cohort, due to the limited number of cases that received each treatment in our dataset. Furthermore, although PAV was widely used rather than PCV in Japan and has been reported to be as effective as PCV in patients with malignant glioma, there are no randomized controlled trials showing the efficacy of PAV (36). Further studies are warranted on the differential effects of temozolomide and PCV (or PAV) on MT.

Our mathematical model revealed that early diagnosis and gross total resection (GTR) followed by intensive adjuvant therapy are the optimal treatment strategies for delaying MT. Ideally, postoperative therapy would not be necessary if all tumor cells, including cells infiltrating the surrounding area, can be removed. Recently, “supratotal” resection, an extended resection beyond MRI-defined abnormalities, has been suggested for the treatment of LGGs (48). Supratotal resection could be added to the optimal treatment strategy. Of note, patients with slight recurrence or incomplete tumor resection, including GTR, should receive intensive adjuvant therapy. Further studies are required to confirm this hypothesis.

Overall, we have developed a model that can accurately estimate the MT-free survival and predict the optimal treatment strategy to minimize the risk for IDHmut-LGGs based on genetic alterations. The methodology of our study can be applied to different types of cancer, which benefits to find an optimal treatment strategy to improve quality of life and survival. Of course, we are proposing a simplified model for the tumor proliferation of IDHmut-LGGs and therefore miss several key characteristics that can be readily explored in future models. An important area of future exploration for this model is understanding the differential effect of temozolomide and PCV (or PAV) on MT and the effect of supratotal resection in IDHmut-LGGs. Further prospective studies are warranted to support the therapeutic strategies suggested in this study.

K. Aoki reports personal fees from Daiichi Sankyo outside the submitted work. M. Takahashi reports grants and non-financial support from Eisai Co., Ltd., outside the submitted work. Y. Narita reports grants and personal fees from Ono Pharmaceutical Company, Eisai, Dainippon-Sumitomo, and Daiichi-Sankyo; grants from Taiho, Ohara, Abbvie, and Stella-Pharma; personal fees from Chugai and Novocure outside the submitted work. Y. Muragaki reports personal fees from Abbvie, Ono Pharmaceutical, Chugai Pharma, Novartis, and Bristol-Myers Squibb; grants and personal fees from Daiichi Sankyo, MSD, Eisai, and Hitachi; and grants from Otsuka during the conduct of the study. H. Haeno reports grants from National Cancer Center Research and Development Fund and SRL and H.U. Group Research Institute during the conduct of the study and grants from Boehringer Ingelheim Japan, AMED, and JSPS outside the submitted work. A. Natsume reports grants from Daiichi-Sankyo and personal fees from Roche outside the submitted work. No disclosures were reported by the other authors.

K. Aoki: Conceptualization, resources, data curation, formal analysis, funding acquisition, validation, investigation, visualization, methodology, writing–original draft, project administration, writing–review and editing. H. Suzuki: Resources, data curation, project administration. T. Yamamoto: Resources, data curation, investigation. K.N. Yamamoto: Methodology. S. Maeda: Resources, data curation. Y. Okuno: Methodology. M. Ranjit: Resources, data curation. K. Motomura: Resources. F. Ohka: Resources. K. Tanahashi: Resources. M. Hirano: Resources. T. Nishikawa: Resources. H. Shimizu: Resources. Y. Kitano: Resources. J. Yamaguchi: Resources. S. Yamazaki: Resources. H. Nakamura: Resources, data curation. M. Takahashi: Resources, data curation. Y. Narita: Resources, data curation. M. Nakada: Resources, data curation. S. Deguchi: Resources, data curation. M. Mizoguchi: Resources, data curation. Y. Momii: Resources, data curation. Y. Muragaki: Resources, data curation. T. Abe: Resources, data curation. J. Akimoto: Resources, data curation. T. Wakabayashi: Resources, supervision. R. Saito: Supervision. S. Ogawa: Supervision, project administration. H. Haeno: Funding acquisition, investigation, methodology, writing–original draft, project administration. A. Natsume: Funding acquisition, investigation, project administration, writing–review and editing.

The authors owe their thanks to the patients. The supercomputing resource was provided by Human Genome Center, the Institute of Medical Science, the University of Tokyo (http://sc.hgc.jp/shirokane.html).

This research was supported by the Grant-in-Aid for Scientific Research (B) 20H03789, Grant-in-Aid for Scientific Research on Innovative Areas 17H06356, and Grant-in-Aid for Young Scientists 20K17960, Grant-in aid from the Japan Brain Foundation, National Cancer Center Research and Development Fund (2020A-7), and a research grant from SRL, and H.U. Group Research Institute.

The publication costs of this article were defrayed in part by the payment of publication fees. Therefore, and solely to indicate this fact, this article is hereby marked “advertisement” in accordance with 18 USC section 1734.