A Cell-Based Model of Extracellular-Matrix-Guided Endothelial Cell Migration During Angiogenesis

Angiogenesis, the formation of new blood vessels sprouting from existing ones, occurs in several situations like wound healing, tissue remodeling, and near growing tumors. Under hypoxic conditions, tumor cells secrete growth factors, including VEGF. VEGF activates endothelial cells (ECs) in nearby vessels, leading to the migration of ECs out of the vessel and the formation of growing sprouts. A key process in angiogenesis is cellular self-organization, and previous modeling studies have identified mechanisms for producing networks and sprouts. Most theoretical studies of cellular self-organization during angiogenesis have ignored the interactions of ECs with the extra-cellular matrix (ECM), the jelly or hard materials that cells live in. Apart from providing structural support to cells, the ECM may play a key role in the coordination of cellular motility during angiogenesis. For example, by modifying the ECM, ECs can affect the motility of other ECs, long after they have left. Here, we present an explorative study of the cellular self-organization resulting from such ECM-coordinated cell migration. We show that a set of biologically-motivated, cell behavioral rules, including chemotaxis, haptotaxis, haptokinesis, and ECM-guided proliferation suffice for forming sprouts and branching vascular trees. Electronic Supplementary Material The online version of this article (doi:10.1007/s11538-013-9826-5) contains supplementary material, which is available to authorized users.


Background
Angiogenesis, the formation of new blood vessels sprouting from existing vessels, is an important process during development, reproduction and tissue repair. However, angiogenesis can also be a pathological process. For example, it is required for tumors to sustain their growth. A lot of research has been done to get more insight in this process, to find the major regulators of angiogenesis and to learn about the mechanisms involved in the initiation, elongation and branching of new blood vessel sprouts. This knowledge could for instance help finding new therapies against cancer [1][2][3].

Scope and Objective
In this thesis we will discuss what is known today about angiogenesis, in particular the process of blood vessel formation induced by growing tumors. We will focus on the chemical and mechanical mechanisms that regulate the directional migration and proliferation of endothelial cells, the cells that form a new growing blood vessel sprout.
In particular we are interested in the interactions between endothelial cells and the extracellular matrix. How do cells respond to the composition of the matrix and in what way does this influence the growth and branching of a new sprout? Can these interactions explain characteristics of sprouting angiogenesis?
We developed a cell based model to find answers to these questions. Cell-based models describe cell behaviors, including cell-cell interactions and the interactions of individual cells with their micro-environment, in terms of simple sets of rules. Thus cell-based models can give new insights in the mechanisms which are important for sprout formation and branching [4]. 2

Organization of rest of thesis
In the following chapter we will give an overview of the processes involved in sprouting angiogenesis. We will focus on the interactions between cells and the extracellular matrix. In chapter 3 we will discuss and compare existing models of angiogenesis. Next we present our cell-based model of angiogenesis and describe the process of translating experimental observations and biological theories into a model.
In the final chapters we will present the simulation results, discuss our model and give suggestions for improvements.

Cell-matrix interactions during angiogenesis 2.1 Hypoxic tumors induce angiogenesis
Small tumors up to a size of ~1 mm can absorb, by simple diffusion, sufficient oxygen and nutrients from their direct environment. When a tumor grows beyond this size tumor cells will lack oxygen and become hypoxic. This turns on the so-called 'angiogenic switch', which leads to an increased expression of several angiogenic factors [5][6].
One of the key angiogenic factors is vascular endothelial growth factor (VEGF). With VEGF we usually mean VEGF-A, which is a member of the VEGF family. VEGF-A has several isoforms, resulting from alternative splicing. A number of these isoforms have a heparin-binding domain and bind to the extracellular matrix (ECM), others are soluble factors [6].
VEGF released by a tumor diffuses into the surrounding tissue, establishing a chemical gradient between the tumor and nearby blood vessels. When it reaches a vessel, it binds to cell surface receptors on endothelial cells (ECs) which form the inner lining of blood vessel walls. This activates the ECs, resulting in increased cell survival, migration and proliferation [7].
ECs activated by VEGF first degrade the basement membrane of the parent vessel and then migrate into the extracellular matrix (ECM) towards the tumor. First, small sprouts are formed by aggregation and migration of ECs that are recruited from the parent vessel. The sprout will further extend when some of the ECs in the sprout wall begin to divide [5,8]. During this process the vessel sprout will form new branches, and these branches can reconnect again, a process called anastomosis. Recently, Fantin and co-workers reported that macrophages can act as 'bridge-cells' in this process, mediating the fusion of the tips of two sprouts [9]. Finally the vascular tree reaches 4 the tumor, providing it with nutrients and oxygen. Once vascularized, a tumor is more likely to become malignant and to spread and metastasize to other parts of the body [10].

Endothelial cells take on different roles
During angiogenesis ECs can have different phenotypes [11][12]. Cells at the tip of the sprout, tip cells, lead the new formed vessel. Their task is to navigate and they actively extend filopodia, which are spiky protrusions, to sense and respond to guidance cues in their environment. They can secrete proteases to degrade the extracellular matrix.
Tip cells do not proliferate (grow and divide).
Stalk cells are cells trailing the tip cells. They are less motile and they barely extend filopodia. Stalk cells proliferate when stimulated with VEGF [13] and they deposit ECM components. Their task is to elongating the stalk, to form lumen and connect to the circulation.
Once the vessel has formed, cells turn into the phalanx phenotype; they become quiescent and mostly stop dividing. The key function of the vessel is then to supply blood and oxygen to tissues [12].
Tip cell selection and induction is regulated with the endothelial DLL4/NOTCH pathway, which is stimulated by VEGF. This signaling pathway inhibits tip cell formation near other tip cells through lateral inhibition [11].

Proliferation is required to reach the tumor
Sprouting is possible without proliferation of stalk cells. However proliferation is necessary to sustain sprouting for a longer period and to grow a large enough sprout that can reach the tumor. [7,10] Although the general idea is that proliferation occurs behind the tip cell, there is no consensus about the exact location of EC mitosis during angiogenesis. Experiments have shown that proliferation can occur some distance behind the sprout tip [10,14], at the base of a new sprout [10,15], and even at the tip of the sprout [8,[15][16]. 5 Several studies suggest that proliferation only occurs when the connection between adjacent cells has been disrupted. Therefore it is possible that during angiogenesis mitosis occurs as a result of gaps between ECs in the new sprout [8].

The extracellular matrix has many roles in angiogenesis
The extracellular matrix is a mesh-like network of macromolecules secreted locally by cells. The main components of the ECM are proteoglycans, fibrous proteins such as collagen and elastin and adhesive proteins like laminin and fibronectin. In vertebrates the ECM constitutes the major part of the connective tissue and it forms the basement membrane (or basal lamina), which is a thin sheet of ECM underlying epithelial cells [17].
The ECM has many roles in angiogenesis. It is essential for EC migration, proliferation and survival, since it provides structural support and chemical cues for cell adhesion and motility. ECM components like collagen I and fibrin are capable of supporting chemotactic migration. The density and spatial distribution of ECM proteins such as fibronectin and collagen can affect the speed and direction of cell migration. Furthermore, ECs are able to secrete and degrade ECM components [18].

Cells migrate by attachment to the ECM
Many studies of cell migration describe in vitro experiments where cells are plated on a dish coated with ECM components. In other experiments, cells are seeded in a three dimensional matrix, in which they show distinct migrating behavior.
The following steps describe how cells move on two-dimensional substrata. First they extend filopodia to sense signals like VEGF gradients in their direct environment.
Then sheet-like extensions, called lamellipodia, will form dynamic attachments to the ECM at the front of the cell. Next the stress fibers within the cell contract, which results in a detachment of the rear of the cell. Finally the cytoskeleton relaxes, adhesive and signaling components are recycled and the cell repeats the cycle. [19].

6
In three-dimensional situations cells behave differently. They have a more elongated and spindle-like shape due to the physical restriction of matrix fibrils. Instead of extending lamellipodia, the cells form pseudopodia following the direction of matrix fibrils. Migrating cells in a 3D matrix have less stress fibers, focal adhesions and spreading. Because the 3D matrix forms a physical barrier around the cells, proteolysis of ECM by cells is necessary for motility [20].
Cells attach to the ECM by reversibly binding transmembrane receptors, mostly integrins, to ECM proteins. The integrin family includes more than 20 members. They bind their extracellular domain to specific ligands such as fibronectin, collagen and laminin and cluster in the membrane to form adhesive contacts called focal adhesions.
During EC migration focal adhesions and stress fibers are aligned in the direction of movement resulting in polarized cells [19].
Cell-ECM adhesions regulate cell migration in two ways. They have an adhesive function, binding the extracellular matrix to the actin cytoskeleton, and a signal transduction function, regulating molecules important for cell motility [20].
VEGF can stimulate cell migration in several manners. It increases the expression and activation of several integrins involved in angiogenesis. Cell-cell adhesions inhibit cell migration and need to be broken down to allow cells to migrate. VEGF can break endothelial cell-cell contacts by disrupting the VE-cadherin/β-catenin complex at adherens junctions [19].

Cell migration is directed by chemotaxis
Chemotaxis is the directional migration of cells in response to gradients of extracellular soluble chemicals. Various cytokines regulate chemotactic migration of ECs, but the three key players are VEGF, bFGF (basic fibroblast growth factor) and angiopoietins [19].
In their study of retinal angiogenesis Gerhardt et al. demonstrated that VEGF independently regulates EC migration of tip cells and proliferation of stalk cells 7 (Gerhardt, Golding et al. 2003). Their experiments showed that endothelial tip cells extend long filopodia in response to VEGF. These filopodia were guided by a VEGF gradient, resulting in the directed migration of the tip cells ( Figure 1). However stalk cell proliferation did depend on the actual concentration of VEGF instead of the gradient. Both functions are mediated by the receptor VEGFR2, but the signals seem to be interpreted differently by the two EC subtypes. Barkefors et al. used a microfluidic chemotaxis chamber (MCC) to study EC migration in response to different VEGF gradients [21]. They observed that a stable gradient of both the isoforms VEGF165 and VEGF121 is sufficient to induce chemotaxis of ECs. Since VEGF121 is unable to interact with several coreceptors, this proved that a stable gradient of VEGF suffices for chemotaxis and interactions between VEGF and coreceptors are not required for a chemotactic response. Furthermore, the authors identified a minimal gradient steepness required for induction of chemotaxis and they showed that chemotaxis is reduced when cells reach the high end of the VEGF gradient. The experiments demonstrated that the shape of the VEGF gradient controls the migratory response of ECs.

Haptotactic migration can play a role in angiogenesis
The local degradation and deposition of matrix proteins by ECs and the heterogeneity of the extracellular matrix can all create local gradients of ECM components which 8 can drive endothelial cell migration, a process called haptotaxis. Haptotaxis of ECs is mainly triggered by the adhesive interactions between ECM components and integrins [19].
While the role of haptotactic migration in angiogenesis has not yet been established [18], in vitro experiments have shown that collagen [22] and fibronectin [23][24] gradients can guide EC migration, and therefore it is plausible that gradients of these components in vivo may lead to haptotaxis as well.
Senger et al. studied the function of two specific integrins in haptotactic migration of ECs in a gradient of immobilized collagen I [22]. Addition of antibodies against these integrins resulted in a significant reduction of directed migration towards collagen.
Smith et al monitored ECs on substrates with linear gradients or uniform concentrations of fibronectin in a highly controlled environment [23]. The experiments demonstrated that the drift speed of ECs increased on fibronectin gradients compared to uniform substrates. In a subsequent study they measured the response of cells to a range of fibronectin gradient slopes [24]. They showed that the cellular drift speed increased linearly with haptotactic gradient slope.

Haptokinesis: Cell sensitivity to ECM concentrations
While haptotaxis is the directional migration of cells up ECM gradients, haptokinesis is the sensitivity of cells to absolute concentrations of ECM components. Several experiments demonstrated that cell speed, spreading and membrane activity show a biphasic dependence on ECM concentrations, both on 2D substrates [25][26][27][28][29][30] as in 3D matrices [31]. We will describe some of these experiments and give possible explanations for the experimental observations. In chapter 4 we will describe our computational single cell experiments with which we tested each of these hypotheses.

Haptokinesis explained with the detachment theory
Palecek et al. measured the mean speed of migrating CHO (Chinese hamster ovary) cells plated on different concentrations of fibronectin and fibrinogen [28]. This speed 9 was maximal at intermediate ECM ligand concentrations regardless of integrin expression level or integrin-ligand binding affinity. At lower ligand levels cells were more rounded and extended more unstable lamellae that couldn't move the cell body.
At high ligand levels cells were very spread and extended lamellae similar to migrating cells, but the cell body didn't move very well.
In another study [27] DiMilla and coworkers showed that the migration behavior of explanation for this behavior, which we shall call the 'detachment theory', is the idea that at low ECM densities, a cell cannot form strong and stable adhesions at the front to generate a traction force, so no movement is possible and the cell spreads poorly.
At high densities a cell cannot detach adhesions from the substrate and therefore the cell will be well spread and immobilized, so again locomotion does not occur.   Figure 2 Maximal (a) cell speed [28] and (b) spreading [29] at intermediate ECM concentrations.

Haptokinesis explained with the receptor saturation model
Although in Palecek's experiments [28] [29]. They measured the projected area, migration speed and traction force at various type I collagen surface densities in a population of fibroblasts. Initially the cell area was an increasing function of surface density, but above a certain concentration the area declined. This threshold collagen density was approximately equal to the cell surface density of integrin molecules.
The 'receptor saturation' model can explain these observations [29]. At low densities, the number of ligands available to a cell is low, and therefore the cell cannot spread effectively. With increasing densities, more integrins can bind to ligands, leading to increased spreading. At the transition point all the integrins on the cell are bound to the substrate. Further spreading is impossible since there are no free integrins available. If the collagen density is increased beyond this point, the saturation of integrin receptors will be possible with a lower level of spreading. So the cell is less spread even though the substrate is more adhesive.

Haptokinesis explained with altered signaling
In the study of Cox et al. cells were plated on different concentrations of fibronectin [26]. The authors reported that membrane activity was maximal at intermediate Their findings suggest that adhesion-dependent signaling is a mechanism to stop cell migration by regulating cell polarity and protrusion via genes of the Rho family. Cox and coworkers imply with their study another mechanism to explain the haptokinetic observations, which we will refer to as 'altered signaling'.

Other mechanisms can influence cell migration
Apart from chemotaxis, haptotaxis and haptokinesis, several other mechanisms can regulate the migration of ECs. For example shear stress mechanically influences the response of ECs to haptotactic and chemotactic signals [19]. Aligned fibers in the ECM can guide cell migration (topographic guidance) and these guiding structures can in turn be remodeled by EC tip cells [20]. Differences in ECM rigidity or stiffness can also direct migration, a process called 'durotaxis' [20,31].

Sprouting endothelial cells must break through their basement membrane and invade
into the extracellular matrix in order to form a new capillary. This process requires proteolytic degradation of the ECM [33].
Endothelial tip cells express matrix metalloproteinases (MMPs) that can break down ECM components. MMPs are a family of more then 20 zinc-dependent enzymes.

12
They can be divided into two structurally different groups, the secreted MMPs and membrane-type MMPs (MT-MMPs). Together the MMPs can degrade all known mammalian extracellular matrix proteins. This includes fibrin, whose breakdown was usually ascribed solely to plasmin, another important proteolytic enzyme [33][34].
Most MMPs, especially MMP-2, MMP-9 and MT1-MMP, appear to be required for angiogenesis. Quiescent ECs produce little or no MMPs, but during wound healing, inflammation and tumor growth MMPs are strongly induced and activated in capillary sprouts [33][34]. Furthermore, several studies have shown that MMP inhibitors could inhibit angiogenesis [34][35]. Activation of proteases can be induced by angiogenic growth factors and inflammatory cytokines. For example MT1-MMP is induced in ECs by several factors, such as VEGF and HGF (hepatocyte growth factor) [33][34].

Other functions of MMPs and their proteolytic products
It was long thought that the only function of MMPs was to degrade ECM components. Recent studies however show that extracellular proteolyses can also regulate endothelial cell function in a more indirect way.

13
Growth factors bound to ECM components can be released by MMPs. For example, Hawinkels et al. reported that MMP-9 can release matrix-bound VEGF, making it more available to VEGF receptors [38]. Furthermore several angiogenic growth factors like VEGF and TGF-b (transforming growth factor beta) require proteolytic processing to become active [33].
Proteolytic fragments of the ECM and other molecules have been reported to show regulatory activity in angiogenesis, either positive or negative. They are often called 'matrikines' [33].

Summary
In this chapter we discussed some of the main mechanisms that are involved in tumor To what extent play chemotaxis, haptotaxis, haptokinesis and proteolysis a role in sprouting angiogenesis? Are they sufficient mechanisms for the growth of new vessel sprouts and the formation of branching vascular structures? Mathematical models of angiogenesis can be a means to answer these questions. In the next chapter we will discuss several existing models of angiogenesis.

The position of modeling in angiogenesis research
Mathematical models of angiogenesis can offer insight into the processes driving angiogenesis, and test or propose new hypotheses. Many models have been developed of angiogenesis over the last 30 years [39]. This was enabled by the boost in available biological data on this topic and the increasing computational capabilities.
Since tumor angiogenesis is a complex process, the models to date only address a part of the aspects involved in angiogenesis. By focusing on specific mechanisms that can influence the formation of capillary sprouts, these models can help to find the necessary conditions that are required for angiogenesis, such as cell proliferation, haptotaxis, haptokinesis or chemotaxis.
A model is of most value if it is able to reproduce experimental observed phenomena, without explicitly prescribing such events. Such a model should be able to produce emergent behavior from lower level rules. For example, many models of sprouting angiogenesis use high-level rules for branching. We will discuss some of them in the next section. In these models branching is not an emerging phenomenon, but a prescribed event.
In this chapter we will discuss three model categories: continuum models, discrete models and hybrid models ( Figure 3). We will describe examples of each category and discuss the advantages and disadvantages of these types of models. We will focus on

Continuum models
Continuum models represent new blood vessels in terms of cell densities, using partial differential equations (PDEs) to describe the average migration and proliferation of cell populations. In these models areas with high densities are blood vessels. Several models make a distinction between the tip and the stalk of a sprout, assuming that tip cells guide the sprout which is formed just behind it.
Although continuum models can provide valuable insight into aspects of angiogenesis, the disadvantage of these models is their use of cell densities and PDEs.
This assumes a large amount of cells to be involved in the process, whereas sprouting angiogenesis concerns a limited number of endothelial cells. As a result continuous models cannot take individual cell interactions and behaviors into account. In addition, these models are not able to predict the actual tree-like vascular structure, because they do not distinct separate sprouts.

A study of chemotaxis and haptotaxis in angiogenesis
The two-dimensional continuum model of Anderson and Chaplain [40] uses PDEs to describe EC, tumor angiogenic factor (TAF) and fibronectin densities. EC migration at the tip of a sprout is influenced by three factors: random motility, (saturated) chemotaxis towards TAF and haptotaxis towards fibronectin. In this model cells do not proliferate. The model incorporates uptake of TAF and fibronectin by ECs.
An important result from simulations was that for the outgrowth of the capillary network a sufficient strong chemotactic response was essential. The interactions between ECs and the ECM were important as well. The uptake of fibronectin and TAF created local gradients that allowed for lateral movement. Without this haptotactic response, cells migrated directly to the tumor.

17
The model contradicts experimental observations on two points: the speed of the vascular front decreases when approaching the tumor, while in reality the speed increases [39]. Second, although the model has not included proliferation, it produces vessels that in some cases reach the tumor, which disagrees with experimental observations as well [10]. Furthermore, this model is not able to capture relevant processes that happen on a smaller scale, such as sprout branching. Therefore the continuum model was converted to a discrete model, applicable at the level of a single cell. This model will be discussed in the next section where we describe discrete models.

Modeling the onset of angiogenesis
Levine et al. [41] presented a one-dimensional continuum model that describes the onset of angiogenesis. It concentrates on the first phase of neovascularization, namely the changes within the existing vessel. The model tries to predict the site of sprout formation.
In the model EC migration is considered to be a diffusive process which can be Simulations demonstrated that if there was enough angiogenic growth factor supplied to the capillary wall, the basal lamina would break down. Inside the fibronectin opening two aggregated peaks of EC concentration arose, forming the lining of the growing new sprout.

Extending the model with EC migration into the ECM
Levine and coworkers [42] coupled this one-dimensional model of the initiation of angiogenesis to a two-dimensional model that describes the migration of ECs into the ECM towards the tumor, using mostly the same cell behaviors and biochemical kinetics from the earlier model.
In this complex model cells proliferate in the ECM in response to the proteolytic enzyme and this proliferation is localized behind the leading tip of the sprout. The complex structure of the ECM is accounted for using a porosity constant.
The model uses a phenomenological approach where several experimental observations, such as the localization of proliferation and the acceleration of the growing sprout near the tumor, are explicitly included in the model. Therefore these observations cannot be predicted or understood with this model, they can only be accounted for [39].

Mechanochemical forces in blood vessel formation
The model presented by Manoussaki [43] considers both mechanical and chemical interactions during vasculogenesis (the formation of the initial vascualar network during development) and angiogenesis and investigates the effects of these mechanochemical forces on blood vessel formation.

Discrete models
Discrete models represent cells as single entities that can behave independently and move, grow and divide given certain prescribed rules. With discrete models one can define cell behaviors and interactions with their local environment in order to show how these can yield complex structures. Since a vessel sprout is normally one or a few cells wide and sprout formation involves stochastic mechanisms, cell-based models are better suited to describe the cellular dynamics during angiogenesis than continuum models, which describe sprouts as cell densities using deterministic differential equations.

A lattice free agent based model
Stokes and Lauffenburger were one of the first to use a discrete model to describe angiogenesis [46]. They presented a two dimensional, lattice free, agent based model, where Simulations demonstrated that the migration rate of ECs mainly determined the rate of vessel outgrowth and that a directional movement, in this model provided by chemotaxis, was required for directed network growth.

Discretizing the continuum model
By discretizing the PDEs in their continuum model and translating them into movement probabilities, Anderson and Chaplain derived a discrete biased random walk model of angiogenesis [40].
The model is based on the assumption that an endothelial cell at the tip of the sprout determines the motion of the whole sprout. In addition sprout branching, anastomosis and cell proliferation were incorporated in this model. The generation of new sprouts occurs at existing sprout tips and only when certain conditions on the sprout age, EC density and available space are met. The probability that a branch is formed is dependent on the local TAF concentration. Proliferation is added to the model by allowing division of cells at the sprout tip when they are above a certain age.
The results from simulations were similar to those from the continuum model, only now realistic networks were formed. Again, without haptotaxis, the sprout grows faster and there is less lateral movement and therefore less branching. The results confirm those from the continuum model that both chemotaxis and haptotaxis are necessary for the formation of vascular networks in tumor induced angiogenesis.
The discrete model was able to reproduce realistic networks structures. It did reproduce anastomosis, the dendritic structure of the capillary network and the formation of the 'brush border', the increased branching density near the tumor [39].

Migration following collagen cues
Yin and coworkers [44] used a novel microfluid device to study and quantify the Furthermore ECs move chemotactically to higher VEGF concentrations. Cells can proliferate if they have only one neighbor.

22
Simulations were able to produce patterns resembling growing and progressively branching vessel sprouts. Yin and colleagues showed that no specific rules on branching and no heterogeneous ECM densities were required to form vascular branching patterns

Hybrid models
Hybrid models combine continuum and discrete models to represent endothelial cells and other components that play a role in angiogenesis.
In the deterministic 3D model of Milde and co-workers [45]

The Cellular Potts Model
The Cellular Potts Model (CPM), also known as the Glazier-Graner-Hogeweg (GGH) model, is a lattice based framework, which has frequently been used to model a broad variety of biological phenomena [47][48][49][50]. It was developed by Glazier and Graner and represents cell dynamics, like interactions between cells and changes in shape, in terms 23 of a generalized energy [50]. It finds the lowest energy of a system using a Monte Carlo simulation technique with Metropolis dynamics.
The CPM has many advantages over other cell-based models. It is a simple model, which is easy to extend. The irregular stochastic cell membrane dynamics, like extensions and retractions of protrusions, are explicitly represented, which allows simulating cell interactions with their environment and mimicking the exploratory behavior of cells.

CPM models of angiogenesis
The CPM is used in several simulations of tumor induced angiogenesis. We will discuss four of them below (Figure 4), focusing on the objective of each study, the rules which where incorporated in the models to describe EC behavior during angiogenesis and some of the results.
In  self-organize the ECs into capillary-like networks, the model was extended with autocrine chemotaxis to a very short-diffusing chemoattractant, as described in [48].
The model neglects proteolytic activity of ECs.
Szabo et al. [54] argued that during sprout formation stalk cells cannot be totally

Summary
Many models of angiogenesis have been developed to get more insight in the mechanisms involved in this process. In this chapter we discussed continuum, discrete 26 and hybrid models. We argued that cell based models are the best choice for modeling angiogenesis, since angiogenic sprouting involves only a few cells and we therefore need to able to define rules on the level of an individual cell. In particular, the CPM can be a good basis for this goal, since it explicitly models the irregularities in cell shape and behavior.
Each of the discussed CPM models used different sets of rules on the level of individual cells to model sprout formation during angiogenesis. This shows that different mechanisms, such as matrix inhomogeneity, chemotaxis to autocrines, contact inhibition of chemotaxis, cell polarity, and preferential adhesion to elongated cells can all play a role in the formation of vascular patterns.
In the next chapter we will present our model of tumor induced angiogenesis, where we focus on the interactions between ECs and the ECM. It is based on the CPM and extended with these cell-matrix interactions. We used a two dimensional Cellular Potts Model (CPM) to implement these rules, we therefore start with a detailed description of the CPM in the following section. Next we will describe how we incorporated our model into the CPM. We give a detailed survey of ways to model haptokinesis and their implementation in the CPM

CPM explained in more detail
In the Cellular Potts Model biological cells are represented as patches of lattice sites.
Each cell has a unique index σ, which is assigned to every lattice site that is occupied by that cell. Furthermore the type of a cell σ is denoted with τ(σ). The extracellular matrix consists of all lattice sites not occupied by cells and is labeled with index σ=0 and type τ=0.
An area constraint penalizes cell shapes deviating too much from their preferred area.
To mimic cell elongation a length constraint is added [48]. The ECM has no area or length constraint. The 'effective energy' is given with the CPM Hamiltonian: where x and x' are neighboring lattice sites; δ x,y = {1 , x=y; 0, x≠y}; a σ is the current area of cell σ, A σ its target area and λ the inelasticity; l σ represents the current length of cell σ, L σ its target length and λ' the strength of the length constraint.
To mimic membrane extensions and retractions, we repeatedly attempt to replace the index σ of a random lattice site x by one of its random neighboring sites x'. 29 We calculate ∆H, the change in total effective energy if we performed the copy and accept the attempt with Boltzmann probability: where T corresponds to the intrinsic motility of the cells. By accepting energetically unfavorable moves, we prevent the system from getting trapped in local energy minima.

Chemotaxis can be incorporated by including an extra reduction in energy for
extensions and retractions towards higher concentrations of a chemoattractant (as described in [49]): where x is the site into which neighbor x' copies its spin, c(x) is the local concentration of chemoattractant at site x, and µ is the strength of the chemotactic response.

Model setup
Our model domain is a rectangular dish (size 500 µm × 700 µm) in which endothelial cells are placed behind a vessel wall situated at the 'bottom' of the dish. The cells can migrate through a gap in the wall into the ECM towards the 'top' of the dish in the direction of a tumor which we assume to be located beyond the top of the dish ( Figure 6). We define only one type of cell, so we do not make a distinction between tip cells and stalk cells. Figure 6 The simulation is initialized with ECs placed behind a vessel wall.
A lattice site represents an area of 2 × 2 µm. We will set the target area A σ = 50 lattice sites corresponding to 200 µm 2 and the target length L σ = 30 µm. We set the adhesion energy at cell-cell borders J CC =40 and the cell-matrix energy J CM = 25, in order to make attachments between cells slightly favorable over cell-matrix bonds. The intrinsic motility parameter T =100.
When doing a copy attempt we select the source site from the twenty first to fourth nearest neighbors, to improve the isotropy. During a Monte Carlo Step (MCS) we carry out N copy attempts where N is the number of lattice sites. We define a high cell-border energy to prevent cells from adhering to the boundaries of the lattice.

Modeling VEGF, MMPs and ECM concentrations
We will include three layers in our model that describe the concentration of VEGF,

VEGF layer
Since the tumor will be relatively large, we assume a constant linear VEGF gradient with equal concentrations at equal y positions throughout the simulation. Therefore we initialize the VEGF layer at the start of the simulation and do not alter it by diffusion, secretion or degradation at later steps (Figure 7).
Given the diffusion coefficient D and the degradation rate ε of VEGF and the concentration at a starting point c(0), we can calculate the analytical solution of the gradient at steady state [56]. If the gradient is caused by diffusion and linear decay, the 1D solution is: Figure 7 The constant VEGF gradient (a) and initial ECM density (b).

MMP layer
We will describe the secretion, decay and diffusion of MMPs with a PDE model: where c E (x,t) and c M (x,t) represent the concentration of ECM and MMPs respectively at site x at time t ; α E , ε E_M are the secretion rate and decay rate of the ECM components.
Since in our model degradation only happens at ECM sites, the non-diffusible ECM components will not decay at lattice sites that are occupied by the cells of a growing sprout. Therefore we set the secretion rate α E = 0, creating an ECM concentration at cell sites being the net result of balanced secretion and decay.
We did try an alternative solution by allowing degradation at cell sites and adding secretion of ECM components. However, to avoid a fast accumulation or decay of ECM components the secretion and decay rates needed to be well-tuned, resulting in a model that was very sensitive to both these parameters. 33 We solve the two PDEs numerically with a finite difference scheme on a similar lattice as the one used in the CPM with ∆x = 2 µm. We set ∆t = 2 seconds and run 15 'diffusion steps' between subsequent MCS. We use fixed boundary conditions. To

Modeling chemotaxis and haptotaxis
Since primarily the extending filopodia of (tip) cells are able to sense and react on chemotactic cues, we consider only extensions of cells into the ECM to contribute to the chemotaxis energy term.
The energy change due to chemotaxis then becomes: where x is the site into which neighbor x' copies its spin, c V (x) is the local concentration of VEGF at site x, and µ is the strength of the chemotactic response.
This corresponds to incorporating contact inhibition of chemotaxis, reflecting VEcadherins suppression of pseudopods, and extension-only chemotaxis as described in [52].

34
In a similar way we implement haptotaxis, the migration towards higher ECM densities. One modification was made however: it would not be realistic for cells to move towards very high densities; therefore we include haptotactic saturation that ensures this restriction: Here c E (x) is the local ECM density at site x, τ is the strength of the haptotactic response and s is a saturation factor.
Note that we consider concentrations of ECM components to be present at or below all lattice sites, so at cell sites as well, since cells are placed in the matrix and we consider a matrix to be present below the cell as well.

Modeling Haptokinesis
In chapter 2 we have described the detachment theory, the receptor saturation model

Detachment theory in the CPM
The detachment theory states that with a low concentration of ECM molecules there will be not enough adhesive forces for the cell to adhere to the matrix. On the other hand with high ECM densities the adhesions will be too strong to release and retract the rear end of the cell. We implemented this theory in the CPM by replacing the   ECM density Figure 10 Testing the detachment model. Cell shape is rounded at low densities and irregular at high ECM densities.

Receptor saturation in the CPM
Next we implemented the receptor saturation theory in our model, which states that cells have a limited amount of integrins available to bind to ECM molecules. When the concentration of ECM molecules is low the cell cannot spread optimally. When the concentration is higher, spreading will increase up to a certain threshold level after which the spreading will decrease since all integrins of the cell can make bonds using less cell-matrix surface.
We implemented this theory into the CPM by summing the ECM densities at all lattice sites of a cell. This total concentration is proportional to the matrix molecules that can bind with integrins. We define a maximum number of bonds, which is related to the maximum number of available integrins, the 'saturation threshold'. The number of bonds of a cell is the minimum of this total ECM concentration and the saturation threshold: Making bonds is considered favorable and therefore we introduce an energy difference term proportional to the difference in number of bonds before and after the copy attempt. This model did show an optimal spreading at intermediate levels of collagen concentration, but did not result in optimal speed at those levels ( Figure 11). On contrary, the speed decreased when spreading was optimal.
However, when we take into account that adding a lattice site to a cell with a large area will have less effect on its center of mass compared to a cell with small area and correct the measured speed by multiplying with the area, we see an increase in speed at intermediate ECM densities Figure 12. However the sudden drop in speed and spreading doesn't agree with experimental results.
(a) (b) Figure 11 Testing the receptor saturation model. Average speed (in µm/min) (a) and cell area (in µm) (b) of cells placed in different uniform ECM densities. η = 200. Error bars represent the standard error of the mean.
39 Figure 12 Testing the receptor saturation model. Speed × area as a function of ECM density. η = 200. Error bars represent the s.e.m. n = 50. We integrated this into the CPM in the following way. When cells try to extend into the ECM, the probability of this step is highest at intermediate collagen levels.

Altered signaling in the CPM
Intermediate levels will in this case result in an energy decrease, low or high levels will result in an increase in energy. The energy difference caused by haptokinesis is therefore given by the following Gaussian function ( Figure 13): where µ is the intermediate ECM density, (c max -c min )/2 and ρ is set to (c max -c min )/5. This results in a bell shaped curve for speed and spreading for low to high collagen concentration levels ( Figure 14), which agrees with experimental results [28]. We therefore selected the altered signaling model to describe haptokinesis in the CPM.

Proliferation
During sprouting angiogenesis stalk cells can proliferate. As described in chapter 2, some studies suggest that this only happens just behind the tip cell, while others argue that stalk cells divide at the base of a sprout. The implementation of any of these 41 options is however complicated. Since VEGF promotes proliferation and following the suggestion that cell division occurs at gaps between ECs [8], we allow proliferation for those cells for which a relatively large part of their surface is in contact with the ECM (and thus in contact with VEGF). The probability of mitosis grows when this proportion increases: where ρ is the ratio (cell-ECM surface/total cell surface) and ρ min is a threshold ratio.
We use Hogeweg's model [58] of cell division by assigning a new index to the grid points on one side of the shortest axis of the dividing cell and giving the daughter cells half the target area and a decreased target length. This polarized division is in agreement with experimental observations of proliferating ECs [59]. By slowly incrementing the target area and target length, cell growth is implemented. In our model we increase the target area of cells every 5 MCS with 2 lattice sites.
Furthermore proliferation is only allowed outside the parent vessel, at a minimum distance of one cell length from the vessel wall.

Summary
Our cell based model of tumor induced sprouting angiogenesis uses the CPM framework to describe the behavior of individual endothelial cells. We added three layers to represent the VEGF, MMP and ECM concentrations and incorporated chemotaxis, haptotaxis, haptokinesis, proteolysis and proliferation as energy difference terms.
The next chapter describes the results of our simulations. We will show that this simple model is able to reproduce realistic vascular sprouts. We use compactness, sprout height and sprout size as quantitative measurements to evaluate the consequences of varying parameters like chemotaxis, haptotaxis and haptokinesis strength or ECM decay rate.

Results
We ran several simulations with parameter values as described in the previous chapter and as indicated in Table 1 Our simulations reproduce vessel sprouts that grow towards the tumor and frequently form branches (Figure 15 and Figure 16). Frequently anastomosis occurs when two 44 sprouts rejoin ( Figure 16 d-g). In some cases branches split off, probably due to a too high pulling force from the leading cells ( Figure 16 g-h). 5000 10000 15000 20000 25000 30000 Figure 15 Example of a growing sprout. Subscripts: number of MCS. All parameters as indicated in Table 1. Figure 16 Examples of growing sprouts after 30000 MCSs. All simulations have parameter settings as indicated in Table 1.

Measuring compactness, height and size of sprout
Although properties like sprout length and size and branching frequency can be evaluated upon visually inspecting the graphical output (pictures and movies) of the simulations, we would prefer to use quantitative measurements that can give us an indication of these properties in order to do a faster and more thorough analysis of the results.
Before doing measurements on the growing sprout, we need to define which cells belong to this sprout. We therefore define the largest 'blob' as the largest set of connected cells and we will consider all cells in this blob that are outside the parent vessel to belong to the sprout. This means that individually migrating cells or branches that have split off are not included in our measurements. We use a fast union-and-find algorithm to determine the largest blob among all cells [60].
We calculate the 'compactness' of the sprout as a measure of branching. We first draw a convex hull around those cells in the largest blob that are located outside the parent vessel. The compactness is the ratio between the total area of these cells and the area of the convex hull [52]. Usually sprouts with high compactness have no or few branches and grow in a direct line towards the tumor. However, extensive branching sprouts can have high compactness as well, since these branches will at the end fill up almost all space.
To measure how fast the new vessel is growing we need to measure the distance of the growing vascular network from the parent vessel. We define the 'height' of the sprout as the largest distance in y direction between any lattice site in the largest blob and the bottom of the dish.
We are also interested in the growth in size of the sprout; we therefore define the 'size' of the sprout as the number of cells of the largest blob. A large size can also be an indication of excessive branching. 46

Sensitivity analysis
In order to find out which of the mechanisms of our model plays a key role in sprouting angiogenesis and to study the sensitivity of the model to certain parameter variations we ran simulations where we varied the parameter of interest and kept all other parameters fixed.

Chemotaxis
We first investigated the role of chemotaxis by varying the chemotaxis strength µ.
Compactness and the size of the sprout are not sensitive to this parameter, but the sprout grows faster with higher values of µ ( Figure 17). This makes sense, as increased chemotaxis will more strongly pull the growing sprout in the direction of the tumor.

Haptotaxis
Next we looked at the strength of haptotaxis. While compactness and height show no major differences, the size of the sprout increases with haptotaxis strength ( Figure   20), producing more complex, branching networks ( Figure 18). 0 300 900 1500 Figure 18 Examples of growing sprouts after 30000 MCSs with different haptotaxis strength (subscripts).
We can explain this as follows. Haptotaxis promotes migration to higher ECM densities, but up to an intermediate ECM concentration. The ECM density near the top of a sprout is higher than in areas at the sides of the sprout tip, therefore haptotaxis can stimulate lateral movement of (cells at) the sprout tip ( Figure 19) and this can induce branching. Figure 19 Lower ECM densities at the sides of the sprout compared to the top. Haptotaxis to those densities can promote lateral movement of the sprout.

Haptokinesis
We next looked at haptokinesis. Here we see differences in compactness, height and size ( Figure 21). If haptokinesis is set to zero, only small sprouts are formed, that do not grow much in size. For small values of haptokinesis sprouts grow very fast towards the tumor, with no or few branches. For larger values of haptokinesis strength the sprout speed decreases and more branches are formed. How can we explain this phenomenon? Without haptokinesis, the only mechanisms playing a role in the directional movement of cells are chemotaxis and haptotaxis. The VEGF gradient will stimulate cells to move towards the tumor, but apparently more is needed to form a growing sprout. When haptokinesis is 'switched on', but with a low strength, a sprout can be formed, and the relatively high chemotaxis strength will pull the sprout fast towards the tumor. With increasing haptokinesis, the migration forward through high ECM densities is more constrained, resulting in an increase in branching and a slower growing sprout ( Figure 22). 50 Figure 22 Examples of growing sprouts after 30000 MCS with different values of chemotaxis and haptokinesis strength.

Degradation of ECM by MMPs
Next we investigated the role of proteolysis, the degradation of ECM components by MMPs. Varying the decay rate of ECM affects the compactness, height and size of the sprout ( Figure 23 and Figure 24).
Increasing the decay rate produces fast growing vascular networks with a lot of branches. Variations in secretion rate of MMPs will have similar effects (results not shown).

Variation in VEGF gradient
We also varied the VEGF gradient, by varying the decay rate of VEGF, to see what the effect of a shallow or steep gradient will be. In all simulations the VEGF concentration in the parent vessel is set to the same value (0.05) to ensure that cells start from the same initial situation ( Figure 26).

Sprout formation without proliferation
When we do not allow cells to proliferate, branching sprouts are formed, however they grow must slower and parts of the sprout split off. The sprout will not reach the tumor ( Figure 27 and Figure 28).

Discussion
We demonstrated that our model, which describes cell-matrix interactions on the level of individual cells, is able to reproduce sprouting and branching behavior during angiogenesis. Most models that were discussed in chapter 3 needed specific branching rules to form vascular trees. Our simulation results show that our model can reproduce branching sprouts and networks, without including such specific rules.
We can conclude from our sensitivity analysis that haptokinesis, the sensitivity of cells to ECM densities, and proteolysis, the degradation of ECM components, seem to play a key role in angiogenesis. Neither chemotactic migration nor haptotaxis alone appear to be sufficient mechanisms to form stable vascular sprouts. Their role appears to be more a supporting one; chemotaxis guides the growing sprout towards the tumor and haptotaxis promotes branch formation by stimulating lateral migration of the sprout.

Comparing results with experimental observations
We can validate our model when we compare our results with experimental observations.

MMPs are essential for angiogenesis
In our model proteolysis is required to form growing sprouts. Without ECM degradation cells cannot invade the matrix. This is in agreement with studies that demonstrate that MMPs are essential for cell migration through 3D collagen matrices [33].

The steepness of the VEGF gradient
The spatial distribution of VEGF has several influences on angiogenesis. In a shallow VEGF gradient, tip cells will extend short undirected filopodia, resulting in a slow and undirected growth of sprouts. Shallow gradients will stimulate proliferation of stalk cells. Steeper VEGF gradients lead to excessive branching of vessels [7].

56
In our simulations with different VEGF gradients we see an increase in the speed of the sprout for steeper gradients. However, as discussed before, this must partly be ascribed to increased proteolysis. If we look at the effect of varying chemotaxis strength, which is similar to varying the VEGF gradient, we again see an increase in the speed of the sprout, but no change in branching frequency.

How fast do sprouts grow?
Many in vitro experiments observe that sprouts grow with a speed of ~ 1 mm in three days. This is about five times faster than our experiments. However other experiments [8] mention 10-12 days for sprouts to form and reach a tumor and this is in agreement with our observations.
In vivo and in vitro experiments show that the speed of a sprout increases when it approaches the tumor [39]. In our simulations sprouts accelerate also when they approach the tumor. This can be ascribed to the increased secretion of MMPs due to higher VEGF concentrations.

Cell migration and proliferation
Experiments show that both EC migration as well as proliferation plays a role in the formation of vessel sprouts. If proliferation is inhibited sprouts can form, but will not reach the tumor [10]. Our simulations agree with this, in the first phase of angiogenesis the growth of the sprout is mainly caused by migrating cells, in a later stage proliferation is responsible for sprout growth. Furthermore, we have shown that without proliferation we can reproduce branching sprouts which remain small in size.

Proteolysis
Although we can question whether the endured secretion of MMPs by all cells and the subsequent degradation around the sprout is realistic, classic descriptions of extracellular proteolytic activity during angiogenesis, talk about capillary sprouts surrounded with a clear space, resulting from the dissolvement of fibrin [34]. In our model this 'empty space' along the stalk of the sprout is one of the main reasons for 57 the formation of small and stable sprouts. Cells are not likely to move away from the sprout towards very low ECM densities.

Brush-border effect
A typical phenomenon in sprouting angiogenesis is the 'brush-border' effect. In the proximity of the tumor the branching of sprouts will increase resulting in a dense vascular tree [39]. In our simulations we see an increase in branching half way the dish, but when the sprouts move further, they accelerate and move to the tumor in a straight line. We assume that this is due to stronger chemotaxis and to the increase in proteolysis caused by higher VEGF concentrations.
The 'brush-border' effect could be caused by increasing stalk cell proliferation and tip cell selection induced by higher VEGF concentrations near the tumor. These phenomena are not incorporated in our model.
According to the hypothesis of Yin and coworkers [44], increased ECM degradation near the tumor will decrease the haptokinetic and haptotactic effects on cell migration.
As a result cells lose the ability to follow the ECM traces left by other cells, which results in increased branching. In our model the ECM density at cell sites is rather constant, and therefore migrating cells can follow ECM cues, even near the tumor.

59
C h a p t e r 7

Future work
Improvements can certainly be made to our model. We will discuss a number of them in this chapter. However, when adding new rules and restrictions to the model, it will loose its simplicity, which is in part one of its charms.

Adding tip and stalk cell behavior
We can make our model more realistic by introducing two cell types: tip cells and stalk cells. Tip cells are more motile; they lead the sprout, navigate by extending filopodia, and invade the ECM by releasing proteases. Stalk cells follow the tip cells, they are less active, proliferate and secrete ECM components. In our current model all cells are sensitive to chemotactic and haptotactic cues, they can all proliferate and secrete MMPs and ECM components. It would make the model more realistic when we divide those tasks between tip cells and stalk cells.

Proliferation induced by VEGF
Proliferation is induced by VEGF, so we could improve the model not only by restricting proliferation to stalk cells, but also by increasing the probability of cell division with higher VEGF concentrations.

Improving model of cell-matrix interactions
Although we intensively studied ways to model the cell-matrix interactions, improvements can certainly be made. For example, in our current model we don't make a distinction between the front and rear of a cell, while in reality migrating cells do have polarity. When cells preferably extend protrusions at their front persistence will be much higher.
We studied three explanations of haptokinesis: the detachment theory, receptor saturation and altered signaling and finally implemented the latter in our model of 60 angiogenesis. However, we think that all three theories (and maybe others) play a role in haptokinesis and it would be interesting to look further into the underlying mechanisms.
In our model degradation of ECM proteins is not regulated in the sense that the secretion of MMPs and the degradation of ECM by MMPs are not related to the ECM density in the direct environment of the cell. However, in reality, cells can finetune proteolysis to prevent excessive break down of the matrix. We could therefore add functionality to the model that inhibits or limits proteolysis when ECM densities are low enough for invasion.
Other mechanisms that influence the migration of cells could be included in the model. For instance shear stress, matrix rigidity and the direction of matrix fibers, can all guide cells when migrating into the matrix. Cells in their turn can remodel the fibers in the ECM leaving cues for following cells.

Improving the representation of the ECM
The ECM is now modeled as a field with initially a uniform concentration of ECM components. In reality the matrix is heterogeneous with irregular concentrations of different molecules. In our model we do not make any distinction between different ECM proteins, we generally refer to ECM components and don't have in particular collagen, fibronectin in mind. Including these distinct components with their specific properties will also make the model more realistic.
We could furthermore add matrix bound factors to the ECM, such as certain VEGF isoforms, which can be released or activated by MMPs. These factors will set up steep local gradients and this will certainly affect cell migration [33]. .