最近在接触脉冲神经网络的基础知识,于是准备先了解几个常见的神经元。首先便是最经典的Hodgkin Huxley模型了。
目录
一、介绍
Hodgkin Huxley模型是Alan Hodgkin和Andrew Huxley在1952年描述的。两人通过测量鱿鱼巨型轴突中动作电位的对不同电流的反应,在没有计算机的情况下拟合出了一个模型。这个模型由一组四个常微分方程组成,虽然描述较为准确,但是计算复杂,不适合大范围仿真。
二、电路图
Hodgkin Huxley模型将细胞建立为如上图所示电路。其中脂质双层被表示为电容;钠、钾离子通道分别用一个电源和电导来表示(图中统一建立为电源和电导);泄露通道也用一个电源和电导来表示;最后产生的离子泵由电流表示;膜电位由电压表示。值得一提的是,其中钠、钾离子通道的电导与时间、电压相关,而泄漏通道的电导与他们不同,是恒定值。
由电路公式我们可以简单得知,流过脂质双层的电流大小为:
又由欧姆定律,对于钠、钾离子通道,其电流都可以表达为:
当然,因为膜内外离子分布不同,细胞外钠离子多,细胞内钾离子多,所以这里的两种电流的方向一定是相反的,二者的数值在符号上一定相反。泄漏电流的公式也与两者一样,所以不再表出。
由此,可以得到流过脂质双层的电流大小表达为:
三、离子电流表征
通过上一小节,我们可以了解到了流过脂质双层的电流大小的表达式。而Hodgkin Huxley通过一系列的实验,开发了一个模型,其中可兴奋细胞的特性由一组四个常微分方程表达。
其中为流经单位面积的电流。在该神经元模型中,每一个离子通道都由多个亚基控制,公式中的n、m、h分别对应钾离子通道亚基激活概率、钠离子通道亚基激活概率、钠离子通道亚基失活概率,因为在鱿鱼轴突中,例如钾离子通道由4个亚基控制,所以公式中要乘上3次方,其他几个乘方也是同理。而代表对应离子通道的最大电导率,乘上对应的概率之后,就得到了对应该时刻的电导率。
在下三个公式中,代表亚基从失活到激活的速率,代表亚基从激活到失活的速率,它们与这一时刻的膜电位相关。乘上对应激活亚基和失活亚基的比例,就可以得到这一时刻的变化量。
通过这四个公式,我们就可以较准确地表达出膜电位的生物特性。而在作者原论文中,给出了和的表达式。
其中,
而关于形如的方程,我们可以得到其解为:
其中,
易知,t=0时,
当然,为了能让这个模型顺利运行起来,我们还需要其他参数。原论文中也给出了作者对这些参数的测量值,如下表所示。分别对应脂质双层的电容(),三种离子通道的静息电位(),三种离子通道的最大电导率(),以及静态膜电位()。
1 | 50 | -77 | -54.4 | 120 | 36 | 0.3 | -65 |
接着我们就可以用代码来仿真这个神经元了。
四、代码实现
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
#参数定义
#平衡电位 单位mV
E_Na = 50.0
E_K = -77.0
E_L = -54.4
V_m = -65.0
#最大电导 单位 mS/cm^2
g_Na = 120.0
g_K = 36.0
g_L = 0.3
#脂质双层电容
C_m = 1
#单位时间
dt = 0.01
#总时间
T = 40
#每一时刻
t = np.arange(0, T, dt)
V = np.zeros(len(t))
n = np.zeros(len(t))
m = np.zeros(len(t))
h = np.zeros(len(t))
#输入电流
I = 0.0
V[0] = V_m
#计算n,m,h初始值
V_tmp = V[0]-V_m
alpha_n = (0.1 - 0.01 * V_tmp) / (np.exp(1 - 0.1 * V_tmp) - 1)
alpha_m = (2.5 - 0.1 * V_tmp) / (np.exp(2.5 - 0.1 * V_tmp) - 1)
alpha_h = 0.07 * np.exp(-V_tmp/20.0)
beta_n = 0.125 * np.exp(-V_tmp / 80.0)
beta_m = 4 * np.exp(-V_tmp / 18.0)
beta_h = 1/(np.exp(3 - 0.1 * V_tmp)+1)
n[0] = alpha_n / (alpha_n + beta_n)
m[0] = alpha_m / (alpha_m + beta_m)
h[0] = alpha_h / (alpha_h + beta_h)
for i in range(1, len(t)):
#计算alpha和beta
V_tmp = V[i-1] - V_m
alpha_n = (0.1 - 0.01 * V_tmp) / (np.exp(1 - 0.1 * V_tmp) - 1)
alpha_m = (2.5 - 0.1 * V_tmp) / (np.exp(2.5 - 0.1 * V_tmp) - 1)
alpha_h = 0.07 * np.exp(-V_tmp/20.0)
beta_n = 0.125 * np.exp(-V_tmp / 80.0)
beta_m = 4 * np.exp(-V_tmp / 18.0)
beta_h = 1/(np.exp(3 - 0.1 * V_tmp)+1)
#计算delta_n, delta_m, delta_h
tau_n = 1.0/(alpha_n + beta_n)
tau_m = 1.0/(alpha_m + beta_m)
tau_h = 1.0/(alpha_h + beta_h)
inf_n = alpha_n *tau_n
inf_m = alpha_m *tau_m
inf_h = alpha_h *tau_h
#更新 n,m,h,使用欧拉一阶近似法
n[i] = n[i-1] + dt * (alpha_n * (1-n[i-1]) - beta_n * n[i-1])
m[i] = m[i-1] + dt * (alpha_m * (1-m[i-1]) - beta_m * m[i-1])
h[i] = h[i-1] + dt * (alpha_h * (1-h[i-1]) - beta_h * h[i-1])
#更新电流
I_Na = g_Na * (m[i] ** 3) * h[i] * (V[i-1] - E_Na)
I_K = g_K * (n[i] ** 4) * (V[i-1] - E_K)
I_L = g_L * (V[i-1] - E_L)
dv = (I - (I_Na + I_K + I_L))*dt/C_m
V[i] = V[i-1] + dv
plt.clf()
plt.subplot(1,1,1)
plt.plot(t,V)
# 膜电位变化
plt.title('Membrane potential')
plt.legend(['V'])
plt.ylim([-100, 100])
接下来我们就可以设置不同大小、时间的电流,来给神经元刺激,观察神经元的膜电位变化情况。
将刺激电流设置为0,得到的图像如下:
可见,膜电位保持在-65mV静息电位。
当给予一个短而大的刺激电流时,膜电位会产生一个很明显的突变,图像如下:
如果持续一段时间地给予较大的电流刺激,膜电位就会在这段时间里持续地产生波形,图像如下:
但当持续一段时间地给予过大的电流刺激时,膜电位在产生一个突变之后便会保持一个比静息电位高的数值不变,直到刺激结束恢复静息电位,图像如下:
当给予神经元一次短而小的刺激时,并不能让膜电位突变到正值;但若短时间间隔内再次用小电流刺激后,膜电位反而会产生一次突变,图像如下:
通过不断改变电流的大小和时间,便可以探究神经元的不同特性了。