- Email: [email protected]
MODELLING AND CONTROL OF HIV DYNAMICS
MODELLING AND CONTROL OF HIV DYNAMICS
Alberto Landi*, Alberto Mazzoldi°, Chiara Andreoni °, Matteo Bianchi°, Andrea Cavallini°, Leonardo Ricotti °, Luca Ceccherini Nelli#, Riccardo Iapoce§ Department of Electrical Systems and Automation, University of Pisa, via Diotisalvi 2, 56126, Pisa, Italy. e-mail: [email protected] ° Interdepartmental Research Center E.Piaggio, University of Pisa, Via Diotisalvi 2, 56126, Pisa, Italy # Department of Virology, Via San Zeno, University of Pisa § U.O. Malattie Infettive, Ospedale Cisanello – Azienda Ospedaliera Pisana, Via Paradisa 2–56100, Pisa, Italy
Abstract: Various models have been considered in literature for modelling HIV infection and evolution. This paper considers a modification of the Wodarz and Nowak model, in order to obtain a simple, but effective mathematical model useful to predict the effects of a therapeutic drug regimen. After presenting several simulations for illustrating the effectiveness of the considered model, an accurate analysis of the results and some new ideas for future studies are discussed. Copyright © 2006 IFAC Keywords: Physiological models, HIV, Biomedical system, Differential equations, Numerical simulations.
1. INTRODUCTION Over the last years HIV dynamical models have been the objects of an intensive research. Nevertheless "Human Immunodeficiency Virus," (HIV) still continues to be not completely modelled. Several aspects of the pathology have been identified and modelled in an effective way, but unfortunately other aspects (especially correlated to therapeutic effects) are under research activity and require more accurate experimental and theoretical evaluations. Dynamic interactions between viral infection and immune system are particularly complex (Wodarz and Nowak, 2000). Although the immune responsiveness is potentially able to attack the virus, HIV infection causes the depletion of the cell Helper T CD4+, which has a primary role in the generation of immune response. Moreover, HIV infection attacks other immune cells, as macrophages, dendrites cell, etc. Therefore, in the initial acute phase of this pathology, there is an immune response, which is sub-optimal for the intrinsic viral activity reducing
ability of the immune systems. This peculiar viral activity contributes to the permanence of the virus and to its very high mutation rate (Vergu, Mallet and Golmard, 2005). Some days after the burst of the viral aggression, there is a high increase of number of virus cells at lymph-nodes, with a peak of viremia. In the following 12 weeks the immune response is completed and the viremia decreases up to low values, which can be also below the measurement threshold (Zurakowski, and Teel, 2006). After the initial phase characterized by a strong reduction, also CD4+ cells come back to acceptable values. It starts a new phase of the infection, with a "clinical latency": in the lymph-nodes and spleen there is a continuous replication of the virus, with destruction of immune cells. In this phase, the immune system is apparently able to control the situation and to oppose its action of virus and of other opportunistic microbes, such that there is no particular evidence of the presence of the virus in the patients. Only after few years, a new acute phase bursts, with re-increasing of viremia and CD4+ decreasing. Primarily viral reservoirs (in which
virus cells are conserved without possibility to be destruct by immune system or by external drugs) cause re-emergence of the virus upon cessation of therapy, even after many years of suppression (D. Finzi et al, 1999). Modelling of reservoirs (Ortega and Landrove, 2000), (Di Mascio et al, 2004), at the moment is one of the more difficult aspects to be implemented in the mathematical models. In this context, the pharmacology therapy offers several interesting results to an increasing of the quality and of expectation of life for the patient. Two main classes of drugs are used to reduce the replication of the virus and to delay the progression of the pathology (see Adams et al. and references therein, 2005). The so-called highly active anti-retroviral therapy (HAART), is a combination therapy including - reverse transcriptase inhibitors, to stop the process of inverse transcription and to avoid to permit to the virus to infect other cells, - protease inhibitors, to stop the production of precursors of viral proteins avoid the assembling the productions of virions by infected cells. Some plain comments on the HAART can be introduced: infected cells have a half-life short (from less than one month to 6 month), but hidden reservoirs of virus contribute to an even slower disease phase (Brandt and Chen, 2001). These long lasting time phases of infection make complete eradication of the virus from the body impossible with current therapies. In addition, genetic modifications of the virus and its ability to change its response to drugs com-plicate the problem. Therefore possible mathematical model have to consider and quantify the 'force' of the virus and its response to the drugs. In this context a mathematical model able to balance the model complexity with a simple description of the viral dynamics is a difficult task to satisfy: low order models are usually too simple to be useful, on the other hand high order models are too complex for simulation purposes and have too many unknown parameters to be identified. Some interesting models have been proposed by Wodarz and Nowak (2000), they include state variables representing both the viral dynamics and the immune responses, in terms of precursors of cytotoxic T-lymphocytes (CTLp), responsible of the development of an immune memory, and effectors of cytotoxic T-lymphocytes (CTLe). This paper has the primary goal to improve the Wodarz and Novak model, in the ambitious attempt to add new information to the results of simulation useful to physicians for comparing different treatment policies such as when to start or to switch therapy. More in detail the main equations of the model are analyzed and discussed in various cases (healthy patients, HIV patients without and with drug treatments, evaluated on short and long times). With respect to other models, a new variable, denoted as “aggressiveness,” is considered for a best evaluation of therapeutic protocols. In order to obtain simulation results coherent with the medical practice, a strict cooperation with clinical researchers expert in HIV therapies was successfully considered.
2. MODELS OF WODARZ AND NOWAK The first model presented [9] considers three state variables (expressed as cell counts in blood per cubic millimetre) inside a whole body model. The model is mathematically described by: ⎧ x = λ − dx − β xv ⎪ ⎨ y = β xv − ay ⎪v = ky − uv ⎩
(1)
The first equation represents the dynamics of the concentration of healthy CD4+ cells (x); λ represents the rate (supposed constant) at which new CD4+ Tcells are generated. The death rate of healthy cells is d. In case of active HIV infection the concentration of healthy cells decreases proportionally to the product βxv, where β represents a coefficient depending to various factors, as the velocity of penetration of virus into the cells, frequency of encounters between uninfected cells free virus, etc. The second equation describes the dynamics of the concentration of infected CD4+ cells (y); β is the rate of infection, a is the death rate of infected cells. The last equation describes the concentration of free virions (v), they are produced by the infected cells at a rate k, u is the death rate of the virions. The therapy is considered with the following policy: reverse transcriptase inhibitors stop the infection in cells still healthy. If the drug efficacy is maximum and equal to 100%, and if the system is at equilibrium before the drug treatment, β is set to zero in the model. If the drug efficacy is imperfect, β is substituted by β* = s β, with s<1. Protease inhibitors need a different modelling, since they reduce infection of new cells but do not block production of viruses from cells already infected and model (1) must be suitably modified (Craig, and Xia, 2005). A more complex model was presented in (Wodarz, and Nowak, 1999) and augmented in (Wodarz and Nowak, 2000): it offers important theoretical insights into immune control of virus based on treatment strategies, which can be viewed as a fast subsystem of the dynamics of HIV infection. The final model of Wodarz and Nowak is mathematically described by: ⎧ x = λ − dx − β xv ⎪ y = β xv − ay − pyz ⎪⎪ ⎨v = ky − uv ⎪ w = cxyw − cqyw − bw ⎪ ⎪⎩ z = cqyw − hz
(2)
Two differential equations are added to (1) for describing the dynamics of precursors of cytotoxic Tlymphocytes CTLp (w), responsible of the development of an immune memory, and effectors of cytotoxic T-lymphocytes CTLe (z). This model can discriminate the trend of the infection as a function
of the rate of viral replication: if the rate is high a successful immune memory cannot establish, on the other hand if the replication rate is slow, a immune memory CTL helps the patient to successfully fight against the infection. A more detailed description of this model can be found both in the original papers and in the description of a modified model (2) in (Zurakowski, and Teel, 2006). 3. MODELS PROPOSED The model developed at our university is developed starting from the model (2), in order to reach the following objectives: the first one (not innovative) was to develop a model able to qualitatively describe experimental data existing in literature and synthesized in Fig. 1. The second one was to manage a model directly devoted for studying and predicting the effects of the therapies. The proposed model is: ⎧ x = λ − dx − rxy ⎪ y = rxy − ay − pyz ⎪⎪ ⎨ w = cxyw − cqyw − bw ⎪ z = cqyw − hz ⎪ ⎪⎩r = r0 − µ1 f1 − µ 2 f 2
v = ky − uv
(5)
can be added to model (3), although in the hypothesis of a constant coefficient k and since usually k>>u, the time response of viremia is almost proportional to the concentration of infected CD4+ cells. An advantage of the proposed model is its direct and simple interpretation in terms of the drug efficacy: it shows limitations due to its simplicity with respect to more complex equations presented in (e.g.) in Zurakowski and Teel (2006), or in Adams et al (2005), but it shows a good agreement with the results of HIV clinical researchers. 4. RESULTS OF SIMULATION TESTS
(3)
The proposed model, based on the final model of Wodarz and Nowak, considers the new state variable r as the intrinsic virulence of the virus; the constant β is substituted with the adimensional state variable r, index of the aggressiveness of the virus. This new equation increases linearly in case of an untreated HIV-infected individual, with a growth rate depending on the constant r0 (higher is r0, higher is the virulence growth rate). This equation has the typical structure: r = r0 − ∑ i µi f i
where µi represents a drug efficacy coefficient for the specific drug fi . In the following simulation tests the initial condition of the state variable r0 is set to 10, the different fi’s (drug concentrations in arbitrary units) are set to 1 (treated patients) or 0 (untreated patients). The equation
(4)
Fig. 1. Clinical behavior of HIV infection. Virus concentration and CD4+ concentration, in the case of an untreated HIV-infected individual (Abbas, Lichtman and Pober, 2005).
All the simulations have been performed using the Simulink environment of Matlab™. Chosen parameters are the same considered in Wodarz and Nowak models, unknown parameters are selected from heuristic considerations, in order to obtain time responses of the state variables physiologically consistent and respecting the findings of medical literature (e.g., the time responses shown in Fig.1). Simulations in different cases have been performed using adimensional values reported in Wodarz (1999): λ =1, d = 0.1, β = 0.5, a = 0.2, p = 1, c = 0.1, b = 0.01, q = 0.5, h = 0.1, s = 0.042). Consider now the simulation tests in the cases of: • untreated HIV-infected individual • treated HIV-infected individual In the following figures cases of untreated (from Fig. 2 to Fig. 6) and poorly treated patients (µ1 = µ2 = 0.1) (from Fig. 7 to Fig. 11) are shown. If coefficients µi increase (e.g., µi = 0.8) therapeutic action becomes more and more efficient (from Fig. 12 to Fig. 16). Modelling simulation hypotheses are the following: -the initial burst of the viral infection is simulated with a short impulse, acting as an external input on y and z variables. - the existence of an additional reservoir of CD4+ cells, produced by the bone marrow, has been considered as an external input added to the variable x (T-lymphocytes), which permits to have a fast recovering of CD4+ cells after the minimum value typical after the first 12 weeks of infection (this 'reservoir effect' allows the physiological rebound observed in the experimental studies). 4.1 Untreated HIV-infected individuals − T CD4+ healthy cells show a rebound typical of the acute infection phase, followed by a constant quasi-
Long time simulation 300
250 v RNA copies [Concentration unit]
homeostatic condition in the latency period. After the latency period the viral burst leads to almost total cell depletion during the last phase of the infection; − viremia shows a peak, corresponding to the initial infection, followed by a reduction during the clinical latency. A long time simulation shows a continuous increasing of viremia, often associated to opportunistic infections; − CTLp cells increase up to a maximum value, followed by a decreasing dynamics; − CTLe cells have the same dynamics of CTLp's, but delayed.
200
150
100
50
0 0
40
60 80 Time units
100
120
140
Fig. 4. Long time simulation of state variables v, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks).
Two combined drugs are considered. − T CD4+ healthy cells after the rebound phase tend to increase as in the physiological case (healthy people) when the therapy inhibits the viral replication; − viremia shows again the initial peak, corresponding to the initial infection, followed by a reduction near zero in case of a perfect (ideal) efficacy of the drugs; − CTLp and CTLe cells increase as in healthy people.
Long time simulation 18 16 w CTLp cells [Concentration unit]
4.2 Treated HIV-infected individuals
Short time simulation
14 12 10 8 6 4 2
14
0 12 x CD4+ cells [Concentration unit]
20
0
20
40
60 80 Time units
100
120
140
Fig. 5. Long time simulation of state variables w, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks).
10
8
6
Long time simulation
4
25
0 0
5
10
15 20 Time units
25
30
35
Fig. 2. Short time simulation of state variables x, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks). Long time simulation 14
x CD4+ cells [Concentration unit]
12
z CTLe cells [Concentration unit]
2
20
15
10
5
0 0
10
20
40
60 80 Time units
100
120
140
Fig. 6. Long time simulation of state variables z, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks).
8
6
4
5. CONCLUDING REMARKS
2
0 0
20
40
60 80 Time units
100
120
140
Fig. 3. Long time simulation of state variables x, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks).
The inclusion of virulence as a new state variable in the Wodarz and Nowak model allows an easy evaluation both in terms of the definition of a variable representing an index of the aggressiveness of the virus and the introduction of the drugs in a simple and direct way.
Long time simulation 25
w CTLp cells [Concentration unit]
A simplified simulation model may be a useful starting tool for the analysis and prediction of a therapy, especially in a complex case as the HAART, where the therapy needs to be optimized for reducing side effects, toxicity and resistance to medications. Future work will be devoted to a further testing phase of the model proposed in cooperation with clinical researchers, along with the attempt to adapt the coefficient of the model to single patients.
20
15
10
5
Short time simulation 14
0 0
8
40
60 80 Time units
100
120
140
6 Long time simulation
4
35
2
30
0 0
5
10
15 20 Time units
25
30
35
Fig. 7. Short time simulation of state variables x, vs. time in case of poorly treated (µ1=µ2=0.1) HIVinfected patient (5 time units = 12 weeks). Long time simulation
25
20
15
10
5
14
0 0
12 x CD4+ cells [Concentration unit]
20
Fig. 10. Long time simulation of state variables w, vs. time in case of poorly treated (µ1=µ2=0.1) HIV-infected patient (5 time units = 12 weeks).
10
z CTLe cells [Concentration unit]
x CD4+ cells [Concentration unit]
12
20
40
60 80 Time units
100
120
140
Fig. 11. Long time simulation of state variables z, vs. time in case of poorly treated (µ1=µ2=0.1) HIVinfected patient (5 time units = 12 weeks).
10
8
6
4
Short time simulation 14
2
20
40
60 80 Time units
100
120
140
Fig. 8. Long time simulation of state variables x, vs. time in case of poorly treated (µ1=µ2=0.1) HIVinfected patient (5 time units = 12 weeks).
Long time simulation
8
6
4
0 0
250
v RNA copies [Concentration unit]
10
2
300
5
10
15 20 Time units
25
30
35
Fig. 12. Short time simulation of state variables x, vs. time in case of treated (µ1=µ2=0.8) HIV-infected patient (5 time units = 12 weeks).
200
150
100
50
0 0
x CD4+ cells [Concentration unit]
12
0 0
20
40
60 80 Time units
100
120
140
Fig. 9. Long time simulation of state variables v, vs. time in case of poorly treated (µ1=µ2=0.1) HIVinfected patient (5 time units = 12 weeks).
4
Long time simulation 1200
4.5
Long time simulation
x 10
4 z CTLe cells [Concentration unit]
x CD4+ cells [Concentration unit]
1000
800
600
400
3.5 3 2.5 2 1.5 1
200 0.5 0 0
20
40
60 80 Time units
100
120
140
Fig. 13. Long time simulation of state variables x, vs. time in case of treated (µ1=µ2=0.8) HIV-infected patient (5 time units = 12 weeks).
0 0
20
40
60 80 Time units
100
120
140
Fig. 16. Long time simulation of state variables z, vs. time in case of treated (µ1=µ2=0.8) HIV-infected patient (5 time units = 12 weeks).
Long time simulation 300
v RNA copies [Concentration unit]
250
REFERENCES
200
150
100
50
0 0
20
40
60 80 Time units
100
120
140
Fig. 14. Long time simulation of state variables v, vs. time in case of treated (µ1=µ2=0.8) HIV-infected patient (5 time units = 12 weeks). 8
w CTLp cells [Concentration unit]
2.5
Long time simulation
x 10
2
1.5
1
0.5
0 0
20
40
60 80 Time units
100
120
140
Fig. 15. Long time simulation of state variables w, vs. time in case of treated (µ1=µ2=0.8) HIVinfected patient (5 time units = 12 weeks).
Abbas, A. K., Lichtman, A. H., Pober, J. S. (2005). Cellular and molecular immunology. Philadelphia: WB Saunders Co., IV Edition. Adams B.M. et al. (2005) HIV dynamics: modeling, data analysis, and optimal treatment protocols. Journal of Computational and Applied Mathematics, 184, 10-49. Brandt, M. E., Chen, G. , (2001). Feedback Control of a Biodynamical Model of HIV-1. IEEE Trans. on Biomedical Engineering, 48, 754,759. Craig, I., Xia, X. (2005). Can HIV/AIDS Be Controlled? IEEE Control Systems Magazine, 80-83. Finzi, D. et al. (1999). Latent infection of CD4+ T cells provides a mechanism for lifelong persistence of HIV-1, even in patients on effective combination therapy. Nat. Med. 5, 512–517. Di Mascio M. et al. (2004). Modeling the long-term control of viremia in HIV-1 infected patients treated with antiretroviral therapy. Mathematical Biosciences 188, 47–62. Nowak, M.A. Bangham , C.R.M. (1996). Population dynamics of immune responses to persistent viruses. Science, 272, 74–79. Ortega, H., Landrove, M. M. (2000). A Model for Continuously Mutant HIV-1. Proc. of the 22nd Annual EMBS International Conference, pp. 1917-1920, July 23-28, Chicago USA. Wodarz, D., Nowak M. A. (1999). Specific therapy regimes could lead to long-term immunological control of HIV. Proc. Natl. Acad. Sci. USA 96, 14464–14469. Wodarz, D., Nowak, M. A. (2000). Mathematical models of HIV pathogenesis and treatment. BioEssays 24, 1178–1187. Vergu, E., Mallet, A. , J.L. Golmard. (2005) A modeling approach to the impact of HIV mutations on the immune system. Computers in Biology and Medicine, 35 1–24. Zurakowski, R., Teel A.R. (2006) A model predictive control based scheduling method for HIV therapy. Journal of Theoretical Biology, 238, 368-382.
Alberto Landi*, Alberto Mazzoldi°, Chiara Andreoni °, Matteo Bianchi°, Andrea Cavallini°, Leonardo Ricotti °, Luca Ceccherini Nelli#, Riccardo Iapoce§ Department of Electrical Systems and Automation, University of Pisa, via Diotisalvi 2, 56126, Pisa, Italy. e-mail: [email protected] ° Interdepartmental Research Center E.Piaggio, University of Pisa, Via Diotisalvi 2, 56126, Pisa, Italy # Department of Virology, Via San Zeno, University of Pisa § U.O. Malattie Infettive, Ospedale Cisanello – Azienda Ospedaliera Pisana, Via Paradisa 2–56100, Pisa, Italy
Abstract: Various models have been considered in literature for modelling HIV infection and evolution. This paper considers a modification of the Wodarz and Nowak model, in order to obtain a simple, but effective mathematical model useful to predict the effects of a therapeutic drug regimen. After presenting several simulations for illustrating the effectiveness of the considered model, an accurate analysis of the results and some new ideas for future studies are discussed. Copyright © 2006 IFAC Keywords: Physiological models, HIV, Biomedical system, Differential equations, Numerical simulations.
1. INTRODUCTION Over the last years HIV dynamical models have been the objects of an intensive research. Nevertheless "Human Immunodeficiency Virus," (HIV) still continues to be not completely modelled. Several aspects of the pathology have been identified and modelled in an effective way, but unfortunately other aspects (especially correlated to therapeutic effects) are under research activity and require more accurate experimental and theoretical evaluations. Dynamic interactions between viral infection and immune system are particularly complex (Wodarz and Nowak, 2000). Although the immune responsiveness is potentially able to attack the virus, HIV infection causes the depletion of the cell Helper T CD4+, which has a primary role in the generation of immune response. Moreover, HIV infection attacks other immune cells, as macrophages, dendrites cell, etc. Therefore, in the initial acute phase of this pathology, there is an immune response, which is sub-optimal for the intrinsic viral activity reducing
ability of the immune systems. This peculiar viral activity contributes to the permanence of the virus and to its very high mutation rate (Vergu, Mallet and Golmard, 2005). Some days after the burst of the viral aggression, there is a high increase of number of virus cells at lymph-nodes, with a peak of viremia. In the following 12 weeks the immune response is completed and the viremia decreases up to low values, which can be also below the measurement threshold (Zurakowski, and Teel, 2006). After the initial phase characterized by a strong reduction, also CD4+ cells come back to acceptable values. It starts a new phase of the infection, with a "clinical latency": in the lymph-nodes and spleen there is a continuous replication of the virus, with destruction of immune cells. In this phase, the immune system is apparently able to control the situation and to oppose its action of virus and of other opportunistic microbes, such that there is no particular evidence of the presence of the virus in the patients. Only after few years, a new acute phase bursts, with re-increasing of viremia and CD4+ decreasing. Primarily viral reservoirs (in which
virus cells are conserved without possibility to be destruct by immune system or by external drugs) cause re-emergence of the virus upon cessation of therapy, even after many years of suppression (D. Finzi et al, 1999). Modelling of reservoirs (Ortega and Landrove, 2000), (Di Mascio et al, 2004), at the moment is one of the more difficult aspects to be implemented in the mathematical models. In this context, the pharmacology therapy offers several interesting results to an increasing of the quality and of expectation of life for the patient. Two main classes of drugs are used to reduce the replication of the virus and to delay the progression of the pathology (see Adams et al. and references therein, 2005). The so-called highly active anti-retroviral therapy (HAART), is a combination therapy including - reverse transcriptase inhibitors, to stop the process of inverse transcription and to avoid to permit to the virus to infect other cells, - protease inhibitors, to stop the production of precursors of viral proteins avoid the assembling the productions of virions by infected cells. Some plain comments on the HAART can be introduced: infected cells have a half-life short (from less than one month to 6 month), but hidden reservoirs of virus contribute to an even slower disease phase (Brandt and Chen, 2001). These long lasting time phases of infection make complete eradication of the virus from the body impossible with current therapies. In addition, genetic modifications of the virus and its ability to change its response to drugs com-plicate the problem. Therefore possible mathematical model have to consider and quantify the 'force' of the virus and its response to the drugs. In this context a mathematical model able to balance the model complexity with a simple description of the viral dynamics is a difficult task to satisfy: low order models are usually too simple to be useful, on the other hand high order models are too complex for simulation purposes and have too many unknown parameters to be identified. Some interesting models have been proposed by Wodarz and Nowak (2000), they include state variables representing both the viral dynamics and the immune responses, in terms of precursors of cytotoxic T-lymphocytes (CTLp), responsible of the development of an immune memory, and effectors of cytotoxic T-lymphocytes (CTLe). This paper has the primary goal to improve the Wodarz and Novak model, in the ambitious attempt to add new information to the results of simulation useful to physicians for comparing different treatment policies such as when to start or to switch therapy. More in detail the main equations of the model are analyzed and discussed in various cases (healthy patients, HIV patients without and with drug treatments, evaluated on short and long times). With respect to other models, a new variable, denoted as “aggressiveness,” is considered for a best evaluation of therapeutic protocols. In order to obtain simulation results coherent with the medical practice, a strict cooperation with clinical researchers expert in HIV therapies was successfully considered.
2. MODELS OF WODARZ AND NOWAK The first model presented [9] considers three state variables (expressed as cell counts in blood per cubic millimetre) inside a whole body model. The model is mathematically described by: ⎧ x = λ − dx − β xv ⎪ ⎨ y = β xv − ay ⎪v = ky − uv ⎩
(1)
The first equation represents the dynamics of the concentration of healthy CD4+ cells (x); λ represents the rate (supposed constant) at which new CD4+ Tcells are generated. The death rate of healthy cells is d. In case of active HIV infection the concentration of healthy cells decreases proportionally to the product βxv, where β represents a coefficient depending to various factors, as the velocity of penetration of virus into the cells, frequency of encounters between uninfected cells free virus, etc. The second equation describes the dynamics of the concentration of infected CD4+ cells (y); β is the rate of infection, a is the death rate of infected cells. The last equation describes the concentration of free virions (v), they are produced by the infected cells at a rate k, u is the death rate of the virions. The therapy is considered with the following policy: reverse transcriptase inhibitors stop the infection in cells still healthy. If the drug efficacy is maximum and equal to 100%, and if the system is at equilibrium before the drug treatment, β is set to zero in the model. If the drug efficacy is imperfect, β is substituted by β* = s β, with s<1. Protease inhibitors need a different modelling, since they reduce infection of new cells but do not block production of viruses from cells already infected and model (1) must be suitably modified (Craig, and Xia, 2005). A more complex model was presented in (Wodarz, and Nowak, 1999) and augmented in (Wodarz and Nowak, 2000): it offers important theoretical insights into immune control of virus based on treatment strategies, which can be viewed as a fast subsystem of the dynamics of HIV infection. The final model of Wodarz and Nowak is mathematically described by: ⎧ x = λ − dx − β xv ⎪ y = β xv − ay − pyz ⎪⎪ ⎨v = ky − uv ⎪ w = cxyw − cqyw − bw ⎪ ⎪⎩ z = cqyw − hz
(2)
Two differential equations are added to (1) for describing the dynamics of precursors of cytotoxic Tlymphocytes CTLp (w), responsible of the development of an immune memory, and effectors of cytotoxic T-lymphocytes CTLe (z). This model can discriminate the trend of the infection as a function
of the rate of viral replication: if the rate is high a successful immune memory cannot establish, on the other hand if the replication rate is slow, a immune memory CTL helps the patient to successfully fight against the infection. A more detailed description of this model can be found both in the original papers and in the description of a modified model (2) in (Zurakowski, and Teel, 2006). 3. MODELS PROPOSED The model developed at our university is developed starting from the model (2), in order to reach the following objectives: the first one (not innovative) was to develop a model able to qualitatively describe experimental data existing in literature and synthesized in Fig. 1. The second one was to manage a model directly devoted for studying and predicting the effects of the therapies. The proposed model is: ⎧ x = λ − dx − rxy ⎪ y = rxy − ay − pyz ⎪⎪ ⎨ w = cxyw − cqyw − bw ⎪ z = cqyw − hz ⎪ ⎪⎩r = r0 − µ1 f1 − µ 2 f 2
v = ky − uv
(5)
can be added to model (3), although in the hypothesis of a constant coefficient k and since usually k>>u, the time response of viremia is almost proportional to the concentration of infected CD4+ cells. An advantage of the proposed model is its direct and simple interpretation in terms of the drug efficacy: it shows limitations due to its simplicity with respect to more complex equations presented in (e.g.) in Zurakowski and Teel (2006), or in Adams et al (2005), but it shows a good agreement with the results of HIV clinical researchers. 4. RESULTS OF SIMULATION TESTS
(3)
The proposed model, based on the final model of Wodarz and Nowak, considers the new state variable r as the intrinsic virulence of the virus; the constant β is substituted with the adimensional state variable r, index of the aggressiveness of the virus. This new equation increases linearly in case of an untreated HIV-infected individual, with a growth rate depending on the constant r0 (higher is r0, higher is the virulence growth rate). This equation has the typical structure: r = r0 − ∑ i µi f i
where µi represents a drug efficacy coefficient for the specific drug fi . In the following simulation tests the initial condition of the state variable r0 is set to 10, the different fi’s (drug concentrations in arbitrary units) are set to 1 (treated patients) or 0 (untreated patients). The equation
(4)
Fig. 1. Clinical behavior of HIV infection. Virus concentration and CD4+ concentration, in the case of an untreated HIV-infected individual (Abbas, Lichtman and Pober, 2005).
All the simulations have been performed using the Simulink environment of Matlab™. Chosen parameters are the same considered in Wodarz and Nowak models, unknown parameters are selected from heuristic considerations, in order to obtain time responses of the state variables physiologically consistent and respecting the findings of medical literature (e.g., the time responses shown in Fig.1). Simulations in different cases have been performed using adimensional values reported in Wodarz (1999): λ =1, d = 0.1, β = 0.5, a = 0.2, p = 1, c = 0.1, b = 0.01, q = 0.5, h = 0.1, s = 0.042). Consider now the simulation tests in the cases of: • untreated HIV-infected individual • treated HIV-infected individual In the following figures cases of untreated (from Fig. 2 to Fig. 6) and poorly treated patients (µ1 = µ2 = 0.1) (from Fig. 7 to Fig. 11) are shown. If coefficients µi increase (e.g., µi = 0.8) therapeutic action becomes more and more efficient (from Fig. 12 to Fig. 16). Modelling simulation hypotheses are the following: -the initial burst of the viral infection is simulated with a short impulse, acting as an external input on y and z variables. - the existence of an additional reservoir of CD4+ cells, produced by the bone marrow, has been considered as an external input added to the variable x (T-lymphocytes), which permits to have a fast recovering of CD4+ cells after the minimum value typical after the first 12 weeks of infection (this 'reservoir effect' allows the physiological rebound observed in the experimental studies). 4.1 Untreated HIV-infected individuals − T CD4+ healthy cells show a rebound typical of the acute infection phase, followed by a constant quasi-
Long time simulation 300
250 v RNA copies [Concentration unit]
homeostatic condition in the latency period. After the latency period the viral burst leads to almost total cell depletion during the last phase of the infection; − viremia shows a peak, corresponding to the initial infection, followed by a reduction during the clinical latency. A long time simulation shows a continuous increasing of viremia, often associated to opportunistic infections; − CTLp cells increase up to a maximum value, followed by a decreasing dynamics; − CTLe cells have the same dynamics of CTLp's, but delayed.
200
150
100
50
0 0
40
60 80 Time units
100
120
140
Fig. 4. Long time simulation of state variables v, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks).
Two combined drugs are considered. − T CD4+ healthy cells after the rebound phase tend to increase as in the physiological case (healthy people) when the therapy inhibits the viral replication; − viremia shows again the initial peak, corresponding to the initial infection, followed by a reduction near zero in case of a perfect (ideal) efficacy of the drugs; − CTLp and CTLe cells increase as in healthy people.
Long time simulation 18 16 w CTLp cells [Concentration unit]
4.2 Treated HIV-infected individuals
Short time simulation
14 12 10 8 6 4 2
14
0 12 x CD4+ cells [Concentration unit]
20
0
20
40
60 80 Time units
100
120
140
Fig. 5. Long time simulation of state variables w, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks).
10
8
6
Long time simulation
4
25
0 0
5
10
15 20 Time units
25
30
35
Fig. 2. Short time simulation of state variables x, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks). Long time simulation 14
x CD4+ cells [Concentration unit]
12
z CTLe cells [Concentration unit]
2
20
15
10
5
0 0
10
20
40
60 80 Time units
100
120
140
Fig. 6. Long time simulation of state variables z, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks).
8
6
4
5. CONCLUDING REMARKS
2
0 0
20
40
60 80 Time units
100
120
140
Fig. 3. Long time simulation of state variables x, vs. time in case of untreated HIV-infected patient (5 time units = 12 weeks).
The inclusion of virulence as a new state variable in the Wodarz and Nowak model allows an easy evaluation both in terms of the definition of a variable representing an index of the aggressiveness of the virus and the introduction of the drugs in a simple and direct way.
Long time simulation 25
w CTLp cells [Concentration unit]
A simplified simulation model may be a useful starting tool for the analysis and prediction of a therapy, especially in a complex case as the HAART, where the therapy needs to be optimized for reducing side effects, toxicity and resistance to medications. Future work will be devoted to a further testing phase of the model proposed in cooperation with clinical researchers, along with the attempt to adapt the coefficient of the model to single patients.
20
15
10
5
Short time simulation 14
0 0
8
40
60 80 Time units
100
120
140
6 Long time simulation
4
35
2
30
0 0
5
10
15 20 Time units
25
30
35
Fig. 7. Short time simulation of state variables x, vs. time in case of poorly treated (µ1=µ2=0.1) HIVinfected patient (5 time units = 12 weeks). Long time simulation
25
20
15
10
5
14
0 0
12 x CD4+ cells [Concentration unit]
20
Fig. 10. Long time simulation of state variables w, vs. time in case of poorly treated (µ1=µ2=0.1) HIV-infected patient (5 time units = 12 weeks).
10
z CTLe cells [Concentration unit]
x CD4+ cells [Concentration unit]
12
20
40
60 80 Time units
100
120
140
Fig. 11. Long time simulation of state variables z, vs. time in case of poorly treated (µ1=µ2=0.1) HIVinfected patient (5 time units = 12 weeks).
10
8
6
4
Short time simulation 14
2
20
40
60 80 Time units
100
120
140
Fig. 8. Long time simulation of state variables x, vs. time in case of poorly treated (µ1=µ2=0.1) HIVinfected patient (5 time units = 12 weeks).
Long time simulation
8
6
4
0 0
250
v RNA copies [Concentration unit]
10
2
300
5
10
15 20 Time units
25
30
35
Fig. 12. Short time simulation of state variables x, vs. time in case of treated (µ1=µ2=0.8) HIV-infected patient (5 time units = 12 weeks).
200
150
100
50
0 0
x CD4+ cells [Concentration unit]
12
0 0
20
40
60 80 Time units
100
120
140
Fig. 9. Long time simulation of state variables v, vs. time in case of poorly treated (µ1=µ2=0.1) HIVinfected patient (5 time units = 12 weeks).
4
Long time simulation 1200
4.5
Long time simulation
x 10
4 z CTLe cells [Concentration unit]
x CD4+ cells [Concentration unit]
1000
800
600
400
3.5 3 2.5 2 1.5 1
200 0.5 0 0
20
40
60 80 Time units
100
120
140
Fig. 13. Long time simulation of state variables x, vs. time in case of treated (µ1=µ2=0.8) HIV-infected patient (5 time units = 12 weeks).
0 0
20
40
60 80 Time units
100
120
140
Fig. 16. Long time simulation of state variables z, vs. time in case of treated (µ1=µ2=0.8) HIV-infected patient (5 time units = 12 weeks).
Long time simulation 300
v RNA copies [Concentration unit]
250
REFERENCES
200
150
100
50
0 0
20
40
60 80 Time units
100
120
140
Fig. 14. Long time simulation of state variables v, vs. time in case of treated (µ1=µ2=0.8) HIV-infected patient (5 time units = 12 weeks). 8
w CTLp cells [Concentration unit]
2.5
Long time simulation
x 10
2
1.5
1
0.5
0 0
20
40
60 80 Time units
100
120
140
Fig. 15. Long time simulation of state variables w, vs. time in case of treated (µ1=µ2=0.8) HIVinfected patient (5 time units = 12 weeks).
Abbas, A. K., Lichtman, A. H., Pober, J. S. (2005). Cellular and molecular immunology. Philadelphia: WB Saunders Co., IV Edition. Adams B.M. et al. (2005) HIV dynamics: modeling, data analysis, and optimal treatment protocols. Journal of Computational and Applied Mathematics, 184, 10-49. Brandt, M. E., Chen, G. , (2001). Feedback Control of a Biodynamical Model of HIV-1. IEEE Trans. on Biomedical Engineering, 48, 754,759. Craig, I., Xia, X. (2005). Can HIV/AIDS Be Controlled? IEEE Control Systems Magazine, 80-83. Finzi, D. et al. (1999). Latent infection of CD4+ T cells provides a mechanism for lifelong persistence of HIV-1, even in patients on effective combination therapy. Nat. Med. 5, 512–517. Di Mascio M. et al. (2004). Modeling the long-term control of viremia in HIV-1 infected patients treated with antiretroviral therapy. Mathematical Biosciences 188, 47–62. Nowak, M.A. Bangham , C.R.M. (1996). Population dynamics of immune responses to persistent viruses. Science, 272, 74–79. Ortega, H., Landrove, M. M. (2000). A Model for Continuously Mutant HIV-1. Proc. of the 22nd Annual EMBS International Conference, pp. 1917-1920, July 23-28, Chicago USA. Wodarz, D., Nowak M. A. (1999). Specific therapy regimes could lead to long-term immunological control of HIV. Proc. Natl. Acad. Sci. USA 96, 14464–14469. Wodarz, D., Nowak, M. A. (2000). Mathematical models of HIV pathogenesis and treatment. BioEssays 24, 1178–1187. Vergu, E., Mallet, A. , J.L. Golmard. (2005) A modeling approach to the impact of HIV mutations on the immune system. Computers in Biology and Medicine, 35 1–24. Zurakowski, R., Teel A.R. (2006) A model predictive control based scheduling method for HIV therapy. Journal of Theoretical Biology, 238, 368-382.