Immersed Boundary Cell model (IBCell).

Principal modeler: Dr. Katarzyna Rejniak.

The IBCell models cells (normal, tumor or stromal) as individual, fully deformable bodies, each consisting of an elastic plasma membrane enclosing the cell cytoplasm (modeled as a viscous incompressible fluid) and the cell cytoskeleton (represented by a mesh of linear springs). All points on the cell membrane play an additional role of either adherens junctions that provide connections and communication between neighboring cells or various cell surface sensors that are used to recognize external chemical cues such as oxygen concentration or extracellular matrix (ECM) density, as well as the presence of other cells in its vicinity.

In the model, initiation and progression of all cell life processes, such as cell growth, division, polarization, migration or death, depend solely on cues sensed by the cell through its membrane sensors. These processes are regulated independently for each individual cell, however, the behavior of each cell can be influenced by actions of other cells, (e.g., growth of a host cell can be suppressed by contact inhibition exerted by its neighbors).

IBCell can be applied to investigate actions of individual cells (such as cell movement or response to treatment) or cell competitive or collaborative behavior in cell populations that may be arranged into topologically different morphologies, such as circular or linear cross-sections through epithelial ducts, self-organized  developmental structures,such as epithelial acini; or unorganized masses of heterogeneous cells, such as invasive finger-like tumor cell cohorts.

IBCell has been used to model natural emergence of an acinar structure from a single cell and to construct charts of normal acini and mutants. By varying thresholds for cell growth and death sensors that define initiation of cell proliferation and apoptosis, respectively, we have determined the parameter space in which either normal hollow acini, mutant filled structures, or small degenerate clusters form.

Several applications of IBCell have already been published, especially in relation to early carcinoma development, including the mapping of the MCF10A non-tumorigenic breast cell line and its specific tumorigenic mutant MCF10A-HER2^[YVMA]. By investigating the 3D space of growth, death and ECM sensor thresholds, we have identified four regions of distinct final acinar morphologies including normal and mutated acini.

A simplified version of the IBCell model--the Point-IBCell model--has been introduced, in which each cell is represented just by a single point corresponding to the center of a cell’s nucleus, without explicitly modeling cell membranes, but preserving cell-cell adhesion relation by introducing a linear spring between each pair of neighboring cells. This model allows a large numbers of cells in both 2D
and 3D space, without the need for expensive computational power.

 

3D space of final acinar morphologies obtained by varying growth, death and ECM sensor thresholds identifies four different regions: normal hollow acini (red), partially or fully filled acini (blue), small degenerate acini (yellow), and non-stabilized growing acini.

 

Cellular Potts Model (CPM).

Principal modelers: Dr. Sandy Anderson.

The CPM model, also know as the Glazier-Grainer-Hogweg model (GGH), is a spatial
grid-based formalism where a cell is defined over a region composed of multiple lattice sites, with constraints acting on its area while experiencing interactions at its boundary.

While a lattice is still used, we included CPM in the Off-Lattice category because the grid does not restrict cell movements with respect to each other. On the contrary, it provides enhanced spatial detail of cell-to-cell margins, which is why we use it. CPM is also a hybrid model, since the tumor microenvironment is represented with continuous fields obeying appropriate partial differential equations (PDEs). CPM model dynamics are based on the free energy minimization principle, and generated by means of Monte Carlo simulations utilizing a Metropolis algorithm. Effectively, this means that cell motion comes about from the overall minimization of the energy of deformation and stretching of the membrane through stochastic fluctuations, in which the global and local forces upon a cell edge are resolved.

We use CPM to extend the outcome of other models to a much larger tumor, i.e., approaching the size of experimentally observable tumors, in the range of a few cubic centimeters (this corresponds to 1-10 billion cells, which on a 2D section reduces to a few million). This is possible because CPM is perhaps the least computationally intensive model we have, which can still capture tumor morphology with single cell resolution, and therefore connect with experimental observations at the histology level.

Additionally, CPM can readily incorporate details of intracellular dynamics (such as signaling pathways via systems of ordinary differential equations (ODEs)). This is one route by which we are able to connect our efforts at the cell scale with experiments at the molecular scale, e.g., performed at other ICBP centers or in the cancer research community at large.

 

Hybrid Discrete-Continuum Model (HDC)/ Hybrid Cellular Automata Model (HCA).

Principal modeler: Dr. Sandy Anderson.

The HDC model of tumor invasion couples a continuum deterministic model of microenvironmental variables (based on a system of reaction-diffusion equations) and a discrete cellular automata-like model of individual tumor cell migration and interaction (based on a biased random-walk model).

In the model, individual cells are defined via collection of phenotypic traits including proliferation, death, cell-cell adhesion, mutation, and production/degradation of microenvironment (mE) specific components that determine how it behaves and interacts with other cells and its mE. The mE is defined by a system of continuous partially differential equations (PDEs) and consists of concentrations of extracellular matrix (ECM), oxygen, and matrix degrading enzymes (MDE). Further, in the model, cells migrate, proliferate, die and potentially mutate depending on their phenotype (we initially define a pool of 100 possible phenotypes) and the mE. If a cell does mutate during division, it is randomly assigned a different phenotype from the pool.

We have previously used the HDC model to investigate the role of both tumor (in terms of a diversity of phenotypes) and mE heterogeneity (in terms of various ECMs and nutrient conditions) using phenotypes based on literature-derived parameter estimates. Recently, we further developed the model to incorporate data from a parameterized collection of genetically-related breast cell lines and examined their sensitivity to different mE constraints (growth factor levels and ECM constraints) without mutation.

By considering different mE conditions (uniform ECM, grainy ECM, and high/ low oxygen), we showed that in “harsh” mEs tumor morphology was distinctly different, exhibiting more fingering and asymmetry. More importantly, these fingered tumors were composed of fewer more aggressive phenotypes, defined in terms of low cell-cell adhesion, short proliferation age, and high migration coefficients. We are further developing the HDC model to incorporate trait distribution data, through five processes:

1. proliferation,
2. inheritance,
3. metabolism,
4. migration,
5. death.

Tumor growth in three different mE: (A) uniform ECM, (B) Grainy ECM and (C) Low nutrient. Upper row: tumor cell distributuions after 3 months of simulated growth, we can see that the three different mE have produced distinct tumour morphologies. Lower row: the relative abundance of the 100 tumor phenotypes as the tumor invaded each of the different mE. We note that there are approximately 6 dominant phenotypes in the uniform tumor, 2 in the grainy and 3 in the low nutrient tumor. In each tumor, one of the phenotypes is the most aggressive and also the most abundant, particularly in B and C.

 

Evolutionary Hybrid Cellular Automata Model (EHCA). The EHCA model uses a cellular automata-based approach that considers cells as simple grid points, similar to HDC. However, unlike the predefined set of traits of HDC, each cell contains a complex neural network that maps genotype to phenotype. Therefore, in the EHCA model phenotypic outcomes emerge naturally as a consequence of evolutionary dynamics, i.e., they are potentially open-ended and certainly not predefined. The grid itself represents the mE (oxygen, glucose and hydrogen ion levels, ECM density).

Cell core traits (e.g., motility, metabolism, proliferation, and death) are controlled by the neural network within each cell, via a feedforward decision mechanism. When a cell divides, the network parameters are copied to its daughter cells with a degree of variation (that can be parameterized), thus introducing heterogeneity and allowing for mE selection. We have previously applied this modeling framework to examine evolutionary dynamics and tumor morphological changes induced by oxygen deprived environments, as well as to investigate the emergence of a glycolytic phenotype in a spatially constrained and oxygen poor mE. Most recently, we used the EHCA model to examine the relationship between proliferation and migration in an evolving tumor. In summary, the EHCA model defines the fitness landscape onto which the tumor progresses, by tracking the feedback between evolving cancer cell population and the mE.

 

Partial Differential Equation Models (PDEs).

Principal modeler: Dr. Glenn Webb, Dr. Sandy Anderson.

PDEs underly the dynamics of most of the microenvironment (mE) factors in the suite of models discussed here, but they can also be used to predict the growth of the tumor itself. They can incorporate global processes with much greater ease than individual-based models.

As a consequence to treating a tumor as a single dynamic entity, PDEs can capture mechanical interactions more readily. They allow for the simulation of tumors on much larger spatial (10 cm) and temporal (decades) scales but are cumbersome at considering many interacting populations an tend to homogenize cell scale dynamics in favor of population scale dynamics.

We have previously developed purely continuous models of tumor invasion that examined the role of migration and the extracellular matrix (ECM) heterogeneity. We also examined utility of different radiotherapy treatments using a PDE approach.

 

A The age-dependent probability of division [phi] (solid blue line) and the corresponding division rate [beta] (dashed green line) are shown. [phi] is given by a shifted [gamma]-distribution with mean 19 h. Observe that the units of a are hours. However, cytostatic effects result in a slowed progression through the cell cycle. B We show the cell-phase-dependent tendency of cells to become nonproliferating [mu](a). Figure from Theor Biol Med Model. 2007; 4: 14.

 

Evolutionary Game Theory (EGT) Models.

Principal modeler: Dr. David Basanta.

EGT provides a conceptual framework for modeling the ecology of many aspects of cancer biology including carcinogenesis, tumor invasion and metastases, and response to therapy. It views cancers and recipient communities of normal tissue as players with strategies. Cellular populations are considered the players in this game. Heritable traits whether genetic or epigenetic are the strategies and can be represented as either continuous or discrete phenotypes. Per capita growth rates are the payoffs. Environmental and cellular interactions determine the selection forces that define the fitness of any cellular phenotype. Importantly, EGT does not make any assumption about the genetic mutations required to obtain the phenotypic strategies, only that those phenotypes are possible and relevant. In the continuous approach, each population is defined by a strategy vector u and a strategy space U which represents its phenotypic potential. This approach uses a Fitness Generating Function. The Darwinian “struggle for existence” is modeled by linking the fitness function to proliferation

In the discrete approach, the interactions between different phenotypes is instead described by a payoff table. This defines how the fitness of each phenotype changes as a result of interaction. To determine the payoff of phenotype 1 interacting with another phenotype, the table is obtained by reading across the 1st row. The same is true for phenotype 2, except we read across the 2nd row (e.g., a type 1 cell interacting with a type 2 will get a payoff fitness of b, conversely that type 2 cell will get a payoff of c. This mechanism can be generalized in games that consist of many more phenotypes.

An important concept in EGT is that of an Evolutionarily Stable Strategy (ESS). A population that has evolved an ESS cannot be invaded by a population with any available strategy. By contrast, a population at an ESS can generally invade and destroy any other population regardless of its strategy. In general, cancer populations will evolve toward an ESS within the tissue they are invading. These evolutionary dynamics of cancer progression and invasion involve competition among tumor subpopulations and between normal and malignant cells.

In previous research, we showed how EGT can be used to study the different evolutionary dynamics in tumors made of two and three phenotypes with different levels of malignancy under a number of scenarios. EGT can also be used to study the interactions of glycolytic cells with exclusively proliferative cells, with relation to the cost of having a glycolytic metabolism. With such approach, it was shown that cancer is an inevitable outcome of the fact that healthy tissues must accommodate new emerging phenotypes as a result of development, and that this can be exploited by tumor cells evolving towards ESS. We further showed that, starting from similar assumptions in terms of phenotypic strategies and the costs and benefits of nutrients and motility, the resulting steady state using EGT is similar to a Cellular Automaton model, without spatial considerations.

In our Center, both EGT approaches will be used to study the role of heterogeneity in driving tumor evolution, as they provide a framework to study dynamics of interactions between cells.

 

A typical CPM cell-lattice configuration showing portions of nine Generalized Cells. Each generalized Cell is a collection of lattice sites (squares) with the same index value. The colors indicate cell types. The number of lattice sites in a cell is its volume and the number of links along its boundary (interfaces with sites containing other indices) is its surface area.