基于皮质神经集群模型的NREM期睡眠节律仿真研究

张瑞 王文华 李强 潘敏

【主持人语】医学数据纷繁复杂，信息量巨大。对各类医学数据进行深度挖掘、建模分析和信息提取，对于进一步理解疾病机制，自动精准诊断和有效治疗有重要意义。本专题集合了医学数据分析的5篇文章，数据类型上既包含了脑电和关节骨振等一维信号，也包含二维甚至三维的冠脉造影（CTA）和磁共振影像数据。从研究方法上看，有适于大数据量自动挖掘的深度学习方法，有针对小样本和先验特征提取的经典机器学习分类方法，也有进行二维形状分析的交互式区域生长法和主动形状模型法。从拟解决的临床问题上看，既包括对癫痫、关节退行性疾病的自动判读，用于心脏CTA影像的冠状动脉自动分割，也包括对睡眠脑电节律的数学建模与仿真，以及可以看到组织生化信息的化学交换饱和转移（CEST）磁共振图像分析工具。期待通过本专题，能激发更多更好的医学数据分析以及学科交叉的研究和应用，服务于我国和全人类的医疗健康事业。

【主持人】宋小磊，西北大学信息科学与技术学院教授，博士生导师，国家千人计划青年基金项目获得者。

摘要：睡眠是人生命中一種重要的生理现象。根据脑电波的不同特征，可以将睡眠分为觉醒期、非快速眼动睡眠期（NREM）和快速眼动睡眠期（REM），其中，NREM期占据整个睡眠过程的主要时长。K-复合波（KCs）和睡眠慢波（SOs）为NREM期的典型节律，其对睡眠具有一定的调节保护作用。该文从计算模型的角度出发，对K-复合波（KCs）和睡眠慢波（SOs）产生的神经生理机制进行解释。首先，考虑皮层兴奋性神经集群与抑制性神经集群之间的相互作用，并融入M型钾电流的相关机制，建立一种改进的皮质神经集群模型;其次，采用鞍结分岔和Hopf分岔理论对模型关键参数进行分析;最后，通过数值实验实现对K复合波和睡眠慢波的仿真模拟。

关键词：：皮质神经集群模型;K-复合波（KCs）;睡眠慢波（SOs）;分岔分析

中图分类号：O29

DOI：10.16152/j.cnki.xdxbzr.2020-04-005开放科学（资源服务）标识码（OSID）：

Simulation of brain rhythms during NREM sleep

based on cortical neural mass model

ZHANG Rui， WANG Wenhua，  LI Qiang， PAN Min

（Medical Big Data Research Center， Northwest University， Xi′an 710127， China）

Abstract： Sleep is an important part of human life. According to the characteristics of EEG signals， sleep is categorized into wakefulness， non-rapid eye movement sleep （NREM） and rapid eye movement sleep （REM）， where the NREM stage is predominant in duration. The typical rhythms in NREM sleep include K-complexes （KCs） and slow oscillations （SOs）， which are thought to be important for the protection of sleep. Therefore， we focus on explaining the neurophysiological mechanisms for the generation of KCs and SOs from the model-based point of view. Firstly， considering the interaction between excitatory and inhibitory cortical clusters， a new cortical neural mass model is developed by incorporating the mechanism of M-type potassium current. Then， using the saddle-node bifurcation and Hopf bifurcation techniques， we perform the bifurcation analysis for several key parameters of the model. Finally， simulations show that the model is able to generate the KCs and SOs.

Key words： cortical neural mass model; K-complexes （KCs）; slow oscillations（SOs）; bifurcation analysis

睡眠约占人类生命1／3左右的时间，是一项非常重要的生命过程［1］。根据脑电节律的不同特征，可将睡眠分为3个时期：觉醒期，非快速眼动睡眠期（non-rapid eye movement sleep，NREM）和快速眼动睡眠期（rapid eye movement sleep，REM），其中，NREM又分为Ⅰ，Ⅱ，Ⅲ期［2］。K-复合波（KCs）和睡眠慢波（SOs）为NREM期的典型节律，其对睡眠具有一定的调节保护作用。本文将从计算模型的角度出发，对K-复合波（KCs）和睡眠慢波（SOs）产生的神经生理机制进行解释。

睡眠慢波主要发生于NREM Ⅲ期（也称深度睡眠期）， 它是一类低频高振幅脑电波，频率为0.5～2 Hz，振幅大于75 μV。 成年人的睡眠慢波仅在睡眠时出现， 但在深度麻醉、 缺氧或大脑损伤时也可出现［3］。 K-复合波主要发生于NREM Ⅱ期， 被认为是一种皮质慢振荡活动［4］， 一般由一个方向向上的负向波连接一个方向向下的正向波组成， 且持续时间≥0.5 s， 最小的峰值不低于100 μV，最大不超过400 μV，通常情况下会和睡眠纺锤波融合在一起，但也可单独出现［5］。有研究认为，在睡眠慢波期，皮层神经元常处于 Up 和 Down 状态的波动之中。其中，Up是激活状态，神经元接收到同步化的兴奋性和抑制性突触传入，兴奋和抑制保持平衡，对应去极化过程;Down是相对沉默的状态，几乎没有突触输入，对应超极化过程。K-复合波对应于Down状态［6］。然而，这些活动状态的发生、切换，以及K-复合波与睡眠慢波产生的内在神经机制，尚不完全清楚，尤其缺乏理论证据。

为了探索大脑电活动产生的内在机制，研究者们基于生物物理基础，提出了一种用于模拟神经系统电活动的宏观数学模型，即神经集群模型（neural mass model，NMM）。NMM通过描述功能相似、结构相近的神经元整体（神经元集群）的某些平均特性来反映其整体放电行为［7］。由于此类模型能够用较为简洁并具有明确生理学含义的变量来反映神经元集群的放电活动，在分析大脑节律的内在机制方面具有重要作用。特别地，NMM也被用于研究睡眠节律产生的内在机制，其中，Weigenand等人通过考虑皮层兴奋性神经元中的钠钾电流，并对其强度随时间的变化关系进行建模，探索了该电流强度的变化对于睡眠典型节律产生的影响［8］。在此基础上，进一步考虑到K-复合波和睡眠慢波的产生与突触输入敏感性和钾离子通道失活等机制相关［1，9］，本文首先通过融合M型钾电流的相关神经生理机制来建立一种改进的皮质神经集群模型;其次，采用鞍结分岔和Hopf分岔理论对模型关键参数进行分析;最后，通过数值实验实现对K-复合波和睡眠慢波的仿真模拟。

1 一种改进的皮质神经集群模型

皮质神经集群模型由皮质上的兴奋性神经集群和抑制性神经集群相互耦合构成［8］， 其结构如图1所示。 其中， 兴奋性神经集群主要接受来自自身的兴奋性反馈以及抑制性神经集群的抑制性输入， 而抑制性神经集群主要接受来自自身的抑制性反馈以及兴奋性神经集群的兴奋性输入。 同时， 两类集群接受来自周围脑皮质区域的外部输入。

在该皮质神经集群模型中，每个集群的状态演变过程主要通过两个算子进行刻画。

2 实验结果及分析

2.1 分岔分析结果

本小节在所建立模型的基础上，采用分岔理论对模型关键参数进行分岔分析，以了解模型对不同参数取值范围的动力学特性，从而确定产生K-复合波和睡眠慢波的参数取值范围。

根据先验知識，选取兴奋性突触增益的倒数σe与M型钾离子通道的电导gM作为分岔参数，并使用XPP软件完成分岔分析。模型其他参数的取值如表1所示，而外部输入n，′n用具有零均值的高斯白噪声函数表示。

1）单参分岔分析

首先，固定兴奋性突触增益σe的不同取值，对电导gM进行单参数分岔分析，结果如图2所示。其中，图2（a），（b），（c）分别对应于σe取值4 mV，4.8 mV，6 mV时的gM分岔分析图。

图2（a）中标号为2的点是Hopf分岔点，标号为3和5的两个点均为鞍结点，其电导gM取值分别为2.839，3.089，1.908。具体分析结果可见表2。

图2（b）中有两个鞍结点3与5以及两个Hopf分岔点2与6，从标号为2的Hopf分岔点出发产生Hopf分岔，在此分岔过程中产生了一个canard爆炸现象，即形成高振幅极限环，并且该极限环经过通过在左侧鞍结同宿轨分岔处消失。标号为2，3，5的点对应的gM值分别为2.897，3.324，3.216。具体分析结果可见表3。

图2（c）中标号为2与3的两个点为Hopf分岔点， 标号为2， 3的点对应的gM值分别为3.692，4.953。具体分析结果可见表4。

2）双参分岔分析

选取兴奋性突触增益的倒数σe与M型钾离子通道的电导gM作为分岔参数，对模型进行双参数分岔分析，得到双参数分岔图如图3所示，分岔曲线主要由两个分支组成，第一个分支是两个鞍结分岔曲线，它们在尖点消失（如图3中的红色曲线），在两条鞍结分岔曲线之间存在3个平衡点，可导致模型呈现双稳态或兴奋态动力学模式，如图2（a）和图2（b）所示。因此，若在此区域中，给予模型一个刺激输入，可引起积极响应，这对于K-复合波和睡眠慢波的产生有重要作用。第二个分支是Hopf分岔曲线（如图3中的蓝色曲线）。在此分岔过程中产生了一个canard爆炸现象，即形成高振幅极限环，并且该极限环通过在左侧鞍结同宿轨分岔处消失。

2.2 数值模拟

首先，选择关键参数σe和gM位于图3区域I内且靠近区域III的边界。在该区域中，模型处于兴奋状态，能够呈现噪声驱动的中等振幅振荡活动，进而诱导K-复合波的产生。图4展示了一个时间长度为30 s的模拟脑电片段。此时σe和gM的取值如表5所示。

其次，选择关键参数σe和gM位于图3区域I内并且靠近区域Ⅱ的边界，在该区域中，模型接近Hopf分岔，能够呈现0.8 Hz的高振幅振荡活动。图5展示了一个时间长度为30 s的模拟脑电片段，此时σe和gM的取值如表5所示。

最后，让关键参数σe和gM分别在一定范围内变化。当σe取值区间为4.55～5，gM取值区间为2.3～2.73时，模型产生K-复合波;当参数取值大于该区间时，模型产生睡眠慢波，即随着σe和gM的取值变化，模型状态由K-复合波转迁至睡眠慢波。由此说明，本文所建立的新的皮质神经集群模型能够模拟NREM Ⅱ期K-复合波到NREM Ⅲ期的睡眠慢波的转迁过程。图6展示了模型关键参数（σe和gM）取值变化对于NREM Ⅱ期K-复合波到NREM Ⅲ期的睡眠慢波的转迁过程的影响。

3 结 语

本文基于Weigenand等人的工作，通过考虑兴奋性神经集群与抑制性神经集群之间的相互作用，并融入M型钾电流的相关生理机制，构建了一种改进的皮质神经集群模型。同时，采用鞍结分岔和Hopf分岔理论对模型关键参数σe和gM进行分析，理解模型在不同参数取值区间所对应的动力学特性，从而进一步确定产生K-复合波和睡眠慢波的参数取值范围。最后，通过数值实验实现对K-复合波和睡眠慢波的仿真模拟。实验结果表明，模型在处于单稳态且接近于Hopf分岔与鞍结分岔时，可诱导K-复合波和睡眠慢波的产生。所得结论为睡眠节律K-复合波和睡眠慢波的产生提供了理论证据。

此外，相关研究表明，丘脑等皮下结构对于参与调控睡眠相关的节律振荡活动具有重要作用［15］。因此，进一步开展针对皮质丘脑等耦合模型的研究工作，对于探索睡眠波及其转迁的相关神经生理机制具有重要价值和意义。

参考文献：

［1］ 高群霞. 基于EEG的睡眠分期与睡眠评估方法研究［D］.广州：华南理工大学， 2015.

［2］ 王菡侨. 有关美国睡眠医学学会睡眠分期的最新判读标准指南解析［J］.诊断学理论与实践，2009，8（6）：575-578.

WANG H Q.Analysis of the latest interpretation standard of sleep stage of American academy of sleep medicine［J］.Journal of Diagnostics Concepts and Practice，2009，8（6）：575-578.

［3］ 肖余粮. 睡眠脑电的特性分析和睡眠分期研究［D］.镇江：江苏大学，2006.

［4］ STERIADE M， AMZICA F. Slow sleep oscillation， rhythmic K-complexes， and their paroxysmal developments［J］. Journal of Sleep Research， 1998， 7（S1）：30-35.

［5］ 彭小虎， 王国锋， 陈欢，等. 睡眠K复合波的功能及其模型［J］.第二军医大学学报， 2014， 35（1）：94-100.

PENG X H， WANG G F， CHEN H， et al. The function and model of sleep K-complex［J］.Academic Journal of Second Military Medical University， 2014， 35（1）：94-100.

［6］ SYDNEY S，HALGREN E， DEHGHANI N，et al. The human K-complex represents an isolated cortical down-state.［J］. Science， 2009，324：1084-1087.

［7］ FABRICE W， PASCAL B， FABRICE B， et al. Computational models of epileptiform activity［J］.Journal of Neuroscience Methods， 2016， 260：233-251.

［8］ WEIGENAND A， SCHELLENBERGER COSTA M， NGO H V V， et al. Characterization of K-complexes and slow wave activity in a neural mass model［J］.Plos Computational Biology， 2014， 10（11）：e1003923.

［9］ STERIADE M， NUNEZ A， AMZICA F. A novel slow （< 1 Hz） oscillation of neocortical neurons in vivo：Depolarizing and hyperpolarizing components［J］. Journal of Neuroscience， 1993， 13（8）：3252-3265.

［10］HASSCLMO M E， GIOCOMO L M. Cholinergic modulation of cortical function［J］. Journal of Molecular Neuroscience， 2006， 30（1/2）：133-135.

［11］ROBINSON P A， RENNIE C J， WRIGHT J J. Propagation and stability of waves of electrical activity in the cerebral cortex［J］. Physical Review E， 1997， 56（1）：826-840.

［12］TUCKWELL， HENRY C. Introduction to Theoretical Neurobiology［M］.Cambridge：Cambridge University Press， 1988.

［13］LILEY D T J， CADUSCH P J， WRIGHT J J. A continuum theory of electro-cortical activity［J］. Neurocomputing， 1999， 26-27：795-800.

［14］壽天德. 神经生物学［M］. 2版.北京：高等教育出版社， 2006：36.

［15］RICH S， ZOCHOWSKI M. Effects of neuromodulation on excitatory-inhibitory neural network dynamics depend on network connectivity structure［J］.Journal of Nonlinear Science，2018（8）：1-24.

（编 辑 李 静）