Cell-scale degradation of peritumoural extracellular matrix fibre network and its role within tissue-scale cancer invasion

TL;DR

This paper develops a multiscale model to study how specific matrix metalloproteinases degrade the extracellular matrix at both tissue and cellular scales, elucidating their role in cancer invasion.

Contribution

It introduces a novel multiscale framework that explicitly models ECM fibre degradation at micro and macro scales during tumor progression.

Findings

MT1-MMP/MMP-2 cascade influences tumor invasion dynamics

ECM fibre degradation occurs at multiple scales

Model provides insights into enzyme-driven tissue remodeling

Abstract

Local cancer invasion of tissue is a complex, multiscale process which plays an essential role in tumour progression. Occurring over many different temporal and spatial scales, the first stage of invasion is the secretion of matrix degrading enzymes (MDEs) by the cancer cells that consequently degrade the surrounding extracellular matrix (ECM). This process is vital for creating space in which the cancer cells can progress and it is driven by the activities of specific matrix metalloproteinases (MMPs). In this paper, we consider the key role of two MMPs by developing further the novel two-part multiscale model introduced in [33] to better relate at micro-scale the two micro-scale activities that were considered there, namely, the micro-dynamics concerning the continuous rearrangement of the naturally oriented ECM fibres within the bulk of the tumour and MDEs proteolytic micro-dynamics that take place in an appropriate cell-scale neighbourhood of the tumour boundary. Focussing primarily on the activities of the membrane-tethered MT1-MMP and the soluble MMP-2 with the fibrous ECM phase, in this work we investigate the MT1-MMP/MMP-2 cascade and its overall effect on tumour progression. To that end, we will propose a new multiscale modelling framework by considering the degradation of the ECM fibres not only to take place at macro-scale in the bulk of the tumour but also explicitly in the micro-scale neighbourhood of the tumour interface as a consequence of the interactions with molecular fluxes of MDEs that exercise their spatial dynamics at the invasive edge of the tumour.

Figures

Tables

Table 1. Table 1 : The parameters in Σ 1 subscript Σ 1 \Sigma_{1}

Parameter	Value	Description	Reference
$D_1$	$3.5\times {10}^{-4}$	diffusion coeff. for cell population $c_1$	Domschke_et_al_2014
$D_2$	$7\times {10}^{-4}$	diffusion coeff. for cell population $c_2$	Domschke_et_al_2014
$D_m$	${10}^{-3}$	diffusion coeff. for MDEs	Estimated
${\mu }_1$	$0.25$	proliferation coeff. for cell population $c_1$	Domschke_et_al_2014
${\mu }_2$	$0.25$	proliferation coeff. for cell population $c_2$	Domschke_et_al_2014
${\gamma }_1$	2	non-fibrous ECM degradation coeff.	Shutt_2018
${\gamma }_2$	1.5	macroscopic fibre degradation coeff.	Peng2016
${\alpha }_1$	$1$	MDE secretion rate of $c_1$	Estimated
${\alpha }_2$	$1.5$	MDE secretion rate of $c_2$	Estimated
$R$	0.15	sensing radius	Shutt_2018
$r$	0.0016	width of micro-fibres	Shutt_2018
$f_{\text{max}}$	0.6360	max. micro-density of fibres	Shutt_2018
$p$	0.15-0.2	percentage of non-fibrous ECM	Estimated
$h$	0.03125	macro-scale spatial discretisation size	Dumitru_et_al_2013
$ϵ$	0.0625	size of micro-domain $ϵY$	Dumitru_et_al_2013
$\delta$	0.03125	size of micro-domain $\delta Y$	Shutt_2018

Equations

$$\delta Y(x)\ni z⟼z-x\in {\mathbb{R}}^N$$

$${\theta }_{{}_{f,\delta Y(x)}}(x,t)=\frac{\underset{\delta Y(x)}{\int }f(z,t)(z-x)dz}{\underset{\delta Y(x)}{\int }f(z,t)dz}.$$

$${\theta }_{{}_f}(x,t)=\frac{1}{\lambda (\delta Y(x))}{\int }_{\delta Y(x)}f(z,t)dz\cdot \frac{{\theta }_{{}_{f,\delta Y(x)}}(x,t)}{||{\theta }_{{}_{f,\delta Y(x)}}(x,t)||}.$$

$$F(x,t):=||{\theta }_{{}_f}(x,t)|{|}_2.$$

$$\text{u}(x,t)=(c(x,t),F(x,t),l(x,t))$$

$$\rho (\text{u}(x,t))={\vartheta }_v(F(x,t)+l(x,t))+{\vartheta }_{c}c(x,t),$$

$$\frac{\partial c}{\partial t}=\nabla \cdot [D_1\nabla c-c\mathcal{A}(x,t,\text{u}(\cdot ,t),{\theta }_f(\cdot ,t))]+{\mu }_{1}c(1-\rho (\text{u}))$$

$$\begin{split}\mathcal{A}(t,x,\textbf{u}(\cdot,t),\theta_{{}_{f}}(\cdot,t))=\frac{1}{R}\int_{B(0,R)}\mathcal{K}(\parallel\!y\!\parallel_{2})&\big{(}n(y)(\textbf{S}_{{}_{cc}}c(x+y,t)+\textbf{S}_{{}_{cl}}l(x+y,t))\\ &+\hat{n}(y)\ \textbf{S}_{{}_{cF}}F(x+y,t)\big{)}(1-\rho(\textbf{u}))^{+}\end{split}$$

$$n(y):=\left\{\begin{array}[]{l}y/||y||_{2}\quad\text{if}\ y\in\mathbf{B}(0,R)\setminus\{(0,0)\},\\[5.69046pt] (0,0)\quad\quad\text{if}\ y=(0,0).\end{array}\right.$$

$$\textbf{S}_{{}_{cc}}(x,t):=\textbf{S}_{{}_{max}}e^{\big{(}{1-\frac{1}{1-(1-l(x,t))^{2}}}\big{)}}.$$

$$\hat{n}(y):=\left\{\begin{array}[]{lll}\frac{y+\theta_{f}(x+y)}{||y+\theta_{f}(x+y)||_{{}_{2}}},&\text{if}&\ (y+\theta_{f}(x+y))\neq(0,0),\\[5.69046pt] (0,0)\in\mathbb{R}^{2},&&\text{otherwise}.\end{array}\right.$$

$$\frac{dF}{dt}$$

$$\frac{dl}{dt}$$

$$\mathcal{F}(x,t):=D_1\nabla c(x,t)-c(x,t)\mathcal{A}(x,t,\text{u}(\cdot ,t),{\theta }_f(\cdot ,t)).$$

$$r(\delta Y(x),t):=\omega (x,t)\mathcal{F}(x,t)+(1-\omega (x,t)){\theta }_f(x,t),$$

$$\omega (x,t)=\frac{c(x,t)}{c(x,t)+F(x,t)}.$$

$$z^{\ast }:=z+{\nu }_{{}_{\delta Y(x)}}(z,t)$$

$${\nu }_{{}_{\delta Y(x)}}(z,t)=\left(x_{\text{dir}}(z)+r(\delta Y(x),t)\right)\cdot \frac{f(z,t)(f_{\max }-f(z,t))}{f^{\ast }+||r(\delta Y(x))-x_{\text{dir}}(z)|{|}_2}\cdot {\chi }_{{}_{\{f(\cdot ,t)>0\}}}$$

$$\displaystyle\begin{split}1.\quad&g_{\epsilon Y}(y,\tau)=\frac{\int\limits_{\textbf{B}(y,\gamma)\cap\Omega(t_{0})}\alpha c(x,t_{0}+\tau)\ \widetilde{F}(x,t_{0}+\tau)\ dx}{\widetilde{F}\cdot\lambda(\textbf{B}(y,\gamma)\cap\Omega(t_{0}))},\quad y\in\epsilon Y\cap\Omega(t_{0}),\\[19.91684pt] 2.\quad&g_{\epsilon Y}(y,\tau)=0,\quad y\in\epsilon Y\setminus\big{(}\Omega(t_{0})+\{y\in Y|\ ||y||_{2}<\gamma\}),\end{split}$$

$$\widetilde{F}\cdot\lambda(G):=\int\limits_{G}\widetilde{F}(x,{\color[rgb]{0,0,0}t_{0}+\tau})\ dx,\quad G\in\Sigma(Y),$$

$$\frac{\partial m}{\partial\tau}=\underbrace{D_{m}{\color[rgb]{0,0,0}(f)}\Delta m}_{\text{diffusion}}\ \ +\underbrace{g_{\epsilon Y}(y,\tau)}_{\text{source term}},\quad y\in\epsilon Y,\quad\tau\in[0,\Delta t]$$

$$D_{m}{\color[rgb]{0,0,0}(f)}=\frac{D}{1+\alpha_{m}\tilde{f}(y,t_{0}+\tau)},$$

$$\begin{array}[]{cl}&\epsilon Y\ni y\longmapsto z(y)\in\bigcup\limits_{x\in Y}\delta Y(x)\\ \textrm{is given by:}\\ z(y)\quad\textrm{ is the set of}&\textrm{the micro-spatial points within up to four $\delta Y(x)(y)$}\\ &\textrm{where}\\ &\textrm{each selected $\delta Y(x)(y)$ is a micro-cube within $\bigcup\limits_{x\in Y}\delta Y(x)$}\\ &\textrm{with the property that this contains $y\in\epsilon Y$}\\ &\textrm{(please see Figure \ref{fig:shapefunction} for an illustration of this situation)}.\end{array}$$

$$\tilde{f}(y,t_{0}+\tau):=\left\{\begin{array}[]{lll}f(z(y),t_{0}+\tau)&if&card(z(y))=1;\\ \frac{1}{card(z(y))}\sum\limits_{\zeta\in z(y)}f(\zeta,t_{0}+\tau)&if&card(z(y))>1.\end{array}\right.$$

$$\begin{array}[]{l}\Phi_{fm}(\!z,\tau\!)\!=\!\!\left\{\!\!\begin{array}[]{lll}\arccos\!\!\left(\!\frac{<\nabla m(z,\tau\!),\nabla f(z,t_{0}+\tau\!)>}{||\nabla m(z,\tau\!)||_{{}_{2}}||\nabla\!f(z,t_{0}+\tau\!)||_{{}_{2}}}\!\right)\chi_{{}_{\delta Y\!(x)\!\cap\!\epsilon Y}}(z),&if&\parallel\!\!\nabla\!f(z,t_{0}\!+\!\tau\!)\!\!\parallel_{{}_{2}}>\!0,\\[5.69046pt] \frac{\pi}{2}\chi_{{}_{\delta Y\!(x)\!\cap\!\epsilon Y}}(z),&if&\parallel\!\!\nabla f(z,t_{0}\!+\!\tau\!)\!\!\parallel_{{}_{2}}=\!0.\end{array}\right.\end{array}$$

$$D_g(z,\tau ):=d_f\cdot \exp \left(1-\frac{1}{1-(1-\Phi (z,\tau ){)}^2}\right).$$

$$\Phi (z,\tau )=\left|\frac{{\Phi }_{fm}(z,\tau )}{{\Phi }_{\text{max}}}\right|,$$

$$\frac{\partial f}{\partial\tau}=-D_{g}(z,\tau)\cdot f(z,{\color[rgb]{0,0,0}t_{0}+}\tau),\qquad z\in\delta Y,\ \tau\in[0,\Delta t].$$

$$\displaystyle\begin{split}\text{{\color[rgb]{0,0,0}Macro-scale}}&\,\,\text{{\color[rgb]{0,0,0} dynamics:}}\\ \frac{\partial c}{\partial t}&=\nabla\cdot[D_{1}\nabla c-c\mathcal{A}(t,x,\textbf{u}(\cdot,t),\theta_{{}_{f}}(\cdot,t))]+\mu_{1}c(1-\rho(\textbf{u})),\\[5.69054pt] \frac{dF}{dt}&=-\gamma_{1}cF,\\[5.69054pt] \frac{dl}{dt}&=-\gamma_{2}cl+\omega(1-\rho(\textbf{u})),\\ \end{split}$$

$$\displaystyle\begin{split}\text{{\color[rgb]{0,0,0}Micro-scale}}&\,\,\text{{\color[rgb]{0,0,0} boundary MMP-2-dynamics interacting with peritumoural micro-fibres:}}\\ \frac{\partial m}{\partial\tau}&=D_{m}{\color[rgb]{0,0,0}(f)}\Delta m+g_{\epsilon Y}(y,\tau),\quad y\in\epsilon Y,\ \tau\in[0,\Delta t],\\[5.69054pt] \frac{\partial f}{\partial\tau}&=-D_{g}\cdot f(z,{\color[rgb]{0,0,0}t_{0}+}\tau),\quad z\in\delta Y,\ \tau\in[0,\Delta t].\\ \end{split}$$

$$\begin{array}{rl}\nabla \cdot [\nabla m{]}_{i,j}^n& =\text{div}[\nabla m{]}_{i,j}^n\\ & \simeq \frac{[m_x{]}_{i+\frac{1}{2},j}^n-[m_x{]}_{i-\frac{1}{2},j}^n}{\Delta x_{ϵ}}+\frac{[m_y{]}_{i,j+\frac{1}{2}}^n-[m_y{]}_{i,j-\frac{1}{2}}^n}{\Delta y_{ϵ}}\end{array}$$

Full text

Focussing primarily on the activities of the membrane-tethered MT1-MMP and the soluble MMP-2 with the fibrous ECM phase, in this work we investigate the MT1-MMP/MMP-2 cascade and its overall effect on tumour progression. To that end, we will propose a new multiscale modelling framework by considering the degradation of the ECM fibres not only to take place at macroscale in the bulk of the tumour but also explicitly in the micro-scale neighbourhood of the tumour interface as a consequence of the interactions with molecular fluxes of MDEs that exercise their spatial dynamics at the invasive edge of the tumour.

††journal: Bulletin of Mathematical Biology

1 Introduction

Cancer cell invasion of tissue is a complex, multiscale process in which gene mutations in healthy cells promote enhanced proliferation and the production of proteolytic enzymes. One of the first steps of the local invasion of tissue is the secretion of matrix degrading enzymes (MDEs) and the consequential degradation of the extracellular matrix (ECM). The ECM is a non-cellular structure that provides not only support to surrounding cells and tissues, but also acts as a platform for cellular communication. This feature is of particular use to cancer cells, which take advantage of the molecular interactions mediated by various ECM components and favourably utilise these as means for achieving invasion of the surrounding tissue.

The ECM is comprised of a variety of secreted proteins, and the main constituent of the ECM is the structural cross-linked collagen type I, a dense network of fibres that give the ECM its rigidity. In order for a tumour to progress, these strong fibres must be broken down and degraded to free space for the cancer cells. One of the first MDEs to interact with the ECM is the membrane tethered MMP, MT1-MMP, or MMP-14. This MMP exhibits strong collagenolytic capabilities in which they are able to cleave the cross-linked collagen type I fibres and break them into shorter, soluble fibres. These soluble fibres are then degraded by the freely diffusible MMP-2, activated through the cleavage of proMMP-2 molecules by MT1-MMP. This MT1-MMP/MMP-2 cascade is highly effective in promoting tumour invasion through ECM degradation, however, there are some disadvantages to both types of MMPs. MT1-MMP molecules can overcome high density regions of collagen type I, particularly the cross-linked fibres, however they do not degrade the collagen, rather cleaving the fibres into smaller fibrils Tam_2004 . On the other hand, MMP-2 cannot degrade the dense cross-linked fibres, but can degrade the smaller, soluble fibrils within the peritumoural region Doren_2015 . Consequently, these two MMPs work in harmony with one another for successful invasion of tissue.

Besides the abnormal rate of MDE secretion, the invasive capabilities of a tumour are strengthened by many other processes, including increased proliferation rates and the ability to adapt cellular adhesion properties. Cell adhesion is an essential process in which cells interact and attach to neighbouring cells through calcium dependent cell-adhesion molecules, known as CAMs, on the cell surface Humphries2006 . Cell-cell adhesion is dependent on specific cell-signalling pathways that are formed between $\text{Ca}^{2+}$ ions and calcium sensing receptors in the ECM Ko_2001 . This calcium dependent cell-cell signalling is regulated by a subfamily of glycoproteins known as E-cadherins that bind with intra-cellular proteins known as catenins, typically $\beta$-catenin, forming the E-cadherin/catenin complex. Any alteration to the function of $\beta$-catenin will result in a loss of the ability of E-cadherin to initiate cell-cell adhesion Wijnhoven2000 . This connection between E-cadherin and $\beta$-catenin and the calcium cell-signalling mechanism was first recorded in colon carcinoma Bhagavathula2007 . Additionally, cells can also bind to the ECM through cell-matrix adhesion Lodish2000 . Mediated by calcium independent CAMs, known as integrins, cell-matrix adhesion enables binding of cells to various components of the ECM, such as collagen and fibronectin, with this contributing to cancer cell migration within the tissue. Moreover, enhanced cell-matrix adhesion paired with a loss in cell-cell adhesion facilitates a quicker spread of the cancer cells into the surrounding tissue CAVALLARO200139 .

Over the past $25$ years there has been an increasing interest in the mathematical modelling of cancer invasion, see for example Andasari2011 ; Anderson2005 ; Anderson2000 ; Chaplain2011 ; Chaplain2006a ; Gerisch2008 ; Peng2016 ; Ramis-Conde_et_al_2000 ; szymanska_08 ; Dumitru_et_al_2013 . Many models of cancer invasion have focussed on the interactions between cancer cells and the extracellular matrix through a variety of different approaches, ranging from continuous and discrete individual cell-based models to the more complex multiscale modelling approach. Numerous processes of cancer invasion have been addressed, including the effects of proliferation and cellular adhesion Bitsouni_2017 ; Chauviere_2007 ; Domschke_et_al_2014 ; Gerisch2008 ; Painter2008 ; Ramis-Conde_et_al_2000 , as well as models developed to describe the migration strategies of cancer cells in tissue networks Chauviere_2007 ; Hillen2006 ; Perumpanani_et_al_1998 . A model describing both the mesenchymal and amoeboid motion of cells through a fibre network Painter2008 concluded that a structured ECM can induce cell aggregation in amoeboid type cells, whereas the actions of both contact-guidance and ECM remodelling are sufficient processes for mesenchymal type cell invasion to occur.

Numerous processes of cancer invasion have been investigated, including the effects of proliferation and cellular adhesion Bitsouni_2017 ; Chauviere_2007 ; Domschke_et_al_2014 ; Gerisch2008 ; Painter2008 ; Ramis-Conde_et_al_2000 as well as the secretion and interaction of proteolytic enzymes, specifically MMPs and uPAs with the ECM Andasari2011 ; Chaplain2005 ; Chaplain2006 ; Peng2016 ; shutt_chapter ; Dumitru_et_al_2013 during tumour invasion. In particular, a single-scale model developed by Deakin_Chaplain_2013 focussed on the roles of two MMP molecules, namely the membrane-bound MT1-MMP and the soluble MMP-2, where the MT1-MMP/MMP-2 cascade was considered, highlighting the importance of MT1-MMP matrix remodelling within collagen-rich environments. These model dynamics are of great interest to us for the current investigation of the peritumoural MMP processes.

Finally, as the invasion process is naturally multiscale, with its dynamics ranging from sub-cellular-, cellular- to tissue-scale, major advances have been witnessed in the multiscale modelling of cancer invasion anderson_07 ; Peng2016 ; Ramis-Conde_et_al_2000 ; shutt_chapter ; Dumitru_et_al_2013 . In particular, advancements have been made in towards two-scale approaches, modelling and appropriately linking the spatio-temporal dynamics occurring at different scales, as first proposed in Dumitru_et_al_2013 and extended upon in Peng2016 ; Shutt_twopop ; Shutt_2018 .

In this paper, we aim to advance the novel two-part multiscale model developed in Shutt_2018 by considering the multiscale contribution of ECM fibres within the invasive edge proteolytic dynamics of MDEs. This will pave the way for the establishment of an additional cell-scale link between the peritumoural micro-scale fibres rearrangements and the micro-scale proteolytic MDEs activities at the tumour boundary, which ultimately extends and complements the two-parts multiscale framework introduced in Shutt_2018 where the linking of its two constituent multiscale systems was so far mediated only through the shared macro-dynamics.

2 The mathematical modelling approach

As this work builds upon the two-part multiscale model introduced and developed in Shutt_2018 that investigates cancer invasion within a heterogeneous ECM, let us start by revisiting the multi-component structure of the ECM and macroscopic tumour dynamics for one cancer cell population considered in Shutt_2018 , while adapting the model to include further integration between the leading edge MDE and fibres rearrangement micro-scale processes.

Since this paper builds on the modelling developed in Shutt_2018 and Dumitru_et_al_2013 , let us start this section by recalling the key details of the framework terminology that we introduced there. Let us denote the support of the invading tumour region as $\Omega(t_{0})$, and assume this evolves in the maximal reference tissue cube $Y\in\mathbb{R}^{N}$ for $N=2,3$, centred at the origin of the space. At any spatio-temporal node $(x,t)\in\Omega(t_{0})\times[0,T]$, we consider the tumour to be comprised of a cancer cell population $c(x,t)$, integrated within a multiphase heterogeneous ECM density denoted by $v(x,t)$. Specifically, the heterogeneous ECM is regarded as comprising of a fibres component and non-fibres soluble component. We denote the tissue-scale (macro-scale) density of non-fibres ECM phase by $l(x,t)$ and we consider to include all non-fibre components of the ECM, i.e., elastin, laminins, fibroblasts, etc. On the other hand, the macro-scale mass density of the fibres ECM phase is denoted by $F(x,t)$ and accounts for all fibrous proteins such as collagen and fibronectin within the matrix.

For completion, in the following the following subsections, we will briefly revisit the key points about the ECM fibres multiscale structure and functionality alongside the tumour macro-dynamics as well as the associated terminology and notations that were introduced and derived in Shutt_2018 .

2.1 The multiscale ECM fibre structure and its contribution to the tissue dynamics

As derived in Shutt_2018 , at any macro-position $x\in\Omega(t_{0})$, the ECM-fibre phase can be represented through a macro-scale vector field, denoted $\theta_{f}(x,t)$, that captures not only the macroscopic distribution of fibres $F(x,t)$, but also their naturally arising macroscopic orientation that is induced by their mass distribution of micro-fibres, denoted $f(z,t)$, over a given micro-domain $\delta Y(x)$ centred at the given macroscopic point $x$ of cell-scale $\delta>0$ (i.e., $\delta Y(x):=x+\delta Y$). An example of micro-fibres $f(z,t)$ patterns over a micro-domain $\delta Y(x)$ (defined in Appendix C), alongside the naturally emerging macroscopic fibre orientation $\theta_{f}(x,t)$ that was derived in Shutt_2018 , are shown schematically in Figure 1.

In brief, while referring to its full derivation presented in Shutt_2018 , the naturally generated revolving barycentral orientation $\theta_{{}_{f,\delta Y(x)}}(x,t)$ associated with $\delta Y(x)$ is given by the Bochner-mean-value of the position vectors function

[TABLE]

with respect to the density measure $f(z,t)\lambda(\cdot)$, where $\lambda(\cdot)$ is the usual Lebesgue measure (see yosida1980 ), and so this is expressed mathematically as:

[TABLE]

Following on, at any spatio-temporal node $(x,t)$, this revolving barycentral orientation $\theta_{{}_{f,\delta Y(x)}}(x,t)$ induces the naturally arising macroscopic fibre orientation vector field representation that is defined as

[TABLE]

in which macroscopic mean-value fibre representation at any $(x,t)$ is given by the Euclidean magnitude of $\theta_{{}_{f}}(x,t)$, namely,

[TABLE]

Thus, the total ECM distributed at any spatio-temporal node $(x,t)$ is therefore given by $v(x,t)=l(x,t)+F(x,t)$.

2.2 Macro-scale dynamics

For the tissue dynamics of a tumour consisting of a single cell population, again respecting the terminology and notations introduced in Shutt_2018 , we denote the macro-scale tumour population, including the surrounding two-phase ECM, by the vector

[TABLE]

with the tumour’s volume of occupied space given as

[TABLE]

where $\vartheta_{v}$ and $\vartheta_{c}$ represent the fractions of physical space occupied by the entire ECM and the cancer cells respectively. Focussing first on the cancer cell population, the spatial movement of the cancer cells is governed by random motility (approximated here by diffusion) and a cell-adhesion process that includes both cell-cell and cell-matrix adhesion, with cell-matrix adhesion accounting for both cell-fibre and cell-non-fibre adhesion. Assuming the cells are subject to a logistic proliferation law, the dynamics of the cancer cell population can be mathematically represented as

[TABLE]

where: $D_{1}$ and $\mu_{1}$ are the non-negative diffusion and proliferation rates, respectively. The non-local flux $\mathcal{A}(x,t,\textbf{u}(\cdot,t),\theta_{f}(\cdot,t))$ describes the spatial movement of the cancer cells due to cellular adhesion properties between each other (cell-cell adhesion) and the surrounding environment (cell-matrix adhesion accounting for both cell-fibre and cell-non-fibre adhesion). The adhesive flux considers the interactions of cancer cells within a sensing radius $R$ with the other cancer cells and two-phase ECM distributed on the sensing region $\textbf{B}(x,R):=x+\{\xi\in\mathbb{R}^{N}\,|\,\parallel\xi\parallel_{{}_{2}}<R\}$, described by the following adhesive flux term

[TABLE]

Whilst the influence of the spatially distributed adhesive interactions on the sensing region $\textbf{B}(x,R)$ is accounted through a radial kernel $\mathcal{K}(\cdot)$ detailed in Appendix A, we explore the strength of the bonds created between the cells distributed at $x$ and the other cells or non-fibre ECM phase distributed at $y$ in the direction of the unit normal

[TABLE]

through the cell-cell and cell-non-fibre ECM adhesion strength coefficients $\textbf{S}_{cc}$ and $\textbf{S}_{cl}$, respectively. We consider $\textbf{S}_{cl}$ to be constant, whilst the coefficient representing cell-cell adhesion, $\textbf{S}_{cc}$, is monotonically dependent on the levels of extracellular $\text{Ca}^{2+}$ ions which enable robust adhesive bonds between neighbouring cells Gu2014 ; Hofer2000 . Therefore we assume $\textbf{S}_{cc}$ is dependent on the underlying non-fibre ECM phase, smoothly ranging from [math] to a $\text{Ca}^{2+}$-saturation level denoted $\textbf{S}_{{}_{max}}$ and is taken as

[TABLE]

The last term in (5) considers the interactions between the cancer cells and the fibre ECM phase distributed within the region $\textbf{B}(x,R)$. The strength of this interaction is proportional to the distribution of macro-fibres at $F(x+y,t)$ and it is the orientation of these fibres, $\theta_{f}(x+y,t)$, that biases the direction in which adhesive bonds are formed in the direction $\hat{n}(\cdot)$

[TABLE]

Finally, we describe the macroscopic dynamics of the ECM involving both the fibre and non-fibres phase on the invading tumour domain $\Omega(t)$, where the fibre and non-fibre ECM phase are subject to an overall degradation by the cancer cells and can be mathematically represented as

[TABLE]

where $\gamma_{1},\gamma_{2}$ are the degradation rates for fibres and non-fibres ECM phases respectively.

2.3 Microscopic fibre rearrangement instigated by the macro-scale cell-flux

As derived in Shutt_2018 , during the macro-dynamics, as the cancer cells invade the surrounding tissue they have the ability to push the fibres in the direction of their spatial flux and rearrange the micro-fibre mass distributions, thereby reorienting the macro-fibres direction. Thus, at time $t$ and at any spatial location $x\in\Omega(t_{0})$, the micro-fibres $f(z,t),\ \forall z\in\delta Y(x)$ undergo a microscopic rearrangement process fuelled by the cancer cell-flux direction given mathematically as

[TABLE]

Hence, the micro-fibres $f(z,t),\ \forall z\in\delta Y(x)$ are acted upon uniformly by the resultant force of the rearrangement flux vector-valued function

[TABLE]

where the spatial flux of the cancer cells $\mathcal{F}(x,t)$ is balanced in a weighted manner by the orientation $\theta_{f}(x,t)$ of the distribution of fibres at $(x,t)$ and the amount of cells exercising spatial transport at $(x,t)$, with a naturally emerging weight given by

[TABLE]

As detailed in Shutt_2018 , under the influence of the force $r(\delta Y(x),t)$, an appropriate level of micro-fibres mass $f(z,t)$ will undergo a spatial movement towards a new position

[TABLE]

where the relocation direction and magnitude is given by

[TABLE]

with $f_{\max}$ representing the maximum amount of micro-fibres at $z\in\delta Y(x)$, and the relative degree of fibre’s occupancy at $z$ being denoted by $f^{*}(z,t):=\frac{f(z,t)}{f_{\text{max}}}$. Finally, the spatial transport of micro-fibres at $z$ will be exercised in accordance to the physical space available at the new position $z^{*}$, which is explored through a movement probability $p_{move}:=\max(0,1-f^{*}(z^{*},t))$, where an amount of micro-fibres $p_{move}f(z,t)$ will be moved to the new position $z^{*}$ while the rest of the micro-fibres will remain at $z$.

3 Novel modelling for the MDE proteolytic micro-scale dynamics at the tumour interface

To incorporate the dynamics of peritumoural micro-scale fibres with the MDE micro-dynamics at the tumour invasive edge, let us start by exploring the influence of the macro-fibre distributions on the emergence of a cell-scale molecular source of MDE at the tumour interface. This will effectively enhance the top-down macro-micro link derived and introduced in Dumitru_et_al_2013 (and which was considered in Shutt_2018 ) for the boundary micro-dynamics and in return will influence the bottom-up feedback link to macro-dynamics. Finally, MDEs micro-dynamics occurring within a cell-scale neighbourhood of the tumour boundary will be explored on an appropriately selected covering bundle of boundary micro-domains $\{\epsilon Y\}_{\epsilon Y\in\mathcal{P}}$, which was introduced and constructed with complete details in Dumitru_et_al_2013 .

The top-down link

As discussed previously, the expansion of the tumour boundary is dependent on the peritumoural degradation of ECM by the matrix-degrading enzymes (MDEs). However, the secretion of MDEs induced by the distribution of cancer cells in the outer proliferating rim is dependent upon the structure of the ECM and in particular upon the fibre distribution. Indeed, in the presence of a high distribution of ECM fibres, the cancer cells exhibit a strong rate of MT1-MMP secretion, which in turn leads to an increase in the activation of proMMP-2 molecules zigrino_2001 . Thus, since in this emerging MT1-MMP/MMP-2 cascade, the secretion rate of MT1-MMP correlates directly to the amount of MMP-2 molecules, we obtain therefore that the presence of the ECM fibres enhances the production of the MMP-2 molecules that are then released at the tumour boundary.

Thus, the MMPs source at each $y\in\epsilon Y\cap\Omega(t_{0})$ arises as a collective contribution of the cancer cells within the outer proliferating rim which are further enhanced by the presence of ECM fibres, and therefore this can be mathematically expressed as

[TABLE]

where $\gamma$ represents the maximum thickness of the outer proliferating rim, $\mathbf{B}(y,\gamma):=\{\xi\in\mathbb{R}^{N}\,|\,\parallel y-\xi\parallel_{{}_{\infty}}\leq\gamma\}$ and $\alpha$ is an MMP secretion rate for the cancer cell population. Furthermore, the function $\widetilde{F}(x,t+\tau):=1+F(x,t)$, $\forall\,t>0$ explores the spatially distributed enhancement of the source of MMPs produced by the cancer cells that is enabled through the presence of fibres zigrino_2001 (i.e., a higher distribution of fibres inducing a greater number of MMP-2 molecules). Finally, $\widetilde{F}\cdot\lambda$ represents the underlying fibre density measure that is defined by

[TABLE]

where $\lambda$ is the standard Lebesgue measure on $\mathbb{R}^{2}$, and $G$ is a nonempty Borel subset set of $Y$, i.e., $G\in\Sigma(Y)$, with $\Sigma(Y)$ representing the Borel $\sigma$-algebra on $Y$.

In the presence of this source, on any micro-domain $\epsilon Y$, a cross-interface MMP-2 diffusive transport takes place, and as the MMP-2 find it easier to exercise their random movement in regions of lower micro-fibres density, the diffusion rate is micro-fibre density dependent, and therefore this spatio-temporal micro-dynamics can be mathematically formulated as

[TABLE]

with the fibre-dependent diffusion coefficient

[TABLE]

where $D$ is the baseline diffusion rate and $\alpha_{m}$ being a “slowing down” constant factor induced by the presence of the micro-fibres density $\tilde{f}(y,t_{0}+\tau)$, which is defined as follows. Considering first the micro-scale set-valued mapping that surveys the number of points $z$ corresponding to up to four fibres micro-domains $\delta Y(x)$ that represent the same micro-position $y$, namely

[TABLE]

we have that

[TABLE]

The bottom-up link

During their micro-dynamics, the MMPs diffuse into the surrounding ECM and it is the pattern of their advancing spatial distribution that control the degradation of the peritumoural ECM captured within each micro-domain $\epsilon Y$. This degradation ultimately leads to the movement of the tumour boundary, whereby a movement direction $\eta_{\epsilon Y}$ and displacement magnitude $\xi_{\epsilon Y}$ are derived from the pattern of ECM degradation in each micro-cube $\epsilon Y$ (and for full derivation we refer the reader to Dumitru_et_al_2013 ). The microscopic movement of the boundary is represented at the macro-scale through the movement from the boundary midpoint $x^{*}_{\epsilon Y}$ to a new spatial position $\widetilde{x^{*}_{\epsilon Y}}$. Thus, although the bottom-up link of the model is much akin to previous works shutt_chapter ; Shutt_2018 , the specific context in which the presence of micro-fibres distribution causes on one enables on one hand an enhanced source of MMP and on the other hand acts as an impediment for their random motility, the MMP micro-dynamics within the peritumoural region of the tumour domain $\Omega(t_{0})$ incorporates now these important aspects, ultimately resulting in an improved estimate for the macroscopic boundary movement characteristics (i.e., encapsulated by the movement direction $\eta_{\epsilon Y}$ and displacement magnitude $\xi_{\epsilon Y}$ ). This crucial micro-scale-induced boundary relocation is then translated back to the macro-scale, resulting in an expanded tumour domain $\Omega(t_{0}+\Delta t)$ on which the multiscale dynamics is continued (as detailed in Dumitru_et_al_2013 ).

3.1 Microscopic fibre degradation

Given the proteolytic properties of the MMPs, the cell-scale cross-interface MMP-2 micro-dynamics will result in a direct peritumoural ECM degradation, whose pattern not only depends on the amount of MMP-2 transported at a given location $y$ within a given micro-domain $\epsilon Y$, but also on the existing ECM micro-fibres spatial distribution at $y$. To address this micro-scales interaction that takes place between the MMP-2 cross-interface transport and micro-fibres within any given micro-domain $\epsilon Y$, we consider that the spatial patterns of micro-fibres distributions that are aligned with the MMP-2 flux $\nabla m(\cdot,t_{0}+\tau)$ suffer less degradation than those that are positioned orthogonal to it. Thus, for any fibre micro-domain $\delta Y(x)$ that has non-empty intersection with $\epsilon Y$, by denoting $\Phi_{fm}(\cdot,\cdot):\delta Y(x)\times[0,\Delta t]\to[0,\pi]$ the function that explores these emerging angles, given by

[TABLE]

where $\chi_{{}_{\delta Y\!(x)\!\cap\!\epsilon Y}}(\cdot)$ is the usual characteristic function of $\delta Y\!(x)\!\cap\!\epsilon Y$.

The strength of the micro-fibre degradation rate is influenced by their degree of alignment with the flux of the MMPs, which is explored by the angle $\Phi_{fm}(\!z,\tau\!)$ at which the flux of the MMP-2 molecules acts upon the mass distribution of micro-fibres $f(z,t_{0}+\tau)$ at each micro-position $z$, and so this can therefore be mathematically formalised as

[TABLE]

where $d_{f}$ is a non-negative degradation constant, and the function $\Phi(\cdot,\cdot)$

[TABLE]

represents the influence of the collision angle relative to the maximum angle $\Phi_{\text{max}}\!:=\!\!\frac{\pi}{2}$ where the highest rate of degradation will occur (i.e., when the flux of MMPs is perpendicular upon the micro-fibres mass). Thus, the micro-fibre dynamics at location $z$ in the $\delta$-sized micro-scale can be mathematically represented as

[TABLE]

4 Summary of model

The tissue- and cell–scale dynamics of the new multiscale moving-boundary model proposed here, can therefore be summarised as:

[TABLE]

The macro- and micro-dynamics summarised in (18) aggregates two interconnected multiscale frameworks that describe the interactions between two independent micro-scale systems and share the same macro-scale tissue level dynamics (summarised in (18a)) that is linked to the two micro-dynamics through two double feedback loops, as established and detailed in Shutt_2018 .

While the macro-scale tissue dynamics describe the evolution of the spatial distribution of cancer cells and both the non-fibres and macro-fibres ECM phase, the micro-scale part of the first multiscale system governs the dynamic rearrangement of fibres. Specifically, the rearrangement of the fibres emerges as a consequence of the acting macro-scale flux of the cancer cell spatial flux upon the mass distribution of micro-fibres distributed within micro-domains $\delta Y(x)$ which are centred at any given macroscopic point $x\in Y$. The redistributed micro-fibres naturally yield a new macroscopic vector field representation of the newly oriented fibres that will lave a cascade influence upon adhesion processes the cancer cell population tissue scale dynamics biasing their migration, as detailed in Shutt_2018 .

The second multiscale system involved in this modelling brings into the picture contribution of the micro-scale proteolytic dynamics occurring in a cell-scale neighbourhood of the tumour interface. Here we take forward the multiscale moving boundary approach initially introduced in Dumitru_et_al_2013 and recapitulated in Shutt_2018 , by considering the influence that the oriented macro-scale ECM fibres has upon the emergence of the micro-scale MMP-2 source and at the same time exploring the explicit cell-sale interaction that takes place between the MMP-2 spatial flux and the mass distribution of micro-fibres at micro-scale which results in fibre degradation, summarised in (18b), notably triggering changes in the macroscopic fibres orientation. Consequently, these leading edge MMPs dynamics instigate a change in the position of the tissue-scale tumour boundary that corresponds to the pattern of peritumoural ECM degradation, this way allowing the macro-dynamics to continue on the newly enlarged tumour region and thus the invasion process continues.

5 Numerical approaches and initial conditions for computations

Expanding on the multiscale moving boundary framework developed in Shutt_2018 building on the model initially introduced in Dumitru_et_al_2013 , we developed a new modelling and computational approach to address specifically the cell-scale peritumoural interaction between the MDE and the micro-fibres mass distributions at micro-scale. Specifically, we explore the link between fibre distribution and MDE density at the tumour interface, in addition to the MDE induced micro-fibre degradation at the cell-scale.

5.1 Brief description of the numerical approach

To address the tumour macro-dynamics, we use the novel predictor-corrector method developed and fully defined in Shutt_2018 , that accounts for the complexity of the cancer dynamics. For this, we consider a uniform spatial mesh of size $h=0.03125$, and we use a combination of central differences and mid-point methods to discretise the local spatial operators, while involving an off-grid approach (introduced and detailed in Shutt_2018 ) for the calculation of the non-local adhesion terms (that enable the adhesive flux) at each spatio-temporal node.

Furthermore, to obtain the microscopic boundary relocation described in Shutt_chap, we explore the top-down and bottom-up link, using a finite difference approach for computing the MDE micro-dynamics occurring on the bundle of boundary micro-domains $\{\epsilon Y\}_{\epsilon Y\in\mathcal{P}}$ over the time interval $[\tau_{0}+\Delta t]$. Hence, while considering a time discretisation into $p$ uniformly distributed time steps, i.e., $\delta t=\frac{\Delta t}{p}$, each micro-cube $\epsilon Y$ is also discretised uniformly, using a spatial mesh of size $h_{\epsilon}$, i.e., $\Delta x_{\epsilon}=\Delta y_{\epsilon}=h_{\epsilon}$. Therefore, to discretise the reaction-diffusion equation (12), we start by addressing the spatially discretised source term induced by the cancer cell and fibre distribution (11) in the way it was described in Dumitru_et_al_2013 . Then, to solve the spatio-temporal dynamics in (12), we develop a similar predictor-corrector method using trapezoidal corrector for the time marching, while for the spatial operators involved we use again a combination of mid-points and central differences, which in this context are given by

[TABLE]

where

[TABLE]

with $n=0,\ldots,p$ and $i=1,\ldots,q$, $j=1,\ldots,q$ representing the indices for the $y$- and $x$-directions, respectively. Finally, important for capturing the degradation of the mass distribution of peritumoural micro-fibres at micro-scale, we use bilinear shape functions on square elements hughes_book to appropriately interpolate the solutions and the associated fluxes of the MMP-2 micro-dynamics in the eventually overlapping regions of their $\epsilon Y$s boundary micro-domains. This way the appropriate MMP-2 flux information is obtained at any micro-position $z$ in any given intersecting $\delta Y(x)$ fibres micro-domains, i.e., within those fibres micro-domains $\delta Y(x)$ in the peritumoural region for which

[TABLE]

as illustrated in Figure 2.

All of the following simulations of the model have been performed on MATLAB.

5.2 Initial conditions used in computational simulations

We consider the same initial conditions from Shutt_2018 , such that we consider the initial cancer cell population $c(x,0)$ to occupy the region $\Omega(0)=\textbf{B}((2,2),0.5)$ positioned in the centre of the tissue cube $Y$, Figure 3,

[TABLE]

where $\psi$ is the standard mollifier detailed in Appendix B that acts within a radius $\gamma<<\frac{\Delta x}{3}$ from $\partial\textbf{B}((2,2),0.5-\gamma)$ to smooth out the characteristic function $\chi_{{}_{\mathbf{B}((2,2),0.5-\gamma)}}$.

We explore the invasion of cancer within two different compositions of the non-fibres ECM phase, a homogeneous distribution given as

[TABLE]

and a heterogeneous distribution, previously explored in Shutt_2018 and defined below

[TABLE]

where

[TABLE]

These initial conditions are shown in Figure 4.

Finally, the fibres ECM phase will be initialised with both a homogeneous and heterogeneous macro-scale distribution. To this end, we first assume a random distribution of five pre-assigned micro-distributions containing patterns of five different micro-fibres (detailed in Appendix C) that are then randomly assigned onto $\delta Y(x):=x+\delta Y$. To calibrate the macro-fibre distributions, we use the initial conditions for the non-fibres ECM phase, for either a homogeneous or heterogeneous distribution, (21) and (22), respectively, taking a percentage $p$ of the density of the non-fibres ECM at each spatio-temporal position $x$. This allows for the control of the maximal height of the micro-fibres in each $\delta Y(x)$, centred at macro-position $x$, so that the resulting macro-fibre distribution $F(x,0)$ represents the percentage $p$, which here will be $15\%,\ p=0.15$ or $20\%,\ p=0.2$, with the latter representing a denser collagen structure. All of the following simulations of the model have been performed on MATLAB.

6 Results

We first explore tumour invasion in the presence of a homogeneous ECM, where the non-fibres phase is the homogeneous distribution (21) with the fibres phase taking the percentage $p=0.15$ of $l(x,0)$. Using the parameter set $\Sigma_{1}$ from Appendix D, the microscopic fibre degradation rate $d_{f}=1$, and the cell adhesion coefficients

[TABLE]

we investigate the effects micro-fibres degradation on tumour progression. Figures 5, 6 and 7 display the evolution of the tumour at stages $\Delta 25t$, $\Delta 50t$ and $\Delta 75t$, respectively. Each figure contains the subfigures: (a) cancer cell population; (b) non-fibres ECM density; (c) macroscopic fibre distribution; vector fields of fibre orientations at two different resolutions, namely: (d) coarsened two-fold; and (f) coarsened four fold; and finally in (e) a $3$D plot of the orientation of the complete ECM distribution $v(x,t)=F(x,t)+l(x,t)$.

After $25\Delta t$ macro-stages, Figure 5, the boundary of the tumour has remained largely unchanged, subfigure 5(a). Both the non-fibres and fibres ECM phase undergo degradation where the cancer cell distribution is highest, subfigures 5(b), 5(c), whilst minor degradation of the fibres has also occurred at the tumour boundary. The masses of micro-fibres have been rearranged such that the macroscopic orientation of fibres remains in line with the general fibre direction, subfigures 5(d), 5(f) which points towards the origin of the space. However, there are some irregularities, particularly visible in subfigure 5(f) which has been magnified 4-fold, where the orientation of the fibres is mixed with some fibres realigned to near perpendicular of their original orientation. The cancer cells are pushing and rearranging the fibres in a direction opposite that of the initial fibre orientation, subfigure 5(c), where a build up of fibre distributions occurs on the top right of the tumour region, situated close to the bulk of the tumour mass. As the tumour expands in size, the fibres are subsequently degraded during each stage of evolution, Figures 6 and 7, illustrating the simulations after $50\Delta t$ and $75\Delta t$, respectively. At stage $50\Delta t$, some cells have detached and formed a region encircling the initial bulk of tumour cells, subfigure 6(a). The boundary of the tumour is expanding into the surrounding tissue, subfigure 6(b), exhibiting a “rippling” effect along the proliferating edge caused by the fibre mediated movement of the boundary. The macroscopic fibre distribution becomes depressed as the cancer cells increase in distribution and thus increase in their degradative behaviour, subfigure 6(c). Moving on to final stage $75\Delta t$, the tumour has largely increased in size and there are dense regions of cell distribution, subfigure 7(a). The cells are gathering at the tumour interface, pulled in this direction by cell-fibre adhesion, following the direction of the fibres observed at the previous interval in subfigures 6(d), 6(f). The source of MMPs induced by the cancer cells are very high in these dense areas, therefore the degradation of fibres is higher, subfigure 7(c), witnessed here by the absence of macroscopic fibre distribution.

To investigate the effects of fibre distribution, whilst keeping the non-fibre ECM phase initially homogeneous (21), we initialise the fibre distribution $F(x,0)$ with $15\%$, or $p=0.15$, of the heterogeneous distribution (22). The simulations at stage $25\Delta t$, Figure 8, indicate the initial distribution of fibres is significant during the evolution of a tumour. The primary bulk of cancer cells have dispersed within the boundary into regions of high cell distribution, subfigure 8(a), these coinciding with low density regions of the fibres ECM phase. The boundary of the tumour is circular with small defects in the direction of the fibre orientation, subfigure 5(f). The fibre magnitude density in subfigure 8(c) has small regions of high fibre distributions, again pushed outwards towards the proliferating edge. Due to the heterogeneity of the fibre distributions, the fibre orientations are subject to an increased degree of realignment, subfigures 8(d), 8(f), attributed to the initially different levels of density at each macro-spatial position. The fibres are realigned towards the higher density regions, where “frenzied” groups of fibre orientations can be observed. Moving on to simulations at stage $50\Delta t$, Figure 9, many of the behaviours previously observed at stage $25\Delta t$ are magnified. The cancer bundle has increased in size and notably spread further into the areas of initially low fibre density, subfigure 9(a) with the tumour boundary also increasing in size and irregularity. The fibres are being rearranged and pushed further outwards towards the tumour boundary, subfigure 9(c), whilst simultaneously undergoing macroscopic degradation at the tissue scale (8), thus resulting is areas of very low to no fibre density. Finally, the simulations in Figure 10 show the evolution of the tumour at stage $75\Delta t$. The cancer cells are forming patterns within the tumour boundary, subfigure 10(a), in the areas of low ECM density, and cells are migrating towards the tumour boundary collecting in high distribution bundles that more rapidly degrade the surrounding ECM. Furthermore, the cells are flooding areas where there is very low to no fibre density, subfigure 10(c), and forming dense bundles of cells. This behaviour is in accordance with the conclusions presented in Shutt_2018 that cancer cells can more freely invade areas of no ECM density and will progress upon these areas first before engulfing the higher density regions.

6.1 Increased collagen density

We proceed by exploring the cancer dynamics within an initial $20\%$ homogeneous fibre distribution, taken as $p=0.2$ of the non-fibres ECM phase $l(x,0)$, and in the presence of the micro-fibre degradation rate $d_{f}=0.5$. At stage $25\Delta t$, Figure 11, the non-fibres ECM phase has been degraded by the cancer cells, subfigure 11(b), with the highest degradation occurring in the regions of highest cell distribution. The tumour boundary at this stage is larger than in previous simulations, this correlating with results in Shutt_twopop ; Shutt_2018 where an initially higher density fibre ECM phase resulted in an accelerated spread of the tumour. Additionally, the macroscopic fibre density, subfigure 11(c), also exhibits different behaviour than in previous simulations, where a similar pattern of fibres was noted in Shutt_twopop ; Shutt_2018 . This pattern of fibres is witnessed because the tumour boundary is expanding faster than the micro-fibres are being rearranged, thus the distributions are not found to build on the proliferating edge. Moving on to stage $50\Delta t$, the cancer cell distribution is increasing and spreading within the tumour region, subfigure 12(a). There is significant non-fibres ECM degradation stretching the entire area of the tumour, subfigure 12(b), with this only becoming more pronounced at later stages, subfigure 13(b).The macroscopic fibre orientations, subfigures 12(d), 12(f), display high levels of reorientation with a similar “frenzied” appearance to figures in Shutt_twopop . At the final stage $75\Delta t$, Figure 13, the cancer cells are dispersed further into the matrix, subfigure 13(a) and noticeably in the pattern of degradation of the non-fibres ECM phase, subfigure 13(b). Overall, we conclude that in the presence of a homogeneous ECM, with a high initial homogeneous fibre distribution, the cancer cells have increased opportunity for adhesion and thus they are able to easier invade the surrounding matrix, resulting in a larger tumour region and a higher level of ECM degradation.

We continue our investigation of cancer invasion within a collagen dense environment by exploring the multiscale model within an initially $20\%$ heterogeneous fibre ECM phase embedded in a homogeneous non-fibres ECM phase. Figures 14-16 display simulation results at stages $25\Delta t$, $50\Delta t$ and final stage $75\Delta t$. When comparing with the results in Figures 11-13, the significant differences when in a homogeneous or heterogeneous fibre environment is the shape and density of the main tumour bulk and the pattern of the proliferating boundary. The tumour boundary in Figure 14 is spreading first to the low density regions of fibres, this behaviour consistent with previous results in Shutt_2018 , where the cancer cells can more easily migrate to close-by areas of low fibre density due to the physical space available, thus the boundary of the tumour becomes lobular as the cells migrate outwards, subfigure 14(c). The main bulk of tumour cells are also exhibiting this behaviour as they have formed a high distribution region of cells in the area of lowest fibre density, subfigure 14(a). The macroscopic fibre orientations have been rearranged to direct the movement of the tumour boundary, witnessed in Figure 15 at the next time stage interval $50\Delta t$. The main body of the tumour has distributed the cells into several small high distribution bundles, subfigure 15(a) particularly in the low density regions of fibres, subfigure 15(c). The microscopic fibres are continuously rearranged and due to a higher rate of fibres degradation where the cancer cell distribution is highest, we observe a very low central region of fibre density with them being both degraded and pushed outwards from the central region of the tumour, subfigures 15(d), 15(f). The orientation of the fibres become increasingly erratic in areas of higher fibre density, witnessed in the protrusions of the tumour boundary, whereby they are orientated in opposing directions. This trait is exaggerated at final stage $75\Delta t$, Figure 16, where the tumour has spread further into the ECM and the protrusions have increased in size. Furthermore, the bulk of cancer cells have become separated and formed distinct regions within the tumour boundary. When comparing the cancer cell distribution and the macro-fibres density, subfigures 16(a) and 16(c) respectively, it can be noted that the cancer cells have bypassed the higher regions of fibres and formed bundles of cells around these areas. This behaviour can also be seen in the homogeneous case, however it is a more prominent feature when in the presence of a heterogeneous fibre ECM phase. In conclusion, in the presence of an initially high fibre ECM density, tumour progression is accelerated and encourages a more aggressively spreading tumour, in both the case of either a homogeneous or heterogeneous fibre distribution. Additionally, when the initial fibre density is high, the boundary of the tumour spreads further away from the main body of the tumour at a rate which the cells cannot maintain. Hence, the main body of the tumour stays centralised, compared to a lower initial fibre density, whereby the cancer cells spread at a consistent rate within the tumour region and stay close to the boundary.

6.2 Invasion patterns of a heterotypic cancer cell population

To conclude our exploration of cancer invasion within a heterogeneous microenvironment, we focus now on the invasion patterns of a heterotypic cell population using the two cell sub-populations macro-dynamics introduced in shutt_chapter and expanded upon in Shutt_twopop . The macro-dynamics of the two cell sub-populations are similar in flavour and can be mathematically expressed as

[TABLE]

where the individual terms retain their exact meaning from Section 2.2, while $M_{c}(\textbf{u},t)$ represents the mutation of cells from cell population $c_{1}$ to cell sub-population $c_{2}$, defined in Shutt_twopop ,

[TABLE]

Here the mutations are considered to be dependent on the underlying fibre densities, with $\kappa$ representing a certain level of ECM beyond which mutations can occur, $H(\cdot)$ being the usual Heaviside function, with $t_{m}$ being the time at which mutations begin. Furthermore, we consider both cancer cells sub-populations to contribute to the source of MMPs on the tumour boundary, hence the microscopic source term (11) will be readdressed as

[TABLE]

where $\alpha_{1},\alpha_{2}$ are the MMP secretion rates for cell populations $c_{1}$ and $c_{2}$, respectively, $\widetilde{F}(x,t_{0}+\tau):=1+F(x,t)$ and $\widetilde{F}\cdot\lambda$ are defined as in Section 3. Using parameter set $\Sigma_{1}$ from Appendix D, the micro-fibre degradation rate $d_{f}=1$ and the cell adhesion matrices

[TABLE]

in following Figures 17 - 22 we present the computational results for the evolution of: (a) the primary cancer cell sub-population; (b) the secondary cell sub-population; (c) the fibre magnitude density; (d) the non-fibres ECM distribution; and of the vector field of orientated fibres at two different resolutions, namely: (e) coarsened twice; and (f) coarsened fourfold. Considering an initially homogeneous fibre ECM phase, when comparing these results directly with the results in Shutt_twopop , we observe no differences between simulations performed in the absence of micro-fibre degradation at the tumour interface at the first interval $25\Delta t$, Figure 17, with the fibre orientations consistent between the results, subfigures 17(e), 17(f). Proceeding to later stages, $50\Delta t$ and $75\Delta t$, Figures 18 and 19, respectively, we begin to witness changes, whereby the boundary of the tumour is visibly smaller than in previous results, thus implying a slower progression of the tumour. The macroscopic fibre orientations are aligned differently in subfigures 18(e), 18(f), where the fibres are directed inwards to the centre of the tumour, confining the tumour to the centre of the domain, subfigure 17(a). Moving on to the final stage, Figure 19, the tumour region is smaller than in previous results in Shutt_twopop . The bulk of the cancer cells stick closely to the tumour boundary and both the fibre and non-fibre ECM phase have undergone a higher level of degradation within the tumour region.

Finally, to further our understanding of the effects of micro-fibre degradation at the tumour interface, we explore the evolution of a heterotypic cell population within a heterogeneous fibre ECM phase, whilst the non-fibre ECM phase remains homogeneous, Figures 20-19. As with an initially homogeneous fibre density, at stage $25\Delta t$, Figure 20, there is very little difference when compared to simulations presented in Shutt_twopop . Moving on to later stages it is obvious that the process of boundary micro-fibre degradation causes a slower rate of tumour progression. The tumour region is considerably smaller in Figures 21 and 22 when compared with previous results Shutt_twopop . The bulk of tumour cells display a similar pattern, however much closer to the tumour boundary, particularly cancer cell population $c_{2}$ where the cells have formed high distribution bundles at the leading edge. Much like in the presence of a homogeneous fibre density, both the fibre and non-fibre ECM phase has been subject to a higher level of degradation in the tumour region. This is attributed to the micro-fibre degradation at the tumour boundary, as degradation continuously occurs as the boundary expands we witness a lower fibre density as the tumour evolves. It can be concluded from these simulations and comparisons with previous results that the process of micro-fibre degradation at the tumour interface is disadvantageous to tumour progression, inhibiting the full invasive capabilities of the tumour. This is due to the lower levels of fibre density at the tumour interface, as degradation of the micro-fibres occurs both inside and within the peritumoural region of the tumour, the cell-fibre adhesion rate is reduced in line with low fibre levels, therefore the cancer cells do not have the same opportunities for adhesion and thus their migration is greatly reduced.

7 Discussion

In this paper we have presented an integrated two-part multiscale model of cancer invasion, which builds on the approach introduced in Shutt_2018 and extends that to capture explicitly the dynamic cell-scale interaction between the MDE boundary micro-dynamics and the peritumoural mass distribution of micro-fibres.

Structured largely similar to the modelling framework introduced in Shutt_2018 , the model proposed here combines two multiscale systems that share the same tissue (macro-) scale dynamics while having separate cell- (micro-) scale processes that are simultaneously connected via two double feedback loops. Specifically, while this new modelling framework inherits completely the multiscale dynamics of naturally oriented ECM fibres (induced by the mass distribution of micro-fibres) occurring on the topological closure of the invading tumour (and including the dynamic rearrangement of fibres under the incidence of the macro-scale flux of cancer cells), this shares its macro-dynamics with a significantly extended multiscale moving-boundary modelling for the proteolytic dynamics at the tumour invasive edge that explicitly considers the interaction with the peritumoural fibres.

This model expands and takes forward both the initial multiscale moving-boundary framework introduced in Dumitru_et_al_2013 and its further development into the two-part multiscale modelling introduced in Shutt_2018 by bringing in and exploring the cell-scale interactions between the cross-interface diffusion of MMPs and the micro-fibre distributions in the peritumoural region, with direct impact upon microscopic peritumoural degradation of micro-fibres that results in a continuously altered macroscopic vector field of oriented ECM fibres at the tumour boundary. Moreover, these altered peritumoural ECM fibres have major relevance within the macroscopic dynamics of the cancer cells as this affects the cell-fibres adhesion properties at the leading edge of the tumour, impacting this way not only the tumour mechanics close to the tumour interface but the entire-tissue scale dynamics of the tumour.

To that end, we first explore mathematically the positive feedback that the macroscopic distribution of ECM fibres close to the tumour interface has upon the emerging cell-scale source of MMP-2 for the cross-interface micro-dynamics that MMP-2 exercise at the invading edge of the tumour. Specifically, in this new formulation we are able to capture the enhanced sources of MMP-2 in a relevant cell-scale neighbourhood which are enabled non-locally by the presence of elevated distributions of ECM fibres within neighbouring active regions from within the outer proliferating rim of the tumour where cancer cells arrive during their macro-dynamics and produce MMPs.

Further, in the presence of the cell-scale MMP-2 source induced by the macro-dynamics, a cross-interface diffusion of MMP-2 occurs at the invasive edge of the tumour. However, as the MMP-2 find it easier to diverge along their gradient directions in regions with lower micro-fibres levels, by accounting for the presence micro-scale mass distribution, we finally obtain that the diffusion rate of this diffusive molecular transport process of MMP-2 naturally depend on the micro-fibres density. Thus, this cell-scale MMP-2 dynamics focuses the cross-interface molecular transport towards the regions of lower mass-distributions of micro-fibres, taking on-the-fly advantage on the potential “micro-fibres valleys” created by the multiscale dynamic rearrangement of fibres induced by the macro-scale flux of cancer cells, which was derived and explored with full details in Shutt_2018 .

Since the MMP-2 cross-interface molecular transport leads to peritumoural micro-fibres degradation at the micro-scale, we explored this degradation explicitly at the cell scale, remarking here at the same time that this complements the previous modelling framework introduced and discussed in Shutt_2018 where the fibres degradation was only considered on the bulk of the tumour at the macro-scale. To that end, we considered the correlation between the micro-fibres degradation and the incidence angle that the MMP-2 molecular flux makes with the regions of significant levels of micro-fibres distributed with any given $\epsilon Y$ from the covering bundle of $\{\epsilon Y\}_{\mathcal{P}}$ boundary micro-domains. This enabled us to derive mathematically a micro-fibres degradation law occurring on each fibre micro-domain $\delta Y(x)$ that has non-empty intersection with at least one of the boundary micro-domains $\epsilon Y$, in which maximum fibre degradation occurs when the angle between the fibres and MMPs flux is perpendicular, while the degradation decreases with increasing alignment of the fibres with the direction of the MMP-2 molecular flux. This suggests that the highly aligned collagen fibrils will act as a pathway for invasion rather than a barrier against it, this being consistent with the biological evidence presented in Provenzano_06 .Finally, this degradation of peritumoural mass distribution micro-fibres at the cell-scale is continuously in time translated back at macro-scale, having a natural and major impact upon on-the-fly changes in the orientation and magnitudes of macro-scale ECM fibres from the peritumoural region.

This new modelling framework has been explored in several scenarios, within both a homogeneous and heterogeneous initial distribution of fibres, and varying the initial ratio of macroscopic fibre distribution in relation to the non-fibres ECM phase. The non-fibres ECM phase was kept as a homogeneous density throughout the paper for the purpose of exploring only the influence of ECM fibres during cancer invasion. These scenarios were explored through randomly allocated distributions of micro-fibres over the fibres micro-domain, as considered already in Shutt_twopop and defined in Appendix C.

We explored the differences between a homogeneous and heterogeneous initial distribution of fibres where we varied the initial percentage of fibres density from $15\%$ to $20\%$, as well as investigating the morphology of a heterotypic cancer cell population whose macroscopic dynamics were developed in Shutt_twopop . We conclude from these simulations that a heterogeneous distribution of fibres induces a more lobular, fingering pattern of the tumour boundary, and an increase in initial fibre density promotes a more aggressive tumour spreading further into the surrounding tissue, a behaviour which is mirrored in the biological experiments performed in Provenzano_2008 . These remarks are in line with previous results in Shutt_2018 ; Shutt_twopop where the same conclusions are drawn regarding the initial condition of the fibres ECM phase. The simulations performed in this paper exhibit an overall larger tumour spread than the simulations in Shutt_2018 , implying that the degradation of fibres at the tumour interface promotes tumour invasion. The final simulations, which consider the invasive behaviour of a heterotypic cell population, suggest that although the cancer cells migrate more easily into low density regions of ECM, at the tumour interface a lack of fibre density is detrimental to the progression of the tumour. The low levels of fibre density inhibit the migration of cancer cells and thus the movement of the tumour boundary by reducing the opportunities for cell-fibre adhesion. In general, we conclude that the invasion of a heterotypic cancer cell population is accelerated in the presence of a high fibres ECM density, however within a low fibres ECM the tumour undergoes slower progression and the bulk of the cancer cells remain closer to the gradually expanding tumour boundary.

Looking forward, to advance with this model we would look to investigate the full MT1-MMP/MMP-2 cascade, namely, to include within the model the tissue-inhibitor of matrix-metalloproteinases-2, TIMP-2, a key molecule required for the activation of pro-MMP-2 Seiki_2003 , considering the role and regulation of these molecules and the resulting effects on peritumoural tissue degradation. Additionally, further exploration of the fibres network and its structure would permit better, more realistic modelling of human tissue, allowing the model to be compared with current biological experiments, for example, with the experiments performed in Provenzano_2008 that investigate the effects of increasing collagen density in the surrounding matrix.

Appendix

Appendix A The radial kernel $\mathcal{K}(\cdot)$

To explore the influence on adhesion-driven migration decreases as the distance from $x+y$ to $x$ within the sensing region $B(x,r)$ increases, the expression of the radial dependent spatial kernel $\mathcal{K}(\cdot)$ appearing in (5) is taken here to be:

[TABLE]

Appendix B The mollifier $\psi_{\gamma}$

The standard mollifier $\psi_{\gamma}:\mathbb{R}^{N}\to\mathbb{R}_{+}$ (which was used also in Shutt_2018 ; Dumitru_et_al_2013 ) is defined as usual, namely

[TABLE]

where $\psi$ is the smooth compact support function given by

[TABLE]

Appendix C Microscopic fibre domains

For the fibres initial conditions,, on the micro-domains $\delta Y(x)$, we consider a family of five distinctive micro-fibres patterns, $\{P^{1}_{i}\}_{i\in J}$, which are defined by the union of paths $P^{1}_{i}=\bigcup\limits_{j=1..5}h^{1}_{i,j}$, which are given as follows.

For the family of fibre paths $P^{1}$, we have:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

For this family of fibre paths $P^{1}$, as described in Shutt_2018 , the micro-scale fibres pattern within each micro-domain $\delta Y(x)$ is given as

[TABLE]

where $\{\psi_{h^{l}_{i,j}}\}_{i,j=1..5}$ are smooth compact support functions of the form

[TABLE]

Here $r>0$ is the width of the micro-fibres and $C_{h^{l}_{i,j}}$ are constants that determine the maximum height of $\psi_{h^{l}_{i,j}}$ along the smooth paths $\{h^{l}_{i,j}\}_{i,j=1..5}$ in $\delta Y(x)$. Finally, $\psi_{\gamma}$ is the standard mollifier defined in Appendix B, with $\gamma_{{}_{0}}=h/16$.

Finally, the initial spatial configuration of the pattern of macroscopic ECM fibre phase is selected according to a randomly generated matrix of labels $A=(a_{i,j})_{i,j=1..n}$ corresponding to the entire $n\times n$ grid discretising $Y$, in which the entries $a_{i,j}$ are allocated values randomly selected from the set of configuration labels $\{1,2,3,4,5\}$ that will dictate the choice of micro-fibres pattern among those described above that will be assigned to the micro-domains $\delta Y(j\Delta x,i\Delta y)$, for all $i,j=1..n$.

Appendix D Table for the parameter set $\Sigma_{1}$

Here we present a table for the parameter set $\Sigma_{1}$.