# Fuzzy Modeling and Control of HIV Infection

Authors: Hassan Zarei; Ali Vahidian Kamyad; Ali Akbar Heydari
Journal: Computational and Mathematical Methods in Medicine
 (2012)

---

## Abstract

The present study proposes a fuzzy mathematical model of HIV infection consisting of a linear fuzzy differential equations (FDEs) system describing the ambiguous immune cells level and the viral load which are due to the intrinsic fuzziness of the immune system’s strength in HIV-infected patients. The immune cells in question are considered CD4+ T-cells and cytotoxic T-lymphocytes (CTLs). The dynamic behavior of the immune cells level and the viral load within the three groups of patients with weak, moderate, and strong immune systems are analyzed and compared. Moreover, the approximate explicit solutions of the proposed model are derived using a fitting-based method. In particular, a fuzzy control function indicating the drug dosage is incorporated into the proposed model and a fuzzy optimal control problem (FOCP) minimizing both the viral load and the drug costs is constructed. An optimality condition is achieved as a fuzzy boundary value problem (FBVP). In addition, the optimal fuzzy control function is completely characterized and a numerical solution for the optimality system is computed.


---

## Body

## 1. Introduction

Usage of fuzzy differential equations is a natural way to model dynamical systems under uncertainty [1]. For example these equations are used to modeling the cell growth and dynamic of population [2], dry friction [3], tumor growth [4], and the phenomenon of nuclear disintegration [5] under uncertainty. In [6], transition from HIV to AIDS (the acquired immunodeficiency syndrome) is described through a mathematical model with fuzzy transference rate correlated with the viral load and CD4+ T-cells level by rule bases. Moreover, in [7], the authors have proposed a methodology combining a macroscopic HIV-positive population model, which is a differential equation system whose transference rate from asymptomatic to symptomatic population is found through a fuzzy rule-based system, with an individual microscopic model to study the evolution of positive HIV population for manifestation of AIDS. In [8], a fuzzy delay differential equation is proposed to model HIV infection, assuming that there exists delay between the infection of a CD4+ T-cell by the virus and the production of new virus particles. In this model, the delay and the clearance rate of HIV particles are fuzzy numbers where correlation between them is restated by rule bases. However, it should be noted that the whole parameters of a model such as the production and clearance rate of viruses and immune cells can be the source of uncertainty. In real world, there are various HIV-infected patients with different strengths of immune system causing uncertainty as to the immune cells level and the viral load during the different stages of the disease. A number of mathematical models have been formulated to describe various aspects of the interaction between HIV and the immune cells. The basic and simple model of HIV infection that contains three state variables: healthy CD4+ T-cells, infected CD4+ T-cells, and viruses, is presented by Perelson et al. [9], and more complicated models containing other parts of the immune system such as the cytotoxic T-lymphocyte and the macrophages are presented in [10] and references therein. None of these models can mirror the mentioned uncertainties proposing a mathematical model with fuzzy parameters which could reflect such ambiguities would be desirable.One of the earliest suggestions to define the concept of differentiability for fuzzy mappings and, in consequence, to study fuzzy differential equations is the Hukuhara derivative [11]. Nevertheless, the solution of fuzzy differential equation interpreted by Hukuhara derivative became fuzzier as time goes by [12]. Hence, the fuzzy solution behaves quite differently from the crisp solution. In order to overcome this difficulty, a more general definition of derivative for fuzzy-number-valued functions, which is called the strongly generalized differentiability, has been introduced and many papers have been published in this field (see, e.g., [1, 2, 5, 13–17] and references therein). In [13] a generalization of the Hukuhara differentiability to the case of interval-valued functions is introduced and the local existence and uniqueness of solutions for the interval differential equations are obtained under this type of differentiability.First-order linear fuzzy differential equations are one of the simplest fuzzy differential equations which may appear in many applications. However, the form of such an equation is very simple, it raises many problems since under different fuzzy differential equation concepts, the behavior of the solutions is different (depending on the interpretation used). This type of equations has been studied by many researchers. For example, the general form of the solutions for the first-order fuzzy differential equations with crisp coefficients under the generalized differentiability concept is presented in [14]. Moreover, an operator method is proposed for solving a class of first- and second-order linear fuzzy differential equations under the assumption of strongly generalized differentiability which is constructed based on their equivalent integral forms [15]. The existence and uniqueness of the solutions has been demonstrated for a first-order linear fuzzy differential equation with impulses subject to boundary value conditions, and the explicit solutions are obtained by calculating the solutions on each level set [16]. The generalized Euler approximation method is applied to solve numerically fuzzy differential equations under generalized differentiability [17]. A number of works in these fields have dealt with the linear fuzzy differential dynamical systems. For instance, in [18], a complex number representation of the α-level sets of the linear first-order fuzzy differential dynamical systems where the initial condition is described by a vector of fuzzy numbers is presented and the solutions are obtained under such representation. Using this approach, a method is proposed to find the solutions of a class of linear differential dynamical systems with fuzzy matrices [19] and the proposed method is extended to provide the solutions of linear matrix differential dynamical systems with fuzzy matrices [20].Classical control system is described by a differential equation. However, uncertainty is inherit in most dynamic systems. The concept of fuzzy optimal control was presented by Komolov et al. [21] in 1979. Since then, many researchers have studied this type of problems. In [22], fuzzy differential equations are generalized to be fuzzy set control differential equations (FSCDEs) and the problem of stability and controllability of FSCDE are presented. Furthermore, some properties of the fuzzy solution for the linear FSCDE as well as the necessary and sufficient optimality conditions for a linear fuzzy time optimal control problem are obtained in [23].In this paper, we model the uncertain behaviors of CD4+ T-cells and CTLs level and the viral load in different patients by a system of linear fuzzy differential equations and analyze the optimal control regarding minimizing both the viral load and drug costs.Following a preliminary introduction, in Section3, a system of linear differential equations with fuzzy parameters describing the ambiguous behaviors of CD4+ T- cells and CTLs level and the HIV viral load in patients with a weak, moderate and strong immune system is introduced. Moreover, a method for finding explicit solutions to the proposed model is introduced in this section. Some authors have used mathematical models for HIV infection in conjunction with control theory to achieve appropriate goals. Although the proposed model is simple, it can be used to investigate the effects of antiretroviral therapy in preventing the HIV progression. Section 4 is devoted to the latter topic. The last section deals with the conclusion.


## 2. Preliminaries

In this section, we give some definitions and introduce the necessary notations which will be used throughout the paper. See, for example, [12].Definition 1.
(fuzzy set,α-level set, and fuzzy number) A fuzzy set ũ in Rn is defined as a set of all pairs (x,μũ(x))∈Rn×[0,1] for some function μũ:Rn→[0,1], which is called the membership function of ũ, and μũ(x) is interpreted as the membership grade of a element x in the fuzzy set ũ. We define [ũ]α={x∈R∣μũ(x)≥α} the α-level set of ũ, with 0<α≤1. For α=0, the support of ũ is defined as [ũ]0=supp(ũ)=cl{x∈R∣μũ(x)>0}, where cl denotes the closure of a subset. A fuzzy set ũ in R is called a fuzzy number in R if(i)
μũ is upper semicontinuous on R,(ii)
ũ is a convex fuzzy set, that is, μũ(λx+(1-λ)y)≥min{μũ(x),μũ(y)}, for all x,y∈R, λ∈[0,1],(iii)
μũ is normal, that is, there exists a unique x0∈R such that μũ(x0)=1,(iv)
the support ofũ is compact.Example 2.
The normal level of CD4+ T-cells in blood, that is: “close to 1000 cells/μL,” can be represented as a fuzzy set ũ with a membership function defined as μũ(x)=exp(-β(x-1000)2) where  β is a positive real number. It is easy to see that [ũ]α=[1000--β-1lnα,1000+-β-1lnα]. The membership function μũ(·) with β=0.0005 and [ũ]α are shown in Figure 1. Note that [ũ]0=(-∞,+∞); hence, ũ is not a fuzzy number.
The set of all fuzzy numbers inR is denoted by F(R). From (i)–(iv), it follows that if ũ belongs to F(R) then the α -level set [ũ]α is a nonempty compact interval for all α∈[0,1]. The notation [ũ]α=[u̲α,u̅α] denotes explicitly the α-level set of ũ. Triangular fuzzy numbers are one of the most commonly used fuzzy numbers. The membership function of a triangular fuzzy number ũ is completely characterized with the peak (or center)m, left width σ≥0,  and right width β≥0 and has the following form:
(1)μũ(x)={1-m-xσ,ifm-σ≤x<m,1-x-mβ,ifm≤x≤m+β,0,otherwise,
and we use the notation ũ=(m,σ,β). A triangular fuzzy number ũ=(m,σ,β) is named a symmetric triangular fuzzy number if its left width  and right width   are equal and we denote it by ũ=(m,σ), for brevity. Some examples of triangular and symmetric triangular fuzzy numbers are depicted in Figures 2 and 3, respectively. We consider a crisp number a∈R as a symmetric triangular fuzzy number a=(a,0).
Letf be a real-valued mapping on Rn. Assume ũj, j=1,…,n are fuzzy numbers in R. Using the extension principle, we can define Ỹ=f(ũ1,…,ũn) as a fuzzy set in R such that
(2)μỸ(y)={sup⁡min⁡{μũ1(x1),…,μũn(xn)∣(x1,…,xn)∈f-1(y)},iff-1(y)≠ϕ,0,iff-1(y)=ϕ.Figure 1
The membership function of the fuzzy set of real numbers “close to 1000” and itsα-level set.Figure 2
The fuzzy numbersũ=(4,4,1) (- - -), ṽ=(3,1,4) (—), max̃{ũ,ṽ} (∘∘), and miñ{ũ,ṽ} (**).Figure 3
The fuzzy numbersũ=(-1,1) (—), ṽ=(2,2/3) (- - -) , ũΘHṽ (⋯), and ũ⊗Zṽ (–·–).Example 3.
Letf(x)=λx be a linear function. Suppose ũ∈F(R), and let Ỹ=f(ũ).Then, using the extension principle, we obtain μỸ(y)=sup{μũ(x)∣y=λx}=μũ(y/λ). Especially, if λ=-1, then we write f(ũ)=Θũ and we have μΘũ(y)=μũ(-y).Theorem 4.
Letf be a real-valued mapping on Rn, and let ũj, j=1,…,n, be fuzzy numbers in R. Then, [f(ũ1,…,ũn)]α=f([ũ1]α,…,[ũn]α) where f(ũ1,…,ũn)  is defined by the extension principle and f([ũ1]α,…,[ũn]α)={f(u1,…,un)∣u1∈[ũ1]α,…,un∈[ũn]α}.Proof (see [12]).
The above theorem defines arithmetic operations of fuzzy numbers in terms of theirα-level sets by [ũ⊕ṽ]α=[u̲α+v̲α,u̅α+v̅α], [ũΘṽ]α=[u̲α-v̅α,u̅α-v̲α], and [ũ⊗ṽ]α=[min{u̲αv̲α,u̲αv̅α,u̅αv̲α,u̅αv̅α}, max{u̲αv̲α,u̲αv̅α,u̅αv̲α,u̅αv̅α}], where ⊕, Θ, and ⊗ denote the addition, minus, and multiplication operators on F(R), respectively. Moreover, [Θũ]α=[-u̅α,-u̲α].Definition 5 (fuzzy max, fuzzy min, fuzzy inequality, and the weighted center of gravity).
Letũ,ṽ∈F(R), and f(x,y)=max{x,y}. Then, the maximum of ũ and ṽ is defined by f(ũ,ṽ) and applying the extension principle. Similarly, setting f(x,y)=min{x,y}, then f(ũ,ṽ) defines the minimum of ũ and ṽ. We denote the maximum and the minimum of ũ and ṽ by max̃{ũ,ṽ} and miñ{ũ,ṽ}, respectively. The notation “≤̃” will be used for the inequality relation between fuzzy numbers ũ and ṽ and is defined as ũ≤̃ṽ⇔mĩn{ũ,ṽ}=ũ. Moreover, the quantity of fuzzy number ũ can be given by its weighted center of gravity (WCOG) defined as u=∫01α(u̅α+u̲α)dα, where the weights are the membership degrees.Example 6.
Two triangular fuzzy numbersũ=(4,4,1), ṽ=(3,1,4) and their maximum and minimum are shown in Figure 2. Moreover, if ũ=(m,σ,β), then the WCOG of ũ is u=∫01α[2m+(1-α)(σ-β)]dα=m+1/6(σ-β).
Applying Theorem4, we get [mĩn(ũ,ṽ)]α=min{[ũ]α,[ṽ]α}=[min{u̲α,v̲α},min{u̅α,v̅α}]; hence, ũ≤̃ṽ⇔u̅α≤v̅α, u̲α≤v̲α, for all α∈[0,1].Definition 7 (H-difference and Z-product).
Letũ,ṽ∈F(R). If there exists z̃∈F(R) such that ũ=ṽ⊕z̃, then z̃ is called the H-difference of ũ and ṽ and it is denoted by ũΘHṽ. Moreover, if there exist z̃∈F(R) such that z̃⊗ṽ-1=ũ, then we call it the Z-product of ũ and ṽ and we denote it by z̃=ũ⊗Zṽ.
Note that the Z-product of fuzzy numbers is a new concept which is introduced in this paper for the first time. It is easy to see that[ũΘHṽ]α=[u̲α-v̲α,u̅α-v̅α] and [ũ⊗Zṽ]α=[u̲αv̲α,u̅αv̅α], if ũ≤̃(0,0) and (0,0)≤̃ṽ.Example 8.
Letũ=(-1,1) and ṽ=(2,2/3). A straightforward calculation shows [ũΘHṽ]α=[-3-1/3(1-α),-3+1/3(1-α)]; hence, ũΘHṽ=(-3,1/3). Furthermore, we obtain [ũ⊗Zṽ]α=[-8/3+2/3α2,-8/3α+2/3α2] which is valid α-level set of a fuzzy number (not a triangular fuzzy number). The fuzzy numbers ũ,ṽ and their H-difference and Z-product are shown in Figure 3, for the sake of clarity.Definition 9 (the strongly generalized differentiability).
Letx̃:(a,b)→F(R) be a fuzzy function and t0∈(a,b). We say that. x̃ is differentiable at t0 if it exists an element x̃̇(t0)∈F(R) such that, for all h>0 sufficiently near to 0,(i)
there arex̃(t0+h)ΘHx̃(t0) and x̃(t0)ΘHx̃(t0-h) and the limits limh→0+(x̃(t0+h)ΘHx̃(t0)/h)=limh→0+(x̃(t0)ΘHx̃(t0-h)/h)=x̃̇(t0) or(ii)
there arex̃(t0-h)ΘHx̃(t0) and x̃(t0)ΘHx̃(t0+h) and the limits limh→0+(x̃(t0)ΘHx̃(t0+h)/-h)=limh→0+(x̃(t0-h)ΘHx̃(t0)/-h)=x̃̇(t0), where the limits are taken in the metric D defined as D(ũ,ṽ)=supα∈[0,1]max{|u̅α-v̅α|,|u̲α-v̲α|}, for all ũ,ṽ∈F(R).Theorem 10.
Ifx̃(t) is differentiable in the first form (i), then x̲α(t) and x̅α(t) are differentiable functions and [x̃̇(t)]α=[ẋ̲α(t),x̅̇α(t)].  If x̃(t) is differentiable in the second form (ii), then x̲α(t) and x̅α(t) are differentiable functions and [x̃̇(t)]α=[x̅̇α(t),ẋ̲α(t)].Proof (see [17]).
Consider the following fuzzy initial value problem:(3)x̃̇(t)=f(t,x̃(t)),x̃(0)=x̃0,
where f:I×F(R)→F(R) is a fuzzy function and x̃0∈F(R). Let [f(t,x̃)]α=[f̲α(t,x̲α,x̅α),fα¯(t,x̲α,x̅α)]. From Theorem 10, if we consider x̃(t) by using the derivative in the first form (i), then the solution of problem (3) is obtained by solving the following system of ordinary differential equations:
(4)ẋ̲α(t)=fα̲(t,x̲α(t),x̅α(t)),x̲α(0)=x̲0α,x̅̇α(t)=fα̅(t,x̲α(t),x̅α(t)),x̅α(0)=x̅0α,
and ensuring that [x̲α(t),x̅α(t)] and [ẋ̲α(t),x̅̇α(t)] are valid α-level sets. Moreover, if we consider x̃(t) by using the derivative in the second form (ii), then the solution of problem (3) is obtained by solving the following system of ordinary differential equations,
(5)x̅̇α(t)=fα̲(t,x̲α(t),x̅α(t)),x̲α(0)=x̲0α,ẋ̲α(t)=fα̅(t,x̲α(t),x̅α(t)),x̅α(0)=x̅0α,
and ensuring that [x̲α(t),x̅α(t)]and [x̅̇α(t),ẋ̲α(t)] are valid α-level sets.
The integral of fuzzy functionx̃(t) using the Riemann integral concept can be defined as follows.Definition 11.
The integral of a fuzzy mappingx̃:[a,b]→F(R) is defined levelwise by
(6)[∫abx̃(t)dt]α=[∫abx̲α(t)dt,∫abx̅α(t)dt].
Note that if x̃:[a,b]→F(R) is continuous in the metric D, then it is integrable, that is, ∫abx̃(t)dt∈F(R).Example 12.
Define the fuzzy mappingx̃:[0,0.5]→F(R) by x̃(t)=(t,sinπt,cosπt). Then, [∫abx̃(t)dt]α=[∫00.5(t-(1-α)sinπt)dt,∫00.5(t+(1-α)cosπt)dt]=[0.125-(1-α)π-1,0.125+(1-α)π-1]; hence, ∫abx̃(t)dt=(0.125,π-1).


## 3. Linear Fuzzy Model of HIV Infection

HIV infection can be characterized as a disease of the immune system, with progressive depletion of defensive cells, resulting in immunosuppression and susceptibility to opportunistic infections. CD4+ T-cells, CTLs, and the virus particles play important roles in HIV infection. CD4+ T-cells are a fundamental component of the human immune response system. These cells can be considered “messengers” or the command centers of the immune system, and they signal other immune cells that an invader is to be fought. The immune response cells, or cytotoxic lymphocytes (CTLs), are the cells that respond to this message and set out to eliminate infection by killing infected cells. HIV can infect a number of cells in the body however, its main target is the CD4+ T-cells. HIV enters these cells by a complex process and begins to replicate, then the new virus particles are released by bursting the infected cells. CD4+ T-cells are generated from sources within the body and are lost either by having finite life span or by bursting during the proliferation of HIV, which leads to a drop in the number of these cells, after infection and an accelerated decrease during the later stages of the disease that signals the onset AIDS. In accordance with experimental findings, too high a level of HIV impairs establishment of a lasting CTL response. This is a delicate task, since CD4+ T-cell population, which plays an essential role in stimulation of immune response, depletes dramatically with raising the HIV load. The rate of CD4+ T-cells depletion varies greatly from patient to patient, depending on the strength or weakness of the immune system. More precisely, a stronger immune system leads to a lower rate of CD4+ T-cells depletion andvice versa. We have a similar argument about the proliferation rate of HIV particles. Therefore, the levels of the immune cells as well as the HIV viral load during the different stages of the disease can be considered as fuzzy quantities. According to these descriptions, the interaction of HIV with the immune system can be modeled by a system of linear differential equations with fuzzy parameters as follows:(7)x̃̇=λ̃Θσ̃⊗x̃Θc̃⊗ṽ,ṽ̇=k̃⊗ṽΘã⊗z̃,z̃̇=h̃⊗x̃Θτ̃⊗ṽ,
where the fuzzy functions x̃(t), z̃(t), and ṽ(t) indicate the level of CD4+ T-cells, CTLs, and the HIV viral load at timet, respectively. Most of the terms in the model have straightforward interpretations as follows.The first equation in (7) represents the dynamics of the concentration of CD4+ T-cells. The CD4+ T-cells are produced from a source, such as the thymus, at a constant rate λ̃. Here, we have assumed that CD4+ T-cells have a finite life-span and die at a rate σ̃ per cell. Therefore, the number of these cells, which are lost due to natural death, is represented through the loss term σ̃⊗x̃ in the first equation. Moreover, the CD4+ T-cell population is lost through infection by a virus particle at a rate of c̃, and so the term c̃⊗ṽ models the rate that free viruses destroy CD4+ T-cells. The second equation in (7) depicts the rate of change in the virus population. An HIV particle uses a host cell to replicate itself and thus proliferates with a growth rate k̃. Thus, the total amount of produced viruses is given by the term k̃⊗ṽ. Infected cells are killed by CTLs, and hence viruses are lost through an immune response. Assuming that a CTL eliminates the virus particles at a rate ã, the number of virus particles eliminated by the immune response is given by the term ã⊗z̃. The third equation in (7) describes the dynamics of CTLs during HIV infection. A CD4+ T-cell stimulates CTLs to proliferate at a rate h̃. Therefore, CD4+ T-cells effect on proliferation of CTLs is expressed by the term h̃⊗x̃. The term τ̃⊗ṽ takes into account loss of CTLs due to increasing the HIV viral load where τ̃ is the rate at which the virus-induced impairment of CD4+ T-cell function occurs.In this paper, a patient with respect to the strength or weakness of its immune system is considered as a patient with the weak, moderate, or strong immune system and is indicated byW, M, or S, respectively. The initial condition of ṽ(0)=ṽ0 varies in different patients. For that reason, a primary response is provoked when the immune system encounters HIV for the first time and, in this stage, a number of viruses depending on the strength or weakness of the immune system are eliminated, and the proposed model describes the changes in the immune cells level and the viral load after this stage which is called the secondary immune response. Therefore, a stronger immune response implies a lower ṽ0 and vice versa. The values of the model parameters and ṽ0 corresponding to patients W, M, and S are shown in Figure 4 as triangular fuzzy numbers. These parameters were chosen to be consistent with biological plausibility. Moreover, we assume that, at time t=0, the level of CD4+ T-cells is normal and there is no CTL-mediated immune response in all patients, that is, x̃(0)=x̃0=(100,0) and z̃(0)=z̃0=(0,0). We must note that x̃(t) denotes CD4+ T-cells level in percentage at time t. The derivative in the second form (ii) leads to solutions with decreasing length of their support which leads us to the conclusion that the uncertainty decreases with the time lapse which is not consistent with real situation. Moreover, the existence of these solutions implies that the initial conditions should be fuzzy. Therefore, we consider only the solutions with the derivative in the first forms (i) which are more consistent with real situation. Consequently, as mentioned in Section 2, the fuzzy model (7) is transformed to the following system of ordinary differential equations (ODEs):
(8)x̅̇α(t)=λ̅α-σ̲αx̲α(t)-c̲αv̲α(t),ẋ̲α(t)=λ̲α-σ̅αx̅α(t)-c̅αv̅α(t),v̅̇α(t)=k̅αv̅α(t)-a̲αz̲α(t),v̲̇α(t)=k̲αv̲α(t)-a̅αz̅α(t),z̅̇α(t)=h̅αx̅α(t)-τ̲αv̲α(t),ż̲α(t)=h̲αx̲α(t)-τ̅αv̅α(t),x̅α(0)=x̅0α,x̲α(0)=x̲0α,v̅α(0)=v̅0α,v̲α(0)=v̲0α,z̅α(0)=z̅0α,z̲α(0)=z̲0α.
For each α∈[0,1], the ODEs (8) are linear; hence, the exact solutions in discrete times are obtained using the ode45 code in MATLAB. However, it would be appropriate to propose explicit solutions as a function of α and t.The next section is devoted to this topic.Figure 4
The values of model parameters corresponding to patientsW (- -), M (⋯), and S (—).


### 3.1. The Approximate Explicit Solutions Based on a Fitting Method

The proposed method is based on the fact that a linear combination of suitable functions ofα and t can generate the best fit to the exact values obtained by the ode45 in the least squares sense. The following discussion shows that these functions can be exponential. The ODEs (8) can be written in a matrix form as:
(9)Ẋα(t)=AαXα(t)+Bα,Xα(0)=X0α.
By the variation of constants formula for ordinary differential equations, the solution of the initial value problem (9) is Xα(t)=eAαtX0α+∫0teAα(t-τ)Bαdτ. Since the six-dimensional matrix Aα depends on α, the calculation of eAαt becomes difficult. But this matrix can be written as Aα=A1+(1-α)A2, where A1 and A2 are α-independent matrices. We have eAαt=eA1te(1-α)A2t+O(t), where O(t) is a function that limt→0O(t)/t=0. Assuming that ℓjs and ςjs,j=1,…,6, are eigenvalues of A1 and A2, respectively, there are invertible matrices P and Q such that A1=PD1P-1 and A2=QD2Q-1, where D1=diag(ℓ1,…,ℓ6) and D2=diag(ς1,…,ς6). Therefore, for small t,(10)eAαt≈P(el1t0⋯00el2t⋯⋮⋮⋮⋱00⋯0el6t)P-1×Q(e(1-α)ζ1t0⋯00e(1-α)ζ2t⋯⋮⋮⋮⋱00⋯0e(1-α)ζ6t)Q-1.
As a result, the fitting functions are chosen as e(ℓi+(1-α)ςj)t, i,j=1,…,6. Therefore, by choosing the numbers ℓj, j=1,…,n, and ςi, i=1,…,m, where n,m∈{1,2,…,6}, an approximate solution can be found in the following form:(11)Sα(t)≈KEα(t),
where Sα(t)=[x̅α(t),x̲α(t),v̅α(t),v̲α(t)]T, Eα(t)=[e(ℓ1+(1-α)ς1)t,…,e(ℓ1+(1-α)ςm)t,…,e(ℓn+(1-α)ς1)t,…,e(ℓn+(1-α)ςm)t]T and K=[kij]4×nm denotes the coefficients matrix that can be found using the lsqnonlin code of the optimization toolbox in MATLAB. Obviously, z̅α(t) and z̲α(t) are obtained from the 5th and the 6th equations in (8) using the approximate values x̅α(t), x̲α(t), v̅α(t), and v̲α(t).


### 3.2. Dynamic Behavior of the Immune Cells Level and the Viral Load in PatientW

Figure5 shows the level of immune cells and the HIV viral load of patient W during the time interval [0, 1800]. The darker color shows the curve with the higher possibility. Figure 5(a) shows that the gradual declines in CD4+ T-cells level correspond to the low possibilities, while the rapid declines and the progression to full blown AIDS after a gradual decay, have the high possibilities of occurrence. Figure 5(b) shows a rapid increase in the viral load during the later stages of the disease. Moreover, we observe that, at each time, a higher viral load corresponds to a higher possibility. Figure 5(c) shows a clear correlation between CTLs level in the blood and HIV progression. As the viral load increases upon initial infection, CTLs increase in order to decrease the virus. But ultimately the level of these cells begins to decrease, which is due to virus-induced impairment of CD4+ T-cell function, with the high possibilities after about the 1450th day. Besides, a lower CTLs level has a higher possibility of occurrence and vice versa. From Figures 5(a) and 5(b), there is an inverse correlation between the HIV viral load and the level of CD4+ T-cells. Following the proposed method in Section 3.1, choose the eigenvalues of the corresponding matrices A1 and A2 as
(12)(ς1,ς2)=(-0.2600×10-3,0.6537×10-11),(l1,l2,l3,l4,l5)=(-0.2801×10-2,0.2801×10-2,-0.2145×10-6,-0.1128×10-6,0.6760×10-2).
Then, we have an approximate explicit solution in the form of (11) where the corresponding coefficients matrix K is(13)K=[0.668410.1635-0.22450.2229-0.710045.2900-0.702845.2924-0.00090.0005-3.140613.97790.1989-0.19740.676343.89730.683243.89940.0004-0.0008-0.01230.0109-0.00020.00020.0044-0.00380.0044-0.0038-0.00000.00030.8273-0.82220.0219-0.0217-0.16240.1608-0.16230.16090.0004-0.0001].CD4+ T-cells level (a), the viral load (b), and CTLs level (c) versus time in patientW.
(a)(b)(c)A comparison between the exact and approximate solutions which is shown in Figure6 confirmed the effectiveness of this approach.Exact (—) and approximate (×) CD4+ T-cells level (a) and the viral load (b) in patientW.
(a)(b)


### 3.3. Dynamic Behavior of the Immune Cells Level and the Viral Load in PatientS

Figure7 shows the changes in the immune cells level and the viral load in patient S. With respect to Figure 7(a), CD4+ T-cells level decreases gradually during the 1800 days from infection with the high possibilities. Moreover, an increase in CD4+ T-cells as well as the progression to AIDS arises with a low possibility. The HIV viral load is low and a lower viral load has a higher possibility of occurring, as shown in Figure 7(b). A high CD4+ T-cell count and a low HIV viral load lead to establishment of a lasting CTL response which is shown in Figure 7(c). A high HIV viral load and a low CD4+ T-cells level impair the immune response where this arises with the low possibilities as shown in Figure 7 as the light curves. The selected numbers ℓis and ςjs and the corresponding coefficients matrix K representing the approximate explicit solutions in the form of (11) are as(14)(ς1,ς2)=(0.2600×10-3,0.6537×10-11),(l1,l2,l3,l4,l5)=(-0.2800×10-2,0.2800×10-2,-0.1159×10-6,0.6240×10-2,-0.2378×10-6),K=[-3.320313.91750.1287-0.12650.766243.92930.0003-0.00060.772343.94120.83039.7266-0.14310.1412-0.806045.5176-0.00080.0005-0.799645.53000.3895-0.39020.0110-0.0109-0.07240.07250.0004-0.0002-0.07240.07260.0076-0.00990.0002-0.0002-0.00510.0063-0.00000.0002-0.00510.0063].CD4+ T-cells level (a), the viral load (b), and CTLs level (c) versus time in patientS.
(a)(b)(c)


### 3.4. Dynamic Behavior of the Immune Cells and the Viral Load in PatientM

Figure8 shows the changes in the immune cells level and the viruses in patient M. From Figure 8(a), the uncertainty of CD4+ T-cells level increases and new possibilities, varying fromm an increment to normal level to rapid progression to full blown AIDS, arise after a gradual decay. But the most possible scenario is between where the level of these cells decreases at a moderate rate. A moderate viral load occurs with a high possibility, as shown in Figure 8(b). A moderate CD4+ T-cells level as well as a moderate viral load implies a moderate CTLs level, as shown in Figure 8(c). Besides, from this figure, a high (low) viral load and a low (high) CD4+ T-cells number decrease (increase) CTLs level, and this happens with a low possibility.CD4+ T-cells level (a), the viral load (b), and CTLs level (c) versus time in patientM.
(a)(b)(c)Here, a representation for solutions is given in the form of (11) where(15)(l1,l2,l3,l4,l5)=(-0.2800×10-2,0.2800×10-2,0.2257×10-6,-0.1143×10-6,0.6500×10-2),(ς1,ς2,ς3)=(0,-0.2600×10-3,0.2600×10-3),K=[45.7977-20.3001-14.92030.8002-0.6250-0.177226.51149.60398.573426.51889.61698.5785-0.0004-0.00060.000625.4227-12.6502-2.18780.4622-0.0933-0.370531.24918.39455.039631.25498.40975.0415-0.00000.0005-0.0008-0.2693-0.30820.5693-0.0159-0.00160.0174-0.00680.0869-0.0772-0.00690.0869-0.0773-0.00030.00010.0005-0.05170.2515-0.2007-0.00350.0074-0.0040-0.0033-0.03630.0406-0.0033-0.03630.0406-0.00020.00040.0000].


### 3.5. A Comparison between the Immune Cells Level and the Viral Load in PatientsW, M, and S

From Figures5(b), 7(b), and 8(b), there is an inverse correlation between the viral load and the immune system strength. CD4+ T-cell and CTL levels in patient S are more than the level of these cells in patient W, as shown in Figures 5 and 7. A high viral load in patient W leads to a virus-induced impairment of CD4+ T-cell function. Therefore, CTLs level in this patient is less than the level of these cells in patient S, as shown in Figures 5(c) and 7(c). A comparison between CTLs level in patients M and S shows that it is possible that CTLs level in patient M be slightly more than the level of these cells in patient S, which can be due to this fact that patient M has more antigens that are required to stimulate CTLs (see Figures 7(b) and 8(b)), and, thus, a higher level of CD4+ T-cells in patient M is possible as shown in Figures 7(a) and 8(a). Figure 9 shows CD4+ T-cell and CTL levels and the viral load of patients W, M, and S on the 1800th day. Here, we perform a comparison between the immune cells level and the HIV viral load of patients based on their weighted center of gravity (WCOG).The fuzzy numbers indicating CD4+ T-cells level (a), the viral load (b), and CTLs level (c) on the 1800th day  in patientsW (- -), M (⋯) and, S (—).
(a)(b)(c)For this end, the WCOG ofx̃(t), ṽ(t), and z̃(t) is denoted by x(t), v(t), and z(t), respectively; hence,(16)x(t)=∫01α(x̅α(t)+x̲α(t))dα,v(t)=∫01α(v̅α(t)+v̲α(t))dα,z(t)=∫01α(z̅α(t)+z̲α(t))dα.
With respect to (11), an approximate explicit formula for the WCOG of x̃(t) and ṽ(t) is given as S(t)≈NF(t), where S(t)=[x(t),v(t)]Tand N is a matrix with two rows that its first and second rows are obtained by summing the first two rows and the last two rows of the corresponding coefficients matrix K, respectively. Moreover,(17)F(t)=∫01αEα(t)dα=[el1tς1t(eς1tς1t-1ς1t-1),…,el1tςmt(eςmtςmt-1ςmt-1),…,elntς1t(eς1tς1t-1ς1t-1),…,elntςmt(eςmtςmt-1ςmt-1)]T.
The WCOG of x̃(t), ṽ(t), and z̃(t) corresponding to patients W, M, and S is shown in Figure 10. From Figures 10(a) and 10(c), CD4+ T-cells and CTLs level are proportional to the strength of patient’s immune system such that a stronger immune system leads to a higher level of these cells. Moreover, there is an inverse correlation between the viral load and the strength of the immune system as shown in Figure 10(b).The WCOG of CD4+ T-cells level (a), the viral load (b) and CTLs level (c)  versus time in patientsW (- -), M (⋯), and S (—).
(a)(b)(c)


## 3.1. The Approximate Explicit Solutions Based on a Fitting Method

The proposed method is based on the fact that a linear combination of suitable functions ofα and t can generate the best fit to the exact values obtained by the ode45 in the least squares sense. The following discussion shows that these functions can be exponential. The ODEs (8) can be written in a matrix form as:
(9)Ẋα(t)=AαXα(t)+Bα,Xα(0)=X0α.
By the variation of constants formula for ordinary differential equations, the solution of the initial value problem (9) is Xα(t)=eAαtX0α+∫0teAα(t-τ)Bαdτ. Since the six-dimensional matrix Aα depends on α, the calculation of eAαt becomes difficult. But this matrix can be written as Aα=A1+(1-α)A2, where A1 and A2 are α-independent matrices. We have eAαt=eA1te(1-α)A2t+O(t), where O(t) is a function that limt→0O(t)/t=0. Assuming that ℓjs and ςjs,j=1,…,6, are eigenvalues of A1 and A2, respectively, there are invertible matrices P and Q such that A1=PD1P-1 and A2=QD2Q-1, where D1=diag(ℓ1,…,ℓ6) and D2=diag(ς1,…,ς6). Therefore, for small t,(10)eAαt≈P(el1t0⋯00el2t⋯⋮⋮⋮⋱00⋯0el6t)P-1×Q(e(1-α)ζ1t0⋯00e(1-α)ζ2t⋯⋮⋮⋮⋱00⋯0e(1-α)ζ6t)Q-1.
As a result, the fitting functions are chosen as e(ℓi+(1-α)ςj)t, i,j=1,…,6. Therefore, by choosing the numbers ℓj, j=1,…,n, and ςi, i=1,…,m, where n,m∈{1,2,…,6}, an approximate solution can be found in the following form:(11)Sα(t)≈KEα(t),
where Sα(t)=[x̅α(t),x̲α(t),v̅α(t),v̲α(t)]T, Eα(t)=[e(ℓ1+(1-α)ς1)t,…,e(ℓ1+(1-α)ςm)t,…,e(ℓn+(1-α)ς1)t,…,e(ℓn+(1-α)ςm)t]T and K=[kij]4×nm denotes the coefficients matrix that can be found using the lsqnonlin code of the optimization toolbox in MATLAB. Obviously, z̅α(t) and z̲α(t) are obtained from the 5th and the 6th equations in (8) using the approximate values x̅α(t), x̲α(t), v̅α(t), and v̲α(t).


## 3.2. Dynamic Behavior of the Immune Cells Level and the Viral Load in PatientW

Figure5 shows the level of immune cells and the HIV viral load of patient W during the time interval [0, 1800]. The darker color shows the curve with the higher possibility. Figure 5(a) shows that the gradual declines in CD4+ T-cells level correspond to the low possibilities, while the rapid declines and the progression to full blown AIDS after a gradual decay, have the high possibilities of occurrence. Figure 5(b) shows a rapid increase in the viral load during the later stages of the disease. Moreover, we observe that, at each time, a higher viral load corresponds to a higher possibility. Figure 5(c) shows a clear correlation between CTLs level in the blood and HIV progression. As the viral load increases upon initial infection, CTLs increase in order to decrease the virus. But ultimately the level of these cells begins to decrease, which is due to virus-induced impairment of CD4+ T-cell function, with the high possibilities after about the 1450th day. Besides, a lower CTLs level has a higher possibility of occurrence and vice versa. From Figures 5(a) and 5(b), there is an inverse correlation between the HIV viral load and the level of CD4+ T-cells. Following the proposed method in Section 3.1, choose the eigenvalues of the corresponding matrices A1 and A2 as
(12)(ς1,ς2)=(-0.2600×10-3,0.6537×10-11),(l1,l2,l3,l4,l5)=(-0.2801×10-2,0.2801×10-2,-0.2145×10-6,-0.1128×10-6,0.6760×10-2).
Then, we have an approximate explicit solution in the form of (11) where the corresponding coefficients matrix K is(13)K=[0.668410.1635-0.22450.2229-0.710045.2900-0.702845.2924-0.00090.0005-3.140613.97790.1989-0.19740.676343.89730.683243.89940.0004-0.0008-0.01230.0109-0.00020.00020.0044-0.00380.0044-0.0038-0.00000.00030.8273-0.82220.0219-0.0217-0.16240.1608-0.16230.16090.0004-0.0001].CD4+ T-cells level (a), the viral load (b), and CTLs level (c) versus time in patientW.
(a)(b)(c)A comparison between the exact and approximate solutions which is shown in Figure6 confirmed the effectiveness of this approach.Exact (—) and approximate (×) CD4+ T-cells level (a) and the viral load (b) in patientW.
(a)(b)


## 3.3. Dynamic Behavior of the Immune Cells Level and the Viral Load in PatientS

Figure7 shows the changes in the immune cells level and the viral load in patient S. With respect to Figure 7(a), CD4+ T-cells level decreases gradually during the 1800 days from infection with the high possibilities. Moreover, an increase in CD4+ T-cells as well as the progression to AIDS arises with a low possibility. The HIV viral load is low and a lower viral load has a higher possibility of occurring, as shown in Figure 7(b). A high CD4+ T-cell count and a low HIV viral load lead to establishment of a lasting CTL response which is shown in Figure 7(c). A high HIV viral load and a low CD4+ T-cells level impair the immune response where this arises with the low possibilities as shown in Figure 7 as the light curves. The selected numbers ℓis and ςjs and the corresponding coefficients matrix K representing the approximate explicit solutions in the form of (11) are as(14)(ς1,ς2)=(0.2600×10-3,0.6537×10-11),(l1,l2,l3,l4,l5)=(-0.2800×10-2,0.2800×10-2,-0.1159×10-6,0.6240×10-2,-0.2378×10-6),K=[-3.320313.91750.1287-0.12650.766243.92930.0003-0.00060.772343.94120.83039.7266-0.14310.1412-0.806045.5176-0.00080.0005-0.799645.53000.3895-0.39020.0110-0.0109-0.07240.07250.0004-0.0002-0.07240.07260.0076-0.00990.0002-0.0002-0.00510.0063-0.00000.0002-0.00510.0063].CD4+ T-cells level (a), the viral load (b), and CTLs level (c) versus time in patientS.
(a)(b)(c)


## 3.4. Dynamic Behavior of the Immune Cells and the Viral Load in PatientM

Figure8 shows the changes in the immune cells level and the viruses in patient M. From Figure 8(a), the uncertainty of CD4+ T-cells level increases and new possibilities, varying fromm an increment to normal level to rapid progression to full blown AIDS, arise after a gradual decay. But the most possible scenario is between where the level of these cells decreases at a moderate rate. A moderate viral load occurs with a high possibility, as shown in Figure 8(b). A moderate CD4+ T-cells level as well as a moderate viral load implies a moderate CTLs level, as shown in Figure 8(c). Besides, from this figure, a high (low) viral load and a low (high) CD4+ T-cells number decrease (increase) CTLs level, and this happens with a low possibility.CD4+ T-cells level (a), the viral load (b), and CTLs level (c) versus time in patientM.
(a)(b)(c)Here, a representation for solutions is given in the form of (11) where(15)(l1,l2,l3,l4,l5)=(-0.2800×10-2,0.2800×10-2,0.2257×10-6,-0.1143×10-6,0.6500×10-2),(ς1,ς2,ς3)=(0,-0.2600×10-3,0.2600×10-3),K=[45.7977-20.3001-14.92030.8002-0.6250-0.177226.51149.60398.573426.51889.61698.5785-0.0004-0.00060.000625.4227-12.6502-2.18780.4622-0.0933-0.370531.24918.39455.039631.25498.40975.0415-0.00000.0005-0.0008-0.2693-0.30820.5693-0.0159-0.00160.0174-0.00680.0869-0.0772-0.00690.0869-0.0773-0.00030.00010.0005-0.05170.2515-0.2007-0.00350.0074-0.0040-0.0033-0.03630.0406-0.0033-0.03630.0406-0.00020.00040.0000].


## 3.5. A Comparison between the Immune Cells Level and the Viral Load in PatientsW, M, and S

From Figures5(b), 7(b), and 8(b), there is an inverse correlation between the viral load and the immune system strength. CD4+ T-cell and CTL levels in patient S are more than the level of these cells in patient W, as shown in Figures 5 and 7. A high viral load in patient W leads to a virus-induced impairment of CD4+ T-cell function. Therefore, CTLs level in this patient is less than the level of these cells in patient S, as shown in Figures 5(c) and 7(c). A comparison between CTLs level in patients M and S shows that it is possible that CTLs level in patient M be slightly more than the level of these cells in patient S, which can be due to this fact that patient M has more antigens that are required to stimulate CTLs (see Figures 7(b) and 8(b)), and, thus, a higher level of CD4+ T-cells in patient M is possible as shown in Figures 7(a) and 8(a). Figure 9 shows CD4+ T-cell and CTL levels and the viral load of patients W, M, and S on the 1800th day. Here, we perform a comparison between the immune cells level and the HIV viral load of patients based on their weighted center of gravity (WCOG).The fuzzy numbers indicating CD4+ T-cells level (a), the viral load (b), and CTLs level (c) on the 1800th day  in patientsW (- -), M (⋯) and, S (—).
(a)(b)(c)For this end, the WCOG ofx̃(t), ṽ(t), and z̃(t) is denoted by x(t), v(t), and z(t), respectively; hence,(16)x(t)=∫01α(x̅α(t)+x̲α(t))dα,v(t)=∫01α(v̅α(t)+v̲α(t))dα,z(t)=∫01α(z̅α(t)+z̲α(t))dα.
With respect to (11), an approximate explicit formula for the WCOG of x̃(t) and ṽ(t) is given as S(t)≈NF(t), where S(t)=[x(t),v(t)]Tand N is a matrix with two rows that its first and second rows are obtained by summing the first two rows and the last two rows of the corresponding coefficients matrix K, respectively. Moreover,(17)F(t)=∫01αEα(t)dα=[el1tς1t(eς1tς1t-1ς1t-1),…,el1tςmt(eςmtςmt-1ςmt-1),…,elntς1t(eς1tς1t-1ς1t-1),…,elntςmt(eςmtςmt-1ςmt-1)]T.
The WCOG of x̃(t), ṽ(t), and z̃(t) corresponding to patients W, M, and S is shown in Figure 10. From Figures 10(a) and 10(c), CD4+ T-cells and CTLs level are proportional to the strength of patient’s immune system such that a stronger immune system leads to a higher level of these cells. Moreover, there is an inverse correlation between the viral load and the strength of the immune system as shown in Figure 10(b).The WCOG of CD4+ T-cells level (a), the viral load (b) and CTLs level (c)  versus time in patientsW (- -), M (⋯), and S (—).
(a)(b)(c)


## 4. Fuzzy Optimal Control Problem

In this section, we formulate a fuzzy optimal control problem that identifies the parameterk̃ in (8), with a function of the fuzzy control variable ũ. In particular, we will replace the parameter k̃ with the function k̃ΘHũ. This choice then identifies the control variable ũ(t) with the rate of inhibition of virus reproduction, which is modeled as a simple function of drug dosage. Therefore, we have the fuzzy set control differential equations (FSCDEs) as(18)x̃̇=λ̃Θσ̃⊗x̃Θc̃⊗ṽ,ṽ̇=(k̃ΘHũ)⊗ṽΘã⊗z̃,z̃̇=h̃⊗x̃Θτ̃⊗ṽ.
This paper aims to propose a drug regimen that minimizes both the viral load and the drug costs. Here, we assume that the cost of the treatment is proportional to ũ2(t) at time t. Therefore, the fuzzy functional J̃(ṽ,ũ)=∫t0tf(w̃⊗ṽ(t)⊕ũ2(t))dt should be minimized, where the positive fuzzy number w̃ is used to set the relative importance between minimizing the viral load and the systemic cost to the body. Let T=[t0,tf], and assume C(T) be the set of all continuous fuzzy number valued functions on T. Assuming that the minimum and the maximum of allowable drug dosage are denoted by fuzzy numbers ũmin and ũmax, respectively, then we are seeking a ũ*∈Ũ such that J̃(ṽ*,ũ*)≤̃J̃(ṽ,ũ), for all ũ∈Ũ, where Ũ={ũ∈C(T):ũmin≤̃ũ(t)≤̃ũmax, for all t∈T}. By using the derivative in the first form (i), the FSCDEs (18) is converted to the following control system:(19)x̅̇α(t)=λ̅α-σ̲αx̲α(t)-c̲αv̲α(t),ẋ̲α(t)=λ̲α-σ̅αx̅α(t)-c̅αv̅α(t),v̅̇α(t)=(k̅α-u̅α(t))v̅α(t)-a̲αz̲α(t),v̲̇α(t)=(k̲α-u̲α(t))v̲α(t)-a̅αz̅α(t),z̅̇α(t)=h̅αx̅α(t)-τ̲αv̲α(t),ż̲α(t)=h̲αx̲α(t)-τ̅αv̅α(t),x̅α(0)=x̅0α,x̲α(0)=x̲0α,v̅α(0)=v̅0α,v̲α(0)=v̲0α,z̅α(0)=z̅0α,z̲α(0)=z̲0α.
With respect to Definition 5, a fuzzy function ũ*∈Ũ is viewed as an optimal solution, if, for each α∈[0,1], the pair (u̲α*,u̅α*)∈U minimizes the functionals [∫t0tf(w̃⊗ṽ(t)⊕ũ2(t))dt]¯α=∫t0tf(w̅αv̅α(t)+u̅α2(t))dt and [∫t0tf(w̃⊗ṽ(t)⊕ũ2(t))dt]̲α=∫t0tf(w̲αv̲α(t)+u̲α2(t))dt simultaneously, where v̅α(t) and v̲α(t) are the solution of ODEs (19) corresponding to control pair (u̲α,u̅α), and U is the set of all measurable control pairs (u̲α,u̅α) that u̲αmin≤u̲α(t)≤u̲αmax,u̅αmin≤u̅α(t)≤u̅αmax, for all t∈T. Therefore, we restrict our attention to optimizing the functional J(v̅α,v̲α,u̅α,u̲α)=∫t0tf(w̲αv̲α(t)+w̅αv̅α(t)+u̲α2(t)+u̅α2(t))dt over the set U. We now proceed to compute candidates for an optimal solution by applying the Pontryagin’s Maximum Principle [24] and begin by defining the Lagrangian to be
(20)Lα(x̲α,x̅α,v̲α,v̅α,z̲α,z̅α,u̲α,u̅α,λ1α,…,λ6α)=w̲αv̲α(t)+w̅αv̅α(t)+u̲α2+u̅α2+λ1α(λ̅α-σ̲αx̲α-c̲αv̲α)+λ2α(λ̲α-σ̅αx̅α-c̅αv̅α)+λ3α((k̅α-u̅α)v̅α-a̲αz̲α)+λ4α((k̲α-u̲α)v̲α-a̅αz̅α)+λ5α(h̅αx̅α-τ̲αv̲α)+λ6α(h̲αx̲α-τ̅αv̅α)-ω11(u̲α-u̲αmin)-ω12(u̲αmax-u̲α)-ω21(u̅α-u̅αmin)-ω22(u̅αmax-u̅α),
where ωij≥0 are the penalty multipliers satisfying ω11(u̲α*-u̲αmin)=ω12(u̲αmax-u̲α*)=0 and ω21(u̅α*-u̅αmin)=ω22(u̅αmax-u̅α*)=0. Thus, the Maximum Principle gives the existence of adjoint variables λjα, j=1,…,6, satisfying
(21)λ̇1α=-∂Lα∂x̅α=σ̅αλ2α-h̅αλ5α,λ̇2α=-∂Lα∂x̲α=σ̲αλ1α-h̲αλ6α,λ̇3α=-∂Lα∂v̅α=-w̅α+c̅αλ2α-(k̅α-u̅α)λ3α+τ̅αλ6α,λ̇4α=-∂Lα∂v̲α=-w̲α+c̲αλ1α-(k̲α-u̲α)λ4α+τ̲αλ5α,λ̇5α=-∂Lα∂z̅α=a̅αλ4α,λ̇6α=-∂Lα∂z̲α=a̲αλ3α,
where λjα(tf)=0, j=1,…,6, are the transversality conditions. The Lagrangian is minimized with respect to u̲α and u̅α at the optimal pair (u̲α*,u̅α*). So the partial derivatives of the Lagrangian with respect to u̲α and u̅α are zero. Since, ∂L/∂u̲α=2u̲α-λ4αv̲α-ω11+ω12=0, we have u̲α=0.5(v̲αλ4α+ω11-ω12). To determine an explicit expression for the optimal control without ω11 and ω12, we consider the following three cases. (i) On the set {t∣u̲αmin<u̲α*(t)<u̲αmax}, we set ω11(t)=ω12(t)=0; hence, u̲α=0.5v̲αλ4α. (ii) On the set {t∣u̲α*(t)=u̲αmax}, we set ω11(t)=0; hence, u̲α=u̲αmax=0.5(v̲αλ4α-ω12) which implies that 0.5v̲αλ4α≥u̲αmax. (iii) On the set {t∣u̲α*(t)=u̲αmin}, we set ω12(t)=0; hence, u̲α=u̲αmin=0.5(v̲αλ4α+ω11) which implies that 0.5v̲αλ4α≤u̲αmin. Combining all the three cases in compact form gives(22)u̲α(t)=max⁡(u̲αmin⁡,min⁡(0.5v̲α(t)λ4α(t),u̲αmax⁡)).
Using similar arguments, we also obtain the following expression for the second optimal control function:(23)u̅α(t)=max⁡(u̅αmin⁡,min⁡(0.5v̅α(t)λ3α(t),u̅αmax⁡)).
We point out that the optimality system consists of the state system (19) with the initial conditions, the adjoint or costate system (21) with the terminal conditions, together with the expressions (22) and (23) for the control functions. We show the optimal controls and corresponding states and costates satisfying (19) and (21) by a star superscript ∗. Obviously, if we assume that the set [λ4α*(t),λ3α*(t)]is a valid α-level set of a positive fuzzy number valued function, say Λ̃2*(t), then, from (22) and (23),the optimal fuzzy control function ũ*∈C(T) can be written as(24)ũ*(t)=max⁡̃(ũmin⁡,min⁡̃(0.5⊗ṽ*(t)⊗Λ̃2*(t),ũmax⁡)).
Moreover, if we assume that the sets [λ2α*(t),λ1α*(t)]and [λ6α*(t),λ5α*(t)]are valid α-level sets of fuzzy number valued functions, say Λ̃1*(t) and Λ̃3*(t), respectively, and the sets [λ̇1α*(t),λ̇2α*(t)], [λ̇3α*(t),λ̇4α*(t)], and [λ̇5α*(t),λ̇6α*(t)]are valid α-level sets, then it is easy to see that (ΘΛ̃2*)⊗Zã and Λ̃3*⊗Zh̃∈F(R). Moreover, we can verify that the fuzzy functions Λ̃1*(t), Λ̃2*(t), Λ̃3*(t), and ũ* satisfy the following system of FDEs, with the terminal conditions Λ̃j(tf)=(0,0), j=1,2,3, and using the derivative in the second form (ii) (see Definition 9):
(25)Λ̃̇1=Θ{(ΘΛ̃1)⊗σ̃⊕Λ̃3⊗Zh̃},Λ̃̇2=Θ{w̃⊕(ΘΛ̃1)⊗c̃⊕Λ̃2⊗(k̃ΘHũ)⊕(ΘΛ̃3)⊗τ̃},Λ̃̇3=Θ{(ΘΛ̃2)⊗Zã}.
Therefore, if we insure that the assumptions mentioned above are satisfied, finding the optimal fuzzy control function can be equivalent to solving a two-point fuzzy boundary value problem (FBVP) which consists of fuzzy system (18) with the initial conditions and the fuzzy system (25) with the final conditions, together with the expression ũ(t)=max̃(ũmin,miñ(0.5⊗ṽ(t)⊗Λ̃2(t),ũmax)) for the control function. Here, we verify that Λ̃1*(t), Λ̃2*(t), Λ̃3*(t), and ũ* satisfy the second equation in (25) and the rest can be verified similarly. Since the functions λ3α* and λ4α* are positive from the two last equations in (21), we conclude that λ5α* and λ6α* are negative functions. Moreover, from the two first equations in (21), we have(26)λ1α*(t)=∫ttf{σ̅ασ̲αh̅αcosh⁡{σ̅ασ̲α(τ-t)}λ5α*(τ)+σ̅ασ̲αh̲αsinh⁡{σ̅ασ̲α(τ-t)}λ6α*(τ)}dτ,λ2α*(t)=∫ttf{σ̲ασ̅αh̅αsinh⁡{σ̅ασ̲α(τ-t)}λ5α*(τ)+h̲αcosh⁡{σ̅ασ̲α(τ-t)}λ6α*(τ)σ̲ασ̅α}dτ,
which implies that the functions λ1α* and λ2α* are negative; hence,(27)Θ{w̃⊕(ΘΛ̃1*)⊗c̃⊕Λ̃2*⊗(k̃ΘHũ*)⊕(ΘΛ̃3*)⊗τ̃}̲α=-{w̃⊕(ΘΛ̃1*)⊗c̃⊕Λ̃2*⊗(k̃ΘHũ*)⊕(ΘΛ̃3*)⊗τ̃}¯α=-{w̅α-λ2α*c̅α+λ3α*(k̅α-u̅α*)-λ6α*τ̅α}=λ̇3α*.
We summarize our results in the following theorem.Theorem 13.
Assume thatΛ̃2*(t) is a positive fuzzy function on T. A candidate for fuzzy optimal control is given by (24), if the fuzzy functions ũ*(t), x̃*(t), z̃*(t), and ṽ*(t) satisfy (18) using the derivative in the first form (i) and ũ*(t), Λ̃1*(t),  Λ̃3*(t), and Λ̃2*(t) satisfy (25) using the derivative in the second form (ii).We have solved the optimality system corresponding to patientW with the fuzzy weight w̃=(0.24×10-5,0.24×10-6) and the fuzzy bounds ũmin=(0,0,0.5×10-3) and ũmax=(0.15×10-2,0.25×10-3,0), during the time interval [0, 1800] by using the gradient method [24]. Numerical results show that the optimal adjoint variables satisfy the assumptions mentioned above. Figure 11 shows the optimal fuzzy adjoint variables, while the optimal fuzzy control function ũ* and its WCOG as a real-valued output u* indicating an optimal drug regimen are depicted in Figure 12. Moreover, the optimal fuzzy states x̃*(t), z̃*(t), and ṽ*(t) indicating the immune cells level and the HIV viral load in presence of treatment are shown in Figure 13. The values of ũ* on the 425th and the 1750th days are, respectively, shown in Figures 14 and 15, for the sake of clarity.Optimal fuzzy adjoint variablesΛ̃1*(a), Λ̃2*(b), and Λ̃3*(c) versus time.
(a)(b)(c)Optimal fuzzy control (a) and its WCOG (b) versus time.
(a)(b)CD4+ T-cells level (a), the viral load (b), and CTLs level (c) in presence of treatment versus time in patientW.
(a)(b)(c)Figure 14
ũmax(- -), 0.5⊗ṽ*(425)⊗Λ̃2*(425) (–·–), ũ*(425)=miñ(0.5⊗ṽ*(425)⊗Λ̃2*(425),ũmax) (—).Figure 15
ũmin (- -), 0.5⊗ṽ*(1750)⊗Λ̃2*(1750) (–·–), ũ*(1750)=max̃(0.5⊗ṽ*(1750)⊗Λ̃2*(1750),ũmin) (—).The proposed treatment by Figure12, which its intensity is decreasing during the time interval [0, 1800], reduces the proliferation rate of viruses considerably. Therefore, the HIV viral load is very low as shown in Figure 13(b). A very low viral load slows the destruction of CD4+ T-cells which are due to contact with virus particles; hence, a CD4+ T-cells level below 70% is not possible in a treated patient as shown in Figure 13(a). From Figure 13(c), a high level of CD4+ T-cells and a low viral load lead to establishment of a lasting CTL-mediated immune response.


## 5. Conclusion

In this paper, we proposed a fuzzy mathematical model of HIV dynamic. Simulation results show that the proposed three-dimensional FDEs can describe the ambiguous immune cells level and the HIV viral load which are due to existing patients with various strength of their immune system. Moreover, we utilized the proposed fuzzy model and studied a fuzzy optimal control problem minimizing both the viral load and drug costs. Using the Pontryagin’s Maximum Principle leads us to the conclusion that the fuzzy optimal control function may not exist in general, but an optimality system containing fuzzy state and fuzzy adjoin equations is derived under certain assumptions. Motivated by these results, we tend to exploit necessary and sufficient conditions to the existence of fuzzy solutions for linear fuzzy optimal control problems. We expect to address these problems in further works.



---

*Source: 893474-2012-03-22.xml*