Equação de Fitzhugh-Nagumo

Grupo: Bernardo Boatini, Murilo Kessler Azambuja e Natália Ferrazzo

O objetivo deste trabalho é implementar e estudar a dinâmica do modelo FitzHungh-Nagumo, e das equações que o compõem, para potenciais de ação em células e tecidos excitáveis. Os métodos computacionais utilizado para resolver os problemas e implementar o modelo foi o FTCS (Forward Time Centered Space) e o método de Crank-Nicolson.

Potencial de Ação em Neurônios

A células vivas são sistemas eletricamente sensíveis, ou seja, podem reagir a estímulos elétricos. Isso se dá devido ao fato de que substâncias carregadas estão naturalmente vinculadas a seus processos internos de interação com o ambiente, principalmente por intermédio de canais iônicos e proteínas transmebrana como, por exemplo, a Bomba de Sódio e Potássio(Bomba Na⁺/K⁺ ATPase)[1].

Naturalmente todas as células vivas possuem um potencial de repouso(PR) elétrico, ou seja, uma diferença de potencial elétrico, em relação ao meio(cerca de 0,1); mantida por um equilíbrio químico de concentração de íons dentro e fora da membrana plasmática.

Existem células que reagem estímulos elétricos apenas reestabelecendo o PR original por transporte passivo(sem gasto de energia) através da membrana, e estas são ditas células não-excitáveis.

Por outro lado, existem células que sob a ação do mesmo estímulo produzem um tipo de resposta bem característica: potencial de ação(PA); um pulso elétrico intenso(capaz de inverter a polarização do Potencial de Membrana) que se propaga ao longo da membrana da célula, sustentado por uma cadeia de transportes ativos(com gasto de energia) e que não decai ao longo do tempo e espaço; a esse tipo de células damos o nome de excitáveis[1].

Os Neurônios são as células excitáveis do tecido nervoso(que constituem o encéfalo e medula espinhal, gânglios e nervos do reino animal) e com já vimos são capazes de gerar PA. Um potencial de ação pode assumir diversos formatos, mas ao longo do axônio(Figura 1) de um neurônio eles tendem a uma curva como a da Figura 2.

**Figura 1 -**Representação de um potencial de ação(vermelho) ao longo de um axônio de neurônio, partindo do soma neural em direção a arvore dentrítica.

**Figura 2 -**Curva de um Potencial de Ação genérico no tempo, em um ponto do axônio de um neurônio.

Olhando para Figura 2 vemos alguns aspectos importantes:

O potencial de ação necessita de um estímulo mínimo(limiar) para ser ativado, abaixo desse valor o estímulo decai como em uma célula não excitável;
Acima desse limiar a célula segue o principio de "Tudo ou Nada", ou seja, assume o valor máximo possivel dentro de sua capacidade, independente do estímulo aplicado;
A etapa de despolarização(crescimento) é brusca e varia mais rapidamente que a repolarização(decaimento);
O período que contém a repolarização e hiperpolarização da membrana é chamado período refratário, e se caracteriza por não permitir que ocorra nenhum disparo até que a membrana atinja o potencial de repouso.

Modelo

Equação de Nagumo

Para iniciar a modelagem do sistema, devemos antes enfatizar três condições básicas que o potencial deve obedecer para que seja um PA:

Deve existir um limiar de voltagem para que um estímulo desencadeie o PA;
Uma vez atingido o limiar, a voltagem deve aumentar até o máximo possível;
Caso o estímulo não atinja o limiar, ele deve desaparecer rapidamente.

Podemos então definir $υ$ como a variável normalizada que fará o papel do potencial elétrico no axônio, sendo que $υ = 0$ é definido como o potencial de repouso, $υ = 1$ é o estímulo máximo e $υ = α$ é o limiar de voltagem $(0 < α < 1)$ . Desta forma podemos escrever as condições acima como

${\begin{cases} S e υ > α \Rightarrow υ \to 1 \\ S e υ < α \Rightarrow υ \to 0 \end{cases}$

Pode-se dizer que estas condições impõem que $v = 0$ e $v = 1$ sejam pontos de equilíbrio estável, enquanto $v = α$ seja um ponto de equilíbrio instável.

Assim, a variação temporal do potencial elétrico na célula pode ser dada por

\frac{\partial υ}{\partial t} = f (υ) (1)

Onde $f (υ)$ é alguma função que faz $v$ satisfazer as condições do potencial.

No modelo de Nagumo, utiliza-se o polinômio de terceiro grau $f_{N} (υ) = - v (v - α) (v - 1)$ .

Substituindo $f_{N} (v)$ em $(1)$ , é fácil notar que $v = 0$ , $v = α$ e $v = 1$ são pontos de equilíbrio de $v$ . Por outro lado, como podemos ver na Figura 3, se $0 < v < α$ , temos que $\frac{\partial υ}{\partial t} < 0,$ logo o valor de $v$ diminui até atingir o ponto de equilíbrio em $v = 0$ . Se $α < v < 1,$ $\frac{\partial υ}{\partial t} > 0,$ fazendo o potencial $v$ crescer até atingir o valor máximo em $v = 1$ . Desta forma, temos que $v = 1$ e $v = 0$ são pontos de equilíbrio estável, enquanto que $v = α$ é um ponto de equilíbrio instável, assim como se esperava pelas condições do potencial.

**Figura 3** — Curva de $f_{N} (v) \times v$ com $α = 0.4$

Em outras palavras, vemos o caráter "tudo ou nada" do potencial quando utilizamos o polinômio $f_{N} (v)$ . Se o valor do estímulo inicial se encontrar entre $0$ e $α$ , o estímulo desaparece, porém se o estímulo inicial estiver entre $α$ e $1$ , o estímulo é amplificado até o valor máximo. Entretanto, neste modelo inicial o PA não é propagado por todo axônio, ele atua apenas localmente, onde ocorreu o estímulo. Além disso, uma vez estimulado o neurônio nunca volta para o seu estado de repouso, permanecendo no valor máximo de potencial.

Equação de Nagumo difusiva

Para resolver o problema do modelo

Método de Crank-Nicolson

(explicar o método --> Natália)

Equação de Nagumo Difusiva 1D

(aplicar o método na equação e testar estabilidade --> Natália)

Método FTCS

(explicar o método --> Murilo)

Equação de Recuperação 1D

(aplicar o método na equação e testar estabilidade --> Murilo))

Modelo FitzHung-Nagumo 2D

(aplicar o método na equação e testar estabilidade --> Bernardo)

O sistema de EDP's em 2 dimensões, assumindo uma difusão isotrópica, é dado por

$\frac{\partial v}{\partial t} = D (\frac{\partial^{2} v}{\partial x^{2}} + \frac{\partial^{2} v}{\partial y^{2}}) + f_{N} (v) - w$

$\frac{\partial w}{\partial t} = ϵ (v - γ w)$

Como a equação de recuperação já foi discretizada, e não depende da dimensão do problema, precisamos apenas aplicar o FTCS na equação difusiva de naumo em 2D, assumindo $Δ x = Δ y$ , temos

$\frac{v_{i, j}^{n + 1} - v_{i, j}^{n}}{Δ t} = \frac{D}{(Δ x)^{2}} [(v_{i - 1, j}^{n} - 2 v_{i, j}^{n} + v_{i + 1, j}^{n}) + (v_{i, j - 1}^{n} - 2 v_{i, j}^{n} + v_{i, j + 1}^{n})] + f_{N i, j}^{n} - w_{i, j}^{n}$

$v_{l, j}^{n + 1} = v_{l, j}^{n} + \frac{D Δ t}{(Δ x)^{2}} [(v_{l - 1, j}^{n} - 2 v_{l, j}^{n} + v_{l + 1, j}^{n}) + (v_{l, j - 1}^{n} - 2 v_{l, j}^{n} + v_{l, j + 1}^{n})] + Δ t f_{N l, j}^{n} - Δ t w_{l, j}^{n}$

$v_{l, j}^{n + 1} = v_{l, j}^{n} + \frac{D Δ t}{(Δ x)^{2}} [v_{l - 1, j}^{n} + v_{l + 1, j}^{n} + v_{l, j - 1}^{n} + v_{l, j + 1}^{n} - 4 v_{l, j}^{n}] + Δ t f_{N l, j}^{n} - Δ t w_{l, j}^{n}$

onde os índices $i$ e $j$ se referem às coordenadas $x$ e $y$ respectivamente; e o índice $n$ referente ao tempo.

Vamos verificar as condições de estabilidade do problema por modos de Fourier:

substituindo $v_{l, j}^{n} = A^{n} e^{i (x q_{x} + y q_{y})} = A^{n} e^{i β}$
usando $k = \frac{D Δ t}{(Δ x)^{2}}$
por simplicidade, como fizemos na seção 3.1, vamos linearizar $f_{N l, j}^{n} \to - α v_{l, j}^{n} = - α A^{n} e^{i β}$
por fim, desconsiderar $Δ t w_{l, j}^{n}$

$A^{n + 1} e^{i β} = A^{n} e^{i β} + k A^{n} [e^{i ((x - Δ x) q_{x} + y q_{y})} + e^{i ((x + Δ x) q_{x} + y q_{y})} + e^{i (x q_{x} + (y - Δ y) q_{y})} + e^{i (x q_{x} + (y + Δ y) q_{y})} - 4 e^{i β}] - Δ t α A^{n} e^{i β}$

$A^{n + 1} e^{i β} = A^{n} e^{i β} + k A^{n} e^{i β} [e^{- i Δ x q_{x}} + e^{i Δ x q_{x}} + e^{- i Δ y q_{y}} + e^{i Δ y q_{y}} - 4] - Δ t α A^{n} e^{i β}$

Dividindo os dois lados da equação por $A^{n} e^{i β}$ temos

$\frac{A^{n + 1}}{A^{n}} = 1 + k [e^{- i Δ x q_{x}} + e^{i Δ x q_{x}} + e^{- i Δ y q_{y}} + e^{i Δ y q_{y}} - 4] - Δ t α$

Sabendo as seguintes identidades trigonométricas:

$e^{- i θ} + e^{i θ} = 2 \cos (θ)$
$\cos (2 θ) = 1 - 2 \sin^{2} (θ)$

Aplicando consecutivamente a equação

$\frac{A^{n + 1}}{A^{n}} = 1 + 2 k [\cos (Δ x q_{x}) + \cos (Δ y q_{y}) - 2] - Δ t α$

$\frac{A^{n + 1}}{A^{n}} = 1 + 2 k [(1 - 2 \sin^{2} (Δ x q_{x} / 2)) + (1 - 2 \sin^{2} (Δ y q_{y} / 2)) - 2] - Δ t α$

$\frac{A^{n + 1}}{A^{n}} = 1 + 4 k [\sin^{2} (Δ x q_{x} / 2) + \sin^{2} (Δ y q_{y} / 2)] - Δ t α$

Sabemos que a condição de estabilidade por modos de Fourier é obtida quando $| \frac{A^{n + 1}}{A^{n}} | < 1$ , portanto

$| 1 - Δ t α + 4 k [\sin^{2} (Δ x q_{x} / 2) + \sin^{2} (Δ y q_{y} / 2)] | < 1$

Como $0 < Δ t < 1$ e $0 < α < 1$ , no pior dos casos os teremos de senos ao quadrado são 1, e a condição de estabilidade fica

$\frac{α Δ t}{8} < k < \frac{2 + α Δ t}{8}$

Resultados

Equação de Nagumo Difusiva 1D

(explicar os testes e os gráficos/animações --> Natália)

Modelo FitzHung-Nagumo 1D

(explicar os testes e os gráficos/animações --> Murilo)

Modelo FitzHung-Nagumo 2D

(explicar os testes e os gráficos/animações --> Bernardo)

Discussão

(Contextualizar resultados --> Bernardo)

Equação de Fitzhugh-Nagumo

Índice

Potencial de Ação em Neurônios

Modelo

Equação de Nagumo

Equação de Nagumo difusiva

Método de Crank-Nicolson

Equação de Nagumo Difusiva 1D

Método FTCS

Equação de Recuperação 1D

Modelo FitzHung-Nagumo 2D

Resultados

Equação de Nagumo Difusiva 1D

Modelo FitzHung-Nagumo 1D

Modelo FitzHung-Nagumo 2D

Discussão

Programas

Referências

Menu de navegação

Equação de Fitzhugh-Nagumo

Potencial de Ação em Neurônios

Modelo

Equação de Nagumo

Equação de Nagumo difusiva

Método de Crank-Nicolson

Equação de Nagumo Difusiva 1D

Método FTCS

Equação de Recuperação 1D

Modelo FitzHung-Nagumo 2D

Resultados

Equação de Nagumo Difusiva 1D

Modelo FitzHung-Nagumo 1D

Modelo FitzHung-Nagumo 2D

Discussão

Programas

Referências

Menu de navegação

Pesquisa