电场作用下HR神经元的分岔分析及参数辨识

doi:10.3785/j.issn.1008-9497.2022.06.007

电场作用下HR神经元的分岔分析及参数辨识

肖冉^,^,, 安新磊^,^,, 祁慧敏, 乔帅

兰州交通大学数理学院，甘肃兰州 730070

Bifurcation analysis and parameter identification of HR neurons under electric field

XIAO Ran^,^,, AN Xinlei^,^,, QI Huimin, QIAO Shuai

School of Mathematics and Physics，Lanzhou Jiaotong University，Lanzhou 730070，China

通讯作者: ORCID： https：//orcid.org/0000-0003-2424-4870，E-mail：anxin1983@163.com.

收稿日期: 2021-08-16

基金资助:

国家自然科学基金资助项目. 11962012

Received: 2021-08-16

作者简介 About authors

肖冉（1996—），ORCID：https：//orcid.org/0000-0002-2468-4895，男，硕士研究生，主要从事非线性动力学研究，E-mail：455235427@qq.com. , E-mail：455235427@qq.com

摘要

详细分析了电场作用下四维Hindmarsh-Rose（HR）神经元模型的分岔模式及放电行为。通过数值仿真得到该神经元模型的多组双参数分岔图、最大Lyapunov指数图、峰峰间期分岔图等，发现该模型在双参数平面上存在倍周期分岔、加周期分岔等模式及“锯齿状”混沌结构。通过构建合适的目标函数，提出了自适应混合粒子群遗传算法，将神经元模型的参数辨识转化为最优化问题。数值仿真结果表明，算法对神经元模型的参数辨识效果较好，能更准确地辨识未知参数，具有一定优越性。

关键词： HR神经元 ; 双参数分岔 ; 混合粒子群遗传算法 ; 参数辨识

Abstract

The bifurcation mode and firing behavior of the four-dimensional Hindmarsh-Rose (HR) neuron model under electric field are analyzed in detail. Several groups of bifurcation diagrams of the neural system with two parameters are derived by numerical simulation, which correspond to the maximum Lyapunov exponent diagram and the ISI bifurcation diagram. It is found that the system has period-doubling bifurcation, period-adding bifurcation and chaotic structure with "jaggedness" on the biparametric plane. Since parameter identification of neuron model is an important part of neuron dynamics analysis, by constructing a suitable objective function, we propose an adaptive hybrid particle swarm genetic algorithm to convert the parameter identification problem of neuron model into an optimization problem. Numerical simulation results show that the proposed algorithm is effective and feasible in the parameter identification of neuron model.

Keywords： Hindmarsh-Rose (HR) neuron ; bifurcation with two parameters ; particle swarm optimization-genetic algorithm ; parameter identification

PDF (3713KB) 元数据多维度评价相关文章导出 EndNote| Ris| Bibtex 收藏本文

本文引用格式

肖冉, 安新磊, 祁慧敏, 乔帅. 电场作用下HR神经元的分岔分析及参数辨识. 浙江大学学报(理学版)[J], 2022, 49(6): 691-697 doi:10.3785/j.issn.1008-9497.2022.06.007

XIAO Ran, AN Xinlei, QI Huimin, QIAO Shuai. Bifurcation analysis and parameter identification of HR neurons under electric field. Journal of Zhejiang University(Science Edition)[J], 2022, 49(6): 691-697 doi:10.3785/j.issn.1008-9497.2022.06.007

生物神经元是机体产生认知功能的系统，在生理活动调节中起主导作用，是最复杂的非线性系统之一。神经元系统的基本结构和功能单位是神经元，神经元频繁且复杂的放电过程具有非线性，这对进一步研究其动力学行为有重要意义^［1］。对生物神经元放电模型的研究可追溯至20世纪50年代，HODGKIN等^［2］通过乌贼轴突触电压钳实验，建立了Hodgkin-Huxley（HH）神经元模型，这是首个关于神经元放电的数学模型。20世纪80年代，HINDMARSH等^［3-4］以丘脑神经元为研究对象，建立了Hindmarsh-Rose（HR）神经元模型。

神经元的放电活动受离子转运、电磁辐射等外界因素的影响，已成为近年来的研究热点之一。考虑外部电磁辐射对神经元放电活动的影响不可忽视，文献［5］采用忆阻器实现了感应电流与膜电位的耦合。文献［6］在考虑神经元膜电位与外部磁场之间的电磁感应效应的基础上，建立了改进的四维HR神经元模型。文献［7］将磁通量引入FitzHugh-Nagumo（FHN）神经元模型用以描述电磁感应效应，发现电磁辐射可使神经元的放电方式和螺旋波发生改变。文献［8］通过引入电场变量，建立了外电场作用下的FHN神经元模型，并通过数值仿真分析了该模型的动力学行为。结果表明，电场变量E可以很好地描述离子的分布和膜电位变化。文献［9］进一步研究了FHN神经元模型的动力学行为，同时证实了此模型中无共存吸引子。文献［10］利用电场效应下改进的HR神经元模型，研究了当跨膜电位与电场进行突触耦合时两个神经元之间的有限时间同步，结果表明，电场耦合强度的变化对神经元间的同步影响较小，改变突触耦合强度则对神经元间的同步影响较大。说明电场可以成为促进信号交换和编码的有效桥接。文献［11］通过引入磁通量，建立了改进的HH神经元模型，探讨了该模型的分岔模式以及在电磁耦合强度下的同步行为。文献［12］研究了HH神经元模型受电磁辐射影响时放电活动的变化，并证明了电磁辐射的影响可以从单个神经元传递至附近的其他神经元。文献［13］采用数值仿真方法研究了HR神经元模型在电磁辐射影响下的放电行为，并通过构造Washout控制器消除系统的隐藏放电行为。文献［14］分析了电磁感应作用下HR神经元模型的混合模式振荡放电特征，同时通过构造Hamilton能量反馈控制器，将该模型从混沌振荡态转为不同周期的簇放电态。文献［15］用相对阈值耦合方法对混沌神经元实施控制，分析了不同弛豫时间以及阈值对混沌神经元动力学行为的影响。

神经元系统参数辨识是进行神经元动力学分析的重要环节。在实际中，神经元系统不可避免地受周期性扰动或外界因素的影响，引起系统多个参数同时变化，因此辨识变化后的参数准确值十分关键。有效的参数辨识可将实验获得的有用数据与理论模型结合，建立更符合实际的神经元模型。

神经元模型的参数辨识方法主要有2种，一种是基于Laypunov稳定性理论以及自适应同步的方法，以原系统为驱动系统，构建对应的响应系统，设计适当的控制器，实现驱动系统与响应系统同步，辨识未知参数^［16-17］。但某些神经元模型的控制器难以构建，且构建控制器需要原系统多个状态变量的时间序列，所需信息量较大，操作困难。另一种是运用现代群智能优化算法的寻优能力，通过构建合适的目标函数，将参数辨识问题转化为最优化问题进行神经元模型参数辨识。该方法在构建目标函数时只需1个状态变量的时间序列，操作容易，是当前主流的一种参数辨识方法。文献［18］设计了2种目标函数，通过实数编码的遗传算法辨识HR神经元模型的未知参数；文献［19］采用改进的遗传算法参数辨识方法，同时辨识了内埋式永磁同步机的4个参数；文献［20］将混沌算法与粒子群算法融合，用自适应混沌粒子群算法辨识了光伏电池模型参数；文献［21］提出了数据驱动的改进型Newton-Raphson法，辨识的电网线路参数准确性较高。

受上述研究启发，考虑细胞膜离子交换时产生的电场可能对神经元系统有影响，本文引入电场变量，构建了一个四维HR神经元模型。通过MATLAB软件得到该模型的双参数分岔图、最大Lyapunov指数图及峰峰间期分岔图等，分析了模型的放电模式及分岔现象。考虑放电过程中模型的多个参数会发生变化，采用混合粒子群遗传算法对四维HR神经元模型进行了参数辨识，并通过数值仿真验证了算法的有效性和可行性。

1　模型描述

基于三维HR神经元模型，由于细胞内外连续离子交换会产生时变电场，时变电场进一步影响细胞内外离子交换，引起膜电位发生变化，故考虑时变电场的影响，建立改进的四维HR神经元模型^［22］，微分方程为

\{\begin{array}{l} \dot{x} = y - a x^{3} + b x^{2} - z + I, \\ \dot{y} = c - d x^{2} - y + k_{0} E, \\ \dot{z} = r [s (x - χ_{0}) - z], \\ \dot{E} = k_{1} y - k_{2} E + E_{e x t}, \end{array}

（1）

其中， $\dot{x} ， \dot{y} ， \dot{z} ， \dot{E}$ 均为状态变量，分别表示膜电压、钾离子和钠离子相关快电流的恢复变量、钙离子相关慢电流的自适应电流以及电场变量； $I$ 为外界刺激电流； $a, b, c, d, r, s, χ_{0}$ 为系统重要动力学参数； $y$ 为快电流，其对感应电场的变化十分敏感； $k_{1} y$ 为快电流与外界电场的变化速率， $k_{1}$ 为变化系数； $k_{0} E$ 表示时变电场的影响；反馈变量 $- k_{2} E$ 表示神经元系统对时变电场的自适应调节所产生的漏电场； $E_{e x t}$ 表示由外电场产生的周期性扰动或随机噪声的影响，取 $E_{e x t} = A s i n (2 π f t)$ ， $A$ 和 $f$ 分别为周期扰动的振幅和频率。在数值仿真中，各参数基准值分别取 $a = 1.0$ ， $b = 3.0$ ， $c = 1.0$ ， $d = 5.0$ ， $r = 0.006$ ， $χ_{0} = - 1.61$ ， $k_{0} = 0.1$ ， $k_{1} = 0.2$ ， $k_{2} = 0.5$ ， $I = 3$ ， $A = 0$ ， $f = 0.01$ 。

2　双参数分岔分析

通过MATLAB软件得到不同双参数组合下系统放电行为的双参数分岔图及其对应的最大Lyapunov指数图（图1）。

图1

新窗口打开| 下载原图ZIP| 生成PPT

图1 双参数分岔图及最大Lyapunov指数图

Fig.1 Bifurcation and maximum Lyapunov exponent of double parameters

以不同的双参数组合为控制变量、其他参数为基准值，系统在双参数平面中的分岔行为如图1（a）（c）（e）所示，其中不同的颜色表示不同的放电态，对应的数字表示放电周期，如1表示尖峰放电态，2表示周期2簇放电态，4表示周期4簇放电态，…，周期为20及20以上的系统放电态和混沌态均用白色表示，对应的最大Lyapunov指数图如图1（b）（d）（f）所示，其中，蓝色表示静息态，黄色表示周期簇放电态，橙红色表示混沌态。

由图1（a）知，系统出现多种周期簇放电态且伴有“锯齿状”的混沌态。根据图1（a）中黑线自左上方至右下方的走势，系统从混沌放电态进入周期3簇放电态，随后经历倍周期分岔进入周期6，12，…簇放电态，周期簇放电态结束后进入混沌放电态，接着出现周期4窗口，再次经历倍周期分岔进入周期8，16，…簇放电态和混沌放电态……。混沌放电态结束后系统进入周期数+1的簇放电态，即伴随混沌的加周期分岔模式。此外，由 $I \in [2.70,2.95]$ ， $r \in [0.002 0,0.006 9]$ 可知，随着 $r$ 的增大，系统逐渐进入周期数+1的簇放电态，且不伴有混沌区域，即无混沌的加周期分岔模式。图1（b）中的橙红色混沌结构与图1（a）的白色区域相对应，说明白色区域是混沌放电态窗口。

在图1（c）中，底部呈现大面积淡蓝色区域，表明系统在此范围内呈尖峰放电态，2个参数的变化对系统的分岔行为几乎无影响。当 $c \in [0.4,2.0]$ ， $d \in [4.4,5.6]$ 时，随着2个参数的同步增大，系统先从周期2簇放电态通过倍周期分岔进入周期4，8，…簇放电态，再进入混沌放电态，接着出现周期3窗口，通过倍周期分岔进入混沌态……这也是伴随混沌的加周期分岔过程。周期≥20的白色窗口呈“锯齿状”，结合图1（d），可知橙红色混沌态结构与图1（c）的白色“锯齿状”结构相吻合。

由图1（e）可知，当 $a \in [0.558,0.695]$ ， $r \in [0.002,0.022]$ 时，系统始终处于大面积尖峰放电态或周期2簇放电态，说明此时2个参数对系统的分岔行为影响较小。当 $a \in [0.035,0.558]$ ， $a \in [0.695,1.150]$ ， $r \in [0.002,0.022]$ 时，随着 $r$ 的减小，左右两侧均伴随与图1（a）相似的分岔行为，即伴随混沌的加周期分岔模式。在此基础上，由图1（f）知，左右两侧均有“锯齿状”混沌结构，其与图1（e）左右两侧的白色区域对应。

按图1（a）中黑线自右下方至左上方分别作系统膜电压的峰峰间期（ISI）分岔图及其对应的最大Lyapunov指数图（图2）。由图2（a）知，随着 $I$ 的增大和 $r$ 的减小，系统膜电压的放电模式依次为混沌放电态→周期3簇放电态→周期6簇放电态→周期12簇放电态→…→混沌放电态→周期4簇放电态→周期8簇放电态→周期16簇放电态→…→混沌放电态……。由图2（b）知，系统反复进入周期簇放电态与混沌态，且混沌态的区间与峰峰间期分岔图相吻合。

图2

新窗口打开| 下载原图ZIP| 生成PPT

图2 系统（1）关于参数 $I$ 和 $r$ 的峰峰间期分岔图及其最大Lyapunov指数图

Fig.2 The ISI bifurcation diagram and the maximum Lyapunov exponent diagram of parameters $I$ and $r$ in system （1）

综上，神经元系统的放电态发生变化是因为受外界刺激参数发生了变化。当系统的多个参数同时发生变化时，通过已有信息对变化后的参数进行辨识具有重要意义。

3　四维HR神经元模型的参数辨识

3.1　参数辨识方法介绍

在智能优化算法中，粒子群算法和遗传算法较为经典，粒子群算法是一种随机搜索算法，将每个解看作一个粒子，每个粒子为使个体距离目标位置最短，需同时靠近个体历史最优位置及种群最优位置。遗传算法采用“物竞天择，适者生存”理念，先选择保留优秀个体，再通过仿遗传过程的交叉和变异操作对解进行搜索，是一种启发式算法。

这2种算法的优缺点较明显。粒子群算法中各粒子迅速向历史个体最优位置和种群最优位置靠近，种群快速聚集，迭代前期效率非常高，能快速逼近最优解，但较易陷入局部最优，出现早熟现象，且难以跳出局部最优。遗传算法的全局性较好，但整体效率不高。本文将2种算法相结合，在粒子群迭代中引入遗传算法的交叉和变异操作，通过交叉和变异操作跳出局部最优，并用交叉公式计算得到新的染色体，从而增加种群多样性，找到全局最优解。

3.2　自适应混合粒子群遗传算法

在粒子群算法中，惯性权重 $ω$ 以及学习因子 $C_{1}$ ， $C_{2}$ 均能在很大程度上影响算法性能。在迭代前期，较大的惯性权重有利于算法跳出局部最优，加强全局寻优性能；在迭代后期，较小的惯性权重有利于算法局部寻优并达到收敛。本文采用自适应惯性权重：

ω = ω_{s} - \frac{(ω_{s} - ω_{e}) g}{G}

，（2）

其中， $ω_{s}$ ， $ω_{e}$ 分别为初始惯性权重与迭代结束时的惯性权重； $g$ 为当前迭代次数； $G$ 为迭代总次数。

学习因子 $C_{1}$ 代表粒子向个体历史最优学习的能力，学习因子 $C_{2}$ 代表粒子向全种群最优学习的能力。在迭代前期，较大的 $C_{1}$ 和较小的 $C_{2}$ 可抑制各粒子快速聚集，使算法具有更好的全局寻优性能；在迭代后期，较小的 $C_{1}$ 和较大的 $C_{2}$ 有利于算法收敛。本文采用自适应学习因子：

\{\begin{array}{l} C_{1} = C_{1 s} + \frac{(C_{1 e} - C_{1 s}) g^{2}}{G^{2}} ， \\ C_{2} = C_{2 s} + \frac{(C_{2 e} - C_{2 s}) g^{2}}{G^{2}} ， \end{array}

（3）

其中， $C_{1 s}$ ， $C_{1 e}$ 分别为 $C_{1}$ 的初始值和最终值； $C_{2 s}$ ， $C_{2 e}$ 分别为 $C_{2}$ 的初始值和最终值。

在更新粒子的速度与位置后，为增加种群多样性，跳出局部最优，引入遗传算法的交叉和变异操作。考虑遗传算法的编码方式，将优化变量作为染色体进行二进制编码，再进行交叉和变异操作。当遇到高维度优化问题或有高寻优精度要求时，二进制编码较冗余，搜索空间过大，难度增加，影响迭代速度。因为待辨识参数均为实数，对辨识精度要求高，且参数辨识属于高维连续优化问题，所以采用实数编码的方式更合适，获得的辨识结果精度更高。

在一般的实数编码遗传算法中，交叉操作是交换2个个体对应的实数染色体，并未产生新的染色体，得到优质个体，因此，多样性并不丰富。本文所用交叉操作公式为

\{\begin{array}{l} n e w_{1} = l \times o l d_{g o o d} + (1 - l) \times o l d_{b a d} ， \\ n e w_{2} = (1 - l) \times o l d_{g o o d} + \frac{l}{d p} \times o l d_{b a d} ， \end{array}

（4）

其中， $l \in [0,1]$ 为随机数，参数 $d p$ 可以调节子代产生的方向，通常， $d p > 1$ 可使子代往好的父代方向移动，而非盲目随机产生个体。一方面，可以得到大量新的染色体，增加种群多样性；另一方面，2个子代的和与2个父代的和较接近，能增加子代的收敛性。变异操作根据变异概率（取0.1）对挑选的个体边界重新赋值。

自适应混合粒子群遗传算法步骤：

步骤1 将种群的位置与速度初始化，将个体最优位置和全局最优粒子的位置初始化。

步骤2 根据式（2）和式（3）更新惯性权重和学习因子。

步骤3 更新整个种群的位置和速度，并对跑出边界的粒子做边界化处理，防止其越界。

步骤4 用式（4）对每对父代粒子做交叉操作，生成新种群。

步骤5 根据变异概率对种群进行变异操作。

步骤6 更新个体最优位置和全局粒子最优位置及最优值。

步骤7 判断是否达到最大迭代次数，若是，则转步骤8结束循环；若否，则转步骤2继续循环。

步骤8 输出全局粒子最优位置和最优值。

4　数值仿真

为体现本文算法在神经元参数辨识中的高精度和优越性，分别采用标准粒子群算法、实值编码遗传算法以及本文算法通过MATLAB软件对四维HR神经元模型进行参数辨识。将 $a$ ， $b$ ， $d$ ， $r$ 作为待辨识参数。目标函数设定为

m i n f = {‖p_{n} (x (t)) - p_{n} (\hat{x} (t))‖}_{2}

，

其中， $p_{n} (x (t))$ 为真实参数状态变量 $x$ 的前200个时间向量， $p_{n} (\hat{x} (t))$ 为待辨识参数状态变量 $x$ 的前200个时间向量。目标函数越小，待辨识参数系统状态变量的信号与实际参数系统状态变量时间序列的信号相差越小，待辨识参数与真实参数越接近，从而将神经元模型的参数辨识问题转化为最优化问题。

4个待辨识参数的取值范围分别为 $a \in [0.2,1.8]$ ， $b \in [1.8,4.0]$ ， $d \in [4.2,6.8]$ ， $r \in [0.001,0.012]$ 。参数真实值取模型描述中的基准值。设置种群数为30，最大迭代次数为150，取 $C_{1 s} = 2.0$ ， $C_{1 e} = 0.5$ ， $C_{2 s} = 0.5$ ， $C_{2 e} = 2.0$ ， $ω_{s} = 0.8$ ， $ω_{e} = 0.4$ ， $d p = 5$ 。3种算法参数设置相同。

3种算法的目标函数曲线如图3所示。可知，本文算法在迭代前期的辨识效果较好，迭代速度快，且前期辨识精度已高于另2种算法的最终辨识精度。当迭代次数为150时，本文算法辨识精度较另2种算法高。

图3

新窗口打开| 下载原图ZIP| 生成PPT

图3 各算法的目标函数曲线

Fig.3 Objective function curve of different algorithms

3种算法的参数辨识结果如表1和图4所示。

表1 3种算法的参数辨识结果

Table 1 Parameter identification results of three algorithms

参数	遗传算法	粒子群算法	本文算法	真实值
a	0.904 4	0.952 0	0.998 7	1
b	2.806 5	2.851 9	2.997 0	3
d	4.745 0	5.173 5	4.999 7	5
r	0.008 4	0.005 5	0.006 0	0.006

新窗口打开| 下载CSV

图4

新窗口打开| 下载原图ZIP| 生成PPT

图4 3种算法的参数辨识曲线

Fig.4 Parameter identification curve of three algorithms

由表1和图4可知，相比于粒子群算法和遗传算法，本文算法参数辨识结果与真实值最接近，误差最小，尤其对最小参数 $r$ 的识别非常精准，这是因为本文算法将遗传算法的交叉和变异操作引入粒子群迭代，跳出了局部最优，找到了更接近真实值的结果。因此，本文算法能更准确地辨识未知参数，解决神经元模型的参数辨识问题，具有优越性。

5　结论

首先，在三维HR神经元模型的基础上引入了电场变量，构建了电场作用下四维HR神经元模型，并用MATLAB软件作了双参数分岔图及其对应的最大Lyapunov指数图、峰峰间期分岔图及其对应的最大Lyapunov指数图，发现模型具有丰富的动力学性质及放电模式。由于神经元模型维数较高，对参数辨识精度要求高，一般的粒子群算法和遗传算法辨识效果欠理想，为此提出了一种自适应混合粒子群遗传算法，并将其应用于神经元模型的参数辨识。数值仿真结果表明，该算法的寻优精度高，迭代速度快，辨识效果好。研究结果对电场下生物学及医学神经元模型的建立以及参数辨识具有一定的参考价值。

http://dx.doi.org/10.3785/j.issn.1008-9497.2022.06.001

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

古华光

神经系统信息处理和异常功能的复杂动力学

［J］. 力学学报， 2017， 49（2）： 410-420. DOI：10. 6052/0459-1879-16-315