Optimal Control of a Delayed HIV Infection Model with Immune Response Using an Efficient Numerical Method

Abstract

We present a delay-differential equation model with optimal control that describes the interactions between human immunodeficiency virus (HIV), CD4⁺ T cells, and cell-mediated immune response. Both the treatment and the intracellular delay are incorporated into the model in order to improve therapies to cure HIV infection. The optimal controls represent the efficiency of drug treatment in inhibiting viral production and preventing new infections. Existence for the optimal control pair is established, Pontryagin’s maximum principle is used to characterize these optimal controls, and the optimality system is derived. For the numerical simulation, we propose a new algorithm based on the forward and backward difference approximation.

1. Introduction

Human immunodeficiency virus (HIV) is a lentivirus that causes acquired immunodeficiency syndrome (AIDS), a condition in humans in which the immune system begins to fail, leading to life-threatening opportunistic infections. There are several other ways the infection can transfer, for example, open wound, saliva, and ulcers.

There are some antiretroviral (ARV) drugs available nowadays which help the immune system in preventing the infection due to HIV even though it is not possible to cure it. Reverse transcriptase inhibitors (RTIs) are one of the chemotherapies which oppose the conversion of RNA of the virus to DNA (reverse transcription), so that the viral population will be minimum and on the other hand the CD4⁺ count remains higher and the host can survive. Another one is the protease inhibitors (PIs) which prevent the production of viruses from the actively infected CD4⁺ T cells.

In the literature, many mathematical models have been developed in order to understand the dynamics of HIV infection [1–6]. In addition, optimal control methods have been applied to the derivation of optimal therapies for this HIV infection [7–13]. All these methods are based on HIV models which ignore the intracellular delay by assuming that the infectious process is instantaneous; that is, in the very moment that the virus enters an uninfected cell, this one starts to produce virus particles, and we know that this is not biologically reasonable. In this paper, we consider the mathematical model for HIV infection with intracellular delay and cell-mediated immune response presented by Zhu and Zou in [6] and we introduce two controls, one simulating effect of RTIs and the other control simulating effect of PIs, incorporating drug efficacy. The intracellular delay represents the time needed for infected cells to produce virions after viral entry.

The paper is organized as follows. Section 2 describes a delayed mathematical model of HIV infection with two control terms. The analysis of optimization problems is presented in Section 3. In Section 4, we present a numerical appropriate method and the simulation corresponding results. Finally, the conclusions are summarized in Section 5.

2. HIV Model with Intracellular Delay

We consider the mathematical model for HIV-1 infection with intracellular delay and cell-mediated immune response presented by Zhu and Zou in [6]. The dynamics of this model are governed by the following equations: 𝑑𝑥(𝑡)(𝑡)𝑥(𝑡),𝑑𝑡=𝑠−𝑑𝑥−𝑘𝑣𝑑𝑦𝑑𝑡=𝑘𝑒^−𝛿𝜏𝑣(𝑡−𝜏)𝑥(𝑡−𝜏)−𝛿𝑦(𝑡)−𝑝𝑦(𝑡)𝑧(𝑡),𝑑𝑣(𝑡)(𝑡),𝑑𝑡=𝑁𝛿𝑦−𝜇𝑣𝑑𝑧(𝑡)𝑧(𝑡)(𝑡),𝑑𝑡=𝑐𝑦−𝑏𝑧 (1) where 𝑥(𝑡), 𝑦(𝑡), 𝑣(𝑡), and 𝑧(𝑡) denote the concentrations of uninfected cells, infected cells, and virus and cytotoxic T lymphocytes (CTLs), respectively.

Susceptible host cells (𝑥) are produced at a rate 𝑠, die at a rate 𝑑, and become infected by virus at a rate 𝑘. Infected cells die at a rate 𝛿 and are killed by the CTLs response at a rate 𝑝. Free virus is produced by infected cells at a rate 𝑁𝛿 and decays at a rate 𝜇 where 𝑁 is the number of free virus produced by infected cells. CTLs expand in response to viral antigen derived from infected cells at a rate 𝑐 and decay in the absence of antigenic stimulation at a rate 𝑏. The intracellular delay, 𝜏, represents the time needed for infected cells to produce virions after viral entry.

We introduce two controls 𝑢₁ and 𝑢₂ which measure the efficiency of reverse transcriptase and protease inhibitors, respectively. Hence, (1) becomes 𝑑𝑥(𝑡)−𝑑𝑡=𝑠−𝑑𝑥1−𝑢₁(𝑡)𝑘𝑣(𝑡)𝑥(𝑡),𝑑𝑦=𝑑𝑡1−𝑢₁(𝑡)𝑘𝑒^−𝛿𝜏𝑣(𝑡−𝜏)𝑥(𝑡−𝜏)−𝛿𝑦(𝑡)−𝑝𝑦(𝑡)𝑧(𝑡),𝑑𝑣=𝑑𝑡1−𝑢₂(𝑡)𝑁𝛿𝑦(𝑡)−𝜇𝑣(𝑡),𝑑𝑧(𝑡)𝑧(𝑡)(𝑡).𝑑𝑡=𝑐𝑦−𝑏𝑧 (2) The control functions, 𝑢₁(𝑡) and 𝑢₂(𝑡), are bounded, Lebesgue integrable functions. The control 𝑢₂(𝑡) represents the efficiency of drug therapy in inhibiting viral production, such that the virion production rate under therapy is (1−𝑢₂(𝑡))𝑁𝛿.

If 𝑢₂=1, the inhibition is 100% effective, whereas if 𝑢₂=0, there is no inhibition.

The control 𝑢₁(𝑡) represents the efficiency of drug therapy in blocking new infection, so that infection rate in the presence of drug is (1−𝑢₁(𝑡))𝛽.

Let 𝐶=𝐶([−𝜏,0],ℝ⁴) be the Banach space of continuous functions mapping the interval [−𝜏,0] into ℝ⁴ with the topology of uniform convergence. It is easy to show that there exists a unique solution (𝑥(𝑡),𝑦(𝑡),𝑣(𝑡),𝑧(𝑡)) of system (2) with initial data (𝑥₀,𝑦₀,𝑣₀,𝑧₀)∈𝐶.

In addition, for biological reasons, we assume that the initial data for system (2) satisfy 𝑥₀(𝑠)≥0,𝑦₀(𝑠)≥0,𝑣₀(𝑠)≥0,𝑧₀(𝑠)[].≥0,𝑠∈−𝜏,0 (3)

3. The Optimal Control Problems

The problem is to maximize the objective functional 𝐽𝑢₁,𝑢₂=^𝑡𝑓₀𝑥(𝑡)(𝑡)−𝐴+𝑧¹2𝑢²₁(𝑡)+𝐴₂2𝑢²₂(𝑡)𝑑𝑡, (4) where the parameters 𝐴₁≥0 and 𝐴₂≥0 are based on the benefits and costs of the treatment. Our target is to maximize the objective functional defined in (4) by increasing the number of the uninfected cells, maximizing immune response by CTLs, decreasing the viral load, and minimizing the cost of treatment. In other words, we are seeking optimal control pair (𝑢^∗₁,𝑢^∗₂) such that 𝐽𝑢^∗₁,𝑢^∗₂=max𝐽𝑢₁,𝑢₂∶𝑢₁,𝑢₂∈𝑈, (5) where 𝑈 is the control set defined by 𝑈=𝑢=𝑢₁,𝑢₂∶𝑢_𝑖measurable,0≤𝑢_𝑖(𝑡)≤1,𝑡∈0,𝑡_𝑓,𝑖=1,2. (6)

3.1. Existence of an Optimal Control Pair

The existence of the optimal control pair can be obtained using a result by Fleming and Rishel in [14] and by Lukes in [15].

Theorem 1. There exists an optimal control pair (𝑢^∗₁,𝑢^∗₂)∈𝑈 such that𝐽𝑢^∗₁,𝑢^∗₂=max_{𝑢1,𝑢2∈𝑈}𝐽𝑢₁,𝑢₂. (7)

Proof. To use an existence result in [14], we must check the following properties.(1)The set of controls and corresponding state variables is nonempty. (2)The control 𝑈 set is convex and closed. (3)The right-hand side of the state system is bounded by a linear function in the state and control variables. (4)The integrand of the objective functional is concave on 𝑈. (5)There exists constants 𝑐₁,𝑐₂>0, and 𝛽>1 such that the integrand 𝐿(𝑥,𝑧,𝑢₁,𝑢₂) of the objective functional satisfies 𝐿𝑥,𝑧,𝑢₁,𝑢₂≤𝑐₂−𝑐₁𝑢₁²+𝑢₂^2𝛽/2. (8)
In order to verify these conditions, we use a result by Lukes in [15] to give the existence of solutions of system (2) with bounded coefficients, which gives condition 1. We note that the solutions are bounded. Our control set satisfies condition 2. Since our state system is bilinear in 𝑢₁,𝑢₂, the right-hand side of system (2) satisfies condition 3, using the boundedness of the solutions. Note that the integrand of our objective functional is concave. Also we have the last condition needed𝐿𝑥,𝑧,𝑢₁,𝑢₂≤𝑐₂−𝑐₁𝑢₁²+𝑢₂², (9) where 𝑐₂ depends on the upper bound on 𝑥 and 𝑧, and 𝑐₁>0 since 𝐴₁,𝐴₂>0. We conclude that there exists an optimal control pair.

3.2. Optimality System

Pontryagin’s minimum Principle with delay given in [16] provides necessary conditions for an optimal control problem. This principle converts (2), (4), and (5) into a problem of maximizing an Hamiltonian, 𝐻, with 𝐻𝑡,𝑥,𝑦,𝑣,𝑧,𝑥_𝜏,𝑣_𝜏,𝑢₁,𝑢₂=𝐴,𝜆₁2𝑢²₁+𝐴₂2𝑢²₂−𝑥−𝑧+𝜆₁𝑠−𝑑𝑥−1−𝑢₁𝑘𝑣𝑥+𝜆₂1−𝑢₁𝑘𝑒^−𝛿𝜏𝑣_𝜏𝑥_𝜏−𝛿𝑦−𝑝𝑦𝑧+𝜆₃𝑁𝛿𝑦−𝜇𝑣+𝜆₄𝑐𝑦𝑧−𝑏𝑧. (10)

By applying Pontryagin’s minimum principle with delay in state [16], we obtain the following theorem.

Theorem 2. Given optimal controls 𝑢^∗₁, 𝑢^∗₂ and solutions 𝑥^∗, 𝑦^∗,𝑣^∗ and 𝑧^∗ of the corresponding state system (2), there exists adjoint variables, 𝜆₁, 𝜆₂, 𝜆₃, and 𝜆₄ satisfying the equations𝜆^′₁(𝑡)=1+𝜆₁(𝑡)𝑑+1−𝑢^∗₁(𝑡)𝑘𝑣^∗(𝑡)+𝜒_{[0,𝑡𝑓−𝜏]}(𝑡)𝜆₂(𝑡+𝜏)𝑢^∗₁(𝑡+𝜏)−1𝑘𝑒^−𝛿𝜏𝑣^∗(𝑡),𝜆^′₂(𝑡)=𝜆₂(𝑡)𝛿−𝜆₃(𝑡)1−𝑢^∗₂(𝑡)𝑁𝛿−𝑐𝑧^∗(𝑡)𝜆₄(𝑡)+𝑝𝑧^∗(𝑡)𝜆₂(𝑡),𝜆^′₃(𝑡)=𝜆₁(𝑡)1−𝑢^∗₁(𝑡)𝑘𝑥^∗(𝑡)+𝜇𝜆₃(𝑡)+𝜒_{[0,𝑡𝑓−𝜏]}(𝑡)𝜆₂(𝑡+𝜏)𝑢^∗₁(𝑡+𝜏)−1𝑘𝑒^−𝛿𝜏𝑥^∗(𝑡),𝜆^′₄(𝑡)=1+𝑝𝑦^∗(𝑡)𝜆₂(𝑡)+𝜆₄(𝑡)𝑏−𝑐𝑦^∗(𝑡), (11) with transversality conditions 𝜆_𝑖𝑡_𝑓=0,𝑖=1,,4. (12) Moreover, the optimal control is given by 𝑢^∗₁(𝑡)𝑘=min1,max0,𝐴₁𝜆₂(𝑡)𝑒^−𝛿𝜏𝑣^∗(𝑡−𝜏)𝑥^∗(𝑡−𝜏)−𝜆₁(𝑡)𝑣^∗(𝑡)𝑥^∗(𝑡𝑢_∗2)(𝑡)1=min1,max0,𝐴₂𝜆₃(𝑡)𝑁𝛿𝑦^∗(𝑡). (13)

Proof. The adjoint equations and transversality conditions can be obtained by using Pontryagin’s minimum principle with delay in state [16] such that 𝜆^′₁(𝑡)=−𝜕𝐻(𝑡)𝜕𝑥−𝜒_{[0,𝑡𝑓−𝜏]}(𝑡)𝜕𝐻𝜕𝑥_𝜏(𝑡+𝜏),𝜆₁𝑡_𝑓𝜆=0,_′2(𝑡)=−𝜕𝐻(𝑡)𝜕𝑦,𝜆₂𝑡_𝑓𝜆=0,_′3(𝑡)=−𝜕𝐻(𝑡)𝜕𝑣−𝜒_{[0,𝑡𝑓−𝜏]}(𝑡)𝜕𝐻𝜕𝑣_𝜏(𝑡+𝜏),𝜆₃𝑡_𝑓𝜆=0,_′4(𝑡)=−𝜕𝐻(𝑡)𝜕𝑧,𝜆₄𝑡_𝑓=0. (14) The optimal control 𝑢^∗₁ and 𝑢^∗₂ can be solved from the optimality conditions 𝜕𝐻𝜕𝑢₁(𝑡)=0,𝜕𝐻𝜕𝑢₂(𝑡)=0. (15) That is, 𝜕𝐻𝜕𝑢₁(𝑡)=𝐴₁𝑢₁(𝑡)+𝑘𝑣(𝑡)𝑥(𝑡)𝜆₁(𝑡)−𝑘𝑒^−𝛿𝜏𝑣(𝑡−𝜏)𝑥(𝑡−𝜏)𝜆₂(𝑡)=0,𝜕𝐻𝜕𝑢₂(𝑡)=𝐴₂𝑢₂(𝑡)−𝑁𝛿𝑦(𝑡)𝜆₃(𝑡)=0. (16) By the bounds in 𝑈 of the controls, it is easy to obtain 𝑢^∗₁ and 𝑢^∗₂ in the form of (13), respectively.

If we substitute 𝑢^∗₁ and 𝑢^∗₂ in the systems (2) and (11), we obtain the following optimality system: 𝑑𝑥^∗𝑑𝑡=𝑠−𝑑𝑥^∗(𝑡)−1−𝑢^∗₁(𝑡)𝑘𝑣^∗(𝑡)𝑥^∗(𝑡),𝑑𝑦^∗=𝑑𝑡1−𝑢^∗₁(𝑡)𝑘𝑒^−𝛿𝜏𝑣^∗(𝑡−𝜏)𝑥^∗(𝑡−𝜏)−𝛿𝑦^∗(𝑡)−𝑝𝑦^∗(𝑡)𝑧^∗(𝑡),𝑑𝑣^∗=𝑑𝑡1−𝑢^∗₂(𝑡)𝑁𝛿𝑦^∗(𝑡)−𝜇𝑣^∗(𝑡),𝑑𝑧^∗𝑑𝑡=𝑐𝑦^∗(𝑡)𝑧^∗(𝑡)−𝑏𝑧^∗(𝑡),𝑑𝜆₁𝑑𝑡=1+𝜆₁(𝑡)𝑑+1−𝑢^∗₁(𝑡)𝑘𝑣^∗(𝑡)+𝜒_{[0,𝑡𝑓−𝜏]}(𝑡)𝜆₂(𝑡+𝜏)𝑢^∗₁(𝑡+𝜏)−1𝑘𝑒^−𝛿𝜏𝑣^∗(𝑡),𝑑𝜆₂𝑑𝑡=𝜆₂(𝑡)𝛿−𝜆₃(𝑡)1−𝑢^∗₂(𝑡)𝑁𝛿−𝑐𝑧^∗(𝑡)𝜆₄(𝑡)+𝑝𝑧^∗𝑡𝜆₂(𝑡),𝑑𝜆₃𝑑𝑡=𝜆₁(𝑡)1−𝑢^∗₁(𝑡)𝑘𝑥^∗(𝑡)+𝜇𝜆₃(𝑡)+𝜒_{[0,𝑡𝑓−𝜏]}(𝑡)𝜆₂(𝑡+𝜏)𝑢^∗₁(𝑡+𝜏)−1𝑘𝑒^−𝛿𝜏𝑥(𝑡),𝑑𝜆₄𝑑𝑡=1+𝑝𝑦^∗(𝑡)𝜆₂(𝑡)+𝜆₄(𝑡)𝑏−𝑐𝑦^∗(𝑡),𝑢^∗₁(𝑡)𝑘=min1,max0,𝐴₁𝜆₂(𝑡)𝑒^−𝛿𝜏𝑣^∗(𝑡−𝜏)𝑥^∗(𝑡−𝜏)−𝜆₁(𝑡)𝑣^∗(𝑡)𝑥^∗(𝑡,𝑢_∗2)(𝑡)1=min1,max0,𝐴₂𝜆₃(𝑡)𝑁𝛿𝑦^∗(𝑡,𝜆_𝑖)𝑡_𝑓=0,𝑖=1,,4. (17)

4. Numerical Simulations

In this section, we give a numerical method to solve the optimality system (17) and present the results.

Let there exists a step size ℎ>0 and integers (𝑛,𝑚)∈ℕ² with 𝜏=𝑚ℎ and 𝑡_𝑓−𝑡₀=𝑛ℎ. For reasons of programming, we consider 𝑚 knots to left of 𝑡₀ and right of 𝑡_𝑓, and we obtain the following partition: Δ=𝑡_−𝑚=−𝜏<⋯<𝑡₋₁<𝑡₀=0<𝑡₁<⋯<𝑡_𝑛=𝑡_𝑓<𝑡_𝑛+1<⋯<𝑡_𝑛+𝑚. (18) Then, we have 𝑡_𝑖=𝑡₀+𝑖ℎ(−𝑚≤𝑖≤𝑛+𝑚). Next, we define the state and adjoint variables 𝑥(𝑡), 𝑦(𝑡), 𝑣(𝑡), 𝑧(𝑡), 𝜆₁(𝑡), 𝜆₂(𝑡), 𝜆₃(𝑡), 𝜆₄(𝑡) and the controls 𝑢₁(𝑡), 𝑢₂(𝑡) in terms of nodal points 𝑥_𝑖, 𝑦_𝑖, 𝑣_𝑖, 𝑧_𝑖, 𝜆^𝑖₁, 𝜆^𝑖₂, 𝜆^𝑖₃, 𝑢^𝑖₁, and 𝑢^𝑖₂.

Now a combination of forward and backward difference approximation, we obtain Algorithm 1.

Step 1: for 𝑖=−𝑚,,0,do∶      𝑥𝑖=𝑥0, 𝑦𝑖=𝑦0, 𝑣𝑖=𝑣0, 𝑧𝑖=𝑧0, 𝑢𝑖1=0, 𝑢𝑖2=0.    end for    for 𝑖=𝑛,,𝑛+𝑚,do∶      𝜆𝑖1=0, 𝜆𝑖2=0, 𝜆𝑖3=0, 𝜆𝑖4=0.    end for Step 2: for 𝑖=0,,𝑛−1,do∶      𝑥𝑖+1=𝑥𝑖+ℎ[𝑠−𝑑𝑥𝑖−(1−𝑢𝑖1)𝑘𝑣𝑖𝑥𝑖],      𝑦𝑖+1=𝑦𝑖+ℎ[𝑘(1−𝑢𝑖1)𝑒−𝛿𝜏𝑣𝑖−𝑚𝑥𝑖−𝑚−𝛿𝑦𝑖−𝑝𝑦𝑖𝑧𝑖],      𝑣𝑖+1=𝑣𝑖+ℎ[(1−𝑢𝑖2)𝑁𝛿𝑦𝑖−𝜇𝑣𝑖],      𝑧𝑖+1=𝑧𝑖+ℎ(𝑐𝑦𝑖𝑧𝑖−𝑏𝑧𝑖),      𝜆1𝑛−𝑖−1=𝜆1𝑛−𝑖−ℎ{1+𝜆1𝑛−𝑖[𝑑−(1−𝑢𝑖1)𝑘𝑣𝑖+1]              +𝜒[0,𝑡𝑓−𝜏](𝑡𝑛−𝑖)𝜆2𝑛−𝑖+𝑚(𝑢1𝑖+𝑚−1)𝑘𝑒−𝛿𝜏𝑣𝑖+1},      𝜆2𝑛−𝑖−1=𝜆2𝑛−𝑖−ℎ[𝜆2𝑛−𝑖𝛿−𝜆3𝑛−𝑖(1−𝑢𝑖2)𝑁𝛿−𝑐𝑧𝑖+1𝜆4𝑛−𝑖+𝑝𝑧𝑖+1𝜆2𝑛−𝑖],      𝜆3𝑛−𝑖−1=𝜆3𝑛−𝑖−ℎ[𝜆1𝑛−𝑖(1−𝑢𝑖1)𝑘𝑥𝑖+1+𝜇𝜆3𝑛−𝑖              +𝜒[0,𝑡𝑓−𝜏](𝑡𝑛−𝑖)𝜆2𝑛+𝑚−𝑖(𝑢1𝑖+𝑚−1)𝑘𝑒−𝛿𝜏𝑥𝑖+1],      𝜆4𝑛−𝑖−1=𝜆4𝑛−𝑖−ℎ[1+𝑝𝜆2𝑛−𝑖𝑦𝑖+1+𝜆4𝑛−𝑖(𝑏−𝑐𝑦𝑖+1)],      𝑅1𝑖+1=(𝑘/𝐴1)(𝜆2𝑛−𝑖−1𝑒−𝛿𝜏𝑣𝑖−𝑚+1𝑥𝑖−𝑚+1−𝜆1𝑛−𝑖−1𝑣𝑖+1𝑥𝑖+1),      𝑅2𝑖+1=(1/𝐴2)𝜆3𝑛−𝑖−1𝑁𝛿𝑦𝑖+1,      𝑢1𝑖+1=min(1,max(𝑅1𝑖+1,0)),      𝑢2𝑖+1=min(1,max(𝑅2𝑖+1,0)).    end for Step 3: for 𝑖=1,,𝑛,write     𝑥∗(𝑡𝑖)=𝑥𝑖, 𝑦∗(𝑡𝑖)=𝑦∗𝑖, 𝑣∗(𝑡𝑖)=𝑣𝑖, 𝑧∗(𝑡𝑖)=𝑧𝑖, 𝑢∗1(𝑡𝑖)=𝑢𝑖1, 𝑢∗2(𝑡𝑖)=𝑢𝑖2.    end for

Algorithm 1 

For this simulation, we use the parameter values given in Table 1.

Table 1 

Parameters, their symbols and values used in simulation.

The graphs from simulating the model, given below, help to compare the uninfected cells, the infected cells, and the viral load before and after the treatments with controls.

Figure 1 shows that after the treatments, the CD4⁺ T population grows significantly which improves the quality of life of the patient.

Figure 1 

Uninfected cells 𝑥 with and without control.

In Figure 2, the number of infected CD4⁺ T cells at the final time 𝑡_𝑓=50 (days) is 0.1096 in the case with control and 1.8439 without control, and the total cases in blocking new infections at the end of the control program is 1.7343, then the efficiency of drug therapy in blocking the new infections is 94%.

Figure 2 

Infected cells 𝑦 with and without control.

Figure 3 shows that after introducing therapy, the viral load declines towards zero. Specifically, the number of free virus at the final time is 0.0250 in the case with control and 162.7342 without control, and the total cases in blocking viral production at the end of the control program is 162.7092, then the efficiency of drug therapy in inhibiting viral production is 99.98%.

Figure 3 

Virus 𝑣 with and without control.

Figure 4 shows that the cell-mediated immune response is always maintained at a positive level and it is never eliminated. We also note that an increase in infection is followed by a corresponding increase in the immune response, which then serves to remove infection by killing off infected cells. Once the infection is low, the immune response is not needed at such high levels and this is why it drops off too. Finally, Figure 5 represents the optimal controls 𝑢^∗₁ and 𝑢^∗₂ in blocking new infection and inhibiting viral production.

Figure 4 

Function 𝑧 with and without control.

(a)

(b)

Figure 5 

The controls 𝑢1 and 𝑢2.

5. Conclusion

The purpose of this work is two-fold. Firstly, we gave a delay mathematical model with two controls that describe HIV infection of CD4⁺ T cells during therapy. Currently, there is no effective therapy for HIV infection and the cost of treatment is beyond reach of many infected patients. Hence, we presented an optimal therapy in order to minimize the cost of treatment, reduce the viral load, and improve immune response. Secondly, we presented an efficient numerical method based on optimal control to identify the best treatment strategy of HIV infection in order to block new infection and prevent viral production by using drug therapy with minimum side effects.

Our numerical results show that the optimal treatment strategies reduce viral load and increase the uninfected CD4⁺ T-cell count after five days of therapy.

Acknowledgment

The authors would like to thank the anonymous referee for his/her valuable comments on the first version of the paper which have led to an improvement in this revised version.