Modelo de Intoxicação por Chumbo em Mamíferos

Prof. Doherty Andrade - www.metodosnumericos.com.br

Contexto Toxicológico

O chumbo (Pb) é um metal pesado que é de difícil eliminação em mamíferos, apresentando meia-vida estimada de 20-30 anos no tecido ósseo. Sua absorção ocorre predominantemente por:

A excreção ocorre principalmente por filtração renal (≈80%) e, em segundo lugar, por excreção capilar, e também sudorípara e biliar.

Abordagem por Sistemas Dinâmicos

Para simplificar o modelo, propomos um modelo matemático compartimental com três pools fisiológicos interconectados:

Compartimento 1 (x₁): Pool sanguíneo
Volume de distribuição: ≈5% massa corporal
Tempo de residência: Horas a dias
Compartimento 2 (x₂): Tecidos moles
Inclui: Fígado, rins, SNC
Dinâmica: Troca rápida com sangue
Compartimento 3 (x₃): Compartimento ósseo
Característica: Armazenamento profundo
Cinética: Troca lenta (constantes de tempo anuais)

Formulação Matemática do Modelo Compartimental

O sistema dinâmico é descrito pelas seguintes equações de balanço de massa. Aqui mantemos a notação original de Yeargers et al.:

\[ \begin{cases} \frac{dx_1}{dt} = -(a_{01} + a_{21} + a_{31})x_1 + a_{12}x_2 + a_{13}x_3 + I_L(t) & \text{(Sangue)} \\ \frac{dx_2}{dt} = a_{21}x_1 - (a_{02} + a_{12})x_2 & \text{(Tecidos moles)} \\ \frac{dx_3}{dt} = a_{31}x_1 - a_{13}x_3 & \text{(Osso)} \end{cases} \]
Parâmetros cinéticos e suas interpretações biológicas
Parâmetro Interpretação Fisiológica Valor Unidade
a01 Taxa de excreção renal 0.021 dia⁻¹
a02 Depuração hepática (fígado + bile) 0.162 dia⁻¹
a12 Transferência sangue → tecidos moles 0.0124 dia⁻¹
a21 Retorno tecidos moles → sangue 0.0111 dia⁻¹
a13 Deposição óssea 0.000035 dia⁻¹
a31 Reabsorção óssea 0.0039 dia⁻¹
IL(t) Taxa de ingestão/exposição 49.3 μg/dia

Convenções Importantes:

  • Notação aij: Fluxo do compartimento j para i (seguindo Yeargers)
  • ai0: Eliminação do compartimento i para o ambiente externo
  • IL(t): Função de entrada no compartimento sanguíneo

Observações:
1. Assume-se a32 = 0 (fluxo desprezível tecido→osso)
2. \(I_L\)(t) pode modelar exposição aguda (pulso) ou crônica (degrau)

Solução Numérica do Modelo Toxicológico

Método de Runge-Kutta 4.5 (RK4.5)

Para resolver numericamente o sistema de equações diferenciais acopladas, empregamos o método de Runge-Kutta 4.5, uma variante adaptativa do método Runge-Kutta. Esta abordagem é particularmente eficaz para sistemas toxicológicos onde diferentes compartimentos apresentam escalas temporais distintas (rápida dinâmica sanguínea vs. lenta acumulação óssea).

A solução formal do sistema \(\mathbf{\dot{X}} = \mathcal{A}\mathbf{X} + \mathbf{B}u(t)\) com condição inicial \(X_0\) é dada por: \[ \mathbf{X}(t) = e^{\mathcal{A}t}\mathbf{X}_0 + \int_0^t e^{\mathcal{A}(t-\tau)}\mathbf{B}u(\tau)d\tau. \] Requer cálculo da exponencial matricial e inversibilidade de \(\mathcal{A}\).

Para resolver o sistema de equações diferenciais ordinárias (EDOs) não homogêneo, como é o caso em questão:

\(\dot{\mathbf{X}}(t) = \mathcal{A} \mathbf{X}(t) + \mathbf{B} f(t)\),

onde:

  • \(\mathbf{X}(t)\) é o vetor de estado,
  • \(\mathcal{A}\) é uma matriz constante,
  • \(\mathbf{B}\) é um vetor constante,
  • \(f(t)\) é uma função de entrada (controle),

a solução geral pode ser obtida combinando a solução da parte homogênea com uma solução particular da equação não homogênea.

A solução geral da equação homogênea é dada por:

\(\mathbf{X}_h(t) = e^{\mathcal{A}t} \mathbf{C}\),

onde:

  • \(e^{\mathcal{A}t}\) é a matriz exponencial de \(\mathcal{A}\),
  • \(\mathbf{C}\) é um vetor constante determinado pelas condições iniciais.

Podemos encontrar uma solução particular usando o método de variação dos parâmetros. A solução particular é:

\(\mathbf{X}_p(t) = e^{\mathcal{A}t}\displaystyle \int_{t_0}^t e^{-\mathcal{A}\tau} \mathbf{B} f(\tau) \, d\tau\).

A solução geral do sistema não homogêneo é a soma da solução homogênea e da solução particular:

\(\mathbf{X}(t) = \mathbf{X}_h(t) + \mathbf{X}_p(t) = e^{\mathcal{A}t} \mathbf{C} + e^{\mathcal{A}t} \displaystyle\int_{t_0}^t e^{-\mathcal{A}\tau} \mathbf{B} f(\tau) \, d\tau\).

Se \(\mathbf{X}(t_0) = \mathbf{X}_0\) é a condição inicial, então:

\(\mathbf{C} = \mathbf{X}_0\),

e a solução formal fica:

\(\mathbf{X}(t) = e^{\mathcal{A}(t - t_0)} \mathbf{X}_0 + \displaystyle\int_{t_0}^t e^{\mathcal{A}(t - \tau)} \mathbf{B} f(\tau) \, d\tau\).

A matriz exponencial pode ser calculada usando:

  • Autovalores e autovetores (diagonalização, forma de Jordan),
  • Série de Taylor (para matrizes simples),
  • Transformada de Laplace (resolvendo \(e^{\mathcal{A}t} = \mathcal{L}^{-1}\{(sI - \mathcal{A})^{-1}\}\)).

A solução formal do sistema é:

\(\mathbf{X}(t) = e^{\mathcal{A}(t - t_0)} \mathbf{X}_0 +\displaystyle \int_{t_0}^t e^{\mathcal{A}(t - \tau)} \mathbf{B} f(\tau) \, d\tau\).

Se \(t_0 = 0\), simplifica-se para:

\(\mathbf{X}(t) = e^{\mathcal{A}t} \mathbf{X}_0 +\displaystyle \int_0^t e^{\mathcal{A}(t - \tau)} \mathbf{B} f(\tau) \, d\tau\).

Observações:

  • Se \(f(t) = 0\) (sistema homogêneo), a solução é apenas \(\mathbf{X}(t) = e^{\mathcal{A}t} \mathbf{X}_0\).
  • Se \(\mathcal{A}\) é diagonalizável, \(e^{\mathcal{A}t}\) pode ser calculada facilmente via \(e^{\mathcal{A}t} = P e^{Dt} P^{-1}\), onde \(D\) é a matriz diagonal dos autovalores.
  • Para entradas específicas (como degrau, seno, exponencial), a integral pode ser resolvida analiticamente.

Esse é o método geral para resolver sistemas lineares não homogêneos com coeficientes constantes. Se \(\mathcal{A}\) ou \(f(t)\) têm formas específicas, simplificações adicionais podem ser aplicadas.

Vantagens do RK4.5

  • Controle automático do passo de integração
  • Precisão adaptativa (tolerância relativa 10-6 no nosso caso)
  • Eficiência computacional para sistemas rígidos

Parâmetros Toxicológicos

Utilizamos os valores empíricos obtidos por Yearges et al. (1996) para as taxas de transferência na determinação da solução numérica apresentada a seguir:

Constantes Cinéticas do Modelo (unidades: dia-1)
Parâmetro Descrição Fisiológica Valor Intervalo Biológico
\(a_{01}\) Eliminação renal 0.021 0.015–0.025
\(a_{12}\) Distribuição sangue→tecidos moles 0.0124 0.010–0.015
\(a_{13}\) Deposição óssea 0.000035 0.00002–0.00005
\(a_{21}\) Retorno tecidos→sangue 0.0111 0.008–0.012
\(a_{02}\) Depuração hepática 0.162 0.15–0.18
\(a_{31}\) Reabsorção óssea 0.0039 0.003–0.005
\(I_L\) Ingestão diária (μg/dia) 49.3 10–100 (exposição ocupacional)

Desempenho Numérico

  • Tempo de simulação (365 dias): 2.8 ms (CPU Intel i7)
  • Erro relativo médio: 4.2×10-7
  • Passos computados: 187 (adaptativos)
49.3 μg/dia
0.0162 dia⁻¹

Interpretação Biológica

Observamos que ocorre estabilização dos valores no sangue e em tecidos moles, mas estes não retornam aos valores iniciais.

Referências

An Introduction to the Mathematics of Biology: With Computer Algebra Models
Yeargers, E. W., Shonkwiler, R. W., & Herod, J. V.
Birkhäuser, 1996.
ISBN: 978-0-8176-3941-8
Numerical Analysis
Burden, R. L., & Faires, J. D.
9th edition, Cengage Learning, 2011.
ISBN: 978-0-538-73351-9