CN113935376A - 一种基于联合约束典型相关分析的大脑功能子区划分方法 - Google Patents

Abstract

本发明公开了一种基于联合约束典型相关分析算法的大脑功能子区划分方法，其步骤包括：1、将收集到的功能磁共振图像进行预处理，提取待划分脑区和参考脑区的时序信号，构建数据集；2、通过本发明提出的联合约束典型相关分析算法，计算出每个子区域对应的权重向量以及参考脑区的权重向量；3、根据子区域对应的权重向量得到待划分脑区中每个体素的隶属度，并根据参考脑区的权重向量得到与该子区域在功能上密切相关的脑区；4、选取每个权重向量中较大权值对应的体素构成一个子区域，作为划分结果。本发明能同时处理异质性较强的多个受试的功能磁共振成像数据，解决了多个体大脑功能子区的联合划分问题，有助于进一步研究大脑的功能。

Description

一种基于联合约束典型相关分析的大脑功能子区划分方法

技术领域

本发明涉及医学信号处理领域，尤其涉及大脑功能区域划分问题。

背景技术

大脑功能连接网络展示了大脑各个区域之间的相互作用，能够提高对大脑的认识，促进对大脑功能的研究和理解。前期的研究表明，大脑连接网络特征可以辅助研究神经系统疾病，例如抑郁症、阿尔茨海默症和帕金森病等。大脑连接网络可以通过脑电图(EEG)、脑磁图(MEG)和功能磁共振成像(functional magnetic resonance imaging，fMRI)等现代神经成像技术进行估计。其中，fMRI因其较高的空间分辨率和无创性而被广泛应用。大脑连接网络构建的关键步骤之一是网络节点的定义。如果在体素水平上选择节点，会导致网络的维数极高，而且容易受到噪声的干扰，使得后续的分析非常困难。更有效的方法是在脑区(region ofinterest，ROI)水平研究大脑连接网络，其中每个ROI对应于网络的一个节点，通常使用区域内所有体素的平均信号代表整个区域的时序活动。

传统的大脑区域划分方法通过解剖特征获得大脑分区，如AAL图谱和Brodmann图谱。虽然它们在描述大脑的结构方面具有重要意义，但却不能保证脑区功能的同质性。此外，这些图谱通常是在少数样本上生成的，无法刻画出个体间的差异，并且将个体数据配准到标准模板时，也会引入一定的偏差，导致后续分析的不准确。

近年来，基于fMRI数据的大脑连接特征进行功能脑区划分的方法有潜力克服传统方法的不足，因此得到了越来越多的关注。由于这类方法直接从数据中学习连接特征并进行区域分割，所以能产生高度同质和功能相似的大脑区域。其中最常见的方法包括聚类算法，如层次聚类、kmeans和谱聚类等，以及盲源分离算法，如独立成分分析(independentcomponent analysis，ICA)、独立向量分析(independent vector analysis，IVA)和典型相关分析(canonical correlation analysis，CCA)等。

然而，目前大部分的方法仍然面临以下问题。首先，在对特殊人群，如衰老和疾病人群进行fMRI分析时，受试者之间的大脑区域活动和功能分布的异质性往往被忽视。有些方法将不同个体的数据配准到标准空间的公共模板上进行处理，而较大的个体差异性会导致后续分析的偏差，不能真正反映个体的特征。另一些方法则在个体水平上分别对每个受试的数据进行处理，但由于个体间脑区形状和大小的差异，个体化的分割会导致功能子区对齐操作的困难，阻碍后续的群组分析，难以寻找群组的共性特征。此外，传统的划分方法要求脑区边界清晰，分区之间没有重叠，称为严格分割(hard-parcellation)。由于在实际情况下，一个体素可能参与不同子区域的功能活动，即不同状态下隶属于不同的功能子区。因此，无明确边界的模糊分割(soft-parcellation)可能更适用于大脑功能区域的定义。在这种情况下，对于每个子区域，所有体素都有一个概率分布，表明该体素在多大程度上属于这个子区域，并依据体素权重的分布，选择功能子区内所包含的体素。但是模糊的边界不利于定义子区的范围，可能为后续的群组分析带来困难。

发明内容

本发明是为了解决上述现有技术存在的不足之处，提出一种基于联合约束典型相关分析(joint constrained CCA，JC-CCA)的大脑功能子区划分方法，以期能在个体水平上对异质性较强的多个个体进行功能子区的联合划分，从而有利于后续的群组分析，为研究大脑的功能提供新的方法。

本发明为达到上述发明目的，采用如下技术方案：

本发明一种基于联合约束典型相关分析的大脑功能子区划分方法的特点是按如下步骤进行：

步骤一：使用信号采集设备获取S个个体样本长度为T的功能磁共振成像信号并进行预处理后，得到预处理后的功能磁共振信号，并根据实验目标将S个个体预处理后的功能磁共振信号中待划分脑区信号记为X₁,X₂,...,X_i,...,X_S，其中，X_i表示第i个个体的待划分脑区信号，且第i个个体的待划分脑区信号X_i的特征维数为P_i；

根据先验知识，针对第i个个体的待划分脑区选择Q个其它大脑区域作为参考脑区，并在预处理后的功能磁共振信号中取每个参考脑区内所有体素的平均信号代表对应参考脑区内的时序活动，从而得到S个个体参考脑区的信号数据集Y₁,Y₂,...,Y_i,...,Y_S，其中，Y_i表示第i个个体的参考脑区的时序活动，且其特征维数为Q，从而构建样本数据集{{X_i,Y_i}|i＝1,2,...,S}；其中，{X_i,Y_i}表示第i个样本的数据；

步骤二：利用联合约束典型相关分析算法对样本数据集{{X_i,Y_i}|i＝1,2,...,S}进行计算，得到对应的K组权重向量

其中，

表示第i个样本{X_i,Y_i}的第k组权重向量；

表示第i个个体待划分脑区信号X_i的第k组权重向量；v^(k)表示每个个体的参考脑区时序活动的第k组权重向量；

步骤三：对第i个个体的第k个子区域，根据第k组权重向量

中每个分量的大小得到对应体素的隶属度；

将第k组权重向量v^(k)中每个分量的大小作为第k个子区域在功能上与对应的参考脑区的相关程度；

步骤四：设置一个阈值，将第k组权重向量

中大于所述阈值的分量所对应的体素划分为第i个个体待划分脑区中的第k个子区域，从而得到第i个个体待划分脑区中的K个有明确边界的子区域，进而得到S个个体待划分脑区中的K个子区域。

本发明所述的基于联合约束典型相关分析的大脑功能子区划分方法，其特点是，步骤二中联合约束典型相关分析算法包括：

步骤2.1：利用式(1)构建损失函数f：

式(1)中，T表示转置，

为第i个个体待划分脑区信号X_i与参考脑区的时序活动Y_i之间的相关系数，

为S个个体的相关系数之和，||·||₂、||·||₁分别为l₂范数和l₁范数，α₁,α₂,β₁,β₂,γ₁,γ₂为六个调节参数，L_xi为包含第i个个体待划分脑区空间信息的拉普拉斯矩阵，L_y为包含每个个体参考脑区空间信息的拉普拉斯矩阵；

步骤2.2：将样本数据集{{X_i,Y_i}|i＝1,2,...,S}标准化为零均值单位方差的数据集，并根据第i个个体待划分脑区空间信息中体素的欧氏距离生成拉普拉斯矩阵L_xi，根据Q个参考脑区空间信息中的中心体素之间的欧氏距离生成拉普拉斯矩阵L_y，使用交叉验证法选择使所述相关系数之和

最大的六个调节参数；

初始化k＝1；

步骤2.3：随机初始化第k组权重向量

并利用式(2)所示的梯度下降法的解，对所述损失函数f进行交替优化迭代，直到损失函数f收敛，从而得到更新后的第k组权重向量并赋值给

式(2)中，

为第i个个体的第k组权重向量

对应的对角矩阵，其对角元素为向量

为第k组权重向量v^(k)的对角矩阵，其对角元为向量

步骤2.4：将k+1赋值给k后，判断k＞K是否成立，若成立，则表示得到K组权重向量

否则，利用式(3)更新数据集{{X_i,Y_i}|i＝1,2,...,S}后，返回步骤2.3顺序执行；

式(3)中，←表示赋值符号，

表示第i个个体待划分脑区信号X_i的第k-1组权重向量；v^(k-1)表示每个个体的参考脑区时序活动的第k-1组权重向量。

与现有的脑区划分方法相比，本发明的优势体现在：

1、本发明步骤二中，将传统CCA算法拓展到多组数据集的情形，并结合领域知识设置合适的约束条件，提出了新的JC-CCA算法。该算法通过约束参考脑区对应的权重向量相同，并加入稀疏和光滑约束，充分整合了每个个体的空间信息，能够更加有效提取待划分区域的个体特征和参考脑区的群体共性。并且能够在个体水平上同时处理多受试的数据，允许各数据集的维数不同，完成个性化的脑区分割，充分考虑了个体的异质性，避免了在公共模板处理数据时需要进行配准而造成偏差。

2、本发明步骤三和步骤四中同时进行了两类分割方法，即通过权重向量得到了模糊分割、严格分割，以及与子区域有功能连接的参考脑区和对应的相关系数，其中相关系数可以反映功能连接的强度。本发明结合了两类分割方法的优势，有助于后续的群组分析，便于探究不同人群脑区分布和功能上的差异，弥补了传统基于连接特征的大脑功能区域划分方法的不足，对研究大脑的功能和结构具有重要意义。

附图说明

图1为本发明方法流程图；

图2为健康年轻人和老年人的结果比较图；

图3为实验结果与临床数据的相关性分析图。

具体实施方式

本实施例中，一种基于联合约束典型相关分析(joint constrained CCA，JC-CCA)的大脑功能子区划分方法是先将收集到的fMRI信号进行预处理，提取待划分脑区的和参考脑区的数据，构建数据集；再通过本发明提出的JC-CCA算法，计算出每个子区域对应的权重向量以及参考脑区的权重向量；然后根据子区域对应的权重向量得到待划分脑区中每个体素的隶属度，并根据参考脑区的权重向量得到与该子区在功能上密切相关的脑区；最后选取每个权重向量中较大权值对应的体素构成一个子区域，作为划分结果。具体的说，如图1所示，按如下步骤进行：

步骤一：实验使用信号采集设备获取28名健康的年轻人(HC组)和30名有认知衰退症状的老年人(PD组)的S＝58个个体，样本长度为T＝230的fMRI信号，并进行刚体矫正、去趋势等预处理后，得到预处理后的功能磁共振信号。实验希望划分的大脑区域为重要的大脑功能区：左壳核(leftputamen)和右壳核(rightputamen)区域，其余的Q＝46个大脑区域将作为参考脑区用来研究和壳核区域的连接关系。由于左右壳核之间存在一定的差异，所以实验将分别对左右壳核进行划分，以下内容以左壳核为例。将S＝58个个体预处理后的fMRI信号中左壳核区域的信号记为X₁,X₂,...,X_i,...,X_S，其中，X_i表示第i个个体的左壳核区域的信号，且第i个个体的左壳核区域的信号X_i的特征维数为P_i；并在预处理后的fMRI信号中取每个参考脑区内所有体素的平均信号代表对应参考脑区内的时序活动，从而得到S＝58个个体参考脑区的信号数据集Y₁,Y₂,...,Y_i,...,Y_S，其中，Y_i表示第i个个体的参考脑区的时序活动，且其特征维数为Q＝46，从而构建样本数据集{{X_i,Y_i}|i＝1,2,...,S}；其中，{X_i,Y_i}表示第i个个体的数据；

步骤二：利用JC-CCA算法对样本数据集{{X_i,Y_i}|i＝1,2,...,S}进行计算，得到K组权重向量

其中，

表示第i个样本{X_i,Y_i}的第k组权重向量；

表示第i个个体左壳核信号X_i的第k组权重向量；v^(k)表示每个个体的参考脑区时序活动的第k组权重向量；

步骤2.1：利用式(1)构建损失函数f：

式(1)中，T表示转置，

为第i个个体待划分脑区信号X_i与参考脑区的时序活动Y_i之间的相关系数，

为S个个体的相关系数之和，||·||₂、||·||₁分别为l₂范数和l₁范数，α₁,α₂,β₁,β₂,γ₁,γ₂为六个调节参数，L_xi为包含第i个个体待划分脑区空间信息的拉普拉斯矩阵，L_y为包含每个个体参考脑区空间信息的拉普拉斯矩阵；具体实施中，调节参数α₁，α₂用来控制解的唯一性，不影响结果，所以均设置为0.1。而参考脑区之间不需要空间光滑性约束，所以设置对应的空间光滑性调节参数γ₂＝0，且不需要计算对应的拉普拉斯矩阵L_y。左壳核区域对应的拉普拉斯矩阵L_xi由第i个个体左壳核中体素之间的欧氏距离生成。

步骤2.2：将样本数据集{{X_i,Y_i}|i＝1,2,...,S}标准化为零均值单位方差的数据集，并根据第i个个体待划分脑区空间信息中体素的欧氏距离生成拉普拉斯矩阵L_xi，根据Q个参考脑区空间信息中的中心体素之间的欧氏距离生成拉普拉斯矩阵L_y，使用5折交叉验证法选择使相关系数之和

最大的3个调节参数β₁,β₂,γ₁，并使用选出的调节参数进行迭代计算；

初始化k＝1；根据相关研究，将左壳核分为3个子区域较为合适，所以设置K＝3。

步骤2.3：随机初始化第k组权重向量

并利用式(2)所示的梯度下降法的解，对损失函数f进行交替优化迭代，直到损失函数f收敛，从而得到更新后的第k组权重向量并赋值给

式(2)中，

为第i个个体的第k组权重向量

对应的对角矩阵，其对角元素为向量

为第k组权重向量v^(k)的对角矩阵，其对角元为向量

步骤2.4：将k+1赋值给k后，判断k＞K是否成立，若成立，则表示得到K组权重向量

否则，利用式(3)更新数据集{{X_i,Y_i}|i＝1,2,...,S}后，返回步骤2.3顺序执行；

式(3)中，←表示赋值符号，

表示第i个个体待划分脑区信号X_i的第k-1组权重向量；v^(k-1)表示每个个体的参考脑区时序活动的第k-1组权重向量。

步骤三：对第i个个体的第k个子区域，根据第k组权重向量

中每个分量的大小得到对应体素的隶属度，权值越大，对应体素属于该子区域的概率就越高，作为每个个体模糊分割的结果；

将第k组权重向量v^(k)中每个分量的大小作为第k个子区在功能上与对应的参考脑区的相关程度，这里由于算法的设计，控制了参考脑区权重向量相同，解决了在个体水平上分割时的对齐问题，同时可以直接计算数据集的相关系数来衡量连接关系的强度，便于后续的群组分析；

步骤四：根据相关研究，选取最大的前30％的体素作为一个子区域较为合适，于是将第k组权重向量

中前30％的分量所对应的体素划分为第i个个体待划分脑区中的第k个子区域，从而得到第i个个体待划分脑区中的K个有明确边界的子区域，进而得到S个个体待划分脑区中的K个子区域，作为严格分割的结果。同时为了进行后续的分析，需要计算每个子区域内体素的同质性，使用子区域内体素fMRI信号的第一主成分所解释的方差比例来衡量。

同样地，对右壳核重复步骤一至四得到右壳核的脑区划分结果。

进一步，为了探究不同人群大脑功能的差异，对两组受试在组水平上进行两样本t检验，分别比较功能连接强度(connectivity strength)和子区内体素的同质性(homogeneity)。得到两组人群有显著差异的结果如图2，H和P分别代表年轻人和老年人，数字代表子区域的序号，L和R分别代表左侧和右侧，‘*’代表检验的P值小于0.05，‘**’代表检验的P值小于0.01。可以得出在左侧第一个子区和右侧第一个子区，年轻组的功能连接强度显著高于老年组；在左侧第一个子区和右侧第三个子区，年轻组的体素同质性显著高于老年组。这样的结果符合预期，体现了老年人大脑中部分体素的功能连接下降且同质性降低，导致大脑功能出现衰退。

在此基础上，结合大脑功能衰退的临床评分(UPDRS-III)数据，分析老年组的结果与大脑功能衰退程度的相关性(ρ)。如图3所示，(a)(b)为左侧和右侧第一个子区域的功能连接强度与UPDRS-III的相关性，(c)(d)为左侧第一个子区域和右侧第三个子区域的同质性与UPDRS-III的相关性，P值均小于0.05。这体现出本发明通过功能连接特性划分脑区，有效地提取出了一种大脑功能衰退的特征，为大脑功能的相关研究提供了有力的支持。

综上所述，本发明能够解决异质性较强的多个体大脑功能子区联合划分问题，既能整合所有个体空间中的信息，同时能有效提取公共特征，避免了配准时可能出现的偏差。这为处理老年人及患病者等异质性较强的群体的fMRI数据提供了一种有效的手段，对促进大脑功能的研究具有重要的意义。

Claims (2)

1.一种基于联合约束典型相关分析的大脑功能子区划分方法，其特征是按如下步骤进行：

步骤一：使用信号采集设备获取S个个体样本长度为T的功能磁共振成像信号并进行预处理后，得到预处理后的功能磁共振信号，并根据实验目标将S个个体预处理后的功能磁共振信号中待划分脑区信号记为X₁,X₂,...,X_i,...,X_S，其中，X_i表示第i个个体的待划分脑区信号，且第i个个体的待划分脑区信号X_i的特征维数为P_i；

根据先验知识，针对第i个个体的待划分脑区选择Q个其它大脑区域作为参考脑区，并在预处理后的功能磁共振信号中取每个参考脑区内所有体素的平均信号代表对应参考脑区内的时序活动，从而得到S个个体参考脑区的信号数据集Y₁,Y₂,...,Y_i,...,Y_S，其中，Y_i表示第i个个体的参考脑区的时序活动，且其特征维数为Q，从而构建样本数据集{{X_i,Y_i}|i＝1,2,...,S}；其中，{X_i,Y_i}表示第i个样本的数据；

步骤二：利用联合约束典型相关分析算法对样本数据集{{X_i,Y_i}|i＝1,2,...,S}进行计算，得到对应的K组权重向量

其中，

表示第i个样本{X_i,Y_i}的第k组权重向量；

表示第i个个体待划分脑区信号X_i的第k组权重向量；v^(k)表示每个个体的参考脑区时序活动的第k组权重向量；

步骤三：对第i个个体的第k个子区域，根据第k组权重向量

中每个分量的大小得到对应体素的隶属度；

将第k组权重向量v^(k)中每个分量的大小作为第k个子区域在功能上与对应的参考脑区的相关程度；

步骤四：设置一个阈值，将第k组权重向量

中大于所述阈值的分量所对应的体素划分为第i个个体待划分脑区中的第k个子区域，从而得到第i个个体待划分脑区中的K个有明确边界的子区域，进而得到S个个体待划分脑区中的K个子区域。

2.权利要求1所述的基于联合约束典型相关分析的大脑功能子区划分方法，其特征是，步骤二中联合约束典型相关分析算法包括：

步骤2.1：利用式(1)构建损失函数f：

式(1)中，T表示转置，

为第i个个体待划分脑区信号X_i与参考脑区的时序活动Y_i之间的相关系数，

为S个个体的相关系数之和，||·||₂、||·||₁分别为l₂范数和l₁范数，α₁,α₂,β₁,β₂,γ₁,γ₂为六个调节参数，L_xi为包含第i个个体待划分脑区空间信息的拉普拉斯矩阵，L_y为包含每个个体参考脑区空间信息的拉普拉斯矩阵；

步骤2.2：将样本数据集{{X_i,Y_i}|i＝1,2,...,S}标准化为零均值单位方差的数据集，并根据第i个个体待划分脑区空间信息中体素的欧氏距离生成拉普拉斯矩阵L_xi，根据Q个参考脑区空间信息中的中心体素之间的欧氏距离生成拉普拉斯矩阵L_y，使用交叉验证法选择使所述相关系数之和

最大的六个调节参数；

初始化k＝1；

步骤2.3：随机初始化第k组权重向量

并利用式(2)所示的梯度下降法的解，对所述损失函数f进行交替优化迭代，直到损失函数f收敛，从而得到更新后的第k组权重向量并赋值给

式(2)中，

为第i个个体的第k组权重向量

对应的对角矩阵，其对角元素为向量

为第k组权重向量v^(k)的对角矩阵，其对角元为向量

步骤2.4：将k+1赋值给k后，判断k＞K是否成立，若成立，则表示得到K组权重向量

否则，利用式(3)更新数据集{{X_i,Y_i}|i＝1,2,...,S}后，返回步骤2.3顺序执行；

式(3)中，←表示赋值符号，

表示第i个个体待划分脑区信号X_i的第k-1组权重向量；v^(k-1)表示每个个体的参考脑区时序活动的第k-1组权重向量。

Images