Applied Mathematics PERGAMON
Applied Mathematics
Letters
16 (2003)
Letters
567-573
www.elsevier.nl/locate/aml
A Two-Phase Model of Solid Tumour Growth H. M.
BYRNE
AND
J. R.
KING
Centre for Mathematical Medicine, School of Mathematical Sciences University of Nottingham, Nottingham NG7 2RD, U.K.
D. L. S. MCELWAIN Centre in Statistical Science and Industrial Mathematics School of Mathematical Sciences Queensland University of Technology, Brisbane, Qld 4001 Australia
L. PREZIOSI Dipartimento di Matematica, Politecnico di Torino Corso Duca degli Abruzzi 24, 10129 Torino, Italy (Received October 2001; accepted March 2002) Communicated
by J. R. Ockendon
Abstract-Many
solid tumour growth models are formulated es systems of parabolic and/or hyperbolic equations. Here an alternative, two-phase theory is developed to describe solid tumour growth. Versions of earlier models are recovered when suitable limits of the new model are taken. We contend that the multiphase approach represents a more general, and natural, modelling framework for studying solid tumour growth than existing theories. @ 2003 Elsevier Science Ltd. All rights reserved.
Keywords-Two-phase
model, Solid tumour growth, Diffusion, Mechanical
stresses.
1. INTRODUCTION The range of mathematical approaches used to describe solid tumour growth is increasing at a dramatic rate. While stochastic and cellular automata based models have been proposed [1,2], most models are deterministic and comprise mixed systems of partial differential equations. For example, in [3], the tumour is treated as an incompressible porous medium in which diffusion of an externally-supplied nutrient leads to nonuniform cell growth and pressure differentials within the tumour. These drive a cellular velocity governed by Darcy’s law. The model is formulated as a free boundary problem, the tumour boundary moving with the local cellular velocity. In [4], reaction-diffusion equations are used to describe the invasion of tumour cells into normal tissue, without accounting for the velocity field driven by cell growth. It is sometimes difficult to identify the assumptions that such models have in common and those which are model-specific. We aim to shed light on these issues by showing how models of The authors are grateful to Queensland University of Technology, the EPSRC, the Nuffield Foundation, and the European Community for providing financial support. They also thank C. Breward, P. Howell, and C. Please for helpful discussions. 08939659/03/$ - see front PII: SO8939659(03)00038-7
matter
@ 2003 Elsevier
Science
Ltd.
All rights
reserved.
VP-et
by d.M-W
566
BYRNE et al.
H. M.
the type proposed in [3,4] may be recovered as limiting cases of a model in which the tumour is viewed as a two-phase material. A further advantage of our formulation is that it provides the starting point for incorporating any number of additional phases, as required to describe more complex aspects of solid tumour growth. For example, interactions between tumour cells, the tissue matrix, and the vasculature could be studied by including appropriately-defined phases. The paper is organised as follows. The two-phase model is formulated in Section 2. A simplified, one-dimensional model is developed in Section 3. In Section 4, we show how models similar to those of [3,4] are recovered when appropriate asymptotic limits are taken. The paper concludes with a summary in Section 5. Useful background to two-phase modelling is provided by [5] and simpler, two-phase models of solid tumour growth are described in [6,7].
2. MODEL
FORMULATION
We now formulate our mathematical model of solid tumour growth. We suppose that the tumour comprises cells and water alone and denote their respective volume fractions by Q and ,f3. As there are only two phases present, a! and ,f3satisfy the no-voids, or excluded volume, relation a:+p=1.
(1)
We associate with each phase a velocity, a pressure, and a volume-fraction-averaged stress tensor: these are denoted (v,,p,, a,) for the cell phase and (v,,,,pi, a,) for the water phase. Our mathematical model is developed by applying mass and momentum balances to the two phases. Mass Balances We treat the cell and water phases as incompressible fluids whose densities are, to leading order, equal. Consequently, mass balances for the cell and water phases may be expressed in terms of the volume fractions, o and p, in the following way:
g + V.(v,a) = I!?,,
gf + V.(vw/3)= SW= -s,.
In equations (2), S, = -S, is the volume conversion rate (i.e., the net cell proliferation rate; we thus view cells as made essentially of water). Adding equations (2) and appealing to (1)) we may rewrite the mass balance equations in the equivalent form 2 Momentum
+ V.(V&)
= s,,
V.(cw, + (1 - (Y)v,) = 0.
(3)
Balances
By assuming that inertial effects are negligible and no external forces act on the system, the momentum conservation laws reduce to force balances. Denoting by F, the force exerted on the cells by the water, by F,, the force exerted by the cells on the water, and recalling that V.(cra,) is the force exerted on the cells by themselves, we deduce that the momentum balances for the cells, the liquid, and the mixture as a whole may be written (using F, = -F,,) 0 = V.(wc)
+ F,,
(4
0 = 0.(/k,,)
+ F,,,
(5)
0 = V.((Yu, + (1 - cr)a,).
(6)
The contribution to the momentum of each species due to phase change (conversion of water to cells) has thus been neglected, it being very small, in practice, compared to F, and F,,. At this stage, our model comprises equations (3), (4), and (6) and is underdetermined. By choosing appropriate constitutive laws for CT, and uwr F,,, and S,, we now close our model, specialising it to describe avascular tumour growth.
Two-Phase
Constitutive
Model
569
Assumptions
We view the cell and water phases as incompressible fluids and assume that, on the timescale of interest, they can be treated as viscous and inviscid fluids, respectively. Hence, we set u, = -p,I + j~c (Vv, + Vvc’) + A, (V.v,) I
and
CTW= -p,I.
(7)
In (7), CL, and A, are the shear and bulk viscosity coefficients of the cell phase, and I is the identity tensor. We assume that p, is related to p, = p as follows:.
In (8), C,(cr), the pressure difference between the two phases may include contributions due to, for example, cell-cell interactions and membrane stress. By definition, we have C,(O) = 0. We assume that F,, consists of contributions die to pressure and drag, so that
Fwc = pV(1 - 0) - k(a) (vu, - v,) ,
(9)
where /c(a) is the drag coefficient. Under the above assumptions, equations (4) and (6) become v. (--PI - c&I + P,cY[vv, + vvc’] + X,a(V.v,)I) -(l - f2)Vp = k(a) (VW - vc),
= 0,
(10) (11)
and equations (3), (lo), and (11) constitute our two-phase model of solid tumour growth. In practice, the local birth and death rates of tumour cells are regulated by diffusible growth factors. Following [8], we assume that the growth rate of the tumour cells is regulated by a single nutrient which is distributed uniformly in the material surrounding the tumour and is consumed by cells as it diffuses toward the tumour centre. Denoting by C(x, t) the nutrient concentration within the tumour, and making a pseudo-steady-state approximation, it follows that C satisfies 0 = V2C - Q,(a, C).
02)
In (12), Q, 1 0 represents dhe rate at which the tumour cells consume nutrients. The pseudosteady-state approximation applies because the timescale of interest is the tumour doubling timescale (- months) which is significantly longer than the nutrient diffusion timescale (- hours); hence, the nutrient distribution rapidly adapts to changes in tumour size [g-11]. It remains to specify S,, C,, Ic, and Qc and, as representative examples, we may take
k(a) = Icocy(1- Q),
.w4
=
0,
0 I
I&la - (Y*y-1 (0 - ff*), (1 - c+
amin
ff
< I
ffmin, ff <
17
for positive constants si (i = 0,. . . ,4), Qa, &I, ko, q, T, 0 < amin < (Y* < 1, and 2,. When specifying C,, cy* denotes a natural cell density, where a! > o* cells move to reduce their stress; where cr < (Y*, they aggregate, if they are not too sparsely populated (a 2 amin). This short-range phenomenon is responsible for the initial formation of multicell spheroids from cell suspensions (see the conclusions in Section 5). We must also prescribe the domain of interest and appropriate boundary and initial conditions: we postpone their discussion until Section 3.
H. M.
570
BYRNE et al.
3. ONE-DIMENSIONAL
MODEL
Henceforth, we assume that the tumour undergoes one-dimensional growth, parallel to the z-axis, that it occupies the region -R(t) 5 z 5 R(t) at time t and is symmetric about its midline z = 0. Under these assumptions, equations (3) and (lo)-(12) reduce to give
2 + &mc) = SC, a ^ -p-cuC,+p,a, zi (
tic >
&,
+ (1 - 0)V,) = 0,
o = ,
0= g
-(l - a)2
= k(a) (VW- WC),
- Qc(a, c>,
(13) (14) (15)
where PC = 2~~ + A,. By definition, the tumour occupies the domain in which (lr > 0 and, consequently, its growth rate is governed by the cell velocity on x = R(t), so
dR - = w,(R,t). dt
(16)
In order to close equations (13)-(16), the following boundary and initial conditions are imposed:
21,= 21,= -8C = 0. 8X
P =,o,
-C,(a)
Pc$
= 0,
c = c,,
Q = m,(x), R = Ro,
at x = 0,
(17)
on x = R(t),
(18)
at t = 0.
(19)
Equations (17) ensure symmetry about x = 0. In (18), C, denotes the nutrient concentration in the medium surrounding the tumour and we assume that C is continuous across x = R(t); equation (18) guarantees continuity of the normal component of the cell and water stress tensors across x = R(t), with the ambient pressure outside the tumour normalised so that p = 0 there and the tumour boundary is stress-free. Finally, equations (19) specify the initial cell distribution within the tumour and its initial size. Writing aR(t) = a(R, t), equations (13), (16), and (18) imply that dWz aR%(oR) --&
Integrating deduce that
= &(~R,&J)
-
pc
.
the second part of (13) and the first part of (14), subject to (17) and (18), we Q V
W=--lVC
and
^
0%
P = clca --ax
&
C.
(20)
Substituting for p and v, in the second part of (14), we deduce that our one-dimensional model of solid tumour growth comprises the first of (13)) (15)) and (21) This system of equations, which is defined on a moving domain, is subject to the boundary and initial conditions for (Y, uc, C, and R stated in (16)-(19). Finally, by integrating (13) over [0, R(t)] and using (18) and (19), we obtain the following global mass balance, which will be used in Section 4:
~~R?mdx
= ~R(f)Scdx.
(22)
TwePhase
4.
Model
571
FURTHER MODEL SIMPLIFICATIONS
In this section, we show how existing models of solid tumour growth may be recovered from our model when appropriate asymptotic limits are considered. Cell Viscosity
Negligible
(Fc -+ 0)
In the limit PC --f 0, (21) may be used to eliminate v, from the first of (13), so
(23) Our model now consists of a quasi-steady reaction-diffusion equation for C coupled to a reactiondiffusion equation for cy. The associated boundary and initial conditions are
a% aa o -@=z= 9 c=c,, cr=a*,
at x = 0, on 2 = R(t), at
t=
(25)
0,
a: = ao(x:), (26) where (Y* is the positive constant such that E,(a) = 0 ( i.e., the natural packing density) and
~=-[(~)~(~~~)I.=,,,,,
withR=Ro,
att=O.
(27)
We remark that the boundary conditions for (Yat x = 0, R guarantee that the original boundary conditions for v, are also satisfied, even though fit 4 0 is a singular limit of (13)-(15). This model is similar to reaction-diffusion models of tumour growth [4,12]. The key differences here are our formulation of the problem as a moving boundary one and the fact that drag and cell-cell interaction effects are responsible for the cell ‘diffusion’ term, rather than random cellular motion (the latter may be expected to be negligible in densely-packed regions). It is possible to show that when viscous effects are neglected the pressures and velocities in each phase are related in the following way:
ah -=ax
k(a)
and
v, = -kov,
o(1 - cY)
a
o
da:
(-1 -QPC
>=
--
k(a) (1 -(.y)2
21,=
(28) Thus, in one-dimensional problems, a Darcy-type law governs the velocities of each phase. We remark also that p, and p, have opposite signs unless (Y = (Y’. Viscosity
and Drag Negligible
&, k/(c&)’
-t 0)
With viscosity negligible, (28) can be used to define vucand v,,, in terms of the cellular stress, 1 (1-g aa v w - ko,(aE,(a))t~
and
UC = - ako/(c&(a))’
%
Assuming that the velocities are finite and a! 2 amin, in the further limit k/(aC,)’ + 0, we have 2 -+ 0. Hence, the low drag (or high stiffness) of the tumour causes Q:to adjust rapidly to a constant value which, from (25), is a = cr*. The mass balance (22) then reduces to R Sc(a*, C) dx. (30) dt=z s0 Finally, we deduce that, in this limit, our model reduces to give (30), with R = & at t = 0, and dR
1
O=g$&(rr’,C), with $
= 0,
atx=O,
and
c=c,,
when x = R(t).
In particular, the nutrient distribution within the tumour, its rate of consumption by the cells, and its impact on the cell proliferation rate dictate the tumour’s evolution. Models of this form were among the earliest models of avascular tumour growth [8,9].
572
H. M. BYRNE
et al.
5. CONCLUSIONS We have used two-phase theory to develop a new mathematical model of avascular tumour growth. Our moving boundary problem, when formulated in a slab geometry, comprises a firstorder PDE for the cell volume fraction and diffusion equations for the cell velocity and nutrient concentration. We have shown how the neglect of viscous effects reduces the model to a system of reaction-diffusion equations and that, in this limit, modified versions of Darcy’s law relate the pressures and velocities of the cell and water phases. As a result, we are able to provide an alternative physical approach to formulating reaction-diffusion models of solid tumour growth (contrast [4,12], f or example). By considering limits in which both viscous effects and drag are negligible, we recover a model which is qualitatively similar to the earliest models of solid tumour growth, see [8-lo], in which spatiotemporal variations in cell volume fraction are neglected and the tumour’s growth rate is assumed to be limited by a single, externally-supplied nutrient. The scope for further model investigation can be illustrated by considering the initial aggregation in vitro of a suspension of tumour cells to form multicell spheroids which subsequently grow in line with the modelling described above. Retaining, for simplicity, the assumption of one-dimensional growth, taking uniform nutrient concentration (Qc negligible), and suppressing the C-dependence of S,, we have the possibility of spatially-uniform population growth
v = 0,
a = ao(t>,
-dao =Sc(ao). dt
However, if @(a) = -(C&(Q)) is positive for CY= ~0, then this state is unstable. Performing a routine linear stability analysis by writing Q - QO+ !I? [A(t)eixz] ,
v - R [V(t)eixz] ,
and linearising in A and V yields
v=
iX@ (ao) A K (00) + b&OX2
ad
1 dA --= A dt
@(ao)X2 KE(Qo)+ Fcao~2
+ s:(ao),
(31)
where K(QO) = k(ae)/(l - a~)~. Th e viscous term serves to regularise the ill-posed (backward diffusion) ,& = 0 problem, moderating the growth rate as X -+ co. We expect from (31) that spatial patterning will occur for & << 1 on a lengthscale x = O(j$2) during a rapid aggregation phase t = O&). A wide range of avenues for further study is suggested by.such considerations. Finally, the two-phase approach seems more physically-baaed than modelling frameworks used previously to study solid tumour growth. In addition to providing an alternative, physical interpretation of the assumptions that may underpin existing models, we believe that the two-phase approach is the most appropriate for generaiisation to include additional phases such as the extracellular matrix and vasculature. Future work will also investigate the implications of using other constitutive assumptions to close the governing mass and momentum balance equations. For example, in [7], the authors assume that in viable regions of the tumour, the cells are uniformly compacted, and therefore, replace the relation p c = p, + C, with the assumption that a = CY*(in our notation). Such generalisations and modifications should permit investigations of the interactions between tumour cells, the surrounding matrix, and the vasculature.
REFERENCES 1. W. Duchting, Tumour growth simulation, Comput. and Graphics 14, 505-508 (1990). 2. M. Baum, M.A.J. Chaplain, A.R.A. Anderson, M. Douek and J.S. Vaidya, Does breast state of chaos?, Eur. J. Cancer 35, 886-891 (1999). 3. H.P. Greenspan, On the growth and stability of cell cultures and solid turnours, J. Theor. (1976).
cancer Biol.
exist
in a
56, 229-242
Two-Phase 4. 5. 6. 7. 8. 9. 10. 11. 12. 13.
Model
573
J.A. Sherratt and M.A. Nowak, Oncogenes, anti-oncogenes and the immune response to cancer, PTOC. Roy. Sot. Land. B 248, 261-271 (1992). D.A. Drew and L.A. Segel, Averaged equations for two-phase flows, Stud. Appl. Math. 50, 205-231 (1971). K.A. Landman and C.P. Please, Tumour dynamics and necrosis: Surface tension and stability, IMA. J. Math. Appl. Med. 18, 131-158 (2001). C. Please, G.J. Pettet and D.L.S. McElwain, A new approach to modelling the formation of necrotic regions in turnouts, Appl. Math. Lett. 11 (3), 89-94 (1998). H.P. Greenspan, Models for the growth of a solid tumour by diiusion, Stud. Appl. Math. 52,317-340 (1972). J.A. Adam, A simpliied mathematical model of tumour growth, Math. Biosci. 81, 229-242 (1986). H.M. Byrne and M.A.J. Chaplain, Free boundary value problems associated with the growth and development of multicellular spheroids, EUT. J. Appl. Math. 8, 639-658 (1997). J.P. Ward and J.R. King, Mathematical modelling of avascular tumour growth, IMA J. Math. Appl. Med. Biol. 14, 39-69 (1997). R.A. Gatenby and E.T. Gawlinski, A reaction-diffusion model of cancer invasion, Cancer Res. 56, 5745-5753 (1996). A. Farina and L. Presiosi, On Darcy’s law for growing porous media, J. Nonlinear Mcch. 87,485-491 (2001).