CN109222965A - 一种同步EEG-fMRI的在线去伪迹方法 - Google Patents

Abstract

本发明公开了一种同步EEG‑fMRI的在线去伪迹方法，包括步骤：1)对脑电信号高通滤波；2)对滤波后的信号上采样，并对脑电信号与MRI设备发出的标记进行同步；3)得到带有残留伪迹的信号S^r；4)信号S^r中有残留的梯度伪迹，构造梯度伪迹模板去噪及PCA去噪，得梯度伪迹消除的信号S_c；5)信号S_c中含有ECG，对ECG进行QRS峰检测；6)将ECG按周期的不同进行K均值聚类；7)对每种EEG信号按时间顺序进行重新整合成m段EEG信号，对这m段EEG信号进行OBS处理；8)按原信号的时间顺序进行重构脑电信号，得除去心电伪迹后的EEG信号。本发明可以在扫描采集fMRI图像的同时采集脑电信号，通过去噪算法对伪迹进行实时的去除，得到信噪比较高的脑电信号。

Description

一种同步EEG-fMRI的在线去伪迹方法

技术领域

本发明涉及脑信号处理的技术领域，尤其是指一种同步EEG-fMRI的在线去伪迹方法。

背景技术

随着脑科学研究的发展，脑电信号(EEG)与功能磁共振成像(fMRI)的结合引起了广泛的关注和研究。脑电信号具有空间分辨率低而时间分辨率高的特点，功能磁共振成像则具有时间分辨率低而空间分辨率高的特点，因此脑电信号与功能磁共振成像的结合，对研究脑的功能和精神类疾病的发病机理以及认知的研究非常重要。

在同步EEG-fMRI的采集脑电信号的实验中，由于EEG信号(脑电信号)是在核磁扫描环境下采集的，在这样的实验环境下主要会引入三种伪迹：

1、脉冲伪迹：它是射频脉冲信号导致的，这类伪迹的频率位于兆赫范围，与EEG有效信号的频率差别很大，故可以设计低通滤波器去除；

2、梯度伪迹：它是扫描时的梯度磁场的变化引起的，这类伪迹是在线去噪过程中要重点考虑的伪迹。由于梯度伪迹的频率范围与EEG信号频率范围相互重叠，所以不能用传统的频率滤波方法来完全去除梯度伪迹；

3、心电伪迹：在高强度的磁场条件下心跳和动脉血流运动带来的伪迹。

EEG信号将淹没在了较大的噪声中，信噪比较低，因此，需要设计一套有效的算法来提高信号的信噪比，消除伪迹对EEG信号的影响。同时，算法的计算量也不能过大，否则由于较高的采样频率(5kHz)，将对在线性能有较大的影响。

目前，同步EEG-fMRI实验过程中很难通过信号去噪算法进行实时去噪，并取得满意的效果。例如中国专利(一种核磁下脑电噪声处理方法，公开号：CN 201310117717)，该发明通过采用差分放大器同时采集脑电信号以及伪迹信号，从硬件层面上进行去噪。中国专利(一种fMRI环境下的脑电装置，公开号：CN 201520599490)，该发明采用屏蔽罩的方式来屏蔽核磁环境下的噪声信号。以上所述发明，都未能通过信号去噪算法进行同步实时去噪并取得满意的效果。通过同步实时去伪迹算法，除了设备上的简便外，伪迹的去除也能达到更好的效果。

发明内容

本发明的目的在于克服现有技术的不足，提出了一种同步EEG-fMRI的在线去伪迹方法，利用此方法可以在扫描采集fMRI图像的同时采集脑电信号，通过去噪算法对伪迹进行实时的去除，得到信噪比较高的脑电信号。

为实现上述目的，本发明所提供的技术方案为：一种同步EEG-fMRI的在线去伪迹方法，主要包括去除梯度伪迹和去除心电伪迹两部分，包括以下步骤：

1)对采集到的脑电信号进行高通滤波，滤除温漂及脉冲伪迹；

2)对步骤1)中滤波后的信号进行上采样，并对脑电信号与MRI设备发出的标记进行同步，即将每个fMRI扫描切片过程中采集的脑电信号与伪迹模板进行滑动互相关运算，运算结果最大的位移即为同步后调整的位置，得到信号S_h；由于脑电信号是随机信号，梯度伪迹则以一个切片扫描时间T为周期，每个周期内该噪声的形状相似，因此能够构造一个滑动窗口来构造噪声模板A_r，其中τ表示模板的时间序号，基于此模板进行初步去噪；

3)得到新的切片梯度伪迹模板后，将混合着伪迹的信号S_h划出一个长度为T*w个点的切片，其中w为采样频率，将噪声模板以最小二乘拟合切片数据，得到拟合参数y_τ，从S_h减去y_τ*A_τ，得到带有残留伪迹的信号S^r；

4)步骤3)得到的信号S^r中依然有残留的梯度伪迹，对其进行主成分分析PCA，再将PCA得到各成分与梯度伪迹模板进行相关操作，将各成分按照相关程度大小排序，取最大的f个成分作为梯度伪迹的最优基β_j,j＝1,2,...,f，将其以最小二乘拟合S^r，得到拟合参数a_j，j＝1,2,...,f，从S^r减去a_j*β_j，经过构造梯度伪迹模板去噪及PCA去噪，最终可得梯度伪迹消除的信号S_c；

5)步骤4)得到的信号S_c中还含有心电伪迹ECG，对ECG进行QRS峰检测，记录下各个相邻R峰之间的间距，即是心电伪迹的周期T_i，i＝1,2,...,n；

6)将ECG按周期T_i，i＝1,2,...,n的不同进行K均值聚类，初始选定m个质心，分为m种ECG，由于ECG在脑电信号EEG上时延T_s，即能够分为m种EEG信号，并记录下ECG中的R峰反映在EEG中的位置R_i，i＝1,2,...,n；

7)对每种EEG信号按时间顺序进行重新整合成m段EEG信号，对这m段EEG信号进行OBS处理；

8)按原信号的时间顺序进行重构脑电信号，最终得除去心电伪迹后的EEG信号。

在步骤2)中，具体处理过程如下：

2.1)以N个切片扫描时间长度N*T为滑动窗口的宽度，在脑电信号上边进行移动，将窗口中的信号以T为间隔进行划分，设采样频率为w，即q＝T*w个采样点为一个切片间隔的数据长度；

2.2)窗口中的一个切片数据段为将N个进行叠加，再除以N取平均，得到梯度伪迹的模板

在步骤4)中，具体处理过程如下：

4.1)在信号S^r中构造一个PCA滑动窗口，将其以一个切片的时间间隔T划分，即T*w为一个数据段，构造一个PCA样本矩阵其中p为PCA样本的个数，q则为一个样本的维数，即q＝T*w，将每个新采集的数据段去平均后存进并删除最早的一个数据段，保持p个样本数；

4.2)对转置得到矩阵之后计算矩阵的特征向量v_i及特征值σ_i，则样本的协方差矩阵的特征向量

4.3)将特征向量u_i与中的最新一个切片数据段作相关运算，将u_i按照相关运算的结果从大到小排序，取最大的f个u_i作为梯度伪迹的最优基β_j，j＝1,2,...,f，将其以最小二乘拟合得到拟合参数a_j，j＝1,2,...,f,即从减去a_j*β_j，最终经过构造梯度伪迹模板去噪及PCA去噪，能够将梯度伪迹消除。

在步骤6)中，具体处理过程如下：

6.1)将ECG按周期T_i，i＝1,2,...,n的不同进行K均值聚类，第一步：初始选定m个质心；第二步：将每个周期的ECG按周期的不同归纳到某个类别中，并更新质心；第三步：循环第二步，直到质心不发生变化或者所有的ECG都已经处理完为止；经过K均值聚类之后，原ECG能够分为m种不同的ECG；

6.2)由于ECG在脑电信号EEG上时延T_s，对步骤6.1)中分类好的ECG按时间序列上时延T_s，即能够分为m种EEG信号，并记录下每段EEG信号中的时间顺序。

在步骤7)中，具体处理过程如下：

7.1)对步骤6)中的m种EEG信号按时间顺序进行重新整合成m段EEG信号，对这m段EEG信号进行步骤7.2)中的OBS处理；

7.2)对于每段EEG信号的R_i，i＝1,2,...,n，取每个R_i，i＝1,2,...,n的中位数R_m，对于每个周期EEG信号的R峰反映在EEG中的位置R_i，i＝1,2,...,n，取的信号段S^l，l＝1,2,...,n，构造一个PCA样本矩，对其进行主成分分析PCA，将各成分按照PCA后的特征值大小排序，取最大的4个成分作为心电伪迹的最优基γ_j，j＝1,2,3,4，将其以最小二乘拟合S^l，l＝1,2,...,n，得到拟合参数b_j，j＝1,2,3,4，从S^l，l＝1,2,...,n减去γ_j*b_j，最终可得去除心电伪迹的信号S_e。

本发明与现有技术相比，具有如下优点与有益效果：

1、可以对EEG-fMRI同步采集信号进行实时去噪，满足了在线实验的数据处理速度要求。

2、无需通过特殊硬件结构来实现去噪功能，只需用普通的核磁共振下脑电放大器即可完成脑电信号的采集，对设备的要求大大降低。

3、实时去噪后的EEG信号基本将梯度伪迹和心电伪迹滤除，相比基于硬件去噪与其他的去噪算法，本处理算法信噪比更高，延迟小，能更好地满足在线实验的信号信噪比需求。

附图说明

图1为心电伪迹的去噪流程图。

具体实施方式

下面结合具体实施例对本发明作进一步说明。

本实施例所提供的同步EEG-fMRI的在线去伪迹方法，可以应用在fMRI及EEG混合脑机接口平台中，包括以下步骤：

1)设计程序框架，程序包含记录放大器采集的信号和信号的去噪处理两部分。

2)信号采集线程首先对需要进行去噪的每个放大器通道都单独创建一个线程，负责每个通道信号的并行处理。

3)将采集到的脑电信号进行高通滤波，滤除温漂及脉冲伪迹，滤波器设计得到参数，之后将滤波器参数导出到程序框架的文件中使用。

4)对步骤3)中滤波后的信号进行上采样并对脑电信号与MRI设备发出的标记进行同步，即将每个fMRI扫描切片过程中采集的脑电信号与伪迹模板进行滑动互相关运算，运算结果最大的位移即为同步后调整的位置，得到信号S_h；由于脑电信号是随机信号，梯度伪迹则以一个切片扫描时间T为周期，每个周期内该伪迹的形状相似，因此可以构造一个滑动窗口来构造伪迹模板A_r，其中τ表示模板的时间序号，基于此模板进行初步去噪，包括以下步骤：

4.1)以N个切片扫描时间长度N*T为滑动窗口的宽度，在脑电信号上边进行移动，将窗口中的信号以T为间隔进行划分，设采样频率为w，即q＝T*w个采样点为一个切片间隔的数据长度；

4.2)窗口中的一个切片数据段为将N个进行叠加，再除以N取平均，得到梯度伪迹的模板

5)得到新的切片梯度伪迹模板后，将混合着伪迹的信号S_h划出一个长度为T*w个点的切片，其中w为采样频率，将伪迹模板以最小二乘拟合切片数据，得到拟合参数y_τ，从S_h减去y_τ*A_τ，得到带有残留伪迹的信号S^r，因为脑电信号是随机信号，可以认为梯度伪迹模板中由于多次叠加取平均脑电信号成分趋近于0，使得减去模板后的信号中脑电成分得到完整保留。

6)对步骤5)得到的信号S^r中依然有残留伪迹，对其进行主成分分析(PCA)，再将PCA得到各成分与梯度伪迹模板进行相关操作，将各成分按照相关程度大小排序，取最大的f个成分作为梯度伪迹的最优基β_j,j＝1,2,...,f，将其以最小二乘拟合S^r,得到拟合参数a_j，j＝1,2,...,f，从减去a_j*β_j，经过构造梯度伪迹模板去噪及PCA去噪，最终可将梯度伪迹消除，具体步骤如下：

6.1)在信号S^r中构造一个PCA滑动窗口，将其以一个切片的时间间隔T划分，即T*w为一个数据段，构造一个PCA样本矩阵其中p为PCA样本的个数，q则为一个样本的维数，即q＝T*w，将每个新采集的数据段去平均后存进并删除最早的一个数据段，保持p个样本数；

6.2)对转置得到矩阵之后计算矩阵的特征向量v_i及特征值σ_i，则样本的协方差矩阵的特征向量

6.3)将特征向量u_i与中的最新一个切片数据段作相关运算，将u_i按照相关运算的结果从大到小排序，取最大的f个u_i作为梯度伪迹的最优基β_j,j＝1,2,...,f，将其以最小二乘拟合得到拟合参数a_j，j＝1,2,...,f,即从减去a_j*β_j，最终经过构造梯度伪迹模板去噪及PCA去噪，可得到将梯度伪迹消除的信号S_c。

7)步骤6)得到的信号S_c中还含有心电伪迹(ECG)，对ECG进行QRS峰检测，记录下各个相邻R峰之间的间距，即是心电伪迹的周期T_i(i＝1,2,...,n)。

8)将ECG按周期T_i(i＝1,2,...,n)的不同进行K均值聚类，初始选定m个质心，分为m种ECG，由于ECG在脑电信号(EEG)上时延T_s(T_s通常为210ms)，即可以分为m种EEG信号，并记录下ECG中的R峰反映在EEG中的位置R_i(i＝1,2,...,n)，具体步骤如下：

8.1)将ECG按周期T_i(i＝1,2,...,n)的不同进行K均值聚类,第一步：初始选定m个质心；第二步：将每个周期的ECG按周期的不同归纳到某个类别中，并更新质心；第三步：循环第二步，直到质心不发生变化或者所有的ECG都已经处理完为止。经过K均值聚类之后，原ECG可以分为m种不同的ECG；

8.2)由于ECG在脑电信号(EEG)上时延T_s(T_s通常为210ms)，对步骤9.1)中分类好的ECG按时间序列上时延T_s，即可以分为m种EEG信号，并记录下每段EEG信号中的时间顺序。

9)对于每种EEG信号按时间顺序进行重新整合成m段EEG信号，对这m段EEG信号进行OBS处理，具体处理过程如下：

9.1)对于步骤8)中的m种EEG信号按时间顺序进行重新整合成m段EEG信号，对这m段EEG信号进行步骤10.2)中的OBS处理；

9.2)对于每段EEG信号的R_i(i＝1,2,...,n)，取每个R_i(i＝1,2,...,n)的中位数R_m，对于每个周期EEG信号的R峰反映在EEG中的位置R_i(i＝1,2,...,n)，取的信号段S^l(l＝1,2,...,n)，构造一个PCA样本矩，对其进行主成分分析(PCA)，将各成分按照PCA后的特征值大小排序，取最大的4个成分作为心电伪迹的最优基γ_j(j＝1,2,3,4)，将其以最小二乘拟合S^l(l＝1,2,...,n)，得到拟合参数b_j(j＝1,2,3,4)，从S^l(l＝1,2,...,n)减去γ_j*b_j，最终可得去除心电伪迹的信号S_e。

10)按原信号的时间顺序进行重构脑电信号，最终得除去心电伪迹后的EEG信号。

心电伪迹的去噪流程具体请参见图1所示。

以上所述实施例只为本发明之较佳实施例，并非以此限制本发明的实施范围，故凡依本发明之形状、原理所作的变化，均应涵盖在本发明的保护范围内。

Claims (5)

1.一种同步EEG-fMRI的在线去伪迹方法，其特征在于：主要包括去除梯度伪迹和去除心电伪迹两部分，包括以下步骤：

1)对采集到的脑电信号进行高通滤波，滤除温漂及脉冲伪迹；

2)对步骤1)中滤波后的信号进行上采样，并对脑电信号与MRI设备发出的标记进行同步，即将每个fMRI扫描切片过程中采集的脑电信号与伪迹模板进行滑动互相关运算，运算结果最大的位移即为同步后调整的位置，得到信号S_h；由于脑电信号是随机信号，梯度伪迹则以一个切片扫描时间T为周期，每个周期内该噪声的形状相似，因此能够构造一个滑动窗口来构造噪声模板A_r，其中τ表示模板的时间序号，基于此模板进行初步去噪；

3)得到新的切片梯度伪迹模板后，将混合着伪迹的信号S_h划出一个长度为T*w个点的切片，其中w为采样频率，将噪声模板以最小二乘拟合切片数据，得到拟合参数y_τ，从S_h减去y_τ*A_τ，得到带有残留伪迹的信号S^r；

4)步骤3)得到的信号S^r中依然有残留的梯度伪迹，对其进行主成分分析PCA，再将PCA得到各成分与梯度伪迹模板进行相关操作，将各成分按照相关程度大小排序，取最大的f个成分作为梯度伪迹的最优基β_j，j＝1,2,...,f，将其以最小二乘拟合S^r，得到拟合参数a_j，j＝1,2,...,f，从S^r减去a_j*β_j，经过构造梯度伪迹模板去噪及PCA去噪，最终可得梯度伪迹消除的信号S_c；

5)步骤4)得到的信号S_c中还含有心电伪迹ECG，对ECG进行QRS峰检测，记录下各个相邻R峰之间的间距，即是心电伪迹的周期T_i，i＝1,2,...,n；

6)将ECG按周期T_i，i＝1,2,...,n的不同进行K均值聚类，初始选定m个质心，分为m种ECG，由于ECG在脑电信号EEG上时延T_s，即能够分为m种EEG信号，并记录下ECG中的R峰反映在EEG中的位置R_i，i＝1,2,...,n；

7)对每种EEG信号按时间顺序进行重新整合成m段EEG信号，对这m段EEG信号进行OBS处理；

8)按原信号的时间顺序进行重构脑电信号，最终得除去心电伪迹后的EEG信号。

2.根据权利要求1所述的一种同步EEG-fMRI的在线去伪迹方法，其特征在于，在步骤2)中，具体处理过程如下：

2.1)以N个切片扫描时间长度N*T为滑动窗口的宽度，在脑电信号上边进行移动，将窗口中的信号以T为间隔进行划分，设采样频率为w，即q＝T*w个采样点为一个切片间隔的数据长度；

2.2)窗口中的一个切片数据段为将N个进行叠加，再除以N取平均，得到梯度伪迹的模板

3.根据权利要求1所述的一种同步EEG-fMRI的在线去伪迹方法，其特征在于，在步骤4)中，具体处理过程如下：

4.1)在信号S^r中构造一个PCA滑动窗口，将其以一个切片的时间间隔T划分，即T*w为一个数据段，构造一个PCA样本矩阵S_p*q，其中p为PCA样本的个数，q则为一个样本的维数，即q＝T*w，将每个新采集的数据段去平均后存进S_p*q，并删除最早的一个数据段，保持p个样本数；

4.2)对S_p*q转置得到矩阵之后计算矩阵的特征向量v_i及特征值σ_i，则样本的协方差矩阵的特征向量

4.3)将特征向量u_i与S_p*q中的最新一个切片数据段作相关运算，将u_i按照相关运算的结果从大到小排序，取最大的f个u_i作为梯度伪迹的最优基β_j，j＝1,2,...,f，将其以最小二乘拟合得到拟合参数a_j，j＝1,2,...,f,即从减去a_j*β_j，最终经过构造梯度伪迹模板去噪及PCA去噪，能够将梯度伪迹消除。

4.根据权利要求1所述的一种同步EEG-fMRI的在线去伪迹方法，其特征在于，在步骤6)中，具体处理过程如下：

6.1)将ECG按周期T_i，i＝1,2,...,n的不同进行K均值聚类，第一步：初始选定m个质心；第二步：将每个周期的ECG按周期的不同归纳到某个类别中，并更新质心；第三步：循环第二步，直到质心不发生变化或者所有的ECG都已经处理完为止；经过K均值聚类之后，原ECG能够分为m种不同的ECG；

6.2)由于ECG在脑电信号EEG上时延T_s，对步骤6.1)中分类好的ECG按时间序列上时延T_s，即能够分为m种EEG信号，并记录下每段EEG信号中的时间顺序。

5.根据权利要求1所述的一种同步EEG-fMRI的在线去伪迹方法，其特征在于，在步骤7)中，具体处理过程如下：

7.1)对步骤6)中的m种EEG信号按时间顺序进行重新整合成m段EEG信号，对这m段EEG信号进行步骤7.2)中的OBS处理；

7.2)对于每段EEG信号的R_i，i＝1,2,...,n，取每个R_i，i＝1,2,...,n的中位数R_m，对于每个周期EEG信号的R峰反映在EEG中的位置R_i，i＝1,2,...,n，取的信号段S^l，l＝1,2,...,n，构造一个PCA样本矩，对其进行主成分分析PCA，将各成分按照PCA后的特征值大小排序，取最大的4个成分作为心电伪迹的最优基γ_j，j＝1,2,3,4，将其以最小二乘拟合S^l，l＝1,2,...,n，得到拟合参数b_j，j＝1,2,3,4，从S^l，l＝1,2,...,n减去γ_j*b_j，最终可得去除心电伪迹的信号S_e。