首先HH神经元模型的电路结构如下图所示:
LIF神经元模型的电路图如下所示
LIF神经元模型就是对HH神经元模型进行一种简化,把所有电阻都模型化成一个电阻R。
利用基尔霍夫电流定律,流入节点的电流之和等于流出节点的电流之和。可以得到
引入分压原理以及电容和电流之间的关系公式
引入时间常数:τm = RC
在这里就出现了du/dt,电压变化量的求法,用来在代码中通过前一时刻电压值求解后一时刻电压值。
上代码!
def run_LIF(pars, Iinj, stop=False):
"""
Simulate the LIF dynamics with external input current
Args:
pars : parameter dictionary
Iinj : input current [pA]. The injected current here can be a value
or an array
stop : boolean. If True, use a current pulse 用一个电流脉冲
Returns:
rec_v : membrane potential
rec_sp : spike times
"""
# Set parameters
V_th, V_reset = pars['V_th'], pars['V_reset'] #初始化阈值和重置电压
tau_m, g_L = pars['tau_m'], pars['g_L'] #泄露电导,和细胞膜时间常数
V_init, E_L = pars['V_init'], pars['E_L'] #初始化电压,泄露平衡电压
dt, range_t = pars['dt'], pars['range_t'] #时间差
Lt = range_t.size #时间长度
tref = pars['tref'] #不应期
# Initialize voltage and current
v = np.zeros(Lt) #记录迭代时刻的电压
v[0] = V_init #初始化电压
Iinj = Iinj * np.ones(Lt) #输入的信号
if stop: # set end of current to 0 if current pulse
Iinj[:int(len(Iinj) / 2) - 1000] = 0
Iinj[int(len(Iinj) / 2) + 1000:] = 0
tr = 0. # the count for refractory duration 不应期
# Simulate the LIF dynamics
rec_spikes = [] # record spike times 发放脉冲的时间
for it in range(Lt - 1):
if tr > 0: # check for refractoriness 如果tr>0,说明不应期时间没结束
v[it] = V_reset
tr = tr - 1
elif v[it] >= V_th: # reset voltage and record spike event 到达阈值,发放脉冲
rec_spikes.append(it) #放入发放时间
v[it] = V_reset
tr = tref / dt #不应期时间
# calculate the increment of the membrane potential 计算细胞膜电势增量
dv = (-(v[it] - E_L) + Iinj[it] / g_L) * (dt / tau_m)
# update the membrane potential 更新细胞膜电势
v[it + 1] = v[it] + dv
rec_spikes = np.array(rec_spikes) * dt
return v, rec_spikes
结果
这个网站有详解,也有完整代码
https://compneuro.neuromatch.io/tutorials/W2D3_BiologicalNeuronModels/student/W2D3_Tutorial1.html#tutorial-slides