CN113792863A - 一种脉冲神经网络建模方法、系统及其应用 - Google Patents

Abstract

本发明涉及一种脉冲神经网络建模方法、系统及其应用，先构建指数型积分放电脉冲神经元模型，然后基于待模拟真实神经网络的空间信息对所有神经元进行选择性连接，得到对应于待模拟真实神经网络的脉冲神经网络，空间信息包括网络结构约束信息、神经元位置信息以及子模块位置信息，并基于固定时间步长对该脉冲神经网络的状态进行计算，进而能更准确的模拟不同刺激下真实神经元的反应，并且能更为精准的建模真实神经网络的连通关系，更为准确的建立用于模拟生物活体神经回路和体外培养神经网络的脉冲神经网络，从而为研究隐藏在生物记忆和学习行为背后的神经信号处理和传输奠定基础。

Description

一种脉冲神经网络建模方法、系统及其应用

技术领域

本发明涉及计算神经科学与人工智能交叉技术领域，特别是涉及一种基于空间信息的脉冲神经网络建模方法、系统及其应用。

背景技术

近年来，脉冲神经网络因其相对于其他人工神经网络更仿生，网络结构更紧凑，训练能耗更小等特点，在人工智能领域得到了广泛研究。最新研究表明，在自动驾驶领域，仅包含数十个人工神经元的脉冲神经网络在决策功能上能达到与传统深度卷积神经网络相近的水平。不仅在人工智能领域展现出巨大的发展前景，脉冲神经网络更为传统和基础的应用—模拟真实的生物神经网络的功能性，也在脑科学和计算神经科学领域取得了长足的进展。目前，人工脉冲神经元已能准确地预测和呈现出接受外界电刺激的单个生物神经元的反应，包括实时的膜电位和脉冲发放时间。在能准确模拟单个生物神经元状态的基础上，脉冲神经网络建模成为了研究生物活体神经回路或体外培养神经网络和解释甚至破解其编解码方式的主流手段。

根据真实测量的生物活体神经回路或体外培养神经网络的重要参数，如神经元间的连通概率，平均静息膜电位，脉冲发放的阈值，突触后电位的强度等，神经科学家们能建模出对应的脉冲神经网络，并研究其自发电活动、信息处理、信号传导以及记忆和学习等功能。而这一切的基础，建立在对生物活体神经回路或体外培养神经网络中神经元间的功能性连接关系的准确再现上。由于难以准确、无伤且高通量地确定参与信息处理的每一神经元与网络中其他神经元的连通关系，目前的主流方法是随机选择数对神经元，确定其连通关系，进而基于统计手段估计整个网络中的平均连通概率，从而基于这一参数生成随机连接的，对应真实神经网络的脉冲神经网络。然而，由于生物活体神经回路和部分体外培养神经网络呈现出模块化和分布式的特点，而非均匀分布的整体神经网络，模块内外连通概率的差异、结构约束对连通概率的影响等均需要得到考虑，导致上述方法所建立得到的脉冲神经网络模拟精度低。

发明内容

本发明的目的是提供一种脉冲神经网络建模方法、系统及其应用，基于空间信息来更为准确的建立用于模拟生物活体神经回路和体外培养神经网络的脉冲神经网络，更为精准的建立网络动态特性与单个神经元的实时状态之间的联系，从而为研究隐藏在生物记忆和学习行为背后的神经信号处理和传输奠定基础。

为实现上述目的，本发明提供了如下方案：

一种脉冲神经网络建模方法，所述建模方法包括：

构建指数型积分放电脉冲神经元模型；所述指数型积分放电脉冲神经元模型用于基于接收到的所有输入电流改变膜电位；

基于待模拟真实神经网络的空间信息对所有神经元进行选择性连接，得到对应于所述待模拟真实神经网络的脉冲神经网络；所述空间信息包括网络结构约束信息、神经元位置信息以及子模块位置信息；每一所述神经元均采用所述指数型积分放电脉冲神经元模型。

一种脉冲神经网络建模系统，所述建模系统包括：

构建模块，用于构建指数型积分放电脉冲神经元模型；所述指数型积分放电脉冲神经元模型用于基于接收到的所有输入电流改变膜电位；

连接模块，用于基于待模拟真实神经网络的空间信息对所有神经元进行选择性连接，得到对应于所述待模拟真实神经网络的脉冲神经网络；所述空间信息包括网络结构约束信息、神经元位置信息以及子模块位置信息；每一所述神经元均采用所述指数型积分放电脉冲神经元模型。

一种脉冲神经网络的网络状态计算方法，基于固定时间步长对利用上述建模方法所得到的脉冲神经网络的网络状态进行计算，所述计算方法包括：

采用向前欧拉法，逐时间步长、逐个神经元计算所有神经元的膜电位和钙离子浓度；

根据研究需要，记录每一所述时间步长每一所述神经元的所述膜电位和所述钙离子浓度，并记录脉冲神经网络中出现动作电位的神经元序号以及出现动作电位的时间。

根据本发明提供的具体实施例，本发明公开了以下技术效果：

本发明用于提供一种脉冲神经网络建模方法、系统及其应用，先构建指数型积分放电脉冲神经元模型，然后基于待模拟真实神经网络的空间信息对所有神经元进行选择性连接，得到对应于待模拟真实神经网络的脉冲神经网络，空间信息包括网络结构约束信息、神经元位置信息以及子模块位置信息，并基于固定时间步长对该脉冲神经网络的状态进行计算。进而本发明基于指数型脉冲放电脉冲神经元模型，具有基于输入电流自适应放电阈值的特点，能更准确的模拟不同刺激下真实神经元的反应。将待模拟真实神经网络中的空间信息纳入考虑范畴，能更为精准地建模真实神经网络的连通关系，进而能更为准确的建立用于模拟生物活体神经回路和体外培养神经网络的脉冲神经网络，更为精准的建立网络动态特性与单个神经元的实时状态之间的联系，从而为研究隐藏在生物记忆和学习行为背后的神经信号处理和传输奠定基础。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明实施例1所提供的建模方法的方法流程图；

图2为本发明实施例1所提供的连接方法的方法流程图；

图3为本发明实施例2所提供的建模系统的系统框图；

图4为本发明实施例3所提供的计算方法的方法流程图；

图5为本发明实施例3所提供的整体流程图；

图6为本发明实施例3所提供的单向连通的2×1模块化神经网络的连通示意图；

图7为本发明实施例3所提供的神经元的放电事件对应的动作电位栅格图；

图8为本发明实施例3所提供的仿真过程中一次网络局部簇放电事件对应的动作电位栅格图；

图9为本发明实施例3所提供的仿真过程中一次网络整体簇放电事件对应的动作电位栅格图；

图10为本发明实施例3所提供的编号为1的神经元的膜电位随时间变化的曲线；

图11为本发明实施例3所提供的编号为1的神经元的胞内钙离子浓度随时间变化的曲线。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明的目的是提供一种脉冲神经网络建模方法、系统及其应用，对待模拟真实神经网络的空间信息进行考虑，即考虑神经元位置信息、网络中的结构约束信息以及子模块的位置信息等，基于空间信息来更为准确的建立用于模拟模块化的生物活体神经回路和体外培养神经网络的脉冲神经网络，进而更为精准的建立网络动态特性与单个神经元的实时状态之间的联系，从而为研究隐藏在生物记忆和学习行为背后的神经信号处理和传输奠定基础。

为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图和具体实施方式对本发明作进一步详细的说明。

实施例1：

本实施例用于提供一种脉冲神经网络建模方法，如图1所示，所述建模方法包括：

S1：构建指数型积分放电脉冲神经元模型；所述指数型积分放电脉冲神经元模型用于基于接收到的所有输入电流改变膜电位；

具体的，分别建立膜电位计算公式和输入电流总和计算公式，得到指数型积分放电脉冲神经元模型。

本实施例的建模方法考虑脉冲神经网络的基本单元为指数型积分放电神经元，其在t时刻的膜电位V(t)由式(1)所述的膜电位计算公式来决定：

式(1)中，τ_m为神经元的膜时间常数；u_rest为静息膜电位；Δ_T为指数型修正项的锐度；u_rh为指数型修正项的相对阈值参数；R_m为膜电阻；I_input(t)为t时刻输入神经元的输入电流总和，其由输入电流总和计算公式计算得到。

所谓指数型积分放电脉冲神经元模型，即神经元实时整合其接收的所有输入电流，得到输入电流总和，根据输入电流总和以改变其膜电位。当在某一时刻t^f，膜电位V(t)突破某一绝对阈值θ_reset时，该神经元模型出现动作电位(即放电)，这时将时刻t^f之后，在绝对不应期时间Δ_abs以内的膜电位设置为低于绝对阈值θ_reset，而高于静息膜电位u_rest的重置电位，重置电位V_reset＝θ_reset-Δ_reset，Δ_reset为重置电位相对于绝对阈值的下调量。

本实施例的神经元模型包括兴奋性神经元和抑制性神经元两种，其中兴奋性神经元能作为突触前神经元和其他神经元形成兴奋性突触，兴奋性突触包含两种以谷氨酸盐为神经递质的突触：AMPA突触和NMDA突触；抑制性神经元能作为突触前神经元和其他神经元形成抑制性突触，抑制性突触包括GABA突触。上述三种突触分别对应包含三种受体：AMPA受体，NMDA受体和GABA受体，使得这三种突触会被AMPA(α-氨基-3-羟基-5-甲基-4-异恶唑丙酸)，NMDA(N-甲基-D-天门冬氨酸)和GABA(γ-氨基丁酸)这三种物质所激活。AMPA突触即为包含α-氨基-3-羟基-5-甲基-4-异恶唑丙酸受体的突触，NMDA突触即为包含N-甲基-D-天门冬氨酸受体的突触，GABA突触即为包含γ-氨基丁酸受体的突触。

此外，考虑到增强神经元模型的真实性和适应性，本实施例在神经元模型的输入项添加了钙离子浓度依赖的钾离子流和相对不应性电流。因此每一神经元作为突触后神经元，均以AMPA突触电流、NMDA突触电流、GABA突触电流、钙离子浓度依赖的钾离子流和相对不应性电流的电流和作为输入电流总和，输入电流总和计算公式如式(2)：

I_input(t)＝I_AMPA(t)+I_NMDA(t)+I_GABA(t)+I_K(Ca)(t)+I_ref(t)； (2)

式(2)中，I_AMPA(t)为t时刻的AMPA突触电流；+I_NMDA(t)为t时刻的NMDA突触电流；I_GABA(t)为t时刻的GABA突触电流；I_K(Ca)(t)为t时刻的钙离子浓度依赖的钾离子流；I_ref(t)为t时刻的相对不应性电流。

所有突触模型均为基于电导变化的模型，即每一次突触前神经元的放电事件，会引起突触后神经元的电导变化，从而产生电流变化。以上三类突触(AMPA突触、NMDA突触、GABA突触)的电导均由通用公式(3)决定：

式(3)中，g(t)为t时刻的突触电导；g_syn为最大突触电导，三类突触的最大突触电导分别为g_AMPA，g_NMDA和g_GABA；τ_r为突触电导的上升时间常数；τ_d为突触电导的衰减时间常数；t^f为突触前神经元产生动作电位事件的时间；t_d为该动作电位事件从发生到传递到突触后神经元的时间间隔。

AMPA突触电流由下式(4)决定：

式(4)中，g_j(t)表示t时刻神经元的第j个突触的突触电导，在式(4)中是指以该神经元作为突触后神经元，与之形成AMPA突触的第j个AMPA突触电导，其根据式(3)计算得到，

则是指所有AMPA突触电导的和；V_E为兴奋性突触的反转电位。

NMDA突触电流由式(5)决定：

相比式(4)，式(5)中的g_j(t)是指以该神经元作为突触后神经元，与之形成NMDA突触的第j个NMDA突触电导，其根据式(3)计算得到，

则是指所有NMDA突触电导的和；式(5)还添加了一项电压依赖的修正项，这是由于NMDA突触在较低膜电位时，会被Mg²⁺离子阻塞，从而具有电压依赖属性，而细胞外液中的Mg²⁺浓度在一次实验中一般保持恒定，因此此处的两个常数c_NMDA,1和c_NMDA,2均为Mg²⁺浓度相关常数。

GABA突触电流由下式(6)决定。

式(6)中，V_I为抑制性突触的反转电位；式(6)中的g_j(t)是指以该神经元作为突触后神经元，与之形成GABA突触的第j个GABA突触电导，其根据式(3)计算得到，

则是指所有GABA突触电导的和。

钙离子浓度依赖的钾离子流的计算公式如下：

I_K(Ca)(t)＝g_K(Ca)c(t)(E_K-V(t))； (7)

式(7)中，g_K(Ca)为钙离子浓度依赖的钾离子流电导；c(t)为t时刻的胞内钙离子浓度；E_K为钾离子流的反转电位。

其中，胞内钙离子浓度由式(8)决定：

式(8)中，c_step为单次动作电位下伴随着电压门控的钙通道进入胞内的钙离子浓度增加量；δ(·)为狄拉克函数；t_k为t时刻之前上一次动作电位发生的时间；

为发生动作电位的总次数；τ_Ca为钙离子浓度衰减的时间常数。需要说明的是，应该从仿真计算的角度来看，逐时间步长计算，当运行到时刻t_a，神经元a产生了动作电位，这样，对于该神经元a，t_a就被记为t_k，t_k是随着仿真进行而被更新的，它永远表示某一神经元发生上一次动作电位的时刻。

相对不应性电流由式(9)决定：

式(9)中，g_ref为相对不应性电流的最大电导；τ_ref为相对不应性电流的衰减时间常数；A_ref(t-t_k)为使相对不应性电流仅考虑上一次动作电位，而不考虑之前动作电位的影响的修正项；

本实施例的脉冲神经网络中的单个神经元采用指数型积分放电脉冲神经元模型，具有基于输入电流自适应放电阈值的特点，能更准确的模拟不同刺激下真实神经元的反应。同时，在计算输入电流总和时，还引入两种兴奋性突触和一种抑制性突触，能够用于研究兴奋-抑制平衡对脉冲神经网络动态特性的影响，还加入对神经元内钙离子浓度状态以及基于钙离子的钾离子流的建模，可用于研究脉冲神经网络中的簇放电行为，以及其伴随的钙离子动态。

S2：基于待模拟真实神经网络的空间信息对所有神经元进行选择性连接，得到对应于所述待模拟真实神经网络的脉冲神经网络；所述空间信息包括网络结构约束信息、神经元位置信息以及子模块位置信息；每一所述神经元均采用所述指数型积分放电脉冲神经元模型。

本实施例考虑在生理环境下，待模拟真实神经网络中的神经元两两形成突触的概率与其空间距离存在关系，因此对于分布均匀的待模拟真实神经网络，可以根据其两两间的欧式距离确定其连通概率以生成神经元网络。此外，对于分布非均匀的模块化的待模拟真实神经网络，同时需要考虑模块内与模块间的不同情况，确定其连通概率。本实施例的方法对以上两种情况具有通用性。

考虑待模拟真实神经网络中含有n个神经元，其中m个为兴奋性神经元，n-m个为抑制性神经元，给神经元编号，其中1到m的序号对应兴奋性神经元，其余对应抑制性神经元。

具体的，如图2所示，S2可以包括：

S21：基于所述待模拟真实神经网络的空间信息计算任意两个所述神经元之间的连通概率；

式(11)中，P_i,j为神经元i与神经元j之间的连通概率；P_base为基准概率，其由突触前后神经元的自身属性决定，即是兴奋性神经元还是抑制性神经元，具体的，突触前后神经元，即突触前神经元和突触后神经元有四种组合方式，兴奋-兴奋，兴奋-抑制，抑制-兴奋和抑制-抑制，这四种组合分别会对应不同的基准概率；d_i,j为神经元i与神经元j之间的欧式距离，其根据神经元位置信息计算得到，所采用的神经元位置信息忽略神经元的尺寸及形貌参数，仅具有空间三维坐标；

(x_i，y_i，z_i)为第i个神经元的空间三维坐标，(x_j，y_j，z_j)为第j个神经元的空间三维坐标；σ_Ed为均匀神经网络中两神经元间连通概率的正态分布的概率常数；P_m为模块化地形信息修正项，其根据子模块位置信息计算得到；P_uni为模块间有向性修正项，其根据网络结构约束信息计算得到。

式(12)中，d_m为神经元i与神经元j所属模块之间的欧氏距离；σ_m为模块间神经元连通概率的修正概率常数。

P_uni为模块间有向性修正项，如由地形约束，导致神经元i与神经元j所属模块为单向连通或者不连通的情况。

当神经元i与神经元j所属模块连通时，P_uni等于1，否则，P_uni等于0。

S22：根据所述连通概率计算二元连通矩阵；

具体的，对于每一连通概率，随机在0-1范围内生成伪随机数；根据伪随机数与连通概率的比较结果，确定每一连通概率对应的二元连通矩阵的元素值，得到二元连通矩阵。

根据神经元i和j的连通概率P_i,j，确定其对应二元连通矩阵的元素M_i,j的方法为，随机在0到1范围内生成伪随机数x_rand，

连通信息存储在n×n二元连通矩阵M中，其中n为网络中神经元的数量。

S23：根据所述二元连通矩阵对所有所述神经元进行选择性连接，得到对应于所述待模拟真实神经网络的脉冲神经网络。

二元连通矩阵的项M_i,j＝1表示神经元i向j形成了突触，M_i,j＝0则表示神经元i向j没有形成了突触，M_j,i＝1则表示神经元j向神经元i形成了突触，进而根据二元连通矩阵中各个元素的值连接相应的神经元，当M_i,j＝1，则按照从神经元i到神经元j的方向对两个神经元进行连接，从而基于空间信息连接神经元形成脉冲神经网络。

本实施例的建模方法将生物活体神经回路或体外培养神经网络等待模拟真实神经网络中的空间信息纳入考虑范畴，包括神经元位置信息、网络结构约束信息以及子模块位置信息等，具有能更为精准地建模真实模块化神经网络的连通关系的优点。

在建立脉冲神经元模型并将其连通形成脉冲神经网络之后，在结构和功能上，该神经网络已经接近待模拟真实神经网络。为进一步模拟真实神经网络的随机性和动态性，在得到对应于待模拟真实神经网络的脉冲神经网络后，本实施例的建模方法还包括对脉冲神经网络增加噪声电流驱动项，即在脉冲神经网络中加入模拟真实神经系统的噪声驱动项。

随机噪声事件由频率为F_noise的泊松过程产生，且仅施加给兴奋性神经元(即编号为1至m的神经元)，对于神经元x，x＝1，2，...，m，平均每1/F_noise秒会出现一次噪声事件。每一次噪声事件会从事件发生起，持续使某一神经元产生噪声电流，噪声电流I_noise由以下公式决定：

式(15)中，I_noise(t)为t时刻的噪声电流；A_noise为噪声电流的幅值；

为第k次噪声事件发生的时间；r_N和τ_N均为噪声事件的时间常数；阿尔法函数α(a,b,c)为：

将噪声电流加入兴奋性神经元模型的输入电流总和，则每一兴奋性神经元的输入电流总和计算公式包括：以AMPA突触电流、NMDA突触电流、GABA突触电流、钙离子浓度依赖的钾离子流、相对不应性电流和噪声电流的电流和作为输入电流总和，因此每一兴奋性神经元的输入电流总和计算公式如下：

I_input(t)＝I_AMPA(t)+I_NMDA(t)+I_GABA(t)+I_K(Ca)(t)+I_ref(t)+I_noise(t)； (17)

本实施例的建模方法基于指数型脉冲放电脉冲神经元模型，具有基于输入电流自适应放电阈值的特点，能更准确的模拟不同刺激下真实神经元的反应。还包含两种兴奋性突触和一种抑制性突触，能用于研究兴奋-抑制平衡对网络动态特性的影响。还加入对神经元内钙离子浓度状态以及基于钙离子的钾离子流的建模，可用于研究神经网络中的簇放电行为，以及其伴随的钙离子动态。将生物活体神经回路或体外培养神经网络中的空间信息纳入考虑范畴，包括神经元位置信息、网络结构约束信息以及子模块位置信息等，具有能更为精准地建模真实模块化神经网络的连通关系的优点。本实施例的建模方法能提升单个神经元的建模精度，扩大脉冲神经网络的应用范围，以及提高脉冲神经网络对于非均匀神经网络建模的适用性。利用该建模方法所建立的脉冲神经网络能够模拟真实的生物神经网络的功能性，研究生物活体神经回路或体外培养神经网络和解释甚至破解其编解码方式，并研究其自发电活动、信息处理、信号传导以及记忆和学习等功能。

实施例2：

本实施例用于提供一种脉冲神经网络建模系统，如图3所示，所述建模系统包括：

构建模块M1，用于构建指数型积分放电脉冲神经元模型；所述指数型积分放电脉冲神经元模型用于基于接收到的所有输入电流改变膜电位；

连接模块M2，用于基于待模拟真实神经网络的空间信息对所有神经元进行选择性连接，得到对应于所述待模拟真实神经网络的脉冲神经网络；所述空间信息包括网络结构约束信息、神经元位置信息以及子模块位置信息；每一所述神经元均采用所述指数型积分放电脉冲神经元模型。

实施例3：

本实施例用于提供一种脉冲神经网络的网络状态计算方法，基于固定时间步长对实施例1所述的脉冲神经网络的网络状态进行计算，如图4所示，所述计算方法包括：

T1：采用向前欧拉法，逐时间步长、逐个神经元计算所有神经元的膜电位和钙离子浓度；

T2：根据研究需要，记录每一所述时间步长每一所述神经元的所述膜电位和所述钙离子浓度，并记录脉冲神经网络中出现动作电位的神经元序号以及出现动作电位的时间。

在利用实施例1所述的建模方法完成完整的脉冲神经网络建模过程后，需要对模型进行时间离散化以便进行数值计算。本计算方法采用固定时间步长δt作为数值计算的最小单位，每增加时间步长δt，计算此时刻的膜电位V(t+δt)和钙离子浓度c(t+δt)，总仿真时间为T。数值计算方面，采用向前欧拉法，逐时间步长，逐个神经元，计算所有神经元的膜电位和钙离子浓度。

具体的，初始时刻t＝0，膜电位为-74mV，钙离子浓度为0mM，动作电位的数量S_count＝0。嵌套使用两个循环函数，外循环为时间，内循环为神经元序号。神经元i在时刻t+δt的膜电位为：

需要说明的是，先依据t时刻的膜电位，利用输入电流总和计算公式得到t时刻的输入电流总和，然后利用式(18)计算t+δt的膜电位。

神经元i在时刻t+δt的钙离子浓度为：

如果在时刻t+δt，神经元膜电位V_i(t+δt)大于绝对阈值θ_reset，则神经元i产生动作电位，动作电位数量S_count＝S_count+1，将此时刻t+δt和放电神经元的编号i记录为矩阵M_spiking的第S_count行，其时刻t+δt的膜电位被设置为一次典型动作电位的峰值E_peak，且从t+2δt到t+Δ_abs的所有时刻，该神经元处于绝对不应期，其膜电位被固定为V_reset，在之后的时间循环中不需要计算。同时，神经元i的钙离子浓度阶跃增加量c_step。此外，该动作电位的产生，同时也将对其作为突触前神经元所连接的突触后神经元产生影响，使突触后神经元的突触电导产生改变。根据连通矩阵M的第i行，确定其所连接的神经元编号，并根据公式(3)计算其改变量。此外，由于根据公式，突触电导的影响从产生起渐进衰减而不完全消失，因此需要保存每一神经元在每一时刻的突触电导，但实际上，在7倍上升和衰减时间步长后，突触电导已经衰减到其最大电导的千分之一。为了优化数值计算的存储空间，仅保存从时刻t起，到t+7(τ_r+τ_d)间的突触电导的数值。

为优化数值计算的内存占用，每一神经元的膜电位和钙离子浓度随时间的变化并未被保存，但仅需要在每一时刻，存储中间变量V_i(t+δt)和c_i(t+δt)，即可实现对任一神经元的膜电位和钙离子浓度随时间的变化的观察。当时间t＝T时，数值计算结束。脉冲网络出现动作电位的神经元序号以及其产生的时间被保存在矩阵M_spiking中，便于后续分析。

利用上述计算方法对脉冲神经网络的网络状态进行计算，除记录上述数值之外，还可以根据研究需要记录其他数值，以利用该脉冲神经网络对所想模拟的待模拟真实神经网络进行研究，模拟真实的生物神经网络的功能性，研究生物活体神经回路或体外培养神经网络和解释甚至破解其编解码方式。

如图5所示，其为模拟真实神经网络的整体流程图，先构建指数型积分放电神经元模型，再基于空间信息连接神经元形成脉冲神经网络，并向脉冲神经网络加入噪声驱动项，以完成脉冲神经网络的建模过程，对待模拟真实神经网络进行模拟。然后基于固定时间步长计算网络状态，以研究真实神经网络自发电活动、信息处理、信号传导以及记忆和学习等功能。

以下，以一具体实例对图5所述的方法进行进一步的说明：

1)构建指数型积分放电神经元模型；

考虑脉冲神经网络的基本单元为指数型积分放电神经元，利用S1的方法建立指数型积分放电脉冲神经元模型。此时，该神经元模型中可取神经元的膜时间常数τ_m为40ms；静息膜电位u_rest为-74mV；指数型修正项的锐度Δ_T＝0.5；指数型修正项的相对阈值参数u_rh＝-56mV；膜电阻R_m＝40MΩ；绝对阈值θ_reset＝-54mV；绝对不应期时间Δ_abs＝5ms；重置电位V_reset＝-64mV；三种突触的最大电导分别为g_AMPA＝5nS，g_NMDA＝25nS和g_GABA＝40nS；AMPA突触，NMDA突触和GABA突触的上升和衰减时间常数分别为：τ_{r_AMPA}＝0.2ms，τ_{d_AMPA}＝5.3ms，τ_{r_NMDA}＝5ms，τ_{d_NMDA}＝100ms，τ_{r_GABA}＝1.2ms和τ_{d_NMDA}＝15ms；兴奋性突触的反转电位V_E＝0mV；Mg²⁺浓度相关常数c_NMDA,1＝0.28和c_NMDA,2＝-0.062；抑制性突触的反转电位V_I＝-75mV；钾离子流的反转电位E_K＝-84mV；单次动作电位伴随着电压门控的钙通道进入胞内的钙离子浓度增加量c_step＝0.1μM；钙离子浓度衰减的时间常数τ_Ca＝2500ms；相对不应性电流的衰减时间常数τ_ref＝2ms。

2)基于空间信息连接神经元形成神经网络；

本实例考虑一种体外培养的模块化神经元网络，其由单向连通的两个模块构成，如图6所示，每个模块均包含25个兴奋性神经元，模块内部神经元之间的连通概率由其欧氏距离决定，距离越近，连通概率越高，最大连通概率为0.5。模块1到模块2构成10个突触连接，模块2向模块1不构成突触连接，因此被称为单向连通。

因此，考虑神经网络中含有n个神经元，其中m个为兴奋性神经元，n＝m＝50，给神经元编号，其中1到25的序号对应模块1中的神经元，26到50的序号对应模块2中的神经元。

利用S2计算二元连通矩阵，此时，基准概率P_base＝0.5；均匀神经网络中两神经元连通概率的正态分布的概率常数σ_Ed＝1；模块1与模块2的欧式距离为d_m＝400μm；模块间神经元连通概率的修正概率常数σ_m＝250μm。

3)加入模拟真实神经系统的噪声电流驱动项；

在建立脉冲神经元模型并将其连通形成脉冲神经网络之后，在结构和功能上，该神经网络已经接近真实神经系统。为进一步模拟真实神经系统的随机性和动态性，采用式(15)-式(17)增加噪声电流驱动项。此时频率F_noise＝0.4Hz；噪声电流的幅值A_noise＝1000pA；噪声事件的时间常数r_N＝30ms、τ_N＝50ms。

4)基于固定时间步长计算网络状态；

完成了完整的脉冲神经网络建模过程后，需要对模型进行时间离散化以便进行数值计算。本方法采用固定时间步长δt作为数值计算的最小单位，取δt＝0.1ms，每增加时间步长δt，利用式(18)和式(19)计算此时刻的膜电位V(t+δt)和钙离子浓度c(t+δt)。总仿真时间为T，取T＝300s。

如果在时刻t+δt，神经元膜电位V_i(t+δt)大于绝对阈值θ_reset，神经元i产生动作电位，动作电位数量S_count＝S_count+1，将此时刻t+δt和放电神经元的编号i记录为矩阵M_spiking的第S_count行。其t+δt时刻的膜电位被设置为一次典型动作电位的峰值E_peak，取E_peak＝40mV，且从t+2δt到t+Δ_abs的所有时刻，该神经元处于绝对不应期，其膜电位被固定为V_reset，在之后的时间循环中不需要计算，Δ_abs＝5ms，V_reset＝-64mV。同时，钙离子浓度阶跃增加量c_step＝0.1mM。此外，该动作电位的产生，同时也将对其作为突触前神经元所连接的突触后神经元产生影响，使其突触电导产生改变。根据连通矩阵M，通过其第i行，确定其所连接的神经元编号，并根据公式(3)计算其改变量。此外，由于根据公式，突触电导的影响从产生起渐进衰减而不完全消失，因此需要保存每一神经元在每一时刻的突触电导，但实际上，在7倍上升和衰减时间步长后，突触电导已经衰减到其最大电导的千分之一。为了优化数值计算的存储空间，仅保存从时刻t起，到t+7(τ_r+τ_d)间的突触电导的数值。

为优化数值计算的内存占用，每一神经元的膜电位和钙离子浓度随时间的变化并未被保存，但仅需要在每一时刻，存储中间变量V_i(t+δt)和c_i(t+δt)，即可实现对任一神经元的膜电位和钙离子浓度随时间的变化的观察。当时间t＝T时，数值计算结束。脉冲网络出现动作电位的神经元序号以及其产生的时间被保存在矩阵M_spiking中，便于后续分析。

根据上述数据，可分别绘制动作事件对应的动作电位栅格图、一次网络局部簇放电事件对应的动作电位栅格图、一次网络整体簇放电事件对应的动作电位栅格图、各个神经元的膜电位和胞内钙离子浓度随时间变化的曲线。如图7所示，其为神经元的动作事件对应的动作电位栅格图；如图8所示，其为仿真过程中一次网络局部簇放电事件对应的动作电位栅格图；如图9所示，其为仿真过程中一次网络整体簇放电事件对应的动作电位栅格图；如图10所示，其为编号为1的神经元的膜电位随时间变化的曲线；如图11所示，其为编号为1的神经元的胞内钙离子浓度随时间变化的曲线。

本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。对于实施例公开的系统而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。

本文中应用了具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处。综上所述，本说明书内容不应理解为对本发明的限制。

Claims (10)

1.一种脉冲神经网络建模方法，其特征在于，所述建模方法包括：

构建指数型积分放电脉冲神经元模型；所述指数型积分放电脉冲神经元模型用于基于接收到的所有输入电流改变膜电位；

基于待模拟真实神经网络的空间信息对所有神经元进行选择性连接，得到对应于所述待模拟真实神经网络的脉冲神经网络；所述空间信息包括网络结构约束信息、神经元位置信息以及子模块位置信息；每一所述神经元均采用所述指数型积分放电脉冲神经元模型。

2.根据权利要求1所述的建模方法，其特征在于，所述构建指数型积分放电脉冲神经元模型具体包括：

分别建立膜电位计算公式和输入电流总和计算公式，得到所述指数型积分放电脉冲神经元模型。

3.根据权利要求2所述的建模方法，其特征在于，所述输入电流总和计算公式具体包括：

以AMPA突触电流、NMDA突触电流、GABA突触电流、钙离子浓度依赖的钾离子流和相对不应性电流的电流和作为输入电流总和。

4.根据权利要求3所述的建模方法，其特征在于，所述钙离子浓度依赖的钾离子流的计算公式包括：

I_K(Ca)(t)＝g_K(Ca)c(t)(E_K-V(t))；

其中，I_K(Ca)(t)为t时刻的钙离子浓度依赖的钾离子流；g_K(Ca)为钙离子浓度依赖的钾离子流电导；c(t)为t时刻的胞内钙离子浓度；E_K为钾离子流的反转电位；V(t)为t时刻的膜电位。

5.根据权利要求1所述的建模方法，其特征在于，所述基于待模拟真实神经网络的空间信息对所有神经元进行选择性连接，得到对应于所述待模拟真实神经网络的脉冲神经网络具体包括：

基于所述待模拟真实神经网络的空间信息计算任意两个所述神经元之间的连通概率；

根据所述连通概率计算二元连通矩阵；

根据所述二元连通矩阵对所有所述神经元进行选择性连接，得到对应于所述待模拟真实神经网络的脉冲神经网络。

6.根据权利要求5所述的建模方法，其特征在于，所述基于所述待模拟真实神经网络的空间信息计算任意两个所述神经元之间的连通概率具体包括：

其中，P_i,j为神经元i与神经元j之间的连通概率；P_base为基准概率；d_i,j为神经元i与神经元j之间的欧式距离，其根据所述神经元位置信息计算得到；σ_Ed为均匀神经网络中两神经元间连通概率的正态分布的概率常数；P_m为模块化地形信息修正项，其根据所述子模块位置信息计算得到；P_uni为模块间有向性修正项，其根据所述网络结构约束信息计算得到。

7.根据权利要求5所述的建模方法，其特征在于，所述根据所述连通概率计算二元连通矩阵具体包括：

对于每一所述连通概率，随机在0-1范围内生成伪随机数；

根据所述伪随机数与所述连通概率的比较结果，确定每一所述连通概率对应的二元连通矩阵的元素值，得到所述二元连通矩阵。

8.根据权利要求1所述的建模方法，其特征在于，在得到对应于所述待模拟真实神经网络的脉冲神经网络后，所述建模方法还包括：对所述脉冲神经网络增加噪声电流驱动项；具体包括：

每一兴奋性神经元的输入电流总和计算公式包括：以AMPA突触电流、NMDA突触电流、GABA突触电流、钙离子浓度依赖的钾离子流、相对不应性电流和噪声电流的电流和作为输入电流总和。

9.一种脉冲神经网络建模系统，其特征在于，所述建模系统包括：

构建模块，用于构建指数型积分放电脉冲神经元模型；所述指数型积分放电脉冲神经元模型用于基于接收到的所有输入电流改变膜电位；

连接模块，用于基于待模拟真实神经网络的空间信息对所有神经元进行选择性连接，得到对应于所述待模拟真实神经网络的脉冲神经网络；所述空间信息包括网络结构约束信息、神经元位置信息以及子模块位置信息；每一所述神经元均采用所述指数型积分放电脉冲神经元模型。

10.一种脉冲神经网络的网络状态计算方法，基于固定时间步长对权利要求1-8任一项所述的脉冲神经网络的网络状态进行计算，其特征在于，所述计算方法包括：

采用向前欧拉法，逐时间步长、逐个神经元计算所有神经元的膜电位和钙离子浓度；

根据研究需要，记录每一所述时间步长每一所述神经元的所述膜电位和所述钙离子浓度，并记录脉冲神经网络中出现动作电位的神经元序号以及出现动作电位的时间。

Images