# Hopf Bifurcation for a Model of HIV Infection ofCD4+ T Cells with Virus Released Delay

Authors: Jun-Yuan Yang; Xiao-Yan Wang; Xue-Zhi Li
Journal: Discrete Dynamics in Nature and Society
 (2011)

---

## Abstract

A viral model of HIV infection ofCD4+ T-cells with virus
released period is formulated, and the effect of this released period on the stability of the equilibria is investigated. It is shown that the introduction of the viral released period can destabilize the system, and the period solution may arise. The direction and stability of the Hopf bifurcation are also discussed. Numerical simulations are presented to illustrate the results.


---

## Body

## 1. Introduction and Model Formulation

In the last decade, many mathematical models have been developed to describe the immunological response to infection with human immunodeficiency virus (HIV) (see [1–11]). Simple HIV models have played a significant role in the development of a better understanding of the disease and the various drug therapy strategies used against it. Perelson et al. in [1] proposed a basic mathematical model to describe spread of HIV. Many other models [12–14] which take the model proposed in [1] as their inspiration have been formulated. Zhou et al. in [5] discussed the following ODE model:
(1.1)dTdt=s-dT+aT(1-TTmax⁡)-βTV+ρI,dIdt=βTV-(δ+ρ)I,dVdt=qI-cV.
In (1.1), T(t),I(t) represent, respectively, the concentration of healthy CD4+ T cells and infected CD4+ T cells at time t, and V(t) represents the concentration of free HIV at time t. These parameters are defined as follows: s is the source of CD4+ T cells precursors, d is the natural death rate of CD4+ T cells, a is their growth rate (thus, a>d in general), and Tmax is their carrying capacity, β is the contact rate between uninfected CD4+ T cells and virus particles, δ is a blanket death term for infected cells, c is the clearance rate constant of virus, q is the lytic death rate for infected cells, and ρ is cure rate from infected cells to healthy cells. Parameters d, a, δ, β, ρ, c, Tmax, q, and s are positive values. They obtained the conditions for which system (1.1) exists an orbitally asymptotically stable periodic solution.Time delays of one type or another have been introduced to describe the time between viral entry into a target cell and the production of new virus particles by many authors (see [6–9]). Culshaw and Ruan in [6] introduced a discrete time delay to the model to describe the time between infection of a CD4+ T-cell and the emission of viral particles on a cellular level. They discussed locally asymptotically stable and obtained existence of the Hopf bifurcation under some conditions. Herz et al. in [3] used a discrete delay to model the intracellular delay that would substantially shorten the estimate for the half-life of free virus.However, almost all ofCD4+ models were discussed and the discrete delay had denoted the time between infection of a CD4+ T-cell and the emission of viral particles on a cellular level. According to reports of CDC government, they believed that there exists a “window period” when the infected cell released virus. The time interval between point of infection and detection of a seroconversion using US FDA-licensed third-generation antibody assays averages 22 days. In primary HIV infection, a localized viral replication (eclipse) takes place first, and lasts for 1–4 weeks [15]. If the amount of virus released does not attend a certain level, a patient infected by infected CD4+ can not be examined. Therefore, we introduce a discrete delay denoted the “window period” in their model. Efficacy of the inhibition of viral replication is imposed upon the virus-host system by nucleoside analogue therapy. Based on the work of Zhou et al. (see [5]), the HIV infection model with delay can be written as follows: (1.2)dTdt=s-dT+aT(1-TTmax⁡)-βTV+ρI,dIdt=βTV-(δ+ρ)I,dVdt=qI(t-τ)-cV.
All the biological meanings of parameters are same as system (1.1). The time delay is introduced in the system describing the dynamics of “window period,” that is, the released virus term of infected CD4+ T cells is changed from ρI to qI(t-τ).The initial conditions of system (1.2) are (1.3)T(0)=φ10,I(θ)=φ2(θ),V(0)=φ30,θ∈[-τ,0],φ10≥0,φ2(0)≥0,φ30≥0.
With a standard argument given in [16], it easy to show that the solution T(t), I(t), V(t) with initial conditions (1.2) exists and is unique for all t≥0.Theorem 1.1.
The solution of system (1.2) is positive and boundary.Proof.
First we proof the positivity of the solution of (1.2). From the first equation of (1.2), we obtain
(1.4)dTdt≥aT(1-TTmax⁡)-dT-βTV.
Solving it, we obtain
(1.5)T(t)≥T(0)exp⁡∫0t[a(1-T(t)Tmax⁡)-βV(t)]dt-dt≤0.
From the second and the third equation, we obtain
(1.6)I(t)≥I(0)exp⁡-(δ+ρ)t≥0,V(t)≥V(0)exp⁡-ct≥0.
Therefore, the solution of (1.2) is positive.
Second, we will prove the boundary of the solution of (1.2). Define d̂=min{d,δ} and L=T+I. Adding the first equation and the second equation of (1.2), it is easy to see(1.7)dLdt≤s+aT(1-T(t)Tmax⁡)-d̂L;
if L≥Tmax, (1.7) become the following:
(1.8)dLdt≤s-d̂L.
It is easy to get L≥s/d̂≐M1. If L≤Tmax, (1.7) become the following:
(1.9)dLdt≤s+aT(1-T(t)Tmax⁡)-d̂L≤s+aL(1L(t)Tmax⁡)-d̂L.
It leads to L≤2max{s/d̂,Tmax}≐M2. Hence, for ε>0 sufficiently small, there is a T1 such that if t>T1, T(t)≤M+ε, I(t)≤M2+ε.
Again, forε>0 sufficiently small, we drive from the third equation of system (1.2) that for t>T1+τ(1.10)dVdt≤q(M2+ε)-cV.
A comparison argument shows that
(1.11)V≤q(M2+ε)c.
Since this is true for arbitrary ε>0, it follows that V≤q(M2)/c. If we set M=max{M2,q(M2)/c}, such that S≤M, I≤M, V≤M. Hence, the solution of system (1.2) is boundary.


## 2. Equilibria Stability and Hopf Bifurcation

First we find all biologically feasible equilibria admitted by the system (1.2) and then study the dynamics of the system around each equilibrium. We introduce the reproduction number of differential-delay model (1.2), which is given by a similar expression (2.1)R0=T0T¯=βqT0c(δ+ρ).
The R0 stands if one virus is introduced in a population of uninfected cells which infect the total number of secondary infectious during their infectious period 1/c(δ+ρ).The equilibria of the system (1.2) are as follows:(i)
uninfected equilibriumE0=(T0,0,0), where T0=(Tmax/2a)(a-d+(a-d)2+(4as/Tmax));(ii)
an infected equilibriumE¯=(T¯,I¯,V¯), which exists if R0>1, where(2.2)T¯=c(δ+ρ)βq,I¯=1δ[s-dT¯+aT¯(1-T¯Tmax⁡)],V¯=qcI¯.
Following the analysis in [5], we can see that, if R0>1, then the infection-free equilibrium E0 is unstable, and incorporation of a delay will not change the instability.Next, we focus on investigating the local stability and Hopf bifurcation of the positive equilibrium of (1.2). To study the local stability of the positive equilibrium E¯=(T¯,I¯,V¯), we consider the linearization of system (1.2) at the point E¯. Let us define (2.3)T(t)=x+T¯,I(t)=y+I¯,V(t)=z+V¯.
The following transcendental characteristic equation is obtained: (2.4)λ3+a1λ2+a2λ+a3=e-λτ(b1λ+b2),
where the coefficients in this equation are expressed as follows: (2.5)a1=c+δ+ρ+d-a+2aT¯Tmax⁡+βV¯,a2=c(δ+ρ)+(c+δ+ρ)(d-a+2aT¯Tmax⁡+βV¯)-ρβV¯,a3=c(δ+ρ)(d-a+2aT¯Tmax⁡+βV¯)-cβρV¯,b1=qβT¯,b2=qβT¯(d-a+2aT¯Tmax⁡).
For τ=0, the characteristic equation (2.4) reduces to the following (2.6)λ3+a1λ2+(a2-b1)λ+a3-b2=0.
By the Routh-Hurwitz Criterion, it follows that all eigenvalues of (2.6) have negative real parts if and only if (2.7)a1>0,a2-b1>0,a3-b2>0,a1(a2-b1)-(a3-b2)>0.
If R0>1, and d-a+(2aT¯/Tmax)>0, (2.8)a1>0,a2-b1=((c+δ+ρ)(d-a+2aT¯Tmax⁡)+(c+δ)βV¯)>0,a3-b2=cδβV¯>0,a1(a2-b1)-(a3-b2)=[c+δ+ρ+d-a+2aT¯Tmax⁡+βV¯]((c+δ+ρ)(d-a+2aT¯Tmax⁡)+(c+δ)βV¯)-cδβV¯>0.
We know that all eigenvalues of (2.6) have negative real parts.Clearly, ifλ=iw with w>0 is a root of (2.4). This is the case if and only if w satisfies the following equation: (2.9)-iw3-a1w2+ia2w+a3=(cos⁡wτ+isinwτ)(ib1w+b2).
Separating the real and imaginary parts, we have the following system that must be satisfied: (2.10)a3-a1w2=b2cos⁡wτ-b1wsinwτ,a2w-w3=b1wcos⁡wτ+b2sinwτ.
We eliminate the trigonometric functions by squaring both sides of each equation above and adding the resulting equations. We obtain the following sixth degree equation for w: (2.11)w6+(a12-2a2)w4+(a22-2a1a3-b12)w2+a32-b22=0.
Letting z=w2, then (2.11) becomes a third order equation in z: (2.12)z3+m1z2+m2z+m3=0,
where we have used the following notation for the coefficients of (2.12): (2.13)m1=a12-2a2,m2=a22-2a1a3-b12,m3=a32-b22.
In order to show that the endemic equilibrium E¯ is locally stable we have to show that (2.12) does not have a positive real solution which might give the square of w, that is, that (2.4) cannot have purely imaginary solutions. The lemma below establishes conditions leading to that result.Lemma 2.1 (see [10]).
For the polynomial equation (2.12) (i)
ifm3<0, then (2.12) has at least one positive root;(ii)
ifm3≥0 and Δ=m12-3m2≤0, then (2.12) has no positive root;(iii)
ifm3≥0 and Δ=m12-3m2>0, then (2.12) has positive root if and only if z1*=(-m1+Δ)/3, and h(z1*)≤0, where h(z)=z3+m1z2+m2z+m3.Lemma 2.2 (see [10]).
For the polynomial equation (2.12) (i)
ifm3≥0 and Δ=m12-3m2≤0, then all roots with positive real parts of (2.4) have the same sum as those of the polynomial equation (2.12) for all τ;(ii)
ifm3<0 or m3≥0, Δ=m12-3m2>0, and h(z1*)≤0, then all roots with positive real parts of (2.4) have the same sum as those of the polynomial equation (2.12) for τ∈[0,τ0).Summarizing the above analysis and noting that(2.14)m3=a32-b22≥0,
we have the following theorem.Theorem 2.3.
Assume that(i)
R0>1;(ii)
d-a+(2aT¯/Tmax)>0;(iii)
m3≥0 and m12-3m2≤0.Then the endemic equilibrium E¯ of (2.4) is absolutely stable, that is, E¯ is asymptotically stable for all values of the delay τ≥0.Now, we turn to the bifurcation analysis. We use the delayτ as bifurcation parameter. We view the solutions of (2.4) as functions of the bifurcation parameter τ. Let λ(τ)=η(τ)+iw(τ) be the eigenvalue of (2.10) such that for some initial value of the bifurcation parameter τ0 we have η(τ0)=0, and w(τ0)=w0 (without loss of generality we may assume w0>0). From (2.10) we have (2.15)τj=1w0arccos(-b1w04+(a2b1-a1b2)w02+a3b2b12w02+b22)+2jπw0,j=0,1,….Also, we can verify that the following transversal condition(2.16)dReλ(τ)dτ|τ=τ0>0
holds. By continuity, the real part of λ(τ) becomes positive when τ>τ0 and the steady state becomes unstable. Moreover, a Hopf bifurcation occurs when τ passes through the critical value τ0 (see [16]).To establish the Hopf bifurcation atτ=τ0, we need to show that dReλ(τ)/dτ|τ=τ0>0. Differentiating (2.4) from both sides with respect to τ, it follows that (2.17)(3λ2+2a1λ+a2)dλdτ=[-τe-λτ(b1λ+b2)+e-λτb1]dλdτ-λe-λτ(b1λ+b2).
This gives (2.18)(dλdτ)-1=3λ2+2a1λ+a2+τe-λτ(b1λ+b2)-e-λτb1-λe-λτ(b1λ+b2)=3λ2+2a1λ+a2-λe-λτ(b1λ+b2)+b1λ(b1λ+b2)-τλ=2λ3+a1λ2-a3-λ2(λ3+a1λ2+a2λ+a3)-b2λ2(b1λ+b2)-τλ.
Thus, (2.19)Sign{d(Reλ)dτ}|λ=iw0=Sign{Re(dλdτ)-1}λ=iw0=Sign{Re[2λ3+a1λ2-a3-λ2(λ3+a1λ2+a2λ+a3)]λ=iw0+Re[-b2λ2(b1λ+b2)]λ=iw0}=Sign{Re[-2w03i-a1w02-a3w02(-w03i-a1w02+a2w0i+a3)]+Re[-b2-w02(b1w0i+b2)]}=Sign{2w06+(a12-2a2)w04-a32w02[(a1w02-a3)2+(w03-a2w0)2]+b22w02[b22+b12w02]}=1w02Sign{2w06+(a12-2a2)w04+b22-a32b22+b12w02}.
Since (2.20)g(z)=2σ3+(a12-2a2)σ2+b22-a32,
thus, (2.21)dg(z)dz=6σ2+2(a12-2a2)σ.
The roots of (2.21) can be expressed as (2.22)σ1=0,σ2=-m13.
Noting that (2.23)m1=a12-2a2=(c+δ+ρ)2+(d-a+2aT¯Tmax⁡)2-2(δ+ρ)+ρβV¯>0,
hence, (2.24)dReλdτ|w=w0,τ=τ0=1w02Sign{2w06+(a12-2a2)w04+b22-a32b22+b12w02}>0.
The above analysis can be summarized into the following theorem.Theorem 2.4.
Suppose thatR0>1. If both  m3≥0 and  m12-3m2≥0 are satisfied, then the endemic equilibrium  E¯ of the delay model (1.2) is asymptotically stable when  τ<τ0 and unstable when  τ>τ0, where
(2.25)τ0=1w0arccos(-b1w04+(a2b1-a1b2)w02+a3b2b12w02+b22),
when τ=τ0, a Hopf bifurcation occurs, that is, a family of periodic solutions bifurcates from E¯ as τ passes through the critical value τ0.In this way, using time delay as a bifurcation parameter, Theorem2.4 indicates that the delay model could exhibit Hopf bifurcation at a certain value τ0 of the delay if the parameters satisfy conditions. They show that a time delay in the infected-to-viral cells transmission term can destabilize the system and periodic solutions can arise through Hopf bifurcation.


## 3. Direction and Stability of the Hopf Bifurcation

In this section, we will study the direction, stability, and the period of the bifurcating periodic solutions. The approach we used here is based on the normal form approach, the center manifold theory, and delay differential equation theory (see [16–19]). Throughout this section, we always assume that system (1.2) undergoes Hopf bifurcation at the positive equilibrium E¯=(S¯,I¯,V¯), for τ=τk, and then iw is corresponding purely imaginary roots of the characteristic equation at the positive equilibrium E¯=(S¯,I¯,V¯).Letting(3.1)T(t)=u1(t)+T¯,I(t)=u2(t)+I¯,V(t)=u3(t)+V¯,xi(t)=ui(τt),τ=τk+μ,
system (1.2) is transformed into an functional differential equation (FDE) in C=C([-1,0],R3) as (3.2)dxdt=Lμ(xt)+f(t,xt),
where x(t)=(x1(t),x2(t),x3(t))T∈R3, Lμ:C→R, and f:R×C→R are given, respectively, by (3.3)Lμ(ϕ)=(τk+μ)(a-d-βV¯-2aT¯Tmax⁡0-βT¯βV¯-δβT¯00-c)(ϕ1(0)ϕ2(0)ϕ3(0))+(τk+μ)(0000000p0)(ϕ1(-1)ϕ2(-1)ϕ3(-1)),(3.4)f(μ,ϕ)=(τk+μ)(-aTmax⁡ϕ12(0)-βϕ1(0)ϕ3(0)βϕ1(0)ϕ3(0)0).
By the Riesz representation theorem, there exists a function η(θ,μ) of bounded variation for θ∈[-1,0), such that (3.5)Lμϕ=∫-10dη(θ,0)ϕ(θ),forϕ∈C.
In fact, we can choose (3.6)η(θ,μ)=(τk+μ)(a-d-βV¯-2aT¯Tmax⁡0-βT¯βV¯-δβT¯00-c)(ϕ1(0)ϕ2(0)ϕ3(0))-(τk+μ)(0000000p0)(ϕ1(-1)ϕ2(-1)ϕ3(-1))δ(θ+1),
where δ is the Dirac delta function. For ϕ∈C([-1,0],R3), define (3.7)Aμ(ϕ)={dϕ(θ)dθ,θ∈[-1,0),dη(θ,s)ϕ(θ),θ=0R(θ,ϕ)={0,θ∈[-1,0),f(θ,ϕ),θ=0.
Then system (3.2) is equivalent to (3.8)ẋt=A(θ)xt+R(θ)xt,
where x(θ)=x(t+θ), for θ∈[-1,0).Forψ∈C1([0,1],(R3)*), define (3.9)A*ψ(s)={-dϕ(s)ds,s∈(0,1],∫-10dηT(t,0)ψ(-t),s=0,
and a bilinear inner product (3.10)〈ψ(s),ϕ(θ)〉=ψ(0)¯ϕ(0)-∫-10∫ξ=0θψ¯(ξ-θ)dη(θ)ϕ(ξ)dξ,
where η(θ)=η(θ,0). Then A(0) and A* are adjoint operators. By the discussion in Section 2, we know that ±iwτk are eigenvalues of A(0). Thus, they are also eigenvalues of A*. We first need to compute the eigenvectors of A(0) and A* corresponding to iwτk and -iwτk, respectively.Suppose thatq(θ)=∫-10(1,α,β)Teiw*τkθ is the eigenvector of A(0) corresponding to iw0τk, then (3.11)A(0)q(θ)=iwτkq(θ),A(0)q(0)=iwτkq(0).
It follows from the definition of A(0) and (3.3), (3.5), and (3.6) that (3.12)τk[(a-d-βV¯-2aT¯Tmax⁡0-βT¯βV¯-δβT¯00-c)+(0000000p0)]q(0)=iwτkq(0).
That is, (3.13)τk(iw-a+d+βV¯+2aT¯Tmax⁡0βT¯-βV¯iw+δ-βT¯0-qe-iwτkiw+c)(1αβ1)=(000).
Thus, we can easily obtain (3.14)q(0)=(1,α,β1)T,
where (3.15)α=-iw+cqe-iwτkiw-a+d+βV¯+(2aT¯/Tmax⁡)βT¯,β1=-iw-a+d+βV¯+(2aT¯/Tmax⁡)βT¯.Suppose thatq*(s)=D(1,α*,β1*)eiwτks, (3.16)τkD(1,α*,β1*)(-iw-a+d+βV¯+2aT¯Tmax⁡0βT¯-βV¯-iw+δ-βT¯0qe-iwτk-iw+c)=(000),
where (3.17)α*=-iw-a+d+βV¯+(2aT¯/Tmax⁡)βV¯,β1*=--iw+δqe-iwτk-iw-a+d+βV¯+(2aT¯/Tmax⁡)βV¯.In order to assume〈q*(θ),q(θ)〉=1, we need to determine the value of D. From (3.10), we have (3.18)〈q*(θ),q(θ)〉=D¯(1,α*,β1*)(1,α,β1)T-∫-10∫ξ=0θD¯(1,α*,β*)e-iwτkθdη(θ)(1,α,β)Teiwτkξdξ=D¯{1+αα*+β1β1*-∫-10∫ξ=0θD¯(1,α*,β*)θeiwτkθdη(θ)(1,α,β1)T}=D¯{1+αα*+β1β1*+τke-iwτk(1,α*,β1*)(0000000q0)(1αβ1)}=D¯{1+αα*+β1β1*+τke-iwτk(0,β1*q,0)(1αβ1)}=D¯{1+αα*+β1β1*+τke-iwτkαβ*q}.
Thus, we can choose D as (3.19)D¯=11+αα*+β1β1*+τke-iwτkαβ1*q.
In the remainder of this section, we use the same notations as in [19]; we first compute the coordinates to describe the center manifold C0 at μ=0. Let xt be the solution of (3.8) when μ=0. Define (3.20)z(t)=〈q*,xt〉,W(t,θ)=xt(θ)-2Re{z(t)q(θ)}.On the center manifoldC0, we have (3.21)W(t,θ)=W(z(t),z¯(t),θ),
where (3.22)W(z(t),z¯(t),θ)=W20(θ)z22+W11(θ)zz¯+W02(θ)z¯22+W30(θ)z36+⋯,z and z¯ are local coordinates for center manifold C0 in the direction of q and q*. Note that W is real if x is real. We only consider real solutions. For solution xt∈C0 of (3.8), since μ=0, we have (3.23)ż(t)=iw0τkz+q*¯(0)f(0,W(z,z¯,0)+2Re{zq(0)})≐iw0τkz+q*¯(0)f0(z,z¯).We rewrite this equation as(3.24)ż(t)=iw0τkz(t)+g(z,z¯)=g20z22+g11zz¯+g02z¯22+g21z2z¯2+⋯,
where (3.25)g(z,z¯)=q*¯(0)f0(z,z¯)=g20z22+g11zz¯+g02z¯22+g21z2z¯2+⋯.
It follows from (3.20) and (3.23) that (3.26)xt(θ)=w(t,θ)+2Re{z(t)q(θ)}=w(t,θ)+z(t)q(θ)+z¯(t)q(θ¯)=w20(θ)z22+w11zz¯+w02z¯22+(1,α,β1)Teiwτkθz+(1,α¯,β1¯)Te-iwτkθz¯+⋯.
It follows together with (3.4) that
(3.27)f0(z,z¯)=f(0,xt)=τk(-aTmax⁡x1t2(0)-βx1t(0)x3t(0)βx1t(0)x3t(0)0)=τk(-aTmax⁡[w20(1)(0)z22+w11(1)(0)zz¯+w02(1)(0)z¯22+z+z¯+⋯]×[w20(1)(0)z22+w11(1)(0)zz¯+w02(1)(0)z¯22+z+z¯+⋯]-β[w20(1)(0)z22+w11(1)(0)zz¯+w02(1)(0)z¯22+z+z¯+⋯]×[w20(3)(0)z22+w11(3)(0)zz¯+w02(3)(0)z¯22+β1z+β1¯z¯+⋯]β[w20(1)(0)z22+w11(1)(0)zz¯+w02(1)(0)z¯22+z+z¯+⋯]×[w20(3)(0)z22+w11(3)(0)zz¯+w02(3)(0)z¯22+β1z+β1¯z¯+⋯]0)=τk(2(-aTmax⁡-ββ1)z22+2(-aTmax⁡-ββ1¯)z¯22-2(-aTmax⁡+βRe(β1))zz¯-[-aTmax⁡(w20(1)(0)+2w11(1)(0)+w20(3)(0)+2w11(3)(0))-β(β1¯w20(1)(0)+2β1w11(1)(0)+w20(3)(0)+2w11(3)(0))aTmax⁡]z2z¯2+⋯2ββ1z22+ββ1¯z¯22+2βReβ1zz¯+[β1¯w20(1)(0)+2β1w11(1)(0)+w20(3)(0)+2w11(3)(0)]z2z¯2+⋯0).
Hence one can obtain
(3.28)g(z,z¯)=g20z22+g11zz¯+g02z¯22+g21z2z¯2+⋯=q*(0)f0(z,z¯)=D¯(1,α¯,β1¯)τk⋅(2(-aTmax⁡-ββ1)z22+2(-aTmax⁡-ββ1¯)z¯22-2(-aTmax⁡+βRe(β1))zz¯-[-aTmax⁡(w20(1)(0)+2w11(1)(0)+w20(3)(0)+2w11(3)(0))-β(β1¯w20(1)(0)+2β1w11(1)(0)+w20(3)(0)+2w11(3)(0))-aTmax⁡]z2z¯2+⋯2ββ1z22+ββ1¯z¯22+2βReβ1zz¯+[β1¯w20(1)(0)+2β1w11(1)(0)+w20(3)(0)+2w11(3)(0)]z2z¯2+⋯0)=D¯τk{2[-aTmax⁡-ββ1]z22-2[aTmax⁡+ββ1¯]z¯22-2[-aTmax⁡+βReβ1]zz¯-[-aTmax⁡(w20(1)(0)+2w11(1)(0)+w20(3)(0)+2w11(3)(0))-β(β1¯w20(1)(0)+2β1w11(1)(0)+w20(3)(0)+2w11(3)(0))-aTmax⁡-aTmax⁡]z2z¯2}+D¯α¯τk{2ββ1z22+ββ1¯z¯22+2βReβ1zz¯+[β1¯w20(1)(0)+2β1w11(1)(0)+w20(3)(0)+2w11(3)(0)]z2z¯2}+⋯.Comparing the coefficients with (3.25), we have (3.29)g20=-D¯τk2(-aTmax⁡+ββ1)+2D¯α¯τkββ1,g11=-2D¯τk(-aTmax⁡+βRe(β1))+2D¯α¯τkβRe(β1),g02=2D¯τk(-aTmax⁡+ββ1¯)+2Dα¯τkββ1¯,g21=-Dτk[-aTmax⁡(w20(1)(0)+2w11(1)(0)+w20(3)(0)+2w11(3)(0))-β(β1¯w20(1)(0)+2β1w11(1)(0)+w20(3)(0)+2w11(3)(0))-aTmax⁡]+D¯α¯β[β1¯w20(1)(0)+2β1w11(1)(0)+w20(3)(0)+2w11(3)(0)].
Since there are W20(θ) and W11(θ) in g21, we still need to compute them (3.30)ẇ=ẋt-żq-z¯̇q¯={Aw-2Re{q*(0)¯f0q(0)},θ∈[-1,0)Aw-2Re{q*(0)¯f0q(0)}+f0,θ=0≐Aw+H(z,z¯,θ),
where (3.31)H(z(t),z¯(t),θ)=H20(θ)z22+H11(θ)zz¯+H02(θ)z¯22+H30(θ)z36+⋯.
From (3.30) and (3.31), we have (3.32)A(0)w(t,θ)-ẇ=-H(z,z¯,θ)=-H20(θ)z22-H11zz¯-H02z¯22+⋯.
In view of (3.32), one can obtain (3.33)A(0)w(t,θ)=A(0)w20z22+A(0)w11zz¯+⋯,ẇ=wzż+wzz¯̇=w20(θ)z(iwτkz+g(z,z¯))+w11(θ){[iwτkz+g(z,z¯)]z¯+z[-iwτkz¯+g¯(z,z¯)]}+⋯=2iwτkw20(θ)z22+⋯.
It follows from (3.32) and (3.33) that (3.34)A(0)w-ẇ=[A(0)-2iwτkI]w20(θ)z22+A(0)w11(θ)zz¯+⋯=-H20z22-H11(θ)zz¯-H02(θ)z¯22-⋯.
Substituting the corresponding series into (3.30) and comparing the coefficients, we obtain (3.35)[A(0)-2iwτkI]w20(θ)=-H20(θ),A(0)w11(θ)=-H11(θ).From (3.30), we know that, for all θ∈[-1,0), (3.36)H(z,z¯,θ)=-2Re{q¯*f0(z,z¯)}=-q¯*f0(z,z¯)q(θ)-q*(θ)f¯0(z,z¯)q¯(θ)=-g(z,z¯)g(θ)-g¯(z,z¯)q¯(θ)=-(g02q(θ)+g¯02q¯(θ))z22-(g11q(θ)+g¯11q¯(θ))zz¯+⋯.
Comparing the coefficients with (3.31) gives that (3.37)H20(θ)=-(g02q(θ)+g¯02q¯(θ)),(3.38)H11(θ)=-(g11q(θ)+g¯11q¯(θ)).
From (3.35) and (3.37) and the definition of A(0), we have (3.39)ẇ20(θ)=2iw20+g20q(θ)+g¯20q¯(θ).
Note that q(θ)=q(0)eiwτkθ, hence (3.40)w20=ig20wτkq(0)eiwτkθ+ig¯203wτkq¯(0)eiwτkθ+e2iwτkθE1.
Similarly, from (3.35) and (3.38) and the definition of A(0), we have (3.41)w11=-ig11wτkq(θ)+ig¯11wτkq¯(θ)+E2,2iwτkw20(θ)=-2g20q(0)-23g20¯q¯(0)+2iwτkE1.
In what follows, we shall seek appropriate E1 and E2. From the definition of A and (3.35), we obtain (3.42)∫-10dη(θ)w20(θ)=2iwτkw20(0)-H20(0),(3.43)∫-10dη(θ)w11(θ)=-H11(0),
where g(θ)=g(0,θ). By (3.37) and (3.38), we have (3.44)H20(0)=-g20q(0)-g¯20q¯(0)-2τk(-aTmax⁡+ββ1ββ10),(3.45)H11(0)=-g11q(0)-g¯11q¯(0)-2τk(-aTmax⁡+ββ1βReβ10).
Substituting (3.42) into (3.44) and noticing that (3.46)(2iwτk-∫-10e2iwτkθdη(θ))E1=2τk(-aTmax⁡+ββ1ββ10)
which leads to (3.47)(2iw-a+d+βV¯+2aT¯Tmax⁡0βT¯-βV¯2iw+δ-βT¯0-qeiwτk2iw+c)(E11E21E31)=2(-aTmax⁡+ββ1βReβ10),
it follows that (3.48)E11=2A|-aTmax⁡+ββ10βT¯ββ12iw+δ-βT¯0-qeiwτk2iw+c|,E21=2A|2iw-a+d+βV¯+2aT¯Tmax⁡-aTmax⁡+ββ1βT¯-βV¯ββ1-βT¯002iw+c|,E31=2A|2iw-a+d+βV¯+2aT¯Tmax⁡0-aTmax⁡+ββ1-βV¯2iw+δββ10-qeiwτk0|,
where (3.49)A=|2iw-a+d+βV¯+2aT¯Tmax⁡0βT¯-βV¯2iw+δ-βT¯0-qeiwτk2iw+c|.
Similarly, substituting (3.43) into (3.45), we can get (3.50)(-a+d+βV¯+2aT¯Tmax⁡0βT¯-βV¯δ-βT¯0-qc)(E12E22E32)=2(-aTmax⁡+ββ1βReβ10),
and hence, (3.51)E12=2B|-aTmax⁡+ββ10βT¯βReβ12iw+δ-βT¯0-qeiwτk2iw+c|,E22=2B|2iw-a+d+βV¯+2aT¯Tmax⁡-aTmax⁡+ββ1βT¯-βV¯βReβ1-βT¯002iw+c|,E32=2B|2iw-a+d+βV¯+2aT¯Tmax⁡0-aTmax⁡+ββ1-βV¯2iw+δβReβ10-qeiwτk0|,
where (3.52)B=|-a+d+βV¯+2aT¯Tmax⁡0βT¯-βV¯δ-βT¯0-qc|.
It follows from (3.29) that g21 can be expressed explicitly. Thus, we can compute the following values: (3.53)c1(0)=i2w0τk(g11g20-2|g112|-g0222),σ2=-Re(c1(0))Re(λk′(τjk)),β2=2Re(c1(0)),T2=-Im(c1(0))+σ2Im(λk′(τjk))wτk,k=0,1,2,….By the result of Hassard et al. [19], we have the following.Theorem 3.1.
In (3.53), the sign of  σ2 determined the direction of Hopf bifurcation: if  σ2>0,(<0), then the Hopf bifurcation is supercritical (subcritical) and the bifurcating periodic solution exists for  τk>τjk(<τjk).  β2 determines the stability of the bifurcating periodic solution: the bifurcating periodic solution is stable (unstable) if  β2<0(>0), and  T2 determines the period of the bifurcating periodic solution: the period increases (decreases) if  T2>0(<0).


## 4. Simulation

In this section, we use numerical simulations to illustrate the theoretical results obtained in previous sections. As an example, we take the parameter values as follows:s=5, a=0.97, d=0.0002, Tmax=1200,δ=0.26, q=120, c=2.4, β=0.00024, τ=0.75, and ρ=0.01. By using the classical implicit format solving the delay differential equations and the method of steps for differential equations, we can solve the numerical solutions of (2.4) via the software package DEDiscover.Simulation of the model in this situation produces stable dynamics as is presented in Figure1. Plots (a)–(c) of Figure 1 show that uninfected cells, infected, cells and virus converge to their equilibrium with the parametric values. They show that the equilibrium E¯ under some conditions (see Theorem 2.3) is asymptotically stable.(a)–(c) show that uninfected cells, infected cells, and virus converge to their equilibrium with parametric values as stated in the text withτ=0.75. They show that the equilibrium is asymptotically stable.
(a)(b)(c)Next, we use a same set of parameter values as those in Figure1, but we vary the value of τ=2.01. Thus the conditions of Theorem 2.4 are satisfied. Then the system (2.4) has an asymptotically stable periodic orbit (See Figure 2). Plots (a)–(c) of Figure 2 are the oscillations of uninfected cells, infected cells, and virus if τ attend a certain level (see Theorem 2.4). Figure 2 shows that there is a periodic solution.(a)–(c) are the oscillations of uninfected cells, infected cells, and virus.
(a)(b)(c)We also find that the infection would always keep stability when the cure rateρ is larger. This can be analyzed from the expression of R0 and the conditions of Theorem 2.3. For example, we know that the oscillations of uninfected cells, infected cells and virus in Figure 3. And if we select ρ1=0.01, ρ2=0.61, and τ=2.01 and the other parameter values are same in Figure 1, then the infection would be stale (see Figure 3). Thus we can claim that the cure rate ρ is a very important parameter. The results show that if we improve the cure rate, we may control the disease.(a)–(c) show the uninfected cells, infected cells and virus withρ1=0.01 and ρ2=0.61. They show the cure rate is a important parameter.
(a)(b)(c)


## 5. Conclusion

An epidemic model of HIV infection ofCD4+ T cells with virus released period is studied. Mathematical analyses of the model equations with regard to dynamic behaviors of equilibria, Hopf bifurcation are analyzed. The basic reproduction number is obtained. In [5], if R0<1, the disease-free equilibrium is globally stable and the disease dies out. If R0>1, a unique endemic equilibrium exists, and it is globally asymptotically stable. In our model, we determine criteria for Hopf bifurcation using the time delay as the bifurcation parameter based on the differential-delay model. We show that positive equilibrium is locally asymptotically stable when time delay is suitably small, while a Hopf bifurcation can occur as the delay increases. Hopf bifurcation has helped us in finding the existence of a region of instability in the neighborhood of a nonzero endemic equilibrium where the population will survive undergoing regular fluctuations. We should discuss the length of the delay which impact on the stability of our model.



---

*Source: 649650-2011-05-05.xml*