Mathematical modeling of solid cancer growth with angiogenesis
- Hyun M Yang^{1}Email author
https://doi.org/10.1186/1742-4682-9-2
© Yang; licensee BioMed Central Ltd. 2012
Received: 28 September 2011
Accepted: 2 February 2012
Published: 2 February 2012
Abstract
Background
Cancer arises when within a single cell multiple malfunctions of control systems occur, which are, broadly, the system that promote cell growth and the system that protect against erratic growth. Additional systems within the cell must be corrupted so that a cancer cell, to form a mass of any real size, produces substances that promote the growth of new blood vessels. Multiple mutations are required before a normal cell can become a cancer cell by corruption of multiple growth-promoting systems.
Methods
We develop a simple mathematical model to describe the solid cancer growth dynamics inducing angiogenesis in the absence of cancer controlling mechanisms.
Results
The initial conditions supplied to the dynamical system consist of a perturbation in form of pulse: The origin of cancer cells from normal cells of an organ of human body. Thresholds of interacting parameters were obtained from the steady states analysis. The existence of two equilibrium points determine the strong dependency of dynamical trajectories on the initial conditions. The thresholds can be used to control cancer.
Conclusions
Cancer can be settled in an organ if the following combination matches: better fitness of cancer cells, decrease in the efficiency of the repairing systems, increase in the capacity of sprouting from existing vascularization, and higher capacity of mounting up new vascularization. However, we show that cancer is rarely induced in organs (or tissues) displaying an efficient (numerically and functionally) reparative or regenerative mechanism.
Background
A total of 562, 875 cancer deaths were recorded in the United States in 2007, and it is estimated that approximately 570, 000 died from cancer in 2010. The overall estimate of approximately 1.53 million new cases does not include carcinoma in situ of any site except urinary bladder, nor does it include basal cell and squamous cell cancers of the skin. Greater than 2 million unreported cases of basal cell and squamous cell skin cancer, approximately 54, 010 cases of breast carcinoma in situ, and 46, 770 cases of melanoma in situ are expected to be newly diagnosed in 2010 [1]. In the world, the estimate for 1990 suggested a total of 8.1 million new cases, divided almost exactly between developed and developing countries, and 5.2 million cancer deaths, about 55% of which occurred in developing countries [2].
Cells of some organs, as the heart, stop proliferation after reaching their size, but others, like skin cells and cells that line body cavities, must proliferate almost continuously to replenish cells that are lost. Cancer arises when within a single cell multiple malfunctions of control systems occur, which are, broadly, the system that promote cell growth and the system that protect against erratic growth (for instance, tumor suppressor gene p53, which name comes from protein of molecular weight 53,000, and RB, from retinoblastoma). Additional systems within the cell must be corrupted so that a cancer cell, to form a mass of any real size, produces substances (such as VEGF - Vascular Endothelial Growth Factor) that promote the growth of new blood vessels. Cellular growth-control systems can be corrupted either when cellular genes are mutated or when proteins produced during a viral infection interfere with control system function. Multiple mutations are required before a normal cell can become a cancer cell by corruption of multiple (about five) growth-promoting systems [3].
With respect to the control systems that protect against cancer, they are classified in two general types: systems that prevent mutations, and systems that deal with mutations once they occurred. For instance, there are mechanisms in cells that can convert toxic by-products of cellular metabolism or carcinogenic substances from the environment into harmless chemicals, preventing DNA damage. In addition, cells have a number of different systems that can detect damaged DNA and fix it, avoiding mutations. However, if the DNA of a cell has been damaged that repair would be impossible, there are systems which trigger the cell to die. Hence, even though there are trillions of cells, only one in three humans will get cancer during his lifetime, and cancer is mainly a disease of the elderly [3].
In order to develop effective treatments, it is important to identify the mechanisms to controlling cancer growth, how they interact, and how they can most easily be manipulated to eradicate (or manage) the disease. Aiming to gain such insight, it is usually necessary to perform large numbers of time-consuming and intricate experiments - but not always. Through the development and solution of mathematical models that describe different aspects of solid tumor growth, applied mathematics has the potential to prevent excessive experimentation and also to provide biologists with complementary and valuable insight into the mechanisms that may control the development of solid tumors [4]. In [5] it is provided a concise history of the study of tumor growth, discussing some of the most influential mathematical models and their relationship to experimental studies, and illustrating how the field of cancer research has evolved due to these interactions between theoretical and experimental approaches.
In this paper we develop a mathematical model to understand the solid cancer growth inducing angiogenesis. We determine the steady states of the model and the stability of the equilibrium points are performed. We simulate the model to assess the appearance of cancer. Treatments aiming the control of cancer are not dealt here (see [6] and [7] for tumor control).
Methods
Our objective is the development of a simple mathematical model to describe the cancer growth dynamics inducing angiogenesis. The model does not include important physiological processes like the diffusion of oxygen into a solid where it is consumed by metabolic processes, the outward diffusion of lactic acid from a solid which produces it by metabolic processes and the diffusion of oxygen away from a blood vessel into a region with an oxygen debt. First, we describe an organ of human body free of cancer and, then, cancer cells appear in this organ.
The normal cells of an organ (the concentration at time t is designed as C) are feed by an existing vasculature formed by endothelial cells (the concentration at time t is designed as E). Adult endothelial cells and cells of an organ of human body are normally quiescent apart from certain developmental processes (e.g., embryogenesis), and proliferation of these cells aims to replenish losses (for instance, cells dying due to aging and wound). The normal cells and surrounding vascular networks (epithelial cells) are governed by the logistic growth, with intrinsic growth rates α_{1} and α_{2}, respectively, and they are under mortality rates μ_{1} and μ_{2}, respectively. Normal cells produce substances to promote the growth of blood vessels in a suitable network to attend their needs. The extension of vascularization (blood vessels network) depends on the size and function of the organ, and conversely, the blood vessels network determines the activity of the organ. In other words, the size of an organ could be determined by the surrounding network of blood vessels, that is, k_{1}(E), as k_{2} depends on C, k_{2}(C), where k_{1} and k_{2} are the carrying capacities of the organ (normal cells) and the extension of blood vessels network (epithelial cells), respectively. In a cancer free state, normal cells (C) and blood vessels (E) are in a steady state, where the amount of normal cells depends on the feeding network, as well as the amount of epithelial cells depends on the size and function of organs that must feed.
Let us suppose that a series of accumulating mutations in normal cells (one single cell may be sufficient) corrupted systems that promote cell growth and that protect against erratic growth. Also systems that control the production of substances which induce new blood vessels are corrupted in the cancer cells. These sequential events initiate a rapid growth of cancer due to angiogenesis. Angiogenesis is the process by which new blood vessels develop from an existing vasculature, through endothelial cell sprouting, proliferation, and fusion [8]. The angiogenesis links the relatively harmless avascular and the potentially fatal vascular growth phases of the tumor [9].
Let us consider what occurs in cancer growth: (1) there is an initial amount of cancer cells (the model does not deal with the origin of cancer cells) appearing in a completely health organ; (2) these cells induce the formation of pre-vascular cells from existing vascular network; (3) after a period of time, the formation of new blood vessels is initiated; and (4) both normal and cancer cells compete for nutrients and space (proteins, oxygen, etc.).
For cancer cells (the concentration at time t is designed as T), the growth, which is likely to be limited by energetic constraints, is limited by the carrying capacity k_{3} and the surrounding new vessels originated by angiogenesis (the concentration at time t is designed as A). The new vessels, represented by angiogenesis cells A, have growth depending on tumor cells T and are limited by available space (k_{4}). The growth of blood vessels (E and A) are constrained by physical restrictions, even they are never subject to either oxygen or nutrient deprivation. We assume the existence of an early stage of angiogenesis, called pre-angiogenesis (the concentration of cells at time t is designed as P). This stage corresponds to the release by tumoral cells of VEGF, which starts to diffuse into the surrounding tissue and approach the endothelial cells of nearby blood vessels. Endothelial cells subsequently respond to the VEGF concentration gradient by forming sprouts. From these sprouts, there occur migration and proliferation toward the tumor, and these new vessels are called as angiogenesis (endothelial cells A).
The tumoral cells induce new vascular network from the existing one to feed cancer cells. The pre-angiogenesis cells are formed by the mass action law [10] at a constant rate γ, the rate of sprouts formation, from which new vascularizations occur resulting in angiogenesis cells after an average period of time δ^{-1}, where δ is transfer rate from pre-angiogenesis to angiogenesis cells. The cancer, pre-angiogenesis and angiogenesis cells are under the mortality rates μ_{3}, μ_{4} and μ_{5}, respectively. The cancer and angiogenesis cells are governed by the logistic growth: intrinsic growth α_{3} and ε, which depend respectively on angiogenesis and cancer cells, and carrying capacities k_{3} and k_{4}, respectively. Cancer cells can also grow logistically using existing vascularization at rate ${\alpha}_{3}^{\prime}$. The rates at which normal and cancer cells compete themselves for resources and space are ${\beta}_{1}^{\prime}$ and ${\beta}_{2}^{\prime}$. Both parameters also take into account effects of the environment.
Pre-angiogenesis cells deserve some words. New blood vessels to nourish cancer cells must sprout from pre-existing blood vessels. However, there is an elapse of time from the appearance of sprouts until the complete formation of new blood vessels that effectively feed cancer cells. By doing this, we take into account a time-delay between the formation of new network of blood vessels (angiogenesis cells A) and tumor growth. Notice that the probability of sproutings becoming angiogenesis cells is δ/(δ + μ_{4}). If we let δ → ∞, this phase is negligible (probability one).
Let us consider the normal cells C taking into account the above descriptions. In the absence of cancer cells T, normal cells obey $\frac{d}{dt}C={\alpha}_{1}C\left[1-\frac{C}{{k}_{1}\left(E\right)}\right]-{\mu}_{1}C$. The cancer cells influence negatively normal cells, and this term is described by $\frac{{\alpha}_{1}C{{\beta}^{\prime}}_{1}T}{{k}_{1}\left(E\right)}$ with negative signal. Due to mutation, some normal cells become cancer cells, which can occur continuously. However, we will assume that the mutations occur instantaneously at a given time (cancer cells grow exponentially in the beginning, hence further mutations, being in small amount, can be disregarded), when an amount of normal cells are transferred to cancer cells compartment. For other types of cells, the dynamical equations can be obtained similarly.
In this model, an amount C_{ m }of normal cells are mutated to cancer cells at time t = 0, which is described by the Dirac delta function δ_{ D }(t): it assumes ∞ at t = 0, and 0, otherwise.
Symbols and definitions
Symbols | Definitions |
---|---|
C | Concentration of normal cells at time t |
E | Concentration of epithelial cells at time t |
T | Concentration of cancer cells at time t |
P | Concentration of pre-angiogenesis cells at time t |
A | Concentration of angiogenesis cells at time t |
α _{1} | Intrinsic growth rate of normal cells |
α _{2} | Intrinsic growth rate of epithelial cells |
α _{3} | Intrinsic growth rate of cancer cells |
ε | Intrinsic growth rate of angiogenesis cells |
μ _{1} | Mortality rate of normal cells |
μ _{2} | Mortality rate of epithelial cells |
μ _{3} | Mortality rate of cancer cells |
μ _{4} | Mortality rate of pre-angiogenesis cells |
μ _{5} | Mortality rate of angiogenesis cells |
k _{1} | Carrying capacity of normal cells |
k _{2} | Carrying capacity of epithelial cells |
k _{3} | Carrying capacity of cancer cells |
k _{4} | Carrying capacity of angiogenesis cells |
δ | Transfer rate from pre-angiogenesis to angiogenesis cells |
γ | Epithelial sprouting rate |
β _{1} | Rate of inhibition of normal cells by cancer cells |
β _{2} | Rate of inhibition of cancer cells by normal cells |
We next obtain the steady states of the system of equations (1), and the corresponding characteristic equations in order to establish the stability of the steady states.
Steady states
For obvious reason, this point is not considered here.
corresponds to the absence of cancer. The trivial equilibrium point exists for α_{1} > μ_{1} and α_{2} > μ_{2}. In this model the normal and epithelial cells are not dependent on each other. One possibility of dependence can be done by the carrying capacities k_{1}(E_{0}) and k_{2}(C_{0}).
but it is always satisfied when $\stackrel{\u0304}{T}<{T}_{C}$, since T_{ g }> T_{ C }.
where $\stackrel{\u0304}{T}={T}_{C}^{c}$, and the equilibrium $\stackrel{\u0304}{T}$ does not depend anymore on the interacting parameter β_{1}. Hence for ${\beta}_{1}\ge {\beta}_{1}^{c}$ we have a constant value $\stackrel{\u0304}{T}=\left({\alpha}_{1}-{\mu}_{1}\right)/{\beta}_{1}^{c}$, which implies in constant values for Ē, $\stackrel{\u0304}{P}$ and Ā.
See more discussion below, in the stability analysis.
The equilibrium point ${\stackrel{\u0304}{Q}}^{c}$ corresponds to the case where the cancer cells displaced normal cells. By inspecting the equation for $\stackrel{\u0304}{C}$, the first in equation (4), the normal cells drop out to zero for a sufficiently higher values of β_{1}. In another words, when cancer cells have very higher fitness than normal cells and overcome the competition for resources, they can lead to the exclusion of the normal cells (of course, the cancer diseased person dies before reaching this equilibrium). Notice that ${\stackrel{\u0304}{Q}}^{c}$ can be avoided as the equilibrium point when k_{3} assumes small values, that is, if k_{3} < min{T_{ E }, T_{ C }}, minimum between T_{ E }and T_{ C }. In this situation, the equilibrium value $\stackrel{\u0304}{T}$ is always smaller than k_{3}, and all coordinates of the equilibrium point are positive (Appendix A).
When β_{2} increases, $\stackrel{\u0304}{T}$ decreases, and at a certain value of β_{2}, say ${\beta}_{2}^{th},\stackrel{\u0304}{T}$ assumes zero value. When $\stackrel{\u0304}{T}=0$, the unique feasible equilibrium point is the trivial ${\stackrel{\u0304}{Q}}_{0}$. When the influence of the normal cells and the environment is higher (${\beta}_{2}\ge {\beta}_{2}^{th}$), then cancer can not be maintained.
Cancer onset depends on the interplay among the initial condition T(0) = C_{ m }and interacting parameters β_{1} and β_{2}. Once cancer cells are originated by mutation from several sources, these cells interact with healthy cells and the surrounding environment. While the parameter β_{1} is practically insensitive in the behavior of the dynamics, however, there is a threshold in the parameter β_{2}, above which cancer can not be established. This result shows that cancer is rarely induced in organs (or tissues) displaying an efficient (numerically and functionally) reparative or regenerative mechanism [12], which is the reason for the incidence of cancer rising exponentially with age [13].
Finally, let us analyze the effects of varying the sprouting rate γ and angiogenesis increasing rate ε.
When γ increases, Ē decreases according to the second equation of (4). Let us assume that Ē reaches zero at a certain value named γ^{ c }, and Ē < 0 sinceafter. When Ē = 0, differently from $\stackrel{\u0304}{C}=0$, the unique feasible equilibrium point is the trivial ${\stackrel{\u0304}{Q}}_{0}$. Hence, we expect that Ē > 0 for γ ≥ 0, because of $\stackrel{\u0304}{T}<{T}_{E}=\left({\alpha}_{2}-{\mu}_{2}\right)/\gamma $ is always obeyed, and there is not a critical value γ^{ c }. When we decrease γ, we observe that equation (5) may not have biologically feasible solutions.
Equation (10) has two positive solutions (not shown) when T_{ g }> k_{3} and ε > μ_{5}/k_{3}, where T_{ g }is given by equation (8), while a unique positive solution occurs when T_{ g }= k_{3}. Let us introduce a threshold of ε, named ε^{ th }. We remark that this threshold value appears only for γ = 0. Hence, when ε > ε^{ th }, with ε^{ th }> μ_{5}/k_{3}, we have two positive solutions, named ${T}_{<}^{\epsilon}$ (and the corresponding equilibrium point ${\stackrel{\u0304}{Q}}_{<}^{p}$) and ${T}_{>}^{\epsilon}$ (and the corresponding equilibrium point ${\stackrel{\u0304}{Q}}_{>}^{p}$), which collapse to one at ε = ε^{ th }, and for ε < ε^{ th }there is not positive solution. The equilibrium point ${\stackrel{\u0304}{Q}}^{p}$ represents the origin of angiogenesis from other sources, not from existing vascular system feeding an organ of the human body. Two facts must occur to appearing of ${\stackrel{\u0304}{Q}}^{p}$: sufficiently higher capacity of cancer cells to originate new vascular system from external sources (ε > ε^{ th }), and higher initial amount, at time 0, of vascular cells $\left(T\left(0\right)>{T}_{<}^{\epsilon}\right)$.
This equation clearly shows that a positive solution $\stackrel{\u0304}{T}$ is feasible for sufficiently higher values of γ, that is, $\gamma >{\gamma}_{0}^{th}$, where ${\gamma}_{0}^{th}$ corresponds to ε = 0. Due to the fact that ε increases as γ^{ th }decreases, ${\gamma}_{0}^{th}$ is the upper bound of the threshold of γ. Even the cancer cells do not promote the proliferation of the sprouted endothelial (pre-angiogenesis) cells (ε = 0), if the capacity of endothelial cell sprouting is higher (γ > γ^{ th }), then cancer cells can be maintained.
There are two mechanisms by which the tumor's vasculature is built: (1) as a tumor grows, it excretes tumor angiogenesis factor, which help activate endothelial cells of nearby blood vessels, initiating angiogenesis; (2) bone marrow derived endothelial progenitor cells are mobilized into the blood, by means of long-range signaling, and they are recruited from blood by the tumor by short-range signaling, resulting in vasculogenesis (the vascular endothelium of the nearby vessels become activated and allows the endothelial progenitor cells to extravase and start a cycle of differentiation/division) [14]. In both cases, the newly stimulated endothelial cells (described by parameter γ) and the recruited endothelial progenitor cells (this phenomenon is not considered in the model; however, the case γ = 0 can in some extent be understood as vasculogenesis) enter a stage of clonal expansion and continue to form blood vessels, which is described by the parameter ε.
and g_{1}(T) + g_{2}(T) is given be equation (6). Notice that f(T) differs from that in equation (6) by the factor δ/(μ_{4} + δ). The equilibrium point ${\stackrel{\u0304}{Q}}^{\infty}$ represents the absence of the intermediate pre-angiogenesis, that is, angiogenesis occurs directly and instantaneously from epithelial cells.
Summarizing, the non-trivial equilibrium points ${\stackrel{\u0304}{Q}}^{*}$ (${\stackrel{\u0304}{Q}}_{<}^{*}$ and ${\stackrel{\u0304}{Q}}_{>}^{*}$, when they exist), ${\stackrel{\u0304}{Q}}^{c}$ (for ${\beta}_{1}\ge {\beta}_{1}^{c}$), and ${\stackrel{\u0304}{Q}}^{p}$ (when γ = 0) are biologically feasible, but ${\stackrel{\u0304}{Q}}^{c}$ is not a real cancer equilibrium point, since $\stackrel{\u0304}{C}=0$. Especially, ${\stackrel{\u0304}{Q}}^{p}$ is a cancer equilibrium point only if vasculogenesis can act alone to generate new blood vessels. The analysis of the model showed the existence of thresholds ${\beta}_{2}^{th}$ and γ^{ th }(and, eventually, ε^{ th }, for γ = 0), and a critical value ${\beta}_{1}^{c}$.
We analyzed the steady states of the model described by the system of ordinary differential equations (1) assuming that all parameters of the linear terms (μ_{1}, μ_{2}, μ_{3}, μ_{4}, μ_{5} and δ) plus the intrinsic growth rates α_{1}, α_{2} and α_{3} and carrying capacities k_{1}, k_{2}, k_{3} and k_{4} are fixed values. But, the parameters regarded to the interaction between different cells (β_{1}, β_{2}, γ and ε) are generally unknown, hence they were allowed to vary to analyze broad range of variation in these parameters.
Local stability analysis
The Jacobian $\stackrel{\u0304}{J}$ was obtained using equilibrium equations given in (4) assuming that $\stackrel{\u0304}{C}\ne 0$. When $\stackrel{\u0304}{C}=0$, the element in the first line and first column is ${a}_{11}=\left({\alpha}_{1}-{\mu}_{1}\right)-{\beta}_{1}\stackrel{\u0304}{T}$, and in the third column, a_{13} = 0.
We present analytical results corresponding to the equilibrium points ${\stackrel{\u0304}{Q}}^{0}$ and ${\stackrel{\u0304}{Q}}^{c}$.
The eigenvalues (from Jacobian ${\stackrel{\u0304}{J}}^{0}$) corresponding to the trivial equilibrium point ${\stackrel{\u0304}{Q}}^{0}$ are: λ_{1} = -μ_{5}, λ_{2} = - (μ_{4} + δ), λ_{3} = - (μ_{3} + β_{2}C_{0}), λ_{4} = - (α_{1} - μ_{1}), and λ_{5} = - (α_{2} - μ_{2}). In the last two eigenvalues we used the coordinates given by equation (3). Hence, the trivial equilibrium ${\stackrel{\u0304}{Q}}^{0}$ is locally asymptotically stable whenever this point exists (α_{1} > μ_{1} and α_{2} > μ_{2}), because all eigenvalues become negative.
where j_{1}, j_{2}, j_{3} and j_{4} were defined in equation (13). In appendix B we show that, for big solution T_{>}, A is an M-matrix, and, hence, all eigenvalues associated to ${\stackrel{\u0304}{J}}_{1}^{c}$ have negative real part. Hence, ${\stackrel{\u0304}{Q}}^{c}$, which is biologically feasible but not real cancer, is locally asymptotically stable only if $\stackrel{\u0304}{T}>{T}_{C}$. Therefore, ${\stackrel{\u0304}{Q}}^{c}$ appears when one of the constraints to ${\stackrel{\u0304}{Q}}^{*}$ be biologically feasible is violated: in the equation (4), instead of $\stackrel{\u0304}{C}<0$, we have $\stackrel{\u0304}{C}=0$ for ${\beta}_{1}\ge {\beta}_{1}^{c}$.
The local stability of the non-trivial equilibrium ${\stackrel{\u0304}{Q}}^{*}$ is performed numerically. However, we stress the fact that the equilibrium ${\stackrel{\u0304}{Q}}^{c}$ is an extension of ${\stackrel{\u0304}{Q}}^{*}$ for ${\beta}_{1}\ge {\beta}_{1}^{c}$. Hence, the big solution T_{>} of equation (5) for ${\beta}_{1}\ge {\beta}_{1}^{c}$ must be stable. For this reason we showed the local stability of ${\stackrel{\u0304}{Q}}^{c}$, even though this is not a real cancer equilibrium.
Results and Discussion
In this section we present numerical simulations of the model, which are discussed. The numerical methods used are bisection (to find zeros of polynomials) and 4^{ th }order Runge-Kutta (to solve system of ordinary differential equations) [16].
Numerical analysis of the model
Values of the parameters
Parameters | Fixed values | Alternative values** | Units |
---|---|---|---|
α_{1} | 0.1 | days ^{-1} | |
α_{2} | 0.1 | days ^{-1} | |
α_{3} | 0.2 | 5.0 | [A]^{-1} × days^{ -1 } |
ε | 0.01* | 0.1 | [T]^{-1} × days^{-1} |
μ _{1} | 0.01 | days ^{-1} | |
μ _{2} | 0.05 | days ^{-1} | |
μ _{3} | 0.05 | 0.005 | days ^{-1} |
μ _{4} | 0.01 | days ^{-1} | |
μ _{5} | 0.01 | days ^{-1} | |
k _{1} | 10 | [C] | |
k _{2} | 20 | [E] | |
k _{3} | 5 | 0.1 | [T] |
k _{4} | 1 | 0.2 | [A] |
δ | 0.1 | days ^{-1} | |
γ | 0.01* | 0.02 | [T]^{-1} × days^{-1} |
β _{1} | 0.01* | [T]^{-1} × days^{-1} | |
β _{2} | 0.01* | [C]^{-1} × days^{-1} |
Eigenvalues
Eigenvalues | Corresponding to${Q}_{<}^{*}$ | Corresponding to${Q}_{>}^{*}$ |
---|---|---|
λ _{1} | +0.028 | -0.227 - i 0.024 |
λ _{2} | -0.089 - i 0.060 | -0.227+ i 0.024 |
λ _{3} | -0.089 + i 0.060 | -0.068 |
λ _{4} | -0.064 | -0.0293 - i 0.0059 |
λ _{5} | -0.047 | -0.0293 + i 0.0059 |
The initial amount of normal cells that suffer mutation C_{ m }is crucial to trigger cancer disease. Figure 1 shows that if this amount is below a critical value, that is, T(0) = C_{ m }< T*, then repairing systems act efficiently and cancer does not settle in an organ. However, if the mutated cells surpass the critical value, the repairing systems do not avoid the onset of cancer.
Next, numerical simulations are performed in order to assess the effects of varying parameters β_{1}, β_{2}, γ and ε. In all dynamical trajectories, remember that unstable solution T_{ < }is very small, but T(0) is higher due to the initial conditions supplied to system (1) correspond to the coordinates of the trivial equilibrium ${\stackrel{\u0304}{Q}}^{0}$, and not those of the unstable ${\stackrel{\u0304}{Q}}_{<}^{*}$. Hence, T(0) is not comparable with T_{ < }, one of the coordinates of the break-point ${\stackrel{\u0304}{Q}}_{<}^{*}$. The cancer cells proliferate above the subclinical threshold of 10^{3} cells and reach 10^{9} cells which is the X-ray detectable threshold [11].
Interaction between normal and cancer cells - β_{1} and β_{2}
Direct competition between normal and cancer cells for resources and space occur in order to grow. But there are many factors that affects both populations, like indirect interaction as excretion released by cells, changes in the environment and control mechanisms. The parameters β_{1} and β_{2} take into account these factors also.
Let us assess the cancer cells affecting negatively in the normal cells, by varying β_{1}.
Let us assess the normal cells affecting negatively in the cancer cells, by varying β_{2}.
As β_{1} increases, the initial number of cancer cells (for instance, originated by mutation) needed to trigger a cancer T(0) = C_{ m }is decreased, but smoothly. Also, the normal cells can be displaced by cancer cells for values of β_{1} higher than its critical value, that is, ${\beta}_{1}>{\beta}_{1}^{c}$. With respect to β_{2}, we observed a threshold for ${\beta}_{2},{\beta}_{2}^{th}$, above which cancer can not be settled. Also, as β_{2} increases, the initial number of cancer cells needed to trigger a cancer C_{ m }increases quickly, avoiding the process of cancer disease. Therefore, cancer can be settled in an organ if the following combination matches: better fitness of cancer cells (β_{1} increases), and decrease in the efficiency of the repairing systems (β_{2} decreases).
Recruitment of existing epithelial cells - γ
New network of blood vessels is created by cancer cells to provide nutrients and oxygen to support their growth. This new network depends on the capacity of originating sprouts from the existing blood vessel network, and is described by the parameter γ, which is varied.
In cancer disease, it is expected that a small number of epithelial cells must be recruited in order to build up new vascularization from the sproutings. For small values of γ but higher than the threshold (γ > γ^{ th }), the dynamical trajectories follow those shown in Figure 1. In Appendix C we illustrate the Hopf bifurcation, which behavior is not compatible with real cancer. Considering number of tumor cells, concentration of growth factor and volume of blood vessels feeding the tumor, Agur et al. [18] showed that Hopf bifurcation can not occur if ordinary differential equations are used. But, Hopf bifurcation can occur if time-delay is encompassed. Our model presented Hopf bifurcation, however only in a range of values of parameter γ which is not compatible with biological findings.
Solid tumors need extra source of resources to attend the quick growth of cancer cells. Hence, cancer can be settled in an organ if the capacity of sprouting from existing vascularization is sufficiently higher (γ > γ^{ th }). However, it must not be so higher in order to avoid the death of the cancer diseased person due to normal cells being displaced by cancer cells quickly.
Capacity of building up new vascularization - ε
The appearance of shunts from existing blood vessels to initiate new vascularizations was described by the parameter γ. New network of blood vessels is created by cancer cells to provide nutrients and oxygen to support their growth. After a period of time δ^{-1}, new vessels are built up from the shunts. The capacity of mounting up new vessels by cancer cells is analyzed by varying ε.
Due to the influx γP, which is a linear term, in the equation for A, the parameter ε that describes the mounting up of new vascularization does not present any special behavior, except by the dependency with the initial conditions. However, when γ = 0, there arises a threshold of ε, called ε^{ th }, above which new vascularizations promoted by cancer cells can occur. The case γ = 0 means that cancer cells are supported exclusively by the new network of blood vessels originated from surrounding tissues for instance, and the pre-existing one maintains its function of nourishing exclusively the normal cells.
Angiogenesis is the process by which new blood vessels develop from an existing vasculature, through endothelial cell sprouting, proliferation, and fusion. Hence, angiogenesis create new vascularization from sprouting originated in existing vasculature, which was called pre-angiogenesis. Due to the endothelial cell sprouting promoted by the vascular endothelial growth factor, cancer cells can grow even in the absence of new vascularization (ε = 0), being the size of cancer cells big with higher capacity of mounting up new vascularization (increasing ε).
Discussion
In foregoing section we have used for k_{3} and k_{4} values comparable to k_{1} and k_{2} in order to enhance the results. In other words, cancer related cells are allowed to grow comparable to the size of normal cells, which is not true. In real world, cancer related cells are found in much smaller size, hence k_{3} and k_{4} must be lower than k_{1} and k_{2}, being the relative size depending on the organ of the body.
The initial amount of cancer cells originated from normal cells by mutations plays an important role in the dynamics of cancer growth. When the interacting parameters β_{1}, γ and ε increase, the small solution T_{<} decreases. This decrease in the initial cancer cells necessary to trigger the cancer is mediated by cancer cells that: (1) inhibit and occupy normal cells habitat (β_{1}); (2) produce higher amount of substances which cause new blood vessels to grow (γ); and (3) construct new blood networks to nourish themselves (ε). In opposite way, when the interacting parameter β_{2} increases, the small solution T_{<} also increases. Hence, when normal cells inhibit the growth of cancer cells by many factors as better fitness, increasing repairing action, and induction of apoptosis, then the onset of cancer is avoided.
Another important aspect in the cancer growth is the appearance of thresholds. The interacting parameters β_{2} and γ present threshold values, respectively, ${\beta}_{2}^{th}$ and γ^{ th }. Usually the systems that control the production of substances that induce the formation of new blood vessels to grow operate normally, which have as consequence that cancer cells are unable to recruit the blood to supply their need to continue to proliferate, and they fade out at this early stage. However, cancer cells may begin to produce substances which cause new blood vessels to grow. This phenomenon is characterized by the threshold of γ. The threshold of β_{2} can be understood as the well functioning repairing mechanisms and the low fitness of cancer cells in comparison with normal cells. The threshold of ε arises only for γ = 0, which mimics cancer cells being nourished only by the new network of blood vessels originated from surrounding tissues by vasculogenesis [14].
There are also critical values for parameters β_{1} and γ. Critical value for β_{1} is a mathematical artefact, because it is meaningless biologically (in general k_{3} is very low in comparison with constraint T_{ E }). With respect to γ, there are two critical values, named ${\gamma}_{1}^{c}$ and ${\gamma}_{2}^{c}$. When γ increases, the epithelial cells E decrease, and damped oscillations appear. When γ approaches to ${\gamma}_{1}^{c}$, the oscillations are less damped, and when surpasses ${\gamma}_{1}^{c}$, regular oscillations occur. However, the amplitude of regular oscillations increases as γ approaches to ${\gamma}_{2}^{c}$, resulting for lowest values of E, T and A reaching zero values (see figures in Appendix C). When the lowest values are incapable to trigger new burst of cancer cells, the oscillations cease and the trivial is the attractor. This occurs when γ surpasses ${\gamma}_{2}^{c}$. Again both ${\gamma}_{1}^{c}$ and ${\gamma}_{2}^{c}$ do not bear any biological meanings, because the cure of cancer is due to the elimination of pre-existing network of blood vessels. The biological meaningless sustained oscillations begin at $\gamma ={\gamma}_{1}^{c}$ (the supercritical Hopf bifurcation) and cease at $\gamma ={\gamma}_{2}^{c}$ (the subcritical Hopf bifurcation).
For instance, considering values of Table 2, the column with **, we observe the same behavior than that found in figures shown in Appendix C as γ increases: (1) damped oscillations around ${\stackrel{\u0304}{Q}}^{*}$ for γ = 0.75; (2) regular oscillations around ${\stackrel{\u0304}{Q}}^{*}$ in the interval [0.76,0.80]; and (3) to the trivial equilibrium ${\stackrel{\u0304}{Q}}^{0}$ for γ = 0.81. In the region of limit cycle, we observed that: (1) C oscillates between (8.9, 9.0) for all γ; E, P and A oscillate between (0, M), where M increases as γ increases; and (3) T oscillates between (m, 0.1), where m decreases as γ increases. When m becomes very small, there is not burst of cancer cells, limit cycle is destroyed, and cancer fades out for higher values of γ. Notice that the amplitude of oscillations of normal cells C is not affected, in opposite way of epithelial cells E which drops to near 0.
Agur et al. [18] showed that ordinary differential equations admit Hopf bifurcation if and only if at least one time-delay is introduced in the tumor growth modeling. They then concluded that an appropriate candidate for describing the cancer growth is the alternative that includes time-delays in the tumor proliferation or angiogenesis process. They also concluded that further mathematical research is warranted for exploring time-delays in the biologically realistic domains in the parameter spaces. Our model, however, showed Hopf bifurcation without time-delay (in fact, there is an elapse of time between pre-angiogenesis and angiogenesis cells), and sustained oscillations occur only in biologically not realistic domains in the parameter space.
When the initial conditions, especially T(0) = C_{ m }, are such that the non-trivial equilibrium is attractor, the cancer cells reach the level T_{>}. C_{ m }increases with increasing β_{1}, decreases with β_{2} and ε. Cancer cells can grow and reach higher levels when they affect negatively normal cells (β_{1}), but reach lower levels when normal cells acts as a barrier against them (β_{2}). When γ varies, T_{>} increases in the initial phase (E decreases), and then decreases (E is practically zero). The parameter γ plays an equivalent role of β_{1}, but, restricted only to E, which decreases it dramatically. For lower γ, there is sufficient number of E to increase A, but epithelial cells are exhausted as γ increases, and A decreases (see equation (4)). Finally, T_{>} decreases with ε, which comes out due to relative higher value of γ. The behavior of ε is strongly dependent on γ due to the influx γP in equation for A: for lower γ (also γ = 0), T_{>} increases with ε (see Figure 12(b)).
Conclusions
Many models [5] [19] have already been proposed to describe cancer growth, and some of those models have explicitly considered the spatial dimension [20–22], which has been shown to play a key role in the understanding of various tumor growth processes. Other models considered computational approach [23, 24]. Spatial and computational modelings were analyzed numerically. However, we developed a zero-dimensional model for the initial stages of tumor angiogenetic growth in order to obtain analytical results.
In this model, various cell species are supposed to satisfy Lotka-Volterra growth laws. Neither metastasis nor geometry of the solid cancer were taken into account. The initial conditions (2), supplied to the dynamical systems (1), describe an impulsive system: in a steady state, a perturbation is introduced at a time t = t_{0} = 0 as a form of pulse (Dirac delta function). This pulse mimics normal cells mutating to cancer cells. Cancer cells then promote the mounting of new network of blood vessels to nourish them after an elapse of time δ^{-1}. We introduced the pre-angiogenesis cells P to include the time delay in order to completing functional angiogenesis.
From the model, we conclude that the dynamical trajectories depend on the initial conditions supplied to the system, and also on interacting parameters. The cumulative effects of mutation is essential to originate a cancer cell. This effect is captured by the initial amount of cancer cells originating from normal cells, denoted by C_{ m }. A sufficient number of cells must suffer mutation in order to a concentration of C_{ m }cells really bear all necessary mutations to become effectively cancer cells. Our model is spatially homogeneous, hence the initial number of cancer cells is C_{ m }× V, where V is the volume of an organ of human body. In extremely favorable environmental and individual conditions, this initial number can be one.
Initially, cancer cells always can grow. But they fade out if they are unable to build up new blood vessels in order to supply their needs. The capacity of inducting new vascularization from existing blood vessel network must be efficient (γ > γ^{ th }), and better fitness (increasing capacity of proliferation and capturing nutrients, decreasing mortality, etc.) of cancer cells (β_{1} increases) in comparison with normal cells. Another aspect of cancer growth is corruption of the repairing systems, and in some extend we can think of that normal cells influencing negatively cancer cells play the role of repairing (fixing mutated DNA and inducting apoptosis). The parameter β_{2} measures the efficient action of repairing system, and the effect of decreasing this value result in a higher level of corruption in the repairing system. Since cancer cells can recruit epithelial cells to form new blood vessels, the capacity of proliferation of new vessels (ε) does not present threshold value.
The parameter γ can be thought of deviation of nutrients and oxygen from normal cells to feed cancer cells. This deviation does not affect normal cells because the carrying capacity of normal cells k_{1} does not depend on the size of network of blood vessels. As γ increases, the deviated network plus the new one mounted by the action of angiogenesis effectors nourish the cancer cells. Hence, γ = 0 means that pre-existing network of blood vessels feeds normal cells, while the new network nourishes cancer cells. In this case, the capacity of cancer cells in promoting new vascularization is essential, which must surpass the threshold value ε^{ th }. We obtained biologically feasible non-trivial equilibrium points, that is, the coordinates of the equilibrium points are positively defined. However, due to the simplifications assumed by model, the range of variations of parameters like β_{1} and γ must be restricted. When ${\beta}_{1}>{\beta}_{1}^{c}$ we have an equilibrium point ${\stackrel{\u0304}{Q}}^{c}$ with $\stackrel{\u0304}{C}=0$, while for ${\gamma}_{1}^{c}<\gamma <{\gamma}_{2}^{c}$, we have limit cycle (sustained oscillations) with E ~ 0, but for $\gamma >{\gamma}_{2}^{c}$, we have abrupt increase of E, and the trivial equilibrium is attained. Both results can not be acceptable for cancer growth description.
Some results obtained here can be understood as vasculogenesis (when γ = 0) [14]. The dynamics of normal and cancer cells are similar than that presented in the model proposed by Nani and Freedman [11]. However, we did not take into account the action of immune system, while they did not take into account the angiogenesis. Agur et al. [18] proposed to examine the occurrence of Hopf bifurcation in the clinical context, that is, to check whether or not one can contain tumor growth by imposing time-delays in the processes of neo-vascularization. Our results showed that Hopf bifurcation occurs in biologically not realistic domains in the parameter space.
In a future work we will analyze a model in which the sizes of the normal and cancer cells are allowed to depend on the overall network of blood vessels: normal and cancer cells compete for nutrients provided by the pre-existing blood vessels, while cancer cells have additional source originated from angiogenesis. If we take these effects into account in the model, maybe Hopf bifurcation can be avoided. There are several ways to improving model given in equation (1). One is the dependency of normal cells with the size, which decreases with increasing γ, of the existing vasculature.
where s_{ i }and τ_{ i }are, for i = 0, ···, n, respectively, the i-th proportion of epithelial cells generating sproutings and the time at which sproutings occur. The proportion s_{ i }can depend on T. Vaccination campaigns against viral infections were analyzed considering age interval vaccination (Heaviside function) [25] or a series of pulses (Dirac function) [26].
Another improvement of the model (1) is the introduction of the immune response. This can be done by introducing lymphocyte cells action as Nani and Freedman dealt with [11].
where q_{ i }and τ_{ i }are, for i = 0, ···, n, the i-th drug administration rate and the time interval during which drug is administrated. To be intermittent, we must have q_{2i}> 0 and q_{2i+1}= 0.
In this paper we obtained threshold values, but we did not deal with the effects of controlling mechanisms, which are left to a further work. Briefly, for instance, two parameters can be used in order to control cancer growth. The thresholds of the parameters γ and β_{2}, named γ^{ th }and ${\beta}_{2}^{th}$, can be managed as follows.
Let us suppose that γ and β_{2} have values above and below their respective thresholds. To control cancer, parameters of the model (including dynamics of controls) must be managed to increase γ^{ th }to surpass γ, and to decrease ${\beta}_{2}^{th}$ to fall below β_{2}. Another way to control cancer is acting on the probability of sproutings becoming angiogenesis cells, in order to decrease δ/(δ + μ_{4}), which can be done by decreasing the rate of transformation from pre-angiogenesis to angiogenesis cells δ and/or increasing the mortality rate μ_{4} of pre-angiogenesis cells.
Appendix A: Non-trivial equilibrium point
and there is one negative solution and two positive solutions. In both cases T_{ g }is always a positive root of g(T), and g_{0} is given by equation (A.1). The existence of other positive solutions depends on the relative position of the roots of g_{1}(T) and the coefficients of the polynomials.
- 1.
We do not have negative solutions, because f(T) < 0 and g(T) > 0, for T < 0, which implies that f(T) - g(T) < 0. At T = 0 we have f(0) - g(0) = -g_{0} < 0.
- 2.
We have f(∞) - g(∞) → -∞, even when g(∞) → +∞, because f(T) is fourth degree polynomial, and g(T), third degree.
- 3.
f(T) - g(T) is a continuous function, a fourth degree polynomial, where there is not solution for T < 0, f(0) - g(0) < 0 and f(∞) - g(∞) < 0. Hence, by the intermediate value theorem, we have 0, 2 or 4 positive solutions in the interval (0, ∞).
- 4.
We are searching biologically feasible equilibrium, hence, according to equation (7), $\stackrel{\u0304}{T}$ must be lower than the lowest root T_{ m }of f(T), that is, $\stackrel{\u0304}{T}\in \left(0,{T}_{m}\right)$, where T_{ m }= min {(α_{2} - μ_{2})/γ, k_{3}}, the minimum between (α_{2} - μ_{2})/γ and k_{3}. Therefore, in the interval (0, T_{ m }) we have at most 2 positive solutions.
When there is not any positive solution for f(T) - g(T) = 0, the unique equilibrium point is the trivial ${\stackrel{\u0304}{Q}}^{0}$. In the case of two positive solutions, the small root (T_{ < }) forms the unstable equilibrium ${\stackrel{\u0304}{Q}}_{<}^{*}$, and the big one (T_{>}) forms a possibly stable (let us for simplicity say that it is stable) equilibrium ${\stackrel{\u0304}{Q}}_{>}^{*}$. Between two stable equilibria ${\stackrel{\u0304}{Q}}_{>}^{*}$ and ${\stackrel{\u0304}{Q}}^{0}$ (this always exists and is stable), we have an unstable point ${\stackrel{\u0304}{Q}}_{<}^{*}$. We call as the 'break-point' [28, 29] the unstable equilibrium point ${\stackrel{\u0304}{Q}}_{<}^{*}$, which separates two attracting regions containing one of the equilibrium points ${\stackrel{\u0304}{Q}}_{>}^{*}$ and ${\stackrel{\u0304}{Q}}^{0}$. In another words, we have a hyper-surface obeying, e.g., $h\left(\stackrel{\u0304}{C}\left({T}_{<}\right),\u0112\left({T}_{<}\right),{T}_{<},\stackrel{\u0304}{P}\left({T}_{<}\right),\u0100\left({T}_{<}\right)\right)=0$, generated by the coordinates of the equilibrium point ${\stackrel{\u0304}{Q}}_{<}^{*}$ such that one of the equilibrium points ${\stackrel{\u0304}{Q}}_{>}^{*}$ and ${\stackrel{\u0304}{Q}}^{0}$ is an attractor depending on the relative position of the initial conditions supplied to the dynamical system (1) with respect to the hyper-surface h [30].
Appendix B: The local stability of the non-trivial equilibrium ${\stackrel{\u0304}{Q}}^{c}$
Let us show that the matrix A, given by equation (14), is an M-matrix [31] for big non-trivial equilibrium ${\stackrel{\u0304}{Q}}_{>}^{c}$, where T_{>} is one of the coordinates.
where I is the identity matrix and ρ is the spectral radius [32].
Or, equivalently:
Proposition 1. A is a non-singular M-matrix if and only if the real part of its eigenvalues is greater than zero.
for i = 1, 2, ..., n.
According to the first part of Proposition 2, the matrix A has positive diagonal elements, see equation (14).
because Ā can be greater than k_{4}. The equilibrium values correspond to the point ${\stackrel{\u0304}{Q}}^{c}$, given by equation (9).
which is valid in the equilibrium. Since φ_{ d }> 0, if we show that φ_{ n }> 0, then there exists a positive number φ.
and we have ${\Phi}_{2}\left(\stackrel{\u0304}{T}\right)\le 0$, for ${\Phi}_{2<}\le \stackrel{\u0304}{T}\ge {\Phi}_{2>}$; otherwise, ${\Phi}_{2}\left(\stackrel{\u0304}{T}\right)>0$. Notice that Φ_{2>} is out of the range of feasibility, once Φ_{2>} >T_{ E }. We know that $\left({\alpha}_{1}-{\mu}_{1}\right)/{\beta}_{1}<\stackrel{\u0304}{T}<{k}_{3}$, but the bigger solution under consideration is ${T}_{>}=\left({\alpha}_{1}-{\mu}_{1}\right)/{\beta}_{1}^{c}={T}_{E}^{c}$, for ${\beta}_{1}\ge {\beta}_{1}^{c}$. It is easy to show that Φ_{2}(k_{3}) < 0; but, ${\Phi}_{2}\left({T}_{E}^{c}\right)\le 0$, if ${T}_{E}^{c}\ge {\Phi}_{2<}$, and ${\Phi}_{2}\left({T}_{E}^{c}\right)>0$, if ${T}_{E}^{c}<{\Phi}_{2<}$. Hence, ${\Phi}_{2}\left(\stackrel{\u0304}{T}\right)\le 0$ if ${T}_{E}^{c}<{\Phi}_{2<}$. Hence, $\Phi \left(\stackrel{\u0304}{T}\right)>0$, for $\stackrel{\u0304}{T}\ge 0$, when ${T}_{E}^{c}\ge {\Phi}_{2<}$, or when ${T}_{E}^{c}\ge {\Phi}_{2<}$ and $\stackrel{\u0304}{T}\ge {\Phi}_{2<}$; and $\Phi \left(\stackrel{\u0304}{T}\right)>0$, when ${T}_{E}^{c}<{\Phi}_{2<}$ and ${T}_{E}^{c}<\stackrel{\u0304}{T}<{\Phi}_{2<}$ for ε > ε_{min}.
Summarizing, φ_{ n }> 0 occurs, in order to have positive number φ:
1. T_{ A }< T_{ E }or γ < (α_{2} - μ_{2})/k_{3} - weak capacity of recruitment of the normal epithelial cells by cancer cells. We have:1.a If ${T}_{E}^{c}\ge {\Phi}_{2<}:{\phi}_{n}>0$ without restriction about T_{>}.
1.b If ${T}_{E}^{c}<\stackrel{\u0304}{T}<{\Phi}_{2<}$: we have two possibilities
1.b.1 If T_{>} > Φ_{2<}: φ_{ n }> 0 without restriction about T_{>}.
- 2.
T_{ A }> T_{ E }or γ > (α_{2} - μ_{2})/k_{3} - strong capacity of coopting normal epithelial cells by cancer cells: φ_{ n }> 0, if ε > ε_{min}, where T_{>} satisfies ε_{min}.
When γ is small, the big equilibrium ${\stackrel{\u0304}{Q}}_{>}^{c}$, with one coordinate T_{>}, is stable without conditions (cases 1.a and 1.b.1). However, for sufficiently higher values of γ, the big equilibrium point ${\stackrel{\u0304}{Q}}_{>}^{c}$ can be unstable (case 2), in which case Hopf bifurcation can occur (see Appendix C).
The bigger roots of ${\Phi}_{1}\left(\stackrel{\u0304}{T}\right)$ and ${\Phi}_{2}\left(\stackrel{\u0304}{T}\right)$ are not biologically feasible, because they are, respectively, higher than k_{3} and (α_{2} - μ_{2})/γ. Let us define ${T}_{m}=\text{min}\left\{\underset{<}{\overset{1}{T}},{T}_{<}^{2}\right\}$, where $\text{min}\left\{\underset{<}{\overset{1}{T}},{T}_{<}^{2}\right\}$ is the minimum between the small roots of, respectively, ${\Phi}_{1}\left(\stackrel{\u0304}{T}\right)$ and ${\Phi}_{2}\left(\stackrel{\u0304}{T}\right)$ given by ${T}_{<}^{1}={k}_{3}\left(1-\sqrt{1-\frac{{\mu}_{3}}{{k}_{4}{\alpha}_{3}}}\right)$, assuming that k_{4} ≥ 1, and ${T}_{<}^{2}={T}_{E}\left(1-\sqrt{1-\frac{{T}_{A}}{{T}_{E}}}\right)$, assuming that T_{ E }≥ T_{ A }. Hence, if $\stackrel{\u0304}{T}<{T}_{m}$, we have $\Phi \left(\stackrel{\u0304}{T}\right)>0$ and there is not a positive number φ. Notice that, when T_{ E }< T_{ A }, we have ${\Phi}_{2}\left(\stackrel{\u0304}{T}\right)>0$, and it is enough to satisfy $\stackrel{\u0304}{T}<{T}_{<}^{1}$.
Appendix C: Hopf bifurcation
Dynamical trajectories are shown in Figure 19 for γ near ${\gamma}_{1}^{c}:{\gamma}_{1}=2.31\times 1{0}^{-2}\lesssim {\gamma}_{c}^{1}$ and ${\gamma}_{2}=2.32\times 1{0}^{-2}\gtrsim 1{0}^{-2}\gtrsim {\gamma}_{1}^{c}$. When T(0) = 0.44, trajectories of both cases go to ${\stackrel{\u0304}{Q}}^{0}$ (trajectories similar to Figure 1(a), not shown). When T(0) = 0.45, dynamical trajectories for γ_{1} go to stable ${\stackrel{\u0304}{Q}}_{>}^{*}$ (a), while for γ_{2}, they oscillate around unstable ${\stackrel{\u0304}{Q}}_{>}^{*}$ (b), in which case limit cycle arises.