Table of Contents Author Guidelines Submit a Manuscript
Computational and Mathematical Methods in Medicine
Volume 2011 (2011), Article ID 452320, 16 pages
http://dx.doi.org/10.1155/2011/452320
Research Article

Numerical Solutions for a Model of Tissue Invasion and Migration of Tumour Cells

1Department of Mathematics and Computer Science, University of Warmia and Mazury, Zolnierska 14, 10-561, Olsztyn, Poland
2Department of Mathematics, Boise State University, 1910 University Drive, Boise, ID 83725, USA

Received 15 January 2010; Accepted 18 July 2010

Copyright © 2011 M. Kolev and B. Zubik-Kowal. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

The goal of this paper is to construct a new algorithm for the numerical simulations of the evolution of tumour invasion and metastasis. By means of mathematical model equations and their numerical solutions we investigate how cancer cells can produce and secrete matrix degradative enzymes, degrade extracellular matrix, and invade due to diffusion and haptotactic migration. For the numerical simulations of the interactions between the tumour cells and the surrounding tissue, we apply numerical approximations, which are spectrally accurate and based on small amounts of grid-points. Our numerical experiments illustrate the metastatic ability of tumour cells.

1. Introduction

The analysis of data obtained from the World Health Organization (WHO) [1] and the UN [2] databases shows that, at present, cancer is and probably will remain to be among the leading causes of death worldwide [35] being surpassed only by cardiovascular diseases. According to the data provided by the WHO, cancer disease is the cause of the death of roughly six million people yearly [1]. This explains the major significance of the fight against the malignant conditions, which includes prevention [6], cure [7, 8], and cancer research.

Tumour development is a very complex multistep process involving many intracellular and extracellular phenomena which are strongly nonlinear and time varying [4, 911]. Genomic changes as well as microenvironmental factors such as the extracellular matrix (ECM), various growth factors, and substrate concentrations have been shown to play a major role in the process of carcinogenesis [12].

Generally, tumours can be classified as benign and malignant. The growth of benign tumours is self-limiting and their cells tend to stay in the same place. Malignant tumours may grow without limitations and their constituent cells are prone to migrate or metastasize to other parts of the organism [1315]. The ability of malignant cancer to invade the local tissue and to spread throughout the organism is their most insidious and dangerous property. Metastasis is the predominant cause of most cancer deaths [14, 16, 17].

The process of metastasis includes angiogenesis and invasion. Tumour angiogenesis (rapid growth of blood vessels near the tumour cells) is induced by a secretion of various growth factors such as vascular endothelial growth factor (VEGF). These vessels facilitate the influx of oxygen and other nutrients needed for the development of the cancer [18]. The process of angiogenesis is followed by invasion and penetration of cancer cells into surrounding tissues and possibly by dissemination of cancer cells through blood vessels. Thus, tumour cells can be carried to a distant site of the body. There they can implant and initiate the development of a secondary tumour [14, 16, 19]. An important role in the process of cancer invasion is performed by matrix degradative enzymes (MDEs) such as metallo-proteases (MMPs). They are produced by tumour cells and digest the ECM, which enables the migration of cancer cells through the tissue [13, 14, 17].

In the last half century, many mathematical models describing the process of tumourigenesis have been the subject of active research. Mathematical and computational methods have contributed to clarifying the factors that are sufficient to explain experimental and clinical data, to defining these factors in precise terms and to suggesting experiments for calculation of these factors [20]. In addition, analyses and simulations of mathematical models have been used for the reduction of the amounts of costly experiments needed for the development of therapies [21, 22]. It is strongly believed that mathematical and computational methods will play a significant role in cancer research in the future. They may improve the understanding of some complicated features and details of tumour evolution as well as be effectively used in clinical laboratories, by means of appropriate model-based decision support systems [4]. We refer the readers to special issues [2327] for more complete bibliography regarding the applications of mathematical and computational methods to cancer research.

Gatenby and Gawlinski present one of the first models of tumour invasion in the papers [28, 29]. Gatenby [28] considers the competition between healthy host cells and modified (tumour) cells and proposes and analyses several models formulated in terms of ordinary differential equations. Gatenby and Gawlinski [29] present a reaction-diffusion model for the investigation of the role of the alteration of the microenvironmental acidity induced by cancer cells for their invasion into the organism. Subsequently, the series of papers, among others [3043], have appeared offering models and detailed analysis of diverse features of cancer invasion. In this paper, we study the continuum models of avascular tumour growth investigated by Chaplain et al. (cf., e.g., [31, 3437]). The first model of this series is proposed in Anderson et al. [31]. The authors consider three major variables involved in the process of cancer invasion, namely, cancer cells, ECM, and MDEs. In order to study in detail mainly the influence of the surrounding tissue on the process of migration of tumour cells, the proliferation of the latter is not included in the continuum model. The authors analyse numerically in one and two dimensions the impact of ECM gradients resulting from the destruction of ECM by MDE and the role of haptotaxis on cancer invasion. An extension of this model is presented in Chaplain and Anderson [34] who consider the role of oxygen as a nutrient for the tumour cells. The authors propose also a new model equation for endogenous inhibitors, such as tissue inhibiting MMPs, that can neutralize MDEs. We include this equation in our model (8), see Section 2 below. The model of Chaplain and Anderson [34] has been further developed by Lolas [37] and Chaplain and Lolas [35, 36] who have considered terms describing chemotaxis, proliferation of cancer cells and reestablishment of the ECM. Lolas [37] examines a variety of continuum models, in particular incorporating the effects of just chemotaxis, and haptotaxis, and their combination, and so forth. One of the conclusions of the author is that the mechanism of chemotaxis without haptotaxis cannot lead to a successful cancer invasion if there is no proliferation of tumour cells and reestablishment of ECM. Further novel ordinary differential equations that describe the cancer cell proliferation and the remodeling of the extracellular matrix re-establishment function allowing the incorporation of the plasminogen activation cycle are included in the model of Chaplain and Lolas [36] that also investigates the role of the uPA system for the cancer invasion. uPA inhibitors and plasmin have also been investigated in the model by Chaplain and Lolas [35]. Clear and detailed description of the biological processes observed during the cancer invasion and metastasis is provided in [31, 3437]. In particular, in these paper, the key stages of the metastatic cascade, the structure and functions of the major constituents of the ECM and the basic representatives of the MDEs participating in the interactions between the healthy and cancer cells are systematically presented on the basis of broad theoretical and experimental bibliography.

In our paper, we propose a different numerical approach than the approach used, for example, in [31, 3437]. The goal of the paper is to obtain numerical results which are based on small amounts of spatial grid points applied to the model equations so that low-dimensional vectors of data are used to make the numerical computations fast. We construct a new algorithm for the systems [31, 34, 36, 37] by using spectrally accurate approximations to the terms that model the tumour cell random motility, the haptotaxis, the MDE diffusion, and the diffusion of the endogenous inhibitors. Since the algorithm computes the solutions with spectral accuracy, it is based on smaller amounts of spatial grid points than the amounts of grid points used for the less accurate finite difference approximations (strategy applied, e.g., in [31, 3437]), which consequently saves computational time. The idea of using small amounts of spatial grid point and saving time for computing one solution for one set of parameters, which has to be repeated many times for many sets, is important, for example, for the numerical experiments carrying the goal of estimating parameter values from laboratory data. This idea is applied in [44] to estimate parameter values of one of the models presented in [36, 37] from the in vivo experimental data [45] developed by using transgenic mouse models. The numerical approach from [44] is based on a different approximation to the haptotactic term than the approximations used in this paper and our numerical schemes are constructed for systems which are various variants and generalizations of the model investigated in [44]. Furthermore, because of considering different variants of boundary conditions the schemes in this paper differ from that of the paper [44].

Additionally to the model presented in [36, 37] and applied in [44], in this paper, we investigate other models, which are presented in [34] or are combinations of the model equations from [34, 36, 37]. Moreover, in [44], the parameter values are evaluated quantitatively from the laboratory data [45] so that the solutions of the model equations correlate with the data. Contrarily to [44], in this paper, we choose the parameter values qualitatively in order to observe and compare solutions computed with different parameters. This comparison allows to analyse the influence of the parameters on the shape of the solutions and we conclude that complicated interactions between tumour cells, ECM, MDEs, and endogenous inhibitors can be directed by choosing the parameter values. Our sequence of numerical simulations is initiated from the solutions obtained with the parameter values chosen in [34] (for comparison) and next we gradually change the values and analyse their influence on the solutions. Animated graphical visualization of the solutions and how they change according to the parameters is helpful in observing the influence of the parameters on the shape of the solutions. The idea of using small amounts of spatial grid points and saving time for computing solutions of the model equations is crucial in the effective utilization of, for example, animated simulations of tumours, which can be used as a predictive and visualized tool in clinical applications. Decreasing the amounts of spatial grid points used for such visualizations saves not only the time of demonstrations but also the computer memory. It is not possible to demonstrate the animated simulations in papers and we only note that they are interesting and help in visualization of the complicated biological processes. Instead of the animated simulations we include snapshots at different stages in time.

The contents of this paper is as follows: the model equations are described in Section 2, the algorithm is introduced in Section 3, the results of numerical experiments and simulations are presented in Section 4, and Section 5 includes our concluding remarks and future research work.

2. Mathematical Model

In this section, we investigate various models of tissue invasion by cancerous cells. In Section 2.1, we investigate the Chaplain and Anderson model [34] focusing on interactions between ECM and cancer tumour and metastatic abilities of cancer cells. In Section 2.2, we investigate further expansions of the model and its different versions with additional terms connected with proliferation of tumour cells, ECM renewal, and different functions modelling the production of MDEs by the tumour cells. Section 2.3 deals with a more general model with an additional equation, which describes evolution of endogenous inhibitors.

2.1. Cell-Matrix Interactions and Cell Migration

In the next section, we construct a numerical scheme for the following model of tissue invasion:𝜕𝑛𝜕𝑡=𝑑𝑛𝜕2𝑛𝜕𝑥2𝜕𝛾𝑛𝜕𝑥𝜕𝑓,𝜕𝑥𝜕𝑓𝜕𝑡=𝜂𝑚𝑓,𝜕𝑚𝜕𝑡=𝑑𝑚𝜕2𝑚𝜕𝑥2+𝛼𝑛𝛽𝑚(1) with the space variable 𝑥 belonging to the scaled domain [0,1] of tissue, and time 𝑡. The model equations (1) describe interactions between tumour cells, MDEs, and ECM. The interacting variables are 𝑛-tumour cell density, 𝑓-ECM density, and 𝑚-MDEs concentration. The system (1) is derived in detail in [34] and is a part of a more general system consisting of (1) with an additional fourth equation for endogenous inhibitor concentration denoted by 𝑢. In [34], it is assumed that the tumour cells, the MDEs, and the inhibitors remain within the space domain and zero-flux boundary conditions are imposed. The fourth equation for the endogenous inhibitor concentration is dropped under the additional assumption that negative effect of the endogenous inhibitors is overcame by the MDEs in an actively invading tumour. This assumption implies that 𝑢=0 and the general system of four equations is reduced to (1). In Section 2.3, we investigate the model with all four equations.

2.2. Migration and Proliferation of Cancer Cells, ECM Renewal, and MDE Production

We also investigate further expansions of the model (1), which are introduced, for example, in [36, 37]. After adding the proliferation term 𝜇1𝑛(1𝑛𝑓) to the right-hand side of the equation governing tumour cell motion (the first equation in (1)) and the ECM renewal term 𝜇2𝑓(1𝑛𝑓) to the right-hand side of the equation for the ECM (the second equation in (1)), we obtain the following model:𝜕𝑛𝜕𝑡=𝑑𝑛𝜕2𝑛𝜕𝑥2𝜕𝛾𝑛𝜕𝑥𝜕𝑓𝜕𝑥+𝜇1𝑛(1𝑛𝑓),𝜕𝑓𝜕𝑡=𝜂𝑚𝑓+𝜇2𝑓(1𝑛𝑓),𝜕𝑚𝜕𝑡=𝑑𝑚𝜕2𝑚𝜕𝑥2+𝛼𝑛𝛽𝑚,(2) where 𝜇1 is the proliferation rate of the tumour cells and 𝜇2 is the growth rate of the ECM.

We also make experiments with the following modification of (2):𝜕𝑛𝜕𝑡=𝑑𝑛𝜕2𝑛𝜕𝑥2𝜕𝛾𝑛𝜕𝑥𝜕𝑓𝜕𝑥+𝜇1𝑛(1𝑛𝑓),𝜕𝑓𝜕𝑡=𝜂𝑚𝑓+𝜇2𝑓(1𝑛𝑓),𝜕𝑚𝜕𝑡=𝑑𝑚𝜕2𝑚𝜕𝑥2+𝛼𝑛(1𝑛)𝛽𝑚,(3) where the MDE production is modeled by 𝛼𝑛(1𝑛). The motivation for choosing such form of the MDE production in [36, 37] follows from experimental observations of polarized expression of MDEs at the invading leading edge of tumour, see, for example, Estreicher et al. [46].

We investigate the model equations (1), (2), and (3) supplemented by the zero-flux boundary conditions𝜕𝑛𝛾𝜕𝑥(0,𝑡)=𝑑𝑛𝑛(0,𝑡)𝜕𝑓𝜕𝑥(0,𝑡),𝜕𝑚𝜕𝑥(0,𝑡)=0(4) at 𝑥=0 and either the Dirichlet conditions𝑛(1,𝑡)=0,𝑚(1,𝑡)=0(5) or the zero-flux boundary condition𝜕𝑛𝛾𝜕𝑥(1,𝑡)=𝑑𝑛𝑛(1,𝑡)𝜕𝑓𝜕𝑥(1,𝑡),𝜕𝑚𝜕𝑥(1,𝑡)=0,(6) at 𝑥=1. As in [34], we assume that the initial tumour is centered at 𝑥=0, the initial MDE concentration is proportional to the initial tumour cell density with 1/2 as the constant of the proportionality, and the MDE has already degraded the ECM, thus we consider the same initial conditions as in [34], which are the following:𝑛(𝑥,0)=exp𝑥2𝜀,𝑓(𝑥,0)=10.5𝑛(𝑥,0),𝑚(𝑥,0)=0.5𝑛(𝑥,0),(7) for 𝑥[0,1]. The parameter values for the model equations are specified in Section 4.

2.3. Production of Endogenous Inhibitors

We additionally consider the general model𝜕𝑛𝜕𝑡=𝑑𝑛𝜕2𝑛𝜕𝑥2𝜕𝛾𝑛𝜕𝑥𝜕𝑓𝜕𝑥+𝜇1𝑛(1𝑛𝑓),𝜕𝑓𝜕𝑡=𝜂𝑚𝑓+𝜇2𝑓(1𝑛𝑓),𝜕𝑚𝜕𝑡=𝑑𝑚𝜕2𝑚𝜕𝑥2+𝛼𝑛𝜃𝑢𝑚𝛽𝑚,𝜕𝑢𝜕𝑡=𝑑𝑢𝜕2𝑢𝜕𝑥2+𝐹(𝑚,𝑓)𝜃𝑢𝑚𝜌𝑢,(8) where the last equation describes evolution of endogenous inhibitors (concentration of which is denoted by 𝑢). This equation is the fourth equation in the model (10.5) proposed by Chaplain and Anderson in [34], where it is assumed that endogenous inhibitors are produced by ECM as a response to the MDEs and the function 𝐹(𝑚,𝑓) models the inhibitor production. The term 𝜃𝑢𝑚 models neutralization of the MDEs and 𝜌𝑢 models decay of the inhibitors. We assume that the initial inhibitor concentration is 𝑢(𝑥,0)=0(9) and impose the zero-flux boundary conditions𝜕𝑢𝜕𝑥(0,𝑡)=𝜕𝑢𝜕𝑥(1,𝑡)=0.(10)

Our goal is to construct a new efficient algorithm for solving the models (1), (2), (3), and (8) and investigate the ability of cancer cells to produce and secrete the MDE, which then degrade the ECM, and allow the cells to start their migration towards healthy parts of the tissue.

3. Construction of Numerical Approximations to Tumour Cells, ECM, and MDEs

In this section, we construct numerical solutions to the model equations (1), (2), and (3) supplemented by the initial conditions (7) and the boundary conditions (4) and (5). For the numerical solutions, we consider the Chebyshev-Gauss-Lobatto points𝑥𝑖=1212cos𝑖𝜋,𝑁+1(11) with 𝑖=0,1,,𝑁+1, in the scaled domain [0,1] of tissue. Our goal is to construct approximations to 𝑛(𝑥𝑖,𝑡), 𝑓(𝑥𝑖,𝑡), and 𝑚(𝑥𝑖,𝑡), for 𝑖=0,1,,𝑁; the values of the solutions at 𝑥𝑁+1 are known from (5).

Let𝑛𝑥𝑛(𝑡)=0𝑛𝑥,𝑡1𝑛𝑥,𝑡𝑁,,𝑡𝑑𝑛(𝑡)=𝑑𝑡𝜕𝑛𝑥𝜕𝑡0,𝑡𝜕𝑛𝑥𝜕𝑡1,𝑡𝜕𝑛𝑥𝜕𝑡𝑁,𝑛,𝑡𝑥(𝑡)=𝜕𝑛𝑥𝜕𝑥0,𝑡𝜕𝑛𝑥𝜕𝑥1,𝑡𝜕𝑛𝑥𝜕x𝑁,𝑡,𝑛𝑥𝑥𝜕(𝑡)=2𝑛𝜕𝑥2𝑥0𝜕,𝑡2𝑛𝜕𝑥2𝑥1𝜕,𝑡2𝑛𝜕𝑥2𝑥𝑁,,𝑡(12) and we use similar notations for 𝑓 and 𝑚. We shall replace the spatial derivatives in (1) by numerical approximations constructed for the vectors 𝑛𝑥𝑥(𝑡), 𝑛𝑥(𝑡), 𝑓𝑥(𝑡), 𝑓𝑥𝑥(𝑡) in the first equation and for 𝑚𝑥𝑥(𝑡) in the third equation. For 𝑛𝑥(𝑡), we apply the following spectrally accurate approximations𝜕𝑛𝑥𝜕𝑥𝑖,𝑡𝑁+1𝑗=0𝑑𝑖,𝑗𝑛𝑥𝑗,,𝑡(13) with 𝑖=0,1,,𝑁+1, where𝑑𝐷=𝑖,𝑗𝑁+1𝑖,𝑗=0(14) is the first-order differentiation matrix based on the points (11), see [47, 48]. We also apply the analogous spectrally accurate approximations𝜕𝑓𝑥𝜕𝑥𝑖,𝑡𝑁+1𝑗=0𝑑𝑖,𝑗𝑓𝑥𝑗,,𝑡𝜕𝑚𝑥𝜕𝑥𝑖,𝑡𝑁+1𝑗=0𝑑𝑖,𝑗𝑚𝑥𝑗,𝑡(15) for 𝑓𝑥(𝑡) and 𝑚𝑥(𝑡), respectively.

Since the exact value of (𝜕𝑛/𝜕𝑥)(𝑥0,𝑡) is given by (4), the approximation (13) is not needed at the first point 𝑥0, that is, for the first component of 𝑛𝑥(𝑡). Therefore, from (13), the first approximation in (15) with 𝑓 and 𝑖=0, and from the boundary conditions (4) and (5) we obtain𝑛𝑥(𝑡)𝐷0(1)𝑛𝛾(𝑡)+𝑑𝑛𝑛𝑥0𝑠,𝑡𝑓0(𝑡)𝑒1,(16) where𝐷0(1)=𝑑0001,0𝑑1,1𝑑1,𝑁𝑑𝑁,0𝑑𝑁,1𝑑𝑁,𝑁,𝑒1=100,𝑠𝑓0(𝑡)=𝑁+1𝑗=0𝑑0,𝑗𝑓𝑥𝑗.,𝑡(17) From (13) and (16) we obtain the following approximation for the second-order derivatives𝑛𝑥𝑥(𝑡)𝐷(1)𝑛𝑥(𝑡)+𝑠𝑛𝑁+1(𝑡)𝑤𝐷(1)𝐷0(1)𝛾𝑛(𝑡)+𝑑𝑛𝑛𝑥0𝑠,𝑡𝑓0(𝑡)𝑒1+𝑠𝑛𝑁+1(𝑡)𝑤,(18) where𝐷(1)=𝑑0,0𝑑0,1𝑑0,𝑁𝑑1,0𝑑1,1𝑑1,𝑁𝑑𝑁,0𝑑𝑁,1𝑑𝑁,𝑁𝑑,𝑤=0,𝑁+1𝑑1,𝑁+1𝑑𝑁,𝑁+1,𝑠𝑛𝑁+1(𝑡)=𝑁+1𝑗=0𝑑𝑁+1,𝑗𝑛𝑥𝑗.,𝑡(19)

We now construct approximations to 𝑓𝑥(𝑡), 𝑚𝑥(𝑡), 𝑓𝑥𝑥(𝑡), and 𝑚𝑥𝑥(𝑡). From the spectrally accurate approximations (15) and from (4) and (5) we obtain𝑓𝑥(𝑡)𝐷(1)𝑥𝑓(𝑡)+𝑓𝑁+1𝑚,𝑡𝑤,(20)𝑥(𝑡)𝐷0(1)𝑚(𝑡).(21) According to (20), we have the following approximation for the second-order derivative of 𝑓𝑓𝑥𝑥(𝑡)𝐷(1)𝐷(1)𝑥𝑓(𝑡)+𝑓𝑁+1𝑤,𝑡+𝑠𝑓𝑁+1(𝑡)𝑤,(22) with the notation𝑠𝑓𝑁+1(𝑡)=𝑁+1𝑗=0𝑑𝑁+1,𝑗𝑓𝑥𝑗.,𝑡(23) From (21) we have𝑚𝑥𝑥(𝑡)𝐷(1)𝐷0(1)𝑚(𝑡)+𝑠𝑚𝑁+1(𝑡)𝑤,(24) with the similar notation for 𝑚𝑠𝑚𝑁+1(𝑡)=𝑁+1𝑗=0𝑑𝑁+1,𝑗𝑚𝑥𝑗.,𝑡(25)

We now replace the spatial derivatives in the model (1) by their corresponding approximations. We apply (18), (16), (20), and (22) to the first equation in (1) and obtain its discrete version written in the following form𝑑𝑛𝑑𝑡(𝑡)=𝑑𝑛𝐷(1)𝐷0(1)𝛾𝑛(𝑡)+𝑑𝑛𝑛𝑥0𝑠,𝑡𝑓0(𝑡)𝑒1+𝑑𝑛𝑠𝑛𝑁+1𝐷(𝑡)𝑤𝛾0(1)𝛾𝑛(𝑡)+𝑑𝑛𝑛𝑥0𝑠,𝑡𝑓0(𝑡)𝑒1𝐷(1)𝑓𝑥(𝑡)+𝑓𝑁+1𝑤𝐷,𝑡+𝑛(𝑡)(1)𝐷(1)𝑥𝑓(𝑡)+𝑓𝑁+1𝑤,𝑡+𝑠𝑓𝑁+1(,𝑡)𝑤(26) where stands for the component-wise multiplication between two vectors. The discrete version of the second equation in (1) is written in the form𝑑𝑓𝑑𝑡(𝑡)=𝜂(𝑚(𝑡)𝑓(𝑡))(27) and from (24) we obtain the following discrete form of the third equation in  (1)𝑑𝑚𝑑𝑡(𝑡)=𝑑𝑚𝐷(1)𝐷0(1)𝑚(𝑡)+𝑑𝑚𝑠𝑚𝑁+1(𝑡)𝑤+𝛼𝑛(𝑡)𝛽𝑚(𝑡).(28) The resulting system (26)–(28) is composed of 3𝑁+3 ordinary differential equations and is a semidiscrete version of (1). Note that since the spatial derivatives in (1) are approximated with the spectral accuracy, much smaller numbers of grid-points 𝑥𝑖 are needed for (26)–(28) than for finite difference schemes, and time integration of the smaller systems is more robust and more efficient than time integration of the finite difference systems.

For the models (2) and (3) supplemented with (4) and the right-hand side boundary condition (6), which is different than (5), we need to apply different approximations than (16), (18), (20), and (22) as they include (5) instead of (6). For this problem, from (13), instead of (16), we obtain𝑛𝑥(𝑡)𝐷(1)00𝑛𝛾(𝑡)+𝑑𝑛𝑤(0,𝑁+1),(29) where𝐷(1)00=𝑑0001,0𝑑1,1𝑑1,𝑁𝑑𝑁,0𝑑𝑁,1d𝑁,𝑁,𝑤000(0,𝑁+1)=𝑠𝑓0𝑥(𝑡)𝑛000𝑠,𝑡𝑓𝑁+1(𝑥𝑡)𝑛𝑁+1.,𝑡(30) Further, from (29), we obtain𝑛𝑥𝑥(𝑡)𝐷(1)𝐷(1)00𝛾𝑛(𝑡)+𝑑𝑛𝑤(0,𝑁+1).(31) As in (13) and (15), for the vector𝜕(𝑡)=𝑛𝑥𝜕𝑥0,𝑡𝜕𝑓𝑥𝜕𝑥0𝜕,𝑡𝑛𝑥𝜕𝑥1,𝑡𝜕𝑓𝑥𝜕𝑥1𝜕,𝑡𝑛𝑥𝜕𝑥𝑁,𝑡𝜕𝑓𝑥𝜕𝑥𝑁,𝑡,(32) of approximations to the haptotactic term we obtain(𝑡)𝐷(1)𝑛(𝑡)𝑓𝑥(𝑡)(33) and instead of (24), from (6), we obtain𝑚𝑥𝑥(𝑡)𝐷(1)00𝑚(𝑡).(34)

From (29), (31), (33), and (34) we obtain the following scheme for the problem (3), (4), and  (6)𝑑𝑛𝑑𝑡(𝑡)=𝐷(1)𝑑𝑛𝐷(1)00𝑛(𝑡)+𝛾𝑤(0,𝑁+1)𝛾𝑛(𝑡)𝑓𝑥(𝑡)+𝜇1𝑛(𝑡)(𝑒𝑛(𝑡)𝑓(𝑡)),𝑑𝑓𝑑𝑡(𝑡)=𝜂(𝑚(𝑡)𝑓(𝑡))+𝜇2𝑓(𝑡)(𝑒𝑛(𝑡)𝑓(𝑡)),𝑑𝑚𝑑𝑡(𝑡)=𝑑𝑚𝐷(1)00𝑚(𝑡)+𝛼𝑛(𝑡)(𝑒𝑛(𝑡))𝛽𝑚(𝑡),(35) where 𝑒 is a vector entries of which are all 1-s. For (2), the component 𝛼𝑛(𝑡)(𝑒𝑛(𝑡)) needs to be replaced by 𝛼𝑛(𝑡) in the last equation of (35). From the boundary conditions (10), we obtain the approximation for the diffusion of the endogenous inhibitors𝑢𝑥𝑥(𝑡)𝐷(1)00𝑢(𝑡)(36) and the semi-discrete version for the last equation in (8) is written in the following form𝑑𝑢𝑑𝑡(𝑡)=𝑑𝑢𝐷(1)00𝑢(𝑡)+𝜉𝑓(𝑡)𝜃𝑢(𝑡)𝑚(𝑡)𝜌𝑢(𝑡),(37) where we assume that the inhibitor production is modelled by 𝐹(𝑚,𝑓)=𝜉𝑓. The semi-discrete equations have to be closed by initial conditions chosen according to (7) and (9).

4. Numerical Experiments

We apply the approximations introduced in Section 3 and begin our series of numerical simulations from (26)–(28), which correspond to model (1). Results of our numerical experiments are presented in Figures 16. We use the parameter values 𝑑𝑛=0.001, 𝑑𝑚=0.001, 𝛼=0.1, 𝛽=0, 𝜀=0.01, and different values of 𝛾 and 𝜂 specified in the captions of the figures.

fig1
Figure 1: Tumour cell migration and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (1), (4), (5), (7) with the parameter values 𝛾=0.005 and 𝜂=10.
fig2
Figure 2: Evolution in time: the tumour cell density (a), ECM density (b), and MDE concentration (c). Solutions to (1), (4), (5), (7) with the parameter values as in Figure 1.
fig3
Figure 3: Tumour cell migration and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (1), (4), (5), (7) with the parameter values 𝛾=0.01 and 𝜂=10.
fig4
Figure 4: Evolution in time: the tumour cell density (a), ECM density (b), and MDE concentration (c). Solutions to (1), (4), (5), (7) with the parameter values as in Figure 3.
fig5
Figure 5: Tumour cell migration and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (1), (4), (5), (7) with the parameter values 𝛾=0.02 and 𝜂=20.
fig6
Figure 6: Evolution in time: the tumour cell density (a), ECM density (b), and MDE concentration (c). Solutions to (1), (4), (5), (7) with the parameter values as in Figure 5.

We apply 𝑁=30, that is 32 grid-points 𝑥𝑖, for the numerical experiments presented in Figures 14. The time of integration of the system (26)–(28) based on 32 grid points is 0.22 sec to compute the numerical solutions presented in Figures 1 and 2 and 0.38 sec to compute the numerical solutions from Figures 3 and 4. For Figures 5 and 6, we apply 𝑁=43, that is 45 grid points, and in this case, the time of integration of (26)–(28) is 0.35 sec.

Figures 1, 3, and 5 show snapshots in time and Figures 2, 4, and 6 show continuous evolution in time of tumour cells, ECM, and MDE, and their interactions for all 𝑥 in the space domain. The numerical results presented in Figures 1 and 2 were obtained with 𝛾=0.005 and 𝜂=10, see also [34, Figure 10.2]. The results from Figures 3 and 4 were obtained with 𝛾=0.01 and 𝜂=10.

Two distinct clusters of tumour cells are seen in Figures 1 and 3 at 𝑡=1 and 𝑡=10. The numerical results show that the new clusters, which are not seen at 𝑡=0 and appear at 𝑡=1 and 𝑡=10, are created at the leading edge of the tumour as a result of the diffusion and haptotactic migration modeled by the two components from the right-hand side of the first equation in (1): random motility 𝑑𝑛Δ2𝑛 and haptotaxis 𝛾Δ(𝑛Δ𝑓), respectively. Since 𝛾 is greater for Figure 3 than for Figure 1, because of larger haptotactic migration in Figure 3 than in Figure 1, the two clusters seen in Figure 3 are more separated from each other than the two clusters in Figure 1. The pictures show the effect of haptotaxis. The small clusters of cells, which break away from the main body of the tumour, illustrate the potential for the cancer cells to degrade the surrounding tissue, migrate, and start the metastatic cascade. The migrations of the small clusters may not be detected during the processes of medical treatments, and even after resections of the main tumours, the new small clusters may initiate recurrences of the disease. A new cluster of tumour cells broken away from the main body of the tumour is also observed in Figures 5 and 6, which present numerical data computed with 𝛾=0.02 and 𝜂=20.

The next part of our numerical experiments concerns the models (2) and (3) supplemented by (4), (6), and (7). The results for model (2) are presented in Figures 7, 8, 11, and 12 and for model (3) in Figures 9, 10, 13, and 14. These experiments start from the initial condition (7) corresponding to the snapshot in time 𝑡=0 in Figure 1. We observe that the small clusters of cancer cells separated from the main tumours are better formed at 𝑡=2 than at 𝑡=1 and as time evolves the haptotactic migration together with the production of new cancer cells spread the shapes of the tumours over the 𝑥-domain. We also observe that the snapshots in time 𝑡=1 in Figures 1, 7, and 8 look similar to each other and the models (1), (2), and (3) give similar results for 𝑡[0,1] although they are supplemented by the different boundary conditions, either (5) or (6), and solved with different parameters 𝜇1,𝜇2{0,0.1,0.5} and 𝛽{0,0.07}. However, these similarities are observed only for 𝑡[0,1] and as time evolves the corresponding solutions of the models (1), (2), and (3) differ from each other. For example, already at 𝑡=2, Figure 8 shows greater production of tumour cells than Figure 7. Moreover, at 𝑡=10 and 𝑡=20, Figure 7 shows weaker MDE production and greater production of the tumour cells and the ECM than in Figure 1 due to the fact that 𝜇1, 𝜇2, and 𝛽 are greater in Figure 7 than in Figure 1. On the other hand, also at 𝑡=10 and 𝑡=20, Figure 8 shows greater MDE production than Figure 1. Although 𝛽=0.07 for Figure 8 and 𝛽=0 for Figure 1, since the MDE production is greater in Figure 8 than in Figure 1, the MDE concentration is greater in Figure 8 than in Figure 1 and consequently the ECM degradation is more progressive in Figure 8 than in Figure 1. It can also be observed that although the parameters 𝑑𝑚, 𝛼, and 𝛽 in the third equation of (1) and (2) are the same for Figures 1, 7, and 8, the MDE curves show different MDE concentrations in all of these figures.

fig7
Figure 7: Tumour cell proliferation, migration, ECM re-establishment, and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (2), (4), (6), (7) with the parameter values 𝛾=0.005, 𝜂=10, 𝜇1=0.1, 𝜇2=0.5, and 𝛽=0.07.
fig8
Figure 8: Tumour cell proliferation, migration, ECM re-establishment, and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (2), (4), (6), (7) with the parameter values 𝛾=0.005, 𝜂=10, 𝜇1=0.5, 𝜇2=0.1, and 𝛽=0.07.
fig9
Figure 9: Tumour cell proliferation, migration, ECM re-establishment, and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (3), (4), (6), (7) with the parameter values 𝛾=0.005, 𝜂=10, 𝜇1=0.1, 𝜇2=0.5, and 𝛽=0.07.
fig10
Figure 10: Tumour cell proliferation, migration, ECM re-establishment, and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (3), (4), (6), (7) with the parameter values 𝛾=0.005, 𝜂=10, 𝜇1=0.5, 𝜇2=0.1, and 𝛽=0.07.
fig11
Figure 11: Tumour cell proliferation, migration, ECM re-establishment, and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (2), (4), (6), (7) with the parameter values 𝛾=0.01, 𝜂=10, 𝜇1=0.1, 𝜇2=0.5, and 𝛽=0.07.
fig12
Figure 12: Tumour cell proliferation, migration, ECM re-establishment, and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (2), (4), (6), (7) with the parameter values 𝛾=0.01, 𝜂=10, 𝜇1=0.5, 𝜇2=0.1, and 𝛽=0.07.
fig13
Figure 13: Tumour cell proliferation, migration, ECM re-establishment, and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (3), (4), (6), (7) with the parameter values 𝛾=0.01, 𝜂=10, 𝜇1=0.1, 𝜇2=0.5, and 𝛽=0.07.
fig14
Figure 14: Tumour cell proliferation, migration, ECM re-establishment, and interactions between the tumour and the surrounding tissue: tumour cell density (solid), ECM density (dashed), and MDE concentration (dashdot). Solutions to (3), (4), (6), (7) with the parameter values 𝛾=0.01, 𝜂=10, 𝜇1=0.5, 𝜇2=0.1, and 𝛽=0.07.

In Figures 710, we observe differences due to the MDE production terms 𝛼𝑛 and 𝛼𝑛(1𝑛) in (2) and (3), respectively. The parameter values for Figure 9 are the same as for Figure 7 and the parameter values for Figure 10 are the same as for Figure 8 but the MDE production is weaker in Figure 9 than in Figure 7 and also weaker in Figure 10 than in Figure 8. The MDE concentrations in Figures 9 and 10 are more uniformly spread out across the 𝑥-domain than in Figures 7 and 8. Furthermore, Figures 9 and 10 illustrate that the system (3) models a decreasing MDE production as the tumour invasion progresses and consequently less MDE concentrations in the regions where the high tumour cells densities are situated (possibly in these regions the tumour cells already finalized their invasion due to lack of space and move to other regions) than in the regions of less advanced stage of invasion. This shows that cancer cells may not need MDEs in the regions where they deal with lack of space. Similar features are observed in Figures 1114, where the clusters separating from the main tumours are more profound than in Figures 710. The same parameter values were used for Figure 11 as for Figure 13 and the same parameter values were used for Figure 12 as for Figure 14. Figures 11 and 12 illustrate numerical solutions to (2), while Figures 13 and 14 illustrate numerical solutions to (3). Comparison of Figures 1114 confirms that the term 𝛼𝑛(1𝑛) in (3) models lower MDE production than the term 𝛼𝑛 in (2) (with the same parameter values). We also observe that since MDE production is lower in Figures 13 and 14 than in Figures 11 and 12, the ECM degradation is smaller in Figures 13 and 14 than in Figures 11 and 12.

Figure 15 shows snapshots in time of the four solutions to the more general model (8) describing the interactions between the tumour cells, ECM, MDEs, and endogenous inhibitors. All four profiles show that, as time evolves, the ECM produces endogenous inhibitors, concentration of which increases in time and their higher concentration is located in the regions where ECM is not yet entirely degraded rather than in the regions where the degradation is already effectively developed. The inhibitor profile shows that the ECM responds to the MDEs by producing the endogenous inhibitors.

fig15
Figure 15: Tumour cell density (solid), ECM density (dashed), MDE concentration (dashdot), and endogenous inhibitors (dotted). Solutions to (8) with (4), (6), (7), (9), (10) and the parameter values 𝛾=0.005, 𝜂=10, 𝜇1=0.5, 𝜇2=0.1, 𝛽=0.07, 𝑑𝑢=𝑑𝑛=𝑑𝑚=0.001, 𝜃=0.05, 𝜉=0.03, and 𝜌=0.07.

5. Concluding Remarks and Future Directions

We have constructed a new numerical algorithm for fast computations of the solutions of the mathematical models proposed by Chaplain et al. in [34, 36, 37], which consist of systems of nonlinear partial differential equations describing interactions between tumour cells, ECM, and MDEs. The algorithm is based on spectrally accurate approximations and small amounts of grid points, which results in ordinary differential systems of small dimensions and fast computations. We have applied the algorithm and presented and compared the numerical simulations with a variety of model equations. The simulations demonstrate that the models describe important features of the interactions between tumour cells and the surrounding tissue, and in particular the initiation of a new colony of cells and metastasis.

Our future research work will address the question for which parameter values and domains the model [34] and the kinetic type model proposed in [49] are equivalent. We will also address numerical methods with spectrally accurate approximations for the models with two-dimensional spatial domain and with different kinds of the function 𝐹(𝑚,𝑓) modelling the inhibitor production [34].

Acknowledgment

The authors wish to express their gratitude to anonymous referees for the useful remarks and comments which led to the improvements in the presentation of the paper.

References

  1. World Health Organization Statistical Information System, “WHO database of mortality andburden estimates,” January 2010, http://www3.who.int/whosis/mort.
  2. UN Database of Population Projections, January 2010, http://esa.un.org/unpp.
  3. H. M. Byrne, M. R. Owen, T. Alarcon, J. Murphy, and P. K. Maini, “Modelling the response of vascular tumours to chemotherapy: a multiscale approach,” Mathematical Models and Methods in Applied Sciences, vol. 16, no. 1, pp. 1219–1241, 2006. View at Publisher · View at Google Scholar
  4. A. d'Onofrio, “Metamodeling tumor-immune system interaction, tumor evasion and immunotherapy,” Mathematical and Computer Modelling, vol. 47, no. 5-6, pp. 614–637, 2008. View at Publisher · View at Google Scholar
  5. M. J. Quinn, A. d'Onofrio, B. Møller et al., “Cancer mortality trends in the EU and acceding countries up to 2015,” Annals of Oncology, vol. 14, no. 7, pp. 1148–1152, 2003. View at Publisher · View at Google Scholar
  6. H. M. Byrne, M. A. J. Chaplain, G. J. Pettet, and D. L. S. McElwain, “An analysis of a mathematical model of trophoblast invasion,” Applied Mathematics Letters, vol. 14, no. 8, pp. 1005–1010, 2001. View at Publisher · View at Google Scholar
  7. V. T. De Vito Jr., J. Hellman, and S. A. Rosenberg, Eds., Cancer: Principles and Practice of Oncology, Lippincott Williams and Wilkins, Philadelphia, Pa, USA, 2005.
  8. M. Pekham, H. M. Pinedo, and U. Veronesi, Oxford Textbook of Oncology, Oxford Medical Publications, Oxford, UK, 1995.
  9. N. Bellomo and A. Bellouquid, “On the onset of non-linearity for diffusion models of binary mixtures of biological materials by asymptotic analysis,” International Journal of Non-Linear Mechanics, vol. 41, no. 2, pp. 281–293, 2006. View at Publisher · View at Google Scholar · View at MathSciNet
  10. A. D'Onofrio, “Tumor-immune system interaction: modeling the tumor-stimulated proliferation of effectors and immunotherapy,” Mathematical Models and Methods in Applied Sciences, vol. 16, no. 8, pp. 1375–1401, 2006. View at Publisher · View at Google Scholar
  11. P. L. Lollini, S. Motta, and F. Pappalardo, “Modeling tumor immunology,” Mathematical Models and Methods in Applied Sciences, vol. 16, no. 1, pp. 1219–1242, 2006. View at Publisher · View at Google Scholar
  12. R. A. Gatenby, T. L. Vincent, and R. J. Gillies, “Evolutionary dynamics in carcinogenesis,” Mathematical Models and Methods in Applied Sciences, vol. 15, no. 11, pp. 1619–1638, 2005. View at Publisher · View at Google Scholar
  13. L. A. Liotta, C. N. Rao, and S. H. Barsky, “Tumor invasion and the extracellular matrix,” Laboratory Investigation, vol. 49, no. 6, pp. 636–649, 1983. View at Google Scholar
  14. W. G. Stetler-Stevenson, S. Aznavoorian, and L. A. Liotta, “Tumor cell interactions with the extracellular matrix during invasion and metastasis,” Annual Review of Cell Biology, vol. 9, pp. 541–573, 1993. View at Google Scholar
  15. T. E. Wheldon, Mathematical Models in Cancer Research, Hilger, Boston, Mass, USA, 1988.
  16. J. L. Holleran, C. J. Miller, N. L. Edgehouse, T. P. Pretlow, and L. A. Culp, “Differential experimental micrometastasis to lung, liver, and bone with lacZ-tagged CWR22R prostate carcinoma cells,” Clinical and Experimental Metastasis, vol. 19, no. 1, pp. 17–24, 2002. View at Publisher · View at Google Scholar
  17. J. A. Lawrence, P. S. Steeg, and N. C. Institute, “Mechanisms of tumor invasion and metastasis,” World Journal of Urology, vol. 14, no. 3, pp. 124–130, 1996. View at Google Scholar
  18. N. V. Mantzaris, S. Webb, and H. G. Othmer, “Mathematical modeling of tumor-induced angiogenesis,” Journal of Mathematical Biology, vol. 49, no. 2, pp. 111–187, 2004. View at Publisher · View at Google Scholar · View at PubMed · View at MathSciNet
  19. J. C. Ulchaker and E. A. Klein, “Biology of metastasis and its clinical implications: renal-cell cancer,” World Journal of Urology, vol. 14, no. 3, pp. 175–181, 1996. View at Google Scholar
  20. B. Asquith and C. R. M. Bangham, “An introduction to lymphocyte and viral dynamics: The power and limitations of mathematical analysis,” Proceedings of the Royal Society B, vol. 270, no. 1525, pp. 1651–1657, 2003. View at Publisher · View at Google Scholar · View at PubMed
  21. N. Bellomo and L. Preziosi, “Modelling and mathematical problems related to tumor evolution and its interaction with the immune system,” Mathematical and Computer Modelling, vol. 32, no. 3-4, pp. 413–452, 2000. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  22. H. M. Byrne, M. A. J. Chaplain, G. J. Pettet, and D. L. S. McElwain, “A mathematical model of trophoblast invasion,” Journal of Theoretical Medicine, vol. 1, pp. 275–286, 1998. View at Google Scholar
  23. N. Bellomo and E. De Angelis, “Modelling and simulation of tumour development, treatmentand control,” Mathematical and Computer Modelling, vol. 37, no. 11, 2003. View at Google Scholar
  24. N. Bellomo and E. De Angelis, “Modelling complex systems in molecular biology and tumordynamics and control,” Mathematical and Computer Modelling, vol. 41, no. 10, 2005. View at Google Scholar
  25. N. Bellomo and P. K. Maini, “Cancer modelling (I),” Mathematical Models and Methods in Applied Sciences, vol. 15, no. 11, 2005. View at Google Scholar
  26. N. Bellomo and P. K. Maini, “Cancer modelling (II),” Mathematical Models and Methods in Applied Sciences, vol. 16, no. 7b, 2006. View at Google Scholar
  27. N. Bellomo and B. Sleeman, “Multiscale cancer modelling,” Computational and Mathematical Methods in Medicine, vol. 20, no. 2-3, 2006. View at Google Scholar
  28. R. A. Gatenby, “Models of tumor-host interaction as competing populations: implications for tumor biology and treatment,” Journal of Theoretical Biology, vol. 176, no. 4, pp. 447–455, 1995. View at Publisher · View at Google Scholar · View at PubMed
  29. R. A. Gatenby and E. T. Gawlinski, “A reaction-diffusion model of cancer invasion,” Cancer Research, vol. 56, no. 24, pp. 5745–5753, 1996. View at Google Scholar
  30. A. R. A. Anderson, “A hybrid mathematical model of solid tumour invasion: the importance of cell adhesion,” Mathematical Medicine and Biology, vol. 16, pp. 1219–1242, 2005. View at Google Scholar
  31. A. R. A. Anderson, M. A. J. Chaplain, E. L. Newman, R. J. C. Steele, and A. M. Thompson, “Mathematical modelling of tumour invasion and metastasis,” Journal of Theoretical Medicine, vol. 2, pp. 129–154, 2000. View at Google Scholar
  32. P. Boyle, P. Autier, H. Bartelink et al., “European code against cancer and scientific justification: third version,” Annals of Oncology, vol. 14, no. 7, pp. 973–1005, 2003. View at Publisher · View at Google Scholar
  33. P. Boyle, A. d'Onofrio, P. Maisonneuve et al., “Measuring progress against cancer in Europe: has the 15% decline targeted for 2000 come about?” Annals of Oncology, vol. 14, no. 8, pp. 1312–1325, 2003. View at Publisher · View at Google Scholar
  34. M. A. J. Chaplain and A. R. A. Anderson, “Mathematical modelling of tissue invasion,” in Cancer Modelling and Simulation, L. Preziosi, Ed., pp. 269–297, Chapman & Hall/CRC, Boca Raton, Fla, USA, 2003. View at Google Scholar
  35. M. A. J. Chaplain and G. Lolas, “Mathematical modelling of cancer cell invasion of tissue: the role of the urokinase plasminogen activation system,” Mathematical Models and Methods in Applied Sciences, vol. 15, no. 11, pp. 1685–1734, 2005. View at Publisher · View at Google Scholar
  36. M. A. J. Chaplain and G. Lolas, “Mathematical modelling of cancer invasion of tissue: dynamic heterogeneity,” Networks and Heterogeneous Media, vol. 1, no. 3, pp. 399–439, 2006. View at Google Scholar
  37. G. Lolas, Mathematical modelling of the urokinase plasminogen activation system and its rolein cancer invasion of tissue, Ph.D. thesis, University of Dundee, 2003.
  38. B. P. Marchant, J. Norbury, and A. J. Perumpanani, “Travelling shocks waves arising in a model of malignant invasion,” SIAM Journal on Applied Mathematics, vol. 60, pp. 463–476, 2000. View at Google Scholar
  39. M. E. Orme and M. A. J. Chaplain, “A mathematical model of vascular tumour growth and invasion,” Mathematical and Computer Modelling, vol. 23, no. 10, pp. 43–60, 1996. View at Google Scholar
  40. A. J. Perumpanani and H. M. Byrne, “Extracellular matrix concentration exerts selection pressure on invasive cells,” European Journal of Cancer, vol. 35, no. 8, pp. 1274–1280, 1999. View at Publisher · View at Google Scholar
  41. A. J. Perumpanani, J. A. Sherratt, J. Norbury, and H. M. Byrne, “Biological inferences from a mathematical model for malignant invasion,” Invasion and Metastasis, vol. 16, no. 4-5, pp. 209–221, 1996. View at Google Scholar
  42. A. J. Perumpanani, D. L. Simmons, A. J. H. Gearing et al., “Extracellular matrix-mediated chemotaxis can impede cell migration,” Proceedings of the Royal Society B, vol. 265, no. 1413, pp. 2347–2352, 1998. View at Publisher · View at Google Scholar
  43. A. J. Perumpanani, J. A. Sherratt, J. Norbury, and H. M. Byrne, “A two parameter family of travelling waves with a singular barrier arising from the modelling of extracellular matrix mediated cellular invasion,” Physica D, vol. 126, no. 3-4, pp. 145–159, 1999. View at Google Scholar
  44. M. Kolev and B. Zubik-Kowal, “Numerical vs. experimental data for prostate tumour growth,” to appear in Journal of Biological Systems.
  45. A. Calvo, N. Xiao, J. Kang et al., “Alterations in gene expression profiles during prostate cancer progression: functional correlations to tumorigenicity and down-regulation of selenoprotein-P in mouse and human tumors,” Cancer Research, vol. 62, no. 18, pp. 5325–5335, 2002. View at Google Scholar
  46. A. Estreicher, J. Muhlhauser, J. -L. Carpentier, L. Orci, and J. -D. Vassalli, “The receptor for urokinase type plasminogen activator polarizes expression of the protease to the leading edge of migrating monocytes and promotes degradation of enzyme inhibitor complexes,” Journal of Cell Biology, vol. 111, no. 2, pp. 783–792, 1990. View at Publisher · View at Google Scholar
  47. C. Canuto, M. Y. Hussaini, A. Quarteroni, and T. A. Zang, Spectral Methods in Fluid Dynamics, Springer, New York, NY, USA, 1988.
  48. L. N. Trefethen, Spectral Methods in Matlab, SIAM, Philadelphia, Pa, USA, 2000.
  49. M. Lachowicz, “Micro and meso scales of description corresponding to a model of tissue invasion by solid tumours,” Mathematical Models and Methods in Applied Sciences, vol. 15, no. 11, pp. 1667–1683, 2005. View at Publisher · View at Google Scholar