Aula 7 - O modelo de Hodgkin-Huxley.ipynb

Kernel: SageMath 9.3

Modelos Bi-estáveis

O modelo de Hodgkin-Huxley

%display typeset

A transmissão de informações no sistema nervorsos dos animais se dá através de picos de potencial elétrico que se propagam ao longo da membrana de células especiais denominadas Neurônios. Estes potenciais elétricos, por sua vez, são chamados de potenciais de ação. Em um dos trabalhos mais célebres da biomatemática do século XX, Hodgkin e Huxley, nos anos 50, propuseram a partir de experimentos revolucionários, um modelo matemático para a transmissão de impulsos elétricos em neurônios.

 

O Potencial de ação é o resultados de uma orquestração delicada de canais iônicos na membrana do neurônio provocando correntes súbitas que se propagam como ondas.

 

O Modelo de Hodgkin e Huxley foi pioneiro em caracterizar este circuito elétrico celular:

A representação esquemática acima é equivalemte ao modelomatemático abaixo. A bi-camada lipídica, é representada como uma capacitância ($C_m$). Canais iônicos acionados por voltagem e canais de "vazamento" iônico, são representados pelas condutâncias $g_n$ (não-linear) e $g_l$ (linear), respectivamente. Os gradientes eletroquímicos que determinam o fluxo dos íons são representados por baterias ($E$). 

 

$$I = C_m\frac{{\mathrm d} V_m}{{\mathrm d} t} + \bar{g}_Kn^4(V_m - V_K) + \bar{g}_{Na}m^3h(V_m - V_{Na}) + \bar{g}_l(V_m - V_l),$$

$$\frac{dn}{dt} = \alpha_n(1 - n) - \beta_n n$$

$$\frac{dm}{dt} = \alpha_m(1 - m) - \beta_m m$$

$$\frac{dh}{dt} = \alpha_h(1 - h) - \beta_h h$$

$n,m,h$ são variáveis adimensionais, associadas com a ativação dos canais de sódio, ativação dos canais de potássio e a inativação dos canais de potássio, respectivamente. Os alfas e betas são taxas associadas com a velocidade dos canais iônicos e sãoafetados pela temperatura.

Vamos implementar o Modelo de Hodgkin e Huxley e examinar a sua dinâmica.

def fun(t,y, p):
    V,n,m,h = y
    VNa = 50
    VK = -77
    VL = -54.4
    gNa = 120
    gK = 36
    gL = .3
    C = 1
    I = p[0]
    an = .01*(V+55)/(1-exp(-(V+55)/10))
    bn = .125*exp(-(V+65)/80)
    am = .1*(V+40)/(1-exp(-(V+40)/10))
    bm = 4*exp(-(V+65)/18)
    ah = .07*exp(-(V+65)/20)
    bh = 1/(1+exp(-(V+35)/10))
    return [
    (I - gNa*h*(V-VNa)*m^3 - gK*(V-VK)*n^4 - gL*(V-VL))/C,
    an*(1-n) - bn*n,
    am*(1-m) - bm*m,
    ah*(1-h) - bh*h
    ]

T = ode_solver()
T.algorithm = "rk8pd"
T.function = fun
t_range = [0,50]
y0 = [-65, .317, .052, .596]

T.ode_solve(t_range,y0,num_points=2000, params=[15])

T.plot_solution(i=0,interpolate=True)

@interact
def integra(I=(0,100,10)):
    T.ode_solve(t_range,y0,num_points=500, params=[I])
    T.plot_solution(i=0,interpolate=True)