Modeling Tumor-Associated Edema in Gliomas during Anti-Angiogenic Therapy and Its Impact on Imageable Tumor

Glioblastoma, the most aggressive form of primary brain tumor, is predominantly assessed with gadolinium-enhanced T1-weighted (T1Gd) and T2-weighted magnetic resonance imaging (MRI). Pixel intensity enhancement on the T1Gd image is understood to correspond to the gadolinium contrast agent leaking from the tumor-induced neovasculature, while hyperintensity on the T2/FLAIR images corresponds with edema and infiltrated tumor cells. None of these modalities directly show tumor cells; rather, they capture abnormalities in the microenvironment caused by the presence of tumor cells. Thus, assessing disease response after treatments impacting the microenvironment remains challenging through the obscuring lens of MR imaging. Anti-angiogenic therapies have been used in the treatment of gliomas with spurious results ranging from no apparent response to significant imaging improvement with the potential for extremely diffuse patterns of tumor recurrence on imaging and autopsy. Anti-angiogenic treatment normalizes the vasculature, effectively decreasing vessel permeability and thus reducing tumor-induced edema, drastically altering T2-weighted MRI. We extend a previously developed mathematical model of glioma growth to explicitly incorporate edema formation allowing us to directly characterize and potentially predict the effects of anti-angiogenics on imageable tumor growth. A comparison of simulated glioma growth and imaging enhancement with and without bevacizumab supports the current understanding that anti-angiogenic treatment can serve as a surrogate for steroids and the clinically driven hypothesis that anti-angiogenic treatment may not have any significant effect on the growth dynamics of the overall tumor cell populations. However, the simulations do illustrate a potentially large impact on the level of edematous extracellular fluid, and thus on what would be imageable on T2/FLAIR MR. Additionally, by evaluating virtual tumors with varying growth kinetics, we see tumors with lower proliferation rates will have the most reduction in swelling from such treatments.

Introduction

Glioblastoma Multiforme (GBM) is a highly aggressive and invasive primary brain tumor. The standard treatment protocol is to surgically remove as much of the tumor as is reasonably safe, followed by a combination of chemotherapy with radiation. Despite aggressive treatment, the prognosis remains poor with a median survival time of 14 months (Stupp et al., 2005). The inability to accurately determine the extent of diffuse tumor cell infiltration of the normal brain affects the ability to assess response to treatment through clinical imaging, confounding clinical progress. Currently, clinicians rely primarily on three magnetic resonance imaging (MRI) modalities to monitor the development of the tumor, the T2 weighting, FLAIR, and T1 weighting with gadolinium contrast enhancement (T1Gd) sequences illustrated in Figure 1. However, it is known that none of these sequences are able to show the entire extent of the malignant cells (Silbergeld and Chicoine, 1997), since the abnormal regions highlighted in the MR images are as dependent on the microenvironment around the disease, particularly the vasculature, as on the tumor cells themselves.

FIGURE 1

Figure 1. Illustration of primary imaging modalities. The T1Gd image will show enhancement where the contrast agent has been able to diffuse into the extracellular space where the blood brain barrier has been compromised due to tumor-induced neo-angiogenesis. The T2-weighted and FLAIR images are associated with edema or swelling; FLAIR is different from T2 in that the signal from the cerebral spinal fluid (CSF) is removed. In the case of GBM, the non-CSF T2/FLAIR enhancement is primarily vasogenic edema, defined as fluid originating from blood vessels that accumulates around cells (Marmarou, 2007). The fluid leaves the vessels due to pressure and osmotic gradients induced by the breakdown in the blood brain barrier.

Increased vasculature is a primary hallmark of GBM, and while angiogenesis is a hallmark of cancer in general, there are some important factors that separate GBMs from other tumors. First, the glioma cells inhabit an organ that is highly vascularized in its native state. Second, glioma cells are diffusely invasive and are known to co-opt the existing vasculature and migrate and grow along the vessels (Holash et al., 1999; Leenders et al., 2002). Nevertheless, GBMs can form hypoxic regions, often leading to regions of necrosis, and thus, downstream of this hypoxic signaling, emit an abnormally large amount of angiogenic factors such as vascular endothelial growth factor (VEGF) for the recruitment of additional vasculature, analogous to observations in solid tumors (Kerbel, 2000). This process results in the vasculature developing abnormally large vessel radial sizes and, unique to the brain, results in a breakdown of the blood brain barrier in the tumor region.

Features Characterizing MR Imaging Observation

In the case of GBM, the enhancing abnormalities on all of the primary MR imaging modalities, T1Gd, T2, and FLAIR primarily result from a compromised blood brain barrier. The T1Gd image signal is enhanced where the contrast agent has been able to leak into the extracellular space through breakdowns in the blood brain barrier due to tumor-induced neo-angiogenesis. The T2-weighted and FLAIR images show edema or swelling; FLAIR is different from T2 in that the signal from the cerebral spinal fluid (CSF) is inverted. In the case of GBM, the non-CSF T2/FLAIR hyperintense signal is primarily vasogenic edema, defined as fluid originating from blood vessels that accumulates around cells (Marmarou, 2007). The fluid leaves the vessels due to pressure and osmotic gradients induced by the breakdown in the blood brain barrier. Thus, whilst changes in any primary imaging modality (T1Gd or T2/FLAIR) are often interpreted as corresponding to changes in tumor cell density, they may be artifacts of MR imaging.

The Role of Anti-Angiogenics in GBM

The concept of anti-angiogenic treatment for cancer has been popular ever since the landmark paper by Folkman (1971) stating that malignant tumors were angiogenesis-dependent and has been used with some success for other solid tumors in combination with chemotherapy (Hurwitz et al., 2004; Sandler et al., 2006). Since a defining hallmark of GBM is increased vasculature through endothelial cell proliferation (Louis et al., 2007), this disease seems like an obvious candidate for vascular targeting treatment. However, the differences between the vasculature in GBMs and other solid tumors produce different treatment effects. In GBMs, one of the effects of anti-angiogenic treatments is to, at least transiently, repair the blood brain barrier and allow the vessels to return to their normal radial size, increasing their efficiency (Jain, 2005; Batchelor et al., 2007) – referred to as vascular normalization. Ostensibly, this improved efficiency of the vasculature is not the desired impact, though it may help in delivery of other therapeutic agents.

More concerning, however, is that this normalization may directly impact the efficacy of the MR imaging. It is possible for glioma patients with enhancing lesions on T1Gd and T2/FLAIR imaging to have decreased enhancement within a day of anti-angiogenic treatment (Batchelor et al., 2007; Norden et al., 2008), as illustrated by patients 1 and 2 in Figure 2, but upon stopping treatment, the imageable lesion is even larger and more disperse than before (Iwamoto et al., 2009). However, responses are varied; a patient may see no deflection in growth but faster growth after treatment, such as the third patient in Figure 2, or see stabilized disease returning to the previous growth patterns after treatment, as illustrated by the fourth patient in Figure 2. These patients were consented to this study with approval by the local institutional review board at either the University of Washington or the University of California, Los Angeles, and their relevant demographic and therapeutic information is given in Table 1. These conundrums have led to two hypotheses: first, anti-angiogenic treatment has minimal cytotoxic effect but does influence the imaging so that the tumor cannot be effectively visualized, and second, that the treatment may be selecting for a more aggressively invasive phenotype (Verhoeff et al., 2009; Keunen et al., 2011).

FIGURE 2

Figure 2. Four patients with varying imageable responses to anti-angiogenic treatment. Treatment period indicated with gray box on radius plots. (A) Patient 1: A 48-year-old male with Grade III glioma is seen to have significant reduction of enhancing lesion during treatment but recurs almost immediately after treatment is stopped, (B) Patient 2: a 55-year-old male with GBM initially responds to treatment but even while being treated the enhancing region is seen to enlarge again, (C) Patient 3: a 61-year-old male with GBM seems to have no response to treatment and the enhancing region seems to grow faster after treatment, and (D) Patient 4: a 66-year-old female with GBM has imaging stabilized during treatment, but enhancing region begins growing again once treatment is stopped.

TABLE 1

Table 1. Demographic and treatment information corresponding to patients in Figure 2.

Previous studies have shown the use of anti-angiogenic drugs, specifically bevacizumab (Avastin), tends to increase progression-free survival and reduce symptoms of recurrent GBMs, but they have failed to consistently show a significant increase in overall survival and there is concern that the measured radiographic responses do not reflect changes in tumor cell counts (Verhoeff et al., 2009; Deming, 2012). Additionally, animal studies have revealed that treatment with anti-angiogenic drugs may be creating an environment favorable for local invasion and metastasis (Ebos et al., 2009; Pàez-Ribes et al., 2009). Even though there is a lack of significant evidence for increased overall survival after bevacizumab, and it is possible that treatment selects for a more invasive phenotype, the increase in quality of life for some of the patients, due to the relief from edema-related symptoms, means anti-angiogenic therapy is an attractive and relevant treatment option. However, the inability to determine a priori which patients will receive more benefit than harm from anti-angiogenic therapy ultimately keeps clinicians wary (Deming, 2012).

In this paper, we aim to illustrate how a previous mathematical model of glioma growth can be extended to explicitly incorporate edema formation allowing us to directly characterize and potentially predict the effects of anti-angiogenics on imageable tumor growth. The ultimate goal of this model is to help the treatment planning process by identifying exactly those patients that would receive the most benefit from anti-angiogenic treatment.

Materials and Methods

The Proliferation-Invasion-Hypoxia-Necrosis-Angiogenesis-Edema Model

Over the last decade we have made a significant effort toward the development of patient specific mathematical models of GBM that are able to capture the growth kinetics of individual patients (Swanson, 1999; Swanson et al., 2000, 2002a,b; Harpold et al., 2007; Rockne et al., 2010). The simplest form of the model, referred to as the Proliferation-Invasion (PI) model is based on patient specific net rates of proliferation and invasion and has been successful in predicting untreated growth rates for individual patients (Harpold et al., 2007) and providing predictions of outcomes following surgical resections (Swanson et al., 2003), chemotherapy (Swanson et al., 2002a,b, 2003), and radiation (Rockne et al., 2010), while also providing insight into glioma ontogeny (Bohman et al., 2010).

A Mathematical Model of the Angiogenic Cascade in Glioblastoma

The Proliferation-Invasion-Hypoxia-Necrosis-Angiogenesis (PIHNA) model first discussed in Swanson et al. (2011) incorporates the angiogenic cascade and characterizes malignant gliomas with relative proportions of well-oxygenated “normoxic” tumor cells, (c), poorly oxygenated hypoxic tumor cells, (h), necrotic cells, (n), and vascular, or endothelial cells, (v), along with a generic population of angiogenic factors, (a) (Swanson et al., 2011). In words, it assumes the level of nutrients present in the local microenvironment, as inferred from the number of vasculature cells, determines whether the present tumor cells will exhibit normoxic or hypoxic phenotypes. That is, if there is a sufficient level of nutrients present, the cells will remain normoxic, but if the nutrient level falls below a given threshold, the cells will become hypoxic. If the nutrients provided by the vasculature fall below an even lower threshold value, the hypoxic cells will undergo necrosis, at a rate of α_h (1/year) and remain in the necrotic cell population. Normoxic tumor cells are allowed to move (invade) and divide while, due to restricted amounts of nutrients, the hypoxic cells are only allowed to move. The hypoxic cells produce a large amount of angiogenic factors which ultimately cause an increase in the number of vasculature cells. The system is described with a mathematical model composed of the following five coupled reaction-diffusion equations:

$\begin{array}{l}\frac{\partial c}{\partial t}=\stackrel{\begin{array}{c}\text{Net\hspace{0.17em}dispersal\hspace{0.17em}of\hspace{0.17em}normoxic}\\ \text{glioma\hspace{0.17em}cells}\end{array}}{\overbrace{\nabla \cdot \text{(}D(x)(1-T)\nabla c)}}+\stackrel{\begin{array}{c}\text{Net\hspace{0.17em}proliferation\hspace{0.17em}of\hspace{0.17em}normoxic}\\ \text{glioma\hspace{0.17em}cells}\end{array}}{\overbrace{\rho c(1-T)}}\\ +\stackrel{\begin{array}{c}\text{Conversion\hspace{0.17em}of\hspace{0.17em}hypoxic}\\ \text{to\hspace{0.17em}normoxic}\end{array}}{\overbrace{\gamma hV}}-\stackrel{\begin{array}{c}\text{Conversion\hspace{0.17em}of\hspace{0.17em}normoxic}\\ \text{to\hspace{0.17em}hypoxic}\end{array}}{\overbrace{\beta c(1-V)}}\\ -\stackrel{\begin{array}{c}\text{Conversion\hspace{0.17em}of\hspace{0.17em}normoxic}\\ \text{to\hspace{0.17em}necrotic}\end{array}}{\overbrace{{\alpha }_{n}nc}}\\ \frac{\partial h}{\partial t}=\stackrel{\begin{array}{c}\text{Dispersal\hspace{0.17em}of\hspace{0.17em}hypoxic}\\ \text{glioma\hspace{0.17em}cells}\end{array}}{\overbrace{\nabla \cdot (D(x)(1-T)\nabla h)}}-\stackrel{\begin{array}{c}\text{Conversion\hspace{0.17em}of\hspace{0.17em}hypoxic}\\ \text{to\hspace{0.17em}normoxic}\end{array}}{\overbrace{\gamma hV}}\\ +\stackrel{\begin{array}{c}\text{Conversion\hspace{0.17em}of\hspace{0.17em}normoxic}\\ \text{to\hspace{0.17em}hypoxic}\end{array}}{\overbrace{\beta c(1-V)}}-\stackrel{\begin{array}{c}\text{Conversion\hspace{0.17em}of\hspace{0.17em}hypoxic}\\ \text{to\hspace{0.17em}necrotic}\end{array}}{\overbrace{({\alpha }_{h}h(1-T)+{\alpha }_{n}nh}}\\ \frac{\partial n}{\partial t}=\stackrel{\text{Conversion\hspace{0.17em}of\hspace{0.17em}hypoxic,\hspace{0.17em}normoxic,\hspace{0.17em}and\hspace{0.17em}vasculature\hspace{0.17em}to\hspace{0.17em}necrotic}}{\overbrace{{\alpha }_{h}h(1-V)+{\alpha }_{n}n(c+h+\nu )}}\\ \frac{\partial \nu }{\partial t}=\stackrel{\begin{array}{c}\text{Dispersal\hspace{0.17em}of}\\ \text{vasculature}\end{array}}{\overbrace{\nabla \cdot (D_{\nu }(x)(1-T)\nabla \nu )}}+\stackrel{\begin{array}{c}\text{Net\hspace{0.17em}proliferation}\\ \text{of\hspace{0.17em}vasculature}\end{array}}{\overbrace{\mu \frac{a}{K_m+a}\nu (1-T)}}\\ -\stackrel{\begin{array}{c}\text{Conversion\hspace{0.17em}of\hspace{0.17em}vasculature}\\ \text{to\hspace{0.17em}necrotic}\end{array}}{\overbrace{{\alpha }_{n}n\nu }}\\ \frac{\partial a}{\partial t}=\stackrel{\begin{array}{c}\text{Net\hspace{0.17em}dispersal}\\ \text{of\hspace{0.17em}angiogenic\hspace{0.17em}factors}\end{array}}{\overbrace{\nabla \cdot (D_a\nabla a)}}+\stackrel{\begin{array}{c}\text{Net\hspace{0.17em}production\hspace{0.17em}of}\\ \text{angiogenic\hspace{0.17em}factors}\end{array}}{\overbrace{{\delta }_{c}c+{\delta }_{h}h}}\\ -\stackrel{\begin{array}{c}\text{Net\hspace{0.17em}consumption\hspace{0.17em}of}\\ \text{angiogenic\hspace{0.17em}factors}\end{array}}{\overbrace{q\mu \frac{a}{K_m+a}\nu (1-T)-\omega a\nu }}-\stackrel{\begin{array}{c}\text{Decay\hspace{0.17em}of}\\ \text{angiogenic\hspace{0.17em}factors}\end{array}}{\overbrace{\lambda a}}.\end{array}\}$ (1)

In these equations, D(x) is the net rate of invasion (mm²/year). Glioma cells migrate faster along myelinated axons in the white matter than in the dense and less structured cortical gray matter. For this reason, we consider the net rate of invasion as piecewise constant, with non-zero values in the gray and white matter, D_g and D_w, respectively, with D_w > D_g, and zero in the regions of cerebral spinal fluid. Additionally, ρ (1/year) is the net proliferation rate of the normoxic cells, γ (1/year) and β (1/year) are the maximum conversion rates between the hypoxic and normoxic cell populations, α_n (1/year) is the rate at which cells undergo necrosis when in contact with necrotic cells (contact necrosis), α_h (1/year) is the rate of conversion of hypoxic cells to necrotic cells when nutrient levels fall too low, D_v (mm²/year) is the rate of dispersal of vasculature cells, estimated from Sherratt and Murray (1990), Levine et al. (2001) μ (1/year) is the vasculature proliferation rate, estimated from Xiu et al. (2006), T = (c + h + n + v)/K (dimensionless), where K is the carrying capacity (cells/mm³), and V = v/(v + c + h) (dimensionless) and is a surrogate for the local vasculature efficiency. Angiogenic factors are produced by both normoxic and hypoxic cells with rates δ_c (1/year) and δ_h (1/year) respectively, with δ_h > δ_c and are consumed by the vasculature for both regular vasculature maintenance [with rate ω (1/year)] and for vasculature proliferation [with rate q (1/year)]. Finally, the angiogenic factors are assumed to decay over time with rate λ (1/year) and disperse with rate D_v (mm²/year). Values for parameters related to the angiogenic factors were derived in part from work done in Levine et al. (2001), Serini et al. (2003), Mac Gabhann and Popel (2004). The reader is referred to (Swanson et al., 2011) for further details.

It is known that GBM tumors are extremely genetically heterogeneous both within a single tumor and between different tumors (Dunn et al., 2012). A large effort has been put forth to identify subtypes of GBMs by their dominating genotype (Verhaak et al., 2010). While the model parameters do not directly try to capture effects of single mutations, it is our belief that the dominating genotypes characterizing subtypes of GBMs ultimately result in different net rates of proliferation and invasion which would be used in our model. For example, pro-neural tumors are more likely to have the IDH-1 mutation and be secondary GBMs. In our model, this would manifest as a low-D, low-ρ tumor which begins as low grade and progresses into higher grade/malignancy.

In this model, there are different cell populations, normoxic, hypoxic, and necrotic competing for space and each with differing phenotypes: normoxic cells proliferating and invading, hypoxic cells only invading, and necrotic cells which are dead and just taking up space. Since each cell population is evolving in space and time, there is an effective spatial heterogeneity of predicted proliferative activity across space and time which could be analogized to heterogeneous Ki67 labeling across glioma specimens. Thus, while this model is attempting to capture the overarching phenotype of different tumors and assumes global constants for individual tumors, spatial heterogeneity in behavior is possible due to regional levels of vasculature which may result in hypoxia and/or necrosis. We also remark the aim of this model is not to predict cell-level behaviors (ex. Ki67), rather, use information obtained from routine imaging to quantify and explain imaging scale behavior and evolution.

The PIHNA model only captures cellular species and angiogenic factors. While these all have an impact on what is ultimately seen on MR imaging, in and of themselves they are not sufficient to describe enhancing regions of T2 and T1Gd MR images. Here we extend the PIHNA model to capture the imaging responses post anti-angiogenic treatment. To achieve this, we add to the model one additional element: edema, (l), to create the merged proliferation-invasion-hypoxia-necrosis-angiogenesis-edema (PIHNA-E) model. A schematic of the six species interactions is shown in Figure 3.

FIGURE 3

Figure 3. Schematic of the PIHNA-E model. The main components of the model are seen in the flow chart: c represents the normoxic glioma cells, h the hypoxic gliomas cells, v the vascular endothelial cells, n the necrotic cells, a, the angiogenic factors, and l the edematous fluid. Depending on the level of oxygen, normoxic, and hypoxic cells will undergo phenotypic switching. If oxygen levels fall too far and are not compensated for by sufficient angiogenesis, the hypoxic cells will undergo necrosis. Additionally, all cells will undergo necrosis if in contact with necrotic cells. Both hypoxic and normoxic cells release angiogenic factors into the extracellular space which recruit additional vasculature to increase the levels of oxygen. The angiogenic factors are removed from the system by interaction with vascular cells or natural decay. The local levels of angiogenic factors are indicative of the local degree of vessel permeability. Edematous liquid exits the vasculature where the permeability, K_trans(a), allows and enters the extracellular space where it diffuses and will be removed at rate d_l.

A Mathematical Model for Tumor-Induced Edema Formation

Generically, edema refers to a swelling phenomenon. While there are different types of cerebral edema, in the case of GBM, it is almost exclusively vasogenic edema which results from fluid and protein leakage from the breakdown of the blood brain barrier (Marmarou, 2007). Over the last few decades, there have been quite a few attempts to model vasogenic edema (Rapoport, 1978; Kumagai, 1986; Nagashima et al., 1990) and its associated interstitial pressure and interstitial fluid velocity (Baxter and Jain, 1989, 1990). These models were primarily based on Starling’s equation which describes fluid exchange between compartments due to pressure and osmotic gradients. These models are very detailed and are generally solved on shorter time scales, i.e., a few days versus months. In our efforts, while we do have interest in the specific mechanisms, we will take a coarser grained approach allowing us to approximate the phenomena over longer time scales relevant to tumor growth kinetics.

To begin, we make the simplifying assumption that the edema is only composed of fluid which has leaked into the extra-cellular space and has not yet been reabsorbed into the system. This fluid is assumed to leak into the extracellular space where the blood brain barrier has been compromised. From the PIHNA model, we can approximate these regions along with the degree of permeability from the local levels of present angiogenic factors. Once the fluid is in the extra-cellular space it moves via diffusion and is reabsorbed into the system at a constant rate. This process is written in the form of a partial differential equation as:

$\stackrel{\text{Change}\text{in}\text{time}\text{of}\text{fluid}}{\overbrace{\frac{\partial l}{\partial t}}}=\stackrel{\text{Fluid diffusion}}{\overbrace{\nabla \cdot \left(D_l\nabla l\right)}}+\stackrel{\text{Leakage}}{\overbrace{K_{\text{trans}}\left(a\right)\cdot \left(l_v-l\right)}}-\stackrel{\text{Drainage}}{\overbrace{{\delta }_{l}l}}.$ (2)

Here l is the concentration of edematous fluid, D_l (mm²/year) is the diffusion rate of the edematous liquid which would be analogous to an Apparent Diffusion Coefficient, ADC, value derived from diffusion-weighted MRI (Moritani, 2009), l_v (fluid/mm³) is the normal level of fluid in the vasculature, δ_l (1/year) is the reabsorption rate, and K_trans (1/year) is the transmission rate [analogous to the value K_trans measured on dynamic contrast enhanced MRI (DCE-MRI)], capturing the permeability surface area product per unit volume of tissue (Tofts, 1991) and is assumed to depend on the level of angiogenic factors, a, present. Homogeneous Neumann boundary conditions are assumed at the boundary of the brain to ensure no fluid leaves the brain.

The dependence of the K_trans coefficient on the angiogenic factors is assumed to take a Michaelis–Menten type I form:

$K_{\text{trans}}\left(a\right)=K_{\text{max}}\frac{a}{a+K_{\text{half}}}$ (3)

to reflect that VEGF (also known as vascular permeability factor, VPF) strongly influences vascular permeability (Bates, 2010). Here K_max (1/year) is the maximum possible value of K_trans, the value of which is calibrated to Grade IV gliomas (Patankar et al., 2005), and K_half (angiogenic factors/mm³) is the concentration of a at which K_trans reaches half of its maximum value. Thus, as the concentration of the angiogenic factors increases, the degree of vessel permeability will also increase until saturated.

Modeling Anti-Angiogenic Treatment

The PIHNA-E model describes the evolution of the tumor and its microenvironment in an untreated context. By understanding the premise of how specific therapies are meant to alter the system, one can also model the effects of various treatments. Here we are interested in anti-angiogenic treatment and, while there are many different types of drugs for this action, will focus on the drug bevacizumab.

Bevacizumab is a drug specifically targeted at the molecule vascular endothelial growth factor A (VEGF A). This particular angiogenic factor stimulates the growth of new vessels by binding with the vascular endothelial growth factor receptor (VEGFR2) on endothelial cells. Bevacizumab inhibits angiogenesis by binding to the free molecules of VEGF A and preventing them from binding to VEGFR2. An unintended consequence of this drug, however, is that beyond preventing the growth of new vessels, it also “normalizes” pre-existing vasculature (Jain, 2005; Verhoeff et al., 2009). That is, once the levels of stimulating angiogenic factors are reduced, the vessels are able to repair their leakiness and return to a normal size – making them more efficient nutrient deliverers. In our model, both of these phenomena can be captured by requiring higher levels of angiogenic factor to be present to have the same level of “action” in the contexts of both vessel proliferation and vessel permeability. Additionally, since the treatment is making the vessels more efficient, the level of vasculature needed for a cell to be normoxic will decrease, which we can capture by modifying the cell conversion rates from hypoxic to normoxic and from normoxic to hypoxic. Treatment is approximated by decreasing the parameter for conversion from normoxic to hypoxic (β) by a factor of 10, increasing the parameter for conversion from hypoxic and normoxic (γ) by a factor of 10, and increasing the required levels of angiogenic factors for inducing vascular growth and vessel permeability by 2 as supported by the studies in Desjardins et al. (2007), Zhang et al. (2009). The treatment modification of β and γ is representative of a dramatic increase in the efficiency of the blood vessels, though exact changes are not available from experimental data.

Simulations of Glioblastoma Growth and Response to Anti-Angiogenic Therapy

For simplicity, we consider in all simulations here a two-dimensional tumor growing on one axial slice of the brain, with the brain geometry defined from the BrainWeb atlas (Cocosco et al., 1997). The brain is primarily composed of three different types of matter, CSF, gray matter, and white matter. Glioblastomas originate in gray or white matter and due to physical barriers will not enter into the regions of CSF. New mass will often deform the barriers, a phenomena called mass effect, and while there are some models that attempt to capture this (Clatz et al., 2005; Mohamed and Davatzikos, 2005; Hogea et al., 2008), here the brain is considered a stationary domain.

In all simulations, the domain is taken to be a slice of human brain embedded in a grid [0, 147] mm × [0, 185] mm and the equations are spatially discretized on a grid with resolution of 1 mm × 1 mm using first order accurate finite volumes. Time integration is done with an operator splitter technique utilizing backward Euler for the diffusion terms and the TR-BDF2 algorithm (Leveque, 2005) for the reaction terms with a time step size of 1 day. The simulations were initiated with a small amount of normoxic cells distributed as

$c_0\left(x,y,t=0\right)=1000*exp\left(-\left(100\left[{\left(x-x_0\right)}^2+{\left(y-y_0\right)}^2\right]\right)\right)$

where (x₀, y₀) = (103, 83). Vasculature is set at 3% of the cell carrying capacity in all gray and white matter based on estimates from Blinkov and Glezer (1968), and all other quantities in the PIHNA-E model are initiated to zero. Unless otherwise stated, parameter values used in simulations for the edema equation are in Table 2, the additional parameter values are taken as specified in Swanson et al. (2011).

TABLE 2

Table 2. Parameter values for the edema equation used in all simulations unless stated otherwise in the text.

Results

Decoupling Imaging Changes from Tumor Response

To highlight the real impact of anti-angiogenic treatment as captured by our model, we chose parameter values that represent a patient with an aggressive GBM (net invasion rates D_w = 53 mm²/year and D_g = 0.53 mm²/year and net proliferation rate ρ = 75 1/year) and simulate tumor growth without treatment (Figure 4) and then compare to tumor growth with treatment (Figure 5). For comparison to what was done in previous work (Swanson et al., 2008a, 2011; Rockne et al., 2010), we refer to the region with total cell density summing to 80% of the cell carrying capacity (K) to correspond to what would enhance on the T1Gd, and started treatment when the T1Gd spherically symmetric equivalent radius was equal to 1 cm and was terminated 100 days later. This is consistent with a typical size of an abnormality seen clinically for consideration of anti-angiogenic treatment. Although treatment length can vary, 100 days is representative of the length of a typical course of treatment with anti-angiogenics in human GBM.

FIGURE 4

Figure 4. Illustrated here is the comparison of the simulated disease burden to what would be imageable on a FLAIR MRI in an untreated context for one set of growth parameters corresponding to an aggressive GBM. Below the plot shows the spherically symmetric equivalent radial growth of the regions containing tumor cells above a threshold and edematous fluid above a threshold. In the untreated context, these lines are nearly identical.

FIGURE 5

Figure 5. Comparison of the simulated disease burden to what would be imageable on a FLAIR MRI in a treated context for the same set of growth parameters as shown in Figure 4. Below the plot shows the spherically symmetric equivalent radial growth of the regions containing tumor cells above a threshold and edematous fluid above a threshold. Once treatment has begun, we see a drop in the levels of edema. Upon termination of the treatment, the edematous volume is seen to once again increase to the same size of the volume of tumorous cells.

Snapshots of the untreated case are shown in Figure 4 with the analogous snapshots corresponding to the same time points of the treated tumor being shown in Figure 5. The top row in both Figures 4 and 5 shows the density of the bulk tumor (the summed density of all the cell populations: normoxic, hypoxic, necrotic, and vasculature). The second row shows what the simulated FLAIR corresponding to the microenvironmental levels of edematous extracellular fluid. These figures also contain radial growth plots showing the equivalent spherical radii for the regions of interest corresponding to the tumor and the edema. The line representing the bulk tumor is calculated from the volume of tissue containing abnormal cells, normoxic, hypoxic, and necrotic, at levels greater than or equal to 16% of the carrying capacity, i.e., a density five times lower than what can be visualized on T1Gd (Swanson et al., 2008b). The edema radius was defined by considering the volume containing edematous fluid above 50% of the fluid level in the capillaries. There is no literature to guide the choice of the cutoffs for fluid volume constituting T2/FL visible edema. Thus, cutoffs were chosen to roughly match clinically observed behavior.

In the untreated case (Figure 4) the radial plot shows the size of the region impacted by edema evolves very similarly to the size of the region occupied by the bulk tumor throughout the entire course of growth. In the treated case (Figure 5), the edema grows at the same rate as the bulk tumor until the treatment begins at which point the edema begins to decline. Edema begins to increase again once the angiogenic factors have been able to accumulate at levels which overcome the impact of the anti-angiogenic drug. Once treatment is terminated the edema levels rise to again occupy a region of the same size as the bulk tumor.

These simulations support the hypothesis that anti-angiogenic treatment may not have a significant effect on the growth dynamics of the overall cell populations, while having a large impact on the level of edematous extracellular fluid and thus on what would be imageable on T2/FLAIR MRI. This is also in agreement with the current understanding that anti-angiogenic treatment serves as a surrogate to steroids for reducing swelling.

Exploring Response Across Tumor Kinetics

The virtual control experiment illustrated in Figures 4 and 5 is only providing insight into tumor/edema response for the case of one set of tumor growth kinetics. However, the range of radiographic response patterns seen clinically is broad, a few examples of which are shown in Figure 2. Previous work has shown patient specific values of net proliferation and invasion range over many orders of magnitude (Tracqui et al., 1995; Harpold et al., 2007). To investigate if the different types and extents of radiographic responses could be explained by different underlying tumor growth kinetics, we simulated tumor growth and the associated edema under treatment for many different combinations of net proliferation rates, ρ, and net invasion rates, D_w and D_g.

For all cases, treatment was started when the tumor reached a simulated 1 cm T1Gd radius for a total of 100 continuous days. Treatment was implemented in the same manner as the first case. Illustrative results are shown in Figure 6 from simulations with ρ = [5, 75, 125] (1/year) and D_w = [5.3, 53] (mm²/year) and D_g = D_w/10 consistent with the observed range hypothesized in human gliomas (Harpold et al., 2007) and since observed (Wang et al., 2009; Rockne et al., 2010).

FIGURE 6

Figure 6. Using the PIHNA-E model, holding all parameters constant except for D and ρ, one can observe many different responses to anti-angiogenic treatment in terms of the levels of edema. These responses vary from complete disappearance of imageable edema to stabilized lower levels to lowered levels of edema that continue to increase. When treatment is terminated, however, edema levels are always seen rise to once again reflect more closely the underlying disease burden. Treatment times are indicated by the gray boxes. These simulations suggest that the majority of imaging responses can be explained by considering how the drug impacts the tumor microenvironment alone without cytotoxic affects. Additionally, they represent a possible mechanism for identifying the patients who will receive significant benefit from the treatment.

From these six scenarios, summarized in Figure 6, we were able to observe a few trends. First, none of the simulations showed a significant change in the bulk tumor growth rate after treatment had begun. However, levels of edema were impacted and by different degrees in each of the simulations. In general, tumors with higher proliferative capacities (higher ρ’s), due to their higher metabolic needs, have larger regions of hypoxia and thus produce greater levels of angiogenic factors. Treatment of these tumors initially reduces the level of edema, however, the tumor quickly produces enough angiogenic factors to continue progressing under imaging. In contrast, the slower growing tumors (low ρ) produce low levels of angiogenic factors and treatment may result in a complete disappearance of the abnormality on imaging. Additionally, higher dispersion rates (high D_w and D_g) reduce the local metabolic needs and the production of angiogenic factors. Thus, the imaging of tumors with high D_w and high D_g improve for a longer time, however, ultimately the tumor does still produce enough angiogenic factors to be seen progressing on imaging. For all cases, once treatment is discontinued, the edema levels quickly rise to again reflect the underlying tumor burden.

What is particularly encouraging is that many of the different behaviors observed for the four patients illustrated in Figure 2 can be connected to different simulation predictions. For example, the first patient’s behavior is analogous to the moderate proliferation rate with a low invasion rate in that the imageable lesion initially decreased, stabilized, but after treatment dramatically increased. Additionally, the second patient can be compared to either the moderate or highly proliferative rate with a high invasion rate where the hyperintensity is seen to decrease at the beginning of treatment, but even while treatment is continuing, start growing again.

Discussion

Anti-angiogenics remain a controversial form of treatment for GBM due to the difficulty in assessing tumor response using MR imaging. The resulting reduction in swelling and related symptoms for a subset of patients keeps it an attractive option despite the lack of evidence of an increase in overall survival and the possibility of the treatment selecting for a more aggressive phenotype (Ebos et al., 2009; Pàez-Ribes et al., 2009; Verhoeff et al., 2009; Keunen et al., 2011). It is unclear how anti-angiogenic treatment impacts the cell phenotypes present and there is not yet a deep enough understanding or a unifying theory to provide explanation for all the different response patterns. Thus, a priori identification of patients who will receive a significant symptom-reduction benefit remains difficult. As there are other side effects and consequences from anti-angiogenic therapy, being able to make this early distinction would help remove the controversial nature of this therapy.

The model developed in this work, built on the PIHNA model for glioma proliferation and invasion (Swanson et al., 2011; Gu et al., 2012), is meant to illustrate a first step toward the creation of a tool for identifying patients who will receive the greatest benefit from anti-angiogenic treatment. It captures the formation of edema caused by leaky vasculature, and is thus able to decouple what would be seen on the T2/FLAIR MRI from the true underlying disease burden. In effect, this would help the clinicians to “turn the light back on” by being able to infer the disease burden that lies beyond what is captured by imaging alone. While this model is not meant to capture individual cell behavior, it does provide a map between overall tumor growth kinetics and treatment response on the imaging/continuum scale.

Many simplifying assumptions have been made in the creation of this model such as ignoring the possible direct impact of anti-angiogenic therapy on cell proliferation, the likely presence of thrombosis (Tehrani et al., 2008), and that the hyperintensity on the T2/FLAIR image is entirely a result of fluid leaking from the vasculature. The model is clearly incomplete, and future modifications of the model will need to consider these other phenomena as well as other possible factors on T2/FLAIR hyperintensity, such as higher cell density and additional cytoplasm. However, even in its current state, it has been able to exhibit many of the types of response patterns observed clinically. It is particularly encouraging because the modeling effects of treatment were held constant and only varying the net dispersal and net proliferation rates, D and ρ, respectively, was sufficient to produce a wide range of imaging responses analogous to what is seen clinically. That is, by modeling the treatment in the exact same way for different values of proliferation and diffusion in the tumor growth model, the visible levels of edema are seen to respond in different ways. In general, the simulations predict edema (swelling) to decrease, supporting the role of these drugs as surrogates for steroids for reduction of symptoms, analogous to the most current understanding (Deming, 2012). Though, the model clearly illustrates that not all patients would receive the same benefit.

Another interesting implication from these simulations is that while different imaging responses to treatment were achieved for the same treatment conditions, in all cases the bulk tumor is seen to progress with little deflection in overall tumor growth rates. This result could be considered evidence against the cytotoxic effects of anti-angiogenic treatment when administered exclusively, also in agreement with current clinical understanding (Verhoeff et al., 2009). While this work is not directly speaking to survival outcomes, we remark that these results highlight the potential for the mathematical model paradigm to serve in evaluating clinical trial outcomes by analyzing relative benefit from anti-angiogenics, especially in the case of low-N trials. As a particular example, such models may have the potential to be applied to patient cohorts for exploring how differential effects of anti-angiogenics on imaging may or may not relate to overall outcomes.

A drawback of this model is the large number of parameters required. In this document, we assumed the primary influential factor were the net rates of invasion and proliferation and thus held all other parameters constant. These other parameters are likely different patient to patient, however, as demonstrated here, changes in the small number of parameters are sufficient to produce a wide variety of imaging responses. While future sensitivity analysis is required, we believe the work is here is evidence that a complicated explanation for the different imaging responses to anti-angiogenic therapy may not be needed. Major next steps of this work will involve developing techniques for obtaining patient specific growth parameters from pre-treatment images which we believe will be successful from previous accomplishments with a simpler model capturing just the proliferation and invasion tumor characteristics PI (Swanson, 1999, 2002; Swanson et al., 2000, 2003, 2004, 2008a; Szeto et al., 2009; Wang et al., 2009; Rockne et al., 2010).

Also, modification of the model to capture the pressure induced from the vasogenic edema and possible herniation would allow for deeper understanding of the steroid-like reduction in swelling. For each new feature eventually added validation tests will be required, however, we believe the results presented here in and of themselves represent a significant step in overcoming clinical imaging restrictions with mathematical models.

Conflict of Interest Statement

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.

Acknowledgments

The authors would like to greatly acknowledge the funding of this work which was sponsored in part by NIH Grants: R01 CA 16437 and U54 CA 143970, along with the NSF Fellowship 1122322.

References

Batchelor, T. T., Sorensen, A. G., di Tomaso, E., Zhang, W. T., Duda, D. G., Cohen, K. S., et al. (2007). AZD2171, a pan-VEGF receptor tyrosine kinase inhibitor, normalizes tumor vasculature and alleviates edema in glioblastoma patients. Cancer Cell 11, 83–95.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Bates, D. O. (2010). Vascular endothelial growth factors and vascular permeability. Cardiovasc. Res. 87, 262–271.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Baxter, L. T., and Jain, R. K. (1989). Transport of fluid and macromolecules in tumors. I. Role of interstitial pressure and convection. Microvasc. Res. 37, 77–104.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Baxter, L. T., and Jain, R. K. (1990). Transport of fluid and macromolecules in tumors. II. Role of heterogeneous perfusion and lymphatics. Microvasc. Res. 40, 246–263.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Blinkov, S. M., and Glezer, I. I. (1968). The Human Brain in Figures and Tables: A Quantitative Handbook. New York: Basic Books.

Bohman, L. E., Swanson, K. R., Moore, J. L., Rockne, R., Mandigo, C., Hankinson, T., et al. (2010). Preoperative MRI characteristics of glioblastoma multiforme: implications for understanding glioma ontogeny. Neurosurgery 67, 1319–1327.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Clatz, O., Sermesant, M., Bondiau, P. Y., Delingette, H., Warfield, S. K., Malandain, G., et al. (2005). Realistic simulation of the 3-D growth of brain tumors in MR images coupling diffusion with biomechanical deformation. IEEE Trans. Med. Imaging 24, 1334–1346.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Cocosco, C., Kollokian, V., Kwan, R. S., and Evans, A. (1997). “Brainweb: online interface to a 3D MRI simulated brain database,” in Neuroimage, Proceedings of the Third International Conference on the Functional Mapping of the Human Brain, Copenhagen.

Deming, S. (2012). Bevacizumab for glioblastoma: symptom relief proven, overall survival benefit being studied. OncoLog 57. Available at: http://www2.mdanderson.org/depts/oncolog/articles/pf/12/6-jun/6-12-2-pf.html

Desjardins, A., Barboriak, D. P., Herndon, J. E. II, Reardon, D. A., Quinn, J. A., Rich, J. N., et al. (2007). Dynamic contrast-enhanced magnetic resonance imaging (DCE-MRI) evaluation in glioblastoma (GBM) patients treated with bevacizumab (BEV) and irinotecan (CPT-11). J. Clin. Oncol. 25(Suppl.), abstr. 2029.

Dunn, G. P., Rinne, M. L., Wykosky, J., Genovese, G., Quayle, S. N., Dunn, I. F., et al. (2012). Emerging insights into the molecular and cellular basis of glioblastoma. Genes Dev. 26, 756–784.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Ebos, J. M., Lee, C. R., Cruz-Munoz, W., Bjarnason, G. A., Christensen, J. G., and Kerbel, R. S. (2009). Accelerated metastasis after short-term treatment with a potent inhibitor of tumor angiogenesis. Cancer Cell 15, 232–239.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Folkman, J. (1971). Tumor angiogenesis: therapeutic implications. N. Engl. J. Med. 285, 1182–1186.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Gu, S., Chakraborty, G., Champley, K., Alessio, A. M., Claridge, J., Rockne, R., et al. (2012). Applying a patient-specific bio-mathematical model of glioma growth to develop virtual [18F]-FMISO-PET images. Math. Med. Biol. 29, 31–48.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Harpold, H. L., Alvord, E. C. Jr., and Swanson, K. R. (2007). The evolution of mathematical modeling of glioma proliferation and invasion. J. Neuropathol. Exp. Neurol. 66, 1–9.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Hogea, C., Davatzikos, C., and Biros, G. (2008). An image-driven parameter estimation problem for a reaction-diffusion glioma growth model with mass effects. J. Math. Biol. 56, 793–825.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Holash, J., Maisonpierre, P. C., Compton, D., Boland, P., Alexander, C. R., Zagzag, D., et al. (1999). Vessel cooption, regression, and growth in tumors mediated by angiopoietins and VEGF. Science 284, 1994–1998.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Hurwitz, H., Fehrenbacher, L., Novotny, W., Cartwright, T., Hainsworth, J., Heim, W., et al. (2004). Bevacizumab plus irinotecan, fluorouracil, and leucovorin for metastatic colorectal cancer. N. Engl. J. Med. 350, 2335–2342.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Iwamoto, F. M., Abrey, L. E., Beal, K., Gutin, P. H., Rosenblum, M. K., Reuter, V. E., et al. (2009). Patterns of relapse and prognosis after bevacizumab failure in recurrent glioblastoma. Neurology 73, 1200–1206.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Jain, R. K. (2005). Normalization of tumor vasculature: an emerging concept in antiangiogenic therapy. Science 307, 58–62.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Kerbel, R. S. (2000). Tumor angiogenesis: past, present and the near future. Carcinogenesis 21, 505–515.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Keunen, O., Johansson, M., Oudin, A., Sanzey, M., Rahim, S. A., Fack, F., et al. (2011). Anti-VEGF treatment reduces blood supply and increases tumor cell invasion in glioblastoma. Proc. Natl. Acad. Sci. U.S.A. 108, 3749–3754.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Kumagai, Y. (1986). Mathematical model of brain edema and optimal control of intracranial pressure. No To Shinkei 38, 885–891.

Pubmed Abstract | Pubmed Full Text

Leenders, W. P., Küsters, B., and de Waal, R. M. (2002). Vessel co-option: how tumors obtain blood supply in the absence of sprouting angiogenesis. Endothelium 9, 83–87.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

LeVeque, R. (2005). Numerical Methods for Conservation Laws, 2nd Edn. Basel: Birkhäuser.

Levine, H. A., Pamuk, S., Sleeman, B. D., and Nilsen-Hamilton, M. (2001). Mathematical modeling of capillary formation and development in tumor angiogenesis: penetration into the stroma. Bull. Math. Biol. 63, 801–863.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Louis, D. N., Ohgaki, H., Wiestler, O. D., Cavenee, W. K., Burger, P. C., Jouvet, A., et al. (2007). The 2007 WHO classification of tumours of the central nervous system. Acta Neuropathol. 114, 97–109.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Mac Gabhann, F., and Popel, A. S. (2004). Model of competitive binding of vascular endothelial growth factor and placental growth factor to VEGF receptors on endothelial cells. Am. J. Physiol. Heart Circ. Physiol. 286, H153–H164.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Marmarou, A. (2007). A review of progress in understanding the pathophysiology and treatment of brain edema. Neurosurg. Focus 22, E1.

CrossRef Full Text

Mohamed, A., and Davatzikos, C. (2005). Finite element modeling of brain tumor mass-effect from 3D medical images. Med. Image Comput. Comput. Assist. Interv. 2005(Pt 13749), 400–408.

Moritani, T. (2009). Diffusion-Weighted MR Imaging of the Brain. New York: Springer.

Nagashima, T., Horwitz, B., and Rapoport, S. I. (1990). A mathematical model for vasogenic brain edema. Adv. Neurol. 52, 317–326.

Pubmed Abstract | Pubmed Full Text

Norden, A. D., Young, G. S., Setayesh, K., Muzikansky, A., Klufas, R., Ross, G. L., et al. (2008). Bevacizumab for recurrent malignant gliomas: efficacy, toxicity, and patterns of recurrence. Neurology 70, 779–787.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Pàez-Ribes, M., Allen, E., Hudock, J., Takeda, T., Okuyama, H., Viñals, F., et al. (2009). Antiangiogenic therapy elicits malignant progression of tumors to increased local invasion and distant metastasis. Cancer Cell 15, 220–231.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Patankar, T. F., Haroon, H. A., Mills, S. J., Balériaux, D., Buckley, D. L., Parker, G. J., et al. (2005). Is volume transfer coefficient (K(trans)) related to histologic grade in human gliomas? AJNR Am. J. Neuroradiol. 26, 2455–2465.

Pubmed Abstract | Pubmed Full Text

Rapoport, S. I. (1978). A mathematical model for vasogenic brain edema. J. Theor. Biol. 74, 439–467.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Rockne, R., Rockhill, J. K., Mrugala, M., Spence, A. M., Kalet, I., Hendrickson, K., et al. (2010). Predicting the efficacy of radiotherapy in individual glioblastoma patients in vivo: a mathematical modeling approach. Phys. Med. Biol. 55, 3271–3285.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Sandler, A., Gray, R., Perry, M. C., Brahmer, J., Schiller, J. H., Dowlati, A., et al. (2006). Paclitaxel-carboplatin alone or with bevacizumab for non-small-cell lung cancer. N. Engl. J. Med. 355, 2542–2550.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Serini, G., Ambrosi, D., Giraudo, E., Gamba, A., Preziosi, L., and Bussolino, F. (2003). Modeling the early stages of vascular network assembly. EMBO J. 22, 1771–1779.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Sherratt, J. A., and Murray, J. D. (1990). Models of epidermal wound-healing. Proc. Biol. Sci. 241, 29–36.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Silbergeld, D. L., and Chicoine, M. R. (1997). Isolation and characterization of human malignant glioma cells from histologically normal brain. J. Neurosurg. 86, 525–531.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Stupp, R., Mason, W. P., van den Bent, M. J., Weller, M., Fisher, B., Taphoorn, M. J., et al. (2005). Radiotherapy plus concomitant and adjuvant temozolomide for glioblastoma. N. Engl. J. Med. 352, 987–991.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Swanson, K. R. (1999). Mathematical Modeling of the Growth and Control of Tumors. Applied Mathematics. Ph.D. thesis, University of Washington, Seattle.

Swanson, K. R. (2002). Serial imaging observations and postmortem examination of an untreated glioblastoma: a traveling wave of glioma growth and invasion. Neuro-oncology 4, 340.

Swanson, K. R., Alvord, E. C. Jr., and Murray, J. D. (2000). A quantitative model for differential motility of gliomas in grey and white matter. Cell Prolif. 33, 317–329.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Swanson, K. R., Alvord, E. C. Jr., and Murray, J. D. (2002a). Quantifying efficacy of chemotherapy of brain tumors with homogeneous and heterogeneous drug delivery. Acta Biotheor. 50, 223–237.

CrossRef Full Text

Swanson, K. R., Alvord, E. C. Jr., and Murray, J. D. (2002b). Virtual brain tumours (gliomas) enhance the reality of medical imaging and highlight inadequacies of current therapy. Br. J. Cancer 86, 14–18.

CrossRef Full Text

Swanson, K. R., Alvord, E. C., Murray J. D. (2004). Dynamics of a model for brain tumors reveals a small window for therapeutic intervention. Discrete Continuous Dyn. Syst. Ser. B 4, 289–295.

Swanson, K. R., Bridge, C., Murray, J. D., and Alvord, E. C. Jr. (2003). Virtual and real brain tumors: using mathematical modeling’ to quantify glioma growth and invasion. J. Neurol. Sci. 216, 1–10.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Swanson, K. R., Harpold, H. L., Peacock, D. L., Rockne, R., Pennington, C., Kilbride, L., et al. (2008a). Velocity of radial expansion of contrast-enhancing gliomas and the effectiveness of radiotherapy in individual patients: a proof of principle. Clin. Oncol. 20, 301–308.

CrossRef Full Text

Swanson, K. R., Rostomily, R. C., and Alvord, E. C. Jr. (2008b). A mathematical modelling tool for predicting survival of individual patients following resection of glioblastoma: a proof of principle. Br. J. Cancer 98, 113–119.

CrossRef Full Text

Swanson, K. R., Rockne, R. C., Claridge, J., Chaplain, M. A., Alvord, E. C. Jr., and Anderson, A. R. (2011). Quantifying the role of angiogenesis in malignant progression of gliomas: in silico modeling integrates imaging and histology. Cancer Res. 71, 7366–7375.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Szeto, M. D., Chakraborty, G., Hadley, J., Rockne, R., Muzi, M., Alvord, E. C. Jr., et al. (2009). Quantitative metrics of net proliferation and invasion link biological aggressiveness assessed by MRI with hypoxia assessed by FMISO-PET in newly diagnosed glioblastomas. Cancer Res. 69, 4502–4509.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Tehrani, M., Friedman, T. M., Olson, J. J., and Brat, D. J. (2008). Intravascular thrombosis in central nervous system malignancies: a potential role in astrocytoma progression to glioblastoma. Brain Pathol. 18, 164–171.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Tofts, P. S. (1991). Measurement of the blood-brain barrier permeability and leakage space using dynamic MR imaging. 1. Fundamental concepts. Magn. Reson. Med. 17, 357–367.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Tracqui, P., Cruywagen, G. C., Woodward, D. E., Bartoo, G. T., Murray, J. D., and Alvord, E. C. Jr. (1995). A mathematical model of glioma growth: the effect of chemotherapy on spatio-temporal growth. Cell Prolif. 28, 17–31.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Verhaak, R. G., Hoadley, K. A., Purdom, E., Wang, V., Qi, Y., Wilkerson, M. D., et al. (2010). Integrated genomic analysis identifies clinically relevant subtypes of glioblastoma characterized by abnormalities in PDGFRA, IDH1, EGFR, and NF1. Cancer Cell 17, 98–110.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Verhoeff, J. J., van Tellingen, O., Claes, A., Stalpers, L. J., van Linde, M. E., Richel, D. J., et al. (2009). Concerns about anti-angiogenic treatment in patients with glioblastoma multiforme. BMC Cancer 9:444. doi:10.1186/1471-2407-9-444

CrossRef Full Text

Wang, C. H., Rockhill, J. K., Mrugala, M., Peacock, D. L., Lai, A., Jusenius, K., et al. (2009). Prognostic significance of growth kinetics in newly diagnosed glioblastomas revealed by combining serial imaging with a novel biomathematical model. Cancer Res. 69, 9133–9140.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Xiu, M., Turner, S. M., Busch, R., Gee, T. A., and Hellerstein, M. K. (2006). Measurement of endothelial cell proliferation rate in vivo using 2H20 labeling: a kinetics biomarker of angiogenesis. FASEB J. 20, A718–A718.

Zhang, W., Kreisl, T. N., Solomon, J., Reynolds, R. C., Glen, D. R., Cox, R. W., et al. (2009). Acute effects of bevacizumab on glioblastoma vascularity assessed with DCE-MRI and relation to patient survival. Proc. Intl. Soc. Mag. Reson. Med. 17, 282.