人类的大脑被认为是宇宙间最为复杂的物体，试图理解正常或者功能紊乱大脑的内部连接模式是现代科学研究中最具挑战的领域。19世纪末期，研究者已经认识到大脑内部的神经元构成了一个强大而又复杂的网络结构^[1-2]，目前研究中脑网络包括结构网络和功能网络，结构网络描述大脑的解剖连接状态，包括不同区域结构形态上的相关和真实的解剖连接；功能网络主要用来描述大脑不同区域结构之间的功能连接和有效的连接^[2]。

功能网络构建的测量手段主要包括脑电图(electroencephalography，EEG)、脑磁图(magnetoencephalogram，MEG)和功能磁共振成像(functional magenetic resonance imaging，fMRI)等^{[1, 3]}。脑功能网络构建的第一步就是测量不同脑区活动间的相关性，目前普遍采用线性技术衡量记录信号之间的关系，包括皮尔森相关、偏相关(时域)及偏相干(频域)等技术^[2-3]。这些线性方法的应用存在一个局限性，它们只能较好地捕捉信号之间的线性关系，却不能很好地检测两个动态系统之间非线性的关系。由于大脑自身的特征，研究者使用非线性的研究方法挖掘大脑信号之间的关系，包括同步和互信息技术^[4-5]。

同步(synchronization)通常被用来表示两个或多个系统随时间变化的相对关系。研究者发现当混沌系统出现耦合时，将伴随出现许多同步现象，包括完全同步(identical synchronization)、相位同步、滞后同步、广义同步以及间歇延迟同步^[5-6]。其中广义同步是研究最广的一种同步方式，它是指在两个动态系统X和Y间，当系统Y的状态响应是通过系统X的函数F来进行驱动时，当函数F连续时，X_i和X_j是系统X中的两个非常接近的点，那么相应的点Y_i和Y_j也是系统Y中两个非常接近的点^[4-5]。广义同步的概念被提出后，相应地出现了基于广义同步理论的几个算法来检测实验数据的相互关系。随着研究的深入，研究者们对这些算法进行不断的改进和完善。文献[4]提出了似然同步算法(synchronization likelihood)，该算法对两个或多个同时记录的时间序列信号之间的相互依存关系给出了一个简单的标准化估计；相比其他算法，似然同步算法更适合于分析脑电等非平稳信号^{[4, 6]}。

本文基于似然同步这种非线性的衡量方法来测量不同区域(不同电极间)间的同步关系，为了检测这种方法应用于非线性系统的表现及效果，首先利用一个混沌系统(是由两个单向耦合Henon映射所构成)对该算法进行实现并仿真。此外，验证正确性后将该算法应用到实际的静息态下的脑电信号。通过上述算法分析静息态下睁眼和闭眼两种不同条件下的脑电信号同步性，从而找出两种不同静息态下的脑电信号同步性上存在的差异，为后续脑电的研究提供一定的参照。

1 似然同步原理

似然同步算法的程序实现的理论依据如下：假设有M个同一时间段记录的信号序列x_{k, i}，其中k表示信道数量(k=1, 2, …, M)，i表示信号记录的时间点(i=1, 2, …, N)，对每一道信号利用时间延迟重建内嵌向量为：

$ {\mathit{\boldsymbol{X}}_{k, i}} = ({x_{k, i}}, {x_{k, i + l}}, {x_{k, i + 2l}}, \cdots , {x_{k, i + (m - 1)l}}) $

(1)

式中，l和m分别表示时间延迟和嵌入的维数。定义任意两个内嵌向量X_{k, i}和X_{k, j}之间的欧氏距离小于ε的概率为：

$ P_{k, i}^\varepsilon = \frac{1}{{2({w_2} - {w_1})}}\sum\limits_{\begin{array}{*{20}{c}} {j = 1}\\ {{w_1} < |i - j| < {w_2}} \end{array}}^N {\theta (\varepsilon - |{\mathit{\boldsymbol{X}}_{k, i}} - {\mathit{\boldsymbol{X}}_{k, j}}|)} $

(2)

式中，|·|操作符是求两个向量之间的欧氏距离；θ(x)是Heaviside阶梯函数，即满足：

$ \theta (x) = \left\{ {\begin{array}{*{20}{c}} 0&{x \le 0}\\ 1&{x > 0} \end{array}} \right. $

(3)

式(2)中w₁和w₂是两个窗，w₁是泰勒修正参数，用来消除由于自相关特性对向量间距离过近造成的影响；w₂窗用来提高时间分辨率，避免信号向量之间的时间间隔太长。当前的研究中，对于窗口w₁和w₂的取值需要满足条件：${w_1} \ll {w_2} \ll N$(≪表示远小于)，对于EEG数据的分析一般取w₁=100，w₂取时间序列长度的1/10，对式(1)中的时间延迟参数l和嵌入的维数m尚没有唯一的取值标准，本文按照研究经验取l=10，m=10。综上可知，式(2)中对于每个k和每个时间点i对应的距离参数ε_{k, i}由概率P_{k, i}^{ε_{k, i}}所决定，通常研究中给定概率$P_{k, i}^{{\varepsilon _{k, i}}} = {p_{{\rm{ref}}}} = 0.05$，假定M道信号中与内嵌向量X_{k, i}之间距离小于ε_{k, i}，并满足窗口条件(w₁ < |i-j| < w₂)的内嵌向量X_{k, j}有H_{i, j}个，即有：

$ {H_{i, j}} = \sum\limits_{k = 1}^M {\theta ({\varepsilon _{k, i}} - \left| {{\mathit{\boldsymbol{X}}_{k, i}} - {\mathit{\boldsymbol{X}}_{k, j}}} \right|)} $

(4)

由实际意义可知H_{i, j}的值介于0和M之间，表示i和j时刻非常接近或者同步的内嵌向量的数量。综上信息，对于k道信号和另外M-1道信号i和j时刻的似然同步值S_{k, i, j}定义如下：

$ \left| {{\mathit{\boldsymbol{X}}_{k, i}} - {\mathit{\boldsymbol{X}}_{k, j}}} \right| < {\varepsilon _{k, i}}{\rm{ : }}{S_{k, i, j}} = \frac{{{H_{i, j}} - 1}}{{M - 1}} $

(5)

$ \left| {{\mathit{\boldsymbol{X}}_{k, i}} - {\mathit{\boldsymbol{X}}_{k, j}}} \right| \ge {\varepsilon _{k, i}}{\rm{ : }}{S_{k, i, j}} = 0 $

(6)

最后对所有时间点j求平均，得到k道信号和其他M-1道信号在i时刻的似然同步值为：

$ {S_{k, i}} = \frac{1}{{2({w_2} - {w_1})}}\sum\limits_{\begin{array}{*{20}{c}} {j = 1}\\ {{w_1} < \left| {i - j} \right| < {w_2}} \end{array}}^N {{S_{k, i, j}}} $

(7)

由式(4)可知此处的似然同步值S_{k, i}表示着当前的k道信号和所有其他道信号的同步强度，S_{k, i}的值介于0~1之间，0表示两者信号完全不同步，1表示两者信号完全同步，通常0和1是不会真正出现在实际信号之间的^[4-5]。对于似然同步方法在脑电EEG数据中应用时，只需要计算电极间信号的同步强度，对上述算法在实际应用时需要进行改动式(4)~式(5)，即只考虑输入的两道信号之间的关系，而无需遍历所有M-1道信号。

2 数据仿真及静息脑电的分析 2.1 似然同步法的实现及仿真效果

为检测本文所使用程序的正确性以及似然同步方法用于非线性数据时的表现，利用一个混沌系统对程序进行仿真。该系统是两个单向耦合Henon映射^[4-5]，包括一个驱动系统X和一个响应系统Y，并且满足Y=F(X)，具体方程如下：

$ \left\{ \begin{array}{l} {x_{i + 1}} = 1.4 - x_i^2 + 0.3{u_i}{\rm{ }}\\ {u_{i + 1}} = {x_i}\\ {y_{i + 1}} = 1.4 - (C{x_i} + (1 - C){y_i}){y_i} + B{v_i}\\ {v_{i + 1}} = {y_i} \end{array} \right. $

(8)

该系统中两者耦合的强度是可变的，由耦合参数C来控制，参数C取值0(无耦合)到1(完全耦合)，式中参数B有两种取值，当B=0.3时，式(8)中的动态系统X和Y是相同的；当B≠0.3时，两者不同。式(8)中的参数x, u, y, v的初始值各取一个0~1之间的随机数，算法中其他参数取值：时间延迟参数l=1，嵌入维数m=10，窗w₁=100，w₂=410，概率参数p_ref=0.05，仿真中将式(8)迭代5 000次，并去掉前1 000个数据点，使用参数B=0.3和B=0.1两种情况比较似然同步结果，使用MATLAB对该系统进行仿真，主要关注似然同步算法的两个方面：

1) 似然同步算法是否能正确表示系统之间的同步性，即当两系统随着耦合系数C的增加，无论B处于哪种情况，两系统的似然同步值是否也相应增大；

2) 检测本文使用的似然同步算法的正确性。仿真过程中耦合参数C由0以0.1的步长递增到1，两个系统X和Y间的似然同步值如图 1所示。由图 1可知，无论是两个相同系统(B=0.3)或是两个不同系统(B=0.1)间的似然同步值，均随着系统间耦合强度的增加而增加。图 2给出了B=0.3时不同耦合系数情况下X和Y信号在1 000~1 100次迭代变化情况。

当系统的耦合系数较低(C=0.4)时，图 2显示此时X和Y信号间幅度的差异较大，对应图 1中两信号的似然同步值约为0.15，耦合系数较高(C=0.8)时，X和Y信号的幅度更为相似，同时两者的似然同步值也接近最大值1。仿真中使用的所有条件出自参考文献[4-5]，仿真的结果验证了本文的似然同步算法程序的正确性，确保本文后期使用该算法对EEG数据分析的有效性。

2.2 真实脑电数据仿真结果

为了把改算法应用到真实的脑电数据，本文截取一段真实的脑电数据，经过滤波提出alpha波成分(8~13 Hz)，再加入不同程度的高斯噪声信号，使得信噪比(SNR)分别为0、5、10、20、30 dB。图 3显示原始信号和5段仿真信号，横坐标表示时间过程，纵坐标表示幅度大小，两者都属于无量纲变量。由图 4可以看出随着信噪比的增加，仿真信号与原始信号越来越相似，因而可利用似然同步方法衡量原始信号和加入噪声后信号的同步性，比较不同信噪比条件下信号之间的同步情况。

仿真中原始脑电信号和加入不同噪声的脑电信号的似然同步值的变化情况如图 4所示。

由于实际脑电测量都会不同程度地引入各类噪声，因此在仿真中加入了不同程度的噪声数据。仿真结果显示，似然同步方法可以较好地衡量任务态EEG信号之间的同步性。

2.3 静息态闭眼与睁眼的脑电数据分析

为了研究正常人在闭眼和睁眼状态下大脑同步性上存在的差异，本文采集了21个健康被试在静息态下的脑电数据，实验范式包含闭眼和睁眼两部分。脑电采集系统采用的是129导的EGI脑电采集系统，采样率是500 Hz，共采集了3 min闭眼和3 min睁眼的脑电数据。经过一些预处理步骤(去坏导、去漂移、平均参考、滤波等)后得到两种静息状态下alpha波(8~13 Hz)的脑电数据。利用似然同步法对预处理后的脑电数据进行每导信号与所有其他128导信号的同步强度计算，平均后得到一个同步值。

图 5显示了两种状态下的分析结果，横坐标表示129导电极的编号，纵坐标表示似然同步值，两者都属于无量纲变量。从图可知，从闭眼到睁眼状态下，大脑的alpha波在几乎所有电极上的同步性都有一个显著的减小，这个发现也与之前EEG文献报道的从闭眼到睁眼状态alpha波去同步化的结果是一致的^[7-9]。

3 结束语

本文结合混沌系统(包含两个单向耦合Henon映射)对似然同步算法进行了仿真实验，并利用真实脑电数据加入噪声方式对似然同步方法在脑电分析应用中的表现和效果进行了测试，测试结果很好地证实了该方法在脑电应用上的有效性。最后将此算法用于分析静息态闭眼和睁眼两种状态下脑电数据的同步性的比较。主要内容及结果如下：

1) 结合混沌系统对似然同步法的仿真

仿真结果发现：①似然同步算法能正确表示系统之间的同步性，即当两系统随着耦合系数C的增加，无论B处于哪种情况，两系统的似然同步值也相应增大。②当耦合系数较低(C=0.4)时，X和Y信号间幅度的差异较大，对应两信号的似然同步值较低(约0.15)；当耦合系数较高(C=0.8)时，X和Y信号的幅度更为相似，同时两者的似然同步值也接近最大值1。这个检测结果表明了似然同步算法的正确性。

2) 真实脑电数据仿真结果

在真实脑电数据加入噪声来对似然同步法在脑电中的应用效果进行测试，结果发现随着信噪比的增加，仿真信号与原始信号越来越相似，因而可利用似然同步方法衡量原始信号和加入噪声后信号的同步性，比较不同信噪比条件下信号之间的同步情况，仿真结果显示似然同步方法可以较好地衡量任务态EEG信号之间的同步性。也很好地证实了该方法在脑电分析中的有效性。

3) 静息态闭眼与睁眼的脑电数据分析

利用上述的似然同步法来比较静息态下从闭眼到睁眼这个过程中脑电信号alpha波的同步性变化。分析结果发现alpha波在几乎所有电极上的同步性都显著地减小。这个发现也与之前EEG文献报道的从闭眼到睁眼状态alpha波去同步化的结果保持一致。

通过仿真实验、真实的脑电信号的分析，将似然同步的非线性分析方法很好地应用到脑电信号的分析，也说明了该方法的可用性及正确性，同时在静息态下从闭眼到睁眼过程中alpha波的同步性显著减小也进一步支持了先前EEG文献报道的alpha波去同步化的结果，为各种认知任务下大脑动力学特性的研究提供一定的帮助^[9-10]。