Discrete Dynamics in Nature and Society

Volume 2013 (2013), Article ID 712829, 9 pages

http://dx.doi.org/10.1155/2013/712829

Research Article

Hepatitis B Virus Dynamics: Modeling, Analysis, and Optimal Treatment Scheduling

A. M. Elaiw,1,2 M. A. Alghamdi,¹ and Shaban Aly2,3

¹Department of Mathematics, Faculty of Science, King Abdulaziz University, P.O. Box 80203, Jeddah 21589, Saudi Arabia

²Department of Mathematics, Faculty of Science, Al-Azhar University (Assiut Branch), Assiut 71511, Egypt

³Department of Mathematics, Faculty of Science, King Khalid University, Abha 9004, Saudi Arabia

Received 8 December 2012; Revised 11 March 2013; Accepted 13 March 2013

Academic Editor: Manuel De la Sen

Copyright © 2013 A. M. Elaiw et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Modeling, analysis, and control of hepatitis B virus (HBV) infection have attracted the interests of mathematicians during the recent years. Several mathematical models exist and adequately explain the HBV dynamics as well as the effect of antiviral drug therapies. However, none of these models can completely exhibit all that is observed clinically and account the full course of infection. Besides model inaccuracies that HBV dynamics models suffer from, some disturbances/uncertainties from different sources may arise in the modeling. In this paper, the HBV dynamics is described by a system of nonlinear ordinary differential equations. The disturbances or uncertainties are modeled in the HBV dynamics model as additive bounded disturbances. The model is incorporated with two types of drug therapies which are used to inhibit viral production and prevent new infections. The model can be considered as nonlinear control system with control input is defined to be dependent on the drug dose and drug efficiency. We developed treatment schedules for HBV infected patients by using multirate model predictive control (MPC). The MPC is applied to the stabilization of the uninfected steady state of the HBV dynamics model. The inherent robustness properties of the MPC against additive disturbances are also shown.

1. Introduction

Hepatitis B virus (HBV) infection is among the most common causes of hepatitis and can result in serious liver diseases such as chronic hepatic insufficiency, hepatocellular carcinoma, and cirrhosis. Over the last decade, much collaborated effort involving biologists and mathematicians has been devoted towards designing mathematical models of HBV dynamics [1–11], human immunodeficiency virus (HIV) [12–20], and hepatitis C virus (HCV) [21]. Mathematical modeling and model analysis of the HBV dynamics are important for exploring possible mechanisms and dynamical behaviors of the viral infection process, estimating key parameter values, and guiding development of efficient antiviral drug therapies. Stability analysis of HBV dynamics models has been studied by many authors (see e.g., [22–32]).

Several drug therapies have been proposed for treating persons with chronic HBV including adefovir dipivoxil, alpha-interferon, lamivudine, pegylated interferon, entecavir, telbivudine, and tenofovir [33]. Hepatitis antiviral drugs prevent replication of HBVs and save the liver from cirrhosis and cancer. During the treatment, the viral load is reduced and consequently the viral replication in liver is decreased [34].

Optimal control theory can be used to optimize the drug doses required in the treatment of HBV infected patients [33, 35, 36]. In this paper, the optimal treatment schedules will be designed using model predictive control (MPC) method. MPC obtains a feedback control by solving a finite horizon optimal control problem at each time instant using the current state of the system as the initial state for the optimization and applying “the first part” of the optimal control [37]. MPC has many advantages which include its facility of handling constraints, its closed-loop stability, and inherent robustness. Moreover, MPC solves optimal control problem on-line for the current state of the plant which is a mathematical programming problem and is much simpler than determining the feedback solution by dynamic programming [37].

Recently, MPC has been applied for determining the optimal treatment schedules for HIV infected patients in [38–41]. As far as the authors know, the application of MPC to the HBV dynamics model has never been studied before. The aim of the present paper is to apply a multirate MPC approach to an HBV dynamics model. The MPC is constructed based on the approximate discrete-time model. The closed-loop stability of the multirate MPC will be shown. The inherent robustness properties of the MPC against disturbances will also be shown. The disturbances may arise from modeling errors, or disturbances arise from immune system fluctuating or immune effect of a coinfection.

2. HBV Infection Model

We will use the mathematical model proposed in [22], incorporating to allow some additive disturbances and taking into account the effect of two types of antiviral drugs

Here , , and denote the concentration of uninfected hepatocytes, infected hepatocytes, and free virions, respectively. Uninfected hepatocytes are assumed to be produced at the constant rate and die at the rate of . Uninfected hepatocytes can also be created by proliferation of existing cells. Here we represent the proliferation by a logistic function in which is the maximum proliferation rate of target cells, and is the concentration at which proliferation shuts off. The rate of infection is given by saturation functional response , where is the infection rate constant characteristic of the infection efficiency and is constant. The death rate of infected hepatocytes is . Free virions are assumed to be produced from infected hepatocytes at the rate of and is the clearance rate of viral particles. Functions , , describe the model uncertainties/disturbances. Two types of drugs are considered; the first is used to prevent the virus from infecting the cell and is represented by the chemotherapy function , and the second is used to prevent the infected cells from producing new viruses and is represented by the function . Here, , are the control inputs, where the parameters and are the efficiencies of the drugs and and are the drug doses (see [42]). All the parameters of the model are supposed to be positive. System (1)–(3) can be considered as a nonlinear control system with as the state vector and is the control input vector. We are now ready to present a study on the basic mathematical properties of system (1)–(3).

We assume that the model uncertainty/disturbances are given by

where is the external disturbances and satisfies the following bound

Now we show that under which conditions the nonnegative orthant is positively invariant for (1)–(3):

This means that the nonnegative orthant is positively invariant, namely, if a trajectory starts in the nonnegative orthant, it remains there.

Proposition 1. Assume that the disturbances satisfy the boundedness condition (5) and . Then there exist positive numbers and for which the compact set

is positively invariant.

Proof. From (1) we have

This means that, even if reaches , it should decrease.

Let , then

where . Hence for all if . It follows that , for all if , where . On the other hand, . Let , then we get that for all , whenever this inequality holds for and .

Note that contains all the biologically relevant states, thus we can restrict the state space of the system to the compact set . Since the drug doses cannot be arbitrarily increased, we then may consider only a compact control constraint set.

Let us compute the steady states of system (1)–(3) under constant controller in the absence of the disturbances that is for ,, , , . System (1)–(3) has two steady states, uninfected steady state and infected steady state , where

It has been shown in [22] that the infected steady state exists if .

The basic reproductive ratio for system (1)–(3) is given by

Theorem 2 (see [22]). For the nominal system (1)–(3) (i.e., , ), if , then is locally asymptotically stable.

Proposition 3. The nominal system (1)–(3) is locally asymptotically controllable to with piecewise constant controllers.

Proof. Let and with , where

Then ; therefore, the corresponding trajectory will tend to as .

Theorem 4 (see [22]). For the nominal system (1)–(3), if , then is globally asymptotically stable.

Assumption 5. There exist and such that

Proposition 6. Suppose that Assumption 5 is valid. Then the nominal system (1)–(3) is globally asymptotically controllable to with piecewise constant controllers.

Proof. Let and such that and satisfy inequality (13), then , and the corresponding trajectory will tend to as .

Theorem 7 (see [22]). Consider the nominal system (1)–(3). Suppose that

(i), (ii) .

Then the positive steady state is locally asymptotically stable.

Define

Theorem 8 (see [22]). Consider the nominal system (1)–(3). Suppose that

(i), (ii) > .

Then the positive steady state is globally asymptotically stable provided that one of the following two assumptions holds:

(i), (ii) .

3. Multirate MPC for Nonlinear Systems

In this section, we outline the multirate MPC design for nonlinear systems in the presence of bounded disturbances and give a review on the results obtained in [43, 44].

The set of real and natural numbers (including zero) are denoted, respectively, by and . The symbol denotes the set of nonnegative real numbers. A continuous function is of class- if , for all and it is strictly increasing. It is of class- if it is of class- and when . A continuous function is of class- if is of class- in for every , it is strictly decreasing in for every and when . Let will be used in .

Consider a continuous-time nonlinear control system with additive disturbances given by

where , , are the state, control input, and disturbances, respectively, is continuous and Lipschitz continuous with respect to in any compact set and , is compact and , is compact and .

The control is taken to be a piecewise constant signal

where is the control sampling period which is fixed.

In this paper we address the problem of state feedback stabilization of (15) under “low measurement rate.” More precisely, we will assume that state measurements can be performed at the time instants , :

where is the measurement sampling period. We assume that for a fixed integer .

For a given function , we use the following notation: , where . We denote the norm . We assume that there exists such that . Let us define

We will assume that there is a compact set containing the origin, that is positively invariant with respect to system (15) for any and any piecewise constant controller . Let denote the solution of (15) with given , , and . Then the exact discrete-time model can be defined as

where .

Let , , and , then the exact -step discrete-time model is given by

We note that the exact discrete-time models (19) and (20) describe, respectively, the behavior of the system at the time instants and ,.

In this work, the construction of multirate MPC is based on the nominal prediction and only small disturbances are allowed. The nominal system of (15) is given by

and its exact discrete-time model is given by

We note that, since is typically nonlinear, in (22) is not known in most cases, therefore the controller design can be carried out by means of the nominal approximate discrete-time model

where is a modeling parameter, which is typically the step size of the underlying numerical method. The applied numerical scheme approximation has to ensure the closeness of the exact models in the following sense.

Assumption 9. There exists an such that

(i), is continuous in both variables, uniformly continuous in , and Lipschitz continuous with respect to in any compact set, uniformly in small ,(ii)there exists a such that

for all , all , and .

Assumption 10. There exists an such that the nominal exact discrete-time model (22) is practically asymptotically controllable from to the origin with piecewise constant controllers for all (see e.g., [43] for the definition).

For the solutions of (19), (22), and (23) with , and we will use the notations ,, and , respectively.

The following problem is to be solved: for given and find a control strategy

using the nominal approximate discrete-time model (23), to practically stabilize the exact discrete-time system (19).

Let with be given. Let (23) be subject to the cost function

where ,,, and are given functions satisfying the following assumptions.

Assumption 11. Let , (i) is continuous, positive definite radially unbounded, and Lipschitz continuous in any compact set,

(ii) is continuous with respect to and , uniformly continuous in small , and Lipschitz continuous in any compact set,

(iii) there exist an and two class- functions and such that the inequality

holds for all and .

Assumption 12. There exist and such that for all there exists a such that inequality

holds for all .

We define the value function, which represents the optimal value of (26) for a given initial condition as follows

Let be the solution of this optimization problem, then the first elements of are applied at the state , that is,

Let denote the minimum of the values generated by Assumptions 9–12. Let and be such numbers that, for each , , .

Theorem 13 (see [43]). If Assumptions 9–12 hold, then

(i)there exist an with , and a constant independent of , such that for all , and ,(ii)there exist constants ,and and functions such that for all , , and ,

Moreover, for all with ,

for all .

Clearly .

Theorem 14 (see [38, 44]). Suppose that Assumptions 9–12 are valid and is chosen such that . Then, there exist ,, and for any there exists an such that for any , and the trajectory of the -step exact discrete-time system

with the -step MPC and satisfies

4. MPC for the HBV Model

In this section we apply the MPC method proposed in Section 3. We will show that, with a suitable choice of and functions and , the assumptions of the previous section can be satisfied. Introduce new variables as ,,. In these new variables the model (1)–(3) takes the form of (15) with

and .

Let the compact set be defined as

where and are given in Proposition 1.

With this definition, satisfies all regularity assumptions, and according to Proposition 1, is positively invariant.

To verify Assumptions 11 and 12, we linearized the nominal system (35) around the origin in case of constant controllers, that is, , with , where is given in Proposition 3. Let be the coefficient matrix of the linearized system and . Then the discrete-time model for the linearized system is given by

Let the sampling period be chosen to be and . The running cost and the terminal cost can be chosen as

where is a positive definite diagonal matrix and is a positive definite symmetric matrix satisfying the Lyapunov equation for the discrete-time system (37)

Using the parameters given in Table 1, we have verified Assumption 5 numerically by fixing one controller, for example, and solve the inequality (13) with respect to . We have found that for , then inequality (13) is satisfied when .

Table 1: The values of the parameters in the HBV dynamics model.

Parameters

Values

From (38)-(39), Assumption 11 is satisfied. Assumption 10 follows from Proposition 6. Assumption 9 is fulfilled as well, if we choose a suitable numerical integration scheme (e.g., the Runge-Kutta formula). To verify Assumption 12, the matrix has been chosen through a series of numerical experiments . It has been verified numerically by solving a constrained minimization problem with several starting points, for which Assumption 12 is satisfied over the whole set . Thus all assumptions of the proposed method can be satisfied with suitable choice of the parameters of the MPC method.

5. Numerical Results

We performed simulation studies using the following parameter values which are listed in Table 1.

We assumed that the state measurements were performed at the instants , with and the horizon length was chosen to be . All computations were carried out by MATLAB. Especially, the optimal control sequence was computed by the fmincon code of the Optimization toolbox. Simulations for the continuous-time system were carried out using ODE45 program in MATLAB for two cases, untreated case and treated case.

(i) Untreated Case. In this case, we simulate the nominal system (1)–(3) (i.e., , ) without treatment (i.e., ). Using the parameters given in Table 1 and choose , we have obtained and , moreover, all the conditions of Theorem 8 are satisfied. Therefore the uninfected steady state is unstable and the infected steady state exists and is globally asymptotically stable. To show the simulation results for the untreated case, we assume that the infection occurs with a certain amount of virus particles . Thus the initial conditions for the untreated case are , , and . From Figures 1, 2, and 3, we can see that, the concentration of uninfected hepatocytes is decaying while the concentrations of the infected hepatocytes, and free viruses are increasing. Also we note that the trajectory tends to the stable infected steady state . Figures 1–3 show also the effect of the saturation parameter on the evolution of uninfected hepatocytes, infected hepatocytes and free viruses, respectively. It can be seen that, as is increased, the evolution of the disease is postponed.

Figure 1: The evolution of uninfected hepatocytes for untreated case.

Figure 2: The evolution of infected hepatocytes for untreated case.

Figure 3: The evolution of free virus for untreated case.

(ii) Treated Case. In this case, the treatment is designed via MPC strategy for the nominal and disturbed HBV infection models. In this case, we assume that the system is in the infected steady state before initiating the treatment, that is, . In this case, we choose the saturation parameter . Thus the initial conditions for the MPC are given by , , and .

The disturbances are simulated by

where the parameters are uniformly distributed random numbers on , and when the system states lie in the interior of the positive orthant . At the boundary of , the lower bound has to be chosen as the following:

to guarantee that the positive orthant is positively invariant.

We apply the MPC to the HBV model with and without disturbances. The following cases are shown in the simulation results:

(1)MPC-(I): , ,(2)MPC-(II) , with , , ,(3)MPC-(III) , with , , .

Figure 4 shows that, when the MPC is applied, the concentration of uninfected hepatocytes is increasing. From Figures 5 and 6, we can see that the concentrations of infected hepatocytes and free virions are decaying during the treatment. From Figures 4–6, we observe that when the MPC strategy is applied in the presence of bounded disturbances, the trajectory of the system tends to a ball around the uninfected steady state and remains there. The drug doses and as functions of the time for MPC-(I), MPC-(II), and MPC-(III) cases are shown in Figures 7 and 8, respectively.

Figure 4: The evolution of uninfected hepatocytes under MPC.

Figure 5: The evolution of infected hepatocytes under MPC.

Figure 6: The evolution of free virus under MPC.

Figure 7: The drug dose of the first type designed by MPC.

Figure 8: The drug dose of the second type designed by MPC.

6. Conclusion

In this paper we have studied a mathematical model describing the HBV dynamics. The model has been incorporated to allow some additive disturbances. Under the effect of two types of drug therapies the HBV dynamics model is considered as a nonlinear control system, where the control input is defined to be dependent on the drug dose and drug efficiency. We have applied MPC method for determining the optimal treatment schedules and for stabilizing the HBV infection system around the uninfected steady state. The inherent robustness properties of the MPC have been established.

Acknowledgments

This project was funded by the Deanship of Scientific Research (DSR), King Abdulaziz University, Jeddah, under Grant no. 408/130/1431. The authors, therefore, acknowledge with thanks DSR technical and financial support. The authors are grateful to the anonymous reviewers for constructive suggestions and valuable comments, which improved the quality of the paper.

References

M. A. Nowak, S. Bonhoeffer, A. M. Hill, R. Boehme, H. C. Thomas, and H. Mcdade, “Viral dynamics in hepatitis B virus infection,” Proceedings of the National Academy of Sciences of the United States of America, vol. 93, no. 9, pp. 4398–4402, 1996.

M. A. Nowak and R. M. May, Virus Dynamic: Mathematical Principles of Immunology and Virology, Oxford University Press, Oxford, UK, 2000.

S. M. Ciupe, R. M. Ribeiro, P. W. Nelson, and A. S. Perelson, “Modeling the mechanisms of acute hepatitis B virus infection,” Journal of Theoretical Biology, vol. 247, no. 1, pp. 23–35, 2007.

S. R. Lewin, R. M. Ribeiro, T. Walters et al., “Analysis of hepatitis B viral load decline under potent therapy: complex decay profiles observed,” Hepatology, vol. 34, no. 5, pp. 1012–1020, 2001.

M. Tsiang, J. F. Rooney, J. J. Toole, and C. S. Gibbs, “Biphasic clearance kinetics of hepatitis B virus from patients during adefovir dipivoxil therapy,” Hepatology, vol. 29, no. 6, pp. 1863–1869, 1999.

L. Min, Y. Su, and Y. Kuang, “Mathematical analysis of a basic virus infection model with application to HBV infection,” The Rocky Mountain Journal of Mathematics, vol. 38, no. 5, pp. 1573–1585, 2008.

C. Long, H. Qi, and S. H. Huang, “Mathematical modeling of cytotoxic lymphocyte-mediated immune response to hepatitis B virus infection,” Journal of Biomedicine and Biotechnology, vol. 2008, Article ID 743690, 9 pages, 2008.

R. Qesmi, J. Wu, J. Wu, and J. M. Heffernan, “Influence of backward bifurcation in a model of hepatitis B and C viruses,” Mathematical Biosciences, vol. 224, no. 2, pp. 118–125, 2010.

K. Wang and W. Wang, “Propagation of HBV with spatial dependence,” Mathematical Biosciences, vol. 210, no. 1, pp. 78–95, 2007.

K. Wang, W. Wang, and S. Song, “Dynamics of an HBV model with diffusion and delay,” Journal of Theoretical Biology, vol. 253, no. 1, pp. 36–44, 2008.

M. A. Nowak and C. R. M. Bangham, “Population dynamics of immune responses to persistent viruses,” Science, vol. 272, no. 5258, pp. 74–79, 1996.

D. Wodarz and M. A. Nowak, “Mathematical models of HIV pathogenesis and treatment,” BioEssays, vol. 24, no. 12, pp. 1178–1187, 2002.

M. A. Nowak, R. Anderson, M. Boerlijst, S. Bonhoeffer, R. May, and A. McMichael, “HIV-1 evolution and disease progression,” Science, vol. 274, pp. 1008–1010, 1996.

M. A. Obaid, “Global analysis of a virus infection model with multitarget cells and distributed intracellular delays,” Life Science Journal, vol. 9, no. 4, pp. 1500–1508, 2012.

A. M. Elaiw, “Global properties of a class of HIV models,” Nonlinear Analysis: Real World Applications, vol. 11, no. 4, pp. 2253–2263, 2010.

A. M. Elaiw and S. A. Azoz, “Global properties of a class of HIV infection models with Beddington-DeAngelis functional response,” Mathematical Methods in the Applied Sciences, vol. 36, no. 4, pp. 383–394, 2013.

A. M. Elaiw, I. A. Hassanien, and S. A. Azoz, “Global stability of HIV infection models with intracellular delays,” Journal of the Korean Mathematical Society, vol. 49, no. 4, pp. 779–794, 2012.

A. M. Elaiw and M. A. Alghamdi, “Global properties of virus dynamics models with multitarget cells and discrete-time delays,” Discrete Dynamics in Nature and Society, vol. 2011, Article ID 201274, 19 pages, 2011.

A. M. Elaiw, “Global properties of a class of virus infection models with multitarget cells,” Nonlinear Dynamics, vol. 69, no. 1-2, pp. 423–435, 2012.

A. M. Elaiw, “Global dynamics of an HIV infection model with two classes of target cells and distributed delays,” Discrete Dynamics in Nature and Society, vol. 2012, Article ID 253703, 13 pages, 2012.

A. U. Neumann, N. P. Lam, H. Dahari et al., “Hepatitis C viral dynamics in vivo and the antiviral efficacy of interferon-α therapy,” Science, vol. 282, no. 5386, pp. 103–107, 1998.

X. Song and A. U. Neumann, “Global stability and periodic solution of the viral dynamics,” Journal of Mathematical Analysis and Applications, vol. 329, no. 1, pp. 281–297, 2007.

L. Wang and R. Xu, “Mathematical analysis of an improved hepatitis B virus model,” International Journal of Biomathematics, vol. 5, no. 1, Article ID 1250006, 18 pages, 2012.

K. Wang, A. Fan, and A. Torres, “Global properties of an improved hepatitis B virus model,” Nonlinear Analysis: Real World Applications, vol. 11, no. 4, pp. 3131–3138, 2010.

J. Pang, J. Cui, and J. Hui, “The importance of immune responses in a model of hepatitis B virus,” Nonlinear Dynamics, vol. 67, no. 1, pp. 723–734, 2012.

W. Shaoli, F. Xinlong, and H. Yinnian, “Global asymptotical properties for a diffused HBV infection model with CTL immune response and nonlinear incidence,” Acta Mathematica Scientia B, vol. 31, no. 5, pp. 1959–1967, 2011.

C. Vargas-De-León, “Stability analysis of a model for HBV infection with cure of infected cells and intracellular delay,” Applied Mathematics and Computation, vol. 219, no. 1, pp. 389–398, 2012.

C. Vargas-De-Leon, “Analysis of a model for the dynamics of hepatitis B with noncytolytic loss of infected cells,” World Journal of Modelling and Simulation, vol. 8, no. 4, pp. 243–259, 2012.

J. Li, K. Wang, and Y. Yang, “Dynamical behaviors of an HBV infection model with logistic hepatocyte growth,” Mathematical and Computer Modelling, vol. 54, no. 1-2, pp. 704–711, 2011.

R. Xu and Z. Ma, “An HBV model with diffusion and time delay,” Journal of Theoretical Biology, vol. 257, no. 3, pp. 499–509, 2009.

S. Zhang and Y. Zhou, “The analysis and application of an HBV model,” Applied Mathematical Modelling, vol. 36, no. 3, pp. 1302–1312, 2012.

X. Tian and R. Xu, “Asymptotic properties of a hepatitis B virus infection model with time delay,” Discrete Dynamics in Nature and Society, vol. 2010, Article ID 182340, 21 pages, 2010.

M. Sheikhan and S. A. Ghoreishi, “Application of covariance matrix adaptation-evolution strategy to optimal control of hepatitis B infection,” Neural Computing and Applications, 2012.

C. L. Lai and M. F. Yuen, “The natural history and treatment of chronic hepatitis B: a critical evaluation of standard treatment criteria and end points,” Annals of Internal Medicine, vol. 147, no. 1, pp. 58–61, 2007.

M. Sheikhan and S. A. Ghoreishi, “Antiviral therapy using a fuzzy controller optimized by modified evolutionary algorithms: a comparative study,” Neural Computing and Applications, 2012.

K. Hattaf, M. Rachik, S. Saadi, and N. Yousfi, “Optimal control of treatment in a basic virus infection model,” Applied Mathematical Sciences, vol. 3, no. 17–20, pp. 949–958, 2009.

D. Q. Mayne, J. B. Rawlings, C. V. Rao, and P. O. M. Scokaert, “Constrained model predictive control: stability and optimality,” Automatica, vol. 36, no. 6, pp. 789–814, 2000.

A. M. Elaiw and X. Xia, “HIV dynamics: analysis and robust multirate MPC-based treatment schedules,” Journal of Mathematical Analysis and Applications, vol. 359, no. 1, pp. 285–301, 2009.

A. M. Elaiw and A. M. Shehata, “Stability and feedback stabilization of HIV infection model with two classes of target cells,” Discrete Dynamics in Nature and Society, vol. 2012, Article ID 963864, 20 pages, 2012.

É. Gyurkovics and A. M. Elaiw, “A stabilizing sampled-data -step receding horizon control with application to a HIV/AIDS model,” Differential Equations and Dynamical Systems, vol. 14, no. 3-4, pp. 323–352, 2006.

R. Zurakowski and A. R. Teel, “A model predictive control based scheduling method for HIV therapy,” Journal of Theoretical Biology, vol. 238, no. 2, pp. 368–382, 2006.

M. A. L. Caetano and T. Yoneyama, “Short and long period optimization of drug doses in the treatment of AIDS,” Anais da Academia Brasileira de Ciencias, vol. 74, no. 3, pp. 379–392, 2002.

É. Gyurkovics and A. M. Elaiw, “Stabilization of sampled-data nonlinear systems by receding horizon control via discrete-time approximations,” Automatica, vol. 40, no. 12, pp. 2017–2028, 2004.

A. M. Elaiw, “Multirate sampling and input-to-state stable receding horizon control for nonlinear systems,” Nonlinear Analysis: Theory, Methods & Applications, vol. 67, no. 5, pp. 1637–1648, 2007.