Mathematical modeling for glioblastoma treatment: scenario generation and validation for clinical patient counseling

Introduction: Glioblastoma (GBM) is an aggressive primary brain tumor. Despite standard treatment, recurrence is common, and patient counseling remains challenging. Mathematical modeling offers a potential strategy to simulate tumor behavior and personalize care. This study evaluates whether a simple reaction-diffusion model can generate realistic scenarios of treatment outcomes for individual patients with recurrent GBM using clinical imaging data.

Methods: We retrospectively analyzed 132 MRI intervals from 46 patients who underwent treatment for recurrent GBM. T1 post-contrast and T2/FLAIR images were co-registered and manually segmented to identify enhancing tumor and edema. Using a systematic parameter sampling design, tumor growth between successive scans was simulated 18 times with a reaction-diffusion equation, the “ASU-Barrow” model, to generate realistic ranges of tumor response to treatment, as evaluated by clinical imaging.

Results: Model-generated scenarios for changes in tumor volumes well approximated the observed ranges in the patient data. In 86% of the imaging intervals, at least one scenario yielded a simulated tumor volume that agreed to within 20% of the observed one (and to within 10% in 65% of the cases). Spatial accuracy was assessed using agreement and containment scores, indicating how well the predicted tumor matched the real one. The best simulations achieved an agreement of 0.52 and a containment score of 0.69. These results suggest that a simple model can generate a realistic range of outcomes, over intervals of two or three months, in a majority of patient cases.

Conclusion: This reaction-diffusion model simulates likely ranges of GBM progression under treatment with reasonable accuracy and modest computational needs and may yield a clinically practical tool to support patient counseling. Incorporating advanced imaging, such as perfusion MRI, may further improve accuracy. With further development, our approach could provide personalized scenarios of treatment outcomes that could aid in patient counseling.

1 Introduction

Since the early twentieth century, mathematical models have been developed to understand tumor growth dynamics and optimize treatment strategies (1). The simplest population growth model is exponential: the population grows by equal ratios in equal intervals of time (e.g., the number of tumor cells doubles every week). The nineteenth century actuary Benjamin Gompertz proposed a modification in which the time intervals lengthen exponentially, so that the growth rate decelerates with time, reflecting resource limitations and, in the case of tumors, increasing cell death rates as the tumor becomes larger (2). The logistic model, also proposed in the 19th century, is a different modification of the basic exponential in which the growth ratio decreases linearly with the population size.

As oncology shifts towards individualized treatment approaches, some research efforts have explored the potential for mathematical modeling to optimize and personalize therapy. Examples include immunotherapy (3) and efforts to optimize regimens of chemotherapy and radiotherapy (4). Credibility assessments of computational models for medical devices and healthcare applications are an active area of research (5, 6). With respect to glioma specifically, Hormuth and collaborators have explored measures of tumor call heterogeneity (7) to predict future tumor growth and the use of imaging and mathematical modeling to predict treatment response (8). (Section 5 contains further discussion.)

Patients who are newly diagnosed with glioblastoma (GBM) typically are treated with maximal safe resection, radiotherapy, and temozolomide chemotherapy (9, 10). There is no standard of care for recurrent tumors (11, 12), but options include chemotherapy, radiation, electric fields, immunotherapy, and reoperation (13, 14). The choice of therapy depends upon tumor response and the patient’s performance status (15).

Nevertheless, there have been no treatment breakthroughs for GBM in the past 20 years, and improvements in patient survival in the interim are due largely to better supportive care (16). Patients with recurrent GBM face a poor prognosis, with tumor-related neurological decline and impaired quality of life as key concerns. Shared decision-making between patients, their families, and clinicians is important (17). Given the difficulties with the current state of the art, new tools are needed to assist with clinical counseling of patients.

We have undertaken this modeling study with a mindset like that of financial planners, who provide personalized estimates of portfolio performance to give clients an approximate idea of future retirement income. They do not attempt to model the global economy or predict the stock market. Instead, they use simple models (e.g., the time value of money) and sample prior rates of return in a statistically reasonable way to generate a range of plausible scenarios, based on a client’s current holdings, savings rate, and anticipated expenses. This procedure cannot account for every potential outcome, nor does it make specific predictions, but it can provide an idea of what to expect under certain hypotheses and serve as an advising tool.

In a similar vein, the objective of this study is to determine whether a simple mathematical model can simulate a personalized, realistic range of potential outcomes over the next 2 to 3 months, corresponding to typical follow-up intervals, for patients who experience a recurrence of their GBM tumors following initial treatment and whose response to future therapy cannot be predicted with certainty. To our knowledge, the potential for mathematical modeling to provide personalized scenarios of GBM progression in individual patients, which clinicians could use to assist with patient counseling, has not previously been assessed.

The biology of GBM is complex, and many details are poorly understood. Factors that may affect tumor progression and response to treatment include genomic abnormalities (18), neovascularization (19), hypoxia (20), immunosuppression (21), GBM stem cells (22), and the brain’s extracellular matrix (23), to name a few. A fundamental difficulty is that there is no established biochemical principle to justify a particular mathematical model of any of these processes and their interactions. Any associated rate equation is speculative at best and introduces coefficients and initial conditions that cannot be measured in living patients. Consequently, the biological correctness of the resulting model cannot be independently validated (24). A complex model is unlikely to be amenable to rigorous analysis, and its output may depend sensitively on parameters that cannot be estimated reliably (25), making it impractical for clinical application.

For these reasons, we focus on a model whose parameter space can be tractably sampled and that can be initialized from (and compared against) imaging data that are collected as part of routine clinical surveillance. We reiterate that we cannot make a specific prediction about a given patient’s clinical course; instead, our goal is to develop a data-driven modeling system, using a simple model and validated against the clinical trajectories of many previous patients, that can generate a personalized range of realistic scenarios of treatment outcomes. As noted above, our focus is on time horizons of 2 to 3 months, corresponding to typical surveillance imaging intervals. In this way, clinicians could use the outputs as a tool to help patients and their families make better-informed decisions about continued treatment.

2 Materials and methods

2.1 Patient population and data preprocessing

Magnetic resonance imaging (MRI) scans were obtained from the Barrow Neurological Institute (BNI) patient archive pursuant to St. Joseph’s Hospital and Medical Center Institutional Review Board (IRB) protocol PHX-19-500-182-20-08. They consist of 357 axial surveillance scans from 75 unique patients, ranging in age from 25 to 78 years (median 62), who previously had undergone maximal safe surgical resection, followed by radiotherapy and, in most cases, temozolomide chemotherapy. The first scan from each patient series is obtained 2 to 12 weeks after initial resection; each patient has at least 2 (and up to 9, median 4) scans that include T1 plus contrast and T2/FLAIR sequences. The median time interval between scans is 57 days (range: 16–121 days).

The imaging workflow proceeds as described in previous work by the authors (26) and consists of the following steps.

1. The patient imaging data consists of DICOM files generated by the MRI scanners. Identifying information is removed, and the data, consisting of axial slices, is assembled into a single three-dimensional image and converted to the NIfTI format (27) using SPM-12 (28); this is done independently for each T1 and T2 modality.

2. SPM-12 also is used to co-register all images to a common time point (usually the first scan in each patient series). The automated algorithms in SPM-12 are used to generate another NIfTI file that contains only the brain (the skull and eyes are stripped out), and the brain domain is further segmented into regions of cerebrospinal fluid (CSF) and white and gray matter. This labeled, patient-specific brain domain serves as the computational domain for each simulation. Altogether, it contains 24 to 30 axial slices of typically 256 × 256 voxels. The horizontal resolution varies from 0.85 to 0.93 mm and the vertical from 7 to 7.5 mm, depending on the number of slices in the original scans.

3. All MRI scans are manually segmented by neuroimaging experts. For this purpose, the co-registered scans are loaded into the 3D-Slicer image processing platform (29), which is used to “paint” the relevant voxels and separate them from anatomically normal tissue according to the neurosurgical judgment of the operator. Tumor voxels are segmented into three categories: necrotic core (hypointense on T1), enhancing tumor (contrast enhancement on T1), and tumor-associated edema (hyperintense on T2-FLAIR). The modeling software uses the segmentations to generate initial conditions, as described below.

For each patient, and at every time point, this procedure produces six or seven anonymized files in NIfTI format: the three-dimensional T1C and T2 images; the derived computational domain; and the manually produced segmentations of enhancing tumor, edema, resection cavity, and, if applicable, necrotic tumor core. All files are named and stored in a consistent way. A separate spreadsheet tracks the number of days between each scan, and a database, maintained at St. Joseph’s Hospital, links the numbered patients to the original scans and medical records.

Scans from the patient database are excluded from further analysis under any of the following circumstances:

● SPM-12 is unable to segment the brain domain (possibly as a result of motion artifacts in the scan), because we cannot construct a patient-specific computational domain;

● no enhancing tumor is apparent on T1 post-contrast imaging, because we cannot compute relative changes in tumor size; or

● the patient undergoes another surgical resection, which is beyond the scope of the model.

Following these exclusions, we have a complete set of preprocessed scans from 46 unique patients for 132 time intervals. The first two preprocessing steps outlined above can be done on a laptop in a few minutes, but the tumor segmentation is performed manually for reliability and thus is the most time-consuming step.

2.2 Mathematical model

Our modeling effort focuses on the gross total expansion (or contraction) of the tumor. Tumor volume is a prognostic factor for overall survival (30–32) and can be readily measured using the 3D-Slicer platform (33) from the surveillance imaging that serves as our data source.

The Fisher-Kolmogorov-Petrovsky-Piscunov (FKPP) equation is a reaction-diffusion equation that was proposed initially to model the spread of invasive species (34, 35). This model also has been proposed to describe the diffusive spread of GBM tumor cells throughout the brain parenchyma (36) and is given by

$\begin{array}{ll}\frac{\partial u}{\partial t}=\mathrm{\nabla }·(D(x)\mathrm{\nabla }u)+\rho (x)u\left(1-\frac{u}{K}\right)& (1)\end{array}$

where $u=u (x,t)$ is the space- and time-dependent tumor cell population, K is the local carrying capacity, and D and ρ are the diffusion and proliferation rates, respectively, which may be taken as constant or allowed to vary by location (e.g., faster diffusion in white matter than gray). Equation 1 admits traveling-wave solutions whose propagation speed depends on the model parameters, but the details also depend on the rates at which glioma cells switch between proliferative and motile phenotypes (37). Equation 1 has been suggested as a way to quantify the effect of surgical resection (38); the effect of chemotherapy (39); and to explain differences in patient survival after similar courses of treatment (40).

The model (1) assumes implicitly that all cells in the modeled population have the same growth potential. One characteristic of GBM tumors, however, is the presence of a hypoxic or necrotic “core” (41). The authors have previously proposed a modification of the FKPP equation (26) to account for this phenomenon, given by the so-called “ASU Barrow” (ASUB) model, which partitions the tumor cells into proliferating (p) and quiescent (q) subpopulations:

$\begin{array}{ll}\frac{\partial p}{\partial t}=\mathrm{\nabla }·(D(x)\mathrm{\nabla }p)+g(p,q,t)-h(p,q,t)& (2)\end{array}$

$\begin{array}{ll}\frac{\partial q}{\partial t}=h(p,q,t).& (3)\end{array}$

The densities of proliferating and quiescent cells are given by the respective space- and time-dependent functions p = p(x, t), and q = q(x, t). The overall dynamics are simple: proliferating cells diffuse at a rate given by ∇ · (D(x)∇p); they grow an a net per capita rate g(p, q, t) and they become quiescent (i.e., die or become hypoxic) at a per capita rate h(p, q, t). The function D(x) defines the rate at which cells infiltrate the brain. Quiescent cells simply accumulate, reflecting high cellularity but no net tumor growth.

In Equations 2, 3, h(p, q, t) reflects the net cell-killing effect of treatment; g(p, q, t) captures treatment resistance and net proliferative tendencies. In the version simulated here, we take

$\begin{array}{ll}g(p,q,t)=\rho (x)p[1-\delta (p+q)]& (4)\end{array}$

$\begin{array}{ll}h(p,q,t)=k(x)p\delta (p+q),& (5)\end{array}$

where ρ(x) and k(x) are the space-dependent maximum growth and quiescence rates, and δ is a monotonically increasing function of the total cell density.

We take δ(x) to be piecewise constant: we fix one value in edematous tissue and a possibly different value in the enhancing rim of the tumor; similarly for k(x). One rationale is that some treatments, such as chemotherapy wafers or localized radiotherapy, may be applied (and have greater effect) on regions corresponding to contrast enhancement. We treat the diffusion rate analogously: D(x) = D_w when x corresponds to a location in a white-matter tract or tumorous region; D(x) = D_w/2 in gray matter; and D(x) = 0 in CSF.

We normalize the maximum cell density to 1 and require δ(0) = 0 and δ(1) = 1. In regions where the cell density is low, 1 − δ is close to 1; therefore, $g\approx \rho p$ and $h\approx 0$, so that the net proliferation rate is approximately exponential. In other words, small cell populations grow at a rate that is roughly proportional to $e^{\rho t}$, which implies that the overall tumor growth rate can be sensitive to the choices of ρ in edematous regions, where tumor cell densities are presumed to be low. One convenient choice for δ is the beta cumulative distribution function B(x; a, b), which increases monotonically from 0 to 1 across the unit interval; we fix a = 3 and b = 1 (26).

2.3 Simulation of the model

The SPM-12 brain segmentation defines the patient-specific computational domain. No-flux boundary conditions are imposed at the interface with the skull and CSF. The model can be integrated efficiently using the IRKC solver by Shampine et al. (42). One 60-day simulation on a typical 256 × 256 × 25 patient brain domain takes about 18 seconds on a single CPU core of a modern personal computer. Multiple patient scenarios can be simulated in parallel on a multicore machine.

As described above, the growth and quiescence functions $\rho (x)$ and $k(x)$ in Equations 4, 5 are piecewise constant: $\rho (x)={\rho }_c$ and $k(x)=k_c$ if x corresponds to a voxel in a region segmented as contrast enhancing, and $\rho (x)={\rho }_e$ and $k(x)=k_e$ in edematous regions. To simplify their interpretation, the values for ρ and k are given in the tables below as doubling and halving times, respectively, in days, assuming pure exponential growth and decay, which occurs when the tumor cell populations are very small or large. (The actual parameters used in the model are $k=(\log 2)/\widehat{k}$ and $\rho =(\log 2)/\widehat{\rho }).$ The net doubling rate of glioblastoma tumor cells is not known, but one study estimates a tumor volume doubling time of 31 days prior to treatment (43) and another, a range from 14 to 49.6 days (44). A study of unresectable gliomas in children estimates tumor halving times of 60 to 78 days in response to radiotherapy and overall doubling times of 48 to 60 days for high-grade gliomas (45). Values of the diffusion rate D(x) also are unknown; one review of the literature has found published values that vary by four orders of magnitude [cf. Table 2 in (46)]. The FKPP Equation 1 predicts that tumors must reach a volume that is proportional to $\sqrt{D/\rho }$ to persist (47), and various studies have attempted to estimate D/ρ from MR imaging [e.g., (48)]. However, this approach is restricted to untreated tumors, as the prediction assumes a uniformly expanding cell population, and D and ρ are not identifiable from the ratio. Given the considerable uncertainties, the values of D, ρ, and k in Table 1 are roughly consistent with the values in the cited studies.

Table 1

Table 1. The 18 parameter sets used to simulate each tumor.

Initial conditions are imputed from the respective tumor segmentations. Chang and co-workers (49) found an approximate linear relationship between MR signal intensity and GBM cell density, with cellularity increasing with intensity in T1-weighted post-contrast imaging and decreasing with T2/FLAIR signal intensity. In our simulations, we ascribe an initial population within the interval $I_c=[p_1^c,p_2^c]$ that increases linearly with T1C image intensity in regions segmented as enhancing tumor and decreases linearly within the interval $I_e=[p_1^e,p_2^e]$ with T2/FLAIR signal in the edema segmentations. Counting the intervals I_c and I_e as one parameter apiece, there are 7 adjustable parameters, which are listed in decreasing order of sensitivity on the simulation results in Table 2.

Table 2

Table 2. Model and initialization parameters for the tumor simulations.

2.4 Experimental design

Of course, we do not know the “true” values for the model parameters and initial conditions, which probably vary considerably among the patients. Instead, our goal is to determine a likely range of values that yield clinically representative results. To make the computations tractable in such a large parameter space, we adopt a Taguchi experimental design, which is a type of Latin hypercube sampling. For simulation purposes, each parameter in Table 2 is assigned one of three “levels,” corresponding to a “low,” “medium,” and “high” value. There are 3⁷ = 2187 possible combinations, but we choose only 18 by prioritizing the most sensitive parameters, ${\widehat{k}}_0$ and ${\widehat{\rho }}_0$, and sampling according to an orthogonal array. The same set of 18 parameters is used to simulate every tumor. There are no stochastic components in the modeling framework presented here.

The orthogonal array consists of columns 2–8 in Table 3 of Kacker et al. (50), which contains 18 rows. Simulations are run with each of the 9 possible pairs of values of the two most sensitive parameters, ${\widehat{k}}_0$ and ${\widehat{\rho }}_0$. For each pair, two different choices of the remaining parameters are made by uniform sampling. The goal is to estimate the variability in the simulations as a function of the parameters in an economical way. Table 1 displays each of the 18 parameter sets that are used to simulate each tumor.

Table 3

Table 3. Observed and simulated chances that the patient’s tumor grows or shrinks (positive or negative percentages, respectively) by selected thresholds until the next scan time.

The simulation code is highly parallelizable, and results for all of the parameter combinations can be obtained in about 80 seconds of wall-clock time on a Macbook Pro laptop for a typical 60-day imaging interval. We have also done 144 simulations according to a sampling design that uses 12 different levels of the parameters and obtained comparable results.

3 Results

Altogether, we evaluated tumor progression across 132 time intervals from 46 unique individuals with recurrent GBM. For each interval, we compare the number of enhancing voxels, V₀, from the initial scan to the number V₁ in the next, and define the relative change as

$\begin{array}{ll}{\text{Δ}}^{\text{obs}}=\frac{V_1-V_0}{V_0}.& (6)\end{array}$

Notice that Δ^obs = −1 if there are no enhancing voxels in the comparison scan; Δ^obs = +1 if the number doubles; and Δ^obs = 0 if the number does not change.

Figure 1A shows a histogram of the distribution of Δ^obs in the clinically observed tumors, which ranges from −0.917 (i.e., reduction of 91.7%) to an expansion by a factor of 21.28; the median change is 0.0377. (The horizontal axis is truncated at Δ^obs = 2 for ease of visualization.) There are 6 cases where 2< Δ^obs< 5 and one case where Δ^obs >21.) In 75% of the cases, the observed changes range from 42.2% shrinkage to an expansion by about 111%.

Figure 1

Figure 1. (A) Histogram of relative change Δ^obs, Equation 6, of voxels in regions segmented as enhancing tumor between the starting and comparison scans over all 132 time intervals. Not shown are 6 cases (5%) where 2.22< Δ^obs< 4.61 and the one case where Δ^obs = 21.29. The median, shown as a vertical dashed line, is 0.0377. (B) Histogram of the relative change, Δ^sim, in all 18 × 132 simulated tumors over the same time interval. Not shown is the upper tail, containing 211 (8.6%) of the simulations, where Δ^sim >2 (the largest value is 74). The median is 0.0153. (C) Histogram of the relative change, R^best, defined in Equation 8 in the simulated tumors that most closely match the number of enhancing voxels in the corresponding comparison scans. Not shown is one simulation in which R^best = 3.15. The median is −0.0049.

The same set of 18 parameters, as listed in Table 1, is used to simulate each patient and imaging interval. In each simulation, we count the number of voxels in which the final tumor cell population, p + q, lies within the interval I_c. For example, in the first parameter set, I_c= [0.16, 0.80], so the voxels in which 0.16 ≤ p + q ≤ 0.80 are counted as enhancing. Likewise, if p + q lies within the interval I_e= [0.012, 0.03], then the corresponding voxel is counted as edematous. (The calculation is similar for the remaining parameters.) Let $V_i^{\text{s}\text{i}\text{m}}$ denote this number for the ith simulation and define, for $i=1,\dots ,18$, the relative change in tumor volume as

$\begin{array}{ll}{\text{Δ}}_i^{\text{sim}}=\frac{V_i^{\text{s}\text{i}\text{m}}-V_0}{V_0},& (7)\end{array}$

where V₀ is the observed number of enhancing voxels at the start of the simulation interval.

Figure 1B shows a histogram of the distribution of Δ^sim for all 18 × 132 = 2376 simulated tumors. The distribution of Δ^sim is roughly consistent with the distribution Δ^obs of observed tumor volume changes, Table 3 shows the likelihood that Δ^obs and Δ^sim falls within selected percentage ranges (negative percentages indicate the corresponding reduction in the number of enhancing voxels at the next scan time).

Each simulation produces 18 realizations (scenarios) of the evolution of the tumor. Of these, let V^* be the volume of the realization that is closest to the volume V^obs of the observed tumor at the next scan time. We define

$\begin{array}{ll}{R}^{\text{best}}=\frac{V^{*}-V^{\text{o}\text{b}\text{s}}}{V^{\text{o}\text{b}\text{s}}}& (8)\end{array}$

where V^obs is the number of voxels in the subsequent scan that have been segmented as enhancing tumor.

Figure 1C shows the distribution of $R^{\text{best}}$ across all 132 simulated time intervals. In 86 (65%) of the cases, the relative error is less than 10%, and in 113 (86%), less than 20%. Furthermore, each of the 18 parameter choices produces a “best” simulated tumor for some patient imaging interval; none of the combinations in Table 1 is redundant. These results, plus the threshold data in Table 3, suggest that the parameter ranges used in these simulations capture most of the likely variability in potential outcomes for typical patients. The relative errors exceed 100% in 4 (3%) of the cases (the largest is 315%, discussed in more detail below).

Figure 2 shows a representative simulation result. In this example, the actual tumor grows from 1629 to 2175 enhancing voxels, i.e., the observed growth, as defined in Equation 6, is Δ^obs = 33.5%. The ensembles produce simulated tumors that range from 1169 to 3495 voxels; the “best” parameter set produces a tumor with 2167 voxels, so the relative error, $\left|R^{\text{best}}\right|$, is approximately 0.37%. Figure 3A shows a histogram of the simulated changes Δ^sim, Equation 7, of the patient’s tumor for each of the 18 parameter sets. The observed change in tumor volume, Δ^obs, is well approximated by several sets of the parameters.

Figure 2

Figure 2. Representative simulation for a typical patient over an interval of 62 days, during which the volume of enhancing tumor grows from 1629 to 2175 voxels (33.5% increase). The axial slices shown encompass the approximate center of mass of the tumor. Top row: Patient scans and manual segmentations. Bottom row: Simulated tumor and computationally derived segmentations. (A) Initial post-contrast T1-weighted MRI. Cyan curves: edema segmentation; red curves: enhancing tumor. (B) As in (A), but for the initial T2-weighted MRI. (C) As in (A), but for the subsequent scan. (D) As in (B), but for the subsequent scan. (E) As in (A), but with the simulation ensemble superimposed. The red curves show the extent of the enhancing region for each of the 18 simulated tumors. (F) As in (E), but superimposed on the computational domain produced by SPM-12. Dark gray: CSF; medium gray: gray matter; light gray: white matter. (G) As in (C). Red curve: segmentation of the simulated tumor that most closely matches the observed number of enhancing voxels. Cyan curve: the corresponding simulated edema segmentation. (H) As in (F), but superimposed on the brain segmentations produced by SPM-12 from the subsequent scan.

Figure 3

Figure 3. (A) Histogram of the simulated change Δ^sim, defined by Equation 7, over all 18 parameter sets used to simulate the tumor for the patient in Figure 2. The vertical dashed line shows the observed change, Δ^obs, as defined by Equation 6. In this case, the volume of enhancing tumor increases by approximately 33%. (B) As in (A), except for the patient simulated in Figure 5. All of the parameter sets significantly overestimate the number of enhancing voxels in the comparison scan, which shrinks by about 78%.

Because the location of the tumor also matters, we define two spatial error measures, as follows. Let R denote the set of voxels in the region segmented as enhancing tumor in the actual patient scan and let S be those of the simulated tumor. We define $\left|S\right|$ and $\left|R\right|$ as the number of voxels in the respective regions. Then $\left|S{\cup }^R\right|$ is the number of voxels occupied by one or the other and $\left|S{\cap }^R\right|$ is the number of voxels in the overlap. We define the containment as

$\begin{array}{ll}C=\frac{\left|S{\cap }^R\right|}{\left|R\right|}.& (9)\end{array}$

and the agreement as

$\begin{array}{ll}A=\frac{\left|S{\cap }^R\right|}{\left|S{\cup }^R\right|}.& (10)\end{array}$

Figure 4 shows an example of the two measures. A perfect model would produce A = C = 1 if the brain domains were exactly co-registered.

Figure 4

Figure 4. Two examples of the agreement (A) defined i Equation 10 and containment (C) defined in Equation 9 measures. (A) If S is twice as large as R but completely overlaps R, then C = 1 but A = 0.5. (B) If S and R are the same size but their overlap is only half of the area of each, then A = C = 0.5.

We compute the agreement measures by superimposing the simulated tumor on the computational domain derived from the subsequent scan. We compare the sets of voxels segmented as enhancing tumor in the subsequent scan with those from the simulated tumors (there are 18 such comparisons, one for each parameter set). The agreement measures range from 0.205 to 0.514 and the containment measures, from 0.408 to 0.799. The parameter set that most closely approximates the actual enhancing tumor volume gives agreement and containment measures of 0.466 and 0.634, respectively, in Figure 2.

Co-registration errors between successive patient scans are inevitable. To help interpret the agreement and containment measures from the simulations, it is useful to compare the corresponding measures between the brain segmentations derived from SPM-12 (cf. panels (f) and (h) in Figure 2). In this case, the respective agreement measures between the regions segmented as gray and white matter are 0.570 and 0.587 and the containment, 0.714 and 0.835.

Figure 5, which is organized identically to Figure 2, shows the results of the simulation with the largest relative error among the 132 that were performed for this study. In this case, the clinically observed volume of enhancing tumor has shrunk considerably, from 2245 to 500 voxels (Δ^obs = −0.773). The 18 simulated tumors (second row) range from 2075 to 7246 enhancing voxels. Figure 3B shows a histogram of all the simulated changes Δ^sim, none of which well approximates the observed one; $\left|R^{\text{best}}\right|$, is approximately 312%. The tumor agreement and containment measures range from approximately 0.02 to 0.18.

Figure 5

Figure 5. Results of the simulation with the largest relative error in enhancing tumor volume. The enhancing volume in the actual patient scans has decreased from 2245 voxels in panel (A) to 500 voxels in panel (C). over an interval of 75 days. The axiel slices shown encompass the approximate center of mass of the tumor. Top row: Patient scans and manual segmentations. Bottom row: Simulated tumor and computationally derived segmentations. (A) Initial post-contrast T1-weighted MRI. (B) As in (A), but for the initial T2-weighted MRI. (C) As in (A), but for the subsequent scan. (D) As in (B), but for the subsequent scan. (E) As in (A), but with the simulation ensemble superimposed. (F) As in (E), but superimposed on the computational domain produced by SPM-12. Dark gray: CSF; medium gray: gray matter; light gray: white matter. (G) As in (C). Red curve: segmentation of the simulated tumor that most closely matches the observed number of enhancing voxels. Cyan curve: the corresponding simulated edema segmentation. (H) As in F(F0, but superimposed on the brain segmentations produced by SPM-12 from the subsequent scan.

4 Discussion

This preliminary study shows the potential utility of simple mathematical models for generating scenarios of treatment response on an individualized patient basis using imaging data that is collected as part of routine clinical surveillance. The experimental design provides reasonably accurate estimates of the likelihood that a patient’s tumor grows or shrinks by a given threshold. At least one choice of parameters from the 18 used to simulate each tumor generates a result whose volume lies within 20% of the observed tumor in 86% of the time intervals studied. As illustrated in Table 3, the simulations match, within a few percentage points, the observed chances that a given patient’s tumor grows or shrinks by specified thresholds by the next scan time. The simulations and visualizations can be run on a laptop computer within a couple of minutes in a clinical setting. This section briefly addresses some related work as well as limitations and potential improvements to the present study.

4.1 Relation to prior work

Many mathematical models have been proposed to approximate the multiscale aspects of the growth of infiltrating gliomas, including glioblastoma (GBM) tumors. Most of them use reaction-diffusion equations (51), similar to the models presented here. Konukoglu et al. (52) proposed a method to estimate “patient-specific” parameters using a sequence of brain images. Hogea et al. (53) proposed a different scheme (with similar objectives) that also attempts to incorporate the tumor mass effect. While not yet standard clinical tools, some models have incorporated angiogenesis and oxygen transport to improve predictions of tumor dynamics (54, 55). Mostaghi-Kashanian et al. (56) proposed a reaction-diffusion model to distinguish visible and invisible tumor regions. Lipková et al. (57) presented a patient-specific model that links glioma growth with pressure distribution in the brain to estimate intracranial hypertension, midline shift, and cognitive impairment.

Several research groups have proposed mathematical modeling frameworks to optimize treatment strategies for GBM. One notable example is the study by Dean et al. (58), which developed a novel radiation therapy schedule based on a mathematical model of cell-state plasticity. Randles et al. (59) applied a computational model of the spatiotemporal dynamics of the perivascular niche to optimize the standard-of-care treatment for GBM. Brüningk et al. (60) proposed a personalized treatment strategy using intermittent radiotherapy based on a mathematical model of tumor growth, radiation response, and a patient-specific evolution of resistance. Our modeling approach is more modest because our data are limited to surveillance MRI.

4.2 Rationale for our scientific approach and potential refinements

The ASU model has three advantages over contemporary machine-learning methods. First, unlike a neural network that may contain millions of parameters that must be fitted, the ASUB model has only seven (Table 2). Second, it is possible to quantify the uncertainties in the output of the ASUB model and prove that it always produces nonnegative solutions (47). (Put another way, the ASUB model cannot “hallucinate,” unlike some machine-learning models.) Finally, only modest computing resources are needed: a full set of 60-day simulations takes less than two minutes on a modern multicore CPU.

The patients in our data set received a variety of treatments, including radiation, chemotherapy, antiangiogenic drugs, and additional surgery. Except for surgery, which reduces tumor cell populations beyond the scope of the model, we have not attempted to account for treatment timelines. The model’s growth and quiescence terms, together with their parameterizations, reflect the sum total of treatment response and resistance. As explained in the introduction, we simply do not have a sufficiently validated understanding of GBM tumor biology or associated patient data to include more detailed mathematical descriptions of the effects of specific treatments. Instead, we lump treatment effects and tumor biology into one set of net diffusion, growth, and quiescence terms. We take the view that the treatment response of recurrent GBM tumors is unpredictable, and, analogously to simulations of portfolio performance for financial planning, have identified ranges of model parameters that yield clinically representative scenarios for patient counseling.

Nevertheless, a future (and larger) study might stratify patients according to selected biomarkers. For example, MGMT (O⁶-methylguanine (O⁶-MeG)-DNA methyltransferase) promoter methylation status (61) and isocitrate dehydrogenase (IDH) mutations (62) can affect treatment outcomes in GBM. Modified ranges of model parameters may yield more clinically relevant results for selected subsets of patients.

4.3 Potential integration with clinical workflows

Once the MR scans from a given patient have been uploaded to the appropriate digital repository, the preprocessing workflow outlined in Section 3.1 can be applied. The ASUB model is not sensitive to small changes in segmentation boundaries; artificial intelligence or other semi-automated method could be used to generate a preliminary segmentation, as long as there are no gross errors in identifying the tumorous regions. The resulting preprocessed data can be downloaded as any other imaging files, and the simulation program itself could be packaged as an “app” on a physician’s laptop. The simulations used to generate outputs like Figure 2 can be run in a couple of minutes.

4.4 Sources of errors and failure analysis

Co-registration errors are unavoidable and, in some cases, significant. One example is in Figure 5. The tumor is simulated using the computational domain in panel (f). However, in panel (g), the cyan curve representing the boundary of the edematous region for the “best” simulated tumor runs through the resection cavity in the subsequent scan. Panel (h) shows the corresponding SPM12-derived brain segmentations. Differences from panel (f) in the size and placement of the ventricles and gray- and white-matter tracts are evident. Measures of tumor agreement and containment must be interpreted in light of such geometric differences in the reconstructed domains.

The effect of co-registration errors can be especially important in the vicinity of the resection cavity, and tumor-related changes to the brain geometry may be a contributing factor. Our model simulations do not account for mass effect, and a future effort may benefit from including one (but would involve significant additional computational complexity).

Although T1 and T2 sequences were used for each patient, no other uniform imaging protocol has been followed. The data used in this research project were collected on different scanners. The strength of the magnetic field varies, and there are differences in the way that operators positioned the patients. Such variations complicate image comparisons, but they also reflect “real world” clinical data. We believe that our approach will be practicable in contemporary routine clinical settings where uniform, prospective scanning protocols may not exist.

The ASUB model and the FKPP model presume a fixed carrying capacity for tumor cells. When initialized with cell densities below the carrying capacity, neither model can simulate a decreasing cell population. When proliferating cells become sufficiently numerous in the ASUB model, they join the quiescent pool. Sufficiently dense regions of quiescent cells are assumed not to enhance on T1 imaging (they are assumed to form a hypointense region). Thus, our existing model parameterizations cannot always simulate a rapid reduction in the number of enhancing voxels. There are 4 cases in which $R^{\text{best}}$ >100%, of which Figure 5 is one, and each involves a significant (> 66%) reduction in the number of enhancing voxels between the initial and subsequent scans. This result may simply reflect an inherent limitation of the ASUB model. Future work (and a larger dataset) will be needed to get a better idea of how often simulation failures like that in Figure 5 are likely to occur and to understand the issues in model initialization and parameterization that lead to such failures.

Our methods for initializing tumor cell populations and performing comparisons with patient scans are simplistic, and the tumor segmentations are based solely on relative pixel intensity. Some patients in our data set were treated with anti-angiogenic agents, which may affect the appearance of tumor on MRI (63) (“pseudoresponse”). Radiotherapy can cause “pseudoprogression,” where dying cells cause imaging enhancement that can be difficult to distinguish from actively growing tumor cells (64, 65). The addition of perfusion MRI may improve our characterization of treatment effects and provide better initializations of the model.

Finally, it may be possible to improve the agreement and containment measures with a more sophisticated model of cell motility than isotropic diffusion. Diffusion tensor imaging (DTI)—either of individual patients or by applying an averaged DTI map over many patients—may be useful in reaction-diffusion models of glioma progression (66). Analogous approaches involve viscous stress tensors (67). In addition, the traveling-wave characteristics of the Fisher and ASU-Barrow models may allow better initialization of tumor cell populations in edematous regions (48, 52).

Data availability statement

The datasets for this article are not publicly available due to concerns regarding participant/patient anonymity. Requests to access the datasets should be directed to the corresponding authors.

Ethics statement

The studies involving humans were approved by St. Joseph’s Hospital and Medical Center Institutional Review Board. The studies were conducted in accordance with the local legislation and institutional requirements. Written informed consent for participation was not required from the participants or the participants’ legal guardians/next of kin in accordance with the national legislation and institutional requirements.

Author contributions

EK: Conceptualization, Software, Funding acquisition, Writing – review & editing, Supervision, Writing – original draft, Methodology, Formal analysis. YX: Writing – review & editing, Formal analysis, Software, Data curation. CC-V: Data curation, Writing – review & editing. DH: Software, Writing – review & editing, Conceptualization, Formal analysis, Methodology. OA-G: Data curation, Writing – review & editing. GG-C: Data curation, Writing – review & editing. TO: Writing – review & editing, Data curation. RD: Resources, Project administration, Writing – review & editing. YK: Methodology, Writing – review & editing, Conceptualization. MP: Writing – review & editing, Supervision, Conceptualization, Funding acquisition.

Funding

The author(s) declare financial support was received for the research and/or publication of this article. This work was funded in part by the Arizona Biomedical Research Center under grant number ADHS16-162514 and by funds from the Newsome Chair for Neurosurgery Research held by Mark C. Preul.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

The author(s) declared that they were an editorial board member of Frontiers, at the time of submission. This had no impact on the peer review process and the final decision.

Generative AI statement

The author(s) declare that no Generative AI was used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

1. Araujo RP and McElwain DLS. A history of the study of solid tumour growth: The contribution of mathematical modeling. Bull Math Biol. (2004) 66:1039–91. doi: 10.1016/j.bulm.2003.11.002

PubMed Abstract | Crossref Full Text | Google Scholar

2. Frenzen CL and Murray JD. A cell kinetics justification for Gompertz’ equation. SIAM J Appl Math. (1986) 46:514–629. doi: 10.1137/0146042

Crossref Full Text | Google Scholar

3. Butner JD, Dogra P, Chung C, Pasqualini R, Arap W, Lowengrub J, et al. Mathematical modeling of cancer immunotherapy for personalized clinical translation. Nat Comput Sci. (2022) 2:785–96. doi: 10.1038/s43588-022-00377-z

PubMed Abstract | Crossref Full Text | Google Scholar

4. Mathur D, Barnett E, Scher HI, and Xavier JB. Optimizing the future: how mathematical models inform treatment schedules for cancer. Trends Cancer. (2022) 8:506–16. doi: 10.1016/j.trecan.2022.02.005

PubMed Abstract | Crossref Full Text | Google Scholar

5. Galappaththige S, Gray RA, Costa CM, Niederer S, and Pathmanathan P. Credibility assessment of patientspecific computational modeling using patient-specific cardiac modeling as an exemplar. PloS Comput Biol. (2022) 18:e1010541. doi: 10.1371/journal.pcbi.1010541

PubMed Abstract | Crossref Full Text | Google Scholar

6. Pathmanathan P, Aycock K, Badal A, Bighamian R, Bodner J, Craven BA, et al. Credibility assessment of in silico clinical trials for medical devices. PloS Comput Biol. (2024) 20:e1012289. doi: 10.1371/journal.pcbi.1012289

PubMed Abstract | Crossref Full Text | Google Scholar

7. Hormuth DA II, Weis JA, Barnes SL, Miga MI, Rericha EC, Quaranta V, et al. A mechanically coupled reaction-diffusion model that incorporates intra-tumoural heterogeneity to predict in vivo glioma growth. J R Soc Interface. (2017) 14:20161010. doi: 10.1098/rsif.2016.1010

PubMed Abstract | Crossref Full Text | Google Scholar

8. Hormuth DA, Al Reghali KA, Elliott AM, Yankeelov TE, and Chung C. Image-based personalization of computational models for predicting response of high-grade glioma to chemoradiation. Sci Rep. (2021) 11:8520. doi: 10.1038/s41598-021-87887-4

PubMed Abstract | Crossref Full Text | Google Scholar

9. Stupp R, Hegi ME, Mason WP, van den Bent MJ, Taphoorn MJ, and Janzer RC. Effects of radiotherapy with concomitant and adjuvant temozolomide versus radiotherapy alone on survival in glioblastoma in a randomized phase III study: 5-year analysis of the EORTC-NCIC trial. Lancet Oncol. (2009) 10:459–66. doi: 10.1016/S1470-2045(09)70025-7

PubMed Abstract | Crossref Full Text | Google Scholar

10. Stupp R, Mason WP, van den Bent MJ, Weller M, Fisher B, and Taphoorn MJ. Radiotherapy plus concomitant and adjuvant temozolomide for glioblastoma. N Engl J Med. (2005) 352:987–96. doi: 10.1056/NEJMoa043330

PubMed Abstract | Crossref Full Text | Google Scholar

11. Weller M, Cloughesy T, Perry JR, and Wick W. Standards of care for treatment of recurrent glioblastoma— are we there yet? Neuro Oncol. (2013) 15:4–27. doi: 10.1093/neuonc/nos273

PubMed Abstract | Crossref Full Text | Google Scholar

12. Tosoni A, Franceschi E, Poggi R, and Brandes AA. Relapsed glioblastoma: Treatment strategies for initial and subsequent recurrences. Curr Treat Options Oncol. (2016) 17:49. doi: 10.1007/s11864-016-0422-4

PubMed Abstract | Crossref Full Text | Google Scholar

13. Kotecha R, Odia Y, Khosla AA, and Ahluwalia MS. Key clinical principles in the management of glioblastoma. JCO Oncol Pract. (2023) 19:180–9. doi: 10.1200/OP.22.00476

PubMed Abstract | Crossref Full Text | Google Scholar

14. Zhao YH, Wang ZF, Pan ZY, Péus D, Delgado-Fernandez J, Pallud J, et al. A meta-analysis of survival outcomes following reoperation in recurrent glioblastoma: Time to consider the timing of reoperation. Front Neurol. (2019) 286. doi: 10.3389/fneur.2019.00028

PubMed Abstract | Crossref Full Text | Google Scholar

15. Pineda E, Domenech M, Hernández A, Comas S, and Balaña C. Recurrent glioblastoma: Ongoing clinical challenges and future prospects. Onco Target Ther. (2023) 16:71–86. doi: 10.2147/OTT.S366371

PubMed Abstract | Crossref Full Text | Google Scholar

16. Thomas-Joulié A, Tran S, El Houari L, Seyve A, Bielle F, Birzu C, et al. Prognosis of glioblastoma patients improves significantly over time interrogating historical controls. Eur J Cancer. (2024) 202:114004. doi: 10.1016/j.ejca.2024.114004

PubMed Abstract | Crossref Full Text | Google Scholar

17. Musella A, DeVitto R, Anthony M, and Mydland DE. The importance of shared decision-making for patients with glioblastoma. Patient Prefer Adher. (2021) 15:2009–16. doi: 10.2147/PPA.S314792

PubMed Abstract | Crossref Full Text | Google Scholar

18. Sareen H, Ma Y, Becker TM, Roberts TL, de Souza P, and Powter B. Molecular biomarkers in glioblastoma: A systematic review and meta-analysis. Int J Mol Sci. (2022) 23:8835. doi: 10.3390/ijms23168835

PubMed Abstract | Crossref Full Text | Google Scholar

19. Pacheco C, Martins C, Monteiro J, Baltazar F, Costa BM, and Sarmento B. Glioblastoma vasculature: From its critical role in tumor survival to relevant in vitro modelling. Front Drug Deliv. (2022) 2:823412. doi: 10.3389/fddev.2022.823412

Crossref Full Text | Google Scholar

20. Feldman L. Hypoxia within the glioblastoma tumor microenvironment: a master saboteur of novel treatments. Front Immunol. (2024) 15:1384249. doi: 10.3389/fimmu.2024.1384249

PubMed Abstract | Crossref Full Text | Google Scholar

21. Jain P, Vashist S, and Panjiyar BK. Navigating the immune challenge in glioblastoma: exploring immunotherapeutic avenues for overcoming immune suppression. Cureus. (2023) 9:e46089. doi: 10.7759/cureus.46089

PubMed Abstract | Crossref Full Text | Google Scholar

22. Lu C, Kang T, Zhang J, Yang K, Liu Y, Song K, et al. Combined targeting of glioblastoma stem cells of different cellular states disrupts Malignant progression. Nat Commun. (2025) 16:2974. doi: 10.1038/s41467-025-58366-5

PubMed Abstract | Crossref Full Text | Google Scholar

23. Mohiuddin E and Wakimoto H. Extracellular matrix in glioblastoma: opportunities for emerging therapeutic approaches. Am J Cancer Res. (2021) 11:3742–54.

PubMed Abstract | Google Scholar

24. Glasser JW and Feng Z. Mechanistic models are hyotheses: A perspective. Math Biosci. (2025) 383:109419. doi: 10.1016/j.mbs.2025.109419

PubMed Abstract | Crossref Full Text | Google Scholar

25. Yankeelov TL, Atuegwu N, Hormuth D, Weis JA, Barnes SL, Miga MI, et al. Clinically relevant modeling of tumor growth and treatment response. Sci Transl Med. (2013) 5:187ps9. doi: 10.1126/scitranslmed.3005686

PubMed Abstract | Crossref Full Text | Google Scholar

26. Harris DC, Magnucci-Jiménez G, Xu Y, Eikenberry SE, Quarles CC, Preul MC, et al. Tracking glioblastoma progression after initial resection with minimal reaction-diffusion models. Math Biosci Eng. (2022) 19:5446–81. doi: 10.3934/mbe.2022256

PubMed Abstract | Crossref Full Text | Google Scholar

27. Neuroimaging Informatics Technology Initiative. (2020). Available online at: nifti.nimh.nih.gov (Accessed 2024-12-18).

Google Scholar

28. of Imaging Neuroscience University College London FILWD (2024). Available online at: www.fil.ion.ucl.ac.uk/spm/software/spm12/ (Accessed 2024-12-18). SPM12.

Google Scholar

29. 3D Slicer image computing platform (2024). Available online at: www.slicer.org (Accessed 2024-12-18).

Google Scholar

30. Ellingson BM, Harris RJ, Woodworth DC, Leu K, Zaw O, Mason WP, et al. Baseline pretreatment contrast enhancing tumor volume including central necrosis is a prognostic factor in recurrent glioblastoma: evidence from single and multicenter trials. Neuro Oncol. (2017) 19:89–98. doi: 10.1093/neuonc/now187

PubMed Abstract | Crossref Full Text | Google Scholar

31. Garrett MD, Yanagihara TK, Yeh R, McKhann GM, Sisti MB, Bruce JN, et al. Monitoring radiation treatment effects in glioblastoma: FLAIR volume as significant predictor of survival. Tomorgraphy. (2017) 3:131–7. doi: 10.18383/j.tom.2017.00009

PubMed Abstract | Crossref Full Text | Google Scholar

32. Johnstad C, Reinertsen I, Thurin E, Dunås T, Bouget D, Sagberg LM, et al. The prognostic importance of glioblastoma size and shape. Acta Neurochir. (2024) 166:450. doi: 10.1007/s00701-024-06351-0

PubMed Abstract | Crossref Full Text | Google Scholar

33. Egger J, Kapur T, Fedorov A, Pieper S, Miller JV, Veeraraghavan H, et al. GBM volumetry using the 3D slicer medical image computing platform. Sci Rep. (2013) 3:1364. doi: 10.1038/srep01364

PubMed Abstract | Crossref Full Text | Google Scholar

34. Fisher RA. The wave of advance of advantageous genes. Ann Eugenics. (1937) 7:355–69. doi: 10.1111/j.1469-1809.1937.tb02153.x

Crossref Full Text | Google Scholar

35. Kolmogorov AN, Petrovsky IG, and Piscunov NS. Investigation of a diffusion equation connected to the growth of materials, and application to a problem in biology. Bull Univ Moscow. (1937) 1:1–26.

Google Scholar

36. Burgess PK, Kulesa PM, Murray JD, and Alvord EC Jr. The interaction of growth rates and diffusion coefficients in a three-dimensional mathematical model of gliomas. Neuropathol Exp Neurol. (1997) 56:704–13. doi: 10.1097/00005072-199706000-00008

PubMed Abstract | Crossref Full Text | Google Scholar

37. Gerlee P and Nelander S. Traveling wave analysis of a mathematical model of glioblastoma growth. Math Biosci. (2016) 276:75–81. doi: 10.1016/j.mbs.2016.03.004

PubMed Abstract | Crossref Full Text | Google Scholar

38. Woodward DE, Cook J, Tracqui P, Cruywagen GC, Murray JD, and Alvord EC Jr. A mathematical model of glioma growth: The effect of extent of surgical resection. Cell Prolif. (1996) 29:269–88. doi: 10.1111/j.1365-2184.1996.tb01580.x

PubMed Abstract | Crossref Full Text | Google Scholar

39. Swanson KR, Alvord EC Jr., and Murray JD. Quantifying efficacy of chemotherapy of brain tumors (gliomas) with homogeneous and heterogeneous drug delivery. Acta Biotheor. (2002) 50:223–37. doi: 10.1023/A:1022644031905

PubMed Abstract | Crossref Full Text | Google Scholar

40. Murray JD. Glioblastoma brain tumors: Estimating the time from brain tumor initiation and resolution of a patient survival anomaly after similar treatment protocols. J Biol Dynam. (2012) 6:118–27. doi: 10.1080/17513758.2012.678392

PubMed Abstract | Crossref Full Text | Google Scholar

41. Eidel O, Burth S, Neumann JO, Kieslich PJ, Sahm F, Jungk C, et al. Tumor infiltration in enhancing and non-enhancing parts of glioblastoma: A correlation with histopathology. PloS One. (2017) 12:e0169292. doi: 10.1371/journal.pone.0169292

PubMed Abstract | Crossref Full Text | Google Scholar

42. Shampine LF, Sommeijer BP, and Verwer JG. IRKC: An IMEX solver for stiff diffusion-reaction PDEs. J Comput Appl Math. (2006) 196:485–97. doi: 10.1016/j.cam.2005.09.014

Crossref Full Text | Google Scholar

43. Feucht D, Haas P, Skardelly M, Behling F, Rieger D, Bombach P, et al. Preoperative growth dynamics of untreated glioblastoma: Description of an exponential growth type, correlating factors, and association with postoperative survival. Neuro Oncol Adv. (2024) 6:vdae053. doi: 10.1093/noajnl/vdae053

PubMed Abstract | Crossref Full Text | Google Scholar

44. Stoyanov GS, Lyutfi E, Georgiev R, Dzhenkov DL, and Kaprelyan A. The rapid development of glioblastoma: A report of two cases. Cureus. (2022) 14:e26319. doi: 10.7759/cureus.26319

PubMed Abstract | Crossref Full Text | Google Scholar

45. Wechsler-Jentzsch K, Witt JH, Friz CR, McCullough DC, and Harisiadis L. Unresectable gliomas in children: tumor-volume response to radiation therapy. J Radiology. (1988) 169:237–41. doi: 10.1148/radiology.169.1.3420265

PubMed Abstract | Crossref Full Text | Google Scholar

46. Klank RI, Rosenfeld SS, and Odde DJ. A Brownian dynamics tumor progression simulator with application to glioblastoma. Converg Sci Phys Oncol. (2018) 4:015001. doi: 10.1088/2057-1739/aa9e6e

PubMed Abstract | Crossref Full Text | Google Scholar

47. Harris DC, He C, Preul MC, Kostelich EJ, and Kuang Y. Critical patch size of a two-population reaction diffusion model describing brain tumor growth. SIAM J Appl Math. (2024) 84:S249–68. doi: 10.1137/22M1509631

Crossref Full Text | Google Scholar

48. Häger W, Lazzzeroni M, Astaraki M, and Toma-Daşu I. CTV delineation for high-grade gliomas: Is there agreement with tumor cell invasion models? Adv Radiat Oncol. (2022) 7:100987. doi: 10.1016/j.adro.2022.100987

PubMed Abstract | Crossref Full Text | Google Scholar

49. Chang PD, Malone HR, Bowden SG, Chow DS, Gill BJA, and Ung TH. A multiparametric model for mapping cellularity in glioblastoma using radiographically localized biopsies. Am J Neuroradiol. (2017) 38:890–8. doi: 10.3174/ajnr.A5112

PubMed Abstract | Crossref Full Text | Google Scholar

50. Kacker RN, Lagergren ES, and Filliben JJ. Taguchi’s orthogonal arrays are classical designs of experiments. J Res Natl Inst Stand Technol. (1991) 96:577–91. doi: 10.6028/jres.096.034

PubMed Abstract | Crossref Full Text | Google Scholar

51. Roniotis A, Marias K, Sakkalis V, and Zervakis M. Diffusive modeling of glioma evolution: A review. J BioMed Sci Eng. (2010) 3:501–8. doi: 10.4236/jbise.2010.35070

Crossref Full Text | Google Scholar

52. Konukoglu E, Clatz O, Menze BH, Stieltjes B, Weber MA, Mandonnet E, et al. Image guided personalization of reaction-diffusion type tumor growth models using modified anisotropic Eikonal equations. IEEE Trans Med Imaging. (2010) 29:77–95. doi: 10.1109/TMI.2009.2026413

PubMed Abstract | Crossref Full Text | Google Scholar

53. Hogea C, Davatzikos C, and Biros G. An image-driven parameter estimation problem for a reactiondiffusion glioma growth model with mass effects. J Math Biol. (2008) 56:793–825. doi: 10.1007/s00285-007-0139-x

PubMed Abstract | Crossref Full Text | Google Scholar

54. Gandolfi A, Franciscis S, d’Onofrio A, Fasano A, and Sinisgalli C. Angiogenesis and vessel co-option in a mathematical model of diffusive tumor growth: The role of chemotaxis. J Theor Biol. (2021) 512:110526. doi: 10.1016/j.jtbi.2020.110526

PubMed Abstract | Crossref Full Text | Google Scholar

55. Jafariandehkordi A and Daneshmehr A. A biomechanical simulation of the brain glioma growth considering the effects of oxygen Bunsen solubility coefficient. Mech Res Commun. (2020) 106:103539. doi: 10.1016/j.mechrescom.2020.103539

Crossref Full Text | Google Scholar

56. Moshtaghi-Kashanian N, Niroomand-Osculi H, and Meghdadi N. Simulating glioblastoma growth consisting both visible and invisible parts of the tumor using a diffusion-reaction model followed by resection and radiotherapy. Acta Neurol Belg. (2020) 120:629–37. doi: 10.1007/s13760-018-0952-6

PubMed Abstract | Crossref Full Text | Google Scholar

57. Lipková J, Menze B, Wiestler B, Koumoutsakos P, and Lowengrub JS. Modelling glioma progression, mass effect and intracranial pressure in patient anatomy. J R Soc Interface. (2022) 19:20210922. doi: 10.1098/rsif.2021.0922

PubMed Abstract | Crossref Full Text | Google Scholar

58. Dean JA, Tanguturi SK, Cagney D, Shin KY, Youssef G, Aizer A, et al. Phase I study of a novel glioblastoma radiation therapy schedule exploiting cell-state plasticity. Neuro-Oncology. (2023) 25:110012. doi: 10.1093/neuonc/noac253

PubMed Abstract | Crossref Full Text | Google Scholar

59. Randles A, Wirsching HG, Dean JA, Cheng YK, Emerson S, Pattwell SS, et al. Computational modelling of perivascular-niche dynamics for the optimization of treatment schedules for glioblastoma. Nat BioMed Eng. (2021) 5:346–59. doi: 10.1038/s41551-021-00710-3

PubMed Abstract | Crossref Full Text | Google Scholar

60. Brüningk SC, Peacock J, Whelan CJ, Yu HHM, Sahebjam S, and Enderling H. Intermittent radiotherapy as alternative treatment for recurrent high grade glioma: A modelling study based on longitudinal tumor measurements. Sci Rep. (2021) 11:20219. doi: 10.1038/s41598-021-99507-2

PubMed Abstract | Crossref Full Text | Google Scholar

61. Butler M, Pongor L, Su YT, Xi L, Raffeld M, Quezado M, et al. MGMT status as a clinical biomarker in glioblastoma. Trends Cancer. (2020) 6:380–91. doi: 10.1016/j.trecan.2020.02.010

PubMed Abstract | Crossref Full Text | Google Scholar

62. Nasany RA and de la Fuente MI. Therapies for IDH-mutant gliomas. Curr Neurol Neurosci Rep. (2023) 23:255–33. doi: 10.1007/s11910-023-01265-3

PubMed Abstract | Crossref Full Text | Google Scholar

63. Arevalo OD, Soto C, Rabiei P, Kamali A, Ballester LY, Esquenazi Y, et al. Assessment of glioblastoma response in the era of bevacizumab: Longstanding and emergent challenges in the imaging evaluation of pseudoresponse. Front Neurol. (2019) 10:460. doi: 10.3389/fneur.2019.00460

PubMed Abstract | Crossref Full Text | Google Scholar

64. Kamiya-Matsuoka C and Gilbert MR. Treating recurrent glioblastoma: an update. CNS Oncol. (2015) 4:91104. doi: 10.2217/cns.14.55

PubMed Abstract | Crossref Full Text | Google Scholar

65. McKenney AS, Weg E, Bale TA, Wild AT, Um H, Fox MJ, et al. Radiomic analysis to predict histopathologically confirmed pseudoprogression in glioblastoma patients. Adv Radiat Oncol. (2023) 8:100916. doi: 10.1016/j.adro.2022.100916

PubMed Abstract | Crossref Full Text | Google Scholar

66. Stretton E, Geremia E, Menze B, Delingette H, and Ayache N. (2013). Importance of patient DTI’s to accurately model glioma growth using the reaction diffusion equation, in: 2013 IEEE 10th International Symposium on Biomedical Imaging (San Francisco, CA: IEEE). pp. 1142–5.

Google Scholar

67. Yuan J and Liu L. Brain glioma growth model using reaction-diffusion euqation with viscous stress tensor on brain MR images. Magn Reson Imaging. (2016) 34:114–9. doi: 10.1016/j.mri.2015.10.012

PubMed Abstract | Crossref Full Text | Google Scholar