Global stability of an HIV infection model with saturated CTL immune response and intracellular delay

: In this paper, we consider an HIV infection model with saturated infection rate, intracellular delay and saturated cytotoxic T lymphocyte (CTL) immune response. By calculation, we obtain immunity-inactivated reproduction number R 0 and immunity-activated reproduction number R 1 . By analyzing the distribution of roots of the corresponding characteristic equations, we study the local stability of an infection-free equilibrium, an immunity-inactivated equilibrium and an immunity-activated equilibrium of the model. By constructing suitable Lyapunov functionals and using LaSalle’s invariance principle, we show that if R 0 < 1, the infection-free equilibrium is globally asymptotically stable; If R 1 < 1 < R 0 , the immunity-inactivated equilibrium is globally asymptotically stable; If R 1 > 1, the immunity-activated equilibrium is globally asymptotically stable. Sensitivity analyses are carried out to show the e ﬀ ects of parameters on the immunity-activated reproduction number R 1 and the viral load.


Introduction
Mathematical modeling and analysis of HIV pathogenesis are crucial in understanding the transmission mechanism of HIV. Many authors are interested in studying the global stability of HIV infection models (see, for example, [1][2][3][4]). In addition, the analysis of mathematical models contributes to developing new antiviral drugs and designing effective therapies for HIV infection (see, for example, [5][6][7]). Since cytotoxic T lymphocytes (CTLs) play a significant role in the within-host anti-HIV defense by attacking infected cells in HIV infection, in 1996, Nowak and Bangham [8] proposed the following HIV model describing the interaction between host cells, free viruses and CTLs: y(t) = βx(t)v(t) − ay(t) − py(t)z(t), v(t) = ky(t) − uv(t), z(t) = cy(t)z(t) − bz(t). (1.1) Here, x(t), y(t), v(t) and z(t) are the concentrations of uninfected cells, infected cells, free viruses and CTLs at time t, respectively. The parameter λ denotes the production rate of uninfected cells. Uninfected cells, infected cells, free viruses and CTLs die at rate dx(t), ay(t), uv(t) and bz(t), respectively. The term βx(t)v(t) is the rate at which uninfected cells are infected by free viruses, and the term py(t)z(t) represents the rate for infected cells to be killed by CTLs. Free viruses are released from infected cells at rate ky(t) and CTL immune response is activated by infected cells at rate cy(t)z(t). We note that system (1.1) assumes that CTL immune response is activated at bilinear rate. However, in [9], De Boer suggested that bilinear rate cannot model several immune responses that are together controlling a chronic infection, and proposed immune response functions with saturation effect. In [10], Jiang and Wang considered saturated immune response function cy(t)z(t)/(h + z(t)) to replace the bilinear rate, here h is a saturation constant.
It is assumed in system (1.1) that as soon as the free viruses enter a target cell, the target cell is immediately infected and new free viruses are produced simultaneously. In [11], a mathematical model with the intracellular phase of the viral life-cycle was first proposed by Herz et al. There is a fixed time delay τ between infection of a cell and production of new free viruses. Some researchers have incorporated the intracellular delay into HIV infection models and investigated the effect of the intracellular delay on HIV infection dynamics (see, for example, [12][13][14][15]).
On the other hand, system (1.1) assumes that the infection rate of virus-to-cell is bilinear. In fact, however, in [16], Ebert et al. observed a nonlinear relationship between parasite dose and infection rate in experiments. Furthermore, in [17], Regoes et al. showed that the infection rate is a sigmoidal function. In [18], Song and Neumann considered the infection rate with saturation effect, βxv p /(1 + αv q ), where p, q and α are positive constants, and investigated the global stability of the viral model with saturated infection rate βxv/(1 + αv).
Inspired by the biological reasons mentioned above, in this paper, we study the effect of saturated infection rate, intracellular delay and saturated CTL immune response on the dynamics of HIV infection. To this end, we consider the following delay differential equations: (1.2) Here, τ denotes the time between viral entry into a cell and production of free viruses; e −mτ is the probability of cellular survival from time t − τ to time t. All parameters of system (1.2) are positive constants.
We assume that the initial condition of system (1.2) satisfies: 2}. According to the basic theory of functional differential equations [19], we can prove that system (1.2) has a unique solution (x(t), y(t), v(t), z(t)) that satisfies the initial condition (1.3). Furthermore, it is easy to show that under the initial condition (1.3), all solutions of system (1.2) are defined on [0, +∞) and are positive for all t ≥ 0.
The main goal of this paper is to make a complete mathematical analysis of system (1.2) and investigate its dynamical behaviors. This paper is organized as follows. In the next section, we calculate the reproduction numbers of system (1.2) and discuss the existence of feasible equilibria. In Section 3, we study the local asymptotic stability of each of feasible equilibria. The global asymptotic stability of each of feasible equilibria is investigated in Section 4. In Section 5, we perform sensitivity analyses to illustrate the effects of parameters on the immunity-activated reproduction number R 1 and the viral load.

Equilibria and reproduction numbers
It is clear that system (1.2) always admits an infection-free equilibrium E 0 (λ/d, 0, 0, 0). Now, we calculate the immunity-inactivated reproduction number of system (1.2). According to the method of next generation matrix proposed by van den Driessche and Watmough [20], we have Sequentially, we obtain that Hence the next generation matrix is We therefore derive the immunity-inactivated reproduction number representing the number of second infected cells produced by a infective cell in a whole susceptible population [20]. It is easy to prove that if R 0 > 1, system (1.2) has a unique immunity-inactivated equilibrium E 1 (x 1 , y 1 , v 1 , 0), here where R 1 is called immunity-activated reproduction number of system (1.2). Besides, we can show that if R 1 > 1, system (1.2) has an immunity-activated equilibrium E 2 (x 2 , y 2 , v 2 , z 2 ), where here,

Local asymptotic stability
In this section, we study the local dynamics of system (1.2).
Proof. By calculation, we have the following characteristic equation of system (1.2) at E 0 : Clearly, Eq. (3.1) has negative real roots ξ = −b and ξ = −d, and other roots are governed by the following equation: Next, we verify that all roots of Eq. (3.3) have negative real parts. By contradiction, let ξ 1 =Reξ 1 +iImξ 1 with Reξ 1 ≥ 0 be a root of Eq. (3.3). It follows that which contradicts (3.3). Hence, we conclude that if R 0 < 1, all roots of Eq. (3.1) have negative real parts. Accordingly, E 0 is locally asymptotically stable. On the other hand, if R 0 > 1, we denote It is easy to see that G(0) = au(1 − R 0 ) < 0 and G(ξ) → +∞ as ξ → +∞. Noting that G(ξ) is a continuous function in respect to ξ, it follows that if R 0 > 1, Eq. (3.1) has a positive real root. Thus, E 0 is unstable.
is locally asymptotically stable.
Proof. By calculation, we have the following characteristic equation of system (1.2) at E 1 : Eq. (3.5) has a negative real root ξ = cy 1 h − b, and the remaining roots are determined by the following equation: We claim that all roots of Eq. (3.7) have negative real parts. Assume the contrary, Eq. (3.7) has a root ξ 2 =Reξ 2 + iImξ 2 with Reξ 2 ≥ 0. It is easy to show that which contradicts (3.7). Consequently, if R 1 < 1 < R 0 , all roots of Eq. (3.5) have negative real parts. Accordingly, E 1 is locally asymptotically stable.
Proof. The characteristic equation of system (1.2) at E 2 is: In the following, we show that all roots of Eq. (3.9) have negative real parts. If Eq. (3.9) has a root ξ 3 with nonnegative real part, we obtain that From the expression of E 2 , we derive that u(a + pz 2 ) = βkx 2 e −mτ 1+αv 2 . Thus, we have (3.10) This leads to a contradiction. Therefore, all roots of Eq. (3.9) have negative real parts if R 1 > 1. Accordingly, E 2 is locally asymptotically stable.

Global asymptotic stability
In this section, we establish the global dynamics of system (1.2). Now, we discuss the boundedness of solutions, and have the following result.

Proof. Let (x(t), y(t), v(t), z(t)) be any solution of system (1.2) satisfying the initial condition (1.3). Denote B(t) = x(t − τ) + e mτ y(t).
Differentiating B(t) with respect to t, we obtain thaṫ which implies that lim sup t→∞ B(t) ≤ λ min{a,d} . Then x(t) and y(t) are ultimately bounded. From the third and the fourth equations of system (1.2), we derive thaṫ and thus lim sup Hence v(t) and z(t) are ultimately bounded. Therefore, the positively invariant set for system (1.2) is given by: Next, we define a function g(x) = x − 1 − ln x, where g(1) = 0 and g(x) attains its minimum at x = 1.

Sensitivity analyses
In this section, we perform sensitivity analyses to show the effects of parameters β, k, u and c on the immunity-activated reproduction number R 1 and the viral load.
We choose the parameter values as follows [21][22][23]:  First, using the methods of Latin Hypercube Sampling and Partial Rank Correlation Coefficients (PRCCs) developed in [24], we investigate the sensitivity of the immunity-activated reproduction number R 1 on the parameters β, k, u and c, which have important effects on HIV infection (see Figure 1). From Figure 1, we see that β, k and c are positively correlated with R 1 while u are negatively correlated with R 1 . In addition, we obtain that decreasing the activation rate of CTL immune response and the infection rate of virus-to-cell infection can more effectively reduce R 1 . Second, we study the sensitivity of the viral load with respect to the parameters β, k, u and c (see Figure 2). In Figure 2, it is easy to see that β, k and c are positively correlated with the viral load while u are negatively correlated with the viral load. Besides, we get that decreasing the release rate of free viruses or increasing the death rate of free viruses can more effectively reduce the viral load, which is helpful in controling the viral infection.