Stability and Hopf bifurcation of an HIV infection model with saturation incidence and two delays

In this paper, the dynamical behaviors of a viral infection model with cytotoxic T-lymphocyte (CTL) immune response, immune response delay and production delay are investigated. The threshold values for virus infection and immune response are established. By means of Lyapunov functionals methods and LaSalle's invariance principle, sufficient conditions for the global stability of the infection-free and CTL-absent equilibria are established. Global stability of the CTL-present infection equilibrium is also studied when there is no immune delay in the model. Furthermore, to deal with the local stability of the CTL-present infection equilibrium in a general case with two delays being positive, we extend an existing geometric method to treat the associated characteristic equation. When the two delays are positive, we show some conditions for Hopf bifurcation at the CTL-present infection equilibrium by using the immune delay as a bifurcation parameter. Numerical simulations are performed in order to illustrate the dynamical behaviors of the model.


1.
Introduction. As it is well known, infectious diseases have been a serious threat to human health, and viral dynamics play an important role in studying these diseases. Several mathematical models have been proposed to describe the viral infection process with immune response such as influenza, HBV, HIV, HTLV, SARS and Ebola (see [5,27,32,17,30,16,21,12,3,29,23,10]). These within-host viral infection models can capture some essential features of the immune system and are able to produce a variety of immune responses which are observed experimentally and clinically. Findings from viral infection models are valuable to help us improve understanding of both diseases and various drug therapy strategies against them.
Because of the importance of the biological meanings, dynamical properties of HIV infection models have been studied by many authors. There has been much interest in mathematical modeling of within-host viral dynamics, some include the immune response directly [2,25,24,31,14,19,9,26], and others do not contain the immune response [22,15,8,28]. During the process of viral infection, immune 2366 HUI MIAO, ZHIDONG TENG AND CHENGJUN KANG response in the host is induced which is initially rapid and nonspecific (natural killer cells, macrophage cells, etc.) and then delayed and specific (cytotoxic T lymphocyte cells, antibody cell). But in most viral infections, cytotoxic T lymphocyte (CTL) immune response cells which attack infected cells play a critical part in antiviral defense. In order to investigate the dynamics of viral infection model with CTL response, Nowak and Bangham [14] constructed a mathematical model describing the basic dynamics of the interactions between the uninfected target cells x(t), productively infected cells y(t), free virus v(t) and CTL response cells z(t) as follows where the parameter s represents the rate at which new target cells are created, d is the death rate of uninfected target cells, β is the infection rate of uninfected cells by virus, a is the death rate of productively infected cells, p represents the killing rate of infected cells by CTL response cells, k is the rate of the virus particles produced by infected cells, u is the virus clearance rate constant, c is the rate at which the CTL response cells are produced, b is the death rate of the CTL response.
It is assumed in model (1) that the infection process is governed by the massaction principle, i.e. that the infection rate per host and per virus is a constant. However, experiments reported in [6] strongly suggested that the infection rate of microparasitic infections is an increasing function of the parasite dose, and is usually sigmoidal in shape (see [18]). In [20], to place the model on more sound biological grounds, Regoes et al. replaced the mass-action infection rate with a dosedependent infection rates. Song et al. in [20], a more general saturated infection rate, βx(t)v(t) θ 1+αv(t) ρ , was suggested, where θ, ρ and α are positive constants. On the other hand, delayed biological systems of one type or another have been considered by a number of authors (see [13]). Delay can arise for various practical reasons in epidemiology. For example, Hethcote et al. in [7] considered a model in which the delay is introduced in the removed class to account for the period of temporary immunity. Culshaw and Ruan [4] considered the time delay between infection of a CD4+ T-cell and the emission of viral particles on a cellular level to investigate the effect of the time delay on the stability of the endemically infected equilibrium. Delay is also used to model the gestation lag, the incubation time for a infectious vector and the time delay in loss of vaccine, etc. These delay-differential equation systems often exhibit much more complicated dynamical behavior than those ordinary differential equation systems since a time delay could cause a stable equilibrium to become unstable and induce oscillations. We note that in model (1) this process was assumed to occur instantaneously: as soon as virus contacts a target cell, the cell begins producing virus. This is not biologically sensible. In reality, there is a time delay between initial virus entry into a cell and subsequent virus production. τ 1 represents the maturation time of the newly produced viruses. Constant m is assumed to be the death rate for new virus during time period [t − τ 1 , t]. Therefore, the possibility of surviving from t − τ 1 to t is e −mτ1 . The immune response plays an important role in eliminating or controlling the disease after human body is infected by virus. Antigenic stimulation generating CTLs may need a period of time τ 2 , i.e. the CTL response at time t may depend on the population of antigen at time t − τ 2 .
Based on the above discussion in this paper, we consider a delayed viral infection model with saturation incidence rate and CTL immune response given by We will study the global dynamics of model (2). We will show that model (2) has three possible equilibria: the infection-free equilibrium, CTL-absent infection equilibrium and the CTL-present infection equilibrium. By using Lyapunov-LaSalle invariance principle, we will prove the global asymptotic stability of these equilibria. By using a geometric method developed by [1], we will discuss the existence of Hopf bifurcation at the CTL-present infection equilibrium by choosing the delay τ 2 as a bifurcation parameter. By numerical simulations, we verify the theoretical results obtained in this paper. The paper is organized as follows. In the next section we deal with some basic properties such as positivity and boundedness of solutions and the existence of equilibria of model (2). In Section 3, we state and prove the main theorems on the local and global stability of equilibria. In Section 4, providing two delays are positive, we establish the results for the Hopf bifurcation of the CTL-present infection equilibrium. In Section 5, some numerical examples are given to illustrate our main results. In the final section, we present a brief conclusion.
be the solution of model (2) satisfying initial condition (3), then x(t), y(t), v(t) and z(t) are positive and ultimately bounded.
For model (2), the basic reproductive ratio of virus which describes the average number of newly infected cells generated from one infected cell at the beginning of the infectious process is given by and the immune response reproductive ratio is given by .
By direct calculation we know that model (2) has the following three possible equilibria. Infection-free equilibrium E 0 = (x 0 , 0, 0, 0) = ( s d , 0, 0, 0). If R 0 > 1, then there is a unique CTL-absent infection equilibrium If R 1 > 1, then there is a unique CTL-present infection equilibrium 3. Stability analysis. In this section, we discuss the local and global stability of three equilibria E 0 , E 1 and E 2 of model (2).
Theorem 3.1. (a) If R 0 < 1, then the infection-free equilibrium E 0 is locally asymptotically stable.
(b) If R 0 > 1, then the equilibrium E 0 is unstable.
Proof. Consider conclusion (a). At the equilibrium E 0 , the characteristic equation of the corresponding linearized system of model (2) is Two roots of (4) are λ 1 = −b and λ 2 = −d. The remaining roots of (4) satisfy the equation If a root λ of (5) has nonnegative real part, then the modulus of the left-hand side of (5) satisfies while the modulus of the right-hand side (5) satisfies This leads to a contradiction. Hence, when R 0 < 1, all roots of (4) have negative real parts. This implies that the equilibrium E 0 is locally asymptotically stable. Next, we consider conclusion (b).
Thus, f (λ) = 0 has at least one positive real root. Therefore, when R 0 > 1, the equilibrium E 0 is unstable.
Theorem 3.2. If R 0 < 1, then the infection-free equilibrium E 0 is globally asymptotically stable.
Proof. Define a Lyapunov functional V 1 (t) as follows Calculating the time derivative of V 1 (t) along any positive solution of model (2) and noticing that x 0 = s d , we can obtain From the third equation of model (2), we easily obtain y = 0. Hence, M = {E 0 }. It follows from LaSalle's invariance principle [11] that the equilibrium E 0 of model (2) is globally asymptotically stable.

Stability of equilibrium
Proof. Consider conclusion (a). At the equilibrium E 1 , the characteristic equation of the corresponding linearized system of model (2) is where

Consider the equation
When τ 1 = 0, (7) becomes Since by the Routh-Hurwitz criterion, all roots of (8) have negative real parts. Let λ = iω with ω > 0 be a purely imaginary root of (7). Separating real and imaginary parts, it follows that Squaring and adding the two equations of (9) yields Let r = ω 2 , then (10) becomes By calculating, we have By the Routh-Hurwitz criterion, equation (11) has no positive roots. This shows that (7) cannot have any purely imaginary root λ = iω. Next, we analyze the transcendental equation When τ 2 = 0, we have Clearly, if R 1 < 1, then λ < 0. Let λ = iω with ω > 0 be a purely imaginary roots of (12), then we have Therefore, we conclude that for any τ 1 ≥ 0 and τ 2 ≥ 0 characteristic equation (6) has no any root with nonnegative real part. By the general theory of the delay differential equations from [11], we have that when R 1 < 1 < R 0 the equilibrium E 1 is locally asymptotically stable.
3.3. Stability of equilibrium E 2 . On the stability analysis of the CTL-present infection equilibrium E 2 , we only consider the special cases: τ 2 = 0, τ 1 ≥ 0.
Proof. At the equilibrium E 2 , the characteristic equation of the corresponding linearized system of model (2) is When τ 2 = 0, (13) becomes (14) has a root λ = m + iω 0 with m ≥ 0 and ω 0 ≥ 0, then we have we obtain that the modulus of the left-hand side of (15) is greater than the modulus of the right-hand side. This leads to a contradiction. Thus, we conclude that for any τ 1 ≥ 0 equation (14) does not have any root with nonnegative real part. Therefore, the equilibrium E 2 is locally asymptotically stable.
Proof. Define a Lyapunov functional V 3 (t) as follows where x * , y * , v * and z * satisfy the following equations Calculating the time derivative of V 3 (t) along the solution of model (2), and using the similar method with the proof of Theorem 3.4, we have From the above calculation, it follows that dV3(t) dt ≤ 0 for all x, y, v, z > 0 and dV3(t) dt = 0 if and only if (x(t), y(t), v(t), z(t)) = (x * , y * , v * , z * ). Therefore, from Lyapunov LaSalle's invariance principle [11], the equilibrium E 2 of model (2) is globally asymptotically stable.
4. Bifurcation analysis. The following discussions focus on the Hopf bifurcation of the equilibrium E 2 when delays τ 1 ≥ 0 and τ 2 > 0. Equation (13) can be rewritten where From Theorem 3.5 we have obtained that when R 1 > 1 and τ 2 = 0, for any τ 1 ≥ 0, all roots of (16) have negative real parts. Now, we fix τ 1 = τ 0 1 ≥ 0 and increase the value of τ 2 from zero to see the possible of appearing Hopf bifurcation. In this case, the characteristic equation (16) becomes where Notice that all roots of (17) have negative real parts when τ 2 = 0, and then stability switch may occur when an imaginary root λ = iω exists in (17) and crosses the imaginary axis as the value of τ 2 increases.
Next, we determine the direction of geometric crossing. Differentiating (17) with respect to τ 2 leads to where the variable (λ, τ 0 1 , τ 2 ) of functions P, Q, T, G and their derivatives is omitted for brevity. The equation (17)  The equality in (23) implies (27) By the calculation, we have and Due to (26), we further obtain that (28) and (29) become into where From (24), (27) and (30), we have On the other hand, S n (τ * 2 ) = 0 implies , and then By (22) and (23) it follows that for τ 2 ∈ I, where θ(τ 2 ) is defined by (26). Then we have where From (34) and noting that From (31), it becomes From (32), (33) and (37), we finally have Applying Theorem 3.5 and the Hopf bifurcation theorem for functional differential equation in [1], we finally can conclude the existence of a Hopf bifurcation as stated in the following theorem.
6. Conclusions. In this paper, we included two delays in a viral infection model to study HIV dynamics. One represents the time needed for the maturation time of the newly produced viruses (production delay), and the other denotes the time needed for the CTL immune response to generate to control viral replication (immune delay). We studied how virus replication delay and immune response delay impact the dynamics. During virus infections, the CTL immune response to HIV infection is observed in the first few weeks of infection, coincident with the initial decline in the plasma viral load [9]. So, we assumed that the production of CTL response cells depends on the infected cells and CTL cells based above important biological meaning. We see that similar assumption also is given in [13,24,31,14,3,29,9]. In most virus infections, CTL cells which attack infected cells play a critical part    (2) when Hopf bifurcation occurs again for τ 2 = 6.8000.
in antiviral defense. In order to obtain a comprehensive view for the CTL immune dynamics in vivo, we investigated the global stability of model (2) by utilizing the method of constructing suitable Lyapunov functionals which are motivated by recent works of Pawelek et al. [17], Huang et al. [8,9], Zhu et al. [31] and Li and Shu [19,12]. By rigorous analysis, we have shown that when R 0 < 1, the infection-free equilibrium is globally asymptotically stable, which means that the viruses are cleared and immune is not active; when R 1 < 1 < R 0 , the CTL-absent infection equilibrium exists and is globally asymptotically stable, which means that the CTL immune response would not be activated and virus infection becomes vanished; when τ 2 = 0 and R 1 > 1, the CTL-present infection equilibrium is globally asymptotically stable when it exists. Actually, these results suggest that introducing the production delay does not change the stability of the model. Furthermore, for the model including production delay and immune delay, Hopf bifurcation analysis at the CTL-present infection equilibrium is an interesting and challenge problem. Compared to earlier studies such as Shu et al. [19] analyzed the global stability for a viral model with nonlinear incidence rate, infinitely distributed delays and CTL response, but the global stability for CTL-present infection equilibrium of model is established under the assumption that there is no immune delay in the CTL response. Our results in the present paper show that the existence of Hopf bifurcation at CTL-present infection equilibrium is investigated by a geometric method. In Li and Shu [12], it is shown that the existence of Hopf bifurcation in a HTLV model only with immune response delay. Our mathematical analysis reveals the sustained oscillations occur when the production delay and immune delay are incorporated simultaneously in the model. These results have enriched our understanding of the effects of production delay in a viral infection processes and its interaction with the immune response. Meanwhile, the analysis helps us to understand how virus replication delay and immune response delay affect the cell-virus-immune system dynamics.
The existence of Hopf bifurcation in our model with two positive time delays leads to emergence of viral oscillations. Furthermore, when τ 1 = τ 0 1 is fixed and τ 2 is very small, the CTL-present infection equilibrium is still stable. Numerical simulations are shown in Fig.1. With the increase of the delay τ 2 , the CTL-present infection equilibrium becomes unstable, and the Hopf bifurcation occurs. The corresponding numerical simulations are shown in Fig.2. Moreover, if we increase τ 2 , then the CTLpresent infection equilibrium becomes stable again, as shown in Fig.3. Furthermore, and along with the increase of τ 2 , the CTL-present infection equilibrium becomes unstable and periodic oscillations appear, as shown in Fig.4. It is shown that an increase of the virus production delay may stabilize the the infection-free equilibrium and CTL-absent infection equilibrium, but the immune delay can destabilize the CTL-present infection equilibrium, leading to stability switches. This may be a new approach for clinicians to design a new drug which may increase the virus production delay and bring the infection under control eventually.
As it is well known, in our body the immune response is made up of both CTL immune response and humoral response. The CTL immune response is that T cells kill the infected cells, the humoral response is that B cells produce an antibody to neutralize the virus. In this paper, we only consider the CTL immune response. In the future, our work will focus on the combined effects of the two kinds of immune response.