CN109325596B

CN109325596B - 一种基于钙浓度的神经突触可塑性计算方法

Info

Publication number: CN109325596B
Application number: CN201811083069.XA
Authority: CN
Inventors: 曹立宏; 王晔
Original assignee: Communication University of China
Current assignee: Communication University of China
Priority date: 2018-09-17
Filing date: 2018-09-17
Publication date: 2024-03-26
Anticipated expiration: 2038-09-17
Also published as: CN109325596A

Abstract

一种基于钙浓度的神经突触可塑性计算方法涉及大脑模拟仿真领域，尤其是针对大脑模拟仿真中脉冲神经网络的神经突触可塑性的计算问题。首先是钙离子浓度的计算，针对需要计算的突触，根据初始时刻t₀的突触前神经元和突触后神经元的膜电位值，以及突触的初始连接权重w₀，分别计算下一时刻t₁的树突和树突棘中的钙离子浓度；其次，根据树突棘中的钙离子浓度，通过和钙离子浓度阈值Ca0_s及Ca1_s对比得到权重变化的方向，根据突触状态标签Tag和可塑性相关蛋白质PRP浓度，计算突触权重的变化，得到t₁时刻的新权重。重复以上过程，计算模拟时间内的突触连接强度。本发明应用于构建类脑神经网络，完成类脑智能所需的学习和记忆过程的仿真，实现通用强人工智能，并应用于智能传媒、医疗等方面。

Description

一种基于钙浓度的神经突触可塑性计算方法

技术领域

本发明涉及大脑模拟仿真领域，尤其是针对大脑模拟仿真中脉冲神经网络的神经突触可塑性的计算问题。

背景技术

大脑中包含数以亿万计的神经元，这些神经元相互连接构成神经网络，神经元连接的部分称为突触。以人脑为例，每个神经元会与几千个至几万个神经元产生连接，因此，大脑中的连接数量十分庞大。感知和认知过程会改变突触连接的强弱，也就是权重。这种突触连接权重变化的现象称为突触可塑性，是大脑学习和记忆的主要方式。如何准确地计算这些突触权重的变化过程是进行大脑模拟仿真的关键。最接近现有技术的神经突触可塑性计算方法为Graupner等人在 Calcium-based plasticity model explains sensitivityof synaptic changes to spike pattern,rate,and dendritic location一文中描述的方法，其问题在于，并没有考虑神经元内蛋白质浓度和树突棘状态对突触可塑性的影响，导致了计算不准确等问题。

发明内容

本发明所要解决的技术问题是在大脑模拟中突触可塑性模型和计算问题。

为了解决这一问题，本发明采用的技术方案分为以下五个步骤：

1.钙离子浓度的计算：针对需要计算的突触，根据初始时刻t₀的突触前神经元和突触后神经元的膜电位值，以及突触的初始连接权重w₀，分别计算下一时刻t₁的树突和树突棘中的钙离子浓度，其中w₀是非负常数，在实际应用中，建议取值为1，原因是在初始时刻做归一化处理；

2.确定权重变化方向(Flag)：根据树突棘中的钙离子浓度，通过和钙离子浓度阈值Ca0_s及Ca1_s对比得到权重变化的方向，其中Ca0_s和Ca1_s均为非负常数，且Ca1_s>Ca0_s，在实际应用中，根据经验，建议Ca0_s＝0.5mM,Ca1_s＝1mM；

3.计算突触状态标签(Tag)：根据权重变化方向及突触状态变化常数α_T，β_T,LTD和β_T,LTP，计算突触标签状态，其中α_T，β_T,LTD和β_T,LTP均为大于零的常数，在实际应用中，根据经验，建议α_T＝0.01,β_T,LTD＝0.5,β_T,LTP＝0.5；

4.计算可塑性相关蛋白质(PRP)浓度：根据树突中的钙离子浓度计算可塑性相关蛋白质浓度；

5.计算突触权重的变化：根据突触状态标签Tag和可塑性相关蛋白质PRP 浓度，计算突触权重的变化，得到t₁时刻的新权重。

重复以上过程，计算模拟时间内的突触连接强度。

本发明的内容可以应用于构建类脑神经网络，完成类脑智能所需的学习和记忆过程的仿真，实现通用强人工智能，并应用于智能传媒、智能医疗、智能制造、智能教育和自动驾驶等方面。

根据本发明的第一个步骤，采用如下方法：

首先，设定突触前神经元和突触后神经元的初始膜电位值以及突触前神经元和突触后神经元之间的初始连接权重w₀。其中突触前神经元和突触后神经元的初始膜电位值为小于零的常数，在实际应用中，建议设置为-70mV，原因是这是神经元的静息态膜电位，突触前神经元和突触后神经元之间的初始连接权重w₀是非负常数，在实际应用中，建议取值为1，原因是在初始时刻做归一化处理。

其次，采用脉冲神经元模型对突触前神经元和突触后神经元进行建模，其中脉冲神经元模型来自于神经科学领域文献，模型的输入为各离子通道的电导、各离子通道的反转电位、输入电流和膜电容，以上输入数据均来自于神经科学领域文献，模型的输出为膜电位值。

再次，根据突触前神经元的脉冲神经元模型，利用欧拉法或龙格库塔法等数值计算方法计算t₁时刻突触前神经元的膜电位。

然后，根据突触动力学模型对突触连接进行建模，其中突触动力学模型来自于神经科学领域文献，模型的输入为t₁时刻突触前神经元膜电位，t₁时刻突触后神经元膜电位，突触动力学时间常数以及突触受体的反转电位，突触动力学时间常数以及突触受体的反转电位均来自于神经科学领域参考文献，模型的输出为突触连接对突触后神经元产生的电流。

最后，根据突触动力学模型及突触后神经元的神经元模型，利用欧拉法或龙格库塔法等数值计算方法计算t₁时刻突触后神经元的膜电位以及突触后神经元树突和树突棘中的钙离子浓度。突触后神经元树突中的钙离子源自位于细胞膜上的钙离子通道，突触后神经元树突棘中的钙离子源自位于树突棘上允许钙离子通过的突触受体。在t₁时刻由钙离子通道进入神经元内的钙离子浓度[Ca²⁺]_Ca-channel可以由包含该离子通道的脉冲神经元模型计算得出，在t₁时刻由树突棘上允许钙离子通过的突触受体进入神经元内的钙离子浓度[Ca²⁺]_receptor可以由该突触的动力学模型计算得出。突触后神经元树突和树突棘中的钙离子可以相互扩散，记t₁时刻树突棘和该树突棘所处的树突中的钙离子浓度分别为 [Ca²⁺]_s和[Ca²⁺]_d。

根据本发明的第二个步骤，采用如下方法：

根据树突棘中的钙离子浓度[Ca²⁺]_s和钙离子浓度阈值Ca0_s及Ca1_s对比得到 t₁时刻该树突棘所处的突触的权重变化方向Flag。其中Ca0_s和Ca1_s均为非负常数，且Ca1_s>Ca0_s，在实际应用中，根据经验，建议Ca0_s＝0.5mM,Ca1_s＝1mM。如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s低于阈值Ca0_s，则Flag为0，该树突棘所处的突触不会发生连接权重变化；如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s处于阈值Ca0_s和阈值Ca1_s之间，则Flag为-1，该树突棘所处的突触连接权重变化方向为减弱；如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s高于阈值Ca1_s，则Flag 为1，该树突棘所处的突触的连接权重变化方向为增加。

根据本发明的第三个方面，采用如下方法：

根据树突棘所处突触权重变化方向及突触状态变化常数α_T，β_T,LTD和β_T,LTP得到t₁时刻该树突棘所处的突触状态标签Tag。如果Tag为0，则该树突棘所处的突触不会发生连接权重变化；如果Tag大于0，则该树突棘所处的突触进入连接权重可以增加的状态；如果Tag小于0，则该树突棘所处的突触进入连接权重可以减弱的状态。Tag的演化是一个与Flag有关的动力学过程。

如果t₁时刻Flag为0，则其中α_T是大于零的常数，在实际应用中，根据经验，建议α_T＝0.01，若Tag以t₀时刻初始值为0开始演化，则Tag关于时间的导数为0，Tag不发生变化；

如果t₁时刻Flag为-1，则其中β_T,LTD是大于零的常数，在实际应用中，根据经验，建议β_T,LTD＝0.5，若Tag以t₀时刻初始值为0开始演化，则Tag关于时间的导数为负，Tag降低，直至降到Tag关于时间的导数为0时，即/>时停止下降，并保持/>

如果t₁时刻Flag为1，则其中β_T,LTP是大于零的常数，，在实际应用中，根据经验，建议β_T,LTP＝0.5，若Tag以t₀时刻初始值为0开始演化，则Tag关于时间的导数为正，Tag升高，直至升至Tag关于时间的导数为0时，即/>时停止上升，并保持/>

根据本发明的第四个步骤，采用如下方法：

根据树突中的钙离子浓度[Ca²⁺]_d和钙离子浓度阈值Ca0_d得到t₁时刻树突中PRP的浓度。PRP浓度的演化是与[Ca²⁺]_d有关的动力学过程。其中Ca0_d为非负常数，在实际应用中，根据经验，建议Ca0_d＝1mM。树突中的钙离子可以引发生成PRP的化学反应，单位时间内PRP的生成率记为PRP_rate。如果t₁时刻树突中的钙离子浓度[Ca²⁺]_d低于阈值Ca0_d，则PRP_rate为0或衰减为0，PRP不发生变化；如果t₁时刻树突中的钙离子浓度[Ca²⁺]_d高于阈值Ca0_d，则PRP_rate按 alpha函数的形式变化，PRP作为PRP_rate对时间的积分，先增加，后持平。

根据本发明的第五个步骤，采用如下方法：

首先，计算突触连接权重系数相关参数y，y的演化是与Tag和PRP浓度有关的动力学过程，其中τ_y是大于零的常数，在实际应用中，根据经验，τ_y＝1ms。由于PRP浓度是非负的，若Tag和PRP浓度中有一个为0，则 y对时间的导数为0，y不会发生变化；若Tag为负，PRP浓度为正，则y对时间的导数为负，y降低；若Tag为正，PRP浓度为正，则y对时间的导数为正，y 升高。

其次，计算突触连接权重系数z，z与y之间的关系满足其中λ，z_l，z_h均为大于零的常数，z_l代表权重系数z 的下限，z_h代表权重系数z的上限，在实际应用中，根据经验，λ＝1,z_l＝0.5,z_h＝ 2。

再次，将突触连接权重系数z乘以初始连接权重w₀，得到t₁时刻的连接权重w₁＝z×w₀。

最后，重复本发明的第一个步骤、第二个步骤、第三个步骤、第四个步骤和第五个步骤，计算整个模拟时间内神经元之间的突触连接权重变化。

本发明的有益效果是大大提高了大脑模拟中突触可塑性的计算精度，并可以应用于构建类脑神经网络，完成类脑智能所需的学习和记忆过程的仿真，实现通用强人工智能，并应用于智能传媒、智能医疗、智能制造、智能教育和自动驾驶等方面。

附图说明

图1为突触后神经元树突和树突棘示意图

图2为权重变化方向Flag和突触态标签Tag的计算流程图

图3为PRP生成率和PRP浓度的计算流程图

图4为突触连接权重系数相关参数y和突触连接权重系数z的计算流程图

具体实施方式

下面将结合本发明中的附图和具体实施方式对本发明的技术方式进行详细描述。

在钙离子浓度的计算中：

然后，根据突触动力学模型对突触连接进行建模，其中突触动力学模型来自于神经科学领域文献，模型的输入为t₁时刻突触前神经元膜电位，t₁时刻突触后神经元膜电位，突触动力学时间常数以及突触受体的反转电位，模型的输出为突触连接对突触后神经元产生的电流。

最后，根据突触动力学模型及突触后神经元的神经元模型，利用欧拉法或龙格库塔法等数值计算方法计算t₁时刻突触后神经元的膜电位以及突触后神经元树突和树突棘中的钙离子浓度。图1为突触后神经元树突和树突棘示意图，突触后神经元树突中的钙离子源自位于细胞膜上的钙离子通道，突触后神经元树突棘中的钙离子源自位于树突棘上允许钙离子通过的突触受体。在t₁时刻由钙离子通道进入神经元内的钙离子浓度[Ca²⁺]_Ca-channel可以由包含该离子通道的离子通道型神经元模型计算得出，在t₁时刻由树突棘上允许钙离子通过的突触受体进入神经元内的钙离子浓度[Ca²⁺]_receptor可以由该受体的动力学模型计算得出。突触后神经元树突和树突棘中的钙离子可以相互扩散，记t₁时刻树突棘和该树突棘所处的树突中的钙离子浓度分别为[Ca²⁺]_s和[Ca²⁺]_d。

在突触的权重变化方向Flag和突触态标签Tag的计算中：

首先，根据树突棘中的钙离子浓度[Ca²⁺]_s和钙离子浓度阈值Ca0_s及Ca1_s对比得到t₁时刻该树突棘所处的突触的权重变化方向Flag。其中Ca0_s和Ca1_s均为非负常数，且Ca1_s>Ca0_s，在实际应用中，根据经验，建议Ca0_s＝0.5mM, Ca1_s＝1mM。如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s低于阈值Ca0_s，则Flag 为0，该树突棘所处的突触不会发生连接权重变化；如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s处于阈值Ca0_s和阈值Ca1_s之间，则Flag为-1，该树突棘所处的突触连接权重变化方向为减弱；如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s高于阈值 Ca1_s，则Flag为1，该树突棘所处的突触的连接权重变化方向为增加。

其次，根据树突棘所处突触权重变化方向及突触状态变化常数α_T，β_T,LTD和β_T,LTP得到t₁时刻该树突棘所处的突触状态标签Tag。如果Tag为0，则该树突棘所处的突触不会发生连接权重变化；如果Tag大于0，则该树突棘所处的突触进入连接权重可以增加的状态；如果Tag小于0，则该树突棘所处的突触进入连接权重可以减弱的状态。Tag的演化是一个与Flag有关的动力学过程。如果t₁时刻 Flag为0，则其中α_T是大于零的常数，在实际应用中，根据经验，建议α_T＝0.01，若Tag以t₀时刻初始值为0开始演化，则Tag关于时间的导数为0，Tag不发生变化；如果t₁时刻Flag为-1，则/> 其中β_T,LTD是大于零的常数，在实际应用中，根据经验，建议β_T,LTD＝0.5，若Tag以t₀时刻初始值为0开始演化，则Tag关于时间的导数为负，Tag降低，直至降到Tag关于时间的导数为0时，即/>时停止下降，并保持/>如果t₁时刻Flag为1，则/> 其中β_T,LTP是大于零的常数，，在实际应用中，根据经验，建议β_T,LTP＝0.5，若Tag以t₀时刻初始值为0开始演化，则Tag关于时间的导数为正，Tag升高，直至升至Tag关于时间的导数为0时，即/>时停止上升，并保持/>

图2为权重变化方向Flag和突触态标签Tag的计算流程图。如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s低于阈值Ca0_s，则Flag为0，若Tag以t₀时刻初始值为0开始演化，则Tag不发生变化，仍然为0；如果树突棘中的钙离子浓度 [Ca²⁺]_s处于阈值Ca0_s和阈值Ca1_s之间，则Flag为-1，若Tag以t₀时刻初始值为 0开始演化，则Tag降低，直至降到时停止下降，并保持/>如果树突棘中的钙离子浓度[Ca²⁺]_s高于阈值Ca1_s，贼Flag 为1，若Tag以t₀时刻初始值为0开始演化，则Tag升高，直至升至/>时停止上升，并保持/>在该计算流程示例中，Ca0_s＝0.5mM, Ca1_s＝1mM,α_T＝0.01,β_T,LTD＝0.5,β_T,LTP＝0.5。

在PRP生成率和PRP的计算中：

根据树突中的钙离子浓度[Ca²⁺]_d和钙离子浓度阈值Ca0_d得到t₁时刻树突中 PRP的浓度。PRP浓度的演化是与[Ca²⁺]_d有关的动力学过程。其中Ca0_d为非负常数，在实际应用中，根据经验，建议Ca0_d＝1mM。树突中的钙离子可以引发生成PRP的化学反应，单位时间内PRP的生成率记为PRP_rate。图3为PRP生成率和PRP的计算流程图。如果t₁时刻树突中的钙离子浓度[Ca²⁺]_d低于阈值Ca0_d，若PRP_rate以t₀时刻初始值为0开始演化，则PRP_rate不发生变化，仍然为0，PRP 是PRP_rate对时间的积分，因此，PRP的浓度为0；如果t₁时刻树突中的钙离子浓度[Ca²⁺]_d高于阈值Ca0_d，若PRP_rate以t₀时刻初始值为0开始演化，则PRP_rate上升，PRP是PRP_rate对时间的积分，因此，PRP的浓度增加。在该计算流程示例中， Ca0_d＝1mM。

在突触连接权重系数相关参数y和突触连接权重系数z的计算中：

其次，计算突触连接权重系数z，z与y之间的关系满足其中λ，z_l，z_h均为大于零的常数，z_l代表权重系数z 的下限，z_h代表权重系数z的上限，在实际应用中，根据经验，λ＝1,z_l＝0.5,z_h＝ 2。图4为突触连接权重系数相关参数y和突触连接权重系数z的计算流程图。由于PRP浓度是非负的，若Tag和PRP浓度中有一个为0，则y不会发生变化；若Tag为负，PRP浓度为正，则y降低；若Tag为正，PRP浓度为正，则y升高。由y可以计算出z。在该计算流程示例中，τ_y＝1ms,λ＝1,z_l＝0.5,z_h＝2。

最后，将突触连接权重系数z乘以初始连接权重w₀，得到t₁时刻的连接权重w₁＝z×w₀。

Claims

1.一种基于钙浓度的神经突触可塑性计算方法，其特征在于，包括以下五个步骤：

1)钙离子浓度的计算：针对需要计算的突触，根据初始时刻t₀的突触前神经元和突触后神经元的膜电位值，以及突触的初始连接权重w₀，分别计算下一时刻t₁的树突和树突棘中的钙离子浓度，其中w₀是非负常数；

2)确定权重变化方向Flag：根据树突棘中的钙离子浓度，通过和钙离子浓度阈值Ca0_s及Ca1_s对比得到权重变化的方向，其中Ca0_s和Ca1_s均为非负常数，且

Ca1_s＞Ca0_s，Ca0_s＝0.5mM，Ca1_s＝1mM；

3)计算突触状态标签Tag：根据权重变化方向及突触状态变化常数α_T，β_T，LTD和β_T，LTP，计算突触状态标签，其中α_T＝0.01，β_T，LTD＝0.5，β_T，LTP＝0.5；

4)计算可塑性相关蛋白质PRP浓度：根据树突中的钙离子浓度计算可塑性相关蛋白质浓度；

5)计算突触权重的变化：根据突触状态标签Tag和可塑性相关蛋白质PRP浓度，计算突触权重的变化，得到t₁时刻的新权重；

重复以上过程，计算模拟时间内的突触连接强度；

步骤3)具体为：

根据树突棘所处突触权重变化方向及突触状态变化常数α_T，β_T，LTD和β_T，LTP得到t₁时刻该树突棘所处的突触状态标签Tag；如果Tag为0，则该树突棘所处的突触不会发生连接权重变化；如果Tag大于0，则该树突棘所处的突触进入连接权重可以增加的状态；如果Tag小于0，则该树突棘所处的突触进入连接权重减弱的状态；Tag的演化是一个与Flag有关的动力学过程；

如果t₁时刻Flag为0，则其中α_T＝0.01，若Tag以t₀时刻初始值为0开始演化，则Tag关于时间的导数为0，Tag不发生变化；

如果t₁时刻Flag为-1，则其中β_T，LTD＝0.5，若Tag以t₀时刻初始值为0开始演化，则Tag关于时间的导数为负，Tag降低，直至降到Tag关于时间的导数为0时，即/>时停止下降，并保持/>

如果t₁时刻Flag为1，则其中β_T，LTP＝0.5，若Tag以t₀时刻初始值为0开始演化，则Tag关于时间的导数为正，Tag升高，直至升至Tag关于时间的导数为0时，即/>时停止上升，并保持/>

步骤4)具体为：

根据树突中的钙离子浓度[Ca²⁺]_d和钙离子浓度阈值Ca0_d得到t₁时刻树突中PRP的浓度；PRP浓度的演化是与[Ca²⁺]_d有关的动力学过程；其中Ca0_d＝1mM；树突中的钙离子可以引发生成PRP的化学反应，单位时间内PRP的生成率记为PRP_rate；如果t₁时刻树突中的钙离子浓度[Ca²⁺]_d低于阈值Cα0_d，则PRP_rate为0或衰减为0，PRP不发生变化；如果t₁时刻树突中的钙离子浓度[Ca²⁺]_d高于阈值Ca0_d，则PRP_rate按alpha函数的形式变化，PRP作为PRP_rate对时间的积分，先增加，后持平；

步骤5)具体为：

首先，计算突触连接权重系数相关参数y，y的演化是与Tag和PRP浓度有关的动力学过程，其中τ_y＝1ms；由于PRP浓度是非负的，若Tag和PRP浓度中有一个为0，则y对时间的导数为0，y不会发生变化；若Tag为负，PRP浓度为正，则y对时间的导数为负，y降低；若Tag为正，PRP浓度为正，则y对时间的导数为正，y升高；

其次，计算突触连接权重系数z，z与y之间的关系满足其中z_l代表权重系数z的下限，z_h代表权重系数z的上限，λ＝1，z₁＝0.5，z_h＝2；

再次，将突触连接权重系数z乘以初始连接权重w₀，得到t₁时刻的连接权重w₁＝z×w₀；最后，重复步骤1)、步骤2)、步骤3)、步骤4)和步骤5)，计算整个模拟时间内神经元之间的突触连接权重变化。

2.根据权利要求1所述的计算方法，其特征在于，步骤1)具体为：首先，设定突触前神经元和突触后神经元的初始膜电位值以及突触前神经元和突触后神经元之间的初始连接权重w₀；其中突触前神经元和突触后神经元的初始膜电位值设置为-70mV，突触前神经元和突触后神经元之间的初始连接权重w₀取值为1；

其次，采用脉冲神经元模型对突触前神经元和突触后神经元进行建模，输出为膜电位值；

再次，根据突触前神经元的脉冲神经元模型，计算t₁时刻突触前神经元的膜电位；

然后，根据突触动力学模型对突触连接进行建模，模型的输入为t₁时刻突触前神经元膜电位，t₁时刻突触后神经元膜电位，突触动力学时间常数以及突触受体的反转电位，模型的输出为突触连接对突触后神经元产生的电流；

最后，根据突触动力学模型及突触后神经元的神经元模型，计算t₁时刻突触后神经元的膜电位以及突触后神经元树突和树突棘中的钙离子浓度；突触后神经元树突中的钙离子源自位于细胞膜上的钙离子通道，突触后神经元树突棘中的钙离子源自位于树突棘上允许钙离子通过的突触受体；在t₁时刻由钙离子通道进入神经元内的钙离子浓度[Ca²⁺]_Ca-channel由包含该离子通道的脉冲神经元模型计算得出，在t₁时刻由树突棘上允许钙离子通过的突触受体进入神经元内的钙离子浓度[Ca²⁺]_receptor由该突触的动力学模型计算得出；记t₁时刻树突棘和该树突棘所处的树突中的钙离子浓度分别为[Ca²⁺]_s和[Ca²⁺]_d。

3.根据权利要求1所述的计算方法，其特征在于，步骤2)具体为：

根据树突棘中的钙离子浓度[Ca²⁺]_s和钙离子浓度阈值Ca0_s及Ca1_s对比得到t₁时刻该树突棘所处的突触的权重变化方向Flag；Ca0_s＝0.5mM，Ca1_s＝1mM；如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s低于阈值Ca0_s，则Flag为0，该树突棘所处的突触不会发生连接权重变化；如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s处于阈值Ca0_s和阈值Ca1_s之间，则Flag为-1，该树突棘所处的突触连接权重变化方向为减弱；如果t₁时刻树突棘中的钙离子浓度[Ca²⁺]_s高于阈值Ca1_s，则Flag为1，该树突棘所处的突触的连接权重变化方向为增加。