May 20, 2018, Sunday
University of Colorado at Boulder Search A to Z Campus Map University of Colorado at Boulder CU 
Search Links

A Reaction-Diffusion Model of Cancer Invasion

From MathBio

Jump to: navigation, search

Article review by Nicole Look, March 2012.

Article: Gatenby, R. A. & Gawlinski, E. T. (1996). A Reaction-Diffusion Model of Cancer Invasion. Cancer Res, 56, 5745-5753.

Executive Summary

The authors sought to understand the underlying mechanisms that govern invasion and destruction of normal tissues by metastatic cancer. They explore the effectiveness of a modified reaction-diffusion model to predict a pH gradient at the tumor-host interface by comparing mathematical analysis to experimental data and clinical observations. The authors initially hypothesize that as neoplastic tissue reverts to more primitive forms of the glycolytic metabolic pathways, acid production increases and diffuses to surrounding healthy tissue. The acidic environment is ideal for the proliferation of tumor cells, but deadly for normal cells. They justify the model’s critical biological parameters from experimental and observational data. Furthermore, they determine tumor wave front velocities that directly correlate to a predicted, and previously unidentified, hypocellular interstitial gap at the tumor-host interface. Most importantly, the model makes no assumptions about the processes of carcinogenesis or central necrosis. It focuses on the interactions between tumor cells and normal cells at the tumor-host interface, because these interactions significantly influence the progression of invasive cancer.

Context/Biological Phenomena

The mechanism for tumor invasion is thought to begin with the increased use of glycolytic metabolism by evolving tumor populations. The neoplastic cells consume almost one order of magnitude more glucose than normal tissues, effectively increasing the production of H+ ions and lowering the interstitial pH by about 0.5 units. While tumor cells thrive in pH levels nearing 6.8, the optimal pH level for normal cells is 7.4. As the environment acidifies below pH levels of 7.1, at the tumor-host interface, normal cells lose viability and die. The progressive loss of normal cells at the tumor-host interface accelerates tumor invasion.


Gatenby [1] was among the first to apply competition theory to the propagation of tumor cells by relating tumor growth to population invasions found in nature, even though competition theory does not pertain to most normal cellular populations. But this enabled the application of well-developed ecological mathematical models to the interaction of tumor and normal cells. By assuming competition theory as a descriptor of tumor invasion, Gatenby presented a competition model that successfully predicted established metabolic changes of tumor cells: increased glycolysis, glucose utilization, and lactic acid production and decreased extracellular pH. The proposed model in this review derives from the transformation-induced metabolic changes of tumor cells, specifically at the tumor-host interface. It has been suggested that such metabolic changes allow tumor cells to overcome the disadvantage of entering a hostile environment by creating a non-viable environment for normal cells.

A. Description of the Model

The model is a system of three coupled reaction-diffusion equations to describe the tumor-host interface based on acid-mediation hypothesis of tumor invasion. Furthermore, it attempts to predict when benign tumors become aggressively invasive tumors in relation to parameters passing through a critical value. The model evaluates the spatial distribution and temporal evolution of three variables:

N_{1}(x,t): the density of normal tissue (cells/(cm^3 )), N_{2}(x,t): the density of neoplastic tissue (cells/(cm^3 )), L(x,t): the excess concentration of H^{{+}} ions (M).

Behavior of normal tissue growth:

Failed to parse (lexing error): \frac{\partial N_{1}}{\partial t}=r_{1} N_{1} (1-\frac{N_{1}}{K_{1}} -\alpha_{12} \frac{N_{2}}{K_{2}} )-d_{1} LN_{1}+\nabla ∙(D_{N1} [N_{2} ]\nabla N_{1} )


where there is (a) logistic growth of N_{{1}} with a growth rate r_{{1}} and carrying capacity K_{{1}}, (b) population competition with a tumor tissue characterized by a Lotka-Volterra competition strength parameter \alpha _{{12}}, (c) interaction of N_{{1}} with excess proton ions that lead to a death rate proportional to L, and (d) cellular diffusion with an N_{{2}}-dependent diffusion coefficient D_{{N1}}[N_{{2}}].

Behavior of neoplastic tissue growth:

Failed to parse (lexing error): \frac{\partial N_{2}}{\partial t}=r_{2} N_{2} (1-\frac{N_{2}}{K_{2}} -\alpha_{21} \frac{N_{1}}{K_{1}} )+\nabla∙(D_{N2} [N_{1} ]\nabla N_{2} )


wherer_{{1}},K_{{2}} are the growth rate and carrying capacity of the neoplastic tissue and \alpha _{{21}} is the competition parameter describing the decrease of neoplastic tissue growth from lack of space and resources. Equations (1) and (2) differ by the lack of a death term in equation (2) due to the excess acid.

Excess H^{{+}} ion concentration is governed by:

{\frac  {\partial L}{\partial t}}=r_{{3}}N_{{2}}-d_{{3}}L+D_{{3}}\nabla ^{{2}}L (3)

where r_{{3}} is the production rate, d_{{3}} is the reabsorption rate and D_{{3}} is the H^{{+}} ion diffusion constant.

B. Assumptions

(1) D_{{N1}}[N_{{2}}]=0→ assumes healthy tissue is well regulated and will not diffuse in space

(2) D_{{N2}}[N_{{1}}]=D_{{2}}(1-{\frac  {N_{{1}}}{K_{{1}}}})→assumes D_{2} is the diffusion constant of neoplastic tissue in the absence of normal tissue. When the healthy tissue reaches its carrying capacity in concentration, the tumor is confined such that the diffusion coefficient is zero. This derives from the assumption that fully healthy tissue (N_{1}=K_{1}) can confine tumors using immunological and non-immunological defense mechanisms. Therefore, only when normal tissue density significantly reduces can tumor tissue spread.

(3) Equation (3) assumes that excess H^{{+}} ions are generated proportional to neoplastic cell density and diffuse chemically. Upon increased pH, there are certain mechanisms such as buffering and large-scale vascular evacuation that are considered by the reabsorption rate constant.

(4) \alpha _{{12}}=\alpha _{{21}}=0→ assumes aggressively invasive cancer such that tumor and normal tissues do not coexist. This implies that the Lotka-Volterra competition does not affect the structure and dynamics of the tumor-host interface.

C. Transformation of Variables

\eta _{{1}}={\frac  {N_{{1}}}{K_{{1}}}},\eta _{{2}}={\frac  {N_{{2}}}{K_{{2}}}},\Lambda ={\frac  {L}{L_{{0}}}}={\frac  {L}{((r_{{3}}K_{{2}})/d_{{3}})}},\tau =r_{{1}}t,\xi ={\sqrt  {{\frac  {r_{{1}}}{D_{{3}}}}}}x (4)

Dimensionless Form of the Equations:

{\frac  {\partial \eta _{{1}}}{\partial \tau }}=\eta _{{1}}(1-\eta _{{1}})-\delta _{{1}}\Lambda \eta _{{1}} (5)

Failed to parse (lexing error): \frac{\partial \eta_{2}}{\partial \tau}=\rho_{2} \eta_{2} (1-\eta_{2} )-\nabla_{\xi}∙[\Delta_{2} (1-\eta_{1} ) \nabla_{\xi} \eta_{2}]


{\frac  {\partial \Lambda }{\partial \tau }}=\delta _{{3}}(\eta _{{2}}-\Lambda )+\nabla _{{\xi }}^{{2}}\Lambda (7)

with dimensionless parameters:

Failed to parse (lexing error): \delta_{1}=\frac{d_{1}}{d_{3}}×\frac{r_{3}}{r_{1}}×K_{2}, \rho_{2}=\frac{r_{2}}{r_{1}} , \Delta_{2}=\frac{D_{2}}{D_{3}} , \delta_{3}=\frac{d_{3}}{r_{1}}


Linear Stability Analysis

Unconditionally unstable fixed points occur at

FP1:\eta _{{1}}^{{*}}=0,\eta _{{2}}^{{*}}=0,\Lambda ^{{*}}=0

FP2:\eta _{{1}}^{{*}}=1,\eta _{{2}}^{{*}}=0,\Lambda ^{{*}}=0.

Fixed point (1) corresponds to absence of all tissues and acids. Fixed point (2) relates to the absence of tumor tissue and acids and the maximum concentration of normal tissue at its carrying capacity. Perturbations from the two unstable fixed points will grow toward the two stable fixed points at

FP3:\eta _{{1}}^{{*}}=1-\delta _{{1}},\eta _{{2}}^{{*}}=1,\Lambda ^{{*}}=1

FP4:\eta _{{1}}^{{*}}=0,\eta _{{2}}^{{*}}=1,\Lambda ^{{*}}=1.

Fixed point (3) corresponds to the coexistence of normal tissue at a decreased level and tumor tissue at its carrying capacity. Fixed point (4) relates to the absence of normal tissue and tumor tissue at its carrying capacity. Interestingly, fixed point (3) is unstable and fixed point (4) is stable when \delta _{{1}}>1. On the other hand, fixed point (3) is stable and fixed point (4) is unstable when \delta _{{1}}<1.

A. Implications of the Parameter \delta _{{1}}

As \delta _{{1}} increases through unity, the tumor’s structure and dynamics changes from FP 3 to FP 4. The crossover occurs when \delta _{{1}}=1, such that \delta _{{1}}<1 corresponds to non-invasive tumors and \delta _{{1}}>1 corresponds to invasive malignant tumors. Given by the transformation of variables (8), changes to d_{{3}},r_{{3}} and K_{{2}} have been clinical shown to increase \delta_{1} through unity.

1. \delta _{{1}} increases linearly with the carrying capacity of the tumor population K_{{2}}. When there is an inadequate supply of nutrients, the carrying capacity decreases and limits tumor growth. In contrast, when there is greater substrate availability, the carrying capacity increases with tumor growth.

2. \delta _{{1}} increases with the production rate of acid by the tumor r_{{3}}. Studies have shown that malignant tumors are more glycolytic than benign tumors, such that increased production of acid correlates to more aggressive and invasive tumor growth.

3. \delta _{{1}} increases as the reabsorption rate of acid, d_{{3}}, decreases. Although not explicitly demonstrated, the relationship areas of scarring and inflammation seem to be more susceptible to tumor growth.

Solution to the Full Mathematical Model

The author used data, consisting of pH measurements at different locations within and near tumors, from Martin and Jain [2] to determine parameters d_{{3}}, r_{{3}} and r_{{0}}. In addition, the data supported their hypothesis, such that a pH gradient extends from the tumor outward and excessH^{{+}} ion concentration diminishes at distances farther from the tumor. The data displayed in Fig. 1 was further compared to the presented model by evaluating equation (3) of excess H^{{+}} ion concentration with the assumptions: the pH levels change slowly, the Laplacian is two-dimensional, the tumor is at carrying capacity, and the tumor-host interface is smooth and sharp. The assumptions lead to

r_{{3}}K_{{2}}-d_{{3}}L+{\frac  {D_{{3}}}{r}}{\frac  {\partial }{\partial r}}(r{\frac  {\partial L}{\partial r}})=0 r<r_{{0}}

-d_{{3}}L+{\frac  {D_{{3}}}{r}}{\frac  {\partial }{\partial r}}(r{\frac  {\partial L}{\partial r}})=0 r\geq r_{{0}} (9)

where r_{{0}} is the tumor radius. Equations (9) are solved as

L(r)=L_{{0}}[1-qr_{{0}}K_{{1}}(qr_{{0}})I_{{0}}(qr)] r<r_{{0}}

L(r)=qr_{{0}}L_{{0}}I_{{1}}(qr_{{0}})K_{{0}}(qr) r\geq r_{{0}} (10)

where I_{{0,1}},K_{{0,1}} are modified Bessel functions, and q={\sqrt  {{\frac  {d_{{3}}}{D_{{3}}}}}},L_{{0}}={\frac  {r_{{3}}K_{{2}}}{d_{{3}}}}. The solutions (10) were fitted to the data in Fig. 1.

Figure 1: Solutions fitted to observed data from Martin and Jain

Assuming a sharp tumor edge and small wave front propagation velocity, the authors approximated the interfacial profile structures in relation to the wave front propagation velocity. Both assumptions are biologically relevant. By assuming the solutions have the form of f(\xi ,t)=f(\xi -c\tau ), where each field has a wave front profile f that propagates in the +\xi direction with velocity c, they can be substituted into the dimensionless equations (5 – 7), to yield a new system of ordinary differential equations, with respect to the co-moving coordinate \zeta =\xi -c\tau :

-c\eta _{{1}}^{{'}}=\eta _{{1}}(1-\eta _{{1}})-\delta _{{1}}\Lambda \eta _{{1}}

-c\eta _{{2}}^{{'}}=\rho _{{2}}\eta _{{2}}(1-\eta _{{2}})+\Delta _{{2}}[(1-\eta _{{1}})\eta ''_{{2}}-\eta '_{{1}}\eta '_{{2}}]

-c\Lambda '=\delta _{{3}}(\eta _{{2}}-\Lambda )+\Lambda '' (11)

The solutions to system (11) are as follows:


where \rho ={\frac  {\delta _{{1}}}{2c{\sqrt  {\delta _{{3}}}}}},q={\sqrt  {\delta _{{3}}}},r={\frac  {1}{c}} and s={\frac  {\delta _{1}-1}{c}}. Furthermore, \gamma (a,x)=\int {e^{{-t}}\rho ^{{-1}}dt} is an incomplete \gamma function.

Finally, the edge position E(\tau ) and the width W(\tau ) of the wave front profiles are given by


The interfacial widths of acid and tumor tissue profiles were established by using equations (12) and (13) to solve the integrals in equations (15) and (16), giving:

W_{{\Lambda }}={\sqrt  {{\frac  {2}{\delta _{{3}}}}}} (17)

W_{{\eta _{{2}}}}={\frac  {\pi c}{{\sqrt  {3}}\rho _{{2}}}} (18)

where c is the speed of the propagating wave front. Because the derivation for the width of normal tissue profile was too complex to yield a simple result, the analytical approximations were found by numerically integrating equations (15) and (16) with f(\xi ,\tau ) replaced by approximate solution for \eta _{{1}}.


The model was simulated with the parameters given in Table 1 for noninvasive tumors (FP 3: \delta _{{1}}=0.5<1) and for invasive tumors (FP 4: \delta _{{1}}=12.5>1). The numerical solutions are displayed in Figure 2, in which the wave fronts propagate from left to right.

Figure 2: numerical solutions with propagating wave fronts from left to right. The top figure displays diminishing normal tissue and the bottom figure displays coexistence of normal tissue and tumor tissue

Figure 2a shows the normal tissue diminishing in the presence of the advancing wave front, just as expected in the linear stability analysis and implications of \delta _{{1}}. As \delta _{{1}} increases beyond unity, a noticeable hypocellular interstitial gap (\xi _{{0}}) occurs between the advancing tumor edge and the retreating healthy tissue. When \delta _{{1}}>>1, the gap is estimated as

\xi _{{0}}={\frac  {log[{\frac  {\delta _{{1}}}{2}}]}{{\sqrt  {\delta _{{3}}}}}}+{\sqrt  {{\frac  {c}{{\sqrt  {\delta _{{3}}}}}}}}\delta _{{1}}>>1 (19).

Based on equation (19), as\delta _{{1}} increases, the interstitial gap becomes large. In addition, when \eta _{{1}}=0, then \delta _{{2}}\neq 0 and vice versa, as demonstrated in equation (6). Figure 3 shows the healthy tissue at the tumor-host interface in relation to\delta _{{1}}, such that the hypocellular interstitial gap continues to increase as \delta _{{1}} gets large while the healthy tissue diminishes.

Figure 3: Normal tissue at tumor-host interface in relation to \delta _{{1}}

Figure 2b demonstrates the coexistence of tumor tissue and normal tissue behind the wave front when \delta _{{1}}<1, where the concentration of healthy tissue is \eta _{{1}}=1-\delta _{{1}}. Furthermore, in this case, the healthy tissue profile overlaps the tumor edge. No such hypocellular interstitial gap exists between the tumor and healthy tissue.

In both cases, the H^{{+}} ion concentration is located on the sharp tumor edge, noticeable at \xi =0 in Figures 2.

A. Model Predictions Correspond to Clinical Observations

Martin and Jain (1994) observed a similar pH gradient from the tumor tissue to the peritumoral normal tissue, as predicted by the above model. In addition, the growth rates of invasive and non-invasive tumor tissue from the model predictions matched clinical observations. Further analysis compared the structure and dynamics of the tumor-host interface in relation to \delta _{{1}} and the critical biological parameters. Model predictions complied with clinical observations in the following regards:

Figure 4: Interstitial gap observed in vivo
Figure 5: Interstitial gap observed in vitro

(1) When 1<\delta _{{1}}<3, the tumor-host interface is comprised of overlapping tumor and host tissue, as shown in Figure 3.

(2) When \delta _{{1}}\approx 4, there is no overlapping interface in the presence of malignant tumors.

(3) When \delta _{{1}}>>4, there will be a significant hypocellular interstitial gap between the advancing tumor edge and retreating host edge. Interestingly, the authors demonstrated the interstitial gap both in vivo (Figure 4) and in vitro (Figure 5), since no previous experiments had observed this phenomenon.

As demonstrated, the interfacial structure and speed of tumor growth depends on \delta _{{1}}. Therefore, the model may be a reliable tool for patient prognosis. There seems to be better prognosis when the overlap between tumor and normal tissue is small rather than when there is a sharp tumor-host interface. And the presence of a hypocellular interstitial gap leads to a worse prognosis than the sharp interface, where the gap negative correlates to prognosis.


Most importantly, the significant correlation between the model predictions and clinical observations implies the simple, modified reaction-diffusion model is complete for understanding invasive cancer growth as the tumor tissue reverts to glycolytic metabolism. The model makes three significant predictions supported by clinical observations:

(1) An acidic pH gradient extends from the tumor cells to the normal cells, such that acidic environments promote the growth of tumor tissue and the death of normal tissue.

(2) Invasive tumor growth develops from a crossover of noninvasive tumor growth. The transition stems from critical parameters like angiogenesis.

(3) Formation of a hypocellular interstitial gap between tumor tissue and normal tissue strongly correlates to tumor growth velocity and patient prognosis.

The authors expect new methods for diagnosing and treating tumors based on this model. They recognize the need for additional clinical and experimental investigations to further support the above predictions. Both prospective and retrospective clinical studies could test these predictions.


  1. Gatenby, R. A. The potential role of transformation-induced metabolic changes in tumor-host interaction. Cancer Res., 55: 4151-4155, 1995.
  2. Martin, G. R., and Jain, R. K. Noninvasive measurement of interstitial pH profiles in normal and neoplastic tissue using fluorescent ration imaging microscopy. Cancer Res., 54: 5670-5674, 1994.