www.biomathforum.org/biomath/index.php/biomath

ORIGINAL ARTICLE

The impact of infected T lymphocyte burst
rate and viral shedding rate on optimal

treatment scheduling in a human
immunodeficiency virus infection

Anuraag Bukkuri
Moffitt Cancer Center

Anuraag.Bukkuri@moffitt.org

Received: 27 April 2020, accepted: 17 August 2020, published: 12 September 2020

Abstract—We consider a mathematical model of
human immunodeficiency virus (HIV) infection dy-
namics of T lymphocyte (T cell), infected T cell,
and viral populations under reverse transcriptase
inhibitor (RTI) and protease inhibitor (PI) treat-
ment. Existence, uniqueness, and characterization
of optimal treatment profiles which minimize total
amount of drug used, viral, and infected T cell
populations, while maximizing levels of T cells are
determined analytically. Numerical optimal control
experiments are also performed to illustrate how
burst rate of infected T cells and shedding rate of
virions impact optimal treatment profiles. Finally, a
sensitivity analysis is performed to detect how model
input parameters contribute to output variance.

Keywords-Optimal Control, HIV, Reverse Tran-
scriptase Inhibitors, Protease Inhibitors

MSC: 92C50, 93C15

I. INTRODUCTION

HIV, also known as the human immunodefi-
ciency virus, is a virus which impairs the immune
system by entering vital cells (e.g. dendritic cells,
microglial cells, CD4+ T cells, and macrophages)
through CD4 and coreceptors on the cell mem-
brane and destroying them [Cunningham et al.
2010; Cenker et al. 2017]. Once the CD4+ levels
drop to a critical level, cell-mediated immunity
is lost and viral load increases, which leads to
the development of acquired immunodeficiency
syndrome (AIDS), the final stage of HIV infec-
tion. When this occurs, due to low immunity,
the patient’s body becomes a breeding ground for
opportunistic infection and AIDS defining cancers
such as Kaposi sarcoma, non-Hodgkin lymphoma
and cervical cancers. Without treatment, average
survival time after HIV infection is 9-11 years
[UNAIDS 2007].

Pioneering research done by scientists in the

Copyright: c© 2020 Bukkuri. This article is distributed under the terms of the Creative Commons Attribution License (CC BY
4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source
are credited.
Citation: Anuraag Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on
optimal treatment scheduling in a human immunodeficiency virus infection, Biomath 9 (2020),
2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 1 of 12

http://www.biomathforum.org/biomath/index.php/biomath
https://creativecommons.org/licenses/by/4.0/
https://creativecommons.org/licenses/by/4.0/
http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

Fig. 1. Cellular Structures and processes involved in cell-to-cell transmission of HIV. Panels A-G show different pathways
for HIV (green dots) cell-cell transmission between donor (green) and target (pink) cells.

1990s determined that direct cell-to-cell contact
is an important, and perhaps the predominant,
contributor to the propagation of HIV within a
host. Specifically, HIV spreads by entering im-
mune cells and subverting intercellular communi-
cation to enable their own viral replication through
structures such as tunneling nanotubes and filopo-
dia. Immunological synapses and phagocytosis of
infected cells are also hijacked by HIV and used as
gateways to infect target cells. Furthermore, HIV
is able to elicit fusion between infected donor and
target cells, forming horrific infected cell clusters
with a very high capacity for viral reproduction
and survival. There are many mechanisms by
which cell-cell transfer of HIV occurs, some of
which are summarized in Figure 1 from [Bracq et
al. 2018]. For more details on how these mecha-
nisms work, refer to [Bracq et al. 2018].

The main treatments for HIV are a class of
drugs called antiretrovirals (ART). Often, these
drugs are given in combination to suppress viral
replication and reduce plasma HIV viral load
in a treatment called highly active antiretroviral
therapy (HAART) [Autran et al. 1997; Komanduri

et al. 1998; Lederman et al. 1998]. Two common
such ARTs often used in HAART are reverse tran-
scriptase inhibitors and protease inhibitors [Arts et
al. 2012].

Reverse transcriptase inhibitors work to ter-
minate transcription, thus inhibiting viral repli-
cation. They do this in one of two ways: nu-
cleoside/nucleotide reverse transcriptase inhibitors
incorporate into the nascent DNA strand during
reverse transcription, whereas non-nucleoside re-
verse transcriptase inhibitors bind to non-catalytic
enzyme sites and is usually mediated through
steric hindrance that impedes structural changes
in HIV reverse transcriptase [Gulnik et al. 1995;
Gotte et al. 2000]. Protease inhibitors are usu-
ally substrate-based inhibitors which are designed
specifically against the viral protease based on its
crystal structure. They act on this viral protease,
inhibiting maturation of new viral particles. By
doing this, they attack the already formed HIV
before the next cycle of infection begins [Michaud
et al. 2012].

Figure 2, reproduced from [Michaud et al.
2012], depicts the HIV life cycle with various

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 2 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

Fig. 2. Depiction of HIV life cycle with entry, protease, RT, and integrase inhibitor drug targets

antiretroviral drug targets to more clearly elucidate
how they work.

Currently, physicians primarily use three mark-
ers to assess the stage of the disease and create
treatment strategies: CD4 T cell count, viral load,
and drug resistance [PENTA 2004]. Through opti-
mal control and sensitivity analyses, we propose
the inclusion of two more criteria into this as-
sessment: viral shedding rate and infected T cell
burst rate. Medically, it is possible to measure
the rate of viral shedding using polymerase chain
reaction (PCR) from collected samples as done
in [Tronstein et al. 2011]. Burst rate of infected
T cells can be measured using standard burst
rate analysis techniques, such as through viral
inhibitors, washouts of infected cells [Dimitrov et
al. 1993; Eckstein et al. 2001], quantitative image
analyses with in situ hybridization, quantitative
competitive RT-PCR of bulk tissue or single cells
at limiting dilution [Chun et al. 1997; Haase et
al. 1996; Hockett et al. 1999; Chen et al. 2007;
Hyman et al. 2009].

In this paper, we first summarize a model of
HIV infection with RTI and PI treatment devel-
oped by [Mobisa et al. 2018]. Then, in section
3, a detailed optimal control analysis will be

performed, including existence, uniqueness, and
analytical characterization results. In section 4,
parameter estimates are provided and several nu-
merical simulations were performed, specifically
analyzing the impact of changes to the burst rate
of infected cells and shedding rate of virions on
optimal treatment profiles. In section 5, a sensitiv-
ity analysis is performed to quantify the amount of
impact each parameter value has on healthy T cell,
infected T cell, and viral population dynamics.

II. MODEL DESCRIPTION

The model for HIV infection dynamics among
T cells (T), infected T cells (U), and free virus
particles (V) under RTI (u) and PI (v) treatment
is taken from a model created by [Mobisa et al.
2018] and is reproduced below:

dT

dt
= rT

(
1 −

T

Tmax

)
− (1 −u)βV T

−(1 −v)αTU −µT
dU

dt
= (1 −u)βV T + (1 −v)αTU −kU

dV

dt
= wkU − cV

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 3 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

In this model, the T cells are produced at rate
r and die naturally at rate µ; the growth of the T
cells is bounded with the enforcement of a carrying
capacity. The healthy T cells become infected by
the free virus and actively infected T cells via the
mass infection terms βV T and αTU, respectively.
The infected T cells die naturally at a rate of k.
Furthermore, these cells produce free viruses, as
captured by the first term in the viral equation,
which are then removed from circulation at a rate
c per virus. Note that the viral decline depends
on the amount and efficacy of the RTI and PI
treatments.

III. OPTIMAL CONTROL ANALYSIS
To determine best treatments strategies for HIV

infection, we formulate our problem as an optimal
control one. We aim to maximize the number of
healthy T cells, while minimizing the number of
infected T cells, free viruses, and drugs u and v
used over a fixed therapy horizon [0,T].

Our control class are piecewise continuous func-
tions defined for all t such that 0 ≤ u(t) ≤ 1 where
u(t) = 1 represents maximal RTI treatment and
u(t) = 0 represents no RTI treatment. Similarly,
v(t) = 0 represents no PI treatment and v(t) = 1
represents maximal PI treatment. We can then
describe the class of admissible controls as

U(t) = u(t),v(t) piecewise continuous s.t.
0 ≤ u(t),v(t) ≤ 1,∀t ∈ [0,T].

Now, we define our objective functional and op-
timal control problem. For a fixed therapy interval
[0,T], maximize the objective functional

J(u,v) =

∫ T
0
aT−bU−fV−

1

2
φ1u

2−
1

2
φ2v

2dt (1)

over all Lebesgue-measurable functions

u : [0,T] → [0,umax], v : [0,T] → [0,vmax],
subject to the above ODE dynamics and initial
conditions of T(0) = 1000, U(0) = 1, V(0) = 0.1.
Note that the objective functional is quadratic,
instead of linear, in u and v to ensure that the
resulting optimal controls can take on a continuum
of values from 0 to 1, rather than being binarily
confined to either 0 or 1.

A. Existence of Optimal Control

Using the theory developed by [Fleming et
al. 1975], we can determine the existence of an
optimal control for our state system. Specifically,
boundedness of solutions of the system for finite
time is required for existence and uniqueness of an
optimal control. Using the technique of supersolu-
tions and keeping in mind that T ≤ Tmax, upper
bounds on the solutions of the state system are
found. Consider the following system:

dT̄

dt
= rT̄

dŪ

dt
= βV̄ Tmax + αŪTmax

dV̄

dt
= wkŪ

Note that the supersolutions T̄ , Ū, and V̄ of
this system are bounded on a finite time interval.
We can also re-write this system in matrix form
as follows:
T̄Ū
V̄




′

=


r 0 00 αTmax βTmax

0 wk 0




T̄Ū
V̄


 , (2)

where ′ = d
dt

. Since this is a linear system in finite
time with bounded coefficients, the supersolutions
T̄, Ū, and V̄ are uniformly bounded. We can now
prove existence of an optimal control as done in
[Bukkuri 2019].

Theorem 1: There exists optimal controls u and
v that maximizes the objective functional J(u,v) if
the following conditions are met:

1) The class of all initial conditions with con-
trols u and v such that u and v are Lebesgue
integrable functions on [0,T] with values in
the admissible control set along with each
state equation being satisfied is not empty

2) The admissible control set is closed and
convex

3) The right hand side of the state system is
continuous, is bounded above by a sum of
the bounded control and the state, and can
be written as a linear function of u and v

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 4 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

with coefficients depending on time and the
state variables

4) The integrand of the functional is concave
on the admissible control set and is bounded
above by c3−c2|u|θ−c1|v|ψ, where c2,c1 >
0, and θ,ψ > 1.

Proof: First, from an existence result in [Lukes
et al. 1982], since our state system has bounded
coefficients and any solutions are bounded on
the finite time interval, we obtain the existence
of the solution of the state system. Second, the
admissible control set is closed and convex, by
definition. For the third condition, the right hand
side of the state system is continuous since each
term with a denominator is nonzero. Moreover,
the system is bilinear in the controls and can be
rewritten as

~f(t, ~X,u,v) = ~γ(t, ~X), (3)

where ~X = (T,U,V ) and ~γ is a vector-valued
function of ~X. Since the solutions are bounded,
we have

|~f(t, ~X,u,v)| ≤

∣∣∣∣∣∣

r 0 00 αTmax βTmax

0 wk 0




TU
V



∣∣∣∣∣∣

≤ C1|~X|, (4)

where C1 depends on the coefficients of the sys-
tem. Also, note that the integrand of J(u,v) is
concave on the admissible control set. The ex-
istence of optimal control follows from the fact
that aT − bU − fV − 1

2
φ1u

2 − 1
2
φ2v

2 ≤ c3 −
c2|u|θ − c1|v|ψ, where c2,c1 > 0, and θ,ψ > 1,
since T(t) ≤ Tmax.

B. Characterization of Optimal Control

Now, we will characterize the optimal control
pair (u,v) using a version of the Pontryagin Maxi-
mum Principle (PMP), as done in [Bukkuri 2019;
Bukkuri 2020]. First, let’s define the Lagrangian

associated with J(u,v):

L = aT − bU −fV −
1

2
φ1u

2 −
1

2
φ2v

2

+ λ1

(
rT
(

1 −
T

Tmax

)
− (1 −u)βV T

− (1 −v)αTU −µT
)

+ λ2((1 −u)βV T + (1 −v)αTU −kU)
+ λ3(wkU − cV )
+ ν1(t)u + ν2(t)(1 −u) + δ1(t)v
+ δ2(t)(1 −v)

Here, the ν and δ terms have been added in as
penalties for non-optimal dosing patterns, i.e.

ν1(t)u = ν2(t)(1−u) = δ1(t)v = δ2(t)(1−v) = 0

for the optimal controls (u∗,v∗).

Theorem 2: Given optimal controls u∗ and v∗
and solutions of the corresponding state system,
there exist adjoint variables λi for i = 1, 2, 3 such
that
dλ1
dt

=−
∂L

∂T
=

−
[
a+λ1

(
Uα(v−1)−r

( T
Tmax
−1
)
−

Tr

Tmax
−µ

+V β(u−1)
)

+λ2(−Uα(v−1))−V β(u−1)
]
,

dλ2
dt

=−
∂L

∂U
=

−[−b+λ1(Tα(v−1))+λ2(−k−Tα(v−1))
+λ3(kw)],

dλ3
dt

=−
∂L

∂V
=

−[−f +λ1(Tβ(u−1))+λ2(−Tβ(u−1))−λ3(c)],

where λi(T) = 0 for i = 1, 2, 3 by the PMP
transversality condition. Moreover, u∗ is given by:

u∗ = min
(
max

(
0,
λ1TV β−λ2TV β

b1

)
, 1
)
, (5)

while v∗ is similarly given by:

v∗ = min
(
max

(
0,
λ1TUα−λ2TUα

b2

)
, 1
)
. (6)

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 5 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

TABLE I
PARAMETERS USED IN NUMERICAL SIMULATIONS

Parameter Meaning Estimated Value

Tmax Carrying Capacity of Healthy T Cells 1500 cells/mm3

r Production Rate of Uninfected T Cells 0.03 cells/day
µ Natural Death Rate of Uninfected T Cells 0.02 cells/day
k Death Rate of Infected T Cells due to Viral Lysis 0.24 cells/day
c Shedding Rate of Virions 2.4/day
β Viral Infection Rate by Free Virions 2.4 ∗ 10−5/mm3
α Cellular Infection Rate 2.4 ∗ 10−5/mm3
w Burst Rate of Infected T Cells Greater than 0
u Input function for RTI Treatment 0-1
v Input function for PI Treatment 0-1

Proof: To maximize the Lagrangian (with re-
spect to the optimal control pair variables), we
differentiate L with respect to u and v. From this,
we get:

∂L

∂u
= −b1u + λ1TV β −λ2TV β + ν1(t) −ν2(t)

(7)

Thus, the representation of u∗ is λ1TV β−λ2TV β
b1

.

∂L

∂v
= −b2v + λ1TUα−λ2TUα + δ1(t) −δ2(t)

(8)
Thus the representation of v∗ is λ1TUα−λ2TUα

b2
.

Since 0 ≤ u ≤ 1 and 0 ≤ v ≤ 1, the explicit
control profiles given in equations 5 and 6 are
determined. Moreover, since both state and adjoint
solutions are L∞-bounded, the right side of the
adjoint and state equations are Lipschitz for those
solutions. This ensures that the solutions to the
optimality system are unique (and thus the optimal
controls are unique), assuming the final time is not
very large. A rigorous proof of such an argument
can be found in [Burden et al. 2004] and [Fister
et al. 1998].

IV. NUMERICAL SIMULATIONS

We now aim to numerically implement our
optimal control solutions to visualize optimal drug
treatment strategies. The optimality system can be
thought of as a two-point boundary value problem,

which was solved using a fourth-order iterative
Runge-Kutta scheme, as done in [Jung et al. 2002].
In this scheme, a forward sweep of the state
equations with initial guesses for u and v was
performed. Then, a backward calculation using the
adjoint equation and an update of the controls
was made as done by [Duda 1997] and [Deininger
et al. 2003]. This was iteratively performed until
convergence was obtained.

The parameter values in Table I represent stan-
dard HIV infection values, were obtained from
[Mobisa et al. 2018], and were used in the fol-
lowing optimal control simulations.

With these parameter values and the aforemen-
tioned initial conditions, several optimal control
profiles were ran. The following objective func-
tional weighting terms were used (derived from
scaling the healthy and infected T cells and virus
populations appropriately): a = 0.00125, b = 0.01,
f = 1, φ1 = φ2 = 0.2. The following initial
conditions were also tested, but did not lead to
any qualitatively different results, so their results
are omitted here: T(0) = 1000, U(0) = 0, V(0) =
0.001; T(0) = 800, U(0) = 10, V(0) = 0.001; T(0)
= 900, U(0) = 5, V(0) = 0.01.

First, we run a control simulation. In this simu-
lation, we let w = 1000 and c = 2.4. As explained
in [Mobisa et al. 2018], w varies greatly and we
shall later explore how optimal treatment protocols
change when the burst rate is much lower or
higher. The value of 2.4 for c is based on viral
infectivity assays which found that HIV-1 strains

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 6 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

Fig. 3. HIV Control Simulation

Fig. 4. Impact of Burst Rate of Infected T Cells (w)

IIIB and RFII which lost half of their infectivity
in 4-6 h at 37◦C [Layne et al. 1989]. However,
it’s clear that different strains of HIV in different
individuals and at different temperatures may lead
to different viral shedding rates. We shall soon see

the impact of lower viral shedding rates on optimal
protocols.

Figure 3 is a depiction of optimal RTI and PI
treatment for our control case. Note that, in the
simulations below, 100 time steps is equivalent to
one day; thus, the following simulations were run
for over 600 days. From this, it’s clear to see that
the optimal treatment calls for full intensity of both
RTI and PI treatment for about 290 days.

Now, let’s take a look at what happens when
we change the burst rate. In Figure 4, the upper
image represents w = 1 and the lower image is
for w = 100, 000.

In the case of the low burst rate, we notice
that the optimal profiles does not change much–the
only change is that it takes slightly longer to reach
the maximal dosage for both RTI and PI treatment.
However, it is clear that the high burst rate has a
significantly different profile. Specifically, we see
the emergence of bang-bang controls for both RTI
and PI treatment; moreover, we notice that the
RTI treatment is given for slightly longer (≈ 20
days) in the beginning than the PI treatment. Both
treatments then have a drug holiday until ≈ day
210, after which both drugs are given at full
potency for ≈ 80 days, before all treatment is
stopped.

Let’s now consider the impact of different shed-
ding rates, for the following values for c: 0.1, 0.2,
0.3, and 5. In Figure 5, the upper three images
represent c values of 0.1, 0.2, and 0.3, respectively,
while the bottom image represents a c value of 5.
In these simulations, we again see a departure from
our control optimal treatment recommendations.
First, consider the high shedding rate value. In this
case, optimal treatment recommends full dosage
treatment of RTI and PI for ≈ 130 days instead of
≈ 280 days, before ceasing all treatment.

Now, consider the lower shedding rates. When
c is 0.1, we notice that the RTI treatment is given
at full dosage for ≈ 190 days, before switching
to rapid bang-bang controls (between full and no
treatment) until ≈ 300 days. A drug vacation is
then given for ≈ 90 days, before giving the RTI

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 7 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

Fig. 5. Impact of Shedding Rate of Virions (c)

treatment at full potential for ≈ 20 more days, then
ceasing treatment. For the PI treatment, we notice
that the treatment is prescribed at full potency for
the first ≈ 80 days, before it gradually decays
to no treatment over the next ≈ 20 days. The
treatment again picks up around day 210, and
switches to rapid bang-bang controls like the RTI
until about day 300. Then we note the presence
of bang-bang controls (mostly remaining at no
treatment, though) during the RTI drug holiday,
before continuing full treatment at ≈ day 380 for
around 30 days, then ceasing all treatment.

When the c value is 0.2, the optimal protocols
change drastically again. The RTI treatment is
prescribed at full potency for the first ≈ 145 days,
before declining to 0 treatment over the next 5
days or so. A drug holiday is then taken for around
the next 250 days, after which the treatment is
gradually increased (over ≈ 90 days) until the
end of the considered therapy horizon. The PI
treatment is given at full potential for the first
≈ 90 days; it then rapidly declines and stays at
no treatment for the rest of the 600 day interval.
Finally, let’s consider what happens when the
shedding rate is 0.3. In this case, the optimal
profile includes an RTI treatment which stays at
full potency for the entire interval. The PI treat-
ment is given at its maximum for ≈ 190 days,
before rapidly declining to no treatment, and then
gradually increasing back to full potency (over an
≈ 80 day period) after about 10 days, and staying
at this maximum level for the rest of the therapy
horizon.

Thus, from the above numerical experiments,
we can clearly see that the standard optimal treat-
ment protocol, which calls for full treatment of
both RTI and PI for ≈ 290 days, can drastically
change depending on each patient’s specific in-
fected T cell burst rate and viral shedding rate.

V. SENSITIVITY ANALYSIS
To further assess the impact of parameter values

on HIV dynamics we perform a sensitivity analysis
on our ODE system using the Sobol-Martinez
method. Other than u and v, which were given
ranges of 0 to 1, all other parameters were given

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 8 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

Fig. 6. Sensitivity Analysis of HIV with RTI and PI Treatments

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 9 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

a range of 10% lower to 10% higher than the
values indicated in the table in Section 5. The
Sobol-Martinez method uses variance decompo-
sition techniques to measure the contributions of
input parameters to output variance. The algorithm
outlined by [Zhang et al. 2015] was implemented
here.

The left panel of each image is the first
order Sobol index: the contribution to the output
variance by the single input alone. The right
panel is the total Sobol index: the contribution
to the output variance caused by a model input,
including the first-order effects as well as all
higher-order interactions.

Note that, in the Figure 6, uther is equivalent to
u and vther is equivalent to v. Also, all parameters
with sensitivity orders less than 0.01 were omitted
from the figures. First, consider the healthy T cells.
In this case, we see three main parameters which
impact the T cell population: the burst rate of the
infected T cells, the RTI treatment, and the viral
shedding rate. The other parameters which had
minor effects on T cell population dynamics are
the natural death rate of healthy T cells and the
viral infection rate by free virions. Thus, the most
effective treatments, solely in terms of the CD4+
cell population, are the RTI treatment and those
that target the burst rate of infected T cells and the
shedding rate of virions. Considering the infected
T cells, we note that almost all the parameters are
effective in changing its population dynamics–thus
the sensitivity analysis cannot be very enlightening
for treatment planning. However, when we con-
sider the viral graph, we see that the the shedding
rate is most significant, followed by the burst rate
of infected cells.

Thus, for the overall desired dynamics, we deem
that further medical research into the development
of drugs which specifically target the shedding
rate of virions and the burst rate of infected
cells is advised, as these are the most effective
ways of minimizing the virus and infected T cell
population, while maximizing the healthy T cell
count. These results are in accordance with our nu-

merical optimal control experiments, which show
great changes in optimal treatment protocols when
shedding and burst rate parameters are changed.

VI. CONCLUSION

In this paper, we performed optimal control
and sensitivity analyses of an HIV model with
reverse transcriptase inhibitor and protease in-
hibitor antiretroviral treatments. Proofs of exis-
tence, uniqueness, and analytical characterization
of the optimal control profile were given. Numer-
ical optimal control simulations were performed
and it was found that burst rate of infected T
cells and viral shedding rate have great impacts on
optimal control profiles, often suggesting intricate
bang-bang controls for the treatment. We hope that
the results of this analysis will help physicians
more effectively assess and treat HIV patients.
Specifically, we hope physicians will measure and
include viral shedding rates and infected T cell
burst rates into treatment consideration. Moreover,
we hope that increased pharmaceutical efforts go
into developing drugs which specifically target
burst and shedding rates. Finally, we hope that
this inspires future work into creation of optimal
treatment profiles for other similar viral infections.

REFERENCES

[] Arts AJ, Hazuda DK (2012). HIV-1 Antiretroviral Drug
Therapy. Cold Spring Harbor Prespectives in Medicine 2,
1-23.

Autran B, Carcelain G, Li TS, Blanc C, Mathez D,
Tubiana R, Katlama C, Debre P, Leibowitch J (1997).
Positive effects of combined antiretroviral therapy on
CD4+ T cell homeostasis and function in advanced HIV
disease. Science 277, 112 –116.

Bracq L, Xie M, Benichou S, Bouchet J (2018).
Mechanisms for Cell-to-Cell Transmission of HIV-1.
Front. in Immuno. 9:260, 1-14.

Bukkuri A (2019). Optimal control analysis of combined
anti-angiogenic and tumor immunotherapy. Open Jour. of
Math. Sci. 3, 349-357.

Bukkuri A (2020). Optimal control analysis of combined
chemotherapy-immunotherapy treatment regimens in a
PKPD cancer evolution model. Biomath 9(1), 1-12.

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 10 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

Burden T, Ernstberger J, Fister R (2004). Optimal Control
Applied to Immunotherapy. Disc. & Cont. Dyn. Sys. 4(1),
135–146.

Cenker JJ, Stultz RD, McDonald D (2017). Brain
Microglial Cells Are Highly Susceptible to HIV-1
Infection and Spread. AIDS Research and Human
Retroviruses 33(11), 1155-1165.

Chen HY, Mascio MD, Perelson AS, Ho DD, Zhang L
(2007). Determination of virus burst size in vivo using
a single-cycle SIV in rhesus macaques. PNAS 104(48),
19079-19084.

Chun TW, Carruth L, Finzi D, Shen X, Digiuseppe JA,
Taylor H, Hermankova M, Chadwick K, Margolick J,
Quinn TC, et al. (1997). Quantification of latent tissue
reservoirs and total body viral load in HIV-1 infection.
Nature 387, 183–188.

Cunningham AL, Donaghy H, Harman AN, Kim M,
Turville SG (2010). Manipulation of dendritic cell
function by viruses. Curr. Op. in Microbiol. 13 (4),
524–529.

Deininger MWN, Druker BJ (2003). Specific targeted
therapy of chronic myelogenous leukemia with imatinib.
Pharmacol. Rev. 55(3), 401-423.

Dimitrov DS, Willey RL, Sato H, Chang LJ, Blumenthal
R, Martin MA (1993). Quantitation of human
immunodeficiency virus type 1 infection kinetics. J.
Virol. 67, 2182–2190.

Duda Z (1997). Numerical solutions to bilinear models
arising in cancer chemotherapy. Nonlinear World 4, 53-72.

Eckstein DA, Penn ML, Korin YD, Scripture-Adams DD,
Zack JA, Kreisberg JF, Roederer M, Sherman MP, Chin
PS, Goldsmith MA (2001). HIV-1 actively replicates in
naive CD4(+) T cells residing within human lymphoid
tissues. Immunity 15, 671–682.

Fister R, Lenhart S, McNally J (1998). Optimizing
Chemotherapy in an HIV , Elec. J.DE., 32, 1–12.

Fleming W, Rishel R (1975). Deterministic and Stochastic
Optimal Control. Springer-Verlag.

Gotte M, Wainberg MA (2000). Biochemical mechanisms
involved in overcoming HIV resistance to nucleoside
inhibitors of reverse transcriptase. Drug Resist. Update 3,
30 –38.

Gulnik SV, Suvorov LI, Liu B, Yu B, Anderson B,
Mitsuya H, Erickson JW (1995). Kinetic characterization
and cross-resistance patterns of HIV-1 protease mutants

selected under drug pressure. Biochem. 34, 9282–9287.

Haase AT, Henry K, Zupancic M, Sedgewick G, Faust
RA, Melroe H, Cavert W, Gebhard K, Staskus K,
Zhang Z-Q, et al. (1996). Quantitative image analysis of
HIV-1 infection in lymphoid tissue. Science 274, 985–989.

Hockett RD, Kilby JM, Derdeyn CA, Saag MS, Sillers
M, Squires K, Chiz S, Nowak MA, Shaw GM, Bucy RP.
(1999). Constant Mean Viral Copy Number per Infected
Cell in Tissues Regardless of High, Low, or Undetectable
Plasma HIV RNA. J. Exp. Med. 189, 1545–1554.

Hyman P, Abedon ST (2009). Practical methods for
determining phage growth parameters. Methods Mol. Biol.
501, 175-202.

Jung E, Lenhart S, Feng Z (2002). Optimal control of
treatments in a two-strain tuberculosis model. Discrete
Contin.Dyn. Syst. Ser. 2, 473-482.

Komanduri KV, Viswanathan MN, Wieder ED, Schmidt
DK, Bredt BM, Jacobson MA, McCune JM (1998).
Restoration of cytomegalovirus-specific CD4+ T-
lymphocyte responses after ganciclovir and highly active
antiretroviral therapy in individuals infected with HIV-1.
Nat. Med. 4, 953–956.

Layne SP, Spouge JL, Dembo M (1989). Quantifying
the infectivity of HIV. Proc. Natl. Acad. Sci. U.S.A. 86,
4644-4648.

Lederman MM, Connick E, Landay A, Kuritzkes DR,
Spritzler J, St Clair M, Kotzin BL, Fox L, Chiozzi
MH, Leonard JM, et al. (1998). Immunologic responses
associated with 12 weeks of combination antiretroviral
therapy consisting of zidovudine, lamivudine, and
ritonavir: Results of AIDS Clinical Trials Group Protocol
315. J. Infect. Dis. 178, 70–79.

Lukes D (1982). Differential Equations: Classical to
Controlled. Math. in Sci. and Engrg. 162.

Michaud V, Bar-Magen T, Turgeon J, Flockhart D,
Desta Z, Wainberg MA (2012). The dual role of
pharmacogenetics in HIV treatment: mutations and
polymorphisms regulating antiretroviral drug resistance
and disposition. Pharmacol. Rev. 64(3), 803-833.

Mobisa B, Lawi GO, Nthiiri JK (2018). Modelling In
Vivo HIV Dynamics under Combined Antiretroviral
Treatment. Jour. App. Math. (2018), 1-11.

PENTA (2004). Highly active antiretroviral therapy started
in infants under 3 months of age: 72-week follow-up for
CD4 cell count, viral load and drug resistance outcome.
AIDS 18(2), 237-245.

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 11 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173


A. Bukkuri, The impact of infected T lymphocyte burst rate and viral shedding rate on optimal ...

Perelson AS, Kirschner DE, Boer RD (1993). Dynamics
of HIV infection of CD4+ T cells. Math. Biosci. 114(1),
81–125.

Tronstein E, Johnston C, Huang M, Selke S, Magaret A,
Warren T, Corey L, Wald A (2011). Genital Shedding
of Herpes Simplex Virus Among Symptomatic and
Asymptomatic Persons With HSV-2 Infection. JAMA

305(14), 1441-1449.

UNAIDS, WHO (2007). 2007 AIDS epidemic update
(PDF). p.10. Retrieved January 2, 2020.

Zhang XY, Trame MN, Lesko LJ, Schmidt S (2015).
Sobol Sensitivity Analysis: A Tool to Guide the
Development and Evaluation of Systems Pharmacology
Models. CPT Pharmacometrics Syst. Pharmacol. 4, 69-79.

Biomath 9 (2020), 2008173, http://dx.doi.org/10.11145/j.biomath.2020.08.173 Page 12 of 12

http://dx.doi.org/10.11145/j.biomath.2020.08.173

	Introduction
	Model Description
	Optimal Control Analysis
	Existence of Optimal Control
	Characterization of Optimal Control

	Numerical Simulations
	Sensitivity Analysis
	Conclusion
	References