Fail to simulate hodgkin-huxley model

Support
1

Description of problem

failed to simulate hodgkin-huxley model using brian2, the voltage behaves weird

Minimal code to reproduce problem

sorry to say that the code has to be long, but I tried my best to annotate the code

from brian2 import *

start_scope()

# parameters
area = 20000 * umetre ** 2
Cm = 1 * ufarad * cm ** -2 * area
g_l = 0.3 * msiemens * cm ** -2 * area
El = -54 * mV
EK = -77 * mV
ENa = 50 * mV
g_Na = 120 * msiemens * cm ** -2 * area
g_k = 0.3 * msiemens * cm ** -2 * area

eqs = Equations('''
dv/dt = (I - (g_Na * m ** 3 * h * (v - ENa)) - (g_k * n ** 4 * (v - EK)) - (g_l * (v - El))) / Cm : volt

dn/dt = alpha_n * (1 - n) - beta_n * n : 1
dm/dt = alpha_m * (1 - m) - beta_m * m : 1
dh/dt = alpha_h * (1 - h) - beta_h * h : 1

alpha_n = (0.01 * (v / mV + 50) / (1 - exp(-(v / mV + 50) / 10))) / ms : Hz
beta_n = (0.125 * exp(-(v / mV + 60) / 80)) / ms : Hz

alpha_m = (0.1 * (v / mV + 35) / (1 - exp(-(v / mV + 35) / 10))) / ms : Hz
beta_m = (4.0 * exp(-0.0556 * (v / mV + 60))) / ms : Hz

alpha_h = (0.07 * exp(-0.05 * (v / mV + 60))) / ms : Hz
beta_h = (1 / (1 + exp(-(0.1) * (v / mV + 30)))) / ms : Hz

I : amp
g_Na : siemens
g_k : siemens
g_l : siemens
''')

# set initial condition
group = NeuronGroup(1, eqs, method='exponential_euler')
group.v = El
group.g_Na = g_Na
group.g_k = g_k
group.g_l = g_l
group.n = 0.35
group.m = 0.15
group.h = 0.8
group.I = '10*nA'

state = StateMonitor(group, True, record=True)

What you have aready tried

the same model works when I use numpy to simulate (i.e. to solve the differential equation)

Expected output (if relevant)

i expect to have something like: (only care about the first two graph)
(it says new users can only put one embedded media item in a post, but you can find it in the github link below)
which is produced by a numpy based hodgkin huxley model available on GitHub - HeMuling/Hodgkin-Huxley-model: an implementation of Hodgkin-Huxley model using python

Actual output (if relevant)

image
image1008×987 47.7 KB

2

Hi @HeMuling. I think you might simply have specified some parameters incorrectly. In particular, g_Na and g_k should not differ by that much. In the code you linked, your default values are g_Na=120, g_K=36, g_L=0.03, while here they are g_Na=120, g_K=0.3, g_L=0.3.

3

Hi @mstimberg, thanks for correcting my stupid mistake. However, I still have a problem: when I set the current to be 0 at the beginning, a spike was wrongly generated, given that there should not be any spike when there is no current. (graph shown below)

image
image1008×987 86.5 KB

I have no idea how this came, and I am grateful if you could help me

1 Like
4

Did you check that the initial values for m, n and h, as well as the membrane voltage match the equilibrium values?

1 Like
5

@schmitts, thanks for replying. I have found out the problem.

The problem is that the initial membrane potential is too high that it might even reach the threshold (which by rough test, is -60 mV). Therefore, the problem can be solved by setting the initial membrane potential to be -65 mV

here is a brief result:

image
image1008×853 53.8 KB

1 Like