磁流作用下网络神经元中簇放电的模式

The Pattern of Bursting in Network Neurons under Magnetic

摘要: 呼吸节律是哺乳动物重要的生理性活动，其产生机制一直受到广泛的关注。关于正常呼吸节律的形成机制，目前主要为起步细胞学说和神经网络学说。其中，神经元网络学说认为依赖于呼吸神经元之间的复杂相互作用。本文，我们主要研究磁流对吸气周期和呼气周期的影响，当呼吸网络中考虑磁流的影响时，网络中的呼吸节律将会表现出许多特性。由于哺乳动物的脑干呼吸中枢模式发生器(CPG)具有许多对生存至关重要的功能，而pre-Bötzinger复合体作为一种CPG，是呼吸节律产生的关键部位。本文研究在磁流作用下，pre-Bötzinger复合体网络系统的节律转迁，并从动力学的角度研究节律转迁的机制。结果表明磁流反馈系数k₀对呼吸网络系统的节律模式和周期都具有重要影响。结果对呼吸节律的研究提供有益的参考。

Abstract: Respiratory rhythm is an important physiological activity in mammals, and its production mecha-nism has been widely concerned. Regarding the normal respiratory rhythm formation mechanism, the main theory is the pacemaker theory and the network theory. Among them, the theory of neu-ronal networks believes that it relies on complex interactions between respiratory neurons. In this paper, the effects of magnetic flow are considered in the respiratory network, the respiratory rhythm in the network will exhibit many properties. As the brainstem respiratory central pattern generator (CPG) in mammals has many functions that are essential for survival, and the pre-Bötzinger complex, as one type of CPG, is a key site of respiratory rhythm generation. In this paper, we focus on the effect of magnetic currents on the inspiratory and expiratory cycles, and when the effect of magnetic currents is taken into account in a respiratory network, the rhythmic transition of the pre-Bötzinger complex network system under the action of magnetic currents is studied, and the mechanism of the rhythmic transition is investigated from a kinetic point of view. The results suggest that magnetic flow feedback coefficient k₀ has an important influence on both the rhythmic pattern and the period of the respiratory network system. The results provide a bene-ficial reference for the study of respiratory rhythm in neural system.

1. 引言

呼吸节律是由位于延髓中起步型神经元的节律性兴奋引起的，已有的研究成果表明，前包钦格复合体是呼吸节律活动产生最主要的部位 [1] [2] 。Lieske等人发现脑干切片中的pre-Bötzinger复合体可以产生两种与正常呼吸极为相似的节律，该节律是多种呼吸模式的基础 [3] 。最近研究表明，在适当的外部电流刺激下，神经元的簇放电节律会发生改变，例如，磁流可以直接影响神经元的动态特性 [4] [5] [6] ，这种节律的控制模式对治疗节律相关的精神类疾病有一定的参考价值 [7] [8] [9] 。马军等人发现，通过电荷和磁通可以分别描述电场与电磁感应的变化，并以此为基础提出了几个神经元控制模型 [10] [11] 。

对于许多物种的脑干和脊髓中枢模式发生器(CPG)中的神经元网络，也都会表现出有规律有节奏的活动。其中，CPG是指产生节律运动活动的神经环路，在许多中枢模式发生器中，这些模式由不同神经元群的顺序激活组成，而这些神经元群通过突触连接相互作用 [12] 。虽然有许多研究已经探索了呼吸系统中pre-Bötzinger复合体神经元放电模式的产生机制 [13] [14] ，并引入磁流和电流条件下对pre-Bötzinger复合体放电模式的产生进行了讨论，但仍有一些关键问题值得进一步研究，如：可放电的神经元嵌入到通过突触相互作用的完整网络中时，对引入磁流之后的放电模式会产生怎样的影响？若将pre-Bötzinger神经元直接嵌入到网络中，它的放电模式是否会影响这些神经元嵌入后的呼吸网络中的放电机制？这些问题还需要进一步探究。呼吸节律的鲁棒性是指在各种条件下保持输出模式，或者在表示这些条件的模型参数或特征的变化之间保持输出模式。而且可以出现一种更强的鲁棒性的形式，其中不仅保留输出模式，而且还会保留产生该模式的动力机制。根据已有的实验表明，在甲壳类动物中，CPG的节律性具有显著的鲁棒性，它会对温度等外在因素以及离子电流等内在因素的变化产生一定的消化节律 [15] [16] [17] 。

本文中，我们运用呼吸CPG的高度简化模型，由代表重要神经元群体 [18] 的少量单位组成。在该模型网络中会产生三阶段呼吸节律，包括吸气神经元，吸气后神经元和呼气神经元 [19] [20] [21] 。我们引入磁流，对呼吸网络中具有放电模式的神经元进行建模，该微分方程支持已建立的簇放电产生机制，而且已在这些神经元中得到很好的实验记录 [22] [23] [24] 。我们发现网络的三阶段呼吸节律在网络中会受到磁流的严重影响。虽然外部刺激对神经元中混合簇放电的产生已有一定的研究，但是电流和磁流等外激励对混合簇放电的影响，其现象和动力学机制都极其复杂，值得进一步探究。

本文结构如下：第二部分给出了pre-Bötzinger复合体兴奋性神经元和呼吸CPG网络的简化常微分方程模型；第三部分，研究磁流调控对神经元放电模式的影响，通过相平面分析，给出了网络系统的节律的产生和转迁的动力学机制；第四部分给出结论和展望。本文中的相平面分析是通过XPPAUT [24] 软件计算得到。

2. 模型介绍

pre-Bötzinger复合体兴奋性神经元和呼吸CPG网络的简化常微分方程模型如方程(1)~(6)。该呼吸CPG网络包括吸气前神经元(pre-I)，早期吸气神经元(early-I)，后期吸气神经元(post-I)，和增强型呼气神经元(aug-E) [20] [22] 。在这个网络中，我们假设每个神经元都能在静息与放电状态之间正常转换。对于编号为第一个神经元的pre-I，引入磁流后的神经元模型的微分方程表达式如下：

$\text{d}{V}_1/\text{d}t=\left(-I_{NaP}\left(V_1,h_1\right)-I_K\left(V_1,n_1\right)-I_L\left(V_1\right)-I_{syn-I}\left(V_1\right)-I_{syn-E}\left(V_1\right)+I_{FBC}\right)/C,$ (1)

$\text{d}{h}_1/\text{d}t=\left(h_{\infty }\left(V_1\right)-h_1\right)/{\tau }_h\left(V_1\right),$ (2)

其中，兴奋性神经元泄漏电流的方程表达式是 $I_L\left(V_1\right)=g_{L,exc}\cdot \left(V_1-E_{Na}\right)$ ，持续钠电流为 $I_{NaP}\left(V_1,h_1\right)=g_{NaP,exc}\cdot m_{\infty }\left(V_1\right)\cdot h_1\cdot \left(V_1-E_{Na}\right)$ ，钾电流为 $I_K\left(V_1,n_1\right)=g_{K,exc}\cdot n_{\infty }^4\left(V_i\right)\cdot \left(V_i-E_K\right)$ ；而模型中的其他神经元分别编号为：第二个神经元为early-I，第三个神经元为 post-I，第四个神经元为aug-E，且具体微分方程表达式如下：

$\text{d}{V}_i/\text{d}t=\left(-I_{NaP,i}\left(V_i,h_i\right)-I_{K,i}\left(V_i,n_i\right)-I_{L,i}\left(V_i\right)-I_{syn-I,i}\left(V_i\right)-I_{syn-E,i}\left(V_i\right)\right)/C,$ (3)

$\text{d}{h}_i/\text{d}t=\left(h_{\infty }\left(V_i\right)-h_i\right)/{\tau }_h\left(V_i\right),$ (4)

其中 $i\in \left\{2,3,4\right\}$ 。 $h_i$ 和 $n_i$ 分别对应持续钠电导和钾电导的门控变量，抑制性神经元泄漏电流的方程表达式为 $I_{L,i}\left(V_i\right)=g_{L,inh}\cdot \left(V_i-E_{Na}\right)$ ，各单元的持续钠电流为 $I_{NaP,i}\left(V_i,h_i\right)=g_{NaP,inh}\cdot m_{\infty }\left(V_i\right)\cdot h_i\cdot \left(V_i-E_{Na}\right)$ ，钾电流为 $I_{K,i}\left(V_i,n_i\right)=g_{K,inh}\cdot p_i\cdot \left(V_i-E_K\right)$ ，对于抑制性单元，具有专用的门控变量p_i由电压相关的激活函数 $\text{d}{p}_i/\text{d}t=\left(d_i{f}_{out,inh}\left(V_i\right)-p_i\right)/{\tau }_{p,i}$ 表示。

此外，动作电位跨此网络所有神经元时需要通过突触传递信号，抑制性突触电流 $I_{syn-E,i}\left(V_i,h_i\right)$ 和兴奋性突触电流 $I_{syn-I,i}\left(V_i,h_i\right)$ 的方程如下：

$I_{syn-I,i}\left(V_i\right)=g_{syn-I}\cdot \left(V_i-E_{syn-I}\right)\cdot {\displaystyle {\sum }_{j\in N_i^{inh}}{b}_{ji}}{f}_{out,inh}\left(V_j\right),$ (5)

$I_{syn-E,i}\left(V_i\right)=g_{syn-E}\cdot \left(V_i-E_{syn-E}\right)\cdot \left[{\displaystyle {\sum }_{j\in N_i^{exc}}{a}_{ji}}{f}_{out,exc}\left(V_j\right)+{\displaystyle {\sum }_{i=1}^4{c}_{1i}}+{\displaystyle {\sum }_{i=1}^4{c}_{2i}}\right],$ (6)

$N_i^{exc}$ 和 $N_i^{inh}$ 分别表示突触到神经元单元i的兴奋性和抑制性单元集，对于pre-I兴奋性单元电流表现为；其中，常数b_ji，a_ji分别表示抑制性和兴奋性连接的强度，而常数c_ji，是对每个单元的驱动强度的参数。磁流激励的表达式是： $I_{FBC}=k_0\left(\alpha +3\beta {\phi }^2\right)\left(v_e-v\right)$ ，k₀为方程中的反馈系数，模型中涉及的其他函数和参数见附录。

3. 主要结果

3.1. I_FBC对呼吸节律周期的影响

我们探索了磁流I_FBC对整个网络节律周期(T)、吸气阶段的持续时间(T_I)和呼气阶段的持续时间(T_E)的影响，即磁流强度的反馈系数k₀对pre-I单元输出周期的影响。由突触连接的呼吸网络放电模式如图1所示，黑色曲线表示pre-I，红色曲线表示early-I，蓝色曲线表示post-I，粉色曲线表示aug-E。这些网络放电模式在图1(a)~(c)显示了整个网络节律周期(T)、吸气阶段的持续时间(T_I)和呼气阶段的持续时间(T_E)，对于组成网络的pre-I，early-I，post-I和aug-E各单元的输出，在pre-I阶段和early-I阶段的活动水平突然剧增之前，pre-I阶段和early-I阶段的输出都会缓慢上升，这与吸气阶段的出现相对应。Pre-I和early-I然后逐渐下降，使得post-I单元放电活动激活，并对这两个单元进行相应的抑制，标志着此时已经切换到呼气期。在这种转换之后，post-I输出下降，对于aug-E的抑制减少，使其输出上升，产生所谓的呼气期，aug-E输出并抑制但不完全抑制post-I。

随着磁流的反馈系数k₀的增加，整个网络节律周期(T)以及吸气(T_I)和呼气(T_E)阶段输出的周期会发生相应的改变。接下来，我们探索反馈系数k₀对包含pre-I单元的的网络中节律的频率以及周期的影响。当k₀分别是0，0.5，1时，pre-I神经元输出的整个网络节律周期(T)与呼气周期(T_E)的持续时间变小，而吸气周期(T_I)的持续时间不变，由突触连接的呼吸网络的簇放电个数随k₀的增加而增加。

Figure 1. Time series of membrane potential V for the entire network rhythm at a variation of $k_0$ . (a) $k_0=0$ ; (b) $k_0=0.5$ ; (c) $k_0=1$

图1. $k_0$ 变化时整个网络节律的膜电位V时间序列。(a) $k_0=0$ ；(b) $k_0=0.5$ ；(c) $k_0=1$

当磁流的反馈系数k₀继续增加时，发现整个网络节律又出现了新的现象，周期(T)以及吸气(T_I)和呼气(T_E)阶段输出的周期与上述时间序列有所不同，如图2所示。当反馈系数k₀分别是1.2，1.22，1.5时，pre-I神经元输出的整个网络节律周期(T)与吸气周期(T_I)的持续时间变长，而的呼气周期(T_E)持续时间不变，由突触连接的呼吸网络的簇放电个数随k₀的增加呈现先增加后减少的变换趋势。

Figure 2. Time series of membrane potential V for the entire network rhythm at a variation of $k_0$ . (a) $k_0=1.2$ ; (b) $k_0=1.22$ ; (c) $k_0=1.25$

图2. $k_0$ 变化时整个网络节律的膜电位V时间序列。(a) $k_0=1.2$ ；(b) $k_0=1.22$ ；(c) $k_0=1.25$

整个网络节律周期T (黑色)、吸气持续时间T_I (蓝色)和呼气持续时间T_E (红色)与在网络节律内反馈系数k₀之间的关系，如图3所示。随着磁流的反馈系数k₀的增加，整个网络节律先是通过调节呼气T_E来控制整个网络的呼吸周期。结果表明，由于磁流在呼吸节律性中起到了抑制作用，即周期T随着k₀的增加而减少。当pre-Bötzinger复合体受到磁流刺激作用后，T_E持续时间缩短，从最大值约3000 ms到最小值约1000 ms，如图3所示；而pre-I单元随着磁流的反馈系数k₀从1.0增加到1.25时，T_E持续时间处于平缓状态，而T_I持续时间开始急剧增加，此时磁流刺激是通过影响吸气T_I的持续时间来控制网络节律周期T的呼吸周期，增加k₀延长吸气时间，导致呼吸节律的周期增加，如图3所示。

Figure 3. The feedback function k₀ varies through the network rhythm period (T) as well as the durations of the inspiratory phase (T_I) and the expiratory phase (T_E)

图3. 反馈系数k₀变化时整个网络节律周期(T)以及吸气阶段的持续时间(T_I)和呼气阶段的持续时间(T_E)

3.2. I_FBC对网络放电节律影响的动力学分析

由前一个神经元释放信号，通过两个神经元接触的结构，后一个神经元的树突前端的轴突接收信号。因此，每个单元的膜电位方程取决于另一个单元的输出。在相平面中可以用可视化的动态模型单元来模拟网络中的每个神经元，用几何形式来表示网络的振荡机制。考虑投影到各二维相平面的网络输出，在这些相平面内，两个变量h_i，V_i的对应的零等值线，如图4所示。其中，此模型输出的pre-I单元的(h₁, V₁)相平面如图4(a)所示，全系统的轨线(黑色曲线)也叠在此相平面上，蓝色实线是指V₁零等值线，右下方的零等值线对应于对pre-I的没有抑制作用的情况，而最左侧的零等值线对应于pre-I接收的最大抑制时的情况，而中间的零等值线是pre-I单元接收到较低的抑制水平，此时pre-I单元开始转至活动阶段。其中，黑色实心圆和方形分别对应pre-I单元抑制最低值和最高值的点。此外，红色曲线为h₁零等值线，与抑制水平无关。

Figure 4. Phase plane analysis of the neurons in the network at feedback coefficient $k_0=0$ . (a) (h₁, V₁); (b) (h₃, V₃); (c) (h4, V4)

图4. 反馈系数 $k_0=0$ 时网络中各单元的相平面分析。(a) (h₁, V₁)平面；(b) (h₃, V₃)平面；(c) (h₄, V₄)平面

另外，对于post-I单元的(h₃, V₃)相平面上，最右侧的V_pre-I零等值线对应于post-I的最小抑制，如图4(a)中所示。图4(b)中的曲线与图4(a)类似的，而沿轨线获得post-I单元的零等值线V₃最大值时(黑色方块)，此时恰好与图4(a)中的时间点相同。接下来，当对post-I单元达到最左侧零等值线数值，即达到对post-I单元的最大抑制时，pre-I也开始转变为低抑制水平(与图4(a)中相同时间点用黑实心圈标记)，如图4(b)所示。

网络中全系统的轨线(黑色曲线)也叠在aug-E单元的相平面(h₄, V₄)上，如图4(c)所示，从右到左，V₄ 零等值线依次对应于对aug-E的最小抑制以及到pre-I单元开始被激活时时对aug-E的抑制，此时，也与其他图相同时间点时用黑实心圈表示。其中，最右侧的零等值线V₁表示为达到峰值时对aug-E的抑制的曲线，而最左侧的零等值线表示为V₃达到峰值时(用黑色方块标记)对aug-E的抑制的曲线值。

通过对网络呼吸节律中的周期内的几个关键时刻绘的相平面分析，对网络放电模式作出了解释。这些轨迹投影通常位于恰当的电压零等值线上。而我们在讨论中省略了early-I单元，这是由于early-I单元缺乏被激活向吸气阶段过渡时所需的自我调节能力。

随着磁流的反馈系数k₀的增加，呼吸网络中的pre-I单元出现不同的放电状态，如图5所示，因此，也出现了非常丰富的神经元网络的动力学机制。不同的磁流度使得网络各单元都会出现放电振荡机制。与图4(a)具有相似的曲线，右下方的零等值线对应于对pre-I的没有抑制作用的情况，而最左侧的零等值线对应于pre-I接收的最大抑制时的情况，而在pre-I单元的相平面中间的零等值线通过h_1I与V₁交点后，是其pre-I单元接收到较小的抑制水平，让其可以开始变成到活动阶段。引入磁流激励作用后，pre-I单元被激活时的零等值线V_pre-I的立方形状呈现单调的状态。Pre-I单元的开始振荡唯一可能起作用的时刻是在被激活开始的时候，这取决于零等值线h₁与V₁交点处。然而，不同反馈系数k₀却使得pre-I单元具有相似的相平面，如图5所示。在磁流作用下，零等值线V₁在激活开始的时刻时是呈现单调状态的，而转变为呼气活动的机制与网络呼吸周期的持续时间对应(图5(a)~(e))。也就是说，即使在缺乏抑制的情况下，零等值线V_pre-I在这一变化中呈现单调状态，使足够强的抑制可以产生较长的静息状态。当对pre-I单元的抑制在一段时间后开始衰减时，I_Nap的失活可以允许pre-I单元不受抑制，这种机制是构成整个网络的呼吸节律的。

当pre-Bötzinger复合体嵌入网络后振荡时，整个网络周期T随着k₀的增加出现两种变化，一种是网络周期T缩短到一定值然后保持稳定，从初始状态4.2秒降到2.1秒左右，然后再增加到6秒；T_E缩短到一定值然后保持稳定，从最大值3.2秒降到1.3秒左右保持稳定；另一种变化是，随着磁流的反馈系数k₀的增加，T_I保持稳定值，当k₀增加到一定值时，T_I突然增大，即T_I从保持稳定1.8秒上升到4.5秒左右。

Figure 5. Phase plane analysis of pre-I neurons in networks. (a) $k_0=0.5$ ; (b) $k_0=1$ ; (c) $k_0=1.2$ ; (d) $k_0=1.22$ ; (e) $k_0=1.25$

图5. 网络中pre-I单元的相平面分析。(a) $k_0=0.5$ ；(b) $k_0=1$ ；(c) $k_0=1.2$ ；(d) $k_0=1.22$ ；(e) $k_0=1.25$

4. 结论

中枢模式发生器(CPG)对生存相关的节律性运动具有许多至关重要的功能，作为一种CPG，pre-Bötzinger复合体对哺乳动物的呼吸节律的输出，特别是在各种代谢和环境条件下持续存在是至关重要的。因此，即使CPG网络中的特定部分表现出输出驱动和神经自我调节的内在节律性，它们的放电模式也可能随着这些水平的变化而改变。所以，内在动力学变化下网络输出和调节的不变性可能是在各种条件下保持CPG稳定性能的关键特征。在呼吸网络的情况下，这种不变性确保了稳定地输出动态节律，这也是如此重要的稳态CPG的关键属性。

本文主要研究了pre-Bötzinger复合体作为CPG的一类，在磁流的作用下网络系统的节律改变及其动力学研究。本文中显示了磁流的反馈系数k₀在保持网络节律性的范围内变化时，周期、吸气持续时间和呼气持续时间与磁流的反馈系数k₀的值之间的关系，即随着磁流的反馈系数k₀的增加，pre-I神经元输出的整个网络节律周期(T)先缩短到一定值后快速增加，呼气周期(T_I)持续时间呈现出由稳定转为快速增加，而呼气周期(T_E)的持续时间减少到一定值后趋于平缓，而且由突触连接的呼吸网络中簇的个数随k₀的增大而增加。

附录

模型中涉及的方程：

$\begin{array}{l}{h}_{\infty }\left(V\right)=1/\left(1+\exp \left[\left(V-{\theta }_h\right)/{\sigma }_h\right]\right),{\tau }_h\left(V\right)=\epsilon /\cosh \left[\left(V-{\theta }_h\right)/\left(2{\sigma }_h\right)\right],\\ m_{\infty }\left(V\right)=1/\left(1+\exp \left[\left(V_i-{\theta }_m\right)/{\sigma }_m\right]\right),n_{\infty }\left(V\right)=1/\left(1+\exp \left[\left(V_i-{\theta }_n\right)/{\sigma }_n\right]\right),\end{array}$

$f_{out,inh}=1/\left(1+\exp \left[\left(V-{\theta }_{out,inh}\right)/{\sigma }_{out,inh}\right]\right),f_{out,pre-I}=1/\left(1+\exp \left[\left(V-{\theta }_{out,pre-I}\right)/{\sigma }_{out,pre-I}\right]\right).$

Table 1. Parameter values used in the model.

表1. 模型中的参数值

NOTES

^*通讯作者。