动态脑磁源定位方法
【技术领域】
[0001] 本发明设及生物信息技术领域，特别是一种动态脑磁源定位方法。
【背景技术】
[0002] 脑功能正常运转不是由单个神经元或单一脑区独立完成，而是依靠神经元交互作 用形成的神经功能连接网络。神经成像技术广泛用于医学技术领域。在现有的神经成像技 术中，核磁共振MR和计算机断层扫描CT等因时间分辨率低，无法捕捉实时的神经异常连 接，而W脑磁图（magnetoenc巧halography，MEG)为代表的脑功能成像技术，因时间分辨率 高，成为辅助医学和认知神经学研究的重要工具。动态脑磁源定位是脑磁功能成像技术的 关键技术问题。
[0003] 动态脑磁源定位研究的是，当采集待测人头部的，例如静息状态下，由脑内发出的 并被脑外传感器检测到的MEG信号，反推脑内神经源信号的位置和强度，进一步地探索脑 神经连接是否异常。根据测量的脑外MEG信号反演脑内神经信号的活动位置和强度被称为 MEG逆问题。
[0004] 根据头脑表观的磁场强度反演定位磁源的空间活动位置时脑磁研究中的一个重 要问题，其本质上是一个非线性优化逆问题，为了简化计算的复杂性，在脑磁源的反演定位 中，常用一线性方法去逼近非线性问题。现有技术中具体是采用脑磁源成像技术去进行脑 磁源定位的，传统最小范数估计法（MNE，minimumnormestimate),是最具代表性的源成像 方法，主要原理是，基于大脑在特定时刻只有局部神经元活动的前提，对欠定线性方程增加 12范数约束，求解一副能量最小的电流密度分布图像，常用的方式为：
[0005] 假设脑外有m个通道的MEG信号，脑内有n个均匀分布的源信号，那么在i时刻， 脑内源信号与MEG信号的关系可W用W下离散化的线性模型表示：
[0006]bi=Axi+e；
[0007] 其中，bi为第i时刻大小为mX1的MEG测量信号；X1为第i时刻脑内源信号，大小 为nXl是第i时刻和H同维度的噪声信号；A为引导场矩阵，代表脑内源信号与MEG测 量信号的映射关系，大小为mXn。当矩阵A已知时，即可由bi求出脑磁逆问题的解Xi。但 矩阵A是病态的，其条件数，即最大特征值与最小特征值之比很大，因此直接求逆不合适， 通常转化为求解最小二次泛函的问题。不难得知，上式最小二乘解对噪声非常敏感，MEG 测量信号中很小的噪声将对解产生很大的扰动，造成无用解。引入Ti化onov正则化技术来 减小噪声对对解的影响。在i时刻，脑磁逆问题求解转化为求解下式的最小值问题：
[0008]
[0009] 等式右边第一项表示测量数据和估计数据的拟合，第二项为正则项，表示解的先 验信息，其中R为约束解空间的正则算子，当m个信号通道噪声均匀一致时，R取单位矩阵 I，入为正则化参数，调节拟合项和正则项在两项之间达到平衡。上式对应的解的形式为：
[0010]
[0011] 由此可见，代价函数估算出来的源信号各个时刻之间是相互独立的。传统M肥方 法仅适用于静态逆问题，因此不适用于研究脑磁动态逆问题。
[0012] 为此，目前现有的一种研究脑磁动态逆问题方法是引入时域平滑约束项构造双正 则化，该方法的不足之处是：
[0013] 1、此方法的前提是基于相邻时刻的信号变化最小的假设，运样的处理虽然能得到 时域上平滑的估计解，但对于突变脑磁源信号（比如20到40毫秒的异常放电现象，时域上 不是平滑的）运样的测量信号是不适用的；
[0014] 2、此方法由于引入了两个正则项和正则化参数，增加了计算难度和计算量；
[0015] 3、此方法需要人为选定时长t，要依据经验才能求解，可重复性差。
【发明内容】
[0016] 针对上述技术问题，本发明中提出了一种动态脑磁源定位方法，该方法不同于传 统最小范数估计算法（minimumnormestimate,MNE)，该方法为用于解决动态脑磁源定位 的方法，其为通过在时域子空间使用最小L2范数解的来进行脑磁源定位方法。
[0017] 该方法的实施可用于在非侵入性研究应用中实施W检测神经元活动和异常神经 元网络的轨迹。例如，可用于多个脑磁源之间动态传递的脑网络。本发明解决了动态脑磁 源不易定位的技术问题，特别是解决了对突变脑磁源信号传递过程不易研究的技术问题。
[0018] 为了实现根据本发明的运些目的和其它优点，提供了一种动态脑磁源定位方法， 包括W下步骤：
[0019] 步骤1)采集待测试人头部的脑MR结构图像；在所述头部上布置m个磁信号采集 传感器，且在某一时段内同步、连续采集t个时刻的磁信号，生成脑磁图MEG信号B;
[0020] 步骤。记录所述头部的参考点位置，并与所述脑MR结构图像上的相对应位置 进行空间配准，形成统一的测量空间，记录每一个所述传感器在所述测量空间中的位置信 息；
[OOW 步骤扣对所述脑MR结构图像进行分割，保留大脑皮层，将大脑皮层区域均匀划分 成n个位置确定的网格，且n远大于m，每个网格上分布有一个等效电流偶极子，即n个脑磁 源信号，源信号在大脑皮层形成的空间为源空间；
[0022] 步骤4)确定所述测量空间和所述源空间的空间转换关系，用引导场矩阵A表示， 根据测得的脑磁图MEG信号B来反演计算源空间中对应的n个脑磁源信号的强度，用矩阵X 表示每一个所述脑磁源信号的位置信息及其每个时刻脑磁源信号的强度信息，数学模型如 下所示：
[0023] B = AX巧 1)
[0024] 其中，E为噪声信号；
[00巧]步骤5)将脑磁图MEG信号B分解到由第一信号Ui构成的第一空间，其表征了MEG信号B的电磁场域；同时将脑磁图MEG信号B分解到由第二信号Ve构成的第二空间，其表 征了MEG信号B的时间序列信息，其中，B=UiV，所述第二信号Ve由t个列向量构成；
[0026] 步骤6)取第二信号Ve的前k个列向量构成第S空间V，其表征了MEG信号B的时 域子空间，其中，k在1~t之间取值；将MEG信号B和源信号X分别投影到所述第S空间 V上，在所述第S空间V中求解MEG信号B和源信号矩阵X的关系式，计算得出源信号矩阵 X，提取源信号X中的位置信息及其强度信息，即完成了动态的脑磁源定位过程。
[0027] 优选的，所述步骤1)中，通过核磁共振仪采集生成所述脑MR结构图像，通过脑磁 仪采集所述MEG信号，所述脑磁仪包括头位置指示器和m个磁信号采集传感器，所述传感器 均匀布置在所述头部头皮上采集相应位置处的磁信号强度，所述头位置指示器用于记录所 述头部的参考点位置并且记录每一个所述传感器在所述头部的位置。
[0028] 优选的，所述步骤2)中，所述头部参考点位置为双侧耳前点和鼻根，将所述头部 的双侧耳前点和鼻根位置与所述脑MR结构图像中的双侧耳前点和鼻根位置进行空间位置 配准。
[0029] 优选的，所述脑磁仪采集的MEG信号为由所述待测试人头部脑内的所述脑磁源产 生并传递到脑外的磁场强度信号。
[0030] 优选的，所述脑磁仪采集到的数据经过去眼电、滤波W及基线校准后得到所述MEG 信号。
[0031] 优选的，所述步骤4)中，采用边界元或有限元方法结合所述测量空间求解正问题 获取所述引导场矩阵A。
[0032] 优选的，所述步骤5)中，对MEG信号B进行奇异值分解：
[0033] B = Ub A bVb" 2)
[0034] 其中，Ub为由t个m维列向量构成的B的列空间，AB为奇异值矩阵，UbAb=Ui， 即所述第一空间，Ve为由t个t维的列向量构成的B的行空间，即所述第二空间。
[0035] 优选的，所述步骤6)中，将MEG信号B和源信号矩阵X投影到所述第S空间V，数 学模型如下：
[0036]BV=AXV巧V
[0037]Btemp=AXtemp+C 3) 阳03引其中，Btcmp为meg信号B在第立空间V上的投影，Xtemp为源信号矩阵X在第立空 间V上的动态脑磁源的解，C为噪声信号，引入Ti化onov正则化求式3)的最小L2范数的 解，即求Xtpmp的估计值
其中I和Ik均为单位矩阵，
正则化参数人由kcurve准则选取。