Abstract and Applied Analysis

Volume 2015, Article ID 354918, 11 pages

http://dx.doi.org/10.1155/2015/354918

## Hopf Bifurcation, Cascade of Period-Doubling, Chaos, and the Possibility of Cure in a 3D Cancer Model

Departamento de Matemática e Computação, Faculdade de Ciências e Tecnologia, Universidade Estadual Paulista (UNESP), 19060-900 Presidente Prudente, SP, Brazil

Received 27 June 2014; Accepted 13 October 2014

Academic Editor: Yongli Song

Copyright © 2015 Marluci Cristina Galindo et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

We study a cancer model given by a three-dimensional system of ordinary differential equations, depending on eight parameters, which describe the interaction among healthy cells, tumour cells, and effector cells of immune system. The model was previously studied in the literature and was shown to have a chaotic attractor. In this paper we study how such a chaotic attractor is formed. More precisely, by varying one of the parameters, we prove that a supercritical Hopf bifurcation occurs, leading to the creation of a stable limit cycle. Then studying the continuation of this limit cycle we numerically found a cascade of period-doubling bifurcations which leads to the formation of the mentioned chaotic attractor. Moreover, analyzing the model dynamics from a biological point of view, we notice the possibility of both the tumour cells and the immune system cells to vanish and only the healthy cells survive, suggesting the possibility of cure, since the interactions with the immune system can eliminate tumour cells.

#### 1. Introduction

Mathematical models describing cancer tumour growth have been extensively studied in the literature in order to understand the mechanism of disease and to predict its future behaviour; see [1] and references therein as a brief review on the subject. One class of these cancer models, well known as Lotka-Volterra systems, is based on nonlinear mathematical systems of interacting species.

The competition for nutrients, the effects of growth factors, and the actions of the immune system stand out among possible interactions.

The immune system has the function of recognizing internal and external threats to the organism, reacting to eliminate, neutralize, or tolerate such threats. The recognition of tumour cells by the immune system can happen in distinct and complementary ways. According to Chammas et al. [2], the lymphocytes of type B participate in the immune response against tumours by producing specific antibodies against tumour antigens. These antibodies bind with specific antigens on tumour cells which facilitate their recognition and destruction by NK cells and phagocytosis by macrophages. The study of different antibodies found in patients with cancer tumour has been useful for determining tumour-associated antigens and pointed out some information for the development of therapeutic antibodies with clinical application. More specifically, several researches seeking the elimination of tumour cells through the immune cells have been developed [3–5]. Despite the increasing number of researches in this area, in [1] the authors pointed out the existence of a gap in the understanding of tumour growth and invasion, requiring an intense and closer collaboration in investigations by different scientists, such as mathematicians, biologists, and experimentalists.

Aiming to contribute to the understanding of these types of model we perform a bifurcation analysis of a cancer tumour growth model of three competing cell populations (tumour, healthy, and immune cells), which was proposed in [6]. The populations of these cells in a given instant of time will be denoted here by , , and , respectively. Analogous to what is done in competition models in population ecology, the mathematical model for the interaction of these cells is given by the following three-dimensional adimensional system of ordinary differential equations (for details on the derivation of the system, see the nice paper [6]),in which the parameters , , , , , , , are all positive.

The equations of system (1) describe the variation rate of tumour cells, healthy cells, and effector cells with the time . In the first equation, a logistic growth of tumour cells is considered, the term represents the negative effects that the healthy cells have on tumour cells, and represents the negative effects that the effector cells of the immune system exert on tumour cells. In the second equation, the healthy cells also grow logistically, with growth rate . The term represents the negative effects that tumour cells exert on healthy cells. In the last equation, describes the stimulation of the immune system by tumour cells through specific antigens. Therefore, the stimulus of the immune system depends directly on the number of tumour cells. The term represents the negative effects that tumour cells exert on the effector cells of the immune system and, finally, refers to the natural death rate of the effector cells. It is important to notice that system (1) is a slight modification of a model of cancer without treatment, proposed in [7], where a first phase space analysis of the model was presented.

In [6], by using the Shilnikov’s theorem [8] and by the calculation of Lyapunov exponents and Lyapunov dimension, the authors have shown that system (1) presents a chaotic attractor for the following fixed set of parameters: , , , , , , , and . In a second paper concerning the dynamical analysis of system (1) [9], the authors described in detail the topological structure of this chaotic attractor. However, in none of these papers the authors show how this attractor is formed nor do they relate the chaotic behaviour of system (1) with biological aspects of the cancer evolution.

In this paper we analyse the local stability of the equilibria of system (1); then, varying two of the parameters involved, namely, the parameters and , and taking for the other parameters the same values considered in [6], we perform a bifurcation analysis of system (1). Through this analysis we prove that, for fixed and varying the parameter , the system presents a supercritical Hopf bifurcation at an equilibrium point for the critical value , which leads to the creation of a small stable limit cycle. Furthermore, numerically studying the continuation of this stable limit cycle by increasing the parameter from the Hopf point, we numerically found the occurrence of a cascade of period-doubling bifurcations which leads to the formation of the chaotic attractor shown to exist in [6]. Trying to analyse this dynamical behaviour from the biological point of view, we observed that, after the occurrence of chaotic dynamics, both the tumour cells and the immune system cells (represented by the coordinates and ) vanish and only the healthy cells remain positive and their amount tends to the carrying capacity . The chaotic dynamics occurs when the parameter varies from the critical value to and the dynamics is not affected by the variation of parameter , which must only be taken different from zero, as will be shown in the calculations ahead. In this way, the bifurcation analysis presented here is performed focusing on the variation of parameter .

The parameter considered in the mentioned bifurcation analysis represents the interaction between the cancer cells and the healthy cells. However, this parameter can also be related to the action of the immune system cells against cancer cells. In fact, we can rewrite the first equation of system (1) asIn this way, for smaller than and fixed, the negative effect of the immune system cells onto the tumour cells is more significant and it decreases as approaches the critical value , when the chaotic dynamics occurs. Although being a very complex matter, this and another tentative analysis relating the dynamical behaviour and the biological meaning of the variables of system (1) will be further considered later on, in Sections 2, 4, and 5.

We believe that the bifurcation analysis and the biological considerations presented here complement the results and the analysis of system (1) performed in [6, 7, 9].

The paper is organized as follows. In this introductory section, we present the studied model, describing the variables and parameters involved, and state the main results obtained. In Section 2, we make a linear analysis of system (1): calculate the equilibrium points, check whether these equilibria are biologically admissible (i.e., if all their coordinates are nonnegative), analyse their local stability, and give an outline of the Hopf bifurcation theorem. In Section 3 we prove the main result of the paper, the occurrence of a supercritical Hopf bifurcation in one of the equilibrium points of system (1), which gives rise to a small stable limit cycle; we also present some numerical simulations which corroborate with the occurrence of this bifurcation. In Section 4 we numerically show that a cascade of period-doubling bifurcations occurs with the limit cycle created at the Hopf bifurcation, leading to the creation of a chaotic attractor (exactly that one shown to exist in [6]). Finally, Section 5 presents some concluding remarks, including a possible interpretation of the dynamical properties of system (1) from a biological point of view.

#### 2. Linear Analysis

Considering in system (1), we obtain the following equilibrium points (for details on the calculation and on the local analysis of the equilibrium points of system (1), see [6]):

These equilibria are called biologically admissible if their three coordinates are greater than or equal to zero. The equilibrium points , , and satisfy such condition.

For the equilibrium point we have the following possibilities:(i)if , the equilibrium coincides with the equilibrium ;(ii)if , coincides with the equilibrium ;(iii)if , , and , is biologically admissible;(iv)if , , and , is also admissible;(v)if none of the previous inequalities are satisfied, is not biologically admissible.

As the equilibria , and , may present negative coordinates, aiming to simplify the analysis of system (1), we will present a study of the feasibility of these equilibria using the values of parameters extracted from [6], that is, , , , , , and , and letting the parameters and vary. For these parameter values we obtain Note that the equilibria , , remain the same and and are not biologically admissible; therefore it will not be considered here. Let us analyse the conditions on the equilibria , , and . (i) is biologically admissible, since we are considering ;(ii)the equilibrium exists if and coincides with the equilibrium for , as seen previously. If , is biologically admissible, which does not occur if ;(iii) is biologically admissible if , since .

##### 2.1. Linear Stability of the Equilibrium Points

The Jacobian matrix of system (1) at a point is given byIt will be used bellow to study the local stability of the admissible equilibria.(i)Equilibrium : the Jacobian matrix applied at has the eigenvaluesAs and , is (unstable) saddle point, which is coherent from the biological point of view, since this point represents the simultaneous extinction of normal cells, immune cells, and tumour cells.(ii)Equilibrium : the eigenvalues of the matrix at this point areThus, the stability of depends on the value of , whereas the other two eigenvalues are negative. If , that is, if , then the equilibrium is a saddle point, therefore unstable; if , is asymptotically stable (a stable node); and if , is nonhyperbolic, since . The equilibrium point itself represents the absence of cancer and immune system cells and the permanence of the healthy cells, which tend to its carrying capacity . In this way, from the biological point of view the instability of (case ) represents the existence of cancer cells and, consequently, the existence of immune system cells that is the persistence of the disease; on the other hand, the stability of represents the absence of disease, which could be possible through a given treatment or by the action of the immune system. As the model considered does not admit treatment and includes the action of immune system cells, we expect that these cells somehow combat the tumour cells, thus enabling the cure, which justifies the possible stability of the equilibrium point (case ).(iii)Equilibrium : the Jacobian matrix at has the eigenvaluesAs , is a saddle point; therefore it is unstable, which is expected from the biological standpoint, since this equilibrium point represents the extinction of healthy and immune cells, remaining only tumour cells.(iv)Equilibrium : the Jacobian matrix calculated at has the eigenvaluesfrom which follows that is a saddle-focus. As , the equilibrium is unstable, which could be possible from the biological point of view, since this equilibrium represents the extinction of healthy cells and the persistence of tumour cells and immune system cells.(v)Equilibrium : the Jacobian matrix at this point has two eigenvalues of the formwith The third eigenvalue is given by Recalling that is biologically admissible if , we have . Hence is a saddle point. The instability of represents the biological fact that, for the set of parameters adopted, tumour cells and healthy cells cannot coexist without the presence of immune cells; that is, the tumour cells actually activate the immune system cells.(vi)Equilibrium : the Jacobian matrix calculated at the point has the following characteristic polynomialwhere , , and . Hence the linear stability of depends on the value of parameter .

We observe that the equilibrium point plays an important role in our analysis, since it represents the coexistence of the three types of cells considered in the model. Thus a detailed analysis of its stability will be presented here. The following result will be used in such analysis (for a proof, see [10]).

Lemma 1. *Consider the cubic polynomial**with real coefficients. *(a)*Then (Routh-Hurwitz Criterion) is stable (i.e., all its roots have negative real part) if and only if*(b)*If and , then has two roots with positive real part and one negative real root.*(c)*If and , then has two complex conjugate roots with zero real part and one negative real root.*

*Applying the previous lemma to the polynomial (13), we can obtain the equilibrium as a candidate for the occurrence of a Hopf bifurcation in system (1). We will prove that this bifurcation indeed occurs at the point in the next section. Before doing this, let us remember some results about the Hopf bifurcation theory.*

*2.2. Outline of the Hopf Bifurcation Theory and the Projection Method*

*This section is a brief review of the projection method described in [8] for the calculation of the first Lyapunov coefficient associated with the well-known Hopf bifurcation, denoted by .*

*Consider the differential equationwhere and are, respectively, vectors representing phase variables and control parameters. Assume that is of class in . Suppose that (16) has an equilibrium point at and, denoting the variable also by , writeas where and, for , and so on for , , , and .*

*Suppose that is an equilibrium point of (16) where the Jacobian matrix has a pair of purely imaginary eigenvalues , , and admits no other eigenvalue with zero real part. Let be the generalized eigenspace of corresponding to . By this it is meant the largest subspace invariant by on which the eigenvalues are .*

*Let be vectors such thatwhere is the transpose of the matrix . Any vector can be represented as , where . The two dimensional center manifold associated with the eigenvalues can be parameterized by the variables and by means of an immersion of the form , where has a Taylor expansion of the formwith and . Substituting this expression into (16) we obtain the following differential equation:where is given by (17). The complex vectors are obtained solving the system of linear equations defined by the coefficients of (22), taking into account the coefficients of , so that system (22), on the chart for a central manifold, is written as follows:with .*

*The first Lyapunov coefficient is defined bywhere and .*

*A Hopf point of system (16) is an equilibrium point where the Jacobian matrix has a pair of purely imaginary eigenvalues , and the other eigenvalue . From the Center Manifold Theorem, at a Hopf point a two-dimensional center manifold is well-defined, and it is invariant under the flow generated by (16) and can be continued with arbitrary high class of differentiability to nearby parameter values (see [8]).*

*A Hopf point is called transversal if the parameter dependent complex eigenvalues cross the imaginary axis with nonzero derivative. In a neighborhood of a transversal Hopf point with the dynamic behaviour of the system (16), reduced to the family of parameter-dependent continuations of the center manifold, is orbitally topologically equivalent to the following complex normal form: where and , , and are real functions having derivatives of arbitrary higher order, which are continuations of , , and the first Lyapunov coefficient at the Hopf point. See [8] for details. When () one family of stable (unstable) periodic orbits can be found on this family of manifolds, shrinking to an equilibrium point at the Hopf point.*

*Based on these results, a numerical algorithm for calculating the first Lyapunov coefficient of a system with and can be obtained. This algorithm has been used to make the calculations of the next section.*

*3. Hopf Bifurcation at the Point *

*By using the Hopf bifurcation theorem stated in [8] and the projection method presented at the same reference and synthesized in the previous section we prove the following result.*

*Theorem 2. Consider system (1) with the following parameter values: , , , , , , and . If , then the equilibrium point is a stable focus of system (1). If is a weak stable focus (or a vague attractor). On the other hand, for , and close to this value, the equilibrium becomes an unstable focus and a small stable limit cycle arises around . In short, a supercritical Hopf bifurcation occurs at the equilibrium point for the critical value .*

*Proof. *Considering the characteristic polynomial of the Jacobian matrix calculated at the point , given by (13), we have (1);(2) if ;(3) if ;(4) if .

The necessary condition for the point being biologically admissible is ; therefore is always positive. If , then and ; then, according to the Routh-Hurwitz Criterion, is a stable focus. On the other hand, if , then and . Hence, from item (b) of Lemma 1 it follows that (13) has two roots with positive real part and a negative real root. Thus is an unstable focus, remaining unstable for or .

If , then and ; hence from item (c) of Lemma 1, the polynomial (13) has two complex conjugate roots with zero real part and a negative real root, which is one of the conditions for the occurrence of the Hopf bifurcation. In order to guarantee that the Hopf bifurcation actually occurs at the equilibrium and determine the stability of the limit cycle originated from it, we have to calculate the first Lyapunov coefficient of system (1) associated with this point. Using the projection method presented in the previous section, we obtained the value for the first Lyapunov coefficient of system (1) at the point : To finish the proof, the transversality condition remains to be checked. Using software MAPLE, we calculate the derivative of the Jacobian matrix with respect to the parameter applying the critical value from which we find Solving the equation we obtain This ends the proof of Theorem 2.

*The numerical simulations presented in Figures 1 and 2 illustrate the results stated in Theorem 2.*