Mechanistic characterization of endothelial sprouting mediated by pro-angiogenic signaling.

OBJECTIVE: We aim to quantitatively characterize the crosstalk between VEGF- and FGF-mediated angiogenic signaling and endothelial sprouting, to gain mechanistic insights and identify novel therapeutic strategies.
METHODS: We constructed an experimentally validated hybrid agent-based mathematical model that characterizes endothelial sprouting driven by FGF- and VEGF-mediated signaling. We predicted the total sprout length, number of sprouts, and average length by the mono- and co-stimulation of FGF and VEGF.
RESULTS: The experimentally fitted and validated model predicts that FGF induces stronger angiogenic responses in the long-term compared with VEGF stimulation. Also, FGF plays a dominant role in the combination effects in endothelial sprouting. Moreover, the model suggests that ERK and Akt pathways and cellular responses contribute differently to the sprouting process. Last, the model predicts that the strategies to modulate endothelial sprouting are context-dependent, and our model can identify potential effective pro- and anti-angiogenic targets under different conditions and study their efficacy.
CONCLUSIONS: The model provides detailed mechanistic insight into VEGF and FGF interactions in sprouting angiogenesis. More broadly, this model can be utilized to identify targets that influence angiogenic signaling leading to endothelial sprouting and to study the effects of pro- and anti-angiogenic therapies.

protein kinase B

average length

extended Fourier amplitude sensitivity test

endothelial nitric oxide synthase

extracellular regulated kinase

fibroblast growth factor

fibroblast growth factor receptor

fibroblast growth factor receptor 1

human umbilical vein cells

neuropilin‐1

number of sprouts

phosphorylated protein kinase B

platelet‐derived growth factor

phosphorylated extracellular regulated kinase

particle swarm optimization

phosphatase 2

total length

thrombospondin‐1

vascular endothelial growth factor

vascular endothelial growth factor receptor 2

INTRODUCTION

Angiogenesis, the formation of new blood capillaries from pre‐existing blood vessels, plays an important role in the survival of tissues, as the vessels carry blood throughout the body to provide oxygen and nutrients required by the resident cells. , Thus, researchers have been putting effort toward targeting angiogenesis as a strategy in many contexts, including in both tissue engineering and cancer treatment. , Specifically, one challenge for the success of synthetic tissues is to ensure sufficient transport of nutrients and gases such as oxygen to the cells. Thus, stimulating new blood vessel formation is an important strategy for the long‐term viability of engineered tissue constructs. On the other hand, inhibiting angiogenesis is an important strategy for cancer treatment, as tumors grow by obtaining nutrients and oxygen from the blood delivered by vessels, and tumor metastasis is facilitated by the blood circulation. , ,

The process of blood vessel formation, particularly for capillaries, is largely initiated and mediated by endothelial cells responding to the local physiological conditions. Many different pro‐angiogenic growth factors, such as fibroblast growth factor (FGF), vascular endothelial growth factor (VEGF), and platelet‐derived growth factor (PDGF), regulate angiogenesis. , These factors promote different cellular processes involving endothelial cell survival, proliferation, migration, and vessel maturation leading to new blood vessel formation. , Strategies to promote or inhibit angiogenesis focus on modulating the effects of these factors to alter the cellular‐level processes they induce, with a focus on endothelial cells.

However, not all approaches to promote or inhibit angiogenesis lead to successful outcomes. For example, clinical trials have shown no effective improvement in FGF‐ or VEGF‐induced angiogenesis. Also, bevacizumab, an anti‐VEGF agent, has limited effect in certain cancer types, and it is no longer approved for the treatment of metastatic breast cancer due to its disappointing results. Thus, there is a need to better understand the mechanism of angiogenesis, specifically the molecular interactions and signaling required for new blood vessel formation and how they affect cellular behaviors, in order to establish more effective therapeutic strategies.

In addition, there is crosstalk between intracellular signaling pathways. The overall response in endothelial angiogenesis is dependent on the integrated signals activated by these growth factor‐mediated pathways to influence cellular decisions, such as proliferation, survival, and migration, and further promote or inhibit vessel formation. However, the integrated effects of more than one factor on intracellular signaling reactions at a detailed level have not been studied in great detail. In pro‐angiogenic strategies, some synergistic effects between FGF and VEGF in endothelial angiogenic activities have been shown. , , However, the mechanism of how they act together quantitatively on regulating molecular and cellular behaviors is still not clear. Also, in the case of inhibiting angiogenesis, tumors often evade the effects of drugs that target a single factor by making use of alternate compensatory pathways to activate signaling species needed for proliferation and migration. For instance, FGFR activation may play a role in the resistance mechanism of anti‐angiogenic drugs, especially anti‐VEGF treatment. , Additionally, experiments show high levels of FGFR1 in tumors that continue to progress, even during anti‐VEGF therapy. Thus, it is necessary to better mechanistically understand the effects of multiple angiogenic factors and their crosstalk in activating angiogenic signaling and further cellular responses.

In this study, we are interested in the angiogenic signals required to initiate vessel growth. Therefore, we focused on the molecular signaling and further cellular responses in the process of endothelial sprouting. We consider sprouting promoted by FGF and VEGF, as they are particularly important in early stages of angiogenesis. , , We are particularly interested in signaling crosstalk between these factors required to initiate vessel growth. FGF and VEGF bind to their receptors and initiate signaling through the mitogen‐activated protein kinase (MAPK) and phosphatidylinositol 3‐kinase/protein kinase B (PI3K/Akt) pathways to phosphorylate extracellular regulated kinase (ERK) and Akt, respectively. ERK and Akt are important signaling species in the angiogenesis process that influence cell proliferation, survival, and migration. Thus, we aim to quantitatively investigate the combination effects of two major pro‐angiogenic factors, FGF and VEGF, on activating MAPK and PI3K/Akt signaling (at the molecular level) and further promoting endothelial sprouting (at the cellular level) in endothelial cells.

Given the complexity of biochemical reactions comprising angiogenesis signaling networks, a better quantitative understanding of the dynamics of these networks is beneficial for current angiogenesis‐based strategies. Computational modeling serves as a powerful tool to investigate molecular and cellular responses systematically. There are many published models that predict molecular , , and cellular , responses mediated by angiogenic factors. However, such models are mostly designed to predict responses upon single agent stimulation. Targeting more than one growth factor and exploring their effects in intracellular signaling and cellular responses in detail deserves more attention. In addition, many models that focus on specific cellular behaviors significantly reduced the intracellular signaling network such that the output signal is simply linearly proportional to the fraction of bound receptors. , , Therefore, we constructed a hybrid agent‐based model to characterize the intracellular signaling interactions of FGF and VEGF and utilized downstream signals (maximum pERK and pAkt) as inputs to describe angiogenic cellular responses in the process of endothelial spheroid sprouting. Overall, we predict the contributions of: (1) cellular processes including cell proliferation, sprout growth, and the chance of forming a sprout; (2) mono‐ and co‐stimulation of angiogenic factors FGF and VEGF; and (3) underlying molecular signals (pERK and pAkt) that modulate endothelial sprouting. Thus, we generate novel cellular‐ and molecular‐level insights related to cell sprouting. In particular, we describe total sprout total length (TL), number of sprouts (NS), and average length (AL), which are the well‐studied metrics in in vitro studies, by the mono‐ and co‐stimulation of FGF and VEGF. The model predicts that the type and concentration of ligand, length of growth factor stimulation, and initial number of cells are important in endothelial sprouting, and FGF is a dominant factor driving the combination effects of FGF and VEGF in endothelial sprouting. Also, the model suggests that proliferation of endothelial cells and the growth of existing sprouts are more important in the sprouting process compared to the effect of the chance of forming a new sprout. In addition, the MAPK and PI3K/Akt pathways regulate the vessel network in different ways. Specifically, the ERK pathway regulates vessel network mainly via regulating cell proliferation and NS, while the Akt pathway mainly affects vessel network via regulating sprout growth. Moreover, the strategies to modulate endothelial sprouting are context‐dependent, and our model can identify potential effective pro‐ and anti‐angiogenic targets under different conditions and study their efficacy.

MATERIALS AND METHODS

Model construction

We constructed a hybrid agent‐based model (Figure 1) that describes cellular responses, including cell proliferation, sprout growth, and the formation of new sprouts. We note that agents are probabilistic discrete cellular responses; however, we do not include spatial effects. These cellular responses are driven by molecular signals, pERK and pAkt, upon the mono‐ and co‐stimulation of FGF and VEGF. The molecularly detailed biochemical reaction network that characterizes the MAPK and PI3K/Akt pathways induced by FGF and VEGF is adapted from our previous work, which is referred to as the ERK‐Akt model in this study. In the ERK‐Akt model, FGF binding to FGFR1 and HSGAG activates FRS2 and then initiates PI3K/Akt and MAPK pathways, and VEGF binding to its receptor, VEGFR2, phosphorylates VEGFR2 and activates PI3K directly. In addition, activated Raf triggers MAPK pathway upon the stimulation by VEGF. This network is implemented as an ordinary differential equation (ODE) model using MATLAB. The main model includes 97 reactions, 99 species and 100 parameters (see detail in ).

FIGURE 1

Endothelial spheroid sprouting process. (A) Activated endothelial cells become tip cells and start to migrate and the stalk cells behind tip cells are proliferative. Finally, the endothelial cells sprout into linear cord‐like structures. (B) The sprouting process involves tip cell migration, stalk cell proliferation, and elongation. (C) Growth factors, FGF and VEGF, binding to their receptors initiate intracellular signaling and regulate cellular responses in all endothelial cells

To link the molecular signals (pERK and pAkt) produced by the angiogenic factors (FGF and VEGF) with the short‐term cellular responses, we made four assumptions: (1) The endothelial cell responses (cell proliferation, sprout growth, and probability of sprouting) are dependent on the maximum pAkt and pERK levels upon the stimulation of FGF and VEGF within 2 h; (2) the intrinsic properties of endothelial cells to grow and form sprouts are stable within 3 days of simulated cell culture, which leads to a constant average cell proliferation rate, sprout growth rate, and the probability of forming a new sprout within 3 days (described in subsequent sections below in detail); (3) daughter cells inherit all the properties from mother cells, specifically the cell proliferation rate, sprout growth rate, and probability of sprouting; and (4) the maximum pERK and pAkt drive endothelial cell proliferation and sprouting following Hill functions: where S refers to either pERK or pAkt, [S] is their maximal concentration, Km is the substrate concentration where the proliferation or sprouting rate is half of its maximum value, V max, and n is the Hill coefficient.

Module of endothelial cell proliferation

The endothelial cell proliferation module is a hybrid agent‐based model that simulates each endothelial cell as one agent that has its own cell proliferation rate and divides based on its own cell doubling time. The pERK and pAkt levels were used to inform the rate of cell proliferation, as has been done in other computational work. The total number of endothelial cells was quantified to account for cell proliferation, which is a net result of survival and proliferation of endothelial cells in response to FGF and/or VEGF stimulation. The ERK‐Akt model predicts the dynamics of the molecular signals pAkt and pERK within 2 h upon ligand stimulation, which are the inputs to calculate the rate of cell proliferation.

Based on the assumptions above, the average cell proliferation rate (rcp) is given by where r is the basal cell proliferation rate when no stimuli are applied, which is not affected by pERK or pAkt levels, k and k are the cell proliferation rate constants, V max_ and V max_ are the maximum rates of cell proliferation driven by maximum pERK and pAkt, respectively. For simplification, we use V and V to represent k × V max_ and k × V max_ , respectively. [max (pERK)] and [max (pAkt)] are maximum pERK and pAkt levels within 2 h upon the ligand simulation, respectively. Km and Km are the maximum pERK and pAkt levels that produce the half maximal of the cell proliferation rates, respectively. n is the Hill coefficient for the cell proliferation rate.

The average cell proliferation rate (r) indicates that the average doubling time for endothelial cells is 1/r. To account for the cell heterogeneity within a cell population, we assigned a cell proliferation rate chosen from a normal distribution with a mean (μ) of the calculated average cell proliferation rate and a standard deviation (σ) to capture 99.7% of the possible values given the range of μ ± 25%μ (i.e., μ ± 3σ) for each cell. The baseline values are provided in Table S1 and supplementary materials.

The total cell number at a particular time T, N tot (T), is given by where cell i from an initial cell population that consists of N int of HUVECs has a cell proliferation rate of r, T is cell culture time, floor (T × r) rounds (T × r) to the nearest integer less than or equal to (T × r).

Note that the mean cell proliferation from experimental data is usually calculated for several replicates. To compare with experimental data, we calculated the average total cell number for ten simulations per condition.

Model of endothelial cell sprouting

Model overview

The endothelial cell sprouting model is also a hybrid agent‐based model that simulates each endothelial cell as one agent that has its own properties and makes its own cellular decisions. The model utilizes a probabilistic approach to model sprouting. The general flow of the model is shown in Figure 2.

FIGURE 2

Flowchart of the endothelial sprouting agent‐based model. The model simulates each endothelial cell as one agent that has its own properties and makes its own cellular decisions in every time step

First, the model checks each cell to see whether it is a tip cell (i.e., if it is a leading cell in a trail of cells) and if so, the cell migrates, leading to sprout elongation with an assigned sprout growth rate. Overall, sprout elongation is due to cell proliferation, migration, and stretching. If the cell is not a tip cell, it has a chance to become a tip cell based on a certain probability: the model generates a random number, and if the random number is greater than the given probability threshold, then this cell becomes a tip cell and starts migrating instead of proliferating. If the random number is not greater than the threshold for the probability of becoming a tip cell, the model then checks whether at this time, the cell is ready to proliferate based on whether the simulated time has reached the assigned cell doubling time. The cell's doubling time is defined by the inverse of its cell proliferation rate (refer to above section on the cell proliferation module). If the check is yes, the cell divides to generate a daughter cell that is assumed to inherit all the properties from the mother cell. If the elapsed simulation time is not enough for the cell to proliferate, it remains quiescent. The model repeats this process at the next time point until the end of the simulation. We update results every hour, which corresponds to the time scale over which cellular responses are usually studied in vitro. , , This model is implemented using MATLAB, and model details are provided in the supplementary materials. Because we only focus on the number of sprouts and sprout lengths, which are well‐studied metrics in in vitro sprouting assays, we do not consider spatial effects for extracellular protein concentrations, cells, or sprouting directions in this study.

Initial state

The model is initialized with a specific number of endothelial cells (N int) in a spheroid. Each of the initial cells are assigned a cell proliferation rate (refer to above section on the cell proliferation module), sprout growth rate, and a probability of forming a new sprout. The average values of the sprout growth rate and probability of forming a new sprout for a population of the endothelial cells are calculated (described in subsequent sections below in detail). The sprout growth rate and a probability of forming a new sprout for each cell were chosen from normal distributions, where the mean and standard deviation for each is based on the calculated values for the population. The mean (μ) is taken as the calculated average probability of forming a new sprout and average sprout growth rate, and a standard deviation (σ) is set to capture 99.7% of the possible values given the range of μ ± 100%μ and μ ± 25%μ (i.e., μ ± 3σ) for each cell, respectively.

During the first iteration, each of the initial cells is assigned a cell proliferation rate, a sprout growth rate, and a probability to become a tip cell, and then all decisions follow the flowchart shown in Figure 2. Each cell can only become a tip cell once, and it stops proliferating if it becomes a tip cell.

Sprouting

Endothelial sprouting is dependent on cell proliferation, migration, and elongation. To quantify endothelial sprouting, we consider the formation of a new sprout and the growth of existing sprouts, which are determined by the probability of forming a new sprout (p) and sprout growth rate (r), respectively. We note that the sprout growth rate here is the net rate of sprout elongation caused by cell proliferation, migration and elongation.

Based on the assumptions made above, the average probability of forming a new sprout p from a cell during a 1‐h period, ∆t, is: where p is the basal probability rate of forming a new sprout when no stimuli are applied, which is not affected by pERK or pAkt levels, k and k are the probability constants, and V max_ and V max_ are the maximum probability rates of sprout formation driven by maximum pERK and pAkt in an hour, respectively. For simplification, we use V and V to represent k × V max and k × V max respectively. Km and Km are the maximum pERK and pAkt levels that produce the half maximal of the probability rates of forming a new sprout, respectively. n is the Hill coefficient for the probability rate of forming a new sprout. The baseline values are provided in Table S1 and supplementary materials.

The average sprout growth rate (r) is where r is the basal sprout growth rate when no stimuli are applied, which is not affected by pERK or pAkt levels; k and k are the sprout growth rate constants; and V max and V max are the maximum rates of sprout length increase driven by maximum pERK and pAkt, respectively. For simplification, we use V and V to represent k× V max_ and k× V max_ , respectively. Km and Km are the maximum pERK and pAkt levels that produce the half maximal of the sprout growth rates, respectively. n is the Hill coefficient for the sprout growth rate. Note that the sprout growth is a net result of sprout formation, growth, degradation, and anastomosis. The baseline values are provided in Table S1 and supplementary materials.

Next, we used these three parameters (r, p, and r) to characterize endothelial cell sprouting: number of sprouts (NS), total sprout length (TL), and average sprout length (AL) in a certain period of time. The number of sprouts is determined by counting the number of tip cells predicted by the model. The total sprout length is the summation of all sprout lengths. The average sprout length is calculated as the total sprout length divided by the number of sprouts.

Constraint on probability of forming a sprout

It has been reported that tip cells activate Notch signaling and prevent neighboring cells from becoming a tip cell. , , Thus, we applied a constraint in our model to account for the effects of lateral inhibition. We adapted a rate constant P max (5 × 10−4 μm−1 h−1), which determines the maximum probability of sprout formation per unit time and vessel length in a rat corneal assay. To match HUVEC data, , we adjusted the P max to 10−3 μm−1 h−1. Thus, P max is used to limit the maximum number of sprouts at every time step.

Model outputs

In vitro spheroid assays are usually three dimensional; however, the TL, NS, and AL from experimental data are usually obtained by quantifying the sprouts from a two‐dimensional image of a focal plane. In order to make comparisons to experimental data, we calculated the TL, NS, and AL on a focal plane by assuming the sprouts are uniformly distributed on a spheroid since our model does not consider spatial effects. Thus, the TL and NS are scaled to the ratio of the number of cells on the focal plane and the total number of cells. Also, the mean values of TL, NS, and AL from experimental data are usually calculated for several randomly selected spheroids per experimental group. To compare with experimental data, we calculated the average TL, NS, and AL on a focal plane for ten simulations per condition.

Sensitivity analysis

Before fitting the model to experimental data, we first performed a sensitivity analysis to identify the parameters that significantly influence the model outputs, using the extended Fourier Amplitude Sensitivity Test (eFAST)  method. Since the initial concentrations and parameters involved in the ERK‐Akt model are adapted from previous work, we used the median fitted values and held them constant during the sensitivity analysis. All remaining model parameters were varied simultaneously within two orders of magnitude above and below the baseline values, where the baseline values are provided in Table S1. In this way, the effects of multiple model inputs on r, r, and r in the cell proliferation module and r, r, r, p, p and p in the sprouting model were computed. Specifically, the eFAST method gives the total sensitivity indices, “Sti,” which can range from 0 to 1, where a higher Sti index indicates the input is more influential to the output. We calculated the Sti values using eFAST for all the same ligand concentrations as the experimental data that were used for model training. The highest Sti value (Sti max) across all of the concentrations was selected to represent the sensitivity index for each parameter.

We also performed eFAST for the trained and validated model to identify potential targets for pro‐ and anti‐angiogenic strategies. All parameters and initial concentrations in the ERK‐Akt model were varied simultaneously within two orders of magnitude above and below the baseline values. The fitted variables were held constant at the median values estimated from model fitting. We calculated the Sti values to quantify how all the variables affected rates of cell proliferation, sprout growth, and the probability of forming a new sprout. Based on the effects in influencing TL, low, intermediate, and high levels of FGF and VEGF (Table 1) were selected as representative ligand concentrations. We calculated the Sti values using eFAST for the nine possible combinations of low, intermediate, and high levels of FGF and VEGF stimulation. Again, the Sti max across all the combinations were compared for all the variables.

TABLE 1

Representative low, intermediate, high levels of FGF and VEGF

	Low (ng/ml)	Intermediate (ng/ml)	High (ng/ml)
FGF	0.03	0.1	10
VEGF	0.1	4	25

Data extraction

Data from published experimental studies , , , were used for parameter fitting and model validation. Experimental data from plots were extracted using the grabit function in MATLAB.

Parameterization

Cell proliferation module

Fitting

The initial concentrations and parameters involved in the ERK‐Akt model are adapted from previous work. Five influential variables with Sti values greater than 0.5 were identified by performing eFAST in the cell proliferation module (Figure S1A). They were estimated against experimental measurements using Particle Swarm Optimization (PSO) to minimize the objective function (the difference between model predictions and experimental data). PSO starts with a population of initial particles (parameter sets). As the particles move around (i.e., as the algorithm explores the parameter space), an objective function is evaluated at each particle location. Particles communicate with one another to determine which has the lowest objective function value. The objective function for each parameter set was used to identify optimal parameter values. Specifically, we used PSO to minimize the weighted sum of squared residuals (WSSR): where V is the ith experimental measurement, V is the ith predicted value at the corresponding time point, and n is the total number of experimental data points. The minimization is subject to θ, the set of upper and lower bounds on each of the free parameters. The bounds were set to be two orders of magnitude above and below the baseline parameter values, which are listed in Table S1.

The cell proliferation module was fitted 200 times using two datasets: the relative proliferation of HUVECs stimulated by 0.03–1 ng/ml FGF and 0.1–1 ng/ml VEGF for 48 h, compared with the reference FGF and VEGF concentration points of 1 ng/ml, respectively (Figure 3A,B). Note that the simulated initial number of cells are the same as experimental data: West et al. and Jih et al. cultured 104 cells and 5000 cells initially, respectively. The relative change of the HUVEC proliferation was calculated as following: where proliferation(c) is the HUVEC proliferation in response to ligand concentration c, and proliferation(cref) is the HUVEC proliferation at a reference concentration point cref.

FIGURE 3

Model comparison to training and validation data for FGF or VEGF stimulation. (A) Relative change of endothelial cell proliferation for 104 cells cultured for 48 h in response to 0.01–30 ng/ml FGF stimulation compared with the reference FGF concentration of 1 ng/ml. (B) Relative change of endothelial cell proliferation for 5000 cells cultured for 48 h by the stimulation of 0.1–10 ng/ml VEGF, compared with the reference VEGF concentration of 1 ng/ml. (C) Total sprout length induced by 0–64 ng/ml VEGF for 24 h cultured with 500‐cell spheroid initially. (D) The fold change of the average length and number of sprouts induced by 25 ng/ml FGF and 25 ng/ml VEGF for 24 h cultured with 400‐cell spheroid initially compared to the control. Circles, squares, and diamonds in Panels A–C are experimental data. , , , Circles in Panel A, squares in Panels A and B, and diamonds in Panel C are experimental data from West et al., Jih et al., and Heiss et al., respectively. The light yellow circles and light blue squares in Panels A, B are experimental data used for model fitting. The orange circles and squares and dark blue squares are experimental data used for model validation. Curves in Panels A, B and C are the mean model predictions of the 21 and 18 best fits, respectively. Shaded regions show standard deviation of the fits. Solid and dashed bars in Panel D are mean ± standard deviation of Liebler et al. data and model predictions, respectively

Validation

We then validated the model with three datasets not used in the fitting. We predicted results for 1–30 ng/ml FGF‐ and 1–10 ng/ml VEGF‐induced HUVEC relative proliferation on Day 2 (using the reference ligand concentration points of 1 ng/ml) and compared with experimental data, specifically the relative proliferation of HUVECs stimulated by 1–30 ng/ml FGF, 0.5–10 ng/ml FGF, and 1–10 ng/ml VEGF for 48 h.

For all three datasets, , we simulated the experimental conditions without any additional model fitting and compared to the experimental measurements. A total of 21 parameter sets with the smallest errors were taken to be the “best” sets based on the model fitting and validation (Figure 3A,B and S2A and Table S2), and the median values were used for the sprouting model.

Sprouting model

Seven influential variables with Sti values greater than 0.75 were identified by performing eFAST in the sprouting model (Figure S1B). Due to the lack of experimental data, we first estimated all the unknown parameters 500 times by fitting the model to experimental observations: 0–64 ng/ml VEGF‐induced total sprout length for 24 h cultured with a 500‐cell spheroid initially. After model training, we validated the model with another dataset not used in the fitting. Specifically, we used the fold change of the average length and number of sprouts induced by 25 ng/ml FGF and 25 ng/ml VEGF for 24 h cultured with a 400‐cell spheroid initially compared to the control for validation. We simulated the experimental conditions without any additional model fitting and compared to the experimental measurements. A total of 15 parameter sets with the smallest errors were taken to be the “best” sets based on the model fitting and validation (Figures S2B and S3 and Table S3). The non‐influential parameters were held constant at the median of the fitted values, and the seven influential variables were estimated using PSO 300 times using data from Heiss et al. (0–64 ng/ml VEGF‐induced total sprout length). We again compared model predictions to the Liebler et al. data without any additional model fitting. A total of 18 parameter sets with the smallest errors were taken to be the “best” sets based on the model fitting and validation (Figure 3C,D and S2C and Table S4) and were used for all model simulations presented below.

Model availability

The MATLAB files with all model equations and scripts for simulating the model are available at: https://github.com/FinleyLabUSC/Endothelial‐cell‐sprouting‐model.

RESULTS

The fitted hybrid agent‐based model captures the main features of FGF‐ and VEGF‐induced endothelial sprouting characteristics

We developed a hybrid agent‐based mathematical model that describes angiogenic cellular responses in the process of endothelial sprouting driven by integrating molecular signals, pERK and pAkt, upon the mono‐ and co‐stimulation of FGF and VEGF (Figure 1). The model focuses on the endothelial proliferation, new sprout formation, and sprout growth, which are assumed to be dependent on the maximum pERK and pAkt levels. The model parameter values are given in Tables S1–S4. The molecular‐detailed biochemical reaction network that characterizes the MAPK and PI3K/Akt pathways induced by FGF and VEGF is adapted from our previous work, and the parameters and initial concentrations are taken from the median of the fitted values. The newly introduced parameters were estimated by fitting the model to experimental data, as described below.

To identify the influential parameters to the model outputs, r, and r, we performed the eFAST (see Methods for more details) and analyzed the maximal Sti (Sti max) value for the newly introduced parameters (Figure S1A). All five parameters (Table S2) were identified as influential and were estimated by fitting the model to experimental measurements , using PSO (see Methods for more details).

The fitted model shows a good agreement with experimental results (Figure 3A,B). It quantitatively captures the main features of FGF‐ and VEGF‐induced endothelial cell proliferation from experimental observations. , In addition to matching data used for fitting, model predictions were compared to experimental data not used in the model fitting , to validate the model (Figure 3A,B). Although the model slightly underestimates the relative proliferation at low VEGF concentration (Figure 3B), the simulated results are consistent with experimental observations and can capture the plateau behavior at high FGF and VEGF concentrations (Figure 3A,B). The weighted errors for 21 best fits are all approximately 4.1 (Table S2). Also, the estimated values of the fitted variables show good consistency (Figure S2A), and the median values were used for the sprouting model.

We again identified seven variables that are influential to r, and r, and p with Sti values greater than 0.75 in the sprouting model using eFAST (Figure S1B). Due to the lack of experimental data, we first estimated all the unknown parameters by fitting the model to experimental observations showing VEGF‐induced total sprout length. We selected 15 “best” fits that showed good match to experimental observations (Figure S3). The non‐influential parameters were then set at the median of fitted values (Figures S1B and S2B) and the influential parameters were estimated 300 times, from which 18 “best” fits were selected (Figure 3C and S2C). The fitting results can capture the main features of VEGF‐induced endothelial total length for 24 h (Figure 3C). We again compared model predictions to independent experimental data without any additional model fitting. It showed good agreement with the fold change of the average length and number of sprouts induced by 25 ng/ml FGF and 25 ng/ml VEGF for 24 h compared with the control (Figure 3D). The weighted errors for 18 best fits are 0.001–0.013 (Table S4). Also, the estimated values of the fitted variables show good consistency (Figure S2C). The ability to predict experimental data not used for estimating the model parameters suggests the model is reliable to make new predictions.

The type and concentration of ligand, length of growth factor stimulation, and initial number of cells impact endothelial sprouting

Mono‐stimulation

We first compared the effects of four inputs: type of ligand (FGF and VEGF), concentration of the ligand (low, intermediate, and high) (Table 1), length of growth factor stimulation (1–3 days), and initial number of cells (250, 500, and 750 cells) on sprouting characteristics (TL, NS, and AL) (Figure 4). Note that AL equals TL/NS, providing the relative change in TL compared to NS. A higher AL suggests that the TL is due to the growth of existing sprouts, whereas a lower AL indicates that the formation of new sprouts plays a more important role in TL. Also, we take low, intermediate, and high levels of FGF and VEGF based on their effects in inducing TL as representative concentrations to study the effects of ligand concentration. Note that for the high levels of FGF and VEGF, the total sprout length is at a plateau level.

FIGURE 4

Predicted sprouting responses stimulated by single agents. Response to FGF stimulation, left panels: Predicted TL (μm) (A), NS (C), and AL (μm) (E) stimulated by low, intermediate, and high levels of FGF. Response to VEGF stimulation, right panels: Predicted TL (μm) (B), NS (D), and AL (μm) (F) stimulated by low, intermediate, and high levels of VEGF. 250‐, 500‐, and 750‐cell spheroid sprouting responses when simulated for 1, 2, and 3 days. Bars are mean model prediction + standard deviation of 18 best fits

Type of ligand (FGF and VEGF)

Generally, FGF induces greater sprouting responses than VEGF in terms of TL (Figure 4A,B) and NS (Figure 4C,D) at the same concentrations. On Day 1, FGF shows slightly higher sprouting responses compared to VEGF, and the differences between FGF‐ and VEGF‐induced TL (Figure 4A,B) and NS (Figure 4C,D) increase with time. This indicates that FGF‐induced TL and NS increase faster than VEGF‐induced TL and NS, which results in greater sprouting responses induced by FGF in the long‐term, compared to VEGF. In addition, FGF‐induced AL is higher than VEGF for all concentrations on Day 1 (Figure 4E,F). On days 2 and 3, FGF showed higher effects in AL at low to intermediate concentrations (Figure 4E), while VEGF showed higher effects in AL at a high concentration (Figure 4F). While increasing FGF and VEGF concentration both lead to greater NS and TL, a higher FGF concentration causes a greater increase in NS compared to TL, which causes a lower AL (Figure 4A, C, and E). The opposite is true for VEGF, where a higher concentration causes a greater increase in TL compared to NS, producing a higher AL (Figure 4B, D, and F). Thus, increasing concentrations of FGF and VEGF affect the average length of sprouts in different ways.

Concentration of the ligands (low, intermediate, and high)

FGF‐ and VEGF‐induced TL, NS, and AL are dose‐dependent (Figure 4). Specifically, FGF‐ and VEGF‐induced TL (Figure 4A,B) and NS (Figure 4C,D) increase with increasing FGF or VEGF concentration. Also, VEGF‐induced AL on Day 1 increases with the increase in VEGF concentration (Figure 4F). Moreover, FGF‐induced AL has a biphasic dose response over all 3 days of stimulation, and VEGF‐induced AL shows a biphasic dose response for days 2–3, respectively (Figure 4E,F). The dose response of the AL for FGF is U‐shaped (Figure 4E), while the dose response for VEGF is an inverted U (Figure 4F). This is caused by the difference in the relative change in TL compared with NS induced by different FGF and VEGF concentrations. As explained above, at high FGF concentrations, the increase in TL is contributed most by the increase in the formation of new sprouts rather than the growth of existing sprouts, while VEGF showed opposite effects.

Length of growth factor stimulation (1–3 days)

TL, NS, and AL increase with the increase in length of cell stimulation (Figure 4). Note that since we are studying responses within 3 days, we assume that anastomosis has not happened yet or is at a minimum level. Otherwise, we would expect that the AL reaches a plateau and might even decrease in a long term.

Initial number of cells (250, 500, and 750 cells)

TL and NS increase with the increase in the initial number of cells, but AL is independent of the initial number of cells (Figure 4). Since the characteristics of endothelial sprouting are qualitatively similar among the groups of different initial number of cells, we take 500 cells as a representative initial number of cells to investigate the effects of FGF and/or VEGF stimulation for the remainder of this study.

Co‐stimulation

We next studied the effects of FGF and VEGF co‐stimulation in endothelial sprouting. Generally, we found that similar to mono‐stimulation, TL, NS, and AL are dose‐dependent for both FGF and VEGF. Specifically, TL and NS increase with the increase in FGF or VEGF concentrations, while AL has a biphasic dose response (Figure 5). Also, TL, NS, and AL increase with the increase in length of cell stimulation (Figure 5). Furthermore, we found that FGF plays a dominant role in the combination effects in endothelial sprouting as co‐stimulation exhibits the features of FGF‐induced sprouting (Figures 4, 5, and S4). Specifically, the effects of co‐stimulation on TL, NS, and AL are more sensitive to FGF concentration change compared to VEGF, which shows no obvious differences between low, intermediate, and high VEGF levels (Figure 5). Also, TL, NS, and AL induced by co‐stimulation are approximately the same level as FGF stimulation alone, as the ratios of combination effects, relative to FGF mono‐stimulation, are approximately equal to one (Figure S4A–C left). In comparison, TL, NS, and AL induced by co‐stimulation are greater than VEGF stimulation alone, as the ratios of combination effects, relative to VEGF mono‐stimulation, are greater than one (Figure S4A–C right).

FIGURE 5

Predicted sprouting responses in response to FGF and VEGF co‐stimulation. Co‐stimulation of FGF‐ and VEGF‐induced TL (μm) on Day 1 (A), Day 2 (D), and Day 3 (G); NS on Day 1 (B), Day 2 (E), and Day 3 (H); and AL (μm) on Day 1 (C), Day 2 (F), and Day 3 (I)

In summary, endothelial sprouting is ligand‐ and dose‐dependent and has different short‐term and long‐term responses. In addition, the initial number of cells is important in the sprouting process. Also, FGF plays a dominant role in the effects of FGF and VEGF co‐stimulation on endothelial sprouting. Moreover, the predicted effects of the co‐stimulation by FGF and VEGF on endothelial sprouting were not significantly greater than FGF mono‐stimulation alone.

The cell proliferation and sprout growth of existing sprouts are predicted to be more important in the sprouting process

To understand the contributions of cellular behaviors in the process of endothelial sprouting, we next investigated the effects of FGF and VEGF on r, and p. NS is a result of the number of cells and the probability of a sprout formation, which are determined by r and p for a certain number of cells initially present. TL is a result of NS and the growth of the sprouts determined by r. And AL equals TL/NS. Thus, understanding how the rate of cell proliferation, rate of sprout growth, and the probability of sprouting depend on growth factor concentration gives insight into the observable sprouting features (number of sprouts, total sprout length, and average sprout length).

We studied cell proliferation, sprout growth, and the probability of sprouting in response to FGF and VEGF mono‐stimulation. We found that FGF and VEGF mono‐stimulation both show sigmoidal dose response curves for cell proliferation (Figure 6A) and sprout growth (Figure 6B). Additionally, FGF induces faster cell proliferation and sprout growth than VEGF at the same concentrations. This is because FGF‐induced r and r are higher than VEGF‐induced r and r, respectively, at the same ligand concentrations (Figure 6A,B). In addition, the large variation in p suggests that the exact value for the probability of forming a new sprout does not significantly affect the model predictions (Figure 6C). This indicates that NS is more dependent on the number of cells compared to the probability of sprout initiation. It also suggests that the number of cells/cell proliferation and the growth of existing sprouts are more important in contributing to the TL over time than the chance of forming a new sprout. Thus, r and r are the main focus in the remainder of this study.

FIGURE 6

Predicted r, r, and p in response to mono‐ and co‐stimulation of FGF and VEGF. Effects of mono‐stimulation of FGF (yellow) or VEGF (blue) on r (A), r (B), and p (C). Effects of co‐stimulation of FGF and VEGF on r (D), r (E), and p (F). Curves in Panels A–C are the mean model predictions of 18 best fits. Shaded regions show standard deviation of the fits

Similar to mono‐stimulation, we next investigated the effects of FGF and VEGF co‐stimulation in the rate of cell proliferation, rate of sprout growth, and the probability of sprouting. We found that r, and p increase as the ligand concentration increases (Figure 6D–F). Moreover, FGF is dominant in the combination effects, as the co‐stimulation exhibits the dose‐dependent features and approximately the same magnitude in inducing r, and p as in response to FGF mono‐stimulation. Specifically, r and r are more sensitive to FGF concentration and relatively independent of the VEGF concentration (Figure 6D–E). Also, the values of r, and p are only slightly higher with FGF and VEGF co‐stimulation, compared to FGF mono‐stimulation; while with co‐stimulation, r is significantly higher, and the values of r and p are also much higher compared to VEGF mono‐stimulation (Figure S5). It suggests that FGF is dominant in promoting cell proliferation, sprout growth, and the probability of sprouting in the combination effects. It is also consistent with the dominant role of FGF in the combination effects in the sprouting characteristics specifically TL, NS, and AL observed in the previous section.

In summary, r and r are more important in endothelial sprouting compared to p. Moreover, FGF plays a dominant role in the combination effects.

The MAPK and PI3K pathways contribute to cell proliferation, sprout growth, and probability of sprouting in different ways

We next studied the contributions of MAPK and PI3K pathways in r, and p to gain insight into the mechanisms of the endothelial sprouting process in response to mono‐ and co‐stimulation of FGF and VEGF.

Cell proliferation

FGF and VEGF stimulate cell proliferation in different ways: FGF‐induced cell proliferation is dominated by pERK, while VEGF‐induced cell proliferation is promoted by pAkt (Figure 7A,B). The dose‐dependent feature of cell proliferation induced by FGF is mostly contributed by ERK phosphorylation, as the rate of cell proliferation influenced by pERK exhibits the dose‐dependent feature of the overall cell proliferation rate, while the rate of cell proliferation influenced by pAkt is independent of the FGF concentration (Figure 7A). In addition, at low FGF concentrations (<0.12 ng/ml), Akt phosphorylation plays a more important role in the cell proliferation rate, compared to its impact at high FGF concentrations (Figure 7A). As the FGF concentration increases, the impact of r increases and eventually surpasses the influence of r (or Akt activation) at higher FGF concentrations (Figure 7A). In contrast, the cell proliferation behavior and dose‐dependent feature induced by VEGF are mostly contributed by Akt phosphorylation, while ERK phosphorylation shows a negligible contribution (Figure 7B). Moreover, the effect of basal cell proliferation minimally contributes to the response induced by FGF or VEGF stimulation (Figure 7A,B).

FIGURE 7

The contributions of MAPK and PI3K/Akt pathways to r and r in response to FGF and VEGF mono‐stimulation. Contributions of pERK (purple), pAkt (green), and basal (gray) for FGF‐induced rates of cell proliferation, r (A) and sprout growth, r (C). Contributions of pERK (purple), pAkt (green), and basal (gray) for VEGF‐induced cell proliferation, r (B) and sprout growth, r (D)

Sprout growth

The MAPK and PI3K pathways contribute differently to FGF‐ and VEGF‐induced sprout growth (Figure 7C,D). First, we found that the basal r (6.87 μm/h) plays a major role in the overall r, accounting for 42% and 63% of the FGF‐ and VEGF‐induced r at their plateau levels, respectively (Figure 7C,D). In addition to the basal r, at low FGF concentrations (<0.03 ng/ml), Akt phosphorylation also plays a substantial role in r (Figure 7C). As the FGF concentration increases, the impact of r increases. Eventually, the contributions of pAkt (height of the green area) and pERK (height of the purple area) to the sprout growth rate plateau at 4.07 μm/h and 5.52 μm/h (Figure 7C), which are the mean values of fitted VpAkt_sg and V , respectively. When VEGF is lower than 3 ng/ml, basal r is dominant in VEGF‐induced r (Figure 7D). As VEGF concentration increases, r increases and plateaus at 4.07 μm/h (V), and r is negligible in VEGF‐induced r (Figure 7D).

We then investigated the contributions of the MAPK and PI3K pathways in r and r upon co‐stimulation with FGF and VEGF. We found that ERK and Akt activation contribute differently to cell proliferation and sprout growth, and FGF plays a dominant role of in combination effects in cell proliferation and sprout growth.

First, r, and r stimulated by a combination of FGF and VEGF mirror the corresponding responses stimulated by FGF alone (Figures 6A,D, 7A and S6A,B). Specifically, the dose‐dependent manner of the cell proliferation rate induced by FGF and VEGF co‐stimulation is mostly impacted by r, while r is independent of the FGF or VEGF concentrations (Figure S6A,B). Also similar to FGF mono‐stimulation, Akt phosphorylation plays a more substantial role in the cell proliferation rate at FGF concentrations lower than 0.12 ng/ml, while r increases as the FGF concentration increases and surpasses the influence of r at higher FGF concentrations (>0.12 ng/ml) (Figure S6A).

The dose‐dependent manner of sprout growth rate induced by FGF and VEGF co‐stimulation is mostly influenced by r (Figure S6C), while r is independent of FGF or VEGF concentrations (Figure S6D). In addition, Akt phosphorylation is more important in r at FGF concentrations lower than 0.03 ng/ml, and the impact of r increases as the FGF concentration increases (Figure S6C,D). These predictions show the same feature as FGF‐induced sprout growth (Figure 7C), indicating the dominant role of FGF on sprout growth.

ERK pathway regulates vessel network mainly via regulating cell proliferation and NS, while Akt pathway mainly affects vessel network via regulating sprout growth

We applied the model to identify the parameters and initial concentrations in the ERK‐Akt model that are influential to the cell proliferation rate and sprout growth rate. This allows us to gain mechanistic insight into how to modulate endothelial sprouting by targeting intracellular signaling pathways. We identified the influential variables by performing eFAST and evaluating the calculated Sti max values. The variables (Table 2) that have Sti max values greater than 0.3 are identified as influential and considered as potential targets for pro‐ and anti‐angiogenic strategies. This analysis predicts that ERK and MEK influence r and r, while Akt and Ptase2 influence r. We note that no parameter was identified as influential to the probability of sprouting, further justifying our primary focus on the rates of cell proliferation and sprout growth. The model predictions suggest that the MAPK pathway is more influential to r, and MAPK and PI3K/Akt pathways are both influential to r.

TABLE 2

Influential parameters affecting rates of cell proliferation (r) and sprout growth (r)

rcp	rsg
ERK MEK	Akt ERK Ptase2 MEK

However, eFAST only tells us that those variables are influential, but the information of how those variables influence model outputs are limited. Therefore, we varied two representative influential variables, ERK and Akt, by 0.1‐ and 10‐fold and predicted the r and r (Figure 8), as well as TL, NS, and AL (Figures 9 and 10) compared to the baseline model predictions. We found that the strategies to modulate endothelial sprouting are context‐dependent, and ERK and Akt pathways regulate vessel network differently.

FIGURE 8

Predicted representative targets for modulating r and r. Predicted r (A) and r (D) from baseline model. Predicted r (B) and r (E) when ERK is varied by 0.1‐ (left) and 10‐fold (right). Predicted r (C) and r (F) when Akt is varied by 0.1‐ (left) and 10‐fold (right)

FIGURE 9

Predicted effects of varying ERK. Predicted TL (A), NS (B), and AL (C) when ERK is varied by 0.1‐ (left) and 10‐fold (right), compared with baseline model predictions (middle) on days 1–3 (i–iii)

FIGURE 10

Predicted effects of varying Akt. Predicted TL (μm) (A), NS (B), and AL (μm) (C) when Akt is varied by 0.1‐ (left) and 10‐fold (right), compared with baseline model predictions (middle) on days 1–3 (i–iii)

Varying ERK

In Figure 8A,B, and S7A, upregulating ERK is predicted to be effective in promoting cell proliferation at low to intermediate but not at high FGF concentrations. In contrast, downregulating ERK is predicted to be more effective in inhibiting cell proliferation at intermediate to high FGF concentrations (Figure 8A,B and S7A). In addition, upregulating ERK is more effective in promoting sprout growth at low FGF concentrations (Figure 8D‐E and S7A). Downregulating ERK is most effective in inhibiting sprout growth at intermediate and high FGF, particularly in combination with low VEGF concentration (Figure 8D,E and S7A).

We next investigated how finitely varying ERK affects endothelial sprouting, specifically, TL, NS, and AL on days 1–3 (Figures 9 and S7C, E, G). The model showed no obvious effects of increasing or decreasing ERK by 10‐fold on TL, NS, and AL on Day 1 (Figures 9Ai, Bi, and Ci, and S7C). In addition, increasing the ERK level is effective in promoting TL at low FGF level and NS at low to intermediate FGF levels but has no obvious effects at high FGF concentration on Day 2 (Figures 9Aii, Bii, and S7E) and Day 3 (Figures 9Aiii, Biii, and S7G). Also, decreasing ERK is effective in inhibiting TL at intermediate FGF level and NS at intermediate and high FGF levels, but this strategy is not very promising at low FGF level on Day 2 (Figures 9Aii, Bii, and S7E) and Day 3 (Figures 9Aiii, Biii, and S7G). Furthermore, the increase in TL is less than the increase in NS, which causes a decrease in AL (Figures 9Cii, Ciii, and S7E, G) since AL equals TL/NS.

Varying Akt

In Figures 8A, C, and S7B, we did not observe obvious effects in r when increasing or decreasing Akt by 10‐fold. However, downregulating Akt is predicted to be effective in inhibiting sprout growth at all FGF and VEGF combinations, while upregulating Akt does not have obvious effects in r (Figures 8D, F, and S7B).

We then investigated how finitely varying Akt affects TL, NS, and AL on days 1–3 (Figures 10 and S7D, F, H). We found that increasing Akt level has no obvious effects on TL, NS, or AL at any of the simulated combinations of FGF and VEGF on days 1–3 (Figures 10 and S7D, F, H). However, decreasing Akt level leads to a decrease in AL and TL at all combinations of FGF and VEGF on days 1–3, but no obvious effects in NS were observed when the Akt level was decreased (Figures 10 and S7D, F, H).

In summary, targeting the ERK pathway is predicted to control the vessel network mainly via regulating cell proliferation and NS, while targeting the Akt pathway mainly influences the vessel network via regulating sprout growth. In addition, the effects of the molecular signaling pathways on sprouting are not obvious after 1 day of stimulation but is more effective in a long‐term (days 2–3). Overall, our model can identify potential effective pro‐ and/or anti‐angiogenic targets and predicts the effects of perturbing those targets under different conditions.

DISCUSSION

We developed a hybrid agent‐based model characterizing the endothelial sprouting process driven by integrating molecular signals, pERK and pAkt, upon the mono‐ and co‐stimulation of two pro‐angiogenic factors FGF and VEGF. The intracellular signaling model of ERK and Akt activation in response to FGF and VEGF stimulation in endothelial cells was adapted from our previous work. The endothelial sprouting process was modeled by assuming the cellular responses (cell proliferation, sprout growth, and the probability of forming a new sprout) are driven by pERK and pAkt, following Hill functions. Unknown parameters were estimated by fitting the model to experimental data. Additionally, we validated the model using a separate set of data.

The fitted model predicts the TL, NS, and AL upon stimulation by FGF and VEGF, alone and in combination, on days 1–3. We particularly focus on TL, NS, and AL because they are metrics examined most often in in vitro studies. , , , , The model predicts that the type and concentration of ligand, length of cell stimulation by the ligand, and initial number of cells are important in endothelial sprouting. Importantly, the model provides quantitative insight into how these inputs influences sprouting. The predicted dose‐dependent sprouting responses (TL induced by FGF and VEGF in Figure 4) are consistent with experimental data, which shows dose‐dependent total tubular length induced by FGF and VEGF. Also, in the same study, 0.1 ng/ml FGF exhibits approximately the same level of increase in proliferation and migration for HUVECs as 25 ng/ml VEGF after 24 h, which agrees with our model prediction that FGF induces greater sprouting responses than VEGF at same concentrations (Figure 4). Furthermore, our model suggests that FGF promotes greater sprouting responses in the long‐term compared with VEGF (Figure 4). In addition, the model predicts that FGF plays a dominant role in the combination effects in endothelial sprouting, and co‐stimulation of FGF and VEGF only slightly increases TL, NS, and AL, compared to FGF simulation alone within 3 days (Figures 4, 5, and S4). Also, the effects of low FGF concentration in combination with high VEGF concentration in TL, NS, and AL show no obvious difference compared to VEGF mono‐stimulation within 3 days (Figures 4, 5, and S4). This prediction is consistent with experimental observations showing no significant increase in HUVEC proliferation in 72 h, migration in 8 h, or total tubular length in 24 h stimulated by the combination of 0.1 ng/ml FGF and 25 ng/ml VEGF compared with their mono‐stimulation. Moreover, our model can supplement experimental endothelial sprouting assays to differentiate and analyze the contributions of cellular and molecular responses during the overall sprouting process including cell proliferation, sprout growth and formation, and activation of the ERK and Akt pathways. Specifically, the model predicts that the ERK pathway regulates vessel network mainly via regulating cell proliferation and NS, while the Akt pathway mainly affects vessel network via regulating sprout growth (Figures 8, 9, 10). These predictions are in line with the literature results reporting that ERK is believed to mainly promote cell proliferation and Akt is more important in cell survival , , , , and migration. , ,

Compared to other models that study cellular behaviors, our mechanistic model considers intracellular signaling and quantitatively analyzes cellular responses driven by integrating molecular signals, pERK and pAkt. Our model explicitly examines how cellular behaviors are driven by pERK and pAkt, which are downstream signals that regulate angiogenic cellular responses. Thus, we can apply our model to mechanistically study the roles of intracellular signaling species in affecting endothelial sprouting. Tong and Yuan constructed a computational model to study vessel growth in rat cornea based on assumptions that cellular responses are only dependent on FGF‐bound FGFR, where the probability of sprout formation and the speed of vessel growth are linearly proportional to the fraction of FGFR occupied by FGF. Our model can complement such models to understand intracellular mechanisms that regulate cellular responses. In another study, Norton and Popel constructed a computational model to study vessel growth in tumors and showed that the proliferation rate has a greater effect on the spread and extent of vascular growth compared to migration rate. Our model predictions agree with these previous modeling works, as we predict that cell proliferation and the number of cells are critical factors that contribute to endothelial sprouting (Figure 4) and that varying ERK level seems to be more influential in TL than varying Akt level (Figures 9 and 10). Our work goes further in that the model considers the intracellular signaling and provides mechanistic insight into the signaling factors driving endothelial sprouting.

Our model can be utilized to study the efficiency of pro‐ or anti‐angiogenic therapies. The model predicts that potential strategies to modulate endothelial sprouting are context dependent (Figures 8, 9, 10), and it can identify potential pro‐ and/or anti‐angiogenic targets under different conditions and study their efficacy. The intracellular angiogenic signals such as pERK and pAkt are believed to play important roles in cellular behaviors, especially cell survival, proliferation, and migration. , , , , , , , However, there is limited quantitative understanding of how integrating the intracellular angiogenic signals affect cellular behaviors, especially in cases where approaches to inhibit angiogenesis have counterintuitive effects on pERK and pAkt. For example, it has been shown that the MEK inhibitor PD0325901 upregulates the PI3K pathway signaling. Our model can investigate such experimental observations and predict the overall cellular responses driven by the integrated molecular signals from pERK and pAkt.

This model can be utilized in combination with other modeling frameworks that predict intracellular signaling to provide more mechanistic insight into certain cellular responses. For example, there are models that study pro‐angiogenic signaling, including sphingosine kinase 1 and calcium responses induced by VEGF, which are the downstream signals of ERK1/2, and also angiopoietin‐Tie signaling in endothelial cells, which has been shown to be important in vessel development, permeability, vascular homeostasis, and remodeling. , , , Other models study anti‐angiogenic signaling, such as TSP1‐CD36 signaling that influences endothelial cell apoptosis and antagonizes VEGF‐induced eNOS signaling, , as nitric oxide is a major vasodilator and is important in angiogenesis and vascular permeability. , Still other models characterize sprouting angiogenesis behaviors by considering how Notch signaling in endothelial cells determines the tip cell and vessel branching as reported in other studies. , , , Our model can be combined with these existing models of angiogenic signaling and applied in various way to provide quantitative insight at both molecular and cellular levels.

We acknowledge some limitations in our model. Although many studies reported that pERK and pAkt play important roles in cell survival, proliferation, and migration, the upstream species in the network (e.g., pVEGFR2, pFGFR1, and PI3K) can activate other pathways that are not included in this work but may also contribute to relevant angiogenic cellular responses. , Also, to ease model construction, we excluded VEGFR1 and neuropilin‐1 (NRP1), although their binding with VEGF contributes to angiogenesis. Moreover, PDGF plays an important role in maturing the vessels ; however, to focus on the initiation of the vessels, we only explored the effects of FGF and VEGF in endothelial sprouting. We can incorporate the contributions of these additional species into the model in future studies. In addition, it has been reported that mechanical and biophysical properties of the ECM affect angiogenesis. , , Excitingly, there is published computational work that focuses on the effects of ECM in angiogenesis. For example, Edgar et al. explored the interactions between ECM and microvessels during angiogenesis, and they then studied the effects of matrix density in sprouting angiogenesis. Similarly, our model can be expanded to consider the effects of the ECM. In addition to endothelial cells, other cell types, including pericytes, , , fibroblasts, , , macrophages, , and smooth muscle cells, , play an important role in angiogenesis. There are computational models that studied the interactions between endothelial cells and pericytes, , fibroblasts, macrophages, , and smooth muscle cells. Our model can be expanded to include these relevant cell types and capture their interactions and further gain quantitative insight into the mechanism of the intracellular signaling and resulting cellular functions. Also, mechanical forces, including tissue stiffness and fluid shear stress, are important in regulating angiogenesis as well. , Our model can be incorporated with other models that studied the effects of mechanical forces. For example, Edgar et al. investigated the effects of ECM stiffness in regulating vascular topology, and Koo and co‐workers studied eNOS regulation induced by fluid shear stress to gain a more comprehensive understanding of angiogenesis process.

We also made some assumptions to simplify our model. First, the molecular interactions usually happen on the order of seconds to minutes. However, it usually takes hours or even days to respond on a cellular level. To bridge the difference in time scale, we assumed that the cellular behaviors in short term are driven by the molecular signals within 2 h, specifically the maximum pERK and pAkt levels. Other computational models have made similar assumptions. For example, Adlung et al. correlated the molecular signals, cyclin D2, cyclin G2, p27, and pS6, which were characterized within 60 min, to markers of cell proliferation mCFU‐E, BaF3‐EpoR, and 32D‐EpoR, which were analyzed after 14–20, 62, and 38 h, respectively. In addition, Tong and Yuan assumed cellular responses are only dependent on the fraction of FGF‐bound FGFR, and that ligated FGFR and the cellular response are linearly related, to investigate FGF‐induced angiogenic dose responses in rat corneal pocket assay. They analyzed migration distance, total and average vessel length, and total number of vessels at times from 0 to 120 h, while the time for FGF‐bound FGFR to reach 90% of steady‐state value was reported to be within 25 min. Also, Padera et al. constructed a mathematical model to study FGF‐mediated cellular response by assuming the signal that drives proliferation of cultured cells from the F32 cell line is only dependent on ligand bound receptor signaling complexes. The model was validated with F32 cell proliferation data, which were evaluated by counting the cell number after 72‐h incubation.

Moreover, our main focus is to study the intrinsic properties of endothelial cells responding to FGF and/or VEGF stimulation. Of course, the experimental setup in in vitro studies may involve environmental factors that influence cell behaviors other than our focus, FGF and VEGF. For example, to maintain sufficient nutrient supply, cell culture media is typically changed every 3 days, which would inevitably induce a change in temperature and air composition, affecting cellular behaviors. In addition, cell proliferation, , , migration, , , and sprouting assays , , are usually conducted in a short period of time following cell culture, typically within 3 days. Thus, to reduce the effects of other possible factors, we only studied short‐term cellular behaviors, specifically 3 days, when nutrients and space are still sufficient. We assumed the cell intrinsic properties, specifically the cell proliferation rate, sprout growth rate, and probability rate of forming a new sprout, remain constant within 3 days considering a relatively stable experimental environment. Other computational studies have also made similar assumptions. For example, Norton and Popel investigated the effects of endothelial proliferation and migration rates on vascular growth by simulating vasculature using various time‐independent proliferation and migration rates for times up to 200 days. In addition, Tong and Yuan studied FGF‐induced angiogenic dose responses in rat corneal pocket assay for up to 5 days.

In addition, we assumed the daughter cells inherit the same cell proliferation rate, sprout growth rate, and the probability rate of forming a new sprout from the mother cells, again considering a relatively stable experimental environment within 3 days. Similar assumptions have also been made by other computational work. For example, Roy and Finley built a multiscale computational model to study pancreatic tumor growth with an assumption that the daughter cell inherits all cell properties and the last intracellular metabolic state of the parent cell. Letort et al. built a multiscale agent‐based model assuming a daughter cell inherits its signaling network state from the mother cell, and they applied the model to study the effects of cell heterogeneity in tumor growth in response to tumor necrosis factor treatment.

Last, we assumed that the maximum pERK and pAkt drive endothelial cell proliferation and sprouting follow Hill functions. There are indeed some computational studies that apply a threshold function to decide if a cell is eligible to migrate or proliferate. , However, HUVECs have been reported to respond to low levels of FGF or VEGF stimulation. For example, Bai et al. showed that FGF concentration as low as 0.1 ng/ml significantly increased HUVECs total tubular length compared to control on Matrigel for 24 h. VEGF has been shown to induce a half‐maximal effect on tubule formation on Laponite substrates at 0.01 μg/ml. Also, Wolfe et al. showed that the EC50 value for VEGF‐induced tube length response is 0.67 ng/ml for HUVECs co‐cultured with normal human dermal fibroblasts. In addition, angiogenic cellular responses, for instance endothelial sprouting , and vessel density, have been shown to be FGF and VEGF dose dependent and will reach a plateau if FGF or VEGF concentration is higher than a certain saturation level. Thus, we decided not to use a threshold function for angiogenic cellular responses. We instead applied a Hill function to account for the observed angiogenic cellular responses, even for low growth factor concentrations. We can expand the model as more mechanistic information becomes available. Despite these limitations, our model provides quantitative insight into angiogenic signaling and cellular responses and can be utilized as a framework for future mechanistic studies.

Perspectives

In conclusion, we developed a mathematical model to characterize endothelial sprouting driven by pERK and pAkt in response to the stimulation of two main pro‐angiogenic factors, FGF and VEGF. The model quantitatively studied FGF‐ and VEGF‐mediated cell proliferation, sprout growth, and formation of new sprouts and provided mechanistic insight into endothelial sprouting. The understanding of the regulation of angiogenesis signals on a molecular scale, and further on a cellular level, can better aid the development of pro‐ and anti‐angiogenic strategies.

Table S1‐S4

Click here for additional data file.

Supplementary Material

Click here for additional data file.