Evolution of cell motility in an individual-based model of tumour growth

Abstract

Tumour invasion is driven by proliferation and importantly migration into the surrounding tissue. Cancer cell motility is also critical in the formation of metastases and is therefore a fundamental issue in cancer research. In this paper we investigate the emergence of cancer cell motility in an evolving tumour population using an individual-based modelling approach. In this model of tumour growth each cell is equipped with a micro-environment response network that determines the behaviour or phenotype of the cell based on the local environment. The response network is modelled using a feed-forward neural network, which is subject to mutations when the cells divide. With this model we have investigated the impact of the micro-environment on the emergence of a motile invasive phenotype. The results show that when a motile phenotype emerges the dynamics of the model are radically changed and we observe faster growing tumours exhibiting diffuse morphologies. Further we observe that the emergence of a motile subclone can occur in a wide range of micro-environmental growth conditions. Iterated simulations showed that in identical growth conditions the evolutionary dynamics either converge to a proliferating or migratory phenotype, which suggests that the introduction of cell motility into the model changes the shape of fitness landscape on which the cancer cell population evolves and that it now contains several local maxima. This could have important implications for cancer treatments which focus on the gene level, as our results show that several distinct genotypes and critically distinct phenotypes can emerge and become dominant in the same micro-environment.

Fig. 1

The layout of the response network which determines the phenotype or behaviour of the cancer cells. The micro-environmental variables are fed into the input layer of the network which then processes the information and produces a response at the output layer. The response depends on the network matrices w and W, which are subject to mutations when the cells divide.

Fig. 2

The life-cycle of the cancer cells represented as a flowchart.

Fig. 3

The life-cycle response of the initial cell as a function of the number of neighbours, oxygen concentration and ECM gradient. The glucose and hydrogen ion concentrations are kept at to their background values.

Fig. 4

Spatial distribution of the cells at t = 30, 60 and 90 (approx. 20, 40 and 60 days) for c₀ = 0.5 and E = 0.1 on a grid of size 200×200. Proliferating cells are shown as red, quiescent cells as green, necrotic cells as yellow, moving cells as cyan, dead cells as blue and empty grid points are white. In this simulation a haptotactic subclone emerges on the left side of the tumour leading to an asymmetric morphology. A supplementary movie of this simulation can be found at: http://www.nbi.dk/~gerlee/Fig4.mov.

Fig. 5

Spatial distribution of the cells at t = 40, 80 and 120 (approx. 27, 53 and 80 days) for c₀ = 0.1 and E = 0.3 on a grid of size 200×200. Proliferating cells are shown as red, quiescent cells as green, necrotic cells as yellow, moving cells as cyan, dead cells as blue and empty grid points are white. Haptotaxis in low oxygen concentration gives rise to a branched tumour morphology although the branches are not as well defined compared to the lower part of the tumour which is dominated by immobile cells. A supplementary movie of this simulation can be found at: http://www.nbi.dk/~gerlee/Fig5.mov

Fig. 6

The time evolution of the (a) invasive distance and (b) number of cells for the simulations shown in fig. 4 and 5.

Fig. 7

The time evolution of the average response vector for (a) the simulation shown in fig. 4 and (b) fig. 5. Panel (c) and (d) show the result of two simulations that occurred in identical micro-environmental conditions as (a) and (b) respectively, but instead evolved towards a proliferating population.

Fig. 8

The average movement potential S₄ as a function of the oxygen concentration c₀ and the matrix density E.

Fig. 9

The probability p_gl that the glycolytic phenotype dominates the population. Note that the general shape of the surface is similar to the results obtained in Gerlee and Anderson (2008) where the cells were immobile, but that p_gl in this case is smaller.

Fig. 10

Schematic of the reseed-simulations. The dominant sublcone (circled) at the end of the first simulation is used as the initial cell in the next step of the simulation. This process is then repeated n times, and the resulting phenotype is usually highly evolved and aggressive.

Fig. 11

The time evolution of the (a) the phylogenetic depth and (b) the proliferation and movement potentials for two different reseed-experiments in the same micro-environment (c₀, E) = (0.1, 0.35). In both cases the dominant subclone was reseeded 10 times, but the smaller number of total time steps in one of the simulations (dashed lines) implies that the growth rate in that case was higher.

Fig. 12

The life-cycle response of three evolved genotypes as a function of the number of neighbours, oxygen concentration and ECM gradient. (a) and (b) utilise the ECM as a switch between proliferation (red) and haptotaxis (cyan), while in (c) the oxygen concentration also influences the choice.