CN103345749A - 一种基于模态融合的大脑网络功能连接偏侧性检测方法 - Google Patents

Abstract

本发明涉及一种基于模态融合的大脑网络功能连接偏侧性检测方法，该方法包括以下步骤：首先，对fMRI图像进行预处理，并提取大脑网络左右半球内脑区的空间三维坐标信息；其次，基于fMRI图像提取的脑区的空间三维坐标信息，从经过预处理后的MEG数据中提取出对应的时间序列信息；然后，利用MEG数据提取的脑区的时间序列，进行脑区之间的功能连接度分析；最后，计算大脑网络脑区之间功能连接的偏侧性指数。本发明所述方法是一种有效的基于磁共振和脑磁图模态融合的大脑网络功能连接偏侧性检测方法，能够比传统的仅利用fMRI图像的检测方法更完整、全面地检测大脑网络功能连接的偏侧性。

Description

一种基于模态融合的大脑网络功能连接偏侧性检测方法

技术领域

本发明属于图像处理领域，具体涉及一种利用脑功能成像技术研究大脑网络功能连接偏侧性的检测方法。尤其涉及利用功能磁共振成像(fMRI)和脑磁图描记术(MEG)模态融合进行大脑网络功能连接偏侧性的检测。

背景技术

随着科学技术的进步，无创性脑功能成像技术发展日新月异，为我们用多种成像技术进行大脑认知研究提供了广阔的舞台。由于成像机理的不同，每一种成像技术都具有相应的优点和不足，因此没有一种成像模态可以适合所有实验研究或临床应用。截至目前，功能磁共振成像技术(functional Magnetic Resonance Imaging，fMRI)以其高空间分辨率、非侵入式等特点在神经疾病诊断治疗方面得到了广泛应用。fMRI技术是通过测量由神经活动引起的脑血流和脑血氧等成分变化而造成的磁共振信号变化来反应人脑内部活动的，其具有非常高的空间分辨率精度(毫米级)。不过，fMRI的时间分辨率取决于生理动力学响应速率。在神经活动中，突触传导为1毫秒级，信息传输是几百毫秒，但血流动力学反应的长潜伏期严重妨碍了血氧水平依赖信号(Blood Oxygenation Level DependentContrast，BOLD)对神经信号的响应。活动皮层BOLD信号的峰值出现在激励开始后的5～8秒，并且回到基线水平需要同样的时间。因此，由于在代谢水平上测量磁共振信号的变化相较神经元活动具有一定的时间延迟，因此，fMRI在时间分辨率上仅能达到秒级精度。近年来，脑磁图描记术(Magnetoencephalography，MEG)逐渐在神经影像领域得到应用。MEG技术可以实时测量脑内神经元之间信号传递时突触后电位变化所产生的磁场信号，因此能够准确测量毫秒级神经元的电活动，具有非常高的时间分辨率。但是由于该技术难以准确定位颅内信号源，其空间分辨率较低。因此，将fMRI与MEG各自的时空分辨率优势结合起来，是我们准确、全面地探索大脑神经元细胞信息加工处理机制的关键。

大脑左右半球存在差异已经在脑解剖学和功能学研究上得到了初步证实，左右半球的差异性研究对深入了解语言认知、视觉认知以及感觉功能的脑内加工机制方面有重要意义。大脑结构和功能偏侧性是有利于大脑实现高效、多任务功能活动的特性之一。然而，以往的关于大脑网络偏侧性的研究基本上仅采用功能磁共振成像技术测定大脑活动性区域，研究结果虽具有较高的空间分辨率，在时间维度上提供的信息却较少。本发明主要考虑利用fMRI空间分辨率的优势，确定大脑fMRI图像上的功能网络脑区的空间三维信息，继而在MEG数据中获得对应位置的高分辨率时间维度信息，基于提取的高分辨率时间维度信息，从而可以实现在高分辨率的时空二个维度下进行大脑网络功能连接偏侧性的精确检测。

发明内容

本发明的目的在于设计一种同时具备高时空分辨率优势的大脑网络偏侧性的检测方法。

为实现上述目的，本发明提出一种基于磁共振和脑磁图模态融合的大脑网络功能连接偏侧性检测方法，包括以下步骤：

采集实验对象静息态下的fMRI数据和MEG数据；

步骤Sa，fMRI图像预处理，并提取大脑网络左右半球内脑区的空间三维坐标信息；

步骤Sb，MEG数据预处理，并基于fMRI图像提取的脑区的空间三维坐标信息，从MEG数据中提取出对应的时间序列信息；

步骤Sc，利用MEG数据中提取的脑区的时间序列，进行脑区之间的功能连接度分析；

步骤Sd，计算大脑网络脑区之间功能连接的偏侧性指数。

本发明所述的基于磁共振和脑磁图模态融合的大脑网络功能连接偏侧性检测方法，能够有效结合fMRI的高空间分辨率优势和MEG的高时间分辨率优势。真实试验数据结果表明，相较于传统的仅利用fMRI单一模态的研究方法，本发明方法在时间维度上更完整地检测出了大脑网络的偏侧性。本发明方法为脑功能机理研究提供了一种全新的数据处理分析框架，在基于多模态融合尤其是fMRI和MEG融合的医学成像领域具有一定的应用价值。

附图说明

图1为本发明方法的流程示意图。

图2是使用本发明方法检测到的大脑默认网络脑区之间功能连接的偏侧性结果。

图3是仅使用fMRI单一模态检测到的大脑默认网络脑区之间功能连接的偏侧性结果。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实例，并参照附图，对本发明进一步详细说明。

本发明主要利用功能磁共振成像技术(fMRI)高空间分辨率的特点，定位fMRI图像上的大脑左右半球内脑区的空间三维坐标信息，继而从脑磁图描技术(MEG)数据中提取出上述脑区对应的高分辨率时间序列，从而可以在时空二维度内对大脑网络偏侧属性进行高分辨率精确检测。

图1是本发明所涉及的一种基于磁共振和脑磁图模态融合的大脑网络偏侧性检测方法的一个具体实施例，其具体涉及联合利用功能磁共振成像(fMRI)和脑磁图描记术(MEG)的大脑网络偏侧性的检测方法。该实施例的具体实施步骤如下：

首先需要采集实验对象静息态下的fMRI图像和MEG数据，作为后续分析方法的基础；

步骤Sa，fMRI图像预处理，并提取大脑网络左右半球内脑区的空间三维坐标信息；

1.对fMRI图像进行预处理

由于在磁共振成像的扫描过程中，被试的头动、图像内层间扫描时间的差别以及外加磁场的不均匀性等各种因素均会影响采集到的fMRI图像的质量，因此我们需要在保留脑功能图像细节的同时，进行适当的预处理，从而提高脑功能图像的信噪比。预处理步骤包括：切片扫描时间对齐、图像序列对齐、配准、标准化、平滑、滤波等。在本实施例中，对fMRI图像依次进行图像序列对齐、配准、标准化、空间平滑滤波。所属技术领域的技术人员也可以根据实际需要采用其他预处理方式，或者上述预处理步骤的任意组合。

2.提取大脑网络左右半球内脑区的空间三维坐标信息

本发明采用基于数据驱动的独立成分分析(Independent ComponentAnalysis，ICA)方法提取大脑各主要功能网络成分。其数学模型如下：设A∈R^N×M为信号混合矩阵，其中，R表示实数矩阵(下同)，N表示全脑体素的个数；x∈R^N×M为脑功能信号观测向量，即我们实际采集到的脑功能实验数据，其中，M表示数据扫描的时间点，c∈R^N×M为由M(N≥M)维统计独立的信号分量构成的信号向量，即人脑内部体素的激活情况。那么我们可以表征混合模型如下：x＝Ac。在脑功能影像学研究中，混合矩阵A是未知的，因此为了得到大脑内部体素的激活情况，即求解矩阵c，我们需要采取合适的处理方法。作为一种被广泛用于解决盲源信号分离的方法，ICA可以在混合矩阵A未知的情况下，寻找一种最优的线性映射w，实现从观测的脑功能信号中提取不能被直接观测的原始信号c(≈y)：y＝w^Tx＝w^TAc。使用基于最小描述长度(Minimum Description Length，MDL)准则估计独立成分数，继而采用基于互信息最大化的算法计算解混矩阵，分离源信号，获得大脑各主要功能网络成分及其脑区的Talairach坐标。

步骤Sb，MEG数据预处理，并基于fMRI图像提取的脑区的空间三维坐标信息，从MEG数据中提取出对应的时间序列信息；

1.MEG数据预处理

由于MEG数据在采集过程中容易受到环境噪声的影响，我们非常有必要对其进行预处理，以提高数据的信噪比；此外，由于MEG技术的高时间分辨率特性，采集的数据包含了大量的时间维度信息。因此，为保证计算精度与运算时间的平衡，预处理步骤主要包括：数据带通滤波，选用合适的下采样频率等。在本实施例中，我们对MEG数据依次进行了0.5Hz-48Hz的带通滤波和频率为300Hz的下采样。所属技术领域的技术人员也可以根据研究目的或者实验条件的实际需要采用其他预处理步骤，或者上述预处理步骤的任意组合。

2.基于fMRI图像提取的脑区的空间三维坐标信息，从MEG数据中提取出对应的时间序列信息

脑磁图(MEG)测量法是一种对大脑完全无侵袭，无损害的以测量脑磁场来达到脑功能区域定位以及评价被定位区域脑功能的健康状态的测量方法。然而，如何求解脑磁逆问题，即利用扫描的磁场数据来推断脑内电流源的分布，成为了脑磁图研究中具有挑战性的问题。已有研究证明了导体外的磁场数据无法唯一确定导体内的电流分布，因此从本质上讲，脑磁逆问题是不适定的，我们必须在满足条件的解集中通过施加一定限制条件找出合理的解。目前，脑磁逆问题的求解方法主要分为两大类：一类基于图像重建技术，其基本思想是将电流分布区域离散成网格，对网络节点处磁源进行重建，重建的磁源参数主要是电流的强度。由于磁场方程中位置参数是非线性的，电流强度参数是线性的，重建时人为地限定了源的位置，因此磁场方程退化为线性方程。另一类为偶极子定位法，该类方法不限定磁源的位置，通常的求解过程是先建立一个目标函数，然后通过调整源的参数使目标函数达到极小，即通过执行一个非线性优化过程来获得磁源参数的解。

本发明的实施例采用基于图像重建技术的最小范数解(MinimumNorm Solution，MNS)求解方法。假设磁场方程组为D＝LS，其中，D为测得的磁场信号值，S代表节点处电流，L为已知的增益矩阵，取决于网格节点及磁场探测点的位置与头的电导率分布。MNS方法的限制条件是找到满足方程组解中范数最小的S^*，即

S^*＝min{||S||}，其中，||·||一般取Frobenius2范数，即||S||＝(S^TS)^1/2，其中^T表示矩阵的转置。具体求解S^*时可以采用基于Moore-Penrose广义逆L⁺的方法，即

S^*＝L⁺D

对于广义逆L⁺的求解，本发明采用基于奇异值分解(Singular ValueDecomposition，SVD)的方法。设L为m×n阶实矩阵(m＜n)，其中，m为MEG扫描通道数，n为网格节点数。则对L进行奇异值分解可得

L＝U∑V^T

其中U＝(u₁，u₂，…，u_m)为m×m阶正交矩阵，V＝(v₁，v₂，…，v_n)为n×n阶正交矩阵，∑为m×n阶对角矩阵，即∑＝diag(δ₁，δ₂，…，δ_P)，其中P＝m，且δ₁≥δ₂≥…≥δ_P≥0，若矩阵L的秩rank(L)＝r，则上述序列中只有前r个值大于零，即δ₁≥δ₂≥…≥δ_r＞δ_r+1＝…＝δ_P＝0

那么奇异值分解可表述为：

其中，u_i(i＝1，2，…，m)称为奇异值分解的左矢量，v_i(i＝1，2，…，n)称为奇异值分解的右矢量。

由此可求得Moore-Penrose广义逆L⁺为：

$L^{+}=\underset{i=1}{\overset{r}{\mathrm{\Σ}}}\frac{1}{{\mathrm{\δ}}_i}{v}_i{u}_i^T$

由此可根据S^*＝L⁺D求出MEG数据所隐含的源空间信息。

通过步骤Sb，我们可根据ICA分析方法获得的局部最大激活值的体素三维Talairach坐标，从而获得该局部脑区的时间序列。

步骤Sc，利用MEG数据提取的脑区的时间序列，进行脑区之间的功能连接度分析；

人脑在完成高级功能时需要多个脑区之间的协同工作，对信息进行加工、融合。脑功能连接的概念最早出现在电生理研究中，它是描述脑区之间协同工作模式的有效手段之一。目前在脑功能影像学研究中，功能连接度分析方法已逐渐受到重视。绝大多数的功能连接分析，无论是基于任务的还是基于静息状态的，都采用相关分析方法来度量脑区之间的功能连接。在本发明方法中，我们采用了偏相关分析方法来度量大脑网络脑区之间的功能连接度。具体计算步骤如下：

首先计算脑区两两之间的皮尔逊相关系数，计算公式如下：

$r_{\mathrm{XY}}=\frac{\underset{i=1}{\overset{n}{\mathrm{\Σ}}}(X_i-\stackrel{\‾}{X})(Y_i-\stackrel{\‾}{Y})}{\sqrt{\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{(X_i-\stackrel{\‾}{X})}^2\sqrt{\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{(Y_i-\stackrel{\‾}{Y})}^2}}}$

其中，X_i表示脑区X在第i个时刻点的激活值，表示脑区X在整个扫描期间激活值的平均。Y_i表示脑区Y在第i个时刻点的激活值，

表示脑区Y在整个扫描期间激活值的平均。r_XY即可表示为脑区X和Y之间的皮尔逊相关系数。依据上式计算，从而得到皮尔逊相关系数矩阵：

$R={\left(r_{\mathrm{XY}}\right)}_{n\×n}=\left(\begin{array}{cccc}{r}_{11}& r_{12}& \·\·\·& r_{1n}\\ r_{21}& r_{22}& \·\·\·& r_{2n}\\ \·& \·& & \·\\ \·& \·& \·\·\·& \·\\ \·& \·& & \·\\ r_{n1}& r_{n2}& \·\·\·& r_{\mathrm{nn}}\end{array}\right)$

对该系数矩阵求逆，得到逆矩阵：

$R^{-1}={\left(c_{\mathrm{XY}}\right)}_{n\×n}=\left(\begin{array}{cccc}{c}_{11}& c_{12}& \·\·\·& c_{1n}\\ c_{21}& c_{22}& \·\·\·& c_{2n}\\ \·& \·& & \·\\ \·& \·& \·\·\·& \·\\ \·& \·& & \·\\ c_{n1}& c_{n2}& \·\·\·& c_{\mathrm{nn}}\end{array}\right)$

从而得到脑区两两之间的偏相关系数：

$p_{\mathrm{XY}}=-c_{\mathrm{XY}}/\sqrt{c_{\mathrm{XX}}{c}_{\mathrm{YY}}}$

由此可得到偏相关矩阵如下：

$P={\left(p_{\mathrm{XY}}\right)}_{n\×n}=\left(\begin{array}{cccc}{p}_{11}& p_{12}& \·\·\·& p_{1n}\\ p_{21}& p_{21}& \·\·\·& p_{2n}\\ \·& \·& & \·\\ \·& \·& \·\·\·& \·\\ \·& \·& & \·\\ p_{n1}& p_{n1}& \·\·\·& p_{\mathrm{nn}}\end{array}\right)$

步骤Sd，计算大脑网络脑区之间功能连接的偏侧性指数。

在计算偏侧性指数之前，我们需要指定种子点(seed region)和目标点(target region)。在这里，我们以A，B为例。A为种子点，B为目标点，左右半球各含一个种子点和目标点，分别表示为L_A，R_A，L_B，R_B，那么偏侧性指数可通过下式计算得到：

$\mathrm{LateralityIndex}=\frac{(L_A{L}_B-R_A{L}_B)-(R_A{R}_B-L_A{R}_B)}{\left|L_A{L}_B\right|+\left|L_A{R}_B\right|+\left|R_A{R}_B\right|+\left|R_A{L}_B\right|}$

其中，L_AL_B表示左侧半球内种子点A和目标点B之间的功能连接度，R_AL_B表示右侧半球种子点A和左侧半球目标点B之间的功能连接度，R_AR_B表示右侧半球内种子点A和目标点B之间的功能连接度，L_AR_B表示左侧半球种子点A和右侧半球目标点B之间的功能连接度，上述四值均可在步骤Sc中生成的偏相关系数矩阵P中找到两两脑区之间的偏相关系数。如果计算得到的LateralityIndex大于0，则我们认为种子点A和目标点B的功能连接是左偏的，反之，如果LateralityIndex小于0，则我们认为种子点A和目标点B的功能连接是右偏的。

本发明所述的基于磁共振和脑磁图模态融合的大脑网络功能连接偏侧性检测方法，可通过分析真实采集的实验数据得以说明：

(1)实验数据采集及分析

在实验中采用真实数据集作测试，一位健康男性被试参与了静息状态下的fMRI扫描和MEG扫描。

由于大脑系统包含功能网络众多，在本例中，我们以目前受到广泛研究关注的默认脑网络为例，所属技术领域的技术人员也可以根据实际需要选择其他脑功能网络，对网络内脑区之间的功能连接偏侧性进行分析研究。

fMRI扫描采用T2*加权梯度回波平面成像序列获取血氧水平依赖(Blood Oxygen Level Dependent，BOLD)数据。针对fMRI图像的预处理，我们采用统计参数图软件包(SPM5，http://www.fil.ion.ucl.ac.uk/spm/)对图像进行处理。对经过预处理后的fMRI图像，我们利用GIFT软件(http://icatb.sourceforge.net/)包含的独立成分分析方法提取出属于大脑默认网络的独立成分，从而获得局部最大激活值的Talairach坐标如下：左脑额内侧回(Left Medial Frontal Gyrus，LMFG，-2，42，-11)，左脑额上回(LeftSuperior Frontal Gyrus，LSFG，-2，56，23)，左脑前扣带回(Left AnteriorCingulate Cortex，LACC，-2，48，-2)，左脑后扣带回(Left Posterior CingulateCortex，LPCC，-2，-40，24)，左脑颞上回(Left Superior Temporal Gyrus，LSTG，-38，17，-19)，左脑角回(Left Angular Gyrus，LAG，-48，-68，37)，左脑顶下小叶(Left Inferior Parietal Lobule，LIPL，-48，-66，38)，左脑颞中回(Left Middle Temporal Gyrus，LMTG，-42，-81，19)，右脑额内侧回(RightMedial Frontal Gyrus，RMFG，4，46，-11)，右脑额上回(Right SuperiorFrontal Gyrus，RSFG，4，59，23)，右脑前扣带回(Right Anterior CingulateCortex，RACC，2，47，1)，右脑后扣带回(Right Posterior Cingulate Cortex，RPCC，2，-44，21)，右脑颞上回(Right Superior Temporal Gyrus，RSTG，38，17，-19)，右脑角回(Right Angular Gyrus，RAG，51，-63，31)，右脑顶下小叶(Right Inferior Parietal Lobule，RIPL，53，-56，38)，右脑颞中回(RightMiddle Temporal Gyrus，RMTG，61，-24，-4)。

脑磁图数据采集的采样频率设置为600Hz，利用加拿大CTF系统公司的CTF-151系统自带的综合三阶梯度的噪声消除技术对MEG信号进行降噪处理，对采集到的数据进行下采样(300Hz)并在此基础上进行脑磁图逆问题求解，利用前述16个极值点坐标提取MEG数据包含的对应脑区的时间序列。随后用本发明所述方法(A)结果与仅用fMRI单一模态方法(B)的大脑网络功能连接偏侧性的检测结果做对比。在计算过程中，我们设定偏侧性阈值为0.1，即如果LateralityIndex值小于-0.1，我们才认为脑区之间的功能连接是右偏的；如果LateralityIndex值大于0.1，我们才认为脑区之间的功能连接是左偏的。在本方法中，阈值的选择并不是一成不变的，所属技术领域的技术人员也可以根据实际需要采用其他阈值。

(2)数据分析结果

利用方法A在fMRI和MEG两种模态数据进行分析的结果如图2中表格所示，利用方法B在fMRI单一模态实验数据上进行分析的结果如图3中表格所示。

以上实验结果说明，本发明所述的基于磁共振和脑磁图模态融合的大脑网络偏侧性检测方法，可以有效地提高脑网络脑区之间的功能连接偏侧性检测性能，在时空二维上更完整、准确地反应大脑网络信息加工处理机制。

以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (7)

1.一种基于模态融合的大脑网络功能连接偏侧性检测方法，其特征在于，包括以下步骤：

采集实验对象静息态下的fMRI数据和MEG数据；

步骤Sa，fMRI图像预处理，并提取大脑网络左右半球内脑区的空间三维坐标信息；

步骤Sb，MEG数据预处理，并基于fMRI图像提取的脑区的空间三维坐标信息，从MEG数据中提取出对应的时间序列信息；

步骤Sc，利用MEG数据中提取的脑区的时间序列，进行脑区之间的功能连接度分析；

步骤Sd，计算大脑网络脑区之间功能连接的偏侧性指数。

2.如权利要求1所述的方法，其特征在于：

在步骤Sa中，采用基于数据驱动的独立成分分析(IndependentComponent Analysis，ICA)方法提取fMRI图像的功能网络脑区的空间三维坐标信息。

3.如权利要求2所述的方法，其特征在于：

所述的ICA方法包括：首先使用基于最小描述长度准则估计独立成分数，继而采用基于互信息最大化的算法计算解混矩阵，分离出源信号，从而得到大脑各主要功能网络成分及其脑区的Talairach坐标。

4.如权利要求1所述的方法，其特征在于：

在步骤Sb中，采用基于图像重建技术的最小范数解方法，对经过预处理的MEG数据进行脑磁逆问题求解，获得大脑网络各脑区的时间序列信息。

5.如权利要求1所述的方法，其特征在于：

在步骤Sc中，利用从MEG数据中提取的大脑网络脑区的时间序列，使用偏相关分析方法进行脑区之间的功能连接度分析。

6.如权利要求1所述的方法，其特征在于：

在步骤Sd中，通过计算各脑区之间的功能连接偏侧性指数，得到脑区之间功能连接的偏侧性分析结果。

7.如权利要求6所述的方法，其特征在于：

采用的偏侧性指数计算方法如下所示：

$\mathrm{LateralityIndex}=\frac{(L_A{L}_B-R_A{L}_B)-(R_A{R}_B-L_A{R}_B)}{\left|L_A{L}_B\right|+\left|L_A{R}_B\right|+\left|R_A{R}_B\right|+\left|R_A{L}_B\right|}$

其中，L_AL_B表示左侧半球内种子点A和目标点B之间的功能连接度，R_AL_B表示右侧半球种子点A和左侧半球目标点B之间的功能连接度，R_AR_B表示右侧半球内种子点A和目标点B之间的功能连接度，L_AR_B表示左侧半球种子点A和右侧半球目标点B之间的功能连接度；如果计算得到的Laterality Index大于0，则我们认为种子点A和目标点B的功能连接是左偏的，反之，如果Laterality Index小于0，则我们认为种子点A和目标点B的功能连接是右偏的。

Images