# Global Dynamics of a Virus-Immune System with Virus-Guided Therapy and Saturation Growth of Virus

Authors: Qian Li; Yanni Xiao
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

We considered a piecewise virus-immune dynamic model to investigate the effectiveness of the HIV virus loads-guided structured treatment interruptions (STIs). To better describe the biological reality, we extended the existing models by taking the carrying capacity of the virus loads into consideration to indicate the saturated growth of virus loads. We initially investigated the sliding dynamics of the proposed model and then obtained the global dynamics of the proposed model. Our main results showed that the system can exhibit very complex and diverse dynamic behaviors including a globally asymptotically stable equilibrium, bistable equilibra, and tristable equilibria, depending on the dynamics of the subsystems and the threshold level. In particular, an interesting result indicated that, with a proper threshold condition, the virus-guided therapy policy can successfully control the virus loads far below its carrying capacity and maintain the activity of the immune system for the case that the effector cells always go to zero without therapy or with continuous therapy. The finding suggested that the optimal strategy should be individual-based due to coexistence of multiple stable steady states, depending on the threshold conditions and the initial levels of viral loads and effector cells of the patients.


---

## Body

## 1. Introduction

Highly active antiretroviral therapy (HAART) plays significant roles in improving survival and reducing morbidity and mortality in HIV patients [1, 2]. Nevertheless, lifelong HAART confronts many complicated problems such as adherence difficulties and the evolution of drug resistance [3–6]. Structured treatment interruptions (STIs), the alternative strategies, have been designed to effectively achieve sustained specific immunity for early therapy in HIV infection. Indeed, STIs are beneficial for some chronically infected individuals who may need to take drugs throughout their lives as well as the reconstruction of patients’ immune system when they are not taking the drugs [7].Since the first HIV patient received the structured antiretroviral therapy in 1999 [8], STIs have been widely investigated by many researchers. Basically, it includes two kinds of treatment regimes: carrying out the therapies at fixed moments or depending on the CD4 cell counts and/or viral loads. For example, some studies designed a scheduled treatment that therapy was initiated only when the CD4 cell counts decreased below 350 cells/μl [9–12]. Aiming to maintain CD4 cell counts higher than 350 cells/μl and HIV virus loads less than 100,000 copies/μl, Ruiz et al. designed an experiment to evaluate the safety of CD4 cell and HIV virus-guided STIs [13]. Many works have been done to compare the STIs with the continuous antiretroviral therapy (ART), but conflicting results have been reported [7, 11, 13–17]. In particular, Maggiolo et al. made a conclusion that the structured treatment can result in the similar benefits for the patients with the continuous therapies [7], while the SMART study group found that STIs may also endanger the patients’ life if their immune responses are not sufficiently stimulated [18]. STI is still an important topic since lifelong ART brings great challenges to all HIV positive individuals. Moreover, as immune therapy arises, combining the ART with immune therapy (such as interleukin-2 (IL-2) treatment) may be a good choice for treating HIV patients [19]. Modelling and quantifying the efficacy of the combination of STIs and immune therapy remains unclear and falls within the scope of this study.Mathematical models are useful tools for investigating the dynamics of HIV-1 infection [20, 21]. Although many studies modelled the continuous therapy for HIV-infected patients [22–24], few attempts studied the STIs through mathematical modelling. In 2000, using a standard population model, Bonhoeffer et al. argued the benefits and risk generated from STIs [25]. In the study [26], Adams et al. analyzed the efficacy of STIs with fixed moments through the optimal control theory. Further, Rosenberg et al. pointed out that it is needed to include drug resistance and CD4 cell counts when making the structured treatment [27]. In 2006, Park et al. developed a mathematical model to investigate the impact of finite times of structured treatment guided by CD4 cell counts on HIV patients [28]. Then, Tang et al. [29] formulated a piecewise system to describe the CD4 cell-guided STIs, quantitatively explored STI strategies, and explained some controversial conclusions from different clinical studies. To further model immune therapy, Tang et al. [30] proposed a piecewise virus-immune dynamic model considering the combined antiretroviral therapy with interleukin-2 (IL-2) treatment. The model is given as follows:(1)x˙=rx-pxy-ε1Φx,y˙=cxy1+ωx-qxy-δy+ε2Φy,with(2)Φ=0,ifHx=x-VT<0,1,ifHx=x-VT>0,where x and y represent the HIV virus loads and the effector cells, respectively. Note that a critical value of virus loads, denoted by VT, is set to trigger the combined therapy. That is, only when the virus loads exceed the threshold, we carry out the antiretroviral therapy and immune therapy (such as the interleukin-2 (IL-2) treatment) simultaneously for the patients, with ε1 being the elimination rate of HIV virus due to antiretroviral therapy and ε2 representing the growth rate of the effector cells due to immune therapy. Here, r is the growth rate of HIV virus which incorporates both their multiplication and death, δ is the death rate of the effector cells, p denotes the rate of binding of the effector cells to the virus loads, q represents the rate of inactivation of the effector cells, and cx/(1+ωx) denotes that effector cell multiplication due to immune response has a maximum value as the HIV virus loads become large sufficiently.In [30], the authors mainly discussed the global dynamics of system (1), and they further considered the efficacy of the effector cell-guided or both the effector cell and the HIV virus-guided piecewise therapies [31, 32]. Their results indicated that HIV virus can be controlled under some conditions with virus or effector cell-guided therapy, but still there are some parameter regions that the threshold policy can not prevent HIV from growing infinitely. As we can see from model (1), they simply assumed the linear growth of HIV virus. However, this is not what the thing looks like. Some clinical facts indicate that there should be “saturation effects” for the growth of HIV virus [33, 34]. The linear growth of HIV virus is unreasonable. Note that including the saturated growth of the HIV virus makes the model more natural but significantly changes the dynamics of the basic immune-virus system. Little is known about whether HIV virus-guided therapy based on this more natural description of its growth can completely control the growth of HIV and falls within the scope of this study.The main purpose of this paper is to propose a piecewise model with the saturated growth for virus loads to describe HIV virus-guided therapy and then evaluate the effectiveness of the threshold policy through analyzing the global dynamics of the proposed model. The remaining parts of this paper are organized as follows. In Section2, we initially propose the model, give some basic definitions of the Filippov system and also briefly conclude the dynamic behaviors of two subsystems. In Section 3, we mainly discuss the sliding mode and sliding dynamics, including the existence of the sliding domain and the pseudoequilibrium. In Section 4, we investigate the global dynamics of the proposed model. Finally, we provide the discussion and some biological conclusions.


## 2. Filippov Model and Preliminaries

In this paper, we extend model (1) by including the saturated growth of HIV virus, and hence we propose the following model:(3)x˙=rx1-ηx-pxy-ε1Φx,y˙=cxy1+ωx-qxy-δy+ε2Φy,with(4)Φ=0,ifHx=x-VT<0,1,ifHx=x-VT>0.Here, 1/η represents the carrying capacity of virus loads.LettingZ=(x,y)T, (5)FR1Z=rx1-ηx-pxy,cxy1+ωx-qxy-δyT,FR2Z=rx1-ηx-pxy-ε1x,cxy1+ωx-qxy-δy+ε2yT.Then system (3) with (4) can be rewritten as the following general planer Filippov system [35–37]:(6)Zt=FR1Z,Z∈R1,FR2Z,Z∈R2.with R1=Z∈R+2∣H(Z)<0 and R2=Z∈R+2∣H(Z)>0. Additionally, we describe the discontinuity boundary Σ separating the two regions R1 and R2 as Σ=Z∈R+2∣H(Z)=0, where H(Z) is a smooth scale function nonvanishing on Σ. For convenience, we denote Filippov system (6) defined in region R1 as system S1 and that defined in region R2 as system S2. Then we let(7)σZ=HZZ,FR1ZHZZ,FR2Z,where 〈·〉 represents the standard scalar product and HZ(Z)=(1,0)T denotes the nonvanishing gradient of H on Σ. In the following, we replace HZ(Z),FRi(Z) with the notation FRiH(Z) for i=1,2.In order to consider the trajectory of the Filippov vector fieldFR1(Z),FR2(Z), through a point Z∈Σ, we distinguish the following regions on Σ according to whether or not the vector field points towards it:(a) Crossing region:(8)Σc=Z∈Σ∣FR1HZ·FR2HZ>0,(b) Sliding region:(9)Σs=Z∈Σ∣FR1HZ>0,FR2HZ<0,(c) Escaping region:(10)Σe=Z∈Σ∣FR1HZ<0,FR2HZ>0.Further, we can define a scalar differential equation on the sliding domain by employing the Filippov method. Let(11)gZ=λFR1+1-λFR2,λ=HZ,FR2HZ,FR2-FR1with 0<λ<1. Then, we have that(12)Z′=gZ,Z∈Σsdescribes the dynamics of system (6) on the sliding segment Σs. In what follows, we give the definitions of two types of equilibria of Filippov system (6).Definition 1.
A pointZ∗ is referred to as a real equilibrium of system (6) if FR1(Z∗)=0,H(Z∗)<0, or FR2(Z∗)=0,H(Z∗)>0. A point Z∗ is called a virtual equilibrium of system (6) if FR1(Z∗)=0,H(Z∗)>0, or FR2(Z∗)=0,H(Z∗)<0. Both the real equilibrium and virtual equilibrium are called regular equilibria.Definition 2.
A pointZ∗ is referred to as a pseudoequilibrium if it is an equilibrium of the sliding mode of system (6); i.e., the equilibrium of the equation Z′=g(Z)=0 satisfies H(Z∗)=0.Before exploring the global dynamics of the Filippov system (3), it is essential to have a clear view of the dynamics of two basic subsystems. Thus, we first briefly make some conclusions on the dynamics of subsystem S1 and subsystem S2.For subsystemS1, there always exists a trivial equilibrium E10=(0,0) and a boundary equilibrium E1η=(x1η,0)=(1/η,0). The positive equilibrium (x∗,y∗) of subsystem S1 satisfies the following equations:(13)r1-ηx∗-py∗=0,cx∗1+ωx∗-qx∗-δ=0.Let(14)fx≐cx-qx+δ1+ωx=-qωx2+c-q-δωx-δ,and Δ1=c-q-δω-2qωδ. It is easy to verify that if Δ1>0 holds true, then f(x)=0 has two positive roots which are given as(15)x11=c-q-δω-c-q-δω2-4qωδ2qω,x12=c-q-δω+c-q-δω2-4qωδ2qω.Particularly, we know f(x)>0 for x∈(x11,x12) and f(x)<0 for x∈(0,x11) or x∈(x12,+∞). Substituting x11 or x12 into the first equation of (13), we get(16)y11=rp1-ηx11,y12=rp1-ηx12.Thus, if(17)c-q-δω∓c-q-δω2-4qωδ2qω<1ηi.e.,x1i<1η,i=1,2holds true, then E1i(x1i,y1i),  i=1,2, are positive equilibria of system S1.Next, we briefly investigate the stability of the equilibria of systemS1. E1η is locally stable if c/η/(1+ω/η)-q/η-δ<0, which is equivalent to Δ1<0 or Δ1>0, y11<0 or Δ1>0, and y12>0. Define a function K(x,y)=1/xy, and calculate(18)L=∂∂xKx,ydxdt+∂∂yKx,ydydt=-rηy,where dx/dt and dy/dt are defined in system S1. Since all the parameters are positive and solutions to system S1 initiating from R+2 are nonnegative, L does not change its sign in the first quadrant. Hence according to Dulac-Bendixson criterion, we get that there are no limit cycles or homoclinic connections for system S1. Therefore, if there is no positive equilibrium, then the boundary equilibrium E1η is globally asymptotically stable because system S1 is bounded. Further, if Δ1>0 and x11<1/η<x12 hold true, it is easy to verify that f(1/η)>0, and hence the boundary equilibrium E1η is unstable. Therefore, under this scenario, positive equilibrium E11 is globally asymptotically stable. When there are two positive equilibria of system S1, the dynamics of system S1 have been shown in [38] that the equilibrium E12 is always an unstable saddle and E11 and E1η are bistable. Then the dynamics of system S1 can be concluded as follows.Proposition 3.
SystemS1 always has a trivial equilibrium E10 and a boundary equilibrium E1η. If Δ1<0 or Δ1>0 and x11>1/η hold true, then there is no positive equilibrium, and the boundary equilibrium E1η is globally asymptotically stable; If Δ1>0 and x11<1/η<x12, there is one and only one positive equilibrium E11 which is globally asymptotically stable; If Δ1>0 and x12<1/η, there are two positive equilibria E11 and E12, and E11 and E1η are bistable.Ifr>ε1 and δ>ε2, system S2 has similar properties to system S1. Denote Δ2=c-q-(δ-ε2)ω-2qω(δ-ε2) and(19)x2i=c-q-δ-ε2ω∓c-q-δ-ε2ω2-4qωδ-ε22qω,y2i=rp1-ηx2i-ε1p,i=1,2.Similarly, if(20)c-q-δ-ε2ω∓c-q-δ-ε2ω2-4qωδ-ε22qω<r-ε1rηi.e.,x2i<r-ε1rη,i=1,2holds true, then E2i(x2i,y2i), i=1,2, are positive equilibria of system S2. Then we can conclude the dynamics of system S2 as follows.Proposition 4.
SystemS2 has a trivial equilibrium E20=(0,0) and a boundary equilibrium E2η=(x2η,0)=((r-ε1)/rη,0). If Δ2<0 or Δ2>0 and x21>(r-ε1)/rη, system S2 has no positive equilibrium, and the boundary equilibrium E2η is globally asymptotically stable; if Δ2>0 and x22<(r-ε1)/rη<x21 hold true, there is one and only one positive equilibrium E21 of system S2 which is globally asymptotically stable; if Δ2>0 and x22<(r-ε1)/rη hold true, there are two positive equilibria E21 and E22, and both E21 and E2η are bistable.


## 3. Sliding Dynamics

According to the definition ofσ(Z) in Section 2, we have (21)σZ=rx1-ηx-pxyrx1-ηx-pxy-ε1x.Let σ(Z)<0, and one yields(22)rp1-ηx-ε1p<y<rp1-ηx.It is reasonable to assume that the therapy should be done before virus growing and approaching its carrying capacity; therefore, in this study, we always assume that VT<1/η holds true. Then, if (1-ε1/r)/η<VT<1/η  i.e.,x2η<VT<x1η, the sliding segment is Σs1=(x,y)∣x=VT,0≤y≤r/p(1-ηVT); if 0<VT<x2η, then the sliding segment is Σs2=x,y∣x=VT,r/p(1-ηVT)-ε1/p≤y≤r/p(1-ηVT).Next, we employ Utkin’s equivalent control method introduced in [39] to examine the sliding dynamics in the region Σs. It follows from H(x)=0 and the first equation of system (3) that(23)dHdt=dxdt=rx1-ηx-pxy-Φε1x=0.Solving (23) with respect to Φ yields (24)Φ=rε11-ηx-pε1y.Then substituting Φ into the second equation of system (3) gives(25)dydt=y-ε2ε1py+ε2ε1r1-ηVT+cVT1+ωVT-qVT-δ≜QsVT,y.Then, the sliding dynamics of system (3) can be described as (25). Let Qs(VT,y)=0, and we obtain two roots(26)y0=0,yp=rp1-ηVT+ε1ε2pcVT1+ωVT-qVT-δ.This indicates that there are two possible pseudoequilibria E0=(VT,0) and Ep=(VT,yp). In order to verify whether they are located at the sliding segment, we consider two situations 0<VT<x2η and x2η<VT<x1η.It is easy to see thatE0 can not be a pseudoequilibrium when 0<VT<x2η, while it is always a pseudoequilibrium for x2η<VT<x1η. Thus, we only need to analyze the existence of Ep. First, we prove the following lemmas.Lemma 5.
IfΔ1>0 and Δ2>0 hold true, then we have x21<x11<x12<x22.Proof.
Let(27)g1z=c-q-ωz-c-q-ωz2-4qωz,withz>0.To ensure the function to be well defined, we need c-q-ωz>2qωz (i.e., Δ1>0). Further, we know(28)c-q-ωz>2qωz⟺c-q-ωz>0,c-q-ωz2-4qωz>0⟺z<c-qω,0<z<c+q-2qcωorz>c+q+2qcω. Definitely, we have c-q/ω<c+q+2qc/ω, and c-q/ω>c+q-2qc/ω hold true. Therefore, the domain for g1(z) is 0,c+q-2qc/ω, which is a continuous interval. Through simple calculation, we obtain(29)g1′z=-ω+ωc+q-ωzc-q-ωz2-4qωz.Then, we have(30)g1′z>0⟺c+q-ωz>c-q-ωz2-4qωz⟺4qc>0.
Because all the parameters are positive,g1′(z)>0 always holds true in its domain. Then, according to Δ1>0, Δ2>0, and δ>δ-ε2>0, we have x21<x11. Similarly, let(31)g2z=c-q-ωz+c-q-ωz2-4qωz.which is also defined in the interval 0<z<c-q-2qc/ω. Note that g2′(z)<0 is equivalent to c+q-ωz>0. Further, it is obvious that c+q-ωz>c-q-ωz>0. This means that g2′(z)<0 always holds true in its domain. Thus, we obtain that x12<x22 because Δ1>0, Δ2>0 and δ>δ-ε1>0. In conclusion, we have x21<x11<x12<x22 if Δ1>0 and Δ2>0. The proof is completed.Lemma 6.
IfΔ1<0 (or Δ2<0) holds true, then we have yp<r/p(1-ηVT)  oryp<r/p(1-ηVT)-ε1/p. If Δ1>0 and Δ2>0, then we have that r/p(1-ηVT)-ε1/p<yp<r/p(1-ηVT) for x21<VT<x11 or x12<VT<x22; yp<r/p(1-ηVT)-ε1/p for 0<VT<x21 or VT>x22; and yp>r/p(1-ηVT) for x11<VT<x12.Proof.
Letf1(x)=cx/1+ωx-qx-δ, and we know that f1(x)<0 if Δ1<0. Therefore, it follows from the formula of yp that when Δ1<0, then yp<r/p(1-ηVT). If Δ1>0, we get f1(x)>0 for x11<x<x12, and f1(x)<0 for 0<x<x11 or x>x12. This means that, when Δ1>0, we get yp>r/p(1-ηVT) for x11<VT<x12, and yp<r/p(1-ηVT) for 0<VT<x11 or VT>x12. Rearranging the formula of yp yields(32)yp=rp1-ηVT-ε1p+ε1ε2pcVT1+ωVT-qVT-δ+ε2.Therefore, applying similar analysis we have that when Δ2<0, then yp<r/p(1-ηVT)-ε1/p; when Δ2>0, then yp>r/p(1-ηVT)-ε1/p for x21<VT<x22, and yp<r/p(1-ηVT)-ε1/p for 0<VT<x21 or VT>x22. This completes the proof.Lemma 7.
IfΔ1>0 and Δ2>0, then we have yp>min{y21,y11} for x21<VT<x11.Proof.
Taking the derivative ofyp with respect to VT gives(33)dypdVT=-rpη+ε1ε2pc1+ωVT2-q.Let dyp/dVT>0, we have -1/ωcε1/(rηε2+qε1)-1/ω<VT<1/ωcε1/(rηε2+qε1)-1/ω.
Thus, if1/ωcε1/(rηε2+qε1)-1/ω>0, yp(VT) is increasing for 0<VT<1/ωcε1/(rηε2+qε1)-1/ω and decreasing for VT>1/ωcε1/(rηε2+qε1)-1/ω; if 1/ωcε1/(rηε2+qε1)-1/ω<0, then yp(VT) is decreasing for VT>0. Therefore, when x11<VT<x21, there may be three cases: (a) yp(VT) is always decreasing as VT increases, which means yp>yp(x11)=y11; (b) yp(VT) is always increasing as VT increases, then yp>yp(x21)=y21; (c) yp(VT) is initially increasing and then decreasing as VT increases, and hence yp>min{yp(x21),yp(x11)}=min{y21,y11}. This completes the proof.It is easy to verify that ifΔ2<0, then we have Δ1<0, which means that there can not be Δ1>0 when Δ2<0. According to the existence of possible equilibria of subsystems S1 or S2, we consider the following cases to discuss the existence of the pseudoequilibrium Ep.Case  A1.There is no positive equilibrium for subsystem  S2. Then, we have three subcases:(i)There is no positive equilibrium for subsystem  S1(denoted by Case  A11). Hence, there could be(34)Δ1<0,Δ2<0;orΔ1<0,Δ2>0,x21>r-ε1rη;orΔ1>0,x11>1η,Δ2>0,x21>r-ε1rη.It follows from Lemma 6 that for Case A11 there is always yp<r/p(1-ηVT)-ε1/p for 0<VT<x2η. This means that there is no pseudoequilibrium for 0<VT<x2η. Then, we mainly consider x2η<VT<x1η. Now, there is r/p(1-ηVT)-ε1/p<0. Therefore, we have yp<r/p(1-ηVT)-ε1/p<0 for the first subcase of CaseA11 (i.e., Δ1<0, Δ2<0), which means Ep is not feasible. For the other two subcases, if x21>x1η, Ep is not feasible as well. However, if x2η<x21<x1η, then the sign of yp is not determined for the whole interval x2η<VT<x1η. We denote(35)H1:rp1-ηVT+ε1ε2pcVT1+ωVT-qVT-δ>0i.e.,yp>0;H2:rp1-ηVT+ε1ε2pcVT1+ωVT-qVT-δ<0i.e.,yp<0.Therefore, Ep is feasible for H1, while it is not feasible for H2.(ii)There is a unique equilibrium for subsystem  S1(denoted by CaseA12); that is, (36)Δ1>0,x11<1η<x12,Δ2>0,x21>r-ε1rη.For such case, there is no pseudoequilibrium for 0<VT<x2η. There are two pseudoequilibria E0 and Ep for H1 and x21<VT<x11, and there is only one pseudoequilibrium E0 for x2η<VT<x21 or x11<VT<x1η or H2 and x21<VT<x11.(iii)There are two equilibria for subsystem  S1(denoted by Case  A13); that is, (37)Δ1>0,x12<1η,Δ2>0,x21>r-ε1rη.Similar to Case A12, there is no pseudoequilibrium for 0<VT<x2η. Ep is feasible for H1, x21<VT<x11, or H1, x12<VT<min{x1η,x22}.Case  A2.There is one positive equilibrium for subsystem  S2. Then, consider the following:(i)There is no positive equilibrium for subsystem  S1(denoted by Case  A21). Then, there should be(38)Δ1<0,Δ2>0,x21<r-ε1rη<x22;orΔ1>0,x11>1η,Δ2>0,x21<r-ε1rη<x22.It follows from Lemma 6 that if Δ1<0  orΔ1>0,y11<0, then we have yp<r/p(1-ηVT) for 0<VT<x1η. Further, because Δ2>0,y22<0<y21, according to Lemma 6, there is no pseudoequilibrium for 0<VT<x21 while Ep is feasible for x21<VT<x2η. Similarly, for x2η<VT<x1η, the sign of yp can not be determined. Thus, when x2η<VT<x1η, there might be either two pseudoequilibria E0 and Ep for H1, or one pseudoequilibrium E0 for H2. Especially, if x22<x1η, then Ep is not feasible for x22<VT<x1η.(ii)There is a unique positive equilibrium for subsystem  S1(denoted by Case  A22); that is,(39)Δ1>0,x11<1η<x12,Δ2>0,x21<r-ε1rη<x22.For such case, it follows from Lemmas 6 and 7 that Ep is feasible for x21<VT<x11. Especially, if x11>x2η, then there are two pseudoequilibria E0 and Ep for x2η<VT<x11, and there is no pseudoequilibrium for 0<VT<x21.(iii)There are two positive equilibria for subsystem  S1(denoted by Case  A23); that is, (40)Δ1>0,x12<1η,Δ2>0,x21<r-ε1rη<x22.Here, we consider three subcases: x12<x2η, x11<x2η<x12, or x11>x2η in terms of the relationship of x2η to x11 and x12. If x12<x2η, then Ep is feasible for x21<VT<x11 or x12<VT<x2η or H1 and x2η<VT<min{x22,x1η}. If x11<x2η<x12, then Ep is feasible for x21<VT<x11 or H1 and x12<VT<min{x22,x1η}. If x11>x2η, then it follows from Lemma 7 that Ep is feasible for x21<VT<x11 or H1 and x12<VT<min{x22,x1η}.Case A3.There are two positive equilibria for subsystemS2. Under this case, there can not be a unique positive equilibrium for subsystem S1 according to Lemma 5. Then, consider the following:(i)There is no positive equilibrium for subsystem  S1(denoted by Case  A31); that is, (41)Δ1<0,Δ2>0,x22<r-ε1rη.It follows from Lemma 6 that Ep is feasible for x21<VT<x22.(ii)There are two positive equilibria for subsystem  S1(denoted by Case  A32); that is,(42)Δ1>0,x12<1η,Δ2>0,x22<r-ε1rη.Similarly, according to Lemma 6, Ep is a pseudoequilibrium for x21<VT<x11 or x12<VT<x22.Having discussed the existence of the pseudoequilibria, we take the derivative ofQs with respect to y to analyze the stability of the pseudoequilibria. Note that (43)dQsdy=-ε2ε1py+ε2ε1r1-ηVT+cVT1+ωVT-qVT-δ-ε2ε1py.Then, we have dQs/dyE0=ε2/ε1r(1-ηVT)+cVT/1+ωVT-qVT-δ=ε2/ε1pyp and dQs/dyEp=-ε2/ε1pyp. Therefore, if H1 holds true, then the pseudoequilibrium Ep is locally stable on the sliding segment when it is feasible, and E0 is unstable; if H2 holds true, then Ep becomes infeasible and the pseudoequilibrium E0 is locally stable on the sliding segment when it exists.


## 4. Global Dynamics of System (3)

In this section, we will investigate the global dynamics of Filippov system (3). Initially, we shall rule out the existence of limit cycle for system (3). Here, we take Case A11 as an example to illustrate that there is not any limit cycles. The detailed proof is given in the appendix. Note that, for the other cases, we can prove the nonexistence of limit cycle for system (3) using the similar methods, and here we omit the detailed proof.Lemma 8.
For CaseA11, there is no limit cycle for Filippov system (3).Corresponding to the discussion of the sliding dynamics, we then consider the following cases to examine the global dynamics of system (3).Case  A1.There is no positive equilibrium for subsystem  S2.(i)There is no positive equilibrium for subsystemS1(i.e., CaseA11). Then, we discuss the global dynamics of the Filippov system by considering the following subsituations:(a) If0<VT<x2η, then we know that E2η is a real equilibrium which is stable while E1η is a virtual equilibrium which is unstable. As we mentioned above, there is no limit cycle for system (3), and moreover there is not any other stable equilibrium, thus the local stability of E2η indicates that E2η is globally asymptotically stable (shown in Figure 2(a)).(b) Ifx2η<VT<x1η and H2 hold true, then we have that E2η becomes a virtual equilibrium. It follows from the discussion in Section 3 that the pseudoequilibrium E0 is stable on the sliding segment Σs1. Similarly, system (3) does not have any limit cycle. Therefore, the local stability of E0 indicates that E0 is globally asymptotically stable (shown in Figure 2(b)).(c) Ifx2η<VT<x1η and H1 hold true, then there are two pseudoequilibria E0 and Ep on the sliding segment Σs1. As we proved in Section 3, Ep is locally stable while E0 is unstable. Similar to subcase (b), we have that Ep is globally asymptotically stable (shown in Figure 2(c)).Based on the above discussions, the global dynamics of system (3) can be concluded as follows.Theorem 9.
For CaseA11, if  0<VT<x2η, then the boundary equilibrium E2η is globally asymptotically stable; if x2η<VT<x1η, then the pseudoequilibrium E0 is globally asymptotically stable for (H2) being satisfied, while the pseudoequilibrium Ep is globally asymptotically stable for (H1) being satisfied.(ii)There is a unique positive equilibrium for subsystemS1(i.e., CaseA12). If 0<VT<x2η, E1η and E11 are virtual equilibria, and E2η is a real equilibrium which is stable. As mentioned above, under this scenario, there is no limit cycle or no pseudoequilibrium for Filippov system (3), and thus E2η is globally stable. If x2η<VT<x21, E2η becomes a virtual equilibrium, and there emerges the stable pseudoequilibrium E0 on the sliding segment Σs1. If x21<VT<x11, the same results hold compared with the case x2η<VT<x21 if H2 holds true, while there are two pseudoequilibria E0 (unstable) and Ep (stable) on Σs1 if H1 holds true. Furthermore, If x11<VT<x1η, E11 becomes a real equilibrium which is locally stable. There is a pseudoequilibrium E0 on Σs1 which is unstable. Then, the global dynamics of Filippov system (3) are concluded as follows.Theorem 10.
For CaseA12, if 0<VT<x2η, then E2η is globally asymptotically stable (shown in Figure 3(a)); if x2η<VT<x21, then E0 is globally asymptotically stable (shown in Figure 3(b)); if x21<VT<x11, then E0 is globally asymptotically stable for (H2) being satisfied and Ep is globally asymptotically stable for (H1) being satisfied (shown in Figure 3(c)); if x11<VT<x1η, then E11 is globally asymptotically stable (shown in Figure 3(d)).(iii)There are two positive equilibria for subsystemS1(i.e., CaseA13). If there is only one steady state (boundary equilibrium or positive equilibrium or pseudoequilibrium), then it is globally stable. We can analyze the dynamics of system (3) under the same scenario using the similar methods. However, there is a different situation for such case. That is, if x12<VT<min{x1η,x22} holds true, then E11 is a real equilibrium; hence it is locally stable. Furthermore, E0 is a locally stable pseudoequilibrium if H2 holds true, but if H1 holds true, E0 lose its stability and Ep becomes a locally stable pseudoequilibrium. In conclusion, if x12<VT<min{x1η,x22}, E11 is bistable with E0(orEp) for H2(orH1) being satisfied. Based on the discussion above, the global dynamics of system (3) for Case A13 can be concluded as follows.Theorem 11.
For CaseA13, if  0<VT<x2η, then E2η is globally asymptotically stable (shown in Figure 4(a)); if x2η<VT<x21, then E0 is globally asymptotically stable (shown in Figure 4(b)); if x21<VT<x11, then E0 (or Ep) is globally asymptotically stable for (H2)  (or (H1)) being satisfied (shown in Figure 4(c) (or Figure 4(d))); if x11<VT<x12, then E11 is globally asymptotically stable (shown in Figure 4(e)); if x12<VT<min{x1η,x22}, then E0 and E11 are bistable for (H2) being satisfied and Ep and E11 are bistable for (H1) being satisfied (shown in Figure 4(f)). In particular, E0 and E11 are always bistable for x22<VT<x1η.Cases A2  (A21,A22,andA23)andA31. All the dynamics of system (3) for four cases are similar to the former three cases, and thus we just omit the proof here. The global dynamics of these cases are summarized in Table 1 with the dynamical behaviors shown in Figures 5–8, respectively.Table 1
Dynamics of system (3) for Cases A11–A32.
Cases Conditions 0 < V T < x 2 η x 2 η < V T < x 1 η A11 Δ 1 < 0 , Δ 2 < 0 E 2 η GAS E 0 GAS Δ 1 < 0 , Δ 2 > 0 , x 21 > x 2 η E p GAS for H1 E 0 GAS for H2 Δ 1 > 0 , x 11 > x 1 η , Δ 2 > 0 , x 21 > x 2 η E p GAS for H1 E 0 GAS for H2  A12 Δ 1 > 0 , x 11 < x 1 η < x 12, E 2 η GAS x 2 η < V T < x 21 x 21 < V T < x 11 x 11 < V T < x 1 η Δ 2 > 0 , x 21 > x 2 η E 0 GAS E p GAS for H1 E 0 GAS for H2 E 11 GAS  A13 Δ 1 > 0 , x 12 < x 1 η, E 2 η GAS x 2 η < V T < x 21 x 21 < V T < x 11 x 11 < V T < x 12 x 12 < V T < x 1 η Δ 2 > 0 , x 21 > x 2 η E 0 GAS E p GAS for H1 E 0 GAS for H2 E 11 GAS E p , E 11 BS for H1 E 0 , E 11 BS for H2  Cases Conditions 0 < V T < x 21 x 21 < V T < x 1 η  A21 Δ 1 < 0 , Δ 2 > 0, E 21 GAS(Figure 5(a)) x 21 < V T < x 2 η x 2 η < V T < x 1 η x 21 < x 2 η < x 22 or E p GAS (Figure 5(b)) E p GAS (Figure 5(c)) for H1 E 0 GAS (Figure 5(d)) for H2 Δ 1 > 0 , x 11 > x 1 η, Δ 2 > 0 , x 21 < x 2 η < x 22  A22 Δ 1 > 0 , x 11 < x 1 η < x 12, E 21 GAS x 21 < V T < x 11 x 11 < V T < x 1 η Δ 2 > 0 , x 21 < x 2 η < x 22 (Figure6(a)) E p GAS (Figure 6(b)) E 11 GAS (Figure 6(c))  A23 Δ 1 > 0 , x 12 < x 1 η, E 21 GAS x 21 < V T < x 11 x 11 < V T < x 12 x 12 < V T < x 1 η Δ 2 > 0 , x 21 < x 2 η < x 22 (Figure7(a)) E p GAS (Figure 7(b)) E 11 GAS (Figure 7(c)) E p , E 11 BS (Figure 7(d)) for H1 E 0 , E 11 BS for H2  A31 Δ 1 < 0 , Δ 2 > 0, E 21 , E 2 η BS x 21 < V T < x 22 x 22 < V T < x 2 η x 2 η < V T < x 1 η x 22 < x 2 η (Figure8(a)) E p , E 2 η BS (Figure 8(b)) E 2 η GAS (Figure 8(c)) E 0 GAS (Figure 8(d))  A32 Δ 1 > 0 , x 12 < x 1 η, E 21 , E 2 η BS x 21 < V T < x 11 x 11 < V T < x 12 x 12 < V T < x 22 x 22 < V T < x 2 η x 2 η < V T < x 1 η Δ 2 > 0 , x 22 < x 2 η E p , E 2 η BS E 11 , E 2 η BS E p , E 11 , E 2 η TS E 11 , E 2 η BS E 0 , E 11 BS Note: “GAS” represents “globally asymptotically stable”; “BS” means bistable; and “TS” is tristable.Case  A32.There are two positive equilibria for both subsystemsS1andS2, respectively. And the horizontal compartments of these equilibria satisfy x21<x11<x12<x22. In terms of the existence and stability of the positive equilibria, the dynamics of the subsystems are similar to the focused situation in study [30]. Thus, based on the relationships of VT with respect to x21, x11, x12, and x22, we can analyze the dynamics of system (3) through the following subcases.If the threshold condition is relatively small, i.e.,0<VT<x21, only E2η, E21, and E22 are real equilibria. It follows from the dynamics of system S2 that both E2η and E21 are stable, and thus they are bistable for the Filippov system (3). If VT increases and satisfies x21<VT<x11, E21 becomes a virtual equilibrium while the locally stable pseudoequilibrium Ep emerges, and thus E2η and Ep are bistable. Further, we can conclude that E11 is bistable with E2η when x11<VT<x12 or x22<VT<x2η, while it is bistable with E0 when x2η<VT<x1η. If x12<VT<x22, E11 and E2η are stable while there also exists pseudoequilibrium Ep which is locally stable as well. Thus, based on the above discussions, we have the following conclusions.Theorem 12.
For CaseA32, if 0<VT<x21, then equilibrium E21 and the boundary equilibrium E2η are bistable (shown in Figure 9(a)); if x21<VT<x11, then Ep and E2η are bistable (shown in Figure 9(b)); if x11<VT<x12, then E11 and E2η are bistable (shown in Figure 9(c)); if x12<VT<x22, then Ep, E11, and E2η are locally tristable (shown in Figure 9(d)); if x22<VT<x2η, E2η and E11 are bistable (shown in Figure 9(e)); if x2η<VT<x1η, then E0 and E11 are bistable (shown in Figure 9(f)).


## 5. Conclusion and Discussion

Based on the basic immune-virus model, we considered the threshold policy depending on virus loads for the HIV-infected patients and obtained a nonsmooth system. Compared with the model in [30], we extend it by taking the saturated growth of virus into consideration, and hence our model is better in line with the biological reality. Meanwhile, our model also includes a piecewise elimination rate of HIV virus and growth rate of the effector cells to indicate that the antiretroviral therapy and interleukin-2 (IL-2) treatment are carried out simultaneously once the HIV virus loads increase above a threshold level.We first briefly concluded the dynamics of two subsystems. Based on the properties of the subsystems, we discussed the sliding mode and sliding dynamics with different threshold conditions. In more detail, the sliding segments areΣs1=(x,y)∣x=VT,0≤y≤r/p(1-ηVT), for x2η<VT<x1η, and Σs2=x,y∣x=VT,r/p(1-ηVT)-ε1/p≤y≤r/p(1-ηVT) for 0<VT<x2η. We note that there is always a pseudoequilibrium E0 on the sliding segment Σs1. And we further discussed the existence of the other pseudoequilibrium Ep by considering eight different cases. Correspondingly, we then investigated the global dynamics of system (3) for various cases. Our main results show that the system can exhibit very rich dynamic behaviors: (a) one of the equilibria is globally asymptotically stable, which can be the pseudoequilibrium E0 or the pseudoequilibrium Ep or the boundary equilibrium E2η of subsystem S2 or the positive equilibrium E11 or E21 of subsystem S1 or S2; (b) two equilibria are bistable, particularly, which can be that the positive equilibrium E11 is bistable with E0 (or Ep or E2η), or the boundary equilibrium E2η is bistable with the pseudoequilibrium Ep (or the positive equilibrium E21); (c) three equilibria are locally tristable; that is, the pseudoequilibrium Ep, the boundary equilibrium E2η, and the positive equilibrium E11 are stable for x12<VT<x22.Global dynamics of our proposed model means that there is the case that the boundary equilibrium is globally stable for both the subsystems. Note that the immune-virus system without any therapy is actually systemS1 (the free system) while the immune-virus system with continuous therapy can be described by system S2 (the control system). Comparing system S2 with system S1, we can find that the continuous therapy can shrink the final density of the virus (i.e., x2η<x1η). It follows that the continuous therapy failed to maintain the activity of the immune system since the effector cell goes to zero finally. However, threshold policy indicated that, with a proper threshold condition, the STIs can not only successfully control virus loads to a previously given level which can be far below its natural carrying capacity (i.e., x1η) but also maintain the activity of the immune system (corresponding to the case that the pseudoequilibrium Ep is globally stable), as Theorem 9 suggested.Moreover, according to the dynamics of subsystemS1 for Case A21 we can see that virus loads will finally approach its carrying capacity if there is no therapy for the patients. With the continuous therapy, the virus loads of the patients can be successfully controlled far below its carrying capacity and their immune system can always maintain the activity; that is, equilibrium E21 is globally asymptotically stable for system S2. For the virus-guided threshold policy, if the threshold level is set to satisfy x21<VT<x2η, then the pseudoequilibrium Ep is globally stable, and Ep will tend to E21 as VT→x21. This means that the threshold policy can obtain the similar effects for the patients compared with the continuous therapy in terms of virus loads being controlled, which is in agreement with the clinical study [7].It is worth mentioning that if we also assume that the two subsystems have two positive equilibria (i.e., CaseA32), our model results, due to the saturated growth in model equations, show that the virus can not grow to infinity. Under this situation, we also found that the dynamics of system (3) are similar to those of the control system (system S2) when the threshold condition is relatively small (i.e., 0<VT<x21); that is, the positive equilibrium E21 and the boundary equilibrium E2η are bistable. And if the threshold is relatively high (i.e., VT>x22), the dynamics of system (3) are similar to the free system S1. Moreover, for middle values of the threshold (i.e., x12<VT<x22), then equilibria Ep, E11, and E2η are locally tristable. Therefore, the bistability and the tristability suggest that the optimal strategy should also be individual-based by taking the initial conditions of the patients and the threshold conditions into consideration.



---

*Source: 4710586-2018-07-09.xml*