# Dynamics Analysis and Simulation of a Modified HIV Infection Model with a Saturated Infection Rate

Authors: Qilin Sun; Lequan Min
Journal: Computational and Mathematical Methods in Medicine
 (2014)

---

## Abstract

This paper studies a modified human immunodeficiency virus (HIV) infection differential equation model with a saturated infection rate. It is proved that if the basic virus reproductive numberR
0 of the model is less than one, then the infection-free equilibrium point of the model is globally asymptotically stable; if R
0 of the model is more than one, then the endemic infection equilibrium point of the model is globally asymptotically stable. Based on the clinical data from HIV drug resistance database of Stanford University, using the proposed model simulates the dynamics of the two groups of patients’ anti-HIV infection treatment. The numerical simulation results are in agreement with the evolutions of the patients’ HIV RNA levels. It can be assumed that if an HIV infected individual’s basic virus reproductive number R
0
<
1 then this person will recover automatically; if an antiretroviral therapy makes an HIV infected individual’s R
0
<
1, this person will be cured eventually; if an antiretroviral therapy fails to suppress an HIV infected individual’s HIV RNA load to be of unpredictable level, the time that the patient’s HIV RNA level has achieved the minimum value may be the starting time that drug resistance has appeared.


---

## Body

## 1. Introduction

The human immunodeficiency virus (HIV) mainly targets a host’sCD
4
+ T cells. Chronic HIV infection causes gradual depletion of the CD
4
+ T cell pool and thus progressively compromises the hosts immune response to opportunistic infections, leading to Acquired Immunodeficiency Syndrome (AIDS) [1].In recent years, there is much work done on HIV infection from different points of view, such as pathology [2], microbiology [3], and mathematics [4–7]. Mathematical models have become essential tools to make assumptions, suggest new experiments, or help easily explain complex processes [8]. The basic mathematical model widely used for studying the dynamics of HIV infection has the following form [4, 9]:
(1)
x
˙
=
λ
-
d
1
x
-
k
1
x
v
,
y
˙
=
k
1
x
v
-
d
2
y
,
v
˙
=
a
y
-
d
3
v
,
where x
(
t
)
,
y
(
t
), and v
(
t
) are the number of uninfected cells, infected cells, and free virus, respectively. Uninfected cells are produced at a constant rate λ, die at rate d
1
x, and become infected at rate k
1
x
v. Infected cells are produced at rate k
1
x
v and die at rate d
2
y. Free virus is produced from infected cells at rate a
y and dies at rate d
3
v.Equation (1) has a basic virus reproductive number R
0
=
a
k
1
λ
/
(
d
1
d
2
d
3
). According to Nowak and Bangham [4], R
0 is defined as the number of newly infected cells arising from any one infected cell; if R
0 is smaller than 1, then in the beginning of the infection, each virus infected cell produces on average less than one newly infected cell. Thus, the infection cannot spread, and the system returns to the uninfected state; if R
0 is larger than 1, then initially each virus infected cell produces on average more than one newly infected cell. The infected cell population will increase, whereas the uninfected cell population will decline and therefore provide less opportunity for the virus to infect new cells.There is a discussion about the process of the HIV RNA transcribing into DNA: when an HIV enters a restingCD
4
+ T cell, the HIV RNA may not be completely reverse transcribed into DNA [10]. A proportion of resting infected cells can revert to the uninfected state before the viral genome is integrated into the genome of the lymphocyte [11].Recently, some mathematical models of HIV infection have been proposed based on the assumption that a fraction of infectedCD
4
+ T cells return to the uninfected class [12–14]. Srivastava and Chandra [13] have considered a model with three populations: uninfected CD
4
+ T cells (
x
), infected CD
4
+ T cells (
y
), and HIV population (
v
). The model has the following form:
(2)
x
˙
=
λ
-
d
1
x
-
k
1
x
v
+
p
y
,
y
˙
=
k
1
x
v
-
(
d
2
+
p
)
y
,
v
˙
=
a
y
-
d
3
v
,
where the meanings of the variables x
(
t
), y
(
t
), and v
(
t
) and the parameters λ
,
d
1
,
k
1
,
d
2
,
a, and d
3 are the same as those given in (1). The term p
y is the rate of infected cells in the latent stage reverting to the uninfected class. Equation (2) also has a basic virus reproductive number R
0
=
a
k
1
λ
/
(
d
3
d
1
(
p
+
d
2
)
). They have proved that if R
0
≤
1, the infection-free steady state of (2) is globally asymptotically stable; if R
0
>
1, the endemic steady state of (2) is globally asymptotically stable [13].In (2), the mass action term k
1
x
v used to model infection of CD
4
+ T cells by free virions is biologically problematic for several reasons. Firstly, since λ
/
d
1 represents the total number of CD
4
+ T cells in the basic virus reproductive number R
0
=
a
k
1
λ
/
(
d
3
d
1
(
p
+
d
2
)
), this causes R
0 to depend upon the total number of CD
4
+ T cells in vivo. This implies the dubious prediction that individuals with more CD
4
+ T cells will be more easily infected than individuals with less CD
4
+ T cells. Secondly, the rate of HIV infection is assumed to be bilinear by the mass action term k
1
x
v. However, the actual incidence rate is probably not linear over the entire range of virus v
(
t
) and uninfected CD
4
+ T cells x
(
t
) [15–17].On biological grounds, during primary HIV infection, the rate of virus infection should be approximately proportionate to the virus loadk
1
v because of a small amount of viral load with respect to a large number of CD
4
+ T cells. However,since the total number of healthy CD
4
+ T cells in vivo is limited, the HIV infection will approach saturation with more and more virus produced. In this case, it is more reasonable to assume that the rate of virus infection should be approximately proportionate to the number of healthy CD
4
+ T cells k
1
x.Based on the argument above, this paper describes an amended model. In this model, we use a saturated infection ratek
1
x
v
/
(
x
+
v
) to replace the mass action term k
1
x
v in (2). Under the formulation of this saturated infection rate, the basic virus reproductive number R
0 is independent of the total number of CD
4
+ T cells. Meanwhile, the actual incidence rate is not linear over the entire range of virus v
(
t
) and uninfected CD
4
+ T cells x
(
t
) any more. The global stabilities of the infection-free state and the endemic infection state of the modified HIV infection model have been discussed. Based on the clinical data from HIV drug resistance database of Stanford University, using the proposed model simulates the dynamics of two groups of patients’ anti-HIV infection treatment, and then make long-term predictions for the two groups’ anti-HIV infection treatment, respectively.The rest of this paper is organized as follows. Section2 introduces a modified model and discusses the boundedness of the solutions of the model. Sections 3 and 4 discuss the global stability of the infection-free state and the endemic infection state of the modified HIV infection model, respectively. Section 5 simulates the dynamics of two groups of patients’ anti-HIV infection treatment. Section 6 summarizes this paper.


## 2. Modified HIV Infection Model

### 2.1. The Modified HIV Infection Model

Based on (2), our modified HIV infection model has the following form:
(3)
x
˙
=
λ
-
d
1
x
-
k
1
x
v
x
+
v
+
p
y
,
y
˙
=
k
1
x
v
x
+
v
-
(
d
2
+
p
)
y
,
v
˙
=
a
y
-
d
3
v
,
where the meanings of the variables x
(
t
), y
(
t
), and v
(
t
) and the parameters λ
,
d
1
,
k
1
,
d
2
,
p
,
a, and d
3 are the same as those given in (2). Equation (3) has two steady states: (1)
The infection-free steady state(4)
Q
1
=
(
x
0
,
0,0
)
represents the virus infection free.Q
1 is called infection-free equilibrium point. Here,
(5)
x
0
=
λ
d
1
.
(2)
The endemic infected steady state(6)
Q
2
=
(
x
¯
,
y
¯
,
v
¯
)
represents persistent virus infection.Q
2 is called endemic infection equilibrium point. Here,
(7)
x
¯
=
λ
R
0
k
1
(
R
0
-
1
)
+
d
1
R
0
-
p
(
R
0
-
1
)
R
0
(
d
3
/
a
)
,
y
¯
=
d
3
a
(
R
0
-
1
)
x
¯
,
v
¯
=
a
y
¯
d
3
=
(
R
0
-
1
)
x
¯
.
Here,(8)
R
0
=
a
k
1
d
3
(
d
2
+
p
)
. Since the total rate of disappearance of infected cells is d
2
+
p, infected cells live on average for time 1
/
(
d
2
+
p
). Each infected cell produces virus at rate a. Thus, each infected cell produces on average a total of a
/
(
d
2
+
p
) viruses. Since virus dies at rate d
3 per virion, each virus survives on average for time 1
/
d
3. During the time 1
/
d
3, each virus infects on average k
1
x
0
/
(
x
0
+
v
0
)
d
3 cells, where x
0 and v
0 are the preinfection target cells’ density and viruses’ density, respectively. Thus, the total number of cells infected by the a
/
(
d
2
+
p
) viruses is a
k
1
x
0
/
(
x
0
+
v
0
)
d
3
(
d
2
+
p
). According to (4) and (5), x
0
=
λ
/
d
1 and v
0
=
0 at the preinfection steady state. Then one can obtain that the total number of cells infected by each infected cell is a
k
1
/
d
3
(
d
2
+
p
). Hence, R
0 is the basic virus reproductive number of (3) which is independent of the total number of the uninfected CD
4
+ T cells.According to (4), (6), and (7), if R
0
≤
1, then Q
1 is the unique infection-free equilibrium point; if R
0
>
1, then, in addition to the infection-free equilibrium point, (3) has another equilibrium point Q
2.


### 2.2. Boundedness of Solutions

It is easy to show that the solutions of (3) with initial conditions x
(
0
)
>
0
,
y
(
0
)
>
0, and v
(
0
)
>
0 have all positive components for t
>
0. Hence, one begins the analysis of (3) by observing the nonnegative octant
(9)
D
=
{
(
x
,
y
,
v
)
∈
R
+
3
:
x
≥
0
,
y
≥
0
,
v
≥
0
}
.According to the first two equations of (3), one can get
(10)
x
˙
+
y
˙
=
λ
-
d
1
x
-
d
2
y
≤
λ
-
d
(
x
+
y
)
,
d
=
min
⁡
(
d
1
,
d
2
)
,
and then
(11)
x
+
y
≤
λ
d
.
So x
(
t
) and y
(
t
) are bounded. From the last equation of (3), it follows that
(12)
v
˙
=
a
y
-
d
3
v
≤
a
λ
d
-
d
3
v
,
and then
(13)
v
≤
a
λ
d
d
3
.
So v
(
t
) are bounded. Hence there is a bounded subset of D:
(14)
Ω
=
{
(
x
,
y
,
v
)
∈
R
+
3
:
0
≤
x
+
y
≤
λ
d
,
0
≤
v
≤
a
λ
d
d
3
}
such that any solution trajectory (
x
(
t
)
,
y
(
t
)
,
v
(
t
)
) of (3) with initial value (
x
(
0
)
,
y
(
0
)
,
v
(
0
)
) in Ω will keep in the subset Ω.According to (7), x
¯
>
0, y
¯
>
0, and v
¯
>
0. One can get that the endemic infection equilibrium point Q
2 exists in the interior of Ω:
(15)
Ω
0
=
{
a
λ
d
d
3
(
x
,
y
,
v
)
∈
Ω
:
0
<
x
,
0
<
y
,
0
<
x
+
y
<
λ
d
,
0
<
v
<
a
λ
d
d
3
}
.
Therefore, the stability of the endemic infection equilibrium point Q
2 only needs to be discussed in Ω
0.


## 2.1. The Modified HIV Infection Model

Based on (2), our modified HIV infection model has the following form:
(3)
x
˙
=
λ
-
d
1
x
-
k
1
x
v
x
+
v
+
p
y
,
y
˙
=
k
1
x
v
x
+
v
-
(
d
2
+
p
)
y
,
v
˙
=
a
y
-
d
3
v
,
where the meanings of the variables x
(
t
), y
(
t
), and v
(
t
) and the parameters λ
,
d
1
,
k
1
,
d
2
,
p
,
a, and d
3 are the same as those given in (2). Equation (3) has two steady states: (1)
The infection-free steady state(4)
Q
1
=
(
x
0
,
0,0
)
represents the virus infection free.Q
1 is called infection-free equilibrium point. Here,
(5)
x
0
=
λ
d
1
.
(2)
The endemic infected steady state(6)
Q
2
=
(
x
¯
,
y
¯
,
v
¯
)
represents persistent virus infection.Q
2 is called endemic infection equilibrium point. Here,
(7)
x
¯
=
λ
R
0
k
1
(
R
0
-
1
)
+
d
1
R
0
-
p
(
R
0
-
1
)
R
0
(
d
3
/
a
)
,
y
¯
=
d
3
a
(
R
0
-
1
)
x
¯
,
v
¯
=
a
y
¯
d
3
=
(
R
0
-
1
)
x
¯
.
Here,(8)
R
0
=
a
k
1
d
3
(
d
2
+
p
)
. Since the total rate of disappearance of infected cells is d
2
+
p, infected cells live on average for time 1
/
(
d
2
+
p
). Each infected cell produces virus at rate a. Thus, each infected cell produces on average a total of a
/
(
d
2
+
p
) viruses. Since virus dies at rate d
3 per virion, each virus survives on average for time 1
/
d
3. During the time 1
/
d
3, each virus infects on average k
1
x
0
/
(
x
0
+
v
0
)
d
3 cells, where x
0 and v
0 are the preinfection target cells’ density and viruses’ density, respectively. Thus, the total number of cells infected by the a
/
(
d
2
+
p
) viruses is a
k
1
x
0
/
(
x
0
+
v
0
)
d
3
(
d
2
+
p
). According to (4) and (5), x
0
=
λ
/
d
1 and v
0
=
0 at the preinfection steady state. Then one can obtain that the total number of cells infected by each infected cell is a
k
1
/
d
3
(
d
2
+
p
). Hence, R
0 is the basic virus reproductive number of (3) which is independent of the total number of the uninfected CD
4
+ T cells.According to (4), (6), and (7), if R
0
≤
1, then Q
1 is the unique infection-free equilibrium point; if R
0
>
1, then, in addition to the infection-free equilibrium point, (3) has another equilibrium point Q
2.


## 2.2. Boundedness of Solutions

It is easy to show that the solutions of (3) with initial conditions x
(
0
)
>
0
,
y
(
0
)
>
0, and v
(
0
)
>
0 have all positive components for t
>
0. Hence, one begins the analysis of (3) by observing the nonnegative octant
(9)
D
=
{
(
x
,
y
,
v
)
∈
R
+
3
:
x
≥
0
,
y
≥
0
,
v
≥
0
}
.According to the first two equations of (3), one can get
(10)
x
˙
+
y
˙
=
λ
-
d
1
x
-
d
2
y
≤
λ
-
d
(
x
+
y
)
,
d
=
min
⁡
(
d
1
,
d
2
)
,
and then
(11)
x
+
y
≤
λ
d
.
So x
(
t
) and y
(
t
) are bounded. From the last equation of (3), it follows that
(12)
v
˙
=
a
y
-
d
3
v
≤
a
λ
d
-
d
3
v
,
and then
(13)
v
≤
a
λ
d
d
3
.
So v
(
t
) are bounded. Hence there is a bounded subset of D:
(14)
Ω
=
{
(
x
,
y
,
v
)
∈
R
+
3
:
0
≤
x
+
y
≤
λ
d
,
0
≤
v
≤
a
λ
d
d
3
}
such that any solution trajectory (
x
(
t
)
,
y
(
t
)
,
v
(
t
)
) of (3) with initial value (
x
(
0
)
,
y
(
0
)
,
v
(
0
)
) in Ω will keep in the subset Ω.According to (7), x
¯
>
0, y
¯
>
0, and v
¯
>
0. One can get that the endemic infection equilibrium point Q
2 exists in the interior of Ω:
(15)
Ω
0
=
{
a
λ
d
d
3
(
x
,
y
,
v
)
∈
Ω
:
0
<
x
,
0
<
y
,
0
<
x
+
y
<
λ
d
,
0
<
v
<
a
λ
d
d
3
}
.
Therefore, the stability of the endemic infection equilibrium point Q
2 only needs to be discussed in Ω
0.


## 3. Stability of the Infection-Free Equilibrium PointQ
1

In this section, we discuss locally asymptotical stability and globally asymptotical stability of the infection-free equilibrium pointQ
1 of (3).


### 3.1. Locally Asymptotical Stability of the Infection-Free Equilibrium Point Q
1

Theorem 1.
IfR
0
=
a
k
1
/
(
d
3
(
d
2
+
p
)
)
<
1, then the infection-free equilibrium point Q
1 of (3) is locally asymptotically stable. If R
0
>
1, then the infection-free equilibrium point Q
1 is unstable.Proof.
The Jacobi matrix of (3) at an arbitrary point is given by
(16)
J
(
x
,
y
,
v
)
=
[
-
d
1
-
a
1
p
-
a
2
a
1
-
p
-
d
2
a
2
0
a
-
d
3
]
,
where a
1
=
k
1
v
2
/
(
x
+
v
)
2 and a
2
=
k
1
x
2
/
(
x
+
v
)
2.
Substituting the equilibrium pointQ
1 into matrix (16) gives
(17)
J
Q
1
=
[
-
d
1
p
-
k
1
0
-
p
-
d
2
k
1
0
a
-
d
3
]
.
The corresponding eigenequation is
(18)
|
λ
E
-
J
Q
1
|
=
|
λ
+
d
1
-
p
k
1
0
λ
+
p
+
d
2
-
k
1
0
-
a
λ
+
d
3
|
=
0
.
Solving
(19)
|
λ
E
-
J
Q
1
|
=
(
λ
+
d
1
)
[
(
λ
+
p
+
d
2
)
(
λ
+
d
3
)
-
a
k
1
]
=
0
gives
(20)
λ
1
=
-
d
1
<
0
,
(21)
(
λ
+
d
2
+
p
)
(
λ
+
d
3
)
-
a
k
1
k
3
=
0
.
Equation (21) can be written as
(22)
λ
2
+
(
d
2
+
p
+
d
3
)
λ
+
d
3
(
p
+
d
2
)
-
a
k
1
=
0
.Solving equation (22) gives
(23)
λ
2
=
-
(
d
2
+
p
+
d
3
)
-
(
d
2
+
p
+
d
3
)
2
-
4
[
d
3
(
p
+
d
2
)
-
a
k
1
]
2
=
-
(
d
2
+
p
+
d
3
)
-
(
d
2
+
p
+
d
3
)
2
-
4
d
3
(
p
+
d
2
)
(
1
-
R
0
)
2
,
λ
3
=
-
(
d
2
+
p
+
d
3
)
+
(
d
2
+
p
+
d
3
)
2
-
4
[
d
3
(
p
+
d
2
)
-
a
k
1
]
2
=
-
(
d
2
+
p
+
d
3
)
+
(
d
2
+
p
+
d
3
)
2
-
4
d
3
(
p
+
d
2
)
(
1
-
R
0
)
2
.IfR
0
<
1, then λ
2
<
0 and λ
3
<
0. Hence the infection-free equilibrium point Q
1 is locally asymptotically stable. If R
0
>
1, then λ
3
>
0 such that the infection-free equilibrium point Q
1 is unstable.


### 3.2. Globally Asymptotical Stability of the Infection-Free Equilibrium Point Q
1

Theorem 2.
IfR
0
<
1, then the infection-free equilibrium point Q
1 of (3) is globally asymptotically stable in Ω.Proof.
Define a global Lyapunov function by(24)
V
1
(
x
,
y
,
v
)
=
y
+
(
d
2
+
p
)
v
a
.
The derivative of V
1
(
x
,
y
,
v
) along the positive solutions of (3) is
(25)
V
˙
1
=
y
˙
+
(
d
2
+
p
)
v
˙
a
=
k
1
x
v
x
+
v
-
d
2
y
-
p
y
+
(
d
2
+
p
)
y
-
d
3
(
d
2
+
p
)
v
a
≤
k
1
v
-
d
3
(
d
2
+
p
)
v
a
=
[
a
k
1
d
3
(
d
2
+
p
)
-
1
]
d
3
(
d
2
+
p
)
a
v
=
(
R
0
-
1
)
d
3
(
d
2
+
p
)
a
v
.
IfR
0
<
1, then V
˙
1
≤
0 holds in Ω. Moreover, V
˙
1
=
0 if and only if v
=
0. Hence, the largest compact invariant set in Ω is
(26)
E
1
=
{
(
x
,
y
,
v
)
∈
Ω
∣
V
˙
1
=
0
}
=
{
(
x
,
y
,
v
)
∈
Ω
∣
v
=
0
}
.
According to the LaSalle’s invariance principle, lim
⁡
t
→
+
∞
v
(
t
)
=
0. Then one can get limit equations:
(27)
x
˙
=
λ
-
d
1
x
+
p
y
,
y
˙
=
-
d
2
y
-
p
y
.
Define a global Lyapunov function by(28)
V
2
(
x
,
y
)
=
x
-
x
0
-
x
0
ln
⁡
x
x
0
+
y
,
where
(29)
λ
=
d
1
x
0
,
Q
1
=
(
x
0
,
0,0
)
.
The derivative of V
2
(
x
,
y
) along the positive solutions of (27) is
(30)
V
˙
2
=
x
˙
+
y
˙
-
x
0
x
x
˙
=
λ
-
d
1
x
+
p
y
-
d
2
y
-
p
y
-
x
0
x
(
λ
-
d
1
x
+
p
y
)
.
Since λ
=
d
1
x
0,
(31)
V
˙
2
=
d
1
x
0
-
d
1
x
-
x
0
x
d
1
x
0
+
d
1
x
0
-
x
0
x
p
y
-
d
2
y
=
d
1
x
0
[
2
-
x
x
0
-
x
0
x
]
-
(
x
0
x
p
+
d
2
)
y
.
Since the arithmetic mean is greater than or equal to the geometric mean, we obtain2
-
(
x
/
x
0
)
-
(
x
0
/
x
)
≤
0.
Therefore,V
˙
2
≤
0 holds in E
1, and V
˙
2
=
0 if and only if x
=
x
0 and y
=
0. There is the largest compact invariant set in E
1:
(32)
E
2
=
{
(
x
,
y
,
v
)
∈
E
1
∣
V
˙
2
=
0
}
=
{
(
x
,
y
,
v
)
∈
E
1
∣
x
=
x
0
,
y
=
0
}
=
{
Q
1
}
.
Hence ifR
0
<
1, all solution paths in Ω approach the infection-free equilibrium point Q
1.


## 3.1. Locally Asymptotical Stability of the Infection-Free Equilibrium Point Q
1

Theorem 1.
IfR
0
=
a
k
1
/
(
d
3
(
d
2
+
p
)
)
<
1, then the infection-free equilibrium point Q
1 of (3) is locally asymptotically stable. If R
0
>
1, then the infection-free equilibrium point Q
1 is unstable.Proof.
The Jacobi matrix of (3) at an arbitrary point is given by
(16)
J
(
x
,
y
,
v
)
=
[
-
d
1
-
a
1
p
-
a
2
a
1
-
p
-
d
2
a
2
0
a
-
d
3
]
,
where a
1
=
k
1
v
2
/
(
x
+
v
)
2 and a
2
=
k
1
x
2
/
(
x
+
v
)
2.
Substituting the equilibrium pointQ
1 into matrix (16) gives
(17)
J
Q
1
=
[
-
d
1
p
-
k
1
0
-
p
-
d
2
k
1
0
a
-
d
3
]
.
The corresponding eigenequation is
(18)
|
λ
E
-
J
Q
1
|
=
|
λ
+
d
1
-
p
k
1
0
λ
+
p
+
d
2
-
k
1
0
-
a
λ
+
d
3
|
=
0
.
Solving
(19)
|
λ
E
-
J
Q
1
|
=
(
λ
+
d
1
)
[
(
λ
+
p
+
d
2
)
(
λ
+
d
3
)
-
a
k
1
]
=
0
gives
(20)
λ
1
=
-
d
1
<
0
,
(21)
(
λ
+
d
2
+
p
)
(
λ
+
d
3
)
-
a
k
1
k
3
=
0
.
Equation (21) can be written as
(22)
λ
2
+
(
d
2
+
p
+
d
3
)
λ
+
d
3
(
p
+
d
2
)
-
a
k
1
=
0
.Solving equation (22) gives
(23)
λ
2
=
-
(
d
2
+
p
+
d
3
)
-
(
d
2
+
p
+
d
3
)
2
-
4
[
d
3
(
p
+
d
2
)
-
a
k
1
]
2
=
-
(
d
2
+
p
+
d
3
)
-
(
d
2
+
p
+
d
3
)
2
-
4
d
3
(
p
+
d
2
)
(
1
-
R
0
)
2
,
λ
3
=
-
(
d
2
+
p
+
d
3
)
+
(
d
2
+
p
+
d
3
)
2
-
4
[
d
3
(
p
+
d
2
)
-
a
k
1
]
2
=
-
(
d
2
+
p
+
d
3
)
+
(
d
2
+
p
+
d
3
)
2
-
4
d
3
(
p
+
d
2
)
(
1
-
R
0
)
2
.IfR
0
<
1, then λ
2
<
0 and λ
3
<
0. Hence the infection-free equilibrium point Q
1 is locally asymptotically stable. If R
0
>
1, then λ
3
>
0 such that the infection-free equilibrium point Q
1 is unstable.


## 3.2. Globally Asymptotical Stability of the Infection-Free Equilibrium Point Q
1

Theorem 2.
IfR
0
<
1, then the infection-free equilibrium point Q
1 of (3) is globally asymptotically stable in Ω.Proof.
Define a global Lyapunov function by(24)
V
1
(
x
,
y
,
v
)
=
y
+
(
d
2
+
p
)
v
a
.
The derivative of V
1
(
x
,
y
,
v
) along the positive solutions of (3) is
(25)
V
˙
1
=
y
˙
+
(
d
2
+
p
)
v
˙
a
=
k
1
x
v
x
+
v
-
d
2
y
-
p
y
+
(
d
2
+
p
)
y
-
d
3
(
d
2
+
p
)
v
a
≤
k
1
v
-
d
3
(
d
2
+
p
)
v
a
=
[
a
k
1
d
3
(
d
2
+
p
)
-
1
]
d
3
(
d
2
+
p
)
a
v
=
(
R
0
-
1
)
d
3
(
d
2
+
p
)
a
v
.
IfR
0
<
1, then V
˙
1
≤
0 holds in Ω. Moreover, V
˙
1
=
0 if and only if v
=
0. Hence, the largest compact invariant set in Ω is
(26)
E
1
=
{
(
x
,
y
,
v
)
∈
Ω
∣
V
˙
1
=
0
}
=
{
(
x
,
y
,
v
)
∈
Ω
∣
v
=
0
}
.
According to the LaSalle’s invariance principle, lim
⁡
t
→
+
∞
v
(
t
)
=
0. Then one can get limit equations:
(27)
x
˙
=
λ
-
d
1
x
+
p
y
,
y
˙
=
-
d
2
y
-
p
y
.
Define a global Lyapunov function by(28)
V
2
(
x
,
y
)
=
x
-
x
0
-
x
0
ln
⁡
x
x
0
+
y
,
where
(29)
λ
=
d
1
x
0
,
Q
1
=
(
x
0
,
0,0
)
.
The derivative of V
2
(
x
,
y
) along the positive solutions of (27) is
(30)
V
˙
2
=
x
˙
+
y
˙
-
x
0
x
x
˙
=
λ
-
d
1
x
+
p
y
-
d
2
y
-
p
y
-
x
0
x
(
λ
-
d
1
x
+
p
y
)
.
Since λ
=
d
1
x
0,
(31)
V
˙
2
=
d
1
x
0
-
d
1
x
-
x
0
x
d
1
x
0
+
d
1
x
0
-
x
0
x
p
y
-
d
2
y
=
d
1
x
0
[
2
-
x
x
0
-
x
0
x
]
-
(
x
0
x
p
+
d
2
)
y
.
Since the arithmetic mean is greater than or equal to the geometric mean, we obtain2
-
(
x
/
x
0
)
-
(
x
0
/
x
)
≤
0.
Therefore,V
˙
2
≤
0 holds in E
1, and V
˙
2
=
0 if and only if x
=
x
0 and y
=
0. There is the largest compact invariant set in E
1:
(32)
E
2
=
{
(
x
,
y
,
v
)
∈
E
1
∣
V
˙
2
=
0
}
=
{
(
x
,
y
,
v
)
∈
E
1
∣
x
=
x
0
,
y
=
0
}
=
{
Q
1
}
.
Hence ifR
0
<
1, all solution paths in Ω approach the infection-free equilibrium point Q
1.


## 4. Stability of the Endemic Infection Equilibrium PointQ
2

In this section, we analyze local asymptotical stability and global asymptotical stability of the endemic infection equilibrium pointQ
2 of (3).


### 4.1. Locally Asymptotical Stability of the Endemic Infection Equilibrium PointQ
2

Theorem 3.
IfR
0
>
1, then the endemic infection equilibrium point Q
2 of (3) is locally asymptotically stable.Proof.
Put the equilibrium pointQ
2 into matrix (16); then one obtains
(33)
J
Q
2
=
[
-
d
1
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
p
-
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
-
d
2
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
0
a
-
d
3
]
,
|
λ
E
-
J
Q
2
|
=
|
λ
+
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
λ
+
p
+
d
2
-
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
0
-
a
λ
+
d
3
|
=
0
.
Solving the eigenequation of the matrix above, here is(34)
λ
3
+
[
d
1
+
d
2
+
d
3
+
p
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
]
︸
a
1
λ
2
+
[
d
3
(
d
2
+
p
)
+
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
(
d
2
+
p
+
d
3
)
-
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
]
︸
a
2
λ
+
d
3
(
d
2
+
p
)
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
+
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
︸
a
3
=
0
.If R
0
>
1, one obtains that(35)
a
1
=
d
1
+
d
2
+
d
3
+
p
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
>
0
,
a
2
=
d
3
(
d
2
+
p
)
+
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
(
d
2
+
p
+
d
3
)
-
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
=
a
k
1
R
0
+
d
1
(
d
2
+
d
3
+
p
)
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
d
3
)
-
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
=
a
k
1
R
0
-
a
k
1
x
¯
2
R
0
2
x
¯
2
+
d
1
(
d
2
+
d
3
+
p
)
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
d
3
)
=
a
k
1
(
R
0
-
1
)
R
0
2
+
d
1
(
d
2
+
d
3
+
p
)
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
d
3
)
>
0
,
a
3
=
d
3
(
d
2
+
p
)
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
+
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
=
d
1
d
3
(
d
2
+
p
)
-
a
d
1
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
+
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
=
d
1
a
k
1
(
R
0
-
1
)
R
0
2
+
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
>
0
,
a
1
a
2
-
a
3
=
(
d
2
+
d
3
+
p
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
×
a
k
1
(
R
0
-
1
)
R
0
2
+
d
1
(
d
2
+
d
3
+
p
)
×
(
d
2
+
d
3
+
p
+
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
+
d
2
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
d
3
+
p
+
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
+
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
3
+
p
+
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
>
0
.
Hence all inequalities of the Routh-Hurwitz criterion are satisfied. Therefore, the endemic infection equilibrium pointQ
2 is locally asymptotically stable.


### 4.2. Globally Asymptotical Stability of the Endemic Infection Equilibrium PointQ
2

In this subsection, we firstly introduce a lemma outlined by Li and Wang [18], and then using this lemma discusses the globally asymptotical stability of the endemic infection equilibrium point Q
2 of (3).The lemma is briefly summarized as follows.Letx
↦
f
(
x
)
∈
R
n be a C
1 function for x in an open set Γ
⊂
R
n. Consider the differential system
(36)
x
˙
=
f
(
x
)
.
Denote by x
(
t
,
x
0
) the solution to (36) such that x
(
0
,
x
0
)
=
x
0. Let x
~ be an equilibrium point of (36). Li and Wang [18] made the following two basic assumptions:(
H
1
)
there exists a compact absorbing setK
⊂
Γ;
(
H
2
)
equation (36) has a unique equilibrium x
~ in Γ.Li and Wang (see Theorem  2.5 in [18]) have given the following lemma.Lemma 4 (see [18]).
Assume that(1)
assumptions (H
1) and (H
2) hold;
(2)
equation (36) satisfies the Poincaré-Bendixson Property;
(3)
for each periodic solutionx
=
p
(
t
) to (36) with p
(
0
)
∈
Γ, the linear system (the second additive compound system)
(37)
w
˙
(
t
)
=
∂
f
[
2
]
∂
x
(
P
(
t
)
)
w
(
t
)
is asymptotically stable, where∂
f
[
2
]
/
∂
x is the second additive compound matrix of the Jacobian matrix ∂
f
/
∂
x;
(4)
(
-
1
)
n
det
(
(
∂
f
/
∂
x
)
(
x
~
)
)
>
0.
Then the unique equilibriumx
~ is globally asymptotically stable in Γ.Now one uses Lemma4 to show the following.Theorem 5.
IfR
0
>
1, then the endemic infection equilibrium point Q
2 of (3) is globally asymptotically stable in Ω
0, where Ω
0 is defined by (15).Proof.
Based on Lemma4, the proof of Theorem 5 has been implemented via the following four steps.
(
1
) For epidemic models and many other biological models where the feasible region is a bounded cone, (
H
1
) is equivalent to the uniform persistence of the system [19]. By (15), Ω
0 is bounded, so it only needs to show the uniform persistence of (3). According to Proposition  3.3 in [20], the necessary and sufficient condition for the uniform persistence of (3) is equivalent to the equilibrium point Q
1 being unstable. Theorem 1 has shown that Q
1 is unstable if R
0
>
1. Therefore, (3) is uniformly persistent if R
0
>
1 so that (
H
1
) holds if R
0
>
1.
Meanwhile,Q
1
=
(
x
0
,
0,0
) by (4), so Q
1 does not exist in Ω
0. Hence, Q
2 is the unique equilibrium point of (3) in Ω
0 so that (
H
2
) holds.
The results above verify the condition(
1
) of Lemma 4.
(
2
) The Jacobian matrix of (3) is
(38)
J
(
x
,
y
,
v
)
=
[
-
d
1
-
a
1
p
-
a
2
a
1
-
p
-
d
2
a
2
0
a
-
d
3
]
,
where a
1
=
(
k
1
v
2
/
(
x
+
v
)
2
) and a
2
=
(
k
1
x
2
/
(
x
+
v
)
2
).
IfH
=
diag
⁡
(
1
,
-
1,1
), then
(39)
H
J
H
=
[
-
d
1
-
k
1
v
2
x
+
v
-
p
-
k
1
x
2
(
x
+
v
)
2
-
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
-
k
1
x
2
(
x
+
v
)
2
0
-
a
-
d
3
]
,
and one can obtain that H
J
H has nonpositive off-diagonal elements in Ω
0. Therefore (3) is competitive in Ω
0. It is known that 3-dimensional competitive systems have the Poincaré-Bendixson Property [21]. Hence, (3) satisfies the Poincaré-Bendixson Property. This verifies condition (
2
) of Lemma 4.
(
3
) Let P
(
t
)
=
(
x
(
t
)
,
y
(
t
)
,
v
(
t
)
) be a periodic solution in Ω
0.
According to [22], if B
=
(
b
i
j
) is a 3
×
3 matrix, then the second additive compound B
[
2
] of B is
(40)
B
[
2
]
=
[
b
11
+
b
22
b
23
-
b
13
b
32
b
11
+
b
33
b
12
-
b
31
b
21
b
22
+
b
33
]
.
The Jacobian matrix of (3) is
(41)
J
(
x
,
y
,
v
)
=
[
-
d
1
-
a
1
p
-
a
2
a
1
-
p
-
d
2
a
2
0
a
-
d
3
]
,
where a
1
=
k
1
v
2
/
(
x
+
v
)
2 and a
2
=
k
1
x
2
/
(
x
+
v
)
2.
And then the second additive compound matrix of the Jacobian matrix of (3) is given by(42)
J
[
2
]
=
[
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
k
1
x
2
(
x
+
v
)
2
k
1
x
2
(
x
+
v
)
2
a
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
d
3
p
0
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
-
d
3
]
, and the second additive compound system of (3) along the periodic solution P
(
t
)
=
(
x
(
t
)
,
y
(
t
)
,
v
(
t
)
) is
(43)
w
˙
1
=
(
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
)
w
1
+
k
1
x
2
(
x
+
v
)
2
w
2
+
k
1
x
2
(
x
+
v
)
2
w
3
,
w
˙
2
=
a
w
1
+
(
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
d
3
)
w
2
+
p
w
3
,
w
˙
3
=
k
1
v
2
(
x
+
v
)
2
w
2
+
(
-
p
-
d
2
-
d
3
)
w
3
.
Define a global Lyapunov function by
(44)
V
(
w
1
,
w
2
,
w
3
,
P
)
=
∥
(
w
1
,
y
(
t
)
v
(
t
)
w
2
,
y
(
t
)
v
(
t
)
w
3
)
∥
,
where ∥
·
∥ is the norm in set D defined by
(45)
∥
(
w
1
,
w
2
,
w
3
)
∥
=
sup
⁡
{
|
w
1
|
,
|
w
2
|
+
|
w
3
|
}
.
Suppose that the solution P
(
t
) is periodic of least period ω
>
0 and that P
(
0
)
∈
Ω
0. According to [23], (3) is uniformly persistent, if there exists a positive constant μ such that
(46)
lim
⁡
inf
⁡
t
→
∞
x
(
t
)
≥
μ
,
lim
⁡
inf
⁡
t
→
∞
y
(
t
)
≥
μ
,
lim
⁡
inf
⁡
t
→
∞
v
(
t
)
≥
μ
.
Step (
1
) has shown that (3) is uniformly persistent if R
0
>
1. Hence, there always exists a positive constant μ which satisfies (46). The orbit of P
(
t
) remains at a positive distance from the boundary of Ω by the uniform persistence, and one can obtain that
(47)
y
(
t
)
≥
μ
,
v
(
t
)
≥
μ
for
large
enough
t
.
Sincev
<
a
λ
/
(
d
d
3
) by (15),
(48)
V
(
w
1
,
w
2
,
w
3
,
P
)
≥
μ
d
d
3
a
λ
∥
(
w
1
,
w
2
,
w
3
)
∥
,
for all (
w
1
,
w
2
,
w
3
)
∈
R
3.
Along a solution(
w
1
,
w
2
,
w
3
) of (43), V
(
w
1
,
w
2
,
w
3
,
P
) becomes
(49)
V
(
w
1
,
w
2
,
w
3
,
P
)
=
sup
⁡
{
|
w
1
|
,
y
(
t
)
v
(
t
)
(
|
w
2
|
+
|
w
3
|
)
}
.
The right-hand derivative of V
(
t
) along the positive solution of (43) is
(50)
D
+
|
w
1
|
≤
(
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
)
|
w
1
|
+
k
1
x
2
(
x
+
v
)
2
(
|
w
2
|
+
|
w
3
|
)
,
D
+
|
w
2
|
≤
a
|
w
1
|
+
(
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
d
3
)
|
w
2
|
+
p
|
w
3
|
D
+
|
w
3
|
≤
k
1
v
2
(
x
+
v
)
2
|
w
2
|
+
(
-
p
-
d
2
-
d
3
)
|
w
3
|
.
Therefore
(51)
D
+
y
v
(
|
w
2
|
+
|
w
3
|
)
=
y
˙
v
-
y
v
˙
v
2
(
|
w
2
|
+
|
w
3
|
)
+
y
v
D
+
(
|
w
2
|
+
|
w
3
|
)
≤
y
v
(
y
˙
y
-
v
˙
v
)
(
|
w
2
|
+
|
w
3
|
)
+
y
v
(
a
|
w
1
|
+
(
-
d
1
-
d
3
)
|
w
2
|
+
(
-
d
2
-
d
3
)
|
w
3
|
)
≤
a
y
v
|
w
1
|
+
y
v
(
|
w
2
|
+
|
w
3
|
)
×
(
y
˙
y
-
v
˙
v
-
d
3
-
min
⁡
(
d
1
,
d
2
)
)
,
(52)
D
+
V
(
t
)
≤
sup
⁡
{
g
1
(
t
)
,
g
2
(
t
)
}
V
(
t
)
,
where
(53)
g
1
(
t
)
=
-
d
1
-
p
-
d
2
-
k
1
v
2
(
x
+
v
)
2
+
k
1
v
x
2
y
(
x
+
v
)
2
≤
y
˙
y
-
d
1
,
g
2
(
t
)
=
a
y
v
+
y
˙
y
-
v
˙
v
-
d
3
-
min
⁡
(
d
1
,
d
2
)
=
y
˙
y
-
min
⁡
(
d
1
,
d
2
)
.
Denote d
=
min
⁡
(
d
1
,
d
2
), and then
(54)
sup
⁡
{
g
1
(
t
)
,
g
2
(
t
)
}
≤
y
˙
y
-
d
.
By (52) and Gronwall’s inequality, one obtains
(55)
V
(
t
)
≤
V
(
0
)
y
(
t
)
e
-
d
t
≤
V
(
0
)
e
-
d
t
λ
d
.
V
(
t
)
→
0 when t
→
∞, and then (
w
1
,
w
2
,
w
3
)
→
0 when t
→
∞ by (48). The second additive compound system is asymptotically stable. This verifies the condition (
3
) of Lemma 4.
(
4
) The Jacobi matrix of (3) at the endemic infection equilibrium Q
2 is
(56)
J
Q
2
=
[
-
d
1
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
p
-
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
-
d
2
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
0
a
-
d
3
]
,
and then
(57)
det
⁡
(
J
Q
2
)
=
|
-
d
1
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
p
-
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
-
d
2
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
0
a
-
d
3
|
=
-
[
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
]
(
d
2
+
p
)
d
3
-
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
4
+
d
3
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
+
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
[
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
]
=
-
d
1
(
d
2
+
p
)
d
3
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
p
)
d
3
-
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
4
+
d
3
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
+
a
d
1
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
+
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
4
=
d
1
[
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
-
(
d
2
+
p
)
d
3
]
-
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
+
d
3
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
=
a
d
1
k
1
[
(
x
¯
x
¯
+
v
¯
)
2
-
1
R
0
]
-
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
.
According to (7), v
¯
=
(
R
0
-
1
)
x
¯, and then one can obtain
(58)
det
⁡
(
J
Q
2
)
=
a
d
1
k
1
[
1
R
0
2
-
1
R
0
]
-
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
=
a
d
1
k
1
(
1
-
R
0
)
R
0
2
-
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
.
Since J
Q
2 is a 3
×
3 matrix, one gets n
=
3. Then
(59)
(
-
1
)
3
det
⁡
(
J
Q
2
)
=
-
a
d
1
k
1
(
1
-
R
0
)
R
0
2
+
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
.
If R
0
>
1, then (
-
1
)
3
det
⁡
(
J
Q
2
)
>
0 holds in Ω
0. This verifies condition (
4
) of Lemma 4.
Hence, ifR
0
>
1, then the endemic infection equilibrium point Q
2 is globally asymptotically stable in Ω
0 by Lemma 4.


## 4.1. Locally Asymptotical Stability of the Endemic Infection Equilibrium PointQ
2

Theorem 3.
IfR
0
>
1, then the endemic infection equilibrium point Q
2 of (3) is locally asymptotically stable.Proof.
Put the equilibrium pointQ
2 into matrix (16); then one obtains
(33)
J
Q
2
=
[
-
d
1
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
p
-
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
-
d
2
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
0
a
-
d
3
]
,
|
λ
E
-
J
Q
2
|
=
|
λ
+
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
λ
+
p
+
d
2
-
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
0
-
a
λ
+
d
3
|
=
0
.
Solving the eigenequation of the matrix above, here is(34)
λ
3
+
[
d
1
+
d
2
+
d
3
+
p
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
]
︸
a
1
λ
2
+
[
d
3
(
d
2
+
p
)
+
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
(
d
2
+
p
+
d
3
)
-
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
]
︸
a
2
λ
+
d
3
(
d
2
+
p
)
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
+
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
︸
a
3
=
0
.If R
0
>
1, one obtains that(35)
a
1
=
d
1
+
d
2
+
d
3
+
p
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
>
0
,
a
2
=
d
3
(
d
2
+
p
)
+
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
(
d
2
+
p
+
d
3
)
-
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
=
a
k
1
R
0
+
d
1
(
d
2
+
d
3
+
p
)
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
d
3
)
-
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
=
a
k
1
R
0
-
a
k
1
x
¯
2
R
0
2
x
¯
2
+
d
1
(
d
2
+
d
3
+
p
)
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
d
3
)
=
a
k
1
(
R
0
-
1
)
R
0
2
+
d
1
(
d
2
+
d
3
+
p
)
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
d
3
)
>
0
,
a
3
=
d
3
(
d
2
+
p
)
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
+
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
(
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
=
d
1
d
3
(
d
2
+
p
)
-
a
d
1
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
+
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
=
d
1
a
k
1
(
R
0
-
1
)
R
0
2
+
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
>
0
,
a
1
a
2
-
a
3
=
(
d
2
+
d
3
+
p
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
×
a
k
1
(
R
0
-
1
)
R
0
2
+
d
1
(
d
2
+
d
3
+
p
)
×
(
d
2
+
d
3
+
p
+
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
+
d
2
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
d
3
+
p
+
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
+
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
3
+
p
+
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
)
>
0
.
Hence all inequalities of the Routh-Hurwitz criterion are satisfied. Therefore, the endemic infection equilibrium pointQ
2 is locally asymptotically stable.


## 4.2. Globally Asymptotical Stability of the Endemic Infection Equilibrium PointQ
2

In this subsection, we firstly introduce a lemma outlined by Li and Wang [18], and then using this lemma discusses the globally asymptotical stability of the endemic infection equilibrium point Q
2 of (3).The lemma is briefly summarized as follows.Letx
↦
f
(
x
)
∈
R
n be a C
1 function for x in an open set Γ
⊂
R
n. Consider the differential system
(36)
x
˙
=
f
(
x
)
.
Denote by x
(
t
,
x
0
) the solution to (36) such that x
(
0
,
x
0
)
=
x
0. Let x
~ be an equilibrium point of (36). Li and Wang [18] made the following two basic assumptions:(
H
1
)
there exists a compact absorbing setK
⊂
Γ;
(
H
2
)
equation (36) has a unique equilibrium x
~ in Γ.Li and Wang (see Theorem  2.5 in [18]) have given the following lemma.Lemma 4 (see [18]).
Assume that(1)
assumptions (H
1) and (H
2) hold;
(2)
equation (36) satisfies the Poincaré-Bendixson Property;
(3)
for each periodic solutionx
=
p
(
t
) to (36) with p
(
0
)
∈
Γ, the linear system (the second additive compound system)
(37)
w
˙
(
t
)
=
∂
f
[
2
]
∂
x
(
P
(
t
)
)
w
(
t
)
is asymptotically stable, where∂
f
[
2
]
/
∂
x is the second additive compound matrix of the Jacobian matrix ∂
f
/
∂
x;
(4)
(
-
1
)
n
det
(
(
∂
f
/
∂
x
)
(
x
~
)
)
>
0.
Then the unique equilibriumx
~ is globally asymptotically stable in Γ.Now one uses Lemma4 to show the following.Theorem 5.
IfR
0
>
1, then the endemic infection equilibrium point Q
2 of (3) is globally asymptotically stable in Ω
0, where Ω
0 is defined by (15).Proof.
Based on Lemma4, the proof of Theorem 5 has been implemented via the following four steps.
(
1
) For epidemic models and many other biological models where the feasible region is a bounded cone, (
H
1
) is equivalent to the uniform persistence of the system [19]. By (15), Ω
0 is bounded, so it only needs to show the uniform persistence of (3). According to Proposition  3.3 in [20], the necessary and sufficient condition for the uniform persistence of (3) is equivalent to the equilibrium point Q
1 being unstable. Theorem 1 has shown that Q
1 is unstable if R
0
>
1. Therefore, (3) is uniformly persistent if R
0
>
1 so that (
H
1
) holds if R
0
>
1.
Meanwhile,Q
1
=
(
x
0
,
0,0
) by (4), so Q
1 does not exist in Ω
0. Hence, Q
2 is the unique equilibrium point of (3) in Ω
0 so that (
H
2
) holds.
The results above verify the condition(
1
) of Lemma 4.
(
2
) The Jacobian matrix of (3) is
(38)
J
(
x
,
y
,
v
)
=
[
-
d
1
-
a
1
p
-
a
2
a
1
-
p
-
d
2
a
2
0
a
-
d
3
]
,
where a
1
=
(
k
1
v
2
/
(
x
+
v
)
2
) and a
2
=
(
k
1
x
2
/
(
x
+
v
)
2
).
IfH
=
diag
⁡
(
1
,
-
1,1
), then
(39)
H
J
H
=
[
-
d
1
-
k
1
v
2
x
+
v
-
p
-
k
1
x
2
(
x
+
v
)
2
-
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
-
k
1
x
2
(
x
+
v
)
2
0
-
a
-
d
3
]
,
and one can obtain that H
J
H has nonpositive off-diagonal elements in Ω
0. Therefore (3) is competitive in Ω
0. It is known that 3-dimensional competitive systems have the Poincaré-Bendixson Property [21]. Hence, (3) satisfies the Poincaré-Bendixson Property. This verifies condition (
2
) of Lemma 4.
(
3
) Let P
(
t
)
=
(
x
(
t
)
,
y
(
t
)
,
v
(
t
)
) be a periodic solution in Ω
0.
According to [22], if B
=
(
b
i
j
) is a 3
×
3 matrix, then the second additive compound B
[
2
] of B is
(40)
B
[
2
]
=
[
b
11
+
b
22
b
23
-
b
13
b
32
b
11
+
b
33
b
12
-
b
31
b
21
b
22
+
b
33
]
.
The Jacobian matrix of (3) is
(41)
J
(
x
,
y
,
v
)
=
[
-
d
1
-
a
1
p
-
a
2
a
1
-
p
-
d
2
a
2
0
a
-
d
3
]
,
where a
1
=
k
1
v
2
/
(
x
+
v
)
2 and a
2
=
k
1
x
2
/
(
x
+
v
)
2.
And then the second additive compound matrix of the Jacobian matrix of (3) is given by(42)
J
[
2
]
=
[
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
k
1
x
2
(
x
+
v
)
2
k
1
x
2
(
x
+
v
)
2
a
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
d
3
p
0
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
-
d
3
]
, and the second additive compound system of (3) along the periodic solution P
(
t
)
=
(
x
(
t
)
,
y
(
t
)
,
v
(
t
)
) is
(43)
w
˙
1
=
(
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
)
w
1
+
k
1
x
2
(
x
+
v
)
2
w
2
+
k
1
x
2
(
x
+
v
)
2
w
3
,
w
˙
2
=
a
w
1
+
(
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
d
3
)
w
2
+
p
w
3
,
w
˙
3
=
k
1
v
2
(
x
+
v
)
2
w
2
+
(
-
p
-
d
2
-
d
3
)
w
3
.
Define a global Lyapunov function by
(44)
V
(
w
1
,
w
2
,
w
3
,
P
)
=
∥
(
w
1
,
y
(
t
)
v
(
t
)
w
2
,
y
(
t
)
v
(
t
)
w
3
)
∥
,
where ∥
·
∥ is the norm in set D defined by
(45)
∥
(
w
1
,
w
2
,
w
3
)
∥
=
sup
⁡
{
|
w
1
|
,
|
w
2
|
+
|
w
3
|
}
.
Suppose that the solution P
(
t
) is periodic of least period ω
>
0 and that P
(
0
)
∈
Ω
0. According to [23], (3) is uniformly persistent, if there exists a positive constant μ such that
(46)
lim
⁡
inf
⁡
t
→
∞
x
(
t
)
≥
μ
,
lim
⁡
inf
⁡
t
→
∞
y
(
t
)
≥
μ
,
lim
⁡
inf
⁡
t
→
∞
v
(
t
)
≥
μ
.
Step (
1
) has shown that (3) is uniformly persistent if R
0
>
1. Hence, there always exists a positive constant μ which satisfies (46). The orbit of P
(
t
) remains at a positive distance from the boundary of Ω by the uniform persistence, and one can obtain that
(47)
y
(
t
)
≥
μ
,
v
(
t
)
≥
μ
for
large
enough
t
.
Sincev
<
a
λ
/
(
d
d
3
) by (15),
(48)
V
(
w
1
,
w
2
,
w
3
,
P
)
≥
μ
d
d
3
a
λ
∥
(
w
1
,
w
2
,
w
3
)
∥
,
for all (
w
1
,
w
2
,
w
3
)
∈
R
3.
Along a solution(
w
1
,
w
2
,
w
3
) of (43), V
(
w
1
,
w
2
,
w
3
,
P
) becomes
(49)
V
(
w
1
,
w
2
,
w
3
,
P
)
=
sup
⁡
{
|
w
1
|
,
y
(
t
)
v
(
t
)
(
|
w
2
|
+
|
w
3
|
)
}
.
The right-hand derivative of V
(
t
) along the positive solution of (43) is
(50)
D
+
|
w
1
|
≤
(
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
p
-
d
2
)
|
w
1
|
+
k
1
x
2
(
x
+
v
)
2
(
|
w
2
|
+
|
w
3
|
)
,
D
+
|
w
2
|
≤
a
|
w
1
|
+
(
-
d
1
-
k
1
v
2
(
x
+
v
)
2
-
d
3
)
|
w
2
|
+
p
|
w
3
|
D
+
|
w
3
|
≤
k
1
v
2
(
x
+
v
)
2
|
w
2
|
+
(
-
p
-
d
2
-
d
3
)
|
w
3
|
.
Therefore
(51)
D
+
y
v
(
|
w
2
|
+
|
w
3
|
)
=
y
˙
v
-
y
v
˙
v
2
(
|
w
2
|
+
|
w
3
|
)
+
y
v
D
+
(
|
w
2
|
+
|
w
3
|
)
≤
y
v
(
y
˙
y
-
v
˙
v
)
(
|
w
2
|
+
|
w
3
|
)
+
y
v
(
a
|
w
1
|
+
(
-
d
1
-
d
3
)
|
w
2
|
+
(
-
d
2
-
d
3
)
|
w
3
|
)
≤
a
y
v
|
w
1
|
+
y
v
(
|
w
2
|
+
|
w
3
|
)
×
(
y
˙
y
-
v
˙
v
-
d
3
-
min
⁡
(
d
1
,
d
2
)
)
,
(52)
D
+
V
(
t
)
≤
sup
⁡
{
g
1
(
t
)
,
g
2
(
t
)
}
V
(
t
)
,
where
(53)
g
1
(
t
)
=
-
d
1
-
p
-
d
2
-
k
1
v
2
(
x
+
v
)
2
+
k
1
v
x
2
y
(
x
+
v
)
2
≤
y
˙
y
-
d
1
,
g
2
(
t
)
=
a
y
v
+
y
˙
y
-
v
˙
v
-
d
3
-
min
⁡
(
d
1
,
d
2
)
=
y
˙
y
-
min
⁡
(
d
1
,
d
2
)
.
Denote d
=
min
⁡
(
d
1
,
d
2
), and then
(54)
sup
⁡
{
g
1
(
t
)
,
g
2
(
t
)
}
≤
y
˙
y
-
d
.
By (52) and Gronwall’s inequality, one obtains
(55)
V
(
t
)
≤
V
(
0
)
y
(
t
)
e
-
d
t
≤
V
(
0
)
e
-
d
t
λ
d
.
V
(
t
)
→
0 when t
→
∞, and then (
w
1
,
w
2
,
w
3
)
→
0 when t
→
∞ by (48). The second additive compound system is asymptotically stable. This verifies the condition (
3
) of Lemma 4.
(
4
) The Jacobi matrix of (3) at the endemic infection equilibrium Q
2 is
(56)
J
Q
2
=
[
-
d
1
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
p
-
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
-
d
2
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
0
a
-
d
3
]
,
and then
(57)
det
⁡
(
J
Q
2
)
=
|
-
d
1
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
p
-
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
-
d
2
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
0
a
-
d
3
|
=
-
[
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
]
(
d
2
+
p
)
d
3
-
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
4
+
d
3
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
+
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
[
d
1
+
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
]
=
-
d
1
(
d
2
+
p
)
d
3
-
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
(
d
2
+
p
)
d
3
-
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
4
+
d
3
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
+
a
d
1
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
+
a
k
1
2
x
¯
2
v
¯
2
(
x
¯
+
v
¯
)
4
=
d
1
[
a
k
1
x
¯
2
(
x
¯
+
v
¯
)
2
-
(
d
2
+
p
)
d
3
]
-
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
-
p
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
+
d
3
p
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
=
a
d
1
k
1
[
(
x
¯
x
¯
+
v
¯
)
2
-
1
R
0
]
-
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
.
According to (7), v
¯
=
(
R
0
-
1
)
x
¯, and then one can obtain
(58)
det
⁡
(
J
Q
2
)
=
a
d
1
k
1
[
1
R
0
2
-
1
R
0
]
-
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
=
a
d
1
k
1
(
1
-
R
0
)
R
0
2
-
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
.
Since J
Q
2 is a 3
×
3 matrix, one gets n
=
3. Then
(59)
(
-
1
)
3
det
⁡
(
J
Q
2
)
=
-
a
d
1
k
1
(
1
-
R
0
)
R
0
2
+
d
2
d
3
k
1
v
¯
2
(
x
¯
+
v
¯
)
2
.
If R
0
>
1, then (
-
1
)
3
det
⁡
(
J
Q
2
)
>
0 holds in Ω
0. This verifies condition (
4
) of Lemma 4.
Hence, ifR
0
>
1, then the endemic infection equilibrium point Q
2 is globally asymptotically stable in Ω
0 by Lemma 4.


## 5. Numerical Simulation

In the first subsection, we determine some parameter values of an anti-HIV infection treatment model based on (3). In the second subsection, using the anti-HIV infection treatment model simulates the dynamics of the Group I’s anti-HIV infection treatment. In the third subsection, using the anti-HIV infection treatment model simulates the dynamics of the Group II’s anti-HIV infection treatment. In the fourth subsection, we make long-term predictions for the two groups’ anti-HIV infection treatment, respectively.


### 5.1. Modeling

Baxter et al. [24] have reported a randomized study of antiretroviral management based on plasma genotypic antiretroviral resistance testing in HIV patients failing therapy, which was enrolled from 14 units of the Terry Beirn Community Programs for Clinical Research on AIDS and the Walter Reed Army Medical Center (see the HIV drug resistance database of Stanford University [25]). These patients were failing virologically on a combination antiretroviral regimen containing protease inhibitors (PI) and nucleoside reverse transcriptase inhibitors (NRTI) [24]. The patients were seen at 4, 8, and 12 weeks. At each follow-up visit, changes in antiretroviral treatment were recorded and the tested items included patients’ plasma CD
4
+ T cells counts and plasma HIV-1 RNA levels by the Chiron 2.0 bDNA assay [24].In the following subsections, we select, from [24, 25], two group patients’ mean uninfected CD
4
+ T cells counts and mean HIV RNA levels to simulate and make long-term predictions for the patients’ treatment outcomes. Group I consists of 15 patients. Group II consists of 13 patients. The two groups of patients received the same PI: ritonavir (RTV) and saquinavir (SQV). Additionally, Group I received NRTI: strvudine (D4T). Group II received NRTI: strvudine (D4T) and dideoxyinosine (DDI) [24, 25].Based on (3), the anti-HIV infection treatment model has the form
(60)
x
˙
=
λ
-
d
1
x
-
(
1
-
m
)
k
1
v
x
x
+
v
+
p
y
,
y
˙
=
(
1
-
m
)
k
1
v
x
x
+
v
-
d
2
y
-
p
y
,
v
˙
=
(
1
-
n
)
a
y
-
d
3
v
,
where m
,
n
(
0
≤
m
,
n
≤
1
) are the efficacy variables of the treatment.The infection-free equilibrium pointQ
1 of (60) is the same as that defined by (4):
(61)
Q
1
=
(
λ
d
1
,
0,0
)
.The endemic infection equilibrium pointQ
2 of (60) is given by
(62)
Q
2
=
(
x
¯
,
y
¯
,
v
¯
)
,
where,
(63)
x
¯
=
λ
R
0
(
1
-
m
)
k
1
(
R
0

-

1
)

+

d
1
R
0

-

p
(
R
0

-

1
)
R
0
(
d
3
/
(
1

-

n
)
a
)
,
y
¯

=

d
3
(
1
-
n
)
a
(
R
0
-
1
)
x
¯
,
v
¯

=

(
1
-
n
)
a
y
¯
d
3

=

(
R
0
-
1
)
x
¯
,
and R
0 is the basic virus reproductive number of (60):
(64)
R
0
=
(
1
-
n
)
(
1
-
m
)
a
k
1
d
3
(
d
2
+
p
)
.Determine the parameter value ranges of (60).(1)
NaiveCD
4
+ T cells decayed with an average half-life of 50 days [26]. Therefore, one obtains
(65)
d
1
=
-
ln
⁡
(
0.5
)
50
.
Henced
1
≈
0.0139. Because the apoptosis of CD
4
+ T cells is raised by HIV infection [27, 28], d
1 should be more than 0.0139 during the simulation.
(2)
Since healthy individuals have an average of 830/μL CD
4
+ T cells [29],
(66)
λ
=
d
1
×
830
≈
0.0139
×
830
=
11.5370
.
(3)
Since the cells that produce the virus are also short-lived, with a half-life of approximately 1.2 days [30], one obtains
(67)
d
2
=
-
ln
⁡
(
0.5
)
1.2
≈
0.5776
.
(4)
Since the half-life of HIV-1 in the plasma appears to be only 1 to 2 days [30], one selects
(68)
d
3
=
-
ln
⁡
(
0.5
)
1.5
≈
0.4621
.
(5)
Because only a small fraction of cells in the eclipse phase will revert to the uninfected state, it assumes thatp
=
0.01 [12]. Hence one obtains
(69)
p
=
0.01
.
(6)
According to reference [31], one can determine the other parameter value ranges as follows:
(70)
k
1
∈
[
2.5
×
1
0
-
5
,
0.5
]
,
a
∈
[
2,1250
]
,
m
∈
[
0,1
]
,
n
∈
[
0,1
]
.In each group, there was one patient whose clinical data was not complete. Therefore we do not conclude the two patients’ clinical data in the following simulations.


### 5.2. The Mean Dynamics Simulation of Group I’s Anti-HIV Infection Treatment

Using the other 14 patients’ clinical data determines the equation parameter values as follows:(71)
k
1
=
6
×
1
0
-
2
,
a
=
37
.The value changes of the parametersd
1, m, and n and the basic virus reproductive number R
0 are listed in Table 1.Table 1
Parameter values andR
0 at different weeks.
Weeks
d
1
m
n
R
0
0~4
0.033
0.53
0.84
0.6148
4~12
0.038
0.45
0.76
1.0792The simulation results of the mean dynamics of anti-HIV infection treatment of Group I are shown in Figure1. During the first 4 weeks, the treatment reduced the basic virus reproductive number R
0 from 8.1759 to 0.6148. Hence the patients’ mean HIV RNA levels decreased rapidly to approach infection-free steady state Q
1 defined by (4) as Theorem 2 predicts.Outcomes of the treatment efficacy of Group I. Circles: the clinical data; solid line: the numerical simulation of (60). (a) Mean uninfected CD
4
+ T cells counts. (b) Mean HIV RNA levels.
(a)
(b)However, after the 4th week, the resistance to antiretroviral drugs appeared. It made the therapy efficacy parametersm and n decrease from 0.53 to 0.45 and 0.84 to 0.76, respectively. Meanwhile the apoptosis of CD
4
+ T cells was raised by HIV more strongly. Hence d
1 rose from 0.033 to 0.038. The suboptimal treatment increased R
0 value of Group I from 0.6148 to 1.0792. As a result, the patients’ mean HIV RNA levels increased slowly to converge to a new infected steady state Q
2 defined by (62) and (63) as Theorem 5 predicts. On the other hand, observe that the mean uninfected CD
4
+ T cell counts of Group I increased rapidly in the first 4 weeks and decreased slowly in the following 8 weeks.


### 5.3. The Mean Dynamics Simulation of Group II’s Anti-HIV Infection Treatment

Using 12 patients’ clinical data determines the equation parameter values as follows:(72)
k
1
=
6
×
1
0
-
2
,
a
=
26
.The value changes of the parametersd
1, m, and n and the basic virus reproductive number R
0 are listed are listed in Table 2.Table 2
Parameter values andR
0 at different weeks.
Weeks
d
1
m
n
R
0
0~4
0.041
0.65
0.68
0.6435
4~12
0.0428
0.55
0.59
1.0600The simulation results of the mean dynamics of anti-HIV infection treatment of Group II are shown in Figure2. During the first 4 weeks, the treatment reduced the basic virus reproductive number R
0 from 5.7452 to 0.6435. Hence the patients’ mean HIV RNA levels decreased rapidly to approach to infection-free steady state Q
1 defined by (4) as Theorem 2 predicts.Outcomes of the treatment efficacy of Group II. Circles: the clinical data; solid line: the numerical simulation of (60). (a) Mean uninfected CD
4
+ T cells counts. (b) Mean HIV RNA levels.
(a)
(b)However, after the 4th week, the resistance to antiretroviral drugs appeared. It made the therapy efficacy parametersm and n decrease from 0.65 to 0.55 and 0.68 to 0.59, respectively. Meanwhile the apoptosis of CD
4
+ T cells was raised by HIV more strongly. Hence d
1 rose from 0.041 to 0.0428. The suboptimal treatment made R
0 value of Group II increase from 0.6435 to 1.0600. As a result, the patients’ mean HIV RNA levels increased slowly to converge to a new infected steady state Q
2 defined by (62) and (63) as Theorem 5 predicts. On the other hand, observe that the mean uninfected CD
4
+ T cell counts of Group II increased rapidly in the first 4 weeks but rose slowly in the following 8 weeks.


### 5.4. The Long-Term Predictions for the Two Groups’ Anti-HIV Infection Treatment

According to 2013 HIV therapy guidelines published by World Health Organization (WHO) [32], viral load is recommended as the preferred monitoring approach to diagnose and confirm antiretroviral treatment failure; treatment failure is defined by a persistently detectable viral load exceeding 1000 copies/mL after at least six months of using antiretroviral drugs. However, HIV RNA levels of the two groups of patients were only tested at 4, 8, and 12 weeks in the study [24]. Therefore, it is necessary to make long-term predictions to detect whether the treatments for the two groups are failure.Assume that after ending the 12 weeks’ antiretroviral treatment testing, the two groups keep receiving the anti-HIV infection treatment for 2 years. During the 2 years, the drug resistance does not become worse and all parameter values do not change. Using the numerical simulation of Equation (60) makes the long-term predictions for the two groups’ anti-HIV infection treatment. The long-term prediction for Group I’ anti-HIV infection treatment is shown in Figure3. Two years’ outcomes of the therapy for Group II are shown in Figure 4. Observe that after finishing the 12 weeks’ antiretroviral treatment testing, the mean uninfected CD
4
+ T cells counts of the two groups both decline a little and finally keep at a level larger than the mean baseline values; the mean HIV RNA levels of the two groups both rise a lot to a level less than the mean baseline values but keep exceeding 1000 copies/mL all the time. The long-term predictions suggest that the treatments for the two groups are failure and better anti-HIV infection therapies should be considered.The long-term prediction for the treatment efficacy of Group I. Circles: the clinical data; solid line: the numerical simulation of (60). (a) Mean uninfected CD
4
+ T cells counts. (b) Mean HIV RNA levels.
(a)
(b)The long-term prediction for the treatment efficacy of Group II. Circles: the clinical data; solid line: the numerical simulation of (60). (a) Mean uninfected CD
4
+ T cells counts. (b) Mean HIV RNA levels.
(a)
(b)


## 5.1. Modeling

Baxter et al. [24] have reported a randomized study of antiretroviral management based on plasma genotypic antiretroviral resistance testing in HIV patients failing therapy, which was enrolled from 14 units of the Terry Beirn Community Programs for Clinical Research on AIDS and the Walter Reed Army Medical Center (see the HIV drug resistance database of Stanford University [25]). These patients were failing virologically on a combination antiretroviral regimen containing protease inhibitors (PI) and nucleoside reverse transcriptase inhibitors (NRTI) [24]. The patients were seen at 4, 8, and 12 weeks. At each follow-up visit, changes in antiretroviral treatment were recorded and the tested items included patients’ plasma CD
4
+ T cells counts and plasma HIV-1 RNA levels by the Chiron 2.0 bDNA assay [24].In the following subsections, we select, from [24, 25], two group patients’ mean uninfected CD
4
+ T cells counts and mean HIV RNA levels to simulate and make long-term predictions for the patients’ treatment outcomes. Group I consists of 15 patients. Group II consists of 13 patients. The two groups of patients received the same PI: ritonavir (RTV) and saquinavir (SQV). Additionally, Group I received NRTI: strvudine (D4T). Group II received NRTI: strvudine (D4T) and dideoxyinosine (DDI) [24, 25].Based on (3), the anti-HIV infection treatment model has the form
(60)
x
˙
=
λ
-
d
1
x
-
(
1
-
m
)
k
1
v
x
x
+
v
+
p
y
,
y
˙
=
(
1
-
m
)
k
1
v
x
x
+
v
-
d
2
y
-
p
y
,
v
˙
=
(
1
-
n
)
a
y
-
d
3
v
,
where m
,
n
(
0
≤
m
,
n
≤
1
) are the efficacy variables of the treatment.The infection-free equilibrium pointQ
1 of (60) is the same as that defined by (4):
(61)
Q
1
=
(
λ
d
1
,
0,0
)
.The endemic infection equilibrium pointQ
2 of (60) is given by
(62)
Q
2
=
(
x
¯
,
y
¯
,
v
¯
)
,
where,
(63)
x
¯
=
λ
R
0
(
1
-
m
)
k
1
(
R
0

-

1
)

+

d
1
R
0

-

p
(
R
0

-

1
)
R
0
(
d
3
/
(
1

-

n
)
a
)
,
y
¯

=

d
3
(
1
-
n
)
a
(
R
0
-
1
)
x
¯
,
v
¯

=

(
1
-
n
)
a
y
¯
d
3

=

(
R
0
-
1
)
x
¯
,
and R
0 is the basic virus reproductive number of (60):
(64)
R
0
=
(
1
-
n
)
(
1
-
m
)
a
k
1
d
3
(
d
2
+
p
)
.Determine the parameter value ranges of (60).(1)
NaiveCD
4
+ T cells decayed with an average half-life of 50 days [26]. Therefore, one obtains
(65)
d
1
=
-
ln
⁡
(
0.5
)
50
.
Henced
1
≈
0.0139. Because the apoptosis of CD
4
+ T cells is raised by HIV infection [27, 28], d
1 should be more than 0.0139 during the simulation.
(2)
Since healthy individuals have an average of 830/μL CD
4
+ T cells [29],
(66)
λ
=
d
1
×
830
≈
0.0139
×
830
=
11.5370
.
(3)
Since the cells that produce the virus are also short-lived, with a half-life of approximately 1.2 days [30], one obtains
(67)
d
2
=
-
ln
⁡
(
0.5
)
1.2
≈
0.5776
.
(4)
Since the half-life of HIV-1 in the plasma appears to be only 1 to 2 days [30], one selects
(68)
d
3
=
-
ln
⁡
(
0.5
)
1.5
≈
0.4621
.
(5)
Because only a small fraction of cells in the eclipse phase will revert to the uninfected state, it assumes thatp
=
0.01 [12]. Hence one obtains
(69)
p
=
0.01
.
(6)
According to reference [31], one can determine the other parameter value ranges as follows:
(70)
k
1
∈
[
2.5
×
1
0
-
5
,
0.5
]
,
a
∈
[
2,1250
]
,
m
∈
[
0,1
]
,
n
∈
[
0,1
]
.In each group, there was one patient whose clinical data was not complete. Therefore we do not conclude the two patients’ clinical data in the following simulations.


## 5.2. The Mean Dynamics Simulation of Group I’s Anti-HIV Infection Treatment

Using the other 14 patients’ clinical data determines the equation parameter values as follows:(71)
k
1
=
6
×
1
0
-
2
,
a
=
37
.The value changes of the parametersd
1, m, and n and the basic virus reproductive number R
0 are listed in Table 1.Table 1
Parameter values andR
0 at different weeks.
Weeks
d
1
m
n
R
0
0~4
0.033
0.53
0.84
0.6148
4~12
0.038
0.45
0.76
1.0792The simulation results of the mean dynamics of anti-HIV infection treatment of Group I are shown in Figure1. During the first 4 weeks, the treatment reduced the basic virus reproductive number R
0 from 8.1759 to 0.6148. Hence the patients’ mean HIV RNA levels decreased rapidly to approach infection-free steady state Q
1 defined by (4) as Theorem 2 predicts.Outcomes of the treatment efficacy of Group I. Circles: the clinical data; solid line: the numerical simulation of (60). (a) Mean uninfected CD
4
+ T cells counts. (b) Mean HIV RNA levels.
(a)
(b)However, after the 4th week, the resistance to antiretroviral drugs appeared. It made the therapy efficacy parametersm and n decrease from 0.53 to 0.45 and 0.84 to 0.76, respectively. Meanwhile the apoptosis of CD
4
+ T cells was raised by HIV more strongly. Hence d
1 rose from 0.033 to 0.038. The suboptimal treatment increased R
0 value of Group I from 0.6148 to 1.0792. As a result, the patients’ mean HIV RNA levels increased slowly to converge to a new infected steady state Q
2 defined by (62) and (63) as Theorem 5 predicts. On the other hand, observe that the mean uninfected CD
4
+ T cell counts of Group I increased rapidly in the first 4 weeks and decreased slowly in the following 8 weeks.


## 5.3. The Mean Dynamics Simulation of Group II’s Anti-HIV Infection Treatment

Using 12 patients’ clinical data determines the equation parameter values as follows:(72)
k
1
=
6
×
1
0
-
2
,
a
=
26
.The value changes of the parametersd
1, m, and n and the basic virus reproductive number R
0 are listed are listed in Table 2.Table 2
Parameter values andR
0 at different weeks.
Weeks
d
1
m
n
R
0
0~4
0.041
0.65
0.68
0.6435
4~12
0.0428
0.55
0.59
1.0600The simulation results of the mean dynamics of anti-HIV infection treatment of Group II are shown in Figure2. During the first 4 weeks, the treatment reduced the basic virus reproductive number R
0 from 5.7452 to 0.6435. Hence the patients’ mean HIV RNA levels decreased rapidly to approach to infection-free steady state Q
1 defined by (4) as Theorem 2 predicts.Outcomes of the treatment efficacy of Group II. Circles: the clinical data; solid line: the numerical simulation of (60). (a) Mean uninfected CD
4
+ T cells counts. (b) Mean HIV RNA levels.
(a)
(b)However, after the 4th week, the resistance to antiretroviral drugs appeared. It made the therapy efficacy parametersm and n decrease from 0.65 to 0.55 and 0.68 to 0.59, respectively. Meanwhile the apoptosis of CD
4
+ T cells was raised by HIV more strongly. Hence d
1 rose from 0.041 to 0.0428. The suboptimal treatment made R
0 value of Group II increase from 0.6435 to 1.0600. As a result, the patients’ mean HIV RNA levels increased slowly to converge to a new infected steady state Q
2 defined by (62) and (63) as Theorem 5 predicts. On the other hand, observe that the mean uninfected CD
4
+ T cell counts of Group II increased rapidly in the first 4 weeks but rose slowly in the following 8 weeks.


## 5.4. The Long-Term Predictions for the Two Groups’ Anti-HIV Infection Treatment

According to 2013 HIV therapy guidelines published by World Health Organization (WHO) [32], viral load is recommended as the preferred monitoring approach to diagnose and confirm antiretroviral treatment failure; treatment failure is defined by a persistently detectable viral load exceeding 1000 copies/mL after at least six months of using antiretroviral drugs. However, HIV RNA levels of the two groups of patients were only tested at 4, 8, and 12 weeks in the study [24]. Therefore, it is necessary to make long-term predictions to detect whether the treatments for the two groups are failure.Assume that after ending the 12 weeks’ antiretroviral treatment testing, the two groups keep receiving the anti-HIV infection treatment for 2 years. During the 2 years, the drug resistance does not become worse and all parameter values do not change. Using the numerical simulation of Equation (60) makes the long-term predictions for the two groups’ anti-HIV infection treatment. The long-term prediction for Group I’ anti-HIV infection treatment is shown in Figure3. Two years’ outcomes of the therapy for Group II are shown in Figure 4. Observe that after finishing the 12 weeks’ antiretroviral treatment testing, the mean uninfected CD
4
+ T cells counts of the two groups both decline a little and finally keep at a level larger than the mean baseline values; the mean HIV RNA levels of the two groups both rise a lot to a level less than the mean baseline values but keep exceeding 1000 copies/mL all the time. The long-term predictions suggest that the treatments for the two groups are failure and better anti-HIV infection therapies should be considered.The long-term prediction for the treatment efficacy of Group I. Circles: the clinical data; solid line: the numerical simulation of (60). (a) Mean uninfected CD
4
+ T cells counts. (b) Mean HIV RNA levels.
(a)
(b)The long-term prediction for the treatment efficacy of Group II. Circles: the clinical data; solid line: the numerical simulation of (60). (a) Mean uninfected CD
4
+ T cells counts. (b) Mean HIV RNA levels.
(a)
(b)


## 6. Conclusion

This paper introduces a modified HIV infection differential equation model with a saturated infection ratek
1
x
v
/
(
x
+
v
) and the proportion of infected cells reverting to the uninfected state.The basic virus reproductive numberR
0 of the model is independent of a patient’s plasma total CD
4
+ T cell counts λ
/
d
1, and the actual incidence rate is not linear over the entire range of virus v
(
t
) and uninfected CD
4
+ T cells x
(
t
) any more. This suggests that our model is more reasonable than the model proposed by [13].The modified model has two equilibrium points: infection-free equilibrium pointQ
1 and endemic infection equilibrium point Q
2. This paper discusses the locally asymptotical stabilities and globally asymptotical stabilities of the two equilibrium points, simulates the dynamics of two group patients’ anti-HIV infection treatment, and makes long-term predictions for the two groups’ anti-HIV infection treatment.The theoretical results suggest the following.(1)
If the basic virus reproductive numberR
0
<
1, then the infection-free equilibrium point Q
1 of (3) is globally asymptotically stable in Ω. This means that if a person with R
0
<
1, the person can recover automatically even if infected with a large amount of HIV; if a treatment makes a patient’s R
0
<
1, the patient will be cured eventually even if infected with a large amount of HIV.
The recent reports on three HIV infected patients have shown that some HIV infected patients may be cured via bone marrow transplants. The Berlin Patient was the first person cured of HIV [33]. After the Berlin Man, two cases reported cured of HIV in Kenya [34]. These reports can make one postulate that most individuals who connect HIV virus will not be infected by it and are not infecting other people. They will recover automatically without any treatment. The fact has not been well recognized since AIDS has been discovered in 1983. Mathematically, such phenomena can be described also via (3) where k
1
=
a
=
0, and thus R
0
=
0.
(2)
If the basic virus reproductive numberR
0
>
1, then the endemic infection equilibrium point Q
2 of (3) is globally asymptotically stable in the interior of Ω. This means that a person with R
0
>
1 will have endemic infection even if infected with only one HIV; if a treatment cannot make a patient’s R
0
<
1, the patient’s HIV RNA in vivo cannot be cleared up eventually.Based on the simulation results, one can propose the following hypotheses.(1)
For a poor HIV treatment response patient, the drug resistance appears when the patient’s HIV RNA level reduces to the first lowest level.
This hypothesis may interpret why the two group patients’ meanCD
4
+ T cells counts rose and mean HIV RNA levels declined rapidly in the first 4 weeks but contrary in the following weeks (see Tables 1 and 2 and Figures 1 and 2).
(2)
Once a patient’s drug resistance appears, the patient’s HIV in vivo promotes the apoptosis ofCD
4
+ T cells more strongly.
This hypothesis may interpret why between 4th week and 8th week the mean HIV RNA levels of the two groups kept a lower level than the first 4 weeks, but Group I’s meanCD
4
+ T cells counts started to decrease slowly and Group II’s mean CD
4
+ T cells counts rose more slowly than before (see Tables 1 and 2 and Figures 1 and 2).
(3)
According to the 2013 HIV therapy guidelines published by WHO [32], treatment failure is defined by a persistently detectable viral load exceeding 1000 copies/mL after at least six months of using antiretroviral drugs. Our long-term numerical simulation predictions suggest that after ending the 12 weeks’ antiretroviral treatment [24], the two group patients’ mean HIV RNA levels keep exceeding 1000 copies/mL during the additional 2 years’ anti-HIV infection treatment. The treatments for the two groups are failure and better anti-HIV infection therapies should be considered. This means that the additional 2 years’ treatments are not able to make patients obtain better outcomes. This may interpret why WHO defines that half a year’s treatment cannot suppress a patient’s HIV level below 1000 copies/mL to be treatment failure.


## Highlights

(1)
This paper introduces a modified HIV infection (anti-HIV infection therapy) differential equation model with a saturated infection rate. The basic virus reproductive numberR
0 of the model is independent of a patient’s plasma total CD
4
+ T cell counts λ
/
d
1. This suggests that our model is more reasonable than the model proposed by Srivastava and Chandra in 2010 [13] whose basic virus reproductive number R
0 is dependent on λ
/
d
1 which follows that the more CD
4
+ T cell counts an individual has, the more easily the individual is infected by HIV.
(2)
This paper proposes and proves two theorems (Theorems2 and 5) on the globally asymptotical stabilities of the infection-free equilibrium point Q
1 and the endemic infection equilibrium point Q
2 of the modified model.
(3)
This paper points out the implications of the two theorems which are ignored by other similar researches on modelling HIV infection (anti-HIV infection therapy):
(a)
A person with the basic virus reproductive numberR
0
<
1 will recover automatically even if the person is infected with a large amount of HIV. If a treatment makes a patient’s R
0
<
1, this patient will be cured eventually even if infected with a large amount of HIV.
(b)
A person withR
0
>
1 will have endemic infection even if the person is infected with only one HIV. If a treatment cannot make a patient’s R
0
<
1, the patient’s HIV RNA in vivo cannot be cleared up eventually.The recent reports on three HIV infected patients show that HIV infected patients may be cured via bone marrow transplants (e.g., see: Berlin Patient, first person cured of HIV, may soon have company, Los Angeles Times, July 27, 2012; after Berlin Man, two reported cured of HIV in Kenya, Africa Review, May 6, 2013).These reports can make one postulate that most individuals who connect HIV virus will not be infected by it and are not infecting other people, who will recover automatically without any treatment. In this case,R
0
=
0 where k
1
=
a
=
0. The fact has not been recognized since AIDS has been discovered in 1983.In the report [35], a small proportion of human immunodeficiency virus type 1 (HIV-1) infected individuals, called elite and viremic controllers, spontaneously control plasma HIV RNA levels to undetectable (elite controller) or <2000 copies/mL (viremic controller) in the absence of antiretroviral therapy.These phenomena can be interpreted by our Theorem2: HIV infected people’s basic virus reproductive number R
0
<
1.(4)
Based on the simulation results, we can propose the following hypotheses:
(a)
for a poor HIV treatment response patient, the drug resistance appears when the patient’s HIV RNA level reduces to the first lowest level;
(b)
Once a patient’s drug resistance appears, the patient’s HIV promotes the apoptosis ofCD
4
+ T cells more strongly;
(c)
According to the 2013 HIV therapy guidelines published by WHO [32], treatment failure is defined by a persistently detectable viral load exceeding 1000 copies/mL after at least six months of using antiretroviral drugs. Our long-term numerical simulation predictions suggest that after ending the 12 weeks’ antiretroviral treatment [24], the two group patients’ mean HIV RNA levels keep exceeding 1000 copies/mL during the additional 2 years’ anti-HIV infection treatment. The treatments for the two groups are failure and better anti-HIV infection therapies should be considered. This means that the additional 2 years’ treatments are not able to make patients obtain better outcomes. This may interpret why WHO defines that half a year’s treatment cannot suppress a patient’s HIV level below 1000 copies/mL to be treatment failure.



---

*Source: 145162-2014-03-23.xml*