Mathematics in Applied Sciences and Engineering https://ojs.lib.uwo.ca/mase
Volume 1, Number 3, September 2020, pp.236-248 https://doi.org/10.5206/mase/10828

GLOBAL DYNAMICS OF A TWO-STRAIN HIV INFECTION MODEL WITH

INTRACELLULAR DELAY

JIN XU

Abstract. In this paper, we formulate a mathematical model to describe the interaction of two

strains of HIV virus and the target cells within a host. The model is in the form of delay differential

equations with two discrete delays to account for the average time for replication for the two strains.

The model dynamic turns to be generically determined by two composite parameters R1 and R2, the

basic reproduction numbers for strain 1 and strain 2, respectively in the absence of the other strain,

in the sense that except for the critical case R1 = R2 > 1, the solutions are proved to converge to the

corresponding equilibrium globally. The method used is Lyapunov functionals.

1. Introduction

It has been realized that mathematical modelling can provide valuable insight into HIV-1 pathogene-

sis. These mathematical models are formulated by using differential equations to explore the mechanisms

and dynamical behaviors of the viral infection process [3, 8, 17, 18, 19]. Such understanding may offer

guidance for developing efficient anti-viral drug therapies [14, 15, 10].

Most existing mathematical models for HIV virus dynamics are by systems of ordinary differential

equations. A standard and classic differential equation model for HIV infection is the following system

of ODEs [16, 14, 18]:



Ṫ = λ−dT −kTV,

Ṫ∗ = kTV −µT∗,

V̇ = pT∗ − cV,

(1.1)

where T(t),T∗(t) and V (t) are the population sizes of uninfected target cells, infected cells and the

free virus particles, respectively, at time t. The assumption is that uninfected cells are generated at

a constant rate, λ, and die at a rate d. Free virus particles infect uninfected target cells at a rate

proportional to the product of their abundances, kTV . The rate constant, k, describes the efficacy of

this process. Infected cells produce free virus particles at a rate proportional to their abundance, pT∗.

Infected cells die at a rate µT∗ either due to the natural death or the action of the virus and free virus

particles are removed from the system at rate cV by the immune system or natural decay. Therefore,

the average life-time of an infected cell, a free virus particle and an uninfected cell are 1/µ, 1/c and

1/d respectively. The model well predicts the primary phase of HIV infection, showing that during

the first weeks of infection there is a peak in viral load with a subsequent decline to a relatively stable

steady-state.

Received by the editors 29 June 2020; accepted 24 August 2020; online first 6 September, 2020.

2000 Mathematics Subject Classification. 34K20, 34K25, 92B05, 92D25.

Key words and phrases. HIV, virus dynamics, intracellular delay, global stability, Lyapunov functional.

236



GLOBAL DYNAMICS OF A TWO-STRAIN HIV INFECTION MODEL WITH INTRACELLULAR DELAY 237

Now, we assume there is an another subtype of virus in within a host which competes with the

original virus for host cell resource. Assuming that super infection is negligible, an ordinary differential

equations can be formulate along the line of (1.1) to describe the interaction between the two subtype

viruses and host cells, as given below.


Ṫ = λ−dT −k1TV1 −k2TV2,

Ṫ1 = k1TV1 −µ1T1,

Ṫ2 = k2TV2 −µ2T2,

V̇1 = p1T1 − c1V1,

V̇2 = p2T2 − c2V2,

(1.2)

where T1(t) denotes the population size of cells productively infected by strain-1 virus, whereas T2(t)

denotes the population size of cells productively infected by strain-2 virus at time t; V1(t) and V2(t)

represent the respective population sizes of subtype-1 and subtype-2 viruses; k1 and k2 represent the rate

constants at which uninfected target cells are infected by subtype-1 and subtype-2 viruses, respectively.

The two subtypes of infected cells are assumed to have two different death rate µ1 and µ2. Once

uninfected target cells are infected by subtype-1 (subtype-2) viruses, new subtype-1 (subtype-2) virus

particles are produced with constant rate p1 (p2). The new subtypes of virus have the respective

clearance rate c1 and c2. All the parameters of the model are assumed to be positive. Here we omit the

super-infection in host cells.

However, in reality, there is a lag between the time target cells are contacted by virus particles and

the time the contacted cells become actively affected meaning that the virions have enter cells and

started producing new virions [23]. This can be explained by the initial phase of the virus life cycle,

which include all stages from viral attachment until the time that the host cell contains the infectious

viral particles in its cytoplasm. To account for this lag, models that include time delays have been

developed and investigated [8, 15, 23]. One distinct feature of delay differential equation models is

that a delay typically destabilizes an stable equilibrium and causes sustained oscillation through Hopf

bifurctions. By rigorously establishing the global dynamics of the two-strain competitive viral model

with intracellular delays, we show that no sustained oscillations are possible in our model.

To incorporate the intracellular phase of the virus life-cycle, we assume that subtype-1 virus and

subtype-2 virus production occur in average, τ1 and τ2 time units later, after the respective virus enter

the host cells. The recruitment of subtype-1 virus producing cells at time t is given by the number of

cells that were newly infected by strain-1 at time t−τ1 and are still alive at time t. In the same way, the
recruitment of subtype-2 virus producing cells at time t is given by the number of cells that were newly

infected by strain-2 at time t− τ2 and are still alive at time t. If we assume two constant death rates
s1 and s2 for infected but not yet virus-producing cells for subtype-1 and subtype-2, the probability of

subtype-1 surviving the time period from t−τ1 to t is e−s1τ1 , the probability of subtype-2 surviving the
time period from t−τ2 to t is e−s2τ2 . The transfer diagram for the transmission of viral infection under
such a scenario is shown in Figure 1. Thus the following delay differential equations model is proposed:



Ṫ = λ−dT(t) −k1T(t)V1(t) −k2T(t)V2(t),

Ṫ1 = k1T(t− τ1)V1(t− τ1)e−s1τ1 −µ1T1(t),

Ṫ2 = k2T(t− τ2)V2(t− τ2)e−s2τ2 −µ2T2(t),

V̇1 = p1T1(t) − c1V1(t),

V̇2 = p2T2(t) − c2V2(t),

(1.3)



238 J. XU

1 −
Infected cells that can

Free subtype-1 virus Infected cells by subtype-1 produce subtype-1

Uninfected cells

d

Free subtype-2 virus Infected cells by subtype-2 Infected cells that can

produce subtype-21 −

Figure 1. Transfer diagram for model (1.3)

Delays have been incorporated into virus dynamics models in [8, 23, 12], but only for single strain

models. Here we consider two strains. Many previous in-host models also considered the effects of

anti-viral drug therapies such as HAART [15, 1, 21], but only local stability were analysed in these

works. We note that by renaming the coeffiicients due to the effect of reverse transcriptase inhibitors

and protease inhibitors, the model in [15, 1, 21] can be transformed into the form of (1.3). Our results

on the global dynamics of model (1.3) can apply to these models with anti-viral therapies, and hence

can rule out the exitence of periodic solutions. This shows novelty of this work and should benefit other

researchers working on similar models.

In the present section we analyse model (1.3) including intracellular delays. We establish global

asymptotic stability of the infected-free, and single-infected by constructing Lyapunov functionals. To

this end, we first establishes the well-posedness of (1.3) in section 3.2. Then we discuss the existence

of equilibria in the feasible region and derive the basic reproductive number R0. It turns out that R1
is a decreasing function of the delay τ1 and R2 is a decreasing function of the delay τ2. These imply

that ignoring the intracellular delays will overestimate the basic reproduction number. We show that

the basic reproductive number R0 generically determines the global dynamics of model (1.3). More

specifically, if R0 ≤ 1, the infection-free equilibrium E0 is globally asymptotically stable, and two
subtype viruses will be cleared; if R0 > 1 and R1 6= R2, the single-infected equilibrium arising from
the greater basic reproduction number is globally asymptotically stable. The proof utilizes a global

Lyapunov functional that is motivated by the work in [11, 12]. The global stability of single-infected

equlibira rule out any possibility of sustained oscillations. In addition, numerical simulations are also

conducted to demonstrate global dynamics of system (1.3).

2. Well-posedness

In the same way as in the previous section, the system (1.3) is biologically acceptable in the sense that

no population goes negative. We expect that starting from non-negative initial values, the corresponding

solution remains non-negative. To proceed, we follow the convention to denote by C1 = C([−τ1, 0],R)
and C2 = C([−τ2, 0],R) the Banach spaces of continuous functions mapping the interval [−τi, 0] into
R, i = 1, 2, with norm ‖Φi‖ = sup−τi≤θ≤0 |Φi(θ)| for Φi ∈ Ci. Let τ = max{τ1,τ2}, denote by



GLOBAL DYNAMICS OF A TWO-STRAIN HIV INFECTION MODEL WITH INTRACELLULAR DELAY 239

C = C([−τ, 0],R) the Banach space of continuous functions mapping the interval [−τ, 0] into R, with
norm ‖Φ‖ = sup−τ≤θ≤0 |Φ(θ)| for θ ∈ C. The nonnegative cone of C,C1 and C2 are defined as
C+ = C([−τ, 0],R+),C+1 = C([−τ1, 0],R+) and C

+
2 = C([−τ2, 0],R+). The initial conditions for

system (1.3) are chosen at t = 0 as ϕ ∈ C+ × R+ × R+ ×C+1 ×C
+
2 . The well-posedness for our delay

differential equation model (1.3) is established by the following theorem.

Theorem 2.1. Under the above initial conditions, all solutions of system (1.3) are positive and ulti-

mately bounded in C ×R×R×C1 ×C2

Proof. First, we prove that T(t) is positive for all t ≥ 0. Assuming the opposite, let t1 > 0 be the first
time such that T(t1) = 0, which means T(t) > 0 as t ∈ [0, t1). Since

Ṫ = λ−dT(t) −k1T(t)V1(t) −k2T(t)V2(t),

we get Ṫ(t1) = λ > 0,and hence T(t) < 0 for t ∈ (t1 − �,t1) where � > 0 is sufficiently small. This
contradicts T(t) > 0 for t ∈ [0, t1). It follows that T(t) > 0 for t > 0. Next, we show V1(t) ≥ 0 for all
t ≥ 0. Assume the opposite and let t2 > 0 be the first time such that V1(t2) = 0. Since

V̇1(t) = p1T1(t) − c1V1(t),

we have V̇1(t2) = p1T1(t2). On the other hand, solving T1(t) by the second equation of (1.3) gives

T1(t2) = (T1(0) +

∫ t2
0

k1T(θ − τ1)V1(θ − τ1)e−s1τ1eµ1θdθ)e−µ1t2 > 0

Hence V̇1(t2) = p1T1(t2) > 0 implying V1(t) is positive for all t ≥ 0.
The positiveness of T(t) and V1(t) and the following formula

T1(t) = (T1(0) +

∫ t
0

k1T(θ − τ1)V1(θ − τ1)e−s1τ1eµ1θdθ)e−µ1t > 0.

in turn leads to the positiveness of T1(t) for all t ≥ 0. Similarly, we can show that V2(t) and T2(t) are
positive for t ≥ 0 under positive initial conditions.

From the first equation of (1.3), we obtain ˙T(t) ≤ λ−dT(t). Hence limsupt→∞T(t) ≤ λd . Adding
the first three equations of (1.3), it follows

(T(t) + T1(t + τ1) + T2(t + τ2))
′ = λ−dT(t) −µ1T1(t + τ1) −µ2T2(t + τ2)

+ k1T(t)V1(t)(e
−s1τ1 − 1) + k2T(t)V2(t)(e−s2τ2 − 1)

≤ λ− r̃(T(t) + T1(t + τ1) + T2(t + τ2))

where r̃ = min{d,µ1,µ2}. Thus, limsupt→∞(T(t) +T1(t+τ1) +T2(t+τ2)) ≤ λr̃ . For any � > 0,∃t
∗ > 0,

such that T(t) + T1(t + τ1) + T2(t + τ2) ≤ λr̃ + � for all t ≥ t
∗. Thus, T(t),T1(t) and T2(t) are all

ultimately bounded by λ
r̃
. The fourth equation of (1.3) implies

V̇1 = p1T1(t) − c1V1(t) ≤ p1(
λ

r̃
+ �) − c1V1(t), t ≥ t∗

This implies limsupt→∞V1 ≤ p1c1 (
λ
r̃

+�). Since � > 0 is arbitrary, we attain limsupt→∞V1(t) ≤ p1λc1r̃ . Sim-
ilarly, we can obtain limsupt→∞V2(t) ≤ p2λc2r̃ . Therefore, T(t),T1(t),T2(t),V1(t) and V2(t) are ultimately
bounded in C ×R×R×C1 ×C2.



240 J. XU

3. Equilibria and basic reproduction numbers

In system (1.3), without infection (T1,T2,V1,V2) = (0, 0, 0, 0), uninfected target cells stabilizes at

the equilibrium T = λ
d

. The basic reproductive number R1 for in-host models [17, 12, 16] measures

the average number virus-producing target cells produced by a single subtype-1 virus-producing target

cell during its entire infectious period in an entirely uninfected target-cell population. As illustrated in

Figure 2, the basic reproduction number R1 for strain-1 is given by

R1 =
p1
µ1
·
k1e
−s1τ1

c1
·
λ

d
. (3.1)

Similarly, the basic reproduction number R2 for strain-2 which is the average number virus-producing

target cells produced by a single subtype-2 virus producing target cell during its entire infectious period

in an entirely uninfected target-cell population is obtained by

R2 =
p2
µ2
·
k2e
−s2τ2

c2
·
λ

d
. (3.2)

When no intracellular delay is considered, τ1 = τ2 = 0, our R1 and R2 reduce to the respective basic

reproduction number for our previous model (3.1) (i.e. (2.21)). If s > 0, R1 and R2 is the decreasing

functions of the delay τ1 and τ2. It shows that the intracellular delays decrease R1 and R2 if cells die

during the delay periods. Thus, ignoring the intracellular delay in a viral model will overestimate the

basic reproduction number.

From our system (1.3) and our result (3.1) (3.2), we define the system basic reproduction number

R0 = max{R1,R2} . (3.3)

1 −
Infected cells that can

Free subtype-1 virus Infected cells by subtype-1 produce subtype-1

Uninfected cells

d

Free subtype-2 virus Infected cells by subtype-2 Infected cells that can

produce subtype-21 −

Burst size Basic reproduction number= ⋅ ⋅
Figure 2. An illustration of the basic reproduction number of model(1.3)

Model system (1.3) always has the infection-free equilibrium E0 =
(
λ
d
, 0, 0, 0, 0

)
. There are two

possible single-infection equilibria E1 =
(
T̄, T̄1, 0, V̄1, 0

)
and E2 =

(
T̃, 0, T̃2, 0, Ṽ2

)
, where

T̄ =
λ

d

1

R1
, T̄1 =

dc1
k1p1

(R1 − 1), V̄1 =
d

k1
(R1 − 1). (3.4)

and

T̃ =
λ

d

1

R2
, T̃2 =

dc2
k2p2

(R2 − 1), Ṽ2 =
d

k2
(R2 − 1). (3.5)

It turns out that the values of R1 and R2 determine the existence of the single-infection equilibria: E1
exists if and only if R1 > 1 and E2 exists if and only if R2 > 1. Obviously, E1 and E2 are biologically

meaningful under the conditions.



GLOBAL DYNAMICS OF A TWO-STRAIN HIV INFECTION MODEL WITH INTRACELLULAR DELAY 241

It is also possible for our model (1.3) to obtain the double-infection equilibrium which means

a equilibrium with all components being positive. Denote such a possible equilibrium by E3 =

(T∗,T∗1 ,T
∗
2 ,V

∗
1 ,V

∗
2 ), then calculation shows that the components in E3 must satisfy



T∗ =
µ1c1e

s1τ1

k1p1
(i.e.

d

λR1
) =

µ2c2e
s2τ2

k2p2
(i.e.

d

λR2
),

T∗1 =
c1V

∗
1

p1
,

T∗2 =
c2V

∗
2

p2
,

d(R1 − 1) = k1V ∗1 + k2V
∗
2 ,

d(R2 − 1) = k1V ∗1 + k2V
∗
2 .

(3.6)

By the last two equation in (3.6), it is clear that E3 exists if and only if

R1 = R2 > 1. (3.7)

If (3.7) holds, there are actually infinitely many co-existence equilibria.

Summarizing the above results, we have the following conclusion. When R0 ≤ 1, E0 is the only
equilibrium; when R1 > 1,R2 ≤ 1, there are E0 and E1; when R2 > 1,R1 ≤ 1, there are E0 and E2;
when R1 > 1 and R2 > 1, in addition to E0,E1 and E2, there are infinitely many co-exitence equilibria

if R1 = R2 > 1. Considering the fact that there are ten model parameters in R1 and R2, the identity

R1 = R2 is unlikely to hold in practice (or infeasible), and hence, E3 will not be considered here in this

thesis.

4. Global Stability of Equilibria

In this section we study the global stability of equilibria by using the Lyapunov functionals.

We apply Lyapunov functionals similar to those recently used by [11, 6, 20]. A useful function is

used to construct our Lyapunov fuctionals:

g(x) = x− ln(x) − 1.

This function attains the global minimum at x = 1, g(1) = 0, and remains positive for all other

postitive values of x. Our Lyapunov functionals take advantage of these properties of g(x). In the

following theorems we show that the equilibria exhibit global stability under some threshold conditions.

Theorem 4.1. If R0 ≤ 1, the infection free-equilibrium E0 is globally asymptotically stable.

Proof. Let T0 =
λ
d

and consider the Lyapunov functional

V (T,T1,T2,V1,V2) = T0g(T(t)/T0) + e
s1τ1T1(t) + e

s2τ2T2(t)

+
µ1
p1
es1τ1V1(t) +

µ2
p2
es2τ2V2(t)

+k1

∫ 0
−τ1

T(t + θ)V1(t + θ) dθ + k2

∫ 0
−τ2

T(t + θ)V2(t + θ) dθ.

Obviously, V (T,T1,T2,V1,V2) is non-negative in the positive cone C
+×R+×R+×C+1 ×C

+
2 and attains

zero at E0. We will show that the derivative of V along the trajectories of our model (1.3) is negatively

defininte. Differentiation gives



242 J. XU

V̇ = ˙T(t) −
T0
T(t)

˙T(t) + es1τ1 ˙T1(t) + e
s2τ2 ˙T2(t) +

µ1
p1
es1τ1 ˙V1(t) +

µ2
p2
es2τ2 ˙V2(t)

+k1T(t)V1(t) −k1T(t− τ1)V1(t− τ1) + k2T(t)V2(t) −k2T(t− τ2)V2(t− τ2)

= λ−dT(t) −k1T(t)V1(t) −k2T(t)V2(t) −
T0
T(t)

(λ−dT −k1T(t)V1(t) −k2T(t)V2(t))

+es1τ1
(
k1T(t− τ1)V1(t− τ1)e−s1τ1 −µ1T1

)
+ es2τ2 (k2T(t− τ2)V2(t− τ2)e−s2τ2 −µ2T2)

+
µ1
p1
es1τ1

(
p1T1(t) − c1V1(t)) +

µ2
p2
es2τ2 (p2T2(t) − c2V2(t)

)
+k1T(t)V1(t) −k1T(t− τ1)V1(t− τ1) + k2T(t)V2(t) −k2T(t− τ2)V2(t− τ2)

After cancelling terms, using T0 =
λ
d

and rearranging terms, we get

V̇ = λ−dT(t) −
T0
T(t)

λ + dT0

+

(
k1T0 −

c1µ1
p1

es1τ1
)
V1(t) +

(
k2T0 −

c2µ2
p2

es2τ2
)
V2(t)

= λ

(
2 −

T(t)

T0
−

T0
T(t)

)
+
c1µ1
p1

es1τ1
(
k1p1λ

µ1c1d
es1τ1 − 1

)
V1(t) +

c2µ2
p2

es2τ2
(
k2p2λ

µ2c2d
es2τ2 − 1

)
V2(t)

= λ

(
2 −

T(t)

T0
−

T0
T(t)

)
+
µ1c1
p1

es1τ1 (R1 − 1)V1(t) +
µ2c2
p2

es2τ2 (R2 − 1)V2(t).

Since the arithmetic mean is greater than or equal to the geometric mean, if R0 = max{R1,R2}≤ 1,
each of the three terms on the right hand side is non-positive. Hence V̇ (T,T1,T2,V1,V2) ≤ 0, and V̇ = 0
if and only if (T,T1,T2,V1,V2) =

(
λ
d
, 0, 0, 0, 0

)
= E0 Therefore, the globally asymptotical stability of

E0 follows from the Lyaunov-LaSalle invariance principle by [7].

When R0 > 1, then E0 becomes unstable and at least one of the E1 and E2 exists. We now investigate

the global stability of these two possible single-strain equilibria.

Theorem 4.2. Assume that E1 exists (i.e. R1 > 1), if R2 < R1, then, E1 is globally asymptotically

stable.

Proof. Define a Lyapunov functional V : C ×R×R×C1 ×C2 → R by

V (T,T1,T2,V1,V2) = T̄g(
T(t)

T̄
) + T̄1e

s1τ1g(
T1(t)

T̄1
) + es2τ2T2(t) +

µ1
p1
V̄1e

s1τ1g(
V1(t)

V̄1
)

+
µ2
p2
es2τ2V2(t) + k1T̄V̄1

∫ 0
−τ1

g(
T(t + θ)V1(t + θ)

T̄V̄1
) dθ

+k2

∫ 0
−τ2

T(t + θ)V2(t + θ) dθ.

By the properties of g(x), the Lyapunov functional V (T,T1,T2,V1,V2) is non-negative in the positive

cone C+ × R+ × R+ ×C+1 ×C
+
2 and attains zero at E1. In order to show V̇ is negatively definite, we

differentiate V (T,T1,T2,V1,V2) along the trajectories of (1.3) to get



GLOBAL DYNAMICS OF A TWO-STRAIN HIV INFECTION MODEL WITH INTRACELLULAR DELAY 243

V̇ = Ṫ(t) +
T̄

T(t)
Ṫ(t) + es1τ1Ṫ1(t) −es1τ1

T̄1
T1(t)

Ṫ1(t) + e
s2τ2Ṫ2(t)

+
µ1
p1
es1τ1V̇1(t) −

µ1
p1
es1τ1

V̄1
V1(t)

V̇1(t) +
µ2
p2
es2τ2V̇2(t)

+k1T̄V̄1
d

dt

∫ 0
−τ1

g(
T(t + θ)V1(t + θ)

T̄V̄1
) dθ

+k2T(t)V2(t) −k2T(t− τ2)V2(t− τ2). (4.1)

Note that

k1T̄V̄1
d

dt

∫ 0
−τ1

g(
T(t + θ)V1(t + θ)

T̄V̄1
) dθ

= k1T̄V̄1

∫ 0
−τ1

d

dt
g(
T(t + θ)V1(t + θ)

T̄V̄1
) dθ

= k1T̄V̄1

∫ 0
−τ1

d

dθ
g(
T(t + θ)V1(t + θ)

T̄V̄1
) dθ

= k1T̄V̄1

(
g(
T(t)V1(t)

T̄V̄1
) −g(

T(t− τ1)V1(t− τ1)
T̄V̄1

)

)
= k1T̄V̄1

(
T(t)V1(t)

T̄V̄1
− ln

T(t)V1(t)

T̄V̄1
−
T(t− τ1)V1(t− τ1)

T̄V̄1
+ ln

T(t− τ1)V1(t− τ1))
T̄V̄1

)
= k1T(t)V1(t) −k1T̄V̄1lnT(t)V1(t) + k1T̄V̄1lnT̄V̄1 −k1T(t− τ1)V1(t− τ1)

+k1T̄V̄1lnT(t− τ1)V1(t− τ1) −k1T̄V̄1lnT̄V̄1
= k1T(t)V1(t) −k1T(t− τ1)V1(t− τ1)

+k1T̄V̄1lnT(t− τ1)V1(t− τ1) −k1T̄V̄1lnT(t)V1(t) (4.2)

Plugging (4.2) and system of (1.3) into equation (4.1), we obtain

V̇ =λ−dT(t) −k1T(t)V1(t) −k2T(t)V2(t) −
T̄

T(t)
λ + dT̄ + k1T̄V1(t) + k2T̄V2(t)

+ k1T(t− τ1)V1(t− τ1) −µ1es1τ1T1(t) − T̄1
k1T(t− τ1)V1(t− τ1)

T1(t)
+ es1τ1µ1T̄1

+ k2T(t− τ2)V2(t− τ2) −µ2es2τ2T2(t) + µ1es1τ1T1(t) −
µ1c1
p1

es1τ1V1(t)

−µ1es1τ1V̄1
T1(t)

V1(t)
+
µ1c1
p1

es1τ1V̄1

+ µ2e
s2τ2T2(t) −

µ2c2
p2

es2τ2V2(t)

+ k1T(t)V1(t) −k1T(t− τ1)V1(t− τ1)
−k1T̄V̄1lnT(t)V1(t) + k1T̄V̄1lnT(t− τ1)V1(t− τ1
+ k2T(t)V2(t) −k2T(t− τ2)V2(t− τ2).

(4.3)

The components of E1 are related by the equilibrium equation, i.e.,



244 J. XU




λ = dT̄ + k1T̄V̄1

k1T̄V̄1 = µ1T̄1e
s1τ1

p1T̄1 = c1V̄1

k1T̄ =
µ1T̄1e

s1τ1

V̄1
=
µ1T̄1e

s1τ1c1
p1T̄1

=
µ1c1
p1

es1τ1.

(4.4)

Making use of these, we can rearrange and simplify the equation (4.3) as

V̇ = dT̄

(
2 −

Tt

T̄
−

T̄

T(t)

)
−
k1T̄

2V̄1
T(t)

+k1T̄V̄1 + k2T̄V2(t) −k1T̄1
T(t− τ1)V1(t− τ1)

T1(t)

+k1T̄V̄1 −
k1T̄V̄1
T̄1

V̄1
T1(t)

V1(t)
+ k1T̄V̄1 −

µ2c2
p2

es2τ2V2(t)

−k1T̄V̄1lnT(t)V1(t) + k1T̄V̄1lnT(t− τ1)V1(t− τ1)

= dT̄

(
2 −

Tt

T̄
−

T̄

T(t)

)
−k1T̄V̄1

(
g(
T̄1T(t− τ1)V1(t− τ1)

T̄V̄1T1(t)
) − ln

T̄1T(t− τ1)V1(t− τ1)
T̄V̄1T1(t)

)
−k1T̄V̄1

(
g(

T̄

T(t)
) − ln

T̄

T(t)
) −k1T̄V̄1(g(

V̄1T1(t)

T̄1V1(t)
) − ln

V̄1T1(t)

T̄1V1(t)

)
−k1T̄V̄1 [lnT(t)V1(t) − lnT(t− τ1)V1(t− τ1)]

+

(
k2T̄ −

µ2c2
p2

es2τ2
)
V2(t)

= dT̄

(
2 −

Tt

T̄
−

T̄

T(t)

)
−k1T̄V̄1g(

T̄1T(t− τ1)V1(t− τ1)
T̄V̄1T1(t)

)

−k1T̄V̄1g(
T̄

T(t)
) −k1T̄V̄1g(

V̄1T1(t)

T̄1V1(t)
)

+
k2λ

d

(
1

R1
−

1

R2

)
V2(t).

Therefore, by our assumptions, V̇ ≤ 0 with equality holding only at E1. From the Lyapunov- LaSalle
inveriance principle [7], the equilibrium E1 is globally asymptotically stable. The proof is completed.

Parallel to Theorem 4.2, we have the following theorem for E2

Theorem 4.3. Assume that E2 exists (i.e. R2 > 1), if R1 < R2, then E2 is globally asymptotically

stable.

Proof. The proof of this theorem is symmetric to that of Theorem 4.2 by considering the following

Lyapunov functional:V : C ×R×R×C1 ×C2 → R



GLOBAL DYNAMICS OF A TWO-STRAIN HIV INFECTION MODEL WITH INTRACELLULAR DELAY 245

V (T,T1,T2,V1,V2) = T̃g(
T(t)

T̃
) + es1τ1T1(t) + T̃2e

s2τ2g(
T2(t)

T̃2
) +

µ1
p1
es1τ1V1(t)

+
µ2
p2
Ṽ2e

s2τ2g(
V2(t)

Ṽ2
) + k1

∫ 0
−τ1

T(t + θ)V1(t + θ) dθ

k2T̃Ṽ2

∫ 0
−τ2

g(
T(t + θ)V2(t + θ)

T̃Ṽ2
) dθ

We omit the details of the proof.

5. Numerical Simulations

In this section, we present some numeric simulations for the DDE model (3.2) to confirm and illustrate

the theoretic results obtained in Section 3.4, which is not significantly different from those for the ODE

model (2.2), except that some plottings are in logarithmic function for better and clearer displays.

First, we chose the following values for the model parameters: λ = 6,d = 1,k1 = 2,p1 = 1,c1 =

3,µ1 = 10,s1 = 2,τ1 = 0.1,k2 = 3,p2 = 2,c2 = 2.5,µ2 = 15,s2 = 1.5,τ2 = 0.15. This give the values

two individual basic reproduction numbers R1 = 0.327 and R2 = 0.767. Three sets of initial values are

used: (I) T(0) = 80,T1(0) = 50,T2(0) = 40,V1(0) = 45,V2(0) = 35; (II) T(0) = 60,T1(0) = 70,T2(0) =

50,V1(0) = 30,V2(0) = 20; (III) T(0) = 50,T1(0) = 60,T2(0) = 30,V1(0) = 20,V2(0) = 45. We used

a base 10 logarithmic scale for target cells population. The corresponding solutions are presented in

Figure 3.

Second, we chose the following values for the model parameters: λ = 6,d = 1,k1 = 5,p1 = 6,c1 =

4,µ1 = 3,s1 = 2,τ1 = 0.1,k2 = 1,p2 = 4,c2 = 3,µ2 = 4,s2 = 1.5,τ2 = 0.15. This give the values

two individual basic reproduction numbers R1 = 12.28 and R2 = 1.597. Three sets of initial values are

used: (I) T(0) = 80,T1(0) = 50,T2(0) = 40,V1(0) = 45,V2(0) = 35; (II) T(0) = 60,T1(0) = 70,T2(0) =

50,V1(0) = 30,V2(0) = 20; (III) T(0) = 50,T1(0) = 60,T2(0) = 30,V1(0) = 20,V2(0) = 45. A base

10 logarithmic scale for target cells population, subtype-1 infected cells and subtype-1 virus cells was

employed in our figures. The corresponding solutions are presented in Figure 4.

Third, we chose the following values for the model parameters: λ = 6,d = 1,k1 = 4,p1 = 8,c1 =

8,µ1 = 5,s1 = 2,τ1 = 0.1,k2 = 3,p2 = 10,c2 = 5,µ2 = 4, ,s2 = 1.5,τ2 = 0.15. This give the values

two individual basic reproduction numbers R1 = 3.93 and R2 = 7.19. Three sets of initial values are

used: (I) T(0) = 80,T1(0) = 50,T2(0) = 40,V1(0) = 45,V2(0) = 35; (II) T(0) = 60,T1(0) = 70,T2(0) =

50,V1(0) = 30,V2(0) = 20; (III) T(0) = 50,T1(0) = 60,T2(0) = 30,V1(0) = 20,V2(0) = 45. A base

10 logarithmic scale for target cells population, subtype-2 infected cells and subtype-2 virus cells was

employed in our figures. The corresponding solutions are presented in Figure 5.



246 J. XU

0 1 2 3 4 5 6 7 8 9 10
10

−5

10
0

10
5

Target Uninfected Cells

time

lo
g1

0 
ta

rg
et

 c
el

ls

0 2 4 6 8 10
0

200

400
Infected Cells by Subtype Virus−1

time

in
fe

ct
ed

 c
el

ls

0 2 4 6 8 10
0

500

1000
Infected Cells by Subtype Virus−2

time

in
fe

ct
ed

 c
el

ls
0 2 4 6 8 10

0

100

200

300

time

vi
ru

s−
1 

po
pu

la
tio

n

Subtype Virus−1 Cells

0 2 4 6 8 10
0

200

400

time

vi
ru

s−
2 

po
pu

la
tio

n

Subtype Virus−2 Cells

Figure 3. R1 < 1 and R2 < 1: viruses of both strains all die out

0 1 2 3 4 5 6 7 8 9 10
10

−5

10
0

10
5

Target Uninfected Cells

time

lo
g1

0 
ta

rg
et

 c
el

ls

0 5 10
10

0

10
2

10
4

Infected Cells by Subtype Virus−1

time

lo
g1

0 
in

fe
ct

ed
 c

el
ls

0 5 10
0

200

400
Infected Cells by Subtype Virus−2

time

in
fe

ct
ed

 c
el

ls

0 5 10
10

0

10
2

10
4

timelo
g1

0 
vi

ru
s−

1 
po

pu
la

tio
n

Subtype Virus−1 Cells

0 5 10
0

100

200

time

vi
ru

s−
2 

po
pu

la
tio

n Subtype Virus−2 Cells

Figure 4. R1 > 1 and R2 < R1: subtype-1 wins the competition



GLOBAL DYNAMICS OF A TWO-STRAIN HIV INFECTION MODEL WITH INTRACELLULAR DELAY 247

0 1 2 3 4 5 6 7 8 9 10
10

−5

10
0

10
5

Target Uninfected Cells

time

lo
g1

0 
ta

rg
et

 c
el

ls

0 5 10
0

500

1000
Infected Cells by Subtype Virus−1

time

in
fe

ct
ed

 c
el

ls

0 5 10
10

−5

10
0

10
5

Infected Cells by Subtype Virus−2

time

lo
g1

0 
in

fe
ct

ed
 c

el
ls

0 5 10
0

200

400

600

time

vi
ru

s−
1 

po
pu

la
tio

n

Subtype Virus−1 Cells

0 5 10
10

0

10
2

10
4

timelo
g1

0 
vi

ru
s−

2 
po

pu
la

tio
n

Subtype Virus−2 Cells

Figure 5. R2 > 1 and R1 < R2: subtype-2 wins the competition

6. Discussion

It is widely recognized that time delays cause sustained oscillations in form of periodic solutions in in-

host models with cell divisions and intracellular delays [2]. It is interesting to explore the dynamics of the

viral load for two strains with intracellular delays both from mathematical and biological perspective [9].

In this paper, we employ a two-strain mathematical model to study the mechanistic basis of the

emergence of the competitive viral strains in host cells. We have carried out complete analysis for

two-strain in-host model with intracellular delays system (1.3). The analysis suggests that the basic

reproductive ratio palys an important role in predicting viral persistence or eradication. The global

dynamics of model (1.3) is rigorously established: if the basic reproduction number R0 ≤ 1, then all
solutions converge to the infection-free equilibirum E0; if R0 > 1, then all positive solutions converge

to the single chronic-infection equilibrium E1 or E2 which is determined by the relative magnitudes of

R1 and R2. The stability results for E0, E1 and E2 are obtained analytically, while the stability of the

co-existence equilibrium E3 via numerical simulations.

The intracellular delays can reduce the basic reproduction number R0 if cell die during the delay

period (3.1) (3.2). As a consequence, ignoring the delay will produce overestimation of R0. Our result

shows that no sustained oscillation regime exists without cell division even in the presence of intracellular

delays. The two-strain HIV model with intracellular delays could provide worthwile information that

potentially could allow the design of efficient individual strategies of HIV treatment.



248 J. XU

References

[1] S. Bonhoeffer, R. M. May, G. M. Shaw and M. A. Nowak, Virus dynamics and drug therapy, Proc.Natl.Acad.Sci.USA,

94 (1997), 6971–6976.

[2] R. V. Culshaw and S. G. Ruan, A delay-differential equation model of HIV infection of CD4+ T-cells, Mathematical

Biosciences, 165 (2000), 27–39.

[3] R. V. Culshaw, S. Ruan and G. Webb, A mathematical model of cell-to-cell spread of HIV-1 that includes a time

delay, Journal of Mathematical Biology, 46 (2003), 425–444.

[4] E. Fanales-Belasio, M. Raimondo, B. Suligoi and S. Butto, HIV virology and pathogenetic mechanisms of infection:

a brief overview, annali dell istituto superiore di sanita, 46, (2010), 5–14.

[5] H. R. Gelderblom, M. Ozel and G. Pauli, Morphogenesis and morphology of HIV structure-function relations,

Archives of Virology, 106 (1989), 1–13.

[6] H. Guo, M. Y. Li and Z. Shuai, A graph-theoretic approach to the method of global Lyapunov functions, Proc.

Amer.Math.Soc, 136 (2008), 2793–2802.

[7] J. K. Hale and L. S. Verduyn, Introduction to Functional Differential Equations, Springer, New York, 1993.

[8] V. Herz, S. Bonhoeffer, R. Anderson, R. May and M. Nowak, Viral dynamics in vivo: limitations on estimations on

intracellular delay and virus decay, Proc.Natl.Acad.Sci.USA, 93 (1996), 7247–7251.

[9] A. Jung, R. Maier, J. P. Vartanian, G. Bocharov, V. Jung, U. Fischer, E. Meese, S. Wain-Hobson and A. Meyerhans,

Mutiply infected spleen cells in HIV patients, Nature, 418 (2002), 144.

[10] T. Kepler and A. Perelson, Drug concentration heterogeneity facilitates the evolution of drug resistance,

Proc.Natl.Acad.Sci.USA, 95 (1998), 11514–11519.

[11] A. Korobeinikov and P. K. Maini, A Lyapunov function and global properties for SIR and SEIR eidemiological

models with nonlinear incidence, Math. Biosci. Eng, 1 (2004), 57–60.

[12] M. Y. Li and H. Shu, Global Dynamics of an In-host Viral Model with Intracellular Delay, Bulletin of Mathematical

Biology, 72 (2010), 1492–1505.

[13] K. Manavi, A review on infection with human immunodeficiency virus, Best Practic and Research Clinical Obsterics

and Gynaecology, 20 (2006), 923–940.

[14] P. Nelson, J. Mittler and A. Perelson, Effect of drug efficacy and the eclipse phase of the viral life cycle on estimates

of HIV-1 viral dynamic parameters, Journal of Acquired Immune Deficiency Syndromes, 26 (2001), 405–412.

[15] P. Nelson, J. Murray and A. Perelson, A model of HIV-1 pathogenesis that includes an intracellular delay, Mathe-

matical Biosciences, 163 (2000), 201–215.

[16] M. A. Nowak and R. M. May, Virus dynamics mathematical principles of immunology and virology, Oxford University

Press, 2000.

[17] A. Perelson, D. Kirschner and R. De Boer, Dynamics of HIV infection of CD4+ T cells, Mathematical Biosciences,

114 (1993), 81–125.

[18] A. Perelson and P. Nelson, Mathematical models of HIV dynamics in vivo, SIAM Review, 41 (1999), 3–44.

[19] A. Perelson, A. Neumann, M. Markowitz, J. Leonard and D. Ho, HIV-1 dynamics in vivo: viron clearance rate,

infected cell life-span, and viral generation time, Science, 271 (1996), 1582–1586.

[20] S. M. A. Rahman and X. Zou, Flu epidmics: A two strain flu model with a single vaccination, Journal of Biological

Dynamics, DOI:10.1080/17513758.2010.510213.

[21] L. Rong, Z. Feng and A. S. Perelson, Emergence of HIV-1 Drug Resistance During Antiretroviral Treatment, Bulletin

of Mathematical Biology, 69 (2007), 2027–2060.

[22] J. H. Strauss and E. G. Strauss, Viruses and Human Disease, Academic Press, 2008.

[23] H. Zhu and X. Zou, Impact of delays in cell infection and virus production on HIV-1 dynamics, Mathematical

Medicine and Biology, 25 (2008), 99–112.

School of Artificial Intelligence, Jianghan University, Wuhan 430056, P. R. China

E-mail address: jxu259@protonmail.com