1. 引言

神经科学是研究人类和动物神经系统的结构和功能的科学。神经元放电的不同时间模式是神经科学的基础，即每个神经元的不同放电模式与神经系统信息的编码和处理密切相关 [1]。例如，周期尖峰放电模式和簇放电模式对信息传输有不同的影响 [2] [3]。因此，根据真实的生物物理环境改进神经元模型并研究其放电模式与人类对大脑活动和神经行为的认知密切相关 [4]。

神经科学的开端来自于是Hodgkin和Huxley在1952年发表的文章，他们提出了一个四维(4-D)系统的自治非线性微分方程来描述细胞膜电流和电压对巨型鱿鱼轴突的依赖性 [5]。Morris和Lecar于1981年在Hodgkin-Huxley (H-H)模型的基础上，通过强调膜的导电性和对神经元离子通道类型的深入分析，得到了Morris-Lecar (M-L)神经元模型 [6]。Rinzel和Ermentrou通过改变参数值得到了系统的稳定解 [7]。通过分岔分析，Kunichika阐明了M-L神经元模型放电模式的转变 [8]。史等人提出了一个分数阶M-L神经元模型，其中分数阶导数用于描述膜的记忆效应和幂律 [9]。Izhikevich等人将M-L神经元方程中的外部刺激电流转化为动力学方程控制的变量，以更好地描述外部刺激带来的复杂动力学现象 [10] [11]。Wang通过调整生理参数的值研究了改进后的M-L神经元模型的峰峰间隔(ISI)分岔结构，并从数值模拟结果观察到系统可以从静止状态转变为周期或混沌爆发状态 [12]。Ranjit Kumar和Upadhyay研究了耦合M-L神经元中耦合强度对神经元放电模式的影响 [13]。Song等人基于改进的M-L神经元模型研究了单个神经元的兴奋性自反馈，该模型只需改变一个参数即可显示三种类型的放电模式 [14]。我们发现现有文献很少从ISI的角度探讨神经元放电模式的变化趋势。

目前对神经元模型的研究和完善过程中，对于传统的神经元模型的研究一般考虑离子通道电流、噪声甚至自突触神经元特殊结构对神经元细胞膜电位的影响。实际上除了神经元细胞膜内外离子浓度不同会影响神经元电活动以外，细胞膜内外带电离子在膜上连续交换会引发离子浓度的复杂波动，从而形成时变电磁场。细胞膜内外的电生理环境会因时变磁场带来的电磁感应效应而发生改变，神经元的电活动也会受到影响。因此研究神经元电活动非常有必要考虑神经元中的电磁感应。膜电位与时变磁场之间的相互作用部分需要建模。在M-L模型中，这种相互作用的部分没有描述出来。因此需要一个变量来描述神经元的电磁感应和磁场效应。蔡少棠教授在1971年提出了忆阻器的概念，它是一种具有记忆电阻特性的电路元件，是一种表示磁通量与电荷微分关系的器件 [15]。随后，惠普实验室在2008年宣布成功开发了忆阻器，验证了蔡教授对第四个基本电路元件的预测 [16]。由于忆阻器的记忆特性及其相对于传统电子器件的优势，它可以广泛应用于许多领域，如高集成度非易失性存储器、模拟神经网络、神经生物突触等。吕等人在研究HR神经元模型时，采取对膜电容施加额外的电流来描述磁电流对膜电位的影响，并利用ISI分岔和时序图进行动态分析 [17] [18]。Bao等人在研究HR神经元系统时，用忆阻器来描述电磁感应现象，然后用分岔和Lyapunov指数方法来研究其隐藏的动力学行为 [19]。

数值模拟是研究系统动力学的有力工具 [20] [21]，然而通常忽略边界条件，也受采样时间、精度和数值算法等影响，并且当神经元系统的维度较高时，计算机数值模拟的实时性会变差。因此需要设计出能够满足实时计算，或者再现系统特性的实用型电路。设计实用型电路是模拟神经元特性的理想方式 [22] [23]。Multisim作为模拟电路仿真软件，可以很好的模拟实际电路元件和电路逻辑，在能包含所有的设计工作的同时也可以节省大量的时间成本和资源。本文采用Multisim软件对设计的电路进行验证。例如，Wu等人利用直流电作为神经元的外部刺激电流，提出了一种简化的二维M-L神经元电子学实现方案，然后用Multisim模拟神经元的电活动 [24]。Behdad等人在M-L神经元模型的实现中使用了离散电子元件，但激活时间常数的电路实现是基于泰勒级数展开的近似函数 [25]。胡等人提出并详细研究了基于M-L神经元模型的二维非线性微分方程的模型的完整电子实现，在没有进行简化，并用PSpice个人仿真程序对所设计的电路进行了验证 [26]。本文设计的整个电路由电容器，电阻器，双极结型晶体管，高精度运算放大器和模拟乘法器的市售分立元件实现。该电路可在无需简化即可在神经元模型中模拟神经元的生物学特性，从而能够在不同条件下重现尖峰放电或簇放电，帮助了解神经元的放电开关机制。

2. 带有忆阻器的M-L神经元模型的建立

蔡少棠教授在1971年通过研究电路物理量之间关系时提出了忆阻器。忆阻器的代表了电荷与磁通的关系，作为非线性的记忆元件，该元件的电阻随着电流的增加或减少而变化。当电流降至零时，忆阻器保持最终电阻值，当接收到反向电流时，其最终状态值将被推回。这些特征与生物突触非常相似。

突触作为神经元之间或神经元与效应细胞之间的连接部分，有传递信息的作用。如图1所示，生物突触由突触前膜、突触后膜和突触间隙三部分组成。轴突是神经冲动信号传递到突触小体的通路。信号到达突触小体后，突触前膜的通透性增强，并产生突触小泡。突触小泡与突触前膜融合后破裂，小泡中的神经递质乙酰胆碱被释放到突触间隙，自由游离到达突触后膜，引起突触后膜的变化。将生物系统中突触模型与忆阻器模型进行比较，二氧化钛中的游离氧原子在突触中起着神经递质的作用，具有传输信息和保存记忆的功能。当电流消失时，自由氧原子保存之前的电信号数据，当电流重新出现时，可以恢复到记忆状态。而突触间的神经递质即记忆信号在信号传递后也不会马上消失。需要注意的是突触信息传递是单向的，方向是从突触前膜向突触后膜释放神经递质进行信息传递，信息需要回送时则需要另一个突触进行。而忆阻器中游离的氧原子的移动是双向的，由电流的方向决定。这些神经元细胞之间能够形成了强大的神经网络，而具有信息传输控制和记忆功能的忆阻器便是称为构建智能神经形态网络的重要组成部分。

对于单个神经元模型，无论是确定性的还是随机性的神经元模型都能产生神经元电活动的主要动力学特性，甚至于考虑噪声后，也可以通过调整分岔参数去诱导神经元的模态跃迁等现象。实际神经元中另外一个需要考虑的重要因素是，细胞膜内外带电离子的跨膜运动会造成离子浓度的波动，也就会引起时变电磁场的产生，从而导致离子进行跨膜运动时膜电位会被感应电流调制。因此选择磁通量描述电磁场的作用。在传统的M-L神经元模型研究中没有考虑膜电位与时变电磁场之间的耦合关系，这种关系需要被建模，因此这样对神经元模型进行完善是非常合理的。

为了描述电磁感应的阈值效应，本文首先考虑一种简单阈值磁控忆阻器，如公式1所示 [19]。

$W\left(\phi \right)=\frac{\text{d}q\left(\phi \right)}{\text{d}\phi }=-m\tanh \left(\phi \right),$(1)

其中， $\phi$是通过忆阻器的磁通量，q是通过忆阻器的电荷量。m代表了电磁感应强度。

为了更好地模拟复杂的动力学行为，Izhikevich等人通过将膜电位方程中的外部刺激电流转化为动力学方程控制的变量 [10]，把经典的二维M-L模型改进为三维M-L模型，可以写成

$\begin{array}{l}C\frac{\text{d}V}{\text{d}t}=-u-g_l\left(V-V_l\right)-g_k\omega \left(V-V_k\right)-g_{Ca}{m}_{\infty }\left(V\right)\left(V-V_{ca}\right),\\ \frac{\text{d}\omega }{\text{d}t}=\lambda \left(V\right)\left({\omega }_{\infty }\left(V\right)-\omega \right),\\ \frac{\text{d}u}{\text{d}t}=a\left(0.2+V\right),\end{array}$(2)

其中

$\{\begin{array}{l}{m}_{\infty }\left(V\right)=0.5\left(1+\tanh \frac{V-v_1}{v_2}\right),\hfill \\ {\omega }_{\infty }\left(V\right)=0.5\left(1+\tanh \frac{V-v_3}{v_4}\right),\hfill \\ \lambda \left(V\right)=\frac{1}{3}\cosh \frac{V-v_3}{2v_4}.\hfill \end{array}$(3)

公式中V代表细胞膜电位， $\omega$为具有恢复变量的慢变电流，C代表了细胞膜电容。 $g_l,g_k,g_{ca}$分别表示泄漏电流通道，钾离子通道和钙离子通道的最大电导。 $V_l,V_k,V_{Ca}$分别表示这些通道的平衡电势。 ${\omega }_{\infty }\left(V\right)\text{,}{m}_{\infty }\left(V\right)$分别是钾通道和钙通道打开概率的稳态值，其中 $v_1,v_2,v_3$和 $v_4$是稳态参数。 $\lambda \left(V\right)$是激活时间常数。u表示慢适应外激励电流。a是一个无量纲的参数，代表神经元的时间比例因子，通常情况下 $a\ll 1$，以确保u的慢变性。

基于以上三维M-L模型神经元模型，我们用阈值磁控忆阻器模拟细胞内外离子交换产生的磁场效应，得到带有忆阻器的四维M-L神经元模型

$\begin{array}{l}C\frac{\text{d}V}{\text{d}t}=-u-m\tanh \left(\phi \right)V-g_l\left(V-V_l\right)-g_k\omega \left(V-V_k\right)-g_{Ca}{m}_{\infty }\left(V\right)\left(V-V_{ca}\right),\\ \frac{\text{d}\omega }{\text{d}t}=\lambda \left(V\right)\left({\omega }_{\infty }\left(V\right)-\omega \right),\\ \frac{\text{d}u}{\text{d}t}=a\left(0.2+V\right),\\ \frac{\text{d}\phi }{\text{d}t}=-V,\end{array}$(4)

其中

$\{\begin{array}{l}{m}_{\infty }\left(V\right)=0.5\left(1+\tanh \frac{V-v_1}{v_2}\right),\hfill \\ {\omega }_{\infty }\left(V\right)=0.5\left(1+\tanh \frac{V-v_3}{v_4}\right),\hfill \\ \lambda \left(V\right)=\frac{1}{3}\cosh \frac{V-v_3}{2v_4}.\hfill \end{array}$(5)

这里 $\phi$代表了磁通量，是膜电位V在时间上的积分。 $m\tanh \left(\phi \right)V$表示当磁通量变化时膜电位的反馈电流，m表示电磁感应强度。物理量纲可以按如下得到验证。

$i\text{'}=\frac{\text{d}q\left(\phi \right)}{\text{d}t}=\frac{\text{d}q\left(\phi \right)}{\text{d}\phi }\frac{\text{d}\phi }{\text{d}t}=W\left(\phi \right)v=-m\tanh \left(\phi \right)v,$(6)

其中v是忆阻器两端可变电压，它跟膜电位V是相同的物理单位。 $i\text{'}$表示流经膜阻器的电流，其物理单位与膜电位在时变磁场中的反馈电流相同。

3. 带有忆阻器的神经元模型放电模式趋势分析

3.1. 膜电容和磁感应强度双参数下神经元放电模式变化趋势

放电模式影响着神经元时变磁场，同时磁场也影响着神经元放电模式。本小节将在双参数平面上探索神经元丰富的放电模式过渡趋势和规律，同时研究时变电磁场会给神经元放电模式带来怎样的变化。

图2(a)~(c)分别绘制了双参数平面下四维M-L神经元模型中放电模式的变化规律。图2(a)中，黄色区域代表周期性放电行为，而蓝色区域代表非周期性放电行为。图2(b)展示了单个神经元在1000秒放电序列中的放电周期数，用颜色柱表示着不同的周期数，随着周期数不断增加，颜色由蓝色逐渐过渡到黄色。对比图2(a)和图2(b)可以看出，当周期数足够大，区域颜色为黄色，正好对应非周期放电行为，即混沌放电模式。图2(c)中颜色柱代表了神经元模型ISIs的标准方差。在统计描述中，方差一般被用来计算每个变量与总体平均值之间的差异，在这里用其描述从尖峰放电到簇放电的趋势或每次簇放电模式中的簇中尖峰数量的增加。在簇放电模式的膜电位时序图中可以发现每个簇之间存在明显的静息状态，这也意味着簇放电的ISIs方差比其他放电模式的大。图2(c)中从左上角到右下角存在明显的波浪线，即随着电磁感应强度m减小，膜电容参数C增大，图上呈现明显的过渡边界线。结合图2(d)，在不考虑混沌时，发现随着m的减小，每次分岔，簇放电的尖峰数量增加1。实质上图2(c)中的过渡边界线是分岔曲线。此外，对比图2(a)和图2(c)图，我们发现，神经元混沌放电模式存在于簇放电和尖峰放电交界处。固定膜电容参数 $C=1.3\mu F/c{m}^2$，图2(d)反映了随磁感应强度m变化的ISIs分岔图，在图上可以清楚地观察到明显的不规则放电区域，在 $m\in \left[0.05,0.079\right]$和 $m\in \left[0.159,0.189\right]$区段，与图2(a)和图2(b)中的混沌完全对应。

为了验证上述规律并观察引入忆阻器后神经元放电模式的变化，分别选取磁感应强度 $m=0.028$, $m=0.125$, $m=0.227$和 $m=0.365$，带有忆阻器的四维M-L神经元模型膜电位时序图如图3所示。 $m=0.028$时，神经元呈现周期4簇放电模式， $m=0.125$对应于周期3簇放电模式， $m=0.227$时呈现周期2簇放电模式， $m=0.365$对应周期2尖峰放电模式。所以当固定膜电容参数 $C=1.3\mu F/c{m}^2$，神经元放电模式随磁感应强度m增加整体表现为：从周期性簇放电模式经过倍周期到混沌放电模式再到周期尖峰放电模式，其中簇放电模式中簇中峰的数量程逐渐减少的趋势。

利用上述规律，观察对比图2(b)和图2(c)，可以得出神经元放电模式在膜电容C和磁感应强度m双参数平面下的放电模式变化趋势。

1) 当 $C\in \left[1.1,1.168\right]$时，带有忆阻器的四维M-L神经元系统随着磁感应强度增加全部呈现简单的尖峰放电模式，并且放电节律几乎没有变化。

2) 当 $C\in \left[1.168,1.434\right]$时，该系统随着磁感应强度增加呈现出复杂的放电模式转变，即由周期性簇放电到混沌放电模式再到尖峰放电，神经元放电节律逐渐变慢。

3) 当 $C\in \left[1.434,2\right]$时，该系统随着磁感应强度增加全部呈现簇放电模式并且簇中尖峰数量呈逐渐减少趋势，神经元放电节律逐渐减慢。

Figure 3. Fixed membrane capacitance $C=1.3\mu F/c{m}^2$, the M-L neuron membrane potential timing diagram

图3. 当膜电容 $C=1.3\mu F/c{m}^2$时神经元膜电位时序图

3.2. 磁感应强度和时间尺度因子双参数下神经元放电模式变化趋势

图4(a)~(c)为双参数平面的周期分布、周期数图、ISIs方差分布图。图4(a)中，黄色区域代表周期性放电模式，蓝色区域代表非周期性放电模式。图4(b)中，颜色柱上的数值表示10,000秒仿真模拟中的神经元放电周期数，随着颜色从蓝色过渡到黄色，放电周期的数量增加，与图4(a)对应。图4(c)，

颜色柱上数值代表了神经元ISIs的标准方差。由于无法直观地观察得到图4(b)和图4(c)中的颜色过渡趋势，因此结合a的ISIs分叉图，即图4(d)来探究神经元的放电模式，其中磁感应强度 $m=0.05$。图中可以观察到在神经元混沌放电区域中，有明显的周期性区域，即当 $a\in \left[0.00328,0.00338\right]$和 $a\in \left[0.005156,0.005384\right]$时，与图4(a)和图4(b)对应，混沌区域也能对应。

在图4(d)中，分别选择 $a=0.001$， $a=0.002$， $a=0.004$， $a=0.008$和 $a=0.01$，神经元放电时序图如图5所示。可以清楚看到，当 $a=0.001$时神经元呈现周期1的尖峰放电模式， $a=0.002$时神经元转变为周期2尖峰放电模式， $a=0.004$时开始进入不规则放电模式即混沌放电模式， $a=0.008$时变为周期四簇放电模式。随着 $a$不断增大，簇中尖峰数量呈逐渐减少的趋势， $a=0.01$时神经元呈现周期2簇放电模式。随着时间尺度因子a的增加，神经元放电模式的整个转换趋势是由周期性尖峰放电经过倍周期分岔演变为不规则放电，再经历逆倍周期分岔演变为周期性簇放电，峰值之间的间隔呈下降趋势，这意味着放电节奏变得更快，而随着频率的增加，神经元放电模式达到稳定所需的时间越短。将放电趋势与图4(c)进行对比，发现尖峰放电模式仍然对应于深蓝区域，即方差较小的区域，而簇放电模式的方差仍显著大于尖峰区域。

(c) 混沌放电模式( $a=0.004$)(d) 周期3簇放电模式( $a=0.008$)

(e) 周期2簇放电模式( $a=0.01$)

Figure 5. The neuron membrane potential timing diagram with fixed magnetic induction intensity $m=0.05$

图5. 神经元膜电位时序图，其中磁感应强度 $m=0.05$

图6. 神经元膜电位时序图，其中时间尺度因子

为了探究忆阻器对M-L神经元系统放电模式的影响，我们固定时间尺度因子，神经元系统随磁感应强度m的变化的ISIs分岔图如图6(a)所示。选择其中几个点画出神经元膜电位时序图，如图6(b)~(e)所示。当时，神经元呈现周期3簇放电模式，在转变为不规则放电模式即混沌放电模式，对应周期2簇放电模式，当时呈现周期2尖峰放电模式。随着磁感应强度m的变大，神经元整体放电模式趋势为从周期性簇放电经历逆倍周期到不规则放电模式再到周期性尖峰放电，其中放电模式尖峰之间的间隔呈变大的趋势，这意味着放电节奏变慢。这变化趋势遵循之前发现的规律，即深蓝色区域对应的区域均为尖峰放电模式，其余区域对应簇放电和不规则放电模式。此外观察图6(b)，发现其中存在一些分割线，分割线两端周期数成倍数关系，与周围区域相比，分割线上的点所对应的放电模式周期数明显增加，分割线上的点所在区域对应到图6(b) ISIs分岔图上为发生倍周期分岔的点。

4. 带有忆阻器的M-L神经元的系统实现

4.1. 双曲正切函数关系电路实现

观察神经元系统数学模型，即公式4和5，可以发现实现，和忆阻器电路部分的关键就是设计一个输入输出是双曲正切函数关系的电路。差分放大器电路由于其传输特性而可以解决该问题。图7为恒流源差分放大器电路的示意图。

图中各部分的电压关系为

(7)

这里和为恒流源差分放大器电路的输入电压，为双端输出的差值，分别是双极结型晶体管和的基极–发射极电压。

晶体管发射极电流和BE结电压之间的关系近似为

(8)

这里表示饱和电流，是基极–发射极电压，是晶体管热电压，在常温下热电压的数值约为26 mV。图中各部分电流关系为

(9)

这里I是流过电流源的电流，和分别是双极结型晶体管的集电极电流。根据公式7，8和9，进一步得到

(10)

(11)

则两端输出电压为

(12)

基于上述方程，根据系统参数确定恒流源、电阻和双端输入的值，可以实现，和忆阻器电路部分。

4.2. 双曲余弦函数关系电路的实现

双曲余弦函数的表达式如下

(13)

在电路中双曲余弦函数很难直接得到，根据双曲余弦函数与指数函数的关系，双曲余弦函数可以描述为输入值相同但符号相反的两个指数函数之和的一半。因此，问题的关键转化为设计一个可以得到指数函数的电路，电路需要满足输入可以是正值也可以是负值的要求。指数函数的电路设计如图8所示。

图中节点电压和电流之间的关系为

(14)

(15)

这里和分别是指数函数电路的输入电压和输出电压。分别是双极结型晶体管和的基极–发射极电压和集电极电流。电路增益电压，整个指数函数电路单元的输入电压和输出电压之间的关系可以描述为

(16)

根据指数函数和双曲余弦函数之间的关系，双曲余弦函数的电路单元如图9所示。

通过对图9中的电路使用标准节点分析技术，可以将输入电压与输出电压的关系描述为

(17)

在电路图9中，和为平衡电阻。当，时，上式可写为

(18)

基于上述方程，基本可以实现双曲余弦函数的电路单元。

4.3. 带有忆阻器的M-L神经元模型的电路实现

基于上述电路单元和四维M-L神经元的数学模型，我们设计了该模型的完整电路实现，如图10所示。忆阻器部分的电路在方框A中实现。和电路单元分别在方框B中和方框C中实现，方框D中实现了的电路功能。恒定电压源分别是钾离子通道，钙离子通道和泄漏电流通道的平衡电位。电容器是神经元的膜电容。通过对图10中的电路使用标准节点分析技术，控制电路行为的状态方程可以描述为

(19)

(20)

(21)

其中是方框B的输出电压，是方框C的输出电压，是方框D的输出电压，通过与公式4对比，可以确定电路中电器元件对应的其余变量

，，，，，，，i对应于方框A忆阻器输出电流。

电容器代表细胞膜电容参数C，当，同步到系统参数即，。当膜电容参数C发生改变时，神经元展示不同的放电模式。当和，神经元电路输出如图11(a)和图11(c)所示，其实验结果与数值仿真结果图11(b)和图11(d)保持一致。

(c) 当时电路实验结果 (d) 当时数值仿真结果

Figure 11. Comparison diagram between experimental results and numerical simulation

图11. 实验结果与数值模拟对比图

5. 总结

离子交换在神经元细胞膜上产生的时变磁场会影响神经元的放电模式。忆阻器是描述电荷和磁通量之间关系的非线性元件。与生物突触相比，二氧化钛中的自由氧原子具有传递信息和保存记忆的功能，就像突触中的神经递质一样。此外，忆阻器具有体积小、能耗低等优点，在理论和工程上都非常适合描述神经元细胞膜离子交换产生的电磁感应和磁场效应。本文提出了一种新的具有阈值电磁感应的四维M-L神经元模型。证实了不同强度的电磁场能引起神经元电活动的不同模式转变。探讨了电磁感应强度、膜电容、标度因子等参数和参数平面下神经元放电模式的变化趋势。利用多参数平面上的ISIs分岔、ISIs放电周期、ISIs方差等方法研究神经元放电模式的变化趋势。最后，基于神经元模型的四维非线性微分方程，设计了该模型的完整电子实现，特别是双曲正切函数和双曲余弦函数的电路单元。该电路模拟了神经元模型中神经元的生物学特性，在没有进行简化条件，从而能够在不同条件下再现尖峰或爆发，从而对神经元的放电转换机制有了清晰的认识。同时，本文提出的神经元电路可以作为研究单个神经元动力学的实验单元。

参考文献