CN113988143A

CN113988143A - 基于光泵磁强计测量脑磁的坏道自动识别方法

Info

Publication number: CN113988143A
Application number: CN202111608900.0A
Authority: CN
Inventors: 宁晓琳; 王若男; 贾智慧; 高阳; 房建成
Original assignee: Hangzhou Innovation Research Institute of Beihang University
Current assignee: Hangzhou Innovation Research Institute of Beihang University
Priority date: 2021-12-27
Filing date: 2021-12-27
Publication date: 2022-01-28
Anticipated expiration: 2041-12-27
Also published as: CN113988143B

Abstract

本发明涉及一种基于光泵磁强计测量脑磁的坏道自动识别方法，方法包括：基于光泵磁强计，获取待分析对象的原始多通道脑磁信号；将原始多通道脑磁信号变换到频域，并计算变换到频域后的多通道脑磁信号的功率谱密度；采用遍历方式，对每个通道脑磁数据的功率谱密度在指定频段内求平均值，并对平均值取对数，获得所有通道的平均功率密度序列；采用箱线图方法对所述平均功率密度序列进行坏道识别，确定坏道占整体脑磁测量通道的比例；并采用孤立森林方法对原始多通道脑磁信号进行脑磁坏道的自动识别并处理，得到去除坏道的多通道脑磁信号。本发明的方法提高了识别的效率，实现坏道识别的自动化处理。

Description

基于光泵磁强计测量脑磁的坏道自动识别方法

技术领域

本发明涉及生物信号处理技术，尤其涉及一种基于光泵磁强计测量脑磁的坏道自动识别方法。

背景技术

脑磁图是一种能够实时、完全无创地测量大脑神经活动的磁源成像技术，通过磁场反应大脑的生理活动和功能连接，相比于脑电和一些其他的测量手段具有更高的时间和空间分辨率。

用于测量脑磁的设备主要有超导量子干涉仪(Superconducting QuantumInterference Device ,SQUID)和小型化的光泵磁强计（Optical Pumping Magnetometer，OPM），由于光泵磁强计便携以及不需要低温冷却的优点，这种传感器在近几年发展快速。

但在利用光泵磁强计在测量多通道脑磁图的过程中，容易受到坏道的影响。脑磁测量数据中出现坏道主要是由于测量过程中光泵磁强计的不稳定以及实验操作不当导致，例如传感器的过热、超过预定的校准因子等。因此对脑磁数据进行降噪之前，必要的一个过程是对多通道脑磁数据的坏道识别剔除；如果在最初不对坏道进行剔除，会导致坏道的数据对后续的预处理、正逆问题求解过程产生影响。

目前对于多通道脑磁图中坏道的在实际的应用中主要存在三种方法：（1）信号空间分离算法SSS以及拓展的rSSS方法进行坏道识别；（2）参考脑电等其他电生理信号，利用数据时域的特征设计分类器进行坏道自动识别，包括通道数据的相关系数、偏差、振幅等；（3）实验人员进行人工识别剔除。信号空间分离算法SSS及其拓展rSSS在预处理数据的过程中主要可以用来去除噪声、识别坏道、校正头部运动，但整体算法较为复杂而且依赖与硬件传感器校准；利用数据时域的特征是通过计算每个通道的时域特征，对比之后去除超过阈值范围的通道便认为是坏道，但脑磁测量过程中实验范式较长，这种方法和一些基于数学统计量的分类器往往起不到很好的效果；而实验人员肉眼识别在脑磁通道较多的情况下则会耗费较多时间，效率也并不高。

发明内容

（一）要解决的技术问题

针对现有技术无法自动化识别坏道的不足以及在实际脑磁测量实验中对坏道识别的需求，本发明提供一种基于光泵磁强计测量脑磁的坏道自动识别方法，能够实现脑磁坏道识别的自动化流程实现。

（二）技术方案

为了达到上述目的，本发明采用的主要技术方案包括：

第一方面，本发明实施例提供一种基于光泵磁强计测量脑磁的坏道自动识别方法，包括：

步骤A：基于光泵磁强计，获取待分析对象的原始多通道脑磁信号；所述光泵磁强计包括：多个用于接触分析对象指定区域的具有自旋交换无弛豫效应的极弱磁测量传感器和连接每一极弱磁测量传感器的处理装置；

步骤B：将原始多通道脑磁信号变换到频域，并计算变换到频域后的多通道脑磁信号的功率谱密度；

步骤C：采用遍历方式，对每个通道脑磁数据的功率谱密度在指定频段内求平均值，并对平均值取对数，获得所有通道的平均功率密度序列

；

步骤D：采用箱线图方法对所述平均功率密度序列

进行坏道识别，确定坏道占整体脑磁测量通道的比例

；

步骤E：根据坏道占整体脑磁测量通道的比例

，采用孤立森林方法对原始多通道脑磁信号进行脑磁坏道的自动识别并处理，得到去除坏道的多通道脑磁信号。

可选地，所述步骤B包括：

所述原始多通道脑磁信号为二维的数据点矩阵，对每一个通道的脑磁数据利用Welch算法变换到频域，并计算频域下每一通道的脑磁数据的功率谱密度，得到多通道脑磁信号的功率谱密度。

可选地，所述步骤C包括：

将每一个通道的功率谱密度在指定频段内取平均值，获得用于表征多有通道脑磁信号平均功率密度的序列

；

对所述序列

取对数，得到所有通道的平均功率密度序列

。

可选地，所述指定频段为对应脑磁测量实验期望信号的频段；

或者，所述指定频段为0-60Hz、0-20Hz、20-40Hz、20-60Hz、40-60Hz或者0-40Hz。

可选地，所述步骤D包括：

子步骤D1：对平均功率密度序列

建立一个箱线图的模型；

子步骤D2：基于所述模型，获取平均功率密度序列

中P1分位数对应的数值Q1，P2分位数对应的数值Q2，P3分位数对应的数值Q3；

子步骤D3：计算箱线图的上须和下须；

子步骤D4：将功率密度序列中大于上须或小于下须的数值定义为是异常的脑磁测量通道即识别为坏道，获取坏道占所有通道的比例

。

可选地，所述P1为25%、P2为50%、P3为75%；

上须=Q3+1.5*（Q3- Q1），下须=Q1-1.5*（Q3- Q1）；

比例

=m/n，其中，m表示坏道数量，n表示具有自旋交换无弛豫效应的极弱磁测量传感器的数量。

可选地，所述步骤E包括：

子步骤E1：从平均功率密度序列

中随机选择N个功率值点作为样本子集，放入树的根节点，构建一个二叉搜索树；N≤n；

子步骤E2：指定功率的大小作为一个特征，在二叉搜索树的当前根节点中随机产生一个分离点p，以此分离点为界，分为小于p的左子树和大于p的右子树；

子步骤E3：在子节点中递归子步骤E2和E3，不断构造新的子节点，直到子节点中只有一个数据或子节点已到达限定高度；

子步骤E4：循环子步骤E1至子步骤E3，直至生成T个孤立树；之后对于每一个通道的功率数据点在T个孤立树中求解平均高度；并对高度做归一化处理，计算每个功率点的异常值分数；

子步骤E5：根据所述比例

指导异常值得分自动选取阈值，获得孤立森林脑磁坏道自动识别的最终结果；

子步骤E6：将识别到的脑磁坏道从原始多通道脑磁信号中删除。

第二方面，本发明实施例还提供一种电子设备，其包括存储器和处理器，所述存储器中存储计算机指令，所述处理器执行所述存储器中存储的计算机指令，具体执行上述第一方面任一所述的基于光泵磁强计测量脑磁的坏道自动识别方法的步骤。

第三方面，本发明实施例还提供一种计算机存储介质，其存储有计算机程序，所述计算机程序被处理器执行时实现如第一方面任一所述的基于光泵磁强计测量脑磁的坏道自动识别方法的步骤。

（三）有益效果

本发明的方法中将测量得到的多通道脑磁数据变换到频域计算功率谱密度，利用正常工作光泵磁强计与坏道在采集数据功率谱上的明显差别进行坏道的识别，提高了识别的效率，且实现了坏道识别的自动化处理，可有效推广使用。

另外，本发明的方法解决了由于使用范式采集时间较长，时域上的幅值、偏差等特征区别不明显的问题，相比于时域的识别方法提高了识别的效率；同时该自动识别方法是脑磁信号预处理的一部分，在导入数据之后最先进行，避免了在坏道之前还需要滤波等其他的一些步骤，处理效率更高。

此外，在本发明实施例中，首先使用箱线图的方法初步识别多通道脑磁数据中的坏道，得到一个预估的坏道占整体测量通道的比例，之后再通过孤立森林的方法进行进一步的坏道识别。通过这样的方法使整个坏道自动识别的流程实现自动化，不需要单独在孤立森林算法中输入坏道的比例用于确定异常值分数的阈值，使整个坏道识别流程实现自动化；代替医疗工作者或者实验人员耗时耗力的人工识别过程。

本发明中应用孤立森林算法实现更为稳定的识别坏道，但是由于孤立森林算法中还需要一个异常通道的比例，用于指导异常值得分的选取阈值。因此考虑先进行箱线图检测预估坏道的异常值比例，之后再进行孤立森林进一步检测坏道，实现无需输入控制参数的坏道自动识别。提高整体坏道自动识别的稳定性和准确性。

附图说明

图1为本发明一实施例提供的基于光泵磁强计测量脑磁的坏道自动识别方法的流程图；

图2为使用箱线图进行初步坏道识别的示意图；

图3为使用孤立森林方法进行坏道识别的实现流程图；

图4为未剔除坏道的平均功率密度序列和剔除坏道的平均功率密度序列的示意图；

图5为本发明一实施例提供的电子设备的示意图。

具体实施方式

为了更好的解释本发明，以便于理解，下面结合附图，通过具体实施方式，对本发明作详细描述。

本发明公开了一种基于光泵磁强计测量脑磁的坏道自动识别方法，涉及脑科学、脑磁图数据采集预处理技术领域。所述方法采集基于光泵磁强计测量的多通道脑磁数据；对采集的脑磁数据进行频域变换得到多通道数据的功率谱密度；在特定频段范围内求每个通道的平均功率之后再取对数，得到一个表征多通道脑磁数据的功率序列，之后进行功率异常值的坏道检测；使用箱线图估计功率序列中坏道的比例，基于坏道的比例使用孤立森林方法进行脑磁坏道的自动识别。本发明提供的脑磁数据坏道识别方法，可以在预处理前进行坏道的识别，使用箱线图估计坏道比例使整个流程实现自动化，提高处理运行效率，代替医疗工作者耗时耗力的人工识别过程。

实施例一

如图1所示，本实施例提供一种基于光泵磁强计测量脑磁的坏道自动识别方法，本实施例的方法的执行主体可为任一电子设备，本实施例的方法可包括下述的步骤：

A. 使用光泵磁强计测量多通道脑磁信号；

举例来说，在进行脑磁测量的过程中使用基于自旋交换无弛豫效应的小型化极弱磁测量传感器，本实施例中进行脑磁听觉诱发响应数据的测量；坏道是指在进行听觉诱发数据测量过程中数据无法正常使用和反映听觉诱发响应特征的脑磁测量通道，造成这种原因主要是由于光泵磁强计的不稳定以及实验操作不当导致，例如传感器的过热、超过预定的校准因子等导致测量。

B. 将测量得到的脑磁数据变换到频域计算多通道脑磁数据的功率谱密度。

本实施例中，测量得到的多通道脑磁数据是一个二维的数据点矩阵，对其中每一个通道利用Welch算法计算功率谱密度，得到每个脑磁数据通道的功率谱密度。

C. 将得到的功率谱平均并取对数。

即，对每个通道的功率谱密度在特定频段内求平均值，之后再对平均值取对数，得到一个多通道的平均功率密度序列

。

图4中示出了未剔除坏道的平均功率密度序列和剔除坏道的平均功率密度序列的示意图，图4中横轴是功率，单位是HZ，一般脑磁有用信号在1-40Hz内，图4中取前60Hz画的图，纵轴就是在某功率下信号能量密度，单位是

，可以看到纵轴单位值较小，取对数便于后续做基本得统计量分析，所以将纵轴变为对数序列，每一条曲线是一个通道的平均功率谱密度，总共是n个通道，N个功率值点是指对纵轴功率对数序列进行采样点划分，右侧表示的是根据纵轴上的偏离群体的异常曲线进行坏道去除。

将步骤B中得到的功率谱在0-60Hz内按通道取平均值，得到一个表征多通道脑磁数据平均功率的序列

，之后再对整个序列取以10 为底的对数得到

。

在本实施例中听觉诱发响应的期望信号频段在60Hz以内，因此特定频段平均选择0-60Hz；但并不局限于0-60Hz，也可以是0-20Hz，0-40Hz等频段，具体根据不同脑磁测量实验期望信号频段设置。

D. 对上述的序列使用箱线图方法进行初步的坏道识别，预估坏道占整体脑磁测量通道的比例

。

如图2所示，对步骤C中得到的序列

建立一个箱线图的模型；计算25%分位数对应的数值Q1、50%分位数对应的数值Q2、75%分位数对应的数值Q3；以及箱线图的上须和下须：上须=Q3+1.5*（Q3-Q1），下须=Q1-1.5*（Q3-Q1）。将在箱线图中小于下须大于上须的功率点视为异常点，初步识别为脑磁测量中的坏道。

本实施例中假设步骤识别到m个坏道，总共用于测量的光泵磁强计有n个；则预估的坏道比例

=m/n。

E. 根据得到的脑磁坏道比例，使用孤立森林方法进行脑磁坏道的自动识别过程。

如图3所示，对上述的功率序列构建二叉树模型，使用孤立森林算法自动识别脑磁坏道。

在一种具体的实现过程中，上述步骤E可包括下述的图中未示出的子步骤：

子步骤E1，从脑磁平均功率序列中随机选择N（2<N≤n）个功率值点作为样本子集，放入树的根节点，构建一个二叉树搜索树；

子步骤E2，指定功率的大小作为一个特征，在当前根节点中随机产生一个分离点p，以此分离点为界，分为小于p的左子树和大于p的右子树；

子步骤E3，在子节点中递归步骤E2和E3，不断构造新的子节点，直到子节点中只有一个数据（无法再继续切割）或子节点已到达限定高度；

子步骤E4，循环E1至E3，直至生成T（如100个）个孤立树。

之后对每一个通道的功率数据点在100个孤立树中求解平均高度；并对高度做归一化处理，计算每个功率点的异常值分数

；h(x)为样本点x在每棵树中的高度，c(N)为给定样本数为N时路径的平均值，用于高度的归一化处理；在本实施例中c(N)=2H （N-1）-2（N-1）/N，

H（N-1）表示调和数，

，

表示欧拉常数。

也就是说，对于每一个通道的功率数据点x在上述的T个孤立树中求解平均高度h (x)；并对高度做归一化处理，之后根据h(x)得到异常值分数：

H（n）是调和数，可以通过欧拉常数来估算；最终异常值S（x，N）得分越高，是坏道的可能性越大，如图3所示处理过程的示意图。

子步骤E5，根据步骤D中得到的脑磁坏道的比例

指导异常值得分自动选取阈值，得到孤立森林脑磁坏道自动识别的最终结果；之后将识别到的脑磁通道从原始数据中删除。

本实施例中的孤立森林可理解为一种无监督学习算法，属于组合决策树家族。在孤立森林（iForest）中，异常点被定义为容易被孤立的离群点，可以将其理解为分布稀疏且离密度高的群体较远的点。在脑磁测量过程中，坏道则对应于由传感器等问题导致测量数据功率偏移高密度的正常通道功率密度区的脑磁通道，因此选用孤立森林的方法进行坏道的自动识别。

上述方法将测量得到的多通道脑磁数据变换到频域计算功率谱密度，利用正常工作光泵磁强计与坏道在采集数据功率谱上的明显差别进行坏道的识别，提高了识别的效率，且实现了坏道识别的自动化处理，可有效推广使用。

在本发明实施例中，首先使用箱线图的方法初步识别多通道脑磁数据中的坏道，得到一个预估的坏道占整体测量通道的比例，之后再通过孤立森林的方法进行进一步的坏道识别。通过这样的方法使整个坏道自动识别的流程实现自动化，不需要单独在孤立森林算法中输入坏道的比例用于确定异常值分数的阈值，使整个坏道识别流程实现自动化；代替医疗工作者或者实验人员耗时耗力的人工识别过程。

本实施例中，将测量得到的数据变换到频域求得特定频段内的平均功率，之后再取对数得到处理后的一个脑磁功率序列。根据坏道和正常状态下传感器工作功率上的区别，利用箱线图和孤立森林方法进行坏道的自动识别。通过实测脑磁听觉诱发以及空桶数据验证表明，该方法可以有效的识别脑磁测量过程中因传感器等硬件问题导致的坏道，提高了剔除坏道数据的效率，方便后续脑磁的预处理、源定位等过程。

实施例二

本实施例还提供一种电子设备，包括：存储器和处理器；所述处理器用于执行所述存储器中存储的计算机程序，以实现执行上述实施例一和实施例二任意所述的坏道自动识别方法的步骤。

具体地，如图5所示，本实施例的电子设备可包括：至少一个处理器51、至少一个存储器52、至少一个网络接口54和/或其他的用户接口53。电子设备中的各个组件通过总线系统55耦合在一起。可理解，总线系统55用于实现这些组件之间的连接通信。总线系统55除包括数据总线之外，还包括电源总线、控制总线和状态信号总线。但是为了清楚说明起见，在图5中将各种总线都标为总线系统55。

本实施例的电子设备可以执行图1至图3任意所示的方法，其中，用户接口53可以包括显示器、键盘或者点击设备(例如，鼠标或者触感板等)。可以理解，本实施例中的存储器52可以是易失性存储器或非易失性存储器，或可包括易失性和非易失性存储器两者。本文描述的存储器52包括任意其它适合类型的存储器。

在一些实施方式中，存储器52存储了如下的元素，可执行单元或者数据结构，或者他们的子集，或者他们的扩展集：操作系统521和应用程序522。其中，操作系统521，包含各种系统程序，例如框架层、核心库层、驱动层等，用于实现各种基础业务以及处理基于硬件的任务。应用程序522，包含各种应用程序，例如媒体播放器(MediaPlayer)、浏览器(Browser)等，用于实现各种应用业务。实现本发明实施例方法的程序可以包含在应用程序522中。

在本发明实施例中，处理器51通过调用存储器52存储的程序或指令，具体的，可以是应用程序522中存储的程序或指令，处理器51用于执行第一方面所提供的方法步骤。

上述本发明实施例揭示的方法可以应用于处理器51中，或者由处理器51实现。处理器51可能是一种集成电路芯片，具有信号的处理能力。在实现过程中，上述方法的各步骤可以通过处理器51中的硬件的集成逻辑电路或者软件形式的指令完成。上述的处理器51可以是通用处理器、数字信号处理器(DigitalSignalProcessor，DSP)、专用集成电路(ApplicationSpecificIntegratedCircuit，ASIC)、现成可编程门阵列(FieldProgrammableGateArray， FPGA)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件。可以实现或者执行本发明实施例中的公开的各方法、步骤及逻辑框图。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等。结合本发明实施例所公开的方法的步骤可以直接体现为硬件译码处理器执行完成，或者用译码处理器中的硬件及软件单元组合执行完成。软件单元可以位于随机存储器，闪存、只读存储器，可编程只读存储器或者电可擦写可编程存储器、寄存器等本领域成熟的存储介质中。该存储介质位于存储器52，处理器51读取存储器52中的信息，结合其硬件完成上述方法的步骤。

另一方面，本发明实施例还提供一种计算机可读存储介质，其用于存储有计算机程序，所述计算机程序被处理器执行时实现上述任意实施例的坏道自动识别方法的步骤。

应当注意的是，在权利要求中，不应将位于括号之间的任何附图标记理解成对权利要求的限制。词语“包含”不排除存在未列在权利要求中的部件或步骤。位于部件之前的词语“一”或“一个”不排除存在多个这样的部件。本发明可以借助于包括有若干不同部件的硬件以及借助于适当编程的计算机来实现。在列举了若干装置的权利要求中，这些装置中的若干个可以是通过同一个硬件来具体体现。词语第一、第二、第三等的使用，仅是为了表述方便，而不表示任何顺序。可将这些词语理解为部件名称的一部分。

此外，需要说明的是，在本说明书的描述中，术语“一个实施例”、“一些实施例”、“实施例”、“示例”、“具体示例”或“一些示例”等的描述，是指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任一个或多个实施例或示例中以合适的方式结合。此外，在不相互矛盾的情况下，本领域的技术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结合和组合。

尽管已描述了本发明的优选实施例，但本领域的技术人员在得知了基本创造性概念后，则可对这些实施例作出另外的变更和修改。所以，权利要求应该解释为包括优选实施例以及落入本发明范围的所有变更和修改。

显然，本领域的技术人员可以对本发明进行各种修改和变型而不脱离本发明的精神和范围。这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也应该包含这些修改和变型在内。