State estimation and optimal long period clinical treatment of HIV seropositive patients

MATHEMATICAL SCIENCES

State estimation and optimal long period clinical treatment of HIV seropositive patients

Juliana M. Grégio^I; Marco A.L. Caetano^II; Takashi Yoneyama^I

^IDivisão de Engenharia Eletrônica, ITA, Praça Mal. Eduardo Gomes, 50, Vila das Acácias 12228-900 São José dos Campos, SP, Brasil

^IIIBMEC-SP, Rua Quatá, 300, Vila Olímpia, 04546-042 São Paulo, SP, Brasil

^{Correspondence to}

ABSTRACT

Optimal control theory provides a very interesting quantitative method that can be used to assist the decision making process in several areas of application, such as engineering, biology, economics and sociology. The main idea is to determine the values of the manipulated variables, such as drug doses, so that some cost function is minimized, subject to physical constraints. In this work, the cost function reflects the number of CD4⁺T cells, viral particles and the drug doses. It is worth noticing that high drug doses are related to more intense side-effects, apart from the impact on the actual cost of the treatment. In a previous paper by the authors, the LQR – Linear Quadratic Regulator approach was proposed for the computation of long period maintenance doses for the drugs, which turns out to be of state feedback form. However, it is not practical to determine all the components of the state vector, due to the fact that infected and uninfected CD4⁺T cells are not microscopically distinguishable. In order to overcome this difficulty, this work proposes the use of Extended Kalman Filter to estimate the state, even though, because of the nonlinear nature of the involved state equations, the separation principle may not be valid. Extensive simulations were then carried out to investigate numerically if the control strategy consisting of the feedback of estimated states yielded satisfactory clinical results.

Key words: AIDS, filtering, drugs, mathematical modeling, optimization.

RESUMO

A teoria de controle ótimo apresenta um método quantitativo muito interessante que pode ajudar no processo de tomada de decisão em algumas áreas de aplicação, tais como engenharia, biologia, economia e sociologia. A principal idéia é determinar os valores das variáveis controladas, tais como doses de medicamentos, onde alguma função-custo é minimizada, sujeito às restrições físicas. Neste trabalho, a função-custo reflete o número de células CD4⁺T, partículas virais e doses de medicamentos. É fato que altas dosagens de medicamentos estão relacionadas à maior intensidade de efeitos colaterais, além do impacto no custo real do tratamento. Num prévio trabalho nosso, foi proposta a abordagem LQR - Regulador Linear Quadrático para o cálculo das doses de manutenção para os medicamentos, as quais dependiam de ser realimentadas pelo estado. Entretanto, a determinação de todos os componentes do vetor de estado não seria prática, devido ao fato de que células infectadas e não infectadas são indistingüíveis no microscópio. Para contornar essa dificuldade, este trabalho propõe o uso do Filtro de Kalman Estendido para estimar o estado, ainda que, devido à natureza não linear das equações de estado envolvidas, o princípio da separação não seja válido. Simulações extensivas foram realizadas para investigar numericamente se a estratégia de controle consistindo da realimentação de estados estimados produz resultados clínicos satisfatórios.

Palavras-chave: AIDS, filtragem, medicamentos, modelagem matemática, otimização.

INTRODUCTION

Several mathematical models that describe the clinical evolution of HIV seropositive patients can be found in the literature (Perelson 1989, Nowak and Bangham 1996, Tan and Wu 1998). This work is based on a simplification of the model that was originally proposed by (Tan and Wu 1998). This model has four differential equations representing the uninfected CD4⁺T cells, latent infected CD4⁺T cells, active infected CD4⁺T cells and free viruses.

One way to combat viruses is to obstruct their multiplication, ideally up to the point of becoming undetectable. In case of AIDS, one can associate reverse transcriptase inhibitors and protease inhibitors. Reserve transcriptase is an enzyme required by the HIV virus to translate its RNA genetic material into DNA, thus allowing its integration to the host cell. The AZT, the DDI, the 3TC, the D4T, and the Abacavir belong to this category. The protease inhibitors affect the enzyme responsible for the maturation of the HIV virus producing defective copies of it itself, incapable to infect new cells. The Saquinavir, the Indinavir, the Ritonavir, the Nelfinavir, and the Amprenavir belong to this group.

Although antiviral drugs can reduce the viral load, they can also cause some side effects such as diarrhea, vomits, nauseas, cutaneous rash, abdominal girth and abdominal fullness, among others (Mittler et al. 1998).

The idea in the present work is to optimize the doses of reverse transcriptase inhibitors and protease inhibitors by using the optimal control theory, with a view on minimizing the side effects of the chemical treatment, while, at the same time, obstructing the advance of virus in the organism. More specifically, the Linear Quadratic Regulator (LQR) is proposed for the determination of the maintenance doses for long period treatments. The solution to the LQR problem is well known and turns out to be of state feedback form. The numbers of uninfected and infected CD4⁺T cells, as well as the viral load, are components of the state vector. Because uninfected and infected CD4⁺T cells are not distinguishable in the microscope, the optimal control law presented in Caetano and Yoneyama 2001 is not readily implementable. In order to obtain estimates of these components of the state, Extended Kalman Filter (EKF)can be used (Maybeck 1982). Because the state equations are nonlinear, the separation principle may not be valid, so that combining the optimal feedback control law and the state estimator is a heuristic procedure. Hence, in order to verify whether the obtained solution is adequate, extensive numerical simulations were carried out.

MATERIALS AND METHODS

THE DYNAMIC MODEL

The mathematical model used in this work is a simplification of a more general one presented by (Tan and Wu 1998) that includes stochastic terms. The model consists of four coupled ordinary differential equations, given by:

where

with x₁(t) = uninfected cells CD4⁺T; x₂(t) = latently infected cells CD4⁺T; x₃(t) = active infected cells CD4⁺T; x₄(t) = free virus HIV; s: rate of generation of x₁ from precursors; r: rate of stimulated growth of x₁; T_max: maximum CD4⁺T cells population level; µ₁: death rate of x₁; µ₂: death rate of x₂; µ₃: death rate of x₃; µ_v: death rate of x₄; k₁: infection rate; k₂: conversion rate from x₂ to x₃; N: number of infectious virus produced by an actively infected CD4⁺T cell and that depends on constants β₁, β₂ and N₀; m₁ and m₂: doses of transcriptase and protease inhibitors, with effects characterized by pairs of constants (k₁₀, ₁) and (k₂₀, ₂), respectively. According to (Tan and Wu 1998), in the presence of antigen and HIV, x₁ cells are stimulated to proliferate to generate new x₁ cells with rate λ(x₁, x₂, x₃). Without the presence of HIV, the rate of generation is S(x₄). In the presence of free HIV (x₄), uninfected cells x₁ can be infected to become x₂ cells and x₃ cells, depending of probability of cells become actively or latently infected with rate ω. The x₂ cells can be activated to become x₃ cells. The activation rate is k₂. The x₃ cells are short lived and will normally be killed upon activation with death rate µ₃. The x₁, x₂ cells and x₄ free virus have finite life and the death rate in this model is µ₁, µ₂ and µ_v respectively. When a x₃ cell dies it will release N(t) free viruses x₄.

For the sake of notational simplicity, let x(t) = [x₁ (t) x₂ (t) x₃ (t) x₄ (t)]^T, m(t) = [m₁(t)m₂(t)]^T and write (1) as

where f : R⁴ × R² × R⁺→ R⁴ is the left side of (1), together with equations (2)-(5).

THE LINEAR QUADRATIC REGULATOR (LQR)

Once the patient's state reaches a clinically satisfactory region after an initial treatment with constant doses of reverse transcriptase and protease inhibitors, the idea is to switch the scheme and use drugs doses determined by solving a LQR problem.

The original state equations (1) are of non-linear type. Hence, a linearization around a nominal state is required. Moreover, the solution x(t) of the original state equations (1) is continuous in time. It is necessary, therefore, to discretize the state equations, adopting a sampling period T, which is compatible with the usual clinical procedures, as the medications are offered at regularly spaced time intervals, rather that continuously in time.

Let x_n(t) denote the nominal state vector, corresponding to a treatment with nominal constant dose of drugs m_n and denote

The linearized state equations are

where the matrices A(x_n, m_n)_4×4 and B(x_n, m_n)_4×2 are given by

and a₁₁ is just a shorthand for

Having obtained the matrices A(x_n,m_n) and B(x_n,m_n), it is straightforward to compute the discretized version A_d(x_n,m_n) and B_d(x_n,m_n), such that perturbations around nominal values {x_n(k),m_n} are now described by a linear difference equation

using a slight abuse of notation by writing Δx(k) in place of Δx(kT) and Δm(k) instead of Δm(kT). The matrices A_d(x_n,m_n) and B_d(x_n,m_n) satisfy

Note that it is possible to write the first few terms of the infinite series expansions for the sampled matrices.

Now, the LQR problem is posed as the determination of the optimal strategy Δm^* (k), that minimizes the performance index J:

where Q is required to be a real symmetric positive semi-definite matrix and R a real symmetric positive definite matrix. The matrix Q is the weight on the states Δx and the matrix R is the weight on the control variable Δm, so that the objective is to keep the actual state x close to the nominal value x_n, i.e. small Δx, with moderate adjustment of drug doses Δm.

The solution of the LQR problem is well known and can be found, for instance, in the classical texts such as (Lewis 1986, Kwakernaak and Sivan 1972, Kirk 1970):

where the feedback gain matrix F is given by

and

THE EXTENDED KALMAN FILTER

Because Δx can not be directly measured, as x₁, x₂ and x₃ are indistinguishable under a microscope, it is necessary to estimate Δx using the Extended Kalman Filter. Consider that the available measurements are

i.e., one can only observe the total number of CD4⁺T cells (x₁+ x₂+ x₃) and count the viral particles x₄. In matrix notation and considering stochastic uncertainty in the observation process (uncertainty in the measurement), the equation for y(k) = [y₁(k) y₂(k)]^T becomes

where v(k) ~ N(0,Q_m) i.i.d. (independent, identically distributed normal random variables).

The Extended Kalman Filter can be implemented by decomposing it in two stages (see, for instance, Maybeck 1982 or Gelb 1980): one for prediction and the other for update of the estimates. Assume that an estimate of the state is available (updated estimated of the state at time k, given the informations up to time k) with uncertainty represented by the covariance matrix P⁺(k) and let (k) be the estimates of x(k) with superscript ( )⁺ denoting updates considering the new observations, i.e. E[ x(k) | y(k) ] and also, analogously ^-(k) = E[ x(k) | y(k - 1) ]. Introduce the notation Δ⁺ (k) = ⁺ (k) - x_n (k) and Δ⁺ (k) = Δx (k) - Δ⁺ (k), so that one can write E[ Δ⁺(k) ] = 0 and E[ Δ (k) Δ^T (k) ] = P⁺ (k), k = 0, 1, 2, ...

In the propagation phase,

⁺ (k) is used as the initial condition to obtain Δ^- (k + 1) which satisfies the difference equation:

where G is a constant matrix to scale the disturbance noise and w is assumed to be N(0, Q_d) i.i.d. Note that the discretization is now carried out around the best estimates of the state, i.e., A_d (⁺ (k), m(k)) in place of A_d (x_n (k), m_n (k)). An analogous procedure is used in the computation of B_d. In what follows, the arguments of A_d and B_d are dropped to avoid cumbersome notation.

For the propagation of P⁺ (k) the following equation is used:

In the update phase, the estimate

^- (k + 1) is improved by aggregating the new information provided by Δy (k + 1) = y(k + 1) - Hx_n (k), yielding ⁺ (k + 1):

The gain K is optimal in the sense that minimizes the estimation error. It serves to make a compromise between the estimate given by the model ^- (k + 1) and the innovation y(k + 1) - H Δ(k + 1) which is the difference between the actual measurement and its expected value (Gelb 1980).

SIMULATION RESULTS AND DISCUSSION

In order to present a realistic example, the data belonging to Patient A in (Pontesilli et al. 1999) was used. The numerical values for the model parameters of Patient A are shown in Table I.

The weight matrices R and Q in the performance index J, were selected by a trial and error method, so that significant reduction in the overall drug doses were obtained in the long period treatment. In principle, an automated iterative procedure could be used to adjust R and Q. However, in this work, the values of R and Q were selected manually by subjectively evaluating the simulated results in terms of reduction in the doses. It is worth mentioning that for any choice of positive definite R and Q, the obtained feedback control law yields stable closed loop system:

The covariance matrices Q_d and Q_m corresponding to the uncertainty in the state w(k) and the uncertainty in the observations v(k) were set to

The choice of Q_d and Q_m were based on the dispersions observed in the data collected from clinical records obtained from "Centro de Referência e Treinamento em DST-AIDS, São Paulo, Brasil".

In the simulations the following scenario was adopted:

1) The initial treatment to bring the patient to a satisfactory health condition consisted of using constant doses of the drugs, namely m₁ = 900 mg and m₂ = 1200 mg. This initial treatment lasted for 224 days. The simulations were started with Patient A in a rather unfavorable condition with x₁(0) = 357, x₂(0) = 10, x₃(0) = 100 and x₄(0) = 133352, reflecting an advanced stage of AIDS. When the patient reaches a satisfactory clinical condition, he begins to use the doses computed by the proposed LQR formulation.

2) In the period t ∈ [225, 500] days, the doses given to the patient are those computed by the LQR problem. In Figure 1, one can notice that the CD4⁺T count (one of the observed signal) is significantly improved by using the Extended Kalman Filter. Therefore, not only the uninfected, infected in latent state and active infected CD4⁺T cells can be distinguished, but also the overall uncertainty is reduced.

In Figure 2, it is seen that the number of viruses is kept at an adequate low level by the doses of drugs computed using the LQR approach. The doses are slightly reduced under LQR control, as one can see in Figure 3. It is work recalling that up to time 224 days, the doses are constant (m₁ = 900 mg and m₂ = 1200 mg) and then, when LQR procedure is started, reductions of about 9.2% of m₁ and 8.0%of m₂ are obtained.

Of course, different choices of the weight matrices Q and R lead to different levels of reduction, so that the user can adjust them, for instance, if one drug is better accepted by the patient, or if it is more readily available.

Although only the case of Patient A is presented in this work, the methodology was tested with other individuals, with similar results. In each case, the main difficulty was to determine the model parameters from the clinical data.

CONCLUSIONS

The objective of this work was to present a quantitative method to assist the medical staff in conceiving clinical treatment schemes for patients with AIDS, where the positive results such as increase in the CD4⁺T count and reduction in the level of viral particles can be balanced with dose dependent side effect of the drugs. More specifically, the method proposed here complements a previous work in Caetano and Yoneyama 2001, where the LQR - Linear Quadratic Regulator approach was used for the computation of long period maintenance doses for the drugs, which turns out to be of state feedback form. However, because is not feasible to determine all the components of the state vector, due to the fact that infected and uninfected CD4⁺T cells are not microscopically distinguishable, a state estimator was developed. Considering the existence of uncertainties in the measurements, the Extended Kalman Filter was introduced to estimate the state, even though the separation principle may not be valid, because the state equations are non-linear. However, extensive numerical simulations that were carried out indicates that the control strategy consisting of the feedback of estimated states yielded satisfactory clinical results, with slight decrease in the drug doses, while maintaining the viral counts at an adequate low level. Although the proposed method would be difficult to be used routinely in practice because the determination of the model parameters for each individual patient is a formidable task, it is well suited for "what if" type analysis, so that medical staff, or even medical students, can easily visualize the effects of changing the relative importance of the terms appearing in the cost function (actually by manipulating the matrices Q and R), contributing to provide insight into the quite complicated dynamics involved in the clinical description of HIV seropositive patients.

ACKNOWLEDGMENTS

The authors are grateful to the Centro de Referência e Treinamento em DST-AIDS in São Paulo, SP, Brazil, in particular to Dra. Mylva Fonsi for her help in providing the actual data used in this work.

Manuscript received on November 5, 2007; accepted for publication on August 4, 2008; presented by DJAIRO G. FIGUEIREDO

• CAETANO MA AND YONEYAMA T. 2001. Short and long period optimization of drug doses in the treatment of AIDS. An Acad Bras Cienc 74: 379392.
• GELB A. 1980. Applied Optimal Estimation. MIT Press.
• KIRK DE. 1970. Optimal Control Theory: An Introduction. Prentice-Hall. New Jersey, USA, 452 p.
• KWAKERNAAK H AND SIVAN R. 1972. Linear Optimal Control Systems. New York, NY: USA, Wiley-Interscience (Eds).
• LEWIS FL. 1986. Optimal Control. J Wiley & Sons, New York, 362 p.
• MAYBECK PS. 1982. Stochastic Models, Estimation, and Control. Academic Press, New York, Vol. 2.
• MITTLER JE, SULZER B, NEUMANN AU AND PERELSON AS. 1998. Influence of Delayed Viral Production on Viral Dynamics in HIV-1 Infected Pacients. Math Biosci 152: 143163.
• NOWAK MA AND BANGHAM CRM. 1996. Population Dynamics of Immune Responses to Persistent Viruses. Science 272: 7419.
• PERELSON AS. 1989. Modeling the interaction of the immune system with HIV (in archive Mathematical and Statistical Approaches to AIDS Epidemiology). In: CASTILLO-CHAVEZ C (Ed), Springer-Verlag, Lect Notes Biomath, New York, USA 83: 350370.
• PONTESILLI O, KERKHOF-GARDE S, PAKKER NG, NOTERMANS DW, ROOS MTL, KLEIN MR, DANNER SA AND MIEDEMA F. 1999. Antigen-specific T-lymphocyte proliferative responses during highly active antiretroviral therapy (HAART) of HIV-1 infection. Immunol Lett 66: 213217.
• TAN WY AND WU H. 1998. Stochastic Modelling of the Dynamics of CD4⁺T-Cell Infection by HIV and Some Monte Carlo Studies. Math Biosci 147: 173203.

Publication Dates

• Publication in this collection
  19 Feb 2009
• Date of issue
  Mar 2009

History

• Received
  05 Nov 2007
• Accepted
  04 Aug 2008