Multiphase Modelling in Cancer Helen Byrne Wolfson Centre for - - PowerPoint PPT Presentation

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Multiphase Modelling in Cancer

Helen Byrne

Wolfson Centre for Mathematical Biology Mathematical Institute University of Oxford

CRM, Barcelona, April 2018

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Outline of Lectures

• Lecture 1: preliminaries
• Lecture 2: two-phase models
• Lecture 3: multiphase models
• Lecture 4: current work and open problems

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Lecture 1

• Biological introduction
• Moving boundary problems
• Single moving boundary
• Multiple moving boundaries
• Application to radiotherapy

Lecture 1 Lecture 2 Lecture 3 Lecture 4

How do (avascular) tumours grow?

Tumour spheroids are routinely used in experiments to characterise tumour cell lines, test drugs and learn more about cancer biology.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

How do (avascular) tumours grow?

The growth and structure of tumour spheroids are (fairly) reproducible. As such, they provide a good model of the earliest stages of avascular tumour growth.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Moving boundary problems

• Within tumour spheroid, cells cluster together and grow as they absorb

vital nutrients that diffuse through the growing tissue.

• To model growth, we need to predict nutrient concentration C within

tumour, since this determines how tumour cells n grow, and will enable us (ultimately) to determine the position of the outer boundary.

• Domain on which we solve for C grows as the tumour grows: moving

boundary problem.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Basic PDE Model

Tumour cell density, n(x, t)

∂n ∂t + ∂ ∂x (vn) advection = S(c, n) net growth rate

Nutrient concentration, c(x, t)

∂c ∂t + ∂ ∂x (vc) = D ∂2c ∂x2 diffusion − F(c, n) nutrient consumption

Outer tumour radius, R(t)

dR dt = v(R(t), t) i.e. tumour radius moves with cells on boundary (kinematic BC)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Cell velocity, v(x, t)

• Need additional equation to determine cell velocity v
• There are multiple options
• Initially, for simplicity, we consider simple phenomenological laws

Option 1: Cells move by random motion v = −µ(n) n ∂n ∂x ⇒ ∂n ∂t = ∂ ∂x

• µ(n)∂n

∂x

• + S(c, n)

⇒ nonlinear diffusion equation for cell density

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Cell velocity, v(x, t)

Option 2: ’no voids’ assumption n ≡ 1 and ∂n ∂t + ∂ ∂x (vn) = S(c, n) ⇒ ∂v ∂x = S(c, 1) ≡ ˜ S(c) ⇒ v(x, t) = x ˜ S(c)dx (by symmetry, v(0, t) = 0) ⇒ dR dt = R(t) ˜ S(c)dx

Comment: integro-differential equation for R(t)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Now make quasi-steady approximation: nutrient diffusion t-scale

• ≪

tumour doubling time-scale

• ≡ t-scale of interest

⇒ 0 ≈ D ∂2c ∂x2 − F(c, n)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Summary of Basic PDE Model

With n = 1, model reduces to coupled system for c(x, t) and R(t): 0 = D ∂2c ∂x2 − F(c, n) dR dt = R(t) ˜ S(c)dx ∂c ∂x = 0 on x = 0, c = c∞ on x = R(t) R = R0 at t = 0 This is a moving boundary problem

Comment 1: model is 1D, Cartesian geometry analogue of model of avascular tumour growth due to HP Greenspan (1972). Comment 2: Constitutive assumption used to specify v(x, t) has profound impact on model structure.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Prescribe nutrient consumption rate and cell proliferation rate by assuming

• cells consume nutrient and proliferate where nutrient levels are

high (c > c∗, say)

• cells do not consume nutrient and die when nutrient levels are

low (c < c∗) ˜ F(c) = λH(c − c∗) ˜ S(c) = s0cH(c − c∗) − s1H(c∗ − c)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Case I: c > c∗ ∀ x ∈ [0, R]. 0 = D ∂2c ∂x2 − λ ∂c ∂x = 0 on x = 0, c = c∞ on x = R(t) ⇒ c(x, t) = c∞ − λ 2D (R2 − x2) Substitute with c(x, t) in equation for R(t) dR dt = R(t) s0cdx = s0 R(t)

• c∞ − λ

2D (R2 − x2)

• dx

⇒ dR dt = s0R

• c∞ − λ

6D R2

Lecture 1 Lecture 2 Lecture 3 Lecture 4

⇒ dR dt = s0R

• c∞ − λ

6D R2

• Comment 1: model reduces to nonlinear ODE for R(t)

(solve via partial fractions) Comment 2: R = 0 is a steady state: is it stable? Comment 3: R =

• 6Dc∞/λ is a non-trivial steady state: is it

physically realistic?

Lecture 1 Lecture 2 Lecture 3 Lecture 4

dR dt = s0R

• c∞ − λ

6D R2

• and

c(x, t) = c∞ − λ 2D (R2 − x2)

Comment 3: Is steady state with R =

• 6Dc∞/λ physically realistic?

R =

• 6D

λ c∞ ⇒ cmin = c(0, t) = c∞ − λ 2D R2 ≡ −2c∞ < 0

In practice, model breaks down before then Recall that we require c > c∗ ∀ x ∈ [0, R]

cmin = c∗ = c∞ − λ 2D R2 ⇔ R = Rmin =

• 2D

λ (c∞ − c∗)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Case II: Two Free Boundaries

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Case II: two free boundaries. 0 = D ∂2c ∂x2 − λH(c − c∗) ∂c ∂x = 0 on x = 0, c = c∞ on x = R(t) c = c∗ on x = R∗(t) c, ∂c ∂x continuous across x = R∗(t) ⇒ c(x, t) =

• c∗ +

λ 2D(x − R∗)2

for R∗ < x < R c∗ for 0 < x < R∗ with c∞ − c∗ = λ 2D (R − R∗)2

Comment: (R − R∗) = width of proliferating rim is fixed

Lecture 1 Lecture 2 Lecture 3 Lecture 4

dR dt = R {s0cH(c − c∗) − s1H(c∗ − c)} dx = R

R∗ s0cdx −

R∗ s1dx where c(x, t) =

• c∗ +

λ 2D(x − R∗)2

for R∗ < x < R c∗ for 0 < x < R∗ ⇒ dR dt = s0(R − R∗)

• c∗ + λ

6D (R − R∗)2

• − s1R∗

where R − R∗ =

• 2D

λ (c∞ − c∗) = Λ, constant

Eliminate R∗ = R − Λ:

⇒ dR dt = s0Λ 3 (c∞ + 2c∗) − s1(R − Λ)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Case I: dR dt = s0R

• c∞ − λ

6D R2

• valid for 0 < R < Rmin =
• 2D

λ (c∞ − c∗)

Case II: (R > Rmin) dR dt = s0Λ 3 (c∞ + 2c∗) − s1(R − Λ) where R∗ = R −

• 2D

λ (c∞ − c∗) ≡ R − Λ and R → Λ

• 1 + s0

3s1 (c∞ + 2c∗)

• as t → ∞

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Comments

• Simple model admits nontrivial steady state
• Explicit expression for width of proliferating rim and R(t)
• Can fit analytical expressions to experimental data to estimate

model parameters

• 1. Initial growth rate of spheroid (0 < R ≪ 1)
• 2. Width of proliferating rim
• 3. Rate of approach to steady state
• 4. Steady state size of tumour spheroid

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Radiotherapy

The Linear-Quadratic Model: Survival fraction, SF(D∗) = exp(−αD∗ − βD∗2) where α and β are tissue properties, and D∗ = dose (Grays). Note: Proportion of cells killed by RT = 1 − SF(D∗).

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Incorporating Radiotherapy: No Necrosis

Suppose 0 < R < Rmin = (2D(c∞ − c∗)/λ)1/2 so that R∗ = 0. Then dR dt = s0R

• c∞ − λR2

6D

• − (1 − SF(D∗))Σn

i=1δ(t − ti)R(t)

where SF(D∗) = exp(−(αD∗ + βD∗2)). Let A(τ) = r 2. Then dA dt = s0A 2

• c∞ − λA

6D

• − (1 − SF(D∗))Σn

i=1δ(t − ti)A(t)

Exercise: calculate A(t). Question: what happens if necrosis present (i.e. R > Rmin and R∗ > 0)?

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Radiotherapy and Necrosis

• Suppose necrosis present: R > Rmin and R∗ > 0.
• Question: how do we adapt model to account for radiotherapy?
• Irradiate at times t = ti (i = 1, 2, . . . , n)
• Denote by R+/−

i

and R∗

i +/− positions of outer tumour radius and

necrotic boundary just before and after radiotherapy.

• Assume radiotherapy targets proliferating cells only so that

R∗+ = R∗− ≡ R∗, say ⇒ width of proliferating rim just after RT

• = SF(D∗)(R−

i

− R∗

i −) < (R− i

− R∗

i −)

⇒ c(R∗) > c∗ ⇔ mismatch in nutrient profile

• n necrotic boundary

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Radiotherapy and Necrosis

Following RT, there may be mismatch in nutrient profile on necrotic boundary. Question: how should we proceed while c(rN) > cN?

Lecture 1 Lecture 2 Lecture 3 Lecture 4

• Question: how should we proceed while c(R∗) > c∗?
• Assume
• proliferating rim (R∗ < x < R) continues to grow
• necrotic core (0 < x < R∗) continues to degrade

dR∗ dt = −s1R∗ ⇒ R∗(t) = R∗(ti)e−s1(t−ti ) d dt (R − R∗) = s0 R

R∗ c dx

⇒ dR dt = s0 R

R∗ c dx − s1R∗

where c(x, t) = B1 0 ≤ x < R∗(t), B2 + A2x +

λ 2D x2

R∗(t) ≤ x ≤ R(t) and B1, B2 and A2 are chosen such that c = c∞ on x = R(t) and c, ∂c ∂x cts across x = R∗(t)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Radiotherapy and Necrosis

After some algebra, we deduce c(x, t) = B1 0 ≤ x < R∗(t), B2 + A2x +

λ 2Dx2

R∗(t) ≤ x ≤ R(t) where dR∗ dt = −s1R∗ ⇒ R∗(t) = R∗(ti)e−s1(t−ti) dR dt = s0 R

R∗

• c∞ − λ

2D (R − x) [(R − R∗) + (x − R∗)]

• dx − s1R∗

until (i) c = c∗ on x = R∗ (ii) next round of radiotherapy Exercise: plot dynamics of R(t) and R∗(t).

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Radiotherapy and Necrosis (3D Radial Symmetry)

Analytical and simulation results showing evolution of necrotic tumour following treatment with single dose of radiotherapy. Clear evidence that nutrient mismatch on necrotic boundary alters tumour dynamics

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Radiotherapy and Necrosis (3D Radial Symmetry)

Simulation results showing response of two tumours to standard RT treatment over period of 6 weeks (2 Grays per day, Mon-Fri). LHS: tumour shrinks during treatment. RHS: tumour shrinks initially, but eventually evolves to periodic solution, with clear evidence of necrotic mismatch following RT.

TD Lewin, H Enderling et al (2018) Bull Math Biol

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Critique

Recall that biological tissues are complex and heterogeneous How can we relax following assumptions in our PDE model?

• Heterogeneity
• In vivo tissues comprise ECM, immune cells, blood vessels, and

fluid

• Geometry: 1D → 2D, 3D
• Tumours are often invasive, spreading in 2D and 3D
• Mechanics:
• Tumour cells may be mechanosensitive

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Outline of Lectures

• Lecture 1: preliminaries
• Lecture 2: two-phase models
• Lecture 3: multiphase models
• Lecture 4: current work and open problems

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Lecture 2

• Heterogeneity: phenomenological multiphase model for

macrophage-based therapies

• Geometry
• Extension to 2D/3D
• Symmetry breaking and invasion
• Towards mechanics: two-phase model

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Heterogeneity

Consider tumour that comprises

• tumour cells, n(x, t)
• Dead or waste material, w(x, t)
• Macrophages, m(x, t)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Aside: Macrophages and Cancer

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Aside: Macrophages and Cancer

Macrophages tend to localise in low oxygen (hypoxic) tumour regions. Can we exploit this to deliver treatment to these tumour regions?

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Governing Equations

Tumour cells, n(x, t)

∂n ∂t + ∂ ∂x (nv) = µn ∂2n ∂x2 + Sn(n, c, w)

Waste material, w(x, t)

∂w ∂t + ∂ ∂x (wv) = µw ∂2w ∂x2 + Sw(n, m, c, w)

Macrophages, m(x, t)

∂m ∂t + ∂ ∂x (mv) = µm ∂2m ∂x2 − χ ∂ ∂x

• m ∂a

∂x

• + Sm(m, a)

where c(x, t) = externally-supplied nutrient and a(x, t) = macrophage chemoattractant, produced by tumour cells under hypoxia

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Nutrient, c(x, t)

0 = Dc ∂2c ∂x2 − (Γ0 + Γ1n + Γ2m)c

Chemoattractant, a(x, t)

0 = Da ∂2a ∂x2 + Sa(c, n) − Γ3a

Velocity, v(x, t)

∂v ∂x = µn ∂2n ∂x2 + µw ∂2w ∂x2 + µm ∂2m ∂x2 − χ ∂ ∂x

• m ∂a

∂x

• + Sn + Sw + Sm

i.e. ’no voids’ assumption (n + w + m = 1)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

If Sn + Sw + Sm = 0 then ∂v ∂x = µn ∂2n ∂x2 + µw ∂2w ∂x2 + µm ∂2m ∂x2 − χ ∂ ∂x

• m ∂a

∂x

• ⇒ v = µn ∂n

∂x + µw ∂w ∂x + µm ∂m ∂x − χ

• m ∂a

∂x

• and dR

dt = v(R(t), t) PDE for n(x, t) becomes nonlinear reaction-advection-diffusion equation: ∂n ∂t = µn ∂ ∂x

• (1 − n)∂n

∂x

• − ∂

∂x

• n
• µw ∂w

∂x + µm ∂n ∂x − χ(m ∂a ∂x )

• +Sn(n, c, w)

where, for example, Sn(n, c, w) = ˜ S(c)nw

• proliferation

− En(c)n death − An(a)nm

• macrophage killing

References: Ward and King (1997); Webb et al (2007); Chen et al (2014)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Numerical Results

50 100 150 200 250 500 1000

Time Tumour

0.2 0.4 0.6 0.8 50 100 150 200 250 500 1000

Time M−φe

0.1 0.2 50 100 150 200 250 500 1000

Time Prodrug

1 2 3 x 10

−3

50 100 150 200 250 500 1000

Distance from centre Time Drug

0.02 0.04 0.06 0.08

Macrophage therapy slows tumour growth

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Numerical Results

50 100 150 200 250 500 1000

Time Tumour

0.2 0.4 0.6 0.8 50 100 150 200 250 500 1000

Time M−φe

0.1 0.2 50 100 150 200 250 500 1000

Time Prodrug

5 x 10

−3

50 100 150 200 250 500 1000

Distance from centre Time Drug

0.05 0.1 0.15

Macrophage therapy controls tumour growth

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Bifurcation Diagram

0.05 0.1 0.15 0.2 0.25 0.3 0.05 0.1 0.15 0.2 0.25

Wave speed, U Surface prodrug (φ∞) 10

1

10

2

10

3

10

4

Saturation size, R∞ 20 40 60 80 100 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9

Distance from centre Tumour cell volume fraction φ∞=0.1 φ∞=0.15 φ∞=0.2

System evolves to steady travelling wave when pro-drug levels are low System evolves to steady state when pro-drug levels high (tumour control)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Bifurcation Diagram

0.14 0.16 0.18 0.2 0.22 0.24 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Wave speed, U Surface prodrug (φ∞) 10

1

10

2

10

3

Saturation size, R∞ 20 40 60 80 100 10 20 30 40 50 60 Time, t Spheroid radius, R

In certain parameter regimes we observe co-existence of travelling waves and steady state solutions

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Next Steps with Macrophages

• Combination therapies (e.g. macrophage-therapy and

radiotherapy)

• Distinguish different macrophage phenotypes (M1,M2, etc)
• Understand/explain why distributions of macrophages in different

tumours vary

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Modelling Combination Therapy

• 3-phase mixture model (tumour cells, macrophages and fluid)
• apply single dose of radiotherapy at t=150
• introduce macrophage-based virotherapy at t=200
• decompose tumour population into uninfected and infected tumour cells

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Modelling Combination Therapy

• Treatment protocol influences outcome

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Modelling Combination Therapy

• Coordination of RT and virotherapy influences outcome

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Critique of Basic Model

Recall that biological tissues are complex and heterogeneous How can we relax following assumptions in our simple PDE model?

• Heterogeneity
• In vivo tissues comprise ECM, immune cells, blood vessels, and

fluid

• Geometry: 1D → 2D, 3D
• Tumours are often invasive, spreading in 2D and 3D
• Mechanics:
• Tumour cells may be mechanosensitive

Lecture 1 Lecture 2 Lecture 3 Lecture 4

2D Model

Velocity v = (v, w) ∂n ∂t + ∂ ∂x (vn) + ∂ ∂y (wn) = S(c, n) 0 = D ∂2c ∂x2 + ∂2c ∂y2

• − F(c, n)

dΓ dt = v.ˆ n (ˆ n = unit outward normal)

Option 1: cells move by random motion v = −µ∇n: diffusion equation for n Option 2: ’no voids’ (n = 1) ⇒ not enough information to close problem! Possible closures: introduce cell pressure, p(x, y, t)

Darcy’s law: v = −k∇p Stoke’s flow: µ∇2v − ∇p = 0

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Darcy’s Law and ’No Voids’

∇.v = S(c, 1) = ˜ S(c) and v = −k∇p ⇒ −k∇2p = ˜ S(c) Model becomes 0 = D∇2c − ˜ F(c) and − k∇2p = ˜ S(c) dΓ dt = v.ˆ n = −k ∂p ∂n with Γ(x, y, t) = 0 = x − R(y, t). Boundary and initial conditions ∂c ∂x = 0 = ∂p ∂x at x = 0 c = c∞ p = p∞ on x = R(y, t) R(y, 0) = R0(y) prescribed

• Darcy: Greenspan (1976); HMB and Chaplain (1997); HMB (1997).
• Stokes flow: Franks et al (2003)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Symmetry Breaking and Tumour Invasion

• Calculate radially symmetric steady state solution
• Linearise about steady state to determine stability to symmetry-breaking

perturbations (HMB and Chaplain, 1997)

• Can perform weakly nonlinear analysis to characterise bifurcations but

the algebra is very messy (see: Byrne, 1999)

• Happy to share details if there is interest

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Critique of Basic Model

Recall that biological tissues are complex and heterogeneous How can we relax following assumptions in our simple PDE model?

• Heterogeneity
• In vivo tissues comprise ECM, immune cells, blood vessels, and

fluid

• Geometry: 1D → 2D, 3D
• Tumours are often invasive, spreading in 2D and 3D
• Mechanics:
• Tumour cells may be mechanosensitive

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Mechanotransduction

Spheroid size influenced by stiffness of gel in which it is cultured

Lecture 1 Lecture 2 Lecture 3 Lecture 4

General model framework (1D)

• View tumour as mixture of cells (n) and water (w)
• Mass balances for n and w:

∂n ∂t + ∂ ∂x (nvn) = Sn, ∂w ∂t + ∂ ∂x (wvw) = −Sn

• Momentum balances for vn and vw:

∂ ∂x (nσn)+knw(vw−vn)+p ∂n ∂x = 0 = ∂ ∂x (wσw)−knw(vw−vn)+p∂w ∂x

• No voids condition (for pressure, p):

n + w = 1

• Close model by specifying Sn, σn and σw

Lecture 1 Lecture 2 Lecture 3 Lecture 4

General model framework

• Assume n and w are isotropic fluids:

σn = −pn, σw = −pw ≡ −p

• Additional constitutive assumption needed to specify pn:
• 1. cells ∼ bags of water:

pn = p + Σ(n)

• 2. cells ∼ incompressible:

n = n∗, constant

Sketch of Σ(φ). Schematic of cell and fluid movt.

• Conservation of total mass and momentum supplies:

nvn + wvw = 0 = ∂ ∂x (nσn + wσw)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

General model framework

• Assume n and w are isotropic fluids:

σn = −pn, σw = −pw ≡ −p

• Additional constitutive assumption needed to specify pn:
• 1. cells ∼ bags of water:

pn = p + Σ(n)

• 2. cells ∼ incompressible:

n = n∗, constant

Sketch of Σ(φ). Schematic of cell and fluid movt.

• Conservation of total mass and momentum supplies:

nvn + wvw = 0 = ∂ ∂x (nσn + wσw)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

General model framework

• Assume n and w are isotropic fluids:

σn = −pn, σw = −pw ≡ −p

• Additional constitutive assumption needed to specify pn:
• 1. cells ∼ bags of water:

pn = p + Σ(n)

• 2. cells ∼ incompressible:

n = n∗, constant

Sketch of Σ(φ). Schematic of cell and fluid movt.

• Conservation of total mass and momentum supplies:

nvn + wvw = 0 = ∂ ∂x (nσn + wσw)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Model simplifications

• I. pn = p + Σ(n).

⇒ ∂n ∂t = 1 k ∂ ∂x

• (1 − n) ∂

∂x (nΣ(n))

• + Sn

i.e. Nonlinear diffusion equation for n with Darcy’s law for fluid motion so that vw = − 1

k px

(cf. Breward et al, 2002; Byrne et al, 2003)

• II. n = n∗, constant.

∂vn ∂x = Sn and vn = −1 k ∂pn ∂x ⇒ −1 k ∂2pn ∂x2 = Sn (Laplace’s Eqn.) Assume outer tumour boundary (x = R(t)) moves with cell velocity dR dt = vn|x=R(t) = R Sndx

(cf. Greenspan, 1972) Retain pn to allow Sn = Sn(pn)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Model simplifications

• I. pn = p + Σ(n).

⇒ ∂n ∂t = 1 k ∂ ∂x

• (1 − n) ∂

∂x (nΣ(n))

• + Sn

i.e. Nonlinear diffusion equation for n with Darcy’s law for fluid motion so that vw = − 1

k px

(cf. Breward et al, 2002; Byrne et al, 2003)

• II. n = n∗, constant.

∂vn ∂x = Sn and vn = −1 k ∂pn ∂x ⇒ −1 k ∂2pn ∂x2 = Sn (Laplace’s Eqn.) Assume outer tumour boundary (x = R(t)) moves with cell velocity dR dt = vn|x=R(t) = R Sndx

(cf. Greenspan, 1972) Retain pn to allow Sn = Sn(pn)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Remarks

• Structure of phenomenological models depends on underlying

constitutive assumptions

• Multiphase framework allows more general constitutive assumptions (eg

viscous or elastic effects) and investigation of biomechanical effects

• Can use multiphase framework to study problems in e.g. tissue

engineering

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Applications/Extensions

• Including chemotaxis within multiphase framework
• Constrained growth (mechanotransduction)
• Ductal carcinoma in situ
• Tumour encapsulation
• Growth of spheroids embedded in gels

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Outline of Lectures

• Lecture 1: preliminaries
• Lecture 2: two-phase models
• Lecture 3: multiphase models
• Lecture 4: current work and open problems

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Lecture 3

Choose two of these applications

• Tumour encapsulation
• Chemotaxis in a multiphase framework
• Model development
• Application: going against the flow (with chemotaxis)
• Mechanical inhibition of spheroid growth

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Tumour Encapsulation

Modelling approach:

• Develop three-phase model for tumour cells (n), extracellular

fluid (w) and ECM (m)

Jackson and Byrne, Math Biosci (2002)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Governing Equations

Mass Balance Equations nt + (nvn)x = Sn, mt + (mvm)x = Sm, wt + (wvw)x = Sw. No-Voids Assumption n + m + w = 1. Momentum Balance Equations (nσn)x + Fnm + Fnw + pnx = 0, (mσm)x − Fnm + Fmw + pmx = 0, (wσw)x − Fnw − Fmw + pwx = 0, where Fij = knm(vj − vi), and σn = −(p + Σn), σm = −(p + Σm), σw = −p.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Governing Equations

Constitutive Assumptions Sn = αnnw −δn and Sm = αmnmw Σn = snn and Σm = smm(1+θn). Note:

• Pressures in tumour and ECM phases increase with their densities
• Tumour cells cause additional increase in ECM pressure (θ > 0)
• ECM only produced if tumour cells present

After some (!) algebra, possible to derive reduced model for n and m: ∂n ∂t = n(1−n−m)−δn+2βn ∂ ∂x

• n(1 − n)∂n

∂x

• −βm ∂

∂x

• n ∂

∂x (m2(1 + θn))

• ,

∂m ∂t = αmn(1 − n − m)+βm ∂ ∂x

• (1 − m) ∂

∂x (m2(1 + θn))

• −2βn ∂

∂x

• mn ∂n

∂x

• ,

where βn = sn αnkL2 , βm = sm αnkL2 , α = αm αn , δ = δn αn .

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Tumour Encapsulation

Simulation results showing how system evolves when active and passive responses are active (α, θ > 0)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Simulation Results

Simulation results showing how system evolves when only the passive response is active (α = 0, θ > 0)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Simulations Results

Simulation results showing how the capsule width changes over time when

• nly a passive response is active, and βn varies.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Simulation Results

Simulation results showing how system evolves when only an active response is active (α > 0, θ = 0).

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Simulation Results

Simulation results showing how protease production may lead to the destruction of the collagen capsule and, thereby, enhance tumour invasion. Tumour cells produce ECM-degrading protease when the pressure they experience exceeds a threshold value.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

The Patlak-Keller-Segel Equations for Chemotaxis

∂n ∂t = ∂ ∂x

• µ(n, a)∂n

∂x

• − ∂

∂x

• χ(n, a)n ∂a

∂x

• + f(n, a)

∂a ∂t = ∂2a ∂x2 + g(n, a)

Question:

• Can we derive PKS equations using multiphase framework?

Ref: Owen and Byrne (2004), J Math Biol

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Mass Balance Equations

nt + (vnn)x = Sn(n, w, a) wt + (vww)x = −Sn(n, w, a) at + (avw)x = Daaxx + Sa(a, n) Note: a distributed in fluid phase

No voids

n + w = 1

Momentum Balance Equations

0 = (nσn)x + Fnw + pnx 0 = (wσw)x − Fnw + pwx

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Constitutive Equations

σn = −p − Λ(n, a) and σw = −p Overall system momentum balance supplies 0 = (p + nΛ)x ⇒ p = −nΛ + P(t) While system mass balance supplies 0 = (nvn + wvw)x ⇒ nvn + wvw = 0 With Fnw = k(n)(vw − vn), momentum balance for fluid supplies 0 = −wpx + k(n)(vn − vw)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

0 = −wpx + k(n)(vn − vw) Substitute for p = −nΛ + P(t) and vw = nvn/(1 − n) to get vn = −(1 − n)2 k(n) (nΛ)x ⇒ ∂n ∂t = ∂ ∂x

• µ(n, a)∂n

∂x

• − ∂

∂x

• χ(n, a)n ∂a

∂x

• + f(n, a)

where µ(n, a) = n(1 − n)2 k(n) ∂ ∂n(nΛ) and χ(n, a) = −n(1 − n)2 k(n) ∂Λ ∂a

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Comment 1: vw = − nvn 1 − n ⇒ cells and water move in opposite directions Comment 2: ∂Λ ∂a < 0 ⇒ pressure in cell phase relieved when cells move up chemical gradients Comment 3: chemokinesis ∂Λ ∂a = 0 ⇒ chemical stimulates chemotaxis (χ = χ(n, a))

and regulates random motion (µ = µ(n, a)

Comment 4: Λ + nΛn < 0 ⇒ ill-posedness

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Linear Stability Analysis

• Neglect cell proliferation and death
• assume cells respond to concentration of chemical in the water,

a/w:

Sn = 0, Sa = n − s1a, Λ = Λ(a/w) = α1 + β1e−θ1a/w = α1 + β1e−θ1a/(1−n) nt = n(1 − n)2 k(n) (nΛ)x

• x

at + an(1 − n) k(n) (nΛ)x

• x

= Daaxx + n − s1a

We linearise about the spatially uniform steady state

n = n0 + ǫn1eiλx+σt and a = n0 s1 + ǫa1eiλx+σt (0 ≪ ǫ)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Pattern formation if ℜ(σ(λ)) > 0 if α1 β1 <

• θ1n0

s1(1 − n0)2 − 1

• e−θ1n0/s1(1−n0) ≤ e−2
• 1 + 4s1

θ1

• Comment 1:

α1 β1 > e−2

• 1 + 4s1

θ1

• ⇒

no patterns Comment 2: α1 β1 < e−2

• 1 + 4s1

θ1

• ⇒

patterns for n0 ∈ (n−, n+) where n± solve α1 β1 =

• θ1n0

s1(1 − n0)2 − 1

• e−θ1n0/s1(1−n0)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Analytical Results

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Analytical Results

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Interstitial Flow Biases Cell Migration

Experimental set-up Experimental results Interstitial flow biases cell migration against the direction of flow

WJ Polacheck et al., PNAS (2011)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Interstitial Flow Biases Cell Migration

Aim

• Identify mechanisms by which flow may enhance cell aggregation
• Are experimental results consistent with cells producing

chemoattractant

• At a constant rate?
• At a rate which depends on drag experienced by cells?

Acknowledgements: CG Bell, JP Whiteley, SL Waters

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Model Development

Schematic diagram of experimental set up

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Model Development

Assumptions

• Two-phase mixture model for fluid and cells
• Mass and momentum balances to determine
• Cell and fluid volume fraction (θn, θw)
• Cell and fluid velocities (uw, un)
• Cells produce chemoattractant which degrades, diffuses and is

advected with the flow

• Neglect cell proliferation
• Additional constitutive assumptions to close model

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Mass balance equations

∂θn ∂t + ∇.(θnun) = 0, ∂θw ∂t + ∇.(θnuw) = 0

Momentum balance equations

= ∇.(θnσn) + Fnw + Fns + p∇θn, = ∇.(θwσw) − Fnw + Fws + p∇θw, where σn = −Pn = −(P + Dn ˜ ψ(θn) − γ ˜ χ(c) and σw = −P. Note: scaffold static but exerts drag force on cells and fluid

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Mass balance equations

∂θn ∂t + ∇.(θnun) = 0, ∂θw ∂t + ∇.(θnuw) = 0

Momentum balance equations

= −θn∇P Mixture pressure − Dn∇ψ(θn)

• Cell pressure

+ γχ(c)θn∇c

• Chemotactic force

− βwnθnθw(un − uw)

• Fluid-cell drag

− βnsθsθnun

• Cell-scaffold drag

0 = −θw∇P

• Mixture pressure

+ βwnθnθw(un − uw)

• Fluid-cell drag

− βwsθsθwuw

• Fluid-scaffold drag

Note: Pn = P + Dn ˜ ψ(θn) − γ ˜ χ(c)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Chemoattractant

∂(θwc) ∂t + ∇.(θwuwc) = D∇.(θw∇c) + κθnθw − λθwc

Constitutive Assumptions

θn + θw = 1

Approach

• 1D, pressure-driven fluid flow across cells
• Perform linear stability analysis of spatially-uniform steady state

θn = θ0 + ǫθ1eiqx+σt, P = Pu

• 1 − x

L

• + ǫP1eiqx+σt, c = κθ0

λ + ǫc1eiqx+σt

• Identify conditions under which flow enhances cell aggregation
• Compare linear stability analysis with numerical solutions of full model

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Results (Linear Stability Analysis)

Flow inhibits cell aggregation when cells produce chemoattractant at a constant rate.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Results (Linear Stability Analysis)

Flow enhances cell aggregation when chemoattractant production is drag-dependent. Production is increasing saturating function of absolute value of drag force between cells and fluid.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Applications of Multiphase Modelling

• Including chemotaxis within multiphase framework
• Constrained growth (mechanotransduction)
• Ductal carcinoma in situ
• Tumour encapsulation
• Growth of spheroids embedded in gels

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Model Development

• Tumour ≡ 2-phase mixture: cells (φ) and water (1 − φ)
• Tumour embedded in deformable gel

Mass Balance Equations φt + ∇.(φUc) = φS(c) (1 − φ)t + ∇.[(1 − φ)Ue] = −φS(c) Note: no voids assumption implicit in these equations Momentum Balance Equations 0 = ∇.(φσc) − 1 k φ(Uc − Ue) + Pe∇φ 0 = ∇.((1 − φ)σe) + 1 k φ(Uc − Ue) + Pe∇(1 − φ) where σc = −(Pe + 1 φPc) and σe = −Pe

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Model Development

Mass Balance Equations φt + ∇.(φUc) = φS(c) (1 − φ)t + ∇.[(1 − φ)Ue] = −φS(c) Momentum Balance Equations (Landman and Please, 2000) − 1 k φ(Uc − Ue)

• drag

− φ∇Pe hydrostatic pressure − ∇Pc

• intercell. pressure

= 0 1 k φ(Uc − Ue) + φ∇Pe − ∇Pe = 0 Sum mass and momentum balance equations φUc = −(1 − φ)Ue and Pc = −Pe + p(t)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Eliminate Pe and Ue from cell momentum balance ⇒ Uc = −k (1 − φ)2 φ ∇Pc Assume cells compacted in viable region Pc > Pe and φ = φ0 Assume cells not compacted in necrotic core Pc = Pe and 0 < φ < φ0 ⇒ Pc − Pe ≥ 0, φ0 − φ ≥ 0, (Pc − Pe)(φ0 − φ) = 0. Note: Analogy with models of fluidised beds Assume radial symmetry, with r = X(t) =

• uter tumour radius

r = L(t) = interface between necrotic and compacted regions

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Model Summary

Compacted/viable region (L(t) < r < X(t)) 0 = 1 r 2 (r 2Cr)r − S(c), φ = 1 1 r 2 (r 2Uc)r = − 1 r 2 (r 2Pcr)r = H(C − α) − ρH(α − C) Necrotic region (0 < r < L(t)) 0 = 1 r 2 (r 2Cr)r − φ S(c), φt = φ[H(C − α) − ρH(α − C)], Pc = p0(t), Uc = 0 Boundary and continuity conditions Uc = −Pcr = Cr = 0

• n r = 0

C = 1, Pc = p(t) = Γ X − σrr(X, t), Pe = 0, dX dt = Uc on r = X(t). [C]+

− = [Cr] = [φ] = [Uc] = [Pcr] = 0 across r = L(t)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Exponential Growth Phase

Non-necrotic growth (L = 0, S(C) = 1) C = 1 − 1 6(X 2 − r 2) Uc = r 3, X(t) = X0et/3, Pc = Γ X − σrr(X, t) − 1 6(X 2 − r 2) Comments:

• Stress in gel does not affect growth during this stage of

exponential growth

• Stress will affect growth when necrosis is initiated
• Model breaks down at t = t1 when

Cmin = C(0, t1) = α ⇔ X = X(t1) =

• 6(1 − α)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Sub-exponential, Non-necrotic Growth Phase

Suppose C = α on r = rc(t) where 0 < rc(t) < X(t). Then, 3X 2 dX dt = X 3 − (1 + ρ)r 3

c

and r 2

c = X 2 − 6(1 − α)

⇒ X → X∞ = 6(1 − α) (1 + ρ)2/3 (1 + ρ)2/3 − 1 as t → ∞ Solution valid if Pc > Pe ∀ r ∈ (0, X) ⇔ −1 6(2 + 3ρ)X 2 + 3(1 + ρ)(1 − α) + (1 + ρ) 3 r 3

c

X 2 + 1 2 Γ X − σrr(X, t)

• > 0

Constraint depends on stress in gel (also require C > 0) Can repeat analysis for growth with necrotic core but algebra more involved Question: What happens in the gel?

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Steady States and Travelling Waves

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Aside: Elasticity Theory

Let Xi = reference coordinate of a particle and xi = spatial coordinate (deformed frame). Then Fi,j = ∂xi ∂Xj ≡ deformation gradient tensor B = F.FT ≡ Cauchy-Green deformation tensor The (stress or) strain invariants I1, I2 and I3 are intrinsic properties of the (stress or) strain tensor, regardless of the frame of reference. They are the coefficients of the characteristic equation that defines the eigenvalues of the (stress or) strain tensor. For example, I1 = tr(B), I2 = (trB)2 − tr(B2), I3 = det(B). The strain energy function W relates the strain energy density to the deformation gradient

Lecture 1 Lecture 2 Lecture 3 Lecture 4

We define a strain energy function W = η

• eβ(¯

I1−3) − 1

• + γ (¯

I3 − 1)2 (¯ I3 − δ)n where ¯ I1 = I−1/2

3

I1, ¯ I2 = I−2/3

3

I2, ¯ I3 = I3. I1 = tr(B), I2 = (trB)2 − tr(B2), I3 = det(B). Comment 1: no stress at no displacement Comment 2: For isotropic, hyperelastic material, W = W(I1, I2) and σ = 2 √I3

• I3

∂W ∂I3 I + ∂W ∂I1 B

• Neglecting inertial terms, the equations of motion are

∇.σ = 0

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Radial Symmetry and Small Displacements

We denote by r and R the spatial and reference coordinates. Then F =  

∂r ∂R r R r R

  , σ =   σrr σθθ σφφ   where σrr = 2

• I3

∂W ∂I3 + 2 √I3 ∂W ∂I1 ∂r ∂R 2 , σθθ = σφφ = 2

• I3

∂W ∂I3 + 2 √I3 ∂W ∂I1 r R 2 , and the equations of motion (∇.σ = 0) supply ∂σrr ∂r + 2 r (σrr − σθθ) = 0 and σθθ = σφφ

Lecture 1 Lecture 2 Lecture 3 Lecture 4

For small displacements, u = r − R ≃ ǫu0 ≪ 1 ⇒ ∂r ∂R = 1 + ∂u ∂R ≃ 1 + ǫ∂u0 ∂R I1 ≈ 3 + 2ǫ ∂u0 ∂R + 2u0 R

• and

I3 ≈ 1 + 2ǫ ∂u0 ∂R + 2u0 R

• and, after some algebra,

σrr ≃ 8ηβ 3 ∂u0 ∂R − u0 R

• +

8γ (1 − δ)n ∂u0 ∂R + 2u0 R

• σθθ ≃ −4ηβ

3 ∂u0 ∂R − u0 R

• +

8γ (1 − δ)n ∂u0 ∂R + 2u0 R

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Finally, the equations of motion supply ∂2u0 ∂R2 + 2 R ∂u0 ∂R − u0 R

• = 0

with u0 → 0 as R → ∞ and u0 = 1 on R = X0 (i.e. on the tumour boundary) ⇒ u0 = X 2 R2 and σrr(X0, t) ≃ −8ηβX 2 R3

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Displacement and Radial Stress in Gel

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Effect of Varying Gel Stiffness

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Steady States and Travelling Waves

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Spheroid Regrowth on Return to Free Suspension

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Spheroid Regrowth on Return to Free Suspension

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Outline of Lectures

• Lecture 1: preliminaries
• Lecture 2: two-phase models
• Lecture 3: multiphase models
• Lecture 4: current work and open problems

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Lecture 4

• Vascular tumour growth
• Open Problems
• Modelling assumptions eg constitutive laws for different phases
• Homogenisation
• (ECM) anisotropy

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Open Problems

Appearance or disappearance of a phase

• Elimination of a tumour
• The emergence of a new, subpopulation of mutant cells; the

arrival of immune cells

Efficient numerical methods

• Extensions to 2D and 3D
• Multiple phases with different properties (eg extracellular fluid,

anisotropic, visco-elastic tissue matrix, compliant blood vessels, tumour cells, immune cells, fibroblasts)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Vascular Tumour Growth

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Multiple Phases and 2D

• View tumour as mixture of tumour, normal cells, fluid and vessels
• Mass and momentum balance for phase fractions and velocities:

∂θi ∂t + ∇.(θiui) = qi and ∇.(θiσi) + Σj=iFij = 0 (i = 1, . . . , 4)

• Assume ’no-voids’ and that system is closed

Σ4

i=1θi = 1

and Σ4

i=1qi = 0

• Diffusion equation for nutrient, c:

Dc∇2c + k5γ(cv − c) + qc = 0 ⇒ mixed system of PDEs solved in 2D using Finite Element Methods

Ref: Hubbard and Byrne (2013)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Typical Simulation Results: Early Times

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Typical Simulation Results: Early Times

Corresponding nutrient profiles

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Typical Simulation Results: Later Times

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Merging of Two Tumours

At t = 0, tissue seeded with 2 tumours of radius 1

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Tumour Growth in a Tapered Domain

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Tissue Engineering

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Tissue Engineering of Articular Cartilage

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Mechanics of Reinforced Hydrogels

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Scaffold Details

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Modelling Approach

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Modelling Approach (ctd)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Microscale Problem

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Macroscale, Effective Parameters

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Macroscale Governing Equations

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Relaxation Experiment

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Initial Simulation Results vs Experiments

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Model Experiment Discrepancies

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Adjusted Results

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Summary

• Used homogenisation theory to derive effective macroscale

properties of periodic, elastic-poroelastic composite = orthotropic material

• Determine dependence of macroscale parameters on microscale

geometry and material properties of constituents

• Discrepancies between initial model simulations and experiments

(due to eg sagging fibres, hydrogel layer) can be corrected for

• Model can be used to predict stress experienced by cells under

mechanical loading

Refs:Visser et al, Nature Comms (2015); MJ Chen et al, arXiv (2017)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

ECM Anisotropy: Model Overview

Deformation changes fibre alignment Fibres guide cell movement Cells exert force on ECM ECM deforms

The influence of cell-matrix interactions and tissue anisotropy on tissue architecture

Ref: Dyson, Green, Whiteley and Byrne (2016)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Model Structure

• Three-phase mixture model: cells, collagen and fluid
• Neglect cell proliferation and death
• View cells as incompressible viscous fluid
• View collagen as incompressible, transversely isotropic viscous

fluid, with preferred direction defined by local fibre alignment

• Evolution equation for fibre direction: fibres advected with flow of

collagen

• Anisotropic drag between collagen fibres and other phases

(depends on fibre orientation)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Typical Simulation Results: Initial Conditions

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Typical Simulation Results (t > 0)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Typical Simulation Results: High Fibre Tension

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Summary

• Lecture 1: preliminaries
• Lecture 2: two-phase models
• Lecture 3: multiphase models
• Lecture 4: current work and open problems

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Open Problems in Multiphase Modelling

• Constitutive assumptions - mechanical properties of the

constituent phases

• Incorporating therapy
• Model comparison - hierarchies of models
• Relationship with multiscale and hybrid models (stochastic

effects): homogenisation

• Emergence and elimination of phases
• Numerical methods
• Model validation

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Background References

• YC Fung (1969). A First Course in Continuum Mechanics. Prentice-Hall
• YC Fung (1990). Biomechanics: Motion, Flow, Stress and Growth.

Springer.

• YC Fung (1993). Biomechanics: mechanical properties of living tissues.

Springer.

• AJM Spencer (1980). Continuum mechanics. Dover.
• O Gonzalez and AM Stuart (2008). A First Course in Continuum
• Mechanics. Cambridge University Press.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Classical Models of Tumour Growth

• HP Greenspan (1972). Models for the growth of a solid tumour

by diffusion. Stud. Appl. Math. 52 317-340.

• HP Greenspan (1976). On the growth and stability of cell

cultures and solid tumours. J. theor. Biol. 56 229-242.

• JP Ward and JR King (1997). Mathematical modelling of

avascular-tumour growth. IMA. J. Math. Appl. Med. 14 39-69.

• G Helmlinger, PA Netti, HC Lichtenbeld, RJ Melder and RK Jain

(1997). Solid stress inhibits the growth of multicellular tumour

• spheroids. Nature Biotech. 15 778-783.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Tumour Models

• HMB and MAJ Chaplain (1997). Free boundary value problems associated with the growth

and development of multicellular spheroids, Eur J Appl Math 8 639-658.

• HMB (1999). Weakly nonlinear analysis of a model of avascular solid tumour growth, J Math

Biol 39 59-89.

• CY Chen, HMBand JR King (2001). The influence of growth-induced stress from the

surrounding medium on the development of multicell spheroids, J Math Biol 43 191-220.

• CJW Breward, HMB and CE Lewis (2002). The role of cell-cell interactions in a two-phase of

solid tumour growth, J Math Biol 45 125-152.

• HMB and L Preziosi (2003). Modelling solid tumour growth using the theory of mixtures,

Math Med Biol 20 341-366.

• HMB, JR King, DLS McElwain and L Preziosi (2003). A two-phase model of solid tumor

growth, Appl Math Lett 16 567-573

• SD Webb, MR Owen, HMB, C Murdoch and CE Lewis (2007). Macrophage-Based

Anti-Cancer Therapy: Modelling Different Modes of Tumour Targeting. Bull Math Biol 69 1747-1776

• ME Hubbard and HMB (2013). Multiphase modelling of vascular tumour growth in two spatial
• dimensions. J Theor Biol. 316 70-89.
• TD Lewin, H Enderling, EG Moros, PK Maini and HMB (2018). The legacy of hypoxia on

tumour responses to fractionated radiotherapy. Bull Math Biol. (DOI: 10.1007/s11538-018-0391-9)

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Chemotaxis and Ductal Carcinoma in Situ

• HMB and MR Owen (2004). A new interpretation of the Keller-Segel

model based on multiphase modelling. J Math Biol 49 604-626.

• SJ Franks, HMB, JR King and CE Lewis (2003). Modelling the early

growth of ductal carcinoma in situ, J Math Biol 47 424-452.

• SJ Franks, HMB, HS Mudhar, JCE Underwood and CE Lewis (2003).

Mathematical modelling of comedo ductal carcinoma in situ of the breast, Math Med Biol 20 277-308.

• SJ Franks, HMB, JCE Underwood, CE Lewis (2005). Biological

inferences from a mathematical model of comedo ductal carcinoma in situ of the breast, J Theor Biol 232 523:543.

Lecture 1 Lecture 2 Lecture 3 Lecture 4

Tissue Engineering Models

• G Lemon, JR King, HM Byrne, OE Jensen and KE Shakesheff (2006). Mathematical

modelling of engineered tissue growth using a multiphase porous flow mixture theory. J Math Biol 52 571-594.

• RD ODea, SL Waters and HMB (2008). A two-fluid model for tissue growth within a dynamic

flow environment. Eur J Appl Math 19 607-634.

• JM Osborne, RD ODea, JP Whiteley, HMB and SL Waters (2010). The influence of

bioreactor geometry and the mechanical environment on engineered tissues. J Biomech Eng 132(5) 051006.

• RD ODea, JM Osborne, AJ El Haj, HMB, SL Waters (2013). The interplay between tissue

growth and scaffold degradation in engineered tissue constructs. J Math Biol 67 (5) 1199-1225

• LAC Chapman, RJ Shipley, JP Whiteley, MJ Ellis, HMB and SL Waters (2014). Optimising

cell aggregate expansion in a hollow fibre bioreactor via mathematical modelling. PLOS One 9 (8) e105813

• RJ Dyson, JEF Green, JP Whiteley and HMB (2015). An investigation of the influence of

extracellular matrix anisotropy and cell-matrix interactions on tissue architecture. J Math Biol 1-35.

• RD ODea, MR Nelson, AJ El-Haj, SL Waters and HMB (2015). A multiscale analysis of

nutrient transport and biological tissue growth in vitro. Math Med Biol 32 (3) 345-366.

• LG Bowden, HMB, PK Maini and DE Moulton (2015). A morphoelastic model for wound
• closure. Biomech Mod Mechanobiol 1-19.