## Abstract

In this paper, we propose a mathematical model for HIV infection with delays in cell infection and virus production. The model examines a viral therapy for controlling infections through recombining HIV with a genetically modified virus. For this model, we derive two biologically insightful quantities (reproduction numbers) *E*_{0} is globally asymptotically stable. If *E*_{s} is globally asymptotically stable. When *E*_{d}, and there exists a constant *R*_{b} such that *E*_{d} is asymptotically stable if

## 1. Introduction

Acquired immune deficiency syndrome (AIDS) was first recognized by the United States Centers for Disease Control and Prevention in 1981 [1]. AIDS is caused by infection with human immunodeficiency virus (HIV), which is transmitted primarily via unprotected sexual contact, contaminated blood transfusions, hypodermic needles and from mother to child during pregnancy, delivery or breastfeeding. With no cure or vaccine in sight, HIV disease continues to be a serious health issue for parts of the world. According to the report [2] by UNAIDS, there were approximately 36.9 (34.3–41.4) million people living with HIV around the world at the end of 2014, and 1.2 (1.0–1.5) million people died from HIV-related causes in 2014. For more information about AIDS, we refer the reader to [1–3] and references therein.

As its name suggests, AIDS causes deficiency of the human immune system, making people infected by HIV more susceptible to common infections that do not usually affect people with a working immune system, because the number of CD4^{+} T cells declines below a centre critical level owing to HIV infection. The entry of HIV into CD4^{+} T cells begins with the interaction between gp120 on its surface with receptor CD4, and coreceptors CCR5 and CXCR4 on the target membrane. Then, the tips of gp41 are inserted into the target membrane, and viral and cellular membranes fuse together [4]. After HIV has bound to the target cell, the HIV RNA and various enzymes are injected into the cell, including reverse transcriptase and integrase. Reverse transcriptase copies the positive single-stranded RNA genome into a complementary DNA molecule, which together with its complement forms a double-stranded viral DNA. The viral DNA is then transported into the cell nucleus, and is integrated into the host cell's genome by another viral enzyme integrase [5]. To produce the virus, the integrated DNA provirus is transcribed into mRNA, and viral proteins are produced by translation of mRNA. Viral genome and proteins assemble together, bud out of the host cell, and a mature HIV virion is produced. During viral replication, CD4^{+} T cells infected with HIV are killed. For the mechanism of CD4^{+} T cell death in HIV infection, see [6,7].

A general model system for HIV infection is described by the following differential equations
*x*(*t*), *y*(*t*), *v*(*t*) is the density of virus-free host cells (mostly corresponding to CD4^{+} T cells), infected cells and the free HIV at time *t*, respectively. In the model (1.1), the healthy host cells are assumed to be produced at rate λ, and to die at rate *d* per cell. Host cells are contacted by the virus at rate *f*(*x*,*v*). It is assumed that it takes an average time *τ*_{1} for the contacting virions to enter cells, which means that the contacted cells become actively affected. Then, after an average time *τ*_{2}, the infected cells start to create and release new virions at rate *k*. The death rate for infected cells and free virus is *a*_{1} and *p*, respectively. The death rate factor for the latent period and the virus production period is *e*^{−a1τ1} and *e*^{−a2τ2}, respectively. Realistically, *a*_{2} may differ from *a*_{1}.

The simplest and earliest form of system (1.1) was derived in [8,9], where *f*(*x*,*v*)=*βxv*, *τ*_{1}=*τ*_{2}=0. The corresponding basic reproduction number *f*(*x*,*v*)=*βxv*, nonlinear incidence rates were also studied, for instance the saturated incidence rate *f*(*x*,*v*)=*βxv*/(1+*b*_{1}*v*) in [15–17], and the Beddington–DeAngelis infection rate *f*(*x*,*v*)=*βxv*/(1+*b*_{0}*x*+*b*_{1}*v*) in [18–21]. For the models mentioned above with these specific nonlinear incidence rates, the corresponding basic reproduction number was identified and its threshold property was discussed.

Based on model (1.1), various HIV models are developed to investigate drug resistance, immune responses and effects of antiretroviral therapy. These researches have contributed to the understanding of HIV biology [22]. We focus on an HIV virotherapy, which is offered by generic engineering to control HIV infections via introducing recombinant virus [23–25]. One application of genetically modified viruses is cell targeting. For instance, the so-called reverse genetics systems can be used to recover rhabdoviruses from cDNA, which makes it possible to genetically engineer rhabdoviruses [23]. In this case, the recombinant virus is targeted to cells infected by HIV, because the recombinant virus is capable of infecting and killing CD4^{+} T cells previously infected by HIV, and does no harm to healthy cells. For example, in [25], a recombinant vesicular stomastitis virus lacking the glycoprotein gene and expressing CD4 and CXCR4 is developed with the property that this virus is unable to infect normal cells, and the cells first infected with HIV-1 are rapidly superinfected with this virus, and killed before high levels of HIV-1 are released.

To model the virotherapy, we add the recombinant virus *w* and double-infected cells *z* into system (1.1), and propose the following model
*αwy*, before turning into double-infected cells. Double-infected cells die at a rate *bz* and release recombinant virus at rate *c* per cell. The death rate of recombinant virus is *q*. We assume that the recombinant virus infection is much faster than HIV infection. So delays are not considered for *z* or *w* in this model.

To the best of our knowledge, only the case *f*(*x*,*v*)=*βxv* has been considered in the literature, such as [26,27]. In [26,27], the corresponding ordinary differential system of (1.2) was studied. Only the structure of the equilibria was analysed in [27], and some numerical simulations were presented. The authors in [26] gave a dynamical analysis on the stability of all three equilibria. The effects of delay *τ*_{1} on the dynamical behaviour were investigated in [28].

In model (1.2), we suppose that the incidence rate depends on *x* and *v*, and is given by a continuous function *f*(*x*,*v*) with continuous derivatives. To be biologically feasible, *f*(*x*,*v*) must satisfy the conditions
*x*, *v*>0. For the biological meaning, the first two conditions (1.3) and (1.4) are obvious. The third one (1.5) means that the incidence rate is a concave function with respect to the number of free HIV. That could be realistic, because when the number of free HIV is so high that any host cells in contact with HIV is virtually certain, the incidence rate will respond more slowly than linearly to the increase in *v*.

Our main goal in this paper is to study the impact of two time delays on lowering the HIV load and increasing the CD4^{+} T cell count. We obtain two reproduction numbers *f*(*x*,*v*).

The rest of this paper is organized as follows. In §2, for system (1.2), we will discuss the well-posedness of the solutions, and the existence of different equilibria. In addition, in order to properly define biologically meaningful equilibria, two reproduction numbers *E*_{0}, single-infection equilibrium *E*_{s} and double-infection equilibrium *E*_{d}. It will be shown that *E*_{0} is globally asymptotically stable for *E*_{s} is globally asymptotically stable for *E*_{d} is asymptotically stable for *R*_{b} is a number larger than 1. A numerical example is presented in §6 to demonstrate the theoretical predictions and to show Hopf bifurcation at *E*_{d}. Finally, conclusion and discussion are drawn in §7.

## 2. Well-posedness and existence of equilibrium points

We introduce the Banach space *θ*∈[−*τ*,0] for any *t*∈[0,*A*], when *A*≥0 and

It is biologically reasonable to consider initial conditions *ϕ*∈*X* for system (1.2). Using the fundamental theory of functional differential equations [29], we have that there is a unique solution

### Theorem 2.1

*Let* *be a solution of system (1.2) satisfying* *. Then,* *is non-negative and uniformly bounded for t≥0. More precisely, we have
**where* *. Furthermore, the solution semiflow* *has a compact global attractor.*

### Proof.

System (1.2) can be rewritten as *ϕ*=(*ϕ*_{1},*ϕ*_{2},…,*ϕ*_{5})∈*X*. It is easy to see that for any *ϕ*∈*X*, *ϕ*_{i}(0)=0 for some *i*, we have *t*≥0 in its maximal interval of existence.

For the boundedness of the solution, we define
*B*(*t*) with respect to time along the solution of (1.2) yields
*x*(*t*), *y*(*t*), *v*(*t*), *z*(*t*) and *w*(*t*) are non-negative. Then, from (2.2) for any *t*≥−*τ*, we have *e*^{−a1τ1}*x*(*t*)≤*B*(*t*), which implies
*y*(*t*), *z*(*t*), *v*(*t*) and *w*(*t*), we get their boundedness listed in (2.1).

Because of the boundedness of the solution, the semiflow *Φ*(*t*) is compact for any *t*>*τ*. On the basis of theorem 3.4.8 in [32], we obtain that *Φ*(*t*) has a compact global attractor in *X*. The proof of the theorem is completed. ▪

Next, we address the basic reproduction number for the model system (1.2) by the next-generation operator approach [33,34]. It is easy to see that system (1.2) has a disease-free equilibrium *E*_{0}=(*x*_{0},0,0,0,0), where *x*_{0}=λ/*d*. Linearizing the system at *E*_{0}, we obtain the following two disease-related equations for variables *y* and *v*
*u*_{1} and *u*_{2} be the number of infected cells and HIV at time *t*=0, respectively. Then, the remaining numbers of infected cells and HIV at time *t* are given by
*M*_{0} is the next infection operator. As usual, the spectral radius of *M*_{0} is called the basic reproduction number *k*/*p*)*e*^{−a1τ1} is the mean number of host cells infected by each virion, and (*e*^{−a2τ2}/*a*_{1})*f*_{v}(*x*_{0},0) is the average number of HIV virions produced by one single-infected cell.

Now, we want to find the other equilibria of system (1.2), beside the disease-free equilibrium point *E*_{0}. From the last three equations in (1.2), we get
*w*=0 or (ii) *y*=(*bq*)/(*cα*).

(i) If *w*=0, then *z*=0, and the first two equations in (1.2) yields

First, we consider the case *x*,*v*) in (2.4) and (2.5), then 0≤*x*<*x*_{0} from (2.4). Further, we can show that *F*(*x*,*v*)<0 for 0≤*x*<*x*_{0} and *v*>0. In fact, from
*F*_{v}(*x*_{0},*v*)<*F*_{v}(*x*_{0},0)<0 for all *v*>0. Noting that *F*(*x*_{0},0)=0, together with *F*_{v}(*x*_{0},*v*)<0, yields *F*(*x*_{0},*v*)<0 for all *v*>0. Finally, we obtain *F*(*x*,*v*)<0 for 0≤*x*<*x*_{0} and *v*>0 owing to *F*_{x}(*x*,*v*)=e^{−a1τ1}*f*_{x}(*x*,*v*)>0 from (1.4). Therefore, when

For the case *x*=*h*(*v*) satisfying (2.5), and *h*′(*v*)≥0. Then that proving the existence of positive solutions (*x*,*v*) for (2.4) and (2.5) is converted to showing that the curve *x*=*h*(*v*) intersects with the straight line *L*:*x*=*x*_{0}−(*a*_{1}*p*)/(*dke*^{−a1τ1−a2τ2})*v* defined by (2.4), as shown in figure 1.

We first show *h*′(*v*)≥0 when *x*=*h*(*v*) exists. Noting that, from the mean value theorem, for any *x* and *v*>0 there exists *v**∈(0,*v*) such that
*x* and *v* satisfying *x*=*h*(*v*). On the other hand, from (2.5), we have
*h*′(*v*)≥0 is straightforward from (2.7).

To prove the existence of function *x*=*h*(*v*) satisfying *F*(*h*(*v*),*v*)=0, we expand *F*(*x*,*v*) for *x*>0 in the form
*e*^{−a1τ1}*f*_{v}(*x*_{0},0)−(*a*_{1}*p*)/(*ke*^{−a2τ2})>0, from (2.8) we have *F*(*x*_{0},*v*)>0 for *v*∈(0,*ε*), *ε*≪1. Note that *F*(0,*v*)<0 for *v*>0. Then given that *F*(*x*,*v*) is an increasing function in *x* from (1.4), for any *v*∈(0,*ε*), there exists only one *x*∈(0,*x*_{0}) satisfying *F*(*x*,*v*)=0. Then, *x*=*h*(*v*) exists for *v*∈(0,*ε*) and *x*∈(0,*x*_{0}).

Because *h*′(*v*)≥0, *x*_{*}<*x*_{0}, because *x*=*h*(*v*)<*x*_{0} for *v*∈(0,*ε*). In order to make sure that the curve *x*=*h*(*v*) can intersect with the straight line *L*, we still need to study the maximal interval of existence for the function *x*=*h*(*v*).

Let *x*_{ε} is infinite, then (0,*ε*) is the largest interval we can find. Otherwise, we have *F*(*x*_{ε},*ε*)=0 from the continuousness of *F*(*x*,*v*) in *x* and *v*. Because (∂*F*/∂*x*)(*x*_{ε},*ε*)=*e*^{−a1τ1}*f*_{x}(*x*_{ε},*ε*)>0 from (1.4), by implicit function theorem, we have *δ*>0 such that function *v*∈(*ε*−*δ*,*ε*+*δ*) exists and satisfies *x*=*h*(*v*) can be extended to interval (0,*ε*+*δ*). Using this process repeatedly, we can find the maximal interval (0,*N*) for *x*=*h*(*v*). Denote *x** or *N* is infinite.

Then, we claim that if *v*_{s},*x*_{s}) in figure 1 for the curves defined by (2.4) and (2.5). Therefore, if *E*_{s}=(*x*_{s},*y*_{s},*v*_{s},0,0) exists, where *y*_{s}=(*p*)/(*ke*^{−a2τ2})*v*_{s}, *x*_{s} and *v*_{s} satisfy equations (2.4) and (2.5).

(ii) Next, we consider the case *w*>0. Let *v*_{d}=(*ke*^{−a2τ2}/*p*)*y*_{d}. The first two equations in (1.2) yield
*g*_{1}(*x*,*w*)=*e*^{−a1τ1}(λ−*dx*)−*a*_{1}*y*_{d}−*αy*_{d}*w*, *g*_{2}(*x*)=λ−*dx*−*f*(*x*,*v*_{d}). From the first equation in (2.9), we see that *w*>0 if and only if *g*_{1}(*x*,0)>0, which yields
*g*_{2}(*x*) is a decreasing function in *x*, and *g*_{2}(0)=λ>0. Then, equations in (2.9) have positive solutions satisfying (2.10) if and only if *M*>0 and *g*_{2}(*M*)<0, that is

Before we simplify the condition (2.11), we shall first show that equilibrium *E*_{s} exists when (2.11) holds. From (1.3) and the mean value theorem in multiple variables, there exists *θ*∈(0,1) such that
*x*_{θ}=λ/*d*−(*θa*_{1}*y*_{d})/(*de*^{−a1τ1}), *v*_{θ}=*θv*_{d}. Because *v*_{d}=(*ke*^{−a2τ2}/*p*)*y*_{d} and *f*(λ/*d*,0)=0, from (2.11) and (2.12), we can obtain *f*_{vx}(*x*,0)≥0, and further with (1.5) we have *f*_{v}(*x*_{θ},*v*_{θ})≤*f*_{v}(*x*_{θ},0)≤*f*_{v}(*x*_{0},0) for *x*_{θ}<*x*_{0} and *v*_{θ}>0. Then, *R*_{1}>1 from (2.13). Hence, inequality (2.11) implies *E*_{s} exists.

On the other hand, (2.11) could be rewritten as
*v*_{d},*M*) satisfying (2.14) is on the straight line *L* and above the curve *x*=*h*(*v*), which means that (2.11) is equivalent to *v*_{d}<*v*_{s}. Because *y*_{d}=(*p*)/(*ke*^{−a2τ2})*v*_{d} and *y*_{s}=(*p*)/(*ke*^{−a2τ2})*v*_{s}, we have *y*_{d}<*y*_{s}, i.e. *M*>0, because

Therefore, if and only if *x*_{d},*w*_{d}) for (*x*,*w*), and system (1.2) has a double-infection equilibrium *E*_{d}=(*x*_{d},*y*_{d},*v*_{d},*z*_{d},*w*_{d}), where *z*_{d}=(*q*/*c*)*w*_{d}.

For the biological meaning of *y*_{d}=(*bq*)/(*cα*). It is easy to see that *c*/*q* is the average number of recombinant virus that a double-infected cell produces, and *αy*_{s}/*b* gives the mean number of double-infected cells caused by each recombinant virus when the number of single-infected cells stabilizes at *y*_{s}. Then,

Obviously, as *x*_{d},*v*_{d}) and (*x*_{s},*v*_{s}) satisfy λ−*dx*−*f*(*x*,*v*)=0 from the first equation in system (1.2), we get *E*_{d} bifurcates from *E*_{s} at

Summarizing the above discussion, we have the following theorem.

### Theorem 2.2

*System (1.2) has three possible biologically meaningful equilibria: disease-free equilibrium E*_{0}*=(x*_{0}*,0,0,0,0), with x*_{0}*=λ/d, single-infection equilibrium E*_{s}*=(x*_{s}*,y*_{s}*,v*_{s}*,0,0), with y*_{s}*=(p)/(ke*^{−a2τ2}*)v*_{s}*, x*_{s} *and v*_{s} *satisfy (2.4) and (2.5), and double-infection equilibrium E*_{d}*=(x*_{d}*,y*_{d}*,v*_{d}*,z*_{d}*,w*_{d}*), with y*_{d}*=(bq)/(cα), v*_{d}*= (ke*^{−a2τ2}*/p)y*_{d}*, z*_{d}*=(q/c)w*_{d}*, x*_{d} *and w*_{d} *satisfy (2.9). More specifically, (i) if* *E*_{0} *is the only equilibrium; (ii) if* *the single-infection equilibrium E*_{s} *exists and (iii) the double-infection equilibrium E*_{d} *exists if and only if*

## 3. Global stability of equilibrium *E*_{0}

Here, we shall study the global stability of the disease-free equilibrium point *E*_{0}. We have the following theorem.

### Theorem 3.1

*If* *the disease-free equilibrium E*_{0}*=(x*_{0}*,0,0,0,0) is globally asymptotically stable, implying that none of the two virus can invade, regardless of the initial load. If* *then E*_{0} *becomes unstable.*

### Proof.

First, we recall that *E*_{0} is the only equilibrium when *E*_{0}, we construct the following Lyapunov function:
*V* _{0} with respect to time *t* along the solution of system (1.2) can be expressed as
*f*_{vx}(*x*,0)≥0, we have *f*_{v}(*x*,0)≥*f*_{v}(*x*_{0},0) if *x*≥*x*_{0}, and *f*_{v}(*x*,0)≤*f*_{v}(*x*_{0},0) if *x*≤*x*_{0}. Then
*f*(*x*,*v*)=*f*_{v}(*x*,*v**)*v*≤*f*_{v}(*x*,0)*v*, 0≤*v**≤*v*. Then
*V*_{0}/d*t*|_{(1.2)}≤0 and the equality holds for *x*=*x*_{0}, *v*=*w*=0. Thus, by LaSalle's invariance principle [35], we conclude that *E*_{0} is globally asymptotically stable.

For the unstability of *E*_{0}, we have the linearized system of (1.2) at *E*_{0} given by
*E*_{0}, it suffices to only consider the zeros of the following function

When *E*_{0} is unstable. ▪

## 4. Global stability of the single-infection equilibrium *E*_{s}

From the analysis given in §2, we know the single-infection equilibrium *E*_{s}=(*x*_{s},*y*_{s},*v*_{s},0,0) exists when *E*_{s}, we have the following persistence result.

### Theorem 4.1

*Let X*_{0}*={ϕ=(ϕ*_{1}*,ϕ*_{2}*,…,ϕ*_{5}*)∈X:ϕ*_{2}*(0)>0 and ϕ*_{3}*(0)>0}, and denote ∂X*_{0}*=X∖X*_{0}*={ϕ∈X:ϕ*_{2}*(0)=0 or ϕ*_{3}*(0)=0}. When* *system (1.2) is uniformly persistent with respect to (X*_{0}*,∂X*_{0}*) in the sense that there exists some η>0 such that*

### Proof.

By the form of system (1.2), it is easy to see that *X*_{0} is positively invariant. We set *M*_{∂}={*ϕ*∈*X*:*Φ*(*t*)*ϕ*∈∂*X*_{0},∀*t*≥0}. Clearly, *M*_{∂}={*ϕ*∈*X*:*ϕ*_{2}(0)= 0,*ϕ*_{3}(0)=0}.

We claim that *W*^{s}(*E*_{0})∩*X*_{0}=∅. Assume that, on the contrary, there exists *ψ*∈*X*_{0} such that *ε*>0, there exists a positive constant *T*_{0}=*T*_{0}(*ε*), such that for *x*(*t*)>*x*_{0}−*ε*, *v*(*t*)<*ε* and *w*(*t*)<*ε* for all *t*≥*T*_{0}. Here, because *ε* small enough such that
*t*≥*T*_{0}+*τ*_{1}, for *t*≥*T*_{0}+*τ*_{1}, from (4.2), we have
*ξ*_{0} is the principal eigenvalue of the following linear cooperative system
*ξ*_{0}>0 from [30], corollary 5.5.2. Let *u*_{s}=(*u*_{y},*u*_{v})^{T} be the positive right eigenvector associated with *ξ*_{0} for system (4.3). We choose *l*>0 small enough such that *lu*_{y}*e*^{ξ0t}≤*y*(*t*,*ψ*), *lu*_{v}*e*^{ξ0t}≤*v*(*t*,*ψ*), for all *t*∈[*T*_{0}+*τ*,*T*_{0}+2*τ*]. Obviously, *le*^{ξ0t}*u*_{s} satisfies (4.3) for all *t*≥*T*_{0}+*τ*. Then by the comparison principle, we get (*y*(*t*,*ψ*),*v*(*t*,*ψ*))^{T}≥*le*^{ξ0t}*u*_{s} for all *t*≥*T*_{0}+*τ*. Because *lu*_{s}>0 and *ξ*_{0}>0, letting

Define a continuous function *ϕ*∈*X*. Then, *p*_{1}(*Φ*(*t*)*ϕ*)>0 if either *p*_{1}(*ϕ*)=0 and *ϕ*∈*X*_{0}, or if *p*_{1}(*ϕ*)>0. Thus, *p*_{1} is a generalize distance function for the solution semiflow *Φ*(*t*) [36]. We obtain that *E*_{0} is a compact and isolated invariant sets in ∂*X*_{0}, and *E*_{0} forms a cycle in ∂*X*_{0}. From the claim above, *E*_{0} is isolated in *X*, and *W*^{s}(*E*_{0})∩*X*_{0}=∅. By [36], theorem 3, it follows that there exists *η*>0 such that *ϕ*∈*X*_{0}, which implies *Φ*(*t*) is uniformly persistent with respect to (*X*_{0},∂*X*_{0}). Thus, we have *ω*(*ϕ*)⊂*X*_{0} for any *ϕ*∈*X*_{0}. ▪

Further, we have the following result of the global stability at *E*_{s}.

### Theorem 4.2

*If* *and* *then the single-infection equilibrium E*_{s} *is globally asymptotically stable, implying that the recombinant virus cannot survive but the pathogen virus can. E*_{s} *becomes unstable when*

### Proof.

We construct the Lyapunov function *V* _{s}=*V* _{1}+*e*^{−a1τ1}*f*(*x*_{s},*v*_{s})*V* _{2}, where
*E*_{s} satisfies the following relations
*V* _{2},
*f*(*x*,*v*_{s})≥*f*(*x*_{s},*v*_{s}) when *x*≥*x*_{s}, and *f*(*x*,*v*_{s})≤*f*(*x*_{s},*v*_{s}) when *x*≤*x*_{s}. Then
*f*(*x*,*v*) is a concave function in *v* from (1.5). Then for any *v*>0 and any *s* in (0,1], *f*(*x*,(1−*s*)0+*sv*)≥(1−*s*)*f*(*x*,0)+*sf*(*x*,*v*), which implies
*f*(*x*,*v*) should satisfy
*a*_{i} and *b*_{i}, because the function *x*>0, and *g*_{s}(*x*)=0 if and only if *x*=1. Thus, we have

Therefore, d*V* _{s}/d*t*|_{(1.2)}≤0 when *x*=*x*_{s}, *y*=*y*_{s}, *v*=*v*_{s} and *w*=0. Then by LaSalle's invariance principle [35], we conclude that *E*_{s} is globally asymptotically stable when

When *E*_{s}, we calculate the linearized system of (1.2) at *E*_{s}, and obtain
*D*_{1}(*ξ*)*D*_{2}(*ξ*)=0, where
*D*_{1}(*ξ*) in *ξ* can be expanded as
*D*_{1}(*ξ*)=0 has two zeros with negative real part if and only if *D*_{1}(*ξ*) has two real roots with different signs. Therefore, *E*_{s} is unstable if

From the proof of theorems 3.1 and 4.2, it is easy to get the following corollary.

### Corollary 4.3

When *1.1*); when *1.1*).

## 5. Stability of the double-infection equilibrium *E*_{d}

The double-infection equilibrium *E*_{d} comes into existence for *E*_{d}, for any quantity *A* involving *τ*_{1} and *τ*_{2} in the paper, we denote by Å the value of *A* when *τ*_{1}=*τ*_{2}=0. We have the following result for the local stability of *E*_{d}.

### Theorem 5.1

*For system (1.2), there exists an R*_{b}*>1 such that the double-infection equilibrium E*_{d} *is asymptotically stable for*

### Proof.

First, we recall that *E*_{d} exists if and only if *E*_{d}=(*x*_{d},*y*_{d},*v*_{d},*z*_{d},*w*_{d}) is
*F*_{d}=*d*+*f*_{x}(*x*_{d},*v*_{d}) and *A*_{w}=*a*_{1}+*αw*_{d}.

By straightforward but tedious algebraic manipulations, we obtain the characteristic equation of (5.1), given by

When *τ*_{1}=*τ*_{2}=0, (5.2) becomes
_{i}, *i*=2,3,4, because *C*_{0}>0 and *C*_{4}>0. But it is not easy to determine them for general *ẘ*_{d} when *ẘ*_{d}=0, we have
*y*_{s}=(*p*)/(*ke*^{−a2τ2})*v*_{s}, from the second equation in (1.2) we have (*a*_{1}*p*)/(*ke*^{−a1τ1−a2τ2})*v*_{s}−*f*(*x*_{s},*v*_{s})=0. Then, *f*(*x*_{s},*v*_{s})≥*f*_{v}(*x*_{s},*v*_{s})*v*_{s} from (1.5) and (2.6). Thus, from (5.4) and (1.4) *i*=2,3,4. Because *C*_{i}, *i*=0,…,4, are meaningful only for *R*_{ε}) around _{i}>0 when *i*=2,3,4. Therefore, all roots of (5.3) have negative real part when

If at least one of *τ*_{i}≠0 for *i*=1,2, it is easy to see *ξ*=0 is not a zero of (5.2) because *A*_{0}>0. Moreover, there are no roots for (5.2) existing as *τ*_{1} and *τ*_{2}, the only possibility that the roots of (5.2) enter into the right half plane is to cross the imaginary axis as *τ*_{1} and *τ*_{2} increase. Suppose a purely imaginary number *ξ*=*iϖ*, (*ϖ*>0), is a root of (5.2). Then substituting *ξ*=*iϖ*, *ϖ*>0 into *D*(*ξ*)=0 yields
*H*(*ϖ*^{2})=*ϖ*^{10}+*h*_{1}*ϖ*^{8}+*h*_{2}*ϖ*^{6}+*h*_{3}*ϖ*^{4}+*h*_{4}*ϖ*^{2}+*h*_{5}=0, with
*h*_{1}>0 and *h*_{5}>0. For *h*_{i}, *i*=2,3,4, we use the continuity argument again, and have *h*_{i}>0, *i*=1,…,5, which implies that *H*(*ϖ*^{2}) does not have any positive real roots. Therefore, combining with the condition 1<*R*_{z}<*R*_{ϵ}, let *E*_{d} is locally asymptotically stable. ▪

Besides the local stability of *E*_{d}, we have the following uniform persistence result with respect to the recombinant virus and double-infected cells.

### Theorem 5.2

*If* *then there is an η>0 such that any solution* *of the system (1.2) with ϕ∈X, ϕ*_{4}*(0)>0 and ϕ*_{5}*(0)>0 satisfies
*

### Proof.

Define *Y* _{0}={*ϕ*∈*X*:*ϕ*_{4}(0)>0 and *ϕ*_{5}(0)>0}. Then, we have ∂*Y* _{0}=*X*∖*Y* _{0}= {*ϕ*∈*X*:*ϕ*_{4}(0)=0 or *ϕ*_{5}(0)=0}. Define *N*_{∂}={*ϕ*∈*X*:*Φ*(*t*)*ϕ*∈∂*Y* _{0}, ∀*t*≥0}.

By a similar argument as that in the proof of theorem 4.1, we can show that when *W*^{s}(*E*_{0})∩*Y* _{0}=∅. We also claim that there exists a *δ*>0, such that any *ϕ*∈*Y* _{0},

Again, assume that on the contrary, there exists *ψ*∈*Y* _{0} such that *ε*>0, there exists a positive constant *T*_{1}=*T*_{1}(*ε*), such that *y*(*t*)>*y*_{s}−*ε* for all *t*≥*T*_{1}. Here, because *ε* small enough such that
*t*≥*T*_{1}, in (1.2) we have *ξ*_{1}>0 is the positive eigenvalue, and *u*_{d}=(*u*_{z},*u*_{w})^{T} be the corresponding positive right eigenvector. We choose *l*>0 small enough such that *lu*_{z}*e*^{ξ1t}≤*z*(*t*), *lu*_{w}*e*^{ξ1t}≤*w*(*t*), for all *t*∈[*T*_{1},*T*_{1}+*τ*]. Obviously, *le*^{ξ1t}*u*_{d} satisfies (5.6) for all *t*≥*T*_{1}. Then by the comparison principle, we get (*z*(*t*),*w*(*t*))^{T}≥*le*^{ξ1t}*u*_{d} for all *t*≥*T*_{1}+*τ*. Because *lu*_{d}>0 and *ξ*_{1}>0, letting *W*^{s}(*E*_{s})∩*Y* _{0}=∅, when

Next, we claim *ϕ*∈*N*_{∂}, i.e. *Φ*(*t*)*ϕ*∈∂*Y* _{0}, we have *z*(*t*,*ϕ*)≡0, or *w*(*t*,*ϕ*)≡0. From the *w* equation in system (1.2), we have *z*(*t*)≡0, or *w*(*t*)≡0. Hence, we have *ω*(*ϕ*)=*ω*_{1}×{(0,0)} for some *ω*_{1}∈*C*([−*τ*,0];*R*^{3}_{+}), and
*Φ*_{1}(*t*) is the solution semiflow associated with system (1.1). From corollary 4.3, we have that *ω*_{1} is either

Define a continuous function *ϕ*∈*X*. Then, *p*_{2} is also a generalize distance function for the solution semiflow *Φ*(*t*). From the proof above, we conclude that any forward orbit of *Φ*(*t*) in *N*_{∂} converges to *E*_{0} or *E*_{s}, that *E*_{0} and *E*_{s} are two isolated invariant sets in *X*, and (*W*^{s}(*E*_{0})∩*W*^{s}(*E*_{1}))∪*Y* _{0}=∅. Moreover, it is easy to see that no subset of {*E*_{0},*E*_{s}} forms a cycle in ∂*Y* _{0}. By [36], theorem 3, it follows that there exist *η*>0 such that *ϕ*∈*Y* _{0}, which implies *Φ*(*t*) is uniformly persistent with respect to (*Y* _{0},∂*Y* _{0}). ▪

## 6. Numerical simulations

In the above discussions, owing to the general form of *f*(*x*,*v*), we cannot obtain the explicit form of *R*_{z}. Consequently, we are not able to either prove the global stability of the third equilibrium point *E*_{d} or determine whether there are other dynamic phenomena around *E*_{d} for *E*_{d}, including the convergence of orbits to *E*_{d} and the existence of Hopf bifurcation for system (1.2).

Because there are more parameters involved in the model if a nonlinear incidence function is used, we cannot find the proper value range for the new parameter in the literature. A non-reliable parameter value could damage the biological interpretation of the numerical simulations. So we chose a bilinear incidence function *f*(*x*,*v*)=*βxv*, where *β* is the constant rate at which a T cell is contacted by the virus, which is widely used in other papers. Then

For computer simulation, we set the parameter values as the following: λ=1 cell mm^{−3}, *a*_{1}=0.5 day^{−1}, *b*=2 day^{−1}, *p*=*q*=3 day^{−1}, *k*=80 vir cell^{−1}, *c*=1800 vir cell^{−1}, see [27,38]. Let *τ*_{1} be the bifurcation parameter.

The disease-free equilibrium *E*_{0} is now given by *E*_{0}=(180,0,0,0,0), which is globally asymptotically stable from theorem 3.1 for *τ*_{1}>*τ*_{s}≈6.76767349288, i.e. *τ*_{1}<*τ*_{s}, *E*_{0} becomes unstable, and the single-infection equilibrium *E*_{s} exists, given by *τ*_{s}>*τ*_{1}>*τ*_{d}≈1.5217799236.

Further decreasing *τ*_{1} to pass through the critical value *τ*_{d} will cause *E*_{s} to lose its stability, and give rise to the double-infection equilibrium,
*D*(*ξ*,*τ*_{1}) at *E*_{d}. Solving *D*(*iη*,*τ*_{1})=*R*(*η*,*τ*_{1})+*iS*(*η*,*τ*_{1})=0 yields (*η*_{h},*τ*_{h})≈(±0.58060139097,0.17431498237). It follows from theorem 5.1 that *E*_{d} is asymptotically stable when *τ*_{d}>*τ*_{1}>*τ*_{h}, where *τ*_{h}≈0.17431498237. The simulations for *τ*_{1}=0.5,0.9,1.2 are shown in figure 2.

Next, we consider possible Hopf bifurcation. The following condition is held
*D*(*ξ*,*τ*_{1})=0 has a pair of purely imaginary roots at *τ*_{1}=*τ*_{h}, whose real parts become positive when *τ*_{1}<*τ*_{h}, implying existence of a Hopf bifurcation. At the critical point, *τ*_{1}=*τ*_{h}, *E*_{d} loses its stability through a Hopf bifurcation, giving rise to limit cycles (figure 3). When *τ*_{1}=*τ*_{h}, we obtain

To sum up, the bifurcation diagram projected on *y*−*τ*_{1} plane is given in figure 4, which shows what impacts the delay *τ*_{1} could have on the dynamics around the equilibria of model (1.2).

## 7. Discussion

In this paper, we propose an HIV model with a general nonlinear incidence rate and two time delays. Two production numbers *E*_{0}=(*x*_{0},0,0,0,0) is globally asymptotically stable. When *E*_{0} becomes unstable, and the single-infection equilibrium *E*_{s}=(*x*_{s},*y*_{s},*v*_{s},0,0) occurs. When *E*_{s} is globally asymptotically stable. At *E*_{s} bifurcates into double-infection equilibrium *E*_{d}=(*x*_{d},*y*_{d},*v*_{d},*z*_{d},*w*_{d}), and *E*_{s} loses its stability for *η*>0 such that *ϕ*∈*X* with *ϕ*_{4}(0)>0 and *ϕ*_{5}(0)>0 when *E*_{d} is asymptotically stable for

From the expression of *τ*_{1} and *τ*_{2} leads to overestimation of the basic reproduction number *y*_{s} and *v*_{s} with respect to *τ*_{1} and *τ*_{2}, respectively. From (2.4) and (2.5), we get the following equation
*j*=1,2, where *y*_{s}=*p*/*a*_{1}*e*^{−a2τ2}*v*_{s}, we have
*τ*_{1} or *τ*_{2} is not included in system (1.2). Similarly, we can easily get d*v*_{d}/d*τ*_{j}<0 and d*y*_{d}/d*τ*_{j}<0 for *j*=1,2.

From the simulations and figure 4, it is easy to see that choosing different values for delays could change the dynamic behaviours, not only quantitatively, but also sometimes qualitatively. So intracellular delays should be included in the modelling of HIV infection. We should mention that some results (theorems 3.1 and 4.2) in this manuscript still hold if the system has no time delay or if function *f* is bilinear. The new dynamics is mainly derived with the introduction of the new variables *z* and *w*, see corollary 4.3. However, when we release some conditions, e.g. (1.5) for the nonlinear incidence function *f*, the system will become much more complicated, multiple steady-state solutions and multistabilities may exist. This is beyond the scope of this manuscript.

Note that systems (1.1) and (1.2) share the same basic reproduction number *R*_{1} is given in (2.13), *E*_{s} is globally asymptotically stable in (1.2), just as *E*_{d} comes into existence. From §2, we see *x*_{d}>*x*_{s}, *y*_{d}<*y*_{s} and *v*_{d}<*v*_{s}, implying that the virotherapy cannot only decrease HIV load and the number of infected cells by HIV, but also increase the healthy CD4^{+} T cell count. So when the recombinant virus can survive, i.e.

Because *E*_{d}. From numerical simulations we can see that there are some phenomena we should pay attention to when *τ*_{1} becomes smaller, *E*_{d}, and the amplitude of oscillations also becomes larger. Furthermore, relative large

On the other hand, we have *y*_{s}. Obviously, *c*, *α*, *b* and *q* in *z* and *w* equations. Although the explicit expression of *y*_{s} cannot obtained from system (1.2) with general nonlinear incidence function *f*(*x*,*v*), we can see that *y*_{s} is totally determined by system (1.1). In other words, *y*_{s} is not affected with the varies of the recombinant virus *w* and double-infected cells *z*. Therefore, increasing *c* and *α*, or decreasing *b* and *q* can make

## Data accessibility

This work does not have any experimental data.

## Authors' contributions

The work was done when Y.T. visited Y.Y. in Memorial University. Y.T. and Y.Y. conceived the mathematical models and wrote the paper. All authors gave final approval for publication.

## Competing interests

We have no competing interests.

## Funding

This work was supported in part by Natural Sciences and Engineering Research Council of Canada (NSERC) grant no. 203786-46310-2000.

## Acknowledgements

The authors greatly appreciate the anonymous referees’ careful reading and valuable comments. Their critical comments and helpful suggestions greatly improve the presentation of the manuscript.

- Received September 4, 2015.
- Accepted January 14, 2016.

- © 2016 The Author(s)

Published by the Royal Society. All rights reserved.