Abstract
This paper deals with stability and Hopf bifurcation analyses of a mathematical model of HIV infection of T-cells. The model is based on a system of delay differential equations with logistic growth term and antiretroviral treatment with a discrete time delay, which plays a main role in changing the stability of each steady state. By fixing the time delay as a bifurcation parameter, we get a limit cycle bifurcation about the infected steady state. We study the effect of the time delay on the stability of the endemically infected equilibrium. We derive explicit formulae to determine the stability and direction of the limit cycles by using center manifold theory and normal form method. Numerical simulations are presented to illustrate the results.
1. Introduction
Since 1980, the human immunodeficiency virus (HIV) or the associated syndrome of opportunistic infections that causes acquired immunodeficiency syndrome (AIDS) has been considered as one of the most serious global public health menaces. When HIV enters the body, its main target is the CD4 lymphocytes, also called CD4 T-cells (including T-cells). When a CD4 cell is infected with HIV, the virus goes through multiple steps to reproduce itself and create many more virus particles. The AIDS term, which is known as the late stage of HIV, covers the range of infections and illnesses which can result from a weakened immune system caused by HIV. Based on the clinical studies, it is known that, for a normal person, the T-cells count is around and for HIV infected patient it gradually decreases to or below, which leads to AIDS. However, this may take several years for the number of CD4 T-cells to reduce to a level where the immune system is weakened [1–6].
Mathematical models, usingdelay differential equations (DDEs), have provided insights in understanding the dynamics of HIV infection. Discrete or continuous time delays have been introduced to the models to describe the time between infection of a T-cell and the emission of viral particles on a cellular level [7–13]. In general, DDEs exhibit much more complicated dynamics than ODEs since the time delay could cause a stable equilibrium to become unstable and cause the populations to fluctuate [14–16]. In studying the viral clearance rates, Perelson et al. [17] assumed that there are two types of delays that occur between the administration of drug and the observed decline in viral load: a pharmacological delay that occurs between the ingestion of drug and its appearance within cells and an intracellular delay that is between initial infection of a cell by HIV and the release of new virion. In this paper, we incorporate an intracellular delay to the model to describe the time between infection of a T-cell and the emission of viral particles on a cellular level [18]. We study the impact of the presence of such time delay on the dynamics of the model.
The outline of the present paper is as follows. In Section 2, we describe the model. In Section 3, we study the qualitative behavior of the model via stability of the steady states and Hopf bifurcation when time delay is considered as a bifurcation parameter. In Section 4, we provide an explicit formula to determine the direction of bifurcating periodic solution by applying center manifold theory and normal form method. We provide some numerical simulations to demonstrate the effectiveness of the analysis in Section 5 and we conclude in Section 6.
2. Description of the Model
Let us start the analysis with some basic models of the dynamics of target (uninfected) cells and infected T-cells by HIV. As a first approximation, the dynamics between HIV and the macrophage population was described by the simplest model of infection dynamics presented in [19–21]. Denoting uninfected cells by and infected cells by and assuming that viruses are transmitted mainly by cell to cell contact, the model is given by The target (uninfected) T-cells are produced at a rate , die at a rate , and become infected by virus at a rate . The infected host cells die at a rate . The basic reproductive ratio of the virus is then given by . If there is no infection or if , there is only trivial equilibrium () with no virus-producing cells. Whereas if , the virus can establish an infection and the system converges to the equilibrium with both uninfected cells and infected cells, .
However, in most viral infections, the CTL response plays a crucial part in antiviral defence by attacking viral infected cells [22, 23]. As the the cytotoxic T-lymphocyte (CTL) immune response is necessary to eliminate or control the viral infection, we incorporated the antiviral CTL immune response into the basic model (1). Therefore, if we add CTL response, which is denoted by , into model (1) (see [19]), then the extended model is Thus, CTLs proliferate in response to antigen at a rate , die at a rate , and lyse infected cells at a rate . We assume that the CTL pool consists of two populations: the precursors and the effectors . In other words, we assume that there are primary and secondary responses to viral infections. Then, the model (2) becomes The infected cells are killed by CTL effector cells at a rate . Upon contact with antigen, CTLp proliferate at a rate and differentiate into effector cells CTLe at a rate . CTL precursors die at a rate , and effectors die at a rate ; see Figure 1.
Since the proliferation of T-cells is density dependent, that is, the rate of proliferation decreases as T-cells increase and reach the carrying capacity, we then extend the above basic viral infection model to include the density dependent growth of the T-cell population (see [24–26]). It is also known that HIV infection leads to low levels of T-cells via three main mechanisms: direct viral killing of infected cells, increased rates of apoptosis in infected cells, and killing of infected T-cells by cytotoxic T-lymphocytes [26]. Hence, it is reasonable to include apoptosis of infected cells. An average of viral particles is produced by infected cells per day. The treatment with single antiviral drug is considered to be failed, so that the combination of antiviral drugs is needed for the better treatment [25]. Therefore, in the below revised model, we combine the antiretroviral drugs, namely, reverse transcriptase inhibitor (RTI) and protease inhibitor (PI) to make the model realistic (see [27–29]). RTIs can block the infection of target T-cells by infectious virus, and PIs cause infected cells to produce noninfectious virus particles. The modified model takes the form The first equation of model (4) represents the rate of change in the count of healthy T-cells that produced at rate and become infected at rate , with the mortality . We assume that the uninfected T-cells proliferate logistically, thus the growth rate is multiplied by the term and this term approaches zero when the total number of T-cells approaches the carrying capacity . The effects of combination of RTI and PI antiviral drugs are represented by the term , where , , represents the effects of RTI and , , represents the effects of PI. The second equation of model (4) denotes the rate of change in the count of infected T-cells. The infected T-cells decay at a rate and denotes apoptosis rate of infected cell; infected cells are killed by CTL effectors at a rate . The third equation of the model denotes the rate of change in the CTLp population; proliferation rate of the CTLp is given by and is proportional to the infected cells ; CTLp die at a rate and differentiate into CTL effectors at a rate . The last equation of the model represents the concentration of CTL effectors, which die at a rate . In reality, the specific immune system is not immediately effective following invasion by a novel pathogen. There may be an explicit time delay between infection and immune initiation and there may be a gradual build-up in immune efficacy during which the immune response develops, before reaching maximal specificity to the pathogen ([8, 30, 31]). In order to make model (4) more realistic, time delay in the immune response should be included in the following model: The range of parameter values of the model are given in Table 1.
We start our analysis by presenting some notations that will be used in the sequel. Let be the Banach space of continuous functions mapping the interval into , where ; the initial conditions are given by where are smooth functions, for all . From the fundamental theory of functional differential equations (see [32, 33]), it is easy to see that the solutions of system (5) with the initial conditions as stated above exist for all and are unique. It can be shown that these solutions exist for all and stay nonnegative. In fact, if , then for all . The same argument is true for the , , and components. Hence, the interior is invariant for system (5).
3. Steady States
We can obtain the steady state values by setting . The steady state value of the infection-free steady sate is given by , while the infected steady state is given by and is given by the following quadratic equation: where , , .
3.1. Stability and Hopf Bifurcation Analysis of Infected Steady State
In order to study full dynamics of model (4) by using time delay as a bifurcation parameter, we need to linearize the model around the steady state and determine the characteristic equation of the Jacobian matrix. The roots of the characteristic equation determine the asymptotic stability and existence of Hopf bifurcation for the model. The characteristic equation of the linearized system is given by which is equivalent to the equation where and Let us consider the following equation: For the nondelayed model (say ), from (10), we have where
Lemma 1. For , the unique nontrivial equilibrium is locally asymptotically stable if the real parts of all the roots of (13) are negative.
Proof. The proof of the above lemma is based on holding the following conditions: , , , and , as proposed by Routh-Hurwitz criterion. We conclude that equilibrium is locally asymptotically stable if and only if all the roots of the characteristic equation (13) have negative real parts which depends on the numerical values of parameters that are shown in the numerical exploration.
3.2. Existence of Hopf Bifurcation
We here study the impact of the time-delay parameter on the stability of HIV infection of T-cells. We deduce criteria that ensure the asymptotic stability of infected steady state , for all . We arrive at the following theorem.
Theorem 2. Necessary and sufficient conditions for the infected equilibrium to be asymptotically stable for all delay are as follows (i)the real parts of all the roots of are negative;(ii)for all and , , where .
Proof. Assume that Lemma 1 is true. Now, for , we have Substituting into (5) and separating the real and imaginary parts of the equations yields After some mathematical manipulations, we obtain the following equations Let From (16), we have where The conditions and of Theorem 2 hold if and only if (19) has no real positive root.
Let ; then (19) takes the form If , then (19) has at least one positive root. In the case when (19) has four positive roots, we have From (16), we have where and ; we choose .
To establish Hopf bifurcation at , we need to show that By differentiating (10) with respect to , we can get It follows that Then, by combining (10), we get Substituting in (27) (where and ) yields where Thus, Notice that By summarizing the above analysis, we arrive at the following theorem.
Theorem 3. The infected equilibrium of the system (5) is asymptotically stable for and it undergoes Hopf bifurcation at .
4. Direction and Stability of Bifurcating Periodic Solutions
In the previous section, we obtained conditions for Hopf bifurcation to occur when , . It is also important to derive explicit formulae from which we can determine the direction, stability, and period of periodic solutions bifurcating around the infected equilibrium at the critical value . We use the cafeteria of normal forms and center manifold proposed by Hassard [34]. We assume that the model (5) undergoes Hopf bifurcation at the infected equilibrium when , , and then are the corresponding purely imaginary roots of the characteristic equation at the infected equilibrium . Assume also that then using Taylors expansion for system (3) at the equilibrium point yields Here, For convenience, let and for . Denote ; has -order continuous derivative. For initial conditions , (33) can be rewritten as where , , and are given, respectively, by is one parameter family of bounded linear operators in and From the discussion in the above section, we know that if , then model (5) undergoes a Hopf bifurcation at the infected equilibrium , and the associated characteristic equation of model (5) has a pair of purely imaginary roots . By Reisz representation, there exists a function of bounded variation for such that In fact, we can choose where is Dirac delta function. Next, for , define Since , (35) can be written as where , . For , the adjoint operator of can be defined as For and , in order to normalize the eigenvalues of operator and adjoint operator , the following bilinear form is defined by where and is complex conjugate of . It can verify that and are adjoint operators with respect to this bilinear form.
We assume that are eigenvalues of and the other eigenvalues have strictly negative real parts. Thus, they are also eigenvalues of . Now we compute the eigenvector of corresponding to the eigenvalue and the eigenvector of corresponding to the eigenvalue . Suppose that is eigenvector of associated with ; then, . It follows from the definition of and (36), (38), and (40) thatSolving (45), we can easily obtain , where Similarly, suppose that the eigenvector of corresponding to is , . By the definition of and (36), (38), and (40), one gets Solving (47), we easily obtain , where In order to assure that , we need to determine the value of . From (44), one gets Let On the center manifold , we have where and are local coordinates of the center manifold in the direction of and , respectively. Note that is real if is real. So we only consider real solutions. From (50), we obtain For the solution of (35), from (41) and (44), since , we have Rewrite (54) as where Substituting (42) and (54) into (50) yields which can be written as where On the center manifold , we have Substituting (52) and (55) into (61), one obtains Substituting (52) and (60) into (59) yields Comparing the coefficients of (62) and (63), one gets Since , then we have Thus, we obtain It is obvious that So also and hence It follows from (54) that where