CN111259849A

CN111259849A - 一种功能近红外光谱成像静息态脑网络检测方法和装置

Info

Publication number: CN111259849A
Application number: CN202010075023.4A
Authority: CN
Inventors: 段炼
Original assignee: Shenzhen University
Current assignee: Shenzhen University
Priority date: 2020-01-22
Filing date: 2020-01-22
Publication date: 2020-06-09
Anticipated expiration: 2040-01-22
Also published as: CN111259849B

Abstract

本发明提供一种功能近红外光谱成像静息态脑网络检测方法和装置。该方法包括：对于个体水平，通过自动确定簇数进行聚类，将功能近红外光谱测量通道划分为不同的簇；对于个体水平的聚类结果进一步划分为不同的聚类，获得组水平的聚类结果；基于组水平的聚类结果中簇的一致性，挑选出构成静息态脑网络的簇集合。本发明能够自动计算簇数并确定显著静息态脑网络地形图，而无需任何人工干预，也不需要使用者具有丰富的经验，即能够提供更客观的检测结果。

Description

一种功能近红外光谱成像静息态脑网络检测方法和装置

技术领域

本发明涉及脑网络检测技术领域，尤其涉及一种功能近红外光谱成像静息态脑网络检测方法和装置。

背景技术

静息态脑网络(RSN)被认为是许多基本和关键的认知功能的基础。广泛的证据也表明RSN障碍与各种神经和精神疾病有关。近年来，利用功能近红外光谱(fNIRS)研究RSN越来越受到重视。fNIRS具有成本低、便携、舒适、对头部运动不敏感等独特优势，非常适合于静息态的研究，尤其适用于婴儿等特定的人群。

在fNIRS研究中，种子点相关法和独立成分分析(ICA)被广泛用于RSN检测。种子点相关法是最经典的RSN检测方法，该方法为了找到一个功能系统的RSN，首先利用结构和/或功能定位信息对该功能系统中的种子区域进行预定义。然后利用相关或回归计算每个通道与种子之间的静息态功能连接(RSFC)。那些在种子区表现出明显RSFC的大脑区域被确定为RSN。虽然种子相关方法得到了广泛的应用，但它也存在一些固有的不足，例如只考虑了大脑其他区域与种子区域之间的关系，而忽略了这些区域之间的相互作用，并且该方法的检测结果通常取决于种子选择。

独立成分分析(ICA)是盲源分离(BSS)算法的一种。它可以将fNIRS信号分解成一系列统计上独立的源分量。然后根据RSN的频率特性和空间分布特征提取RSN分量。作为一种数据驱动的方法，基于ICA的fNIRS-RSN检测不依赖于种子选择，可以将噪声从神经活动信号中分离出来。然而，ICA的成分选择过程通常是在人的监督下或依赖于人工的选择，这需要使用者具有丰富的经验，可能会对结果产生主观影响。

因此，需要对现有技术改进，提供不依赖先验信息和人为干预的fNIRS-RSN检测新方法，以提高实用价值。

发明内容

本发明的目的在于克服上述现有技术的缺陷，提供一种功能近红外光谱成像静息态脑网络检测方法和装置，能够在缺乏先验信息和人为干预的情况下有效检测静息态脑网络。

根据本发明的第一方面，提供一种功能近红外光谱成像静息态脑网络检测方法。该方法包括以下步骤：

对于个体水平，通过自动确定簇数进行聚类，将功能近红外光谱测量通道划分为不同的簇；

对于个体水平的聚类结果进一步划分为不同的聚类，获得组水平的聚类结果；

基于组水平的聚类结果中簇的一致性，挑选出构成静息态脑网络的簇集合。

在一个实施例中，对于个体水平，通过自动确定簇数进行聚类，将功能近红外光谱测量通道划分为不同的簇包括：

在个体水平，对于M个受试者中的第i个参与者，计算每对通道之间的Pearson相关系数，并构建静息态功能连接矩阵；

将所有通道划分为K_i个簇

并将每个聚类

表示为一个地形图：

在一个实施例中，对于个体水平的聚类结果进一步划分为不同的聚类，获得组水平的聚类结果包括：

通过计算{K_i}(i＝1,2,…,M)的众数来确定组水平簇数K_group；

根据每个聚类对应的地形图之间的相似度，将所有个体水平的地形图

划分为K_group个簇{ξ_i}(i＝1,2,…,K_group)。

在一个实施例中，基于组水平聚类结果中簇的一致性，挑选出构成静息态脑网络的簇集合包括：

计算每个ξ_i中所有个体水平地形图的平均得到群组水平一致性地形图：

检查ξ_i中的各个地形图在所有受试者中的一致性；

对于包含显著结果通道的ξ_i，将ξ_i中的显著通道集合构成一个静息态脑网络。

在一个实施例中，使用二项分布B(M，0.5)对ξ_i中每一个通道上的所有个体水平地形图进行检验，以检查ξ_i中的各地形图在所有受试者中的一致性。

在一个实施例中，簇数K_i是采用Ncut自动估计获得。

根据本发明的第二方面提供一种功能近红外光谱成像静息态脑网络检测装置。该装置包括：第一聚类单元，用于对于个体水平，通过自动确定簇数进行聚类，将功能近红外光谱测量通道划分为不同的簇；第二聚类单元，用于对于个体水平的聚类结果进一步划分为不同的聚类，获得组水平的聚类结果；静息态脑网络检测单元：用于基于组水平的聚类结果中簇的一致性，挑选出构成静息态脑网络的簇集合。

与现有技术相比，本发明的优点在于：能够自动计算簇数并确定显著静息态脑网络地形图，而无需任何人工干预，也不需要使用者具有丰富的经验；不会因种子选择或手动成分选择过程而产生偏差，能够提供更客观的检测结果。

附图说明

以下附图仅对本发明作示意性的说明和解释，并不用于限定本发明的范围，其中：

图1是根据本发明一个实施例的功能近红外光谱成像静息态脑网络检测方法的流程图；

图2是根据本发明一个实施例的模拟仿真的不同信噪比影响示意图；

图3是根据本发明一个实施例不同数量RSN的模拟数据的验证结果；

图4是根据本发明一个实施例的不同RSN体积的验证结果；

图5是不平衡的RSN体积条件的验证结果；

图6是根据本发明一个实施例的感觉运动RSN成分的组水平一致性图；

图7是根据本发明一个实施例的ROC曲线分析结果；

附图中，Goodness of fit-拟合优度；Simulated RSN-仿真的RSN；ClusterNumber-簇数量；Group Consistency Map-群组水平一致性地形图；Clustering-basedgroup consistency map-基于聚类的群组水平一致性地形图；Seed-based correlationt-map-种子相关t-map；Sensorimotor template-感觉运动模板；Seed correlation-种子相关；Spectral clustering-谱聚类；Specificity(False positive rate)-特异性(假阳性率)；Sensitivity(True positive rate)-灵敏度(真阳性率)；ROC curve-受试者工作特征曲线。

具体实施方式

为了使本发明的目的、技术方案、设计方法及优点更加清楚明了，以下结合附图通过具体实施例对本发明进一步详细说明。应当理解，此处所描述的具体实施例仅用于解释本发明，并不用于限定本发明。

在本文示出和讨论的所有例子中，任何具体值应被解释为仅仅是示例性的，而不是作为限制。因此，示例性实施例的其它例子可以具有不同的值。

对于相关领域普通技术人员已知的技术、方法和设备可能不作详细讨论，但在适当情况下，所述技术、方法和设备应当被视为说明书的一部分。

本发明根据不同区域之间静息态的神经活动相似性将大脑划分成簇，而在同一簇中划分的大脑区域可能构成一个静息态脑网络(RSN)。为此，提出了一种数据驱动的基于谱聚类的fNIRS-RSN检测的方法。谱聚类是一种基于谱图论的将数据划分为不同簇的前沿算法。与k-means等传统聚类算法相比，谱聚类不需要假设数据的任何分布或原型。

参见图1所示，本发明提供的功能近红外光谱成像静息态脑网络检测方法包括：步骤S110，对于个体水平，通过自动确定簇数进行聚类，将功能近红外光谱测量通道划分为不同的簇；步骤S120，对于个体水平的聚类结果进一步划分为不同的聚类，获得组水平的聚类结果；步骤S130，基于组水平的聚类结果中簇的一致性，挑选出构成静息态脑网络的簇集合。简言之，在本发明中，首先，在个体层次上，采用自动确定簇数的谱聚类方法，将所有的fNIRS测量通道划分为不同的簇；其次，将个体层次的聚类结果进一步划分为不同的聚类；最终具有高组水平一致性的簇被确定为静息态脑网络(RSN)簇。通过这种方式，能够在缺乏先验信息和人为干预的情况下检测RSN。在下文中，本发明提供的方法也称为基于谱聚类的功能近红外光谱成像静息态脑网络检测方法或简称为基于谱聚类的方法。

为清楚理解本发明的思想，下文将依次介绍谱聚类和Ncut算法、自动簇数目估计、静息态网络检测等。

1)、关于谱聚类和Ncut算法

谱聚类是一种基于图分割的算法。在本发明实施例中，对于静息的功能近红外光谱(fNIRS)数据，可以通过将n个测量通道作为图G的顶点：V＝{v₁,v₂,…,v_n}，每对通道(i,j)之间的静息功能连接(RSFC)作为边：E＝{e_ij|i,j＝1,2,…,n}。该图G可以由一个n×n权重矩阵W表示，其元素w_ij是通道i和j之间的RSFC值，当i＝j时w_ij＝0。图的度矩阵D是对角矩阵，其对角向量为d＝[d₁,d₂,…,d_n]^T，其中

G的拉普拉斯算子由L＝D-W给出。

Ncut是一种用于谱聚类的图分割方法。具体地，令(α，β)为G的一个划分(即α∪β＝V和

)，Ncut方法通过最小化损失函数来确定该划分：

其中：

在公式(1)中，

和

分别是集合α和β的体积。通过Ncut方法进行图分割，能够同时最小化簇之间的解离并最大化集群内的关联，因此可以避免分离出仅包含几个点的小集群。

为划分(α，β)引入指标向量x，其元素满足

则Ncut可以写成：

上述最小化方程式(4)等同于最小化瑞利熵：

在y为分段常数且y^Td＝0的条件下，可以通过求解一个标准的特征值问题来实现。即求解：

其中

是归一化的拉普拉斯算子，λ表示特征值。在二分类情况下，通过划分

的最小非平凡特征向量来完成分类。对于k(k为大于2的整数)分类，首先令

为由

的k个最小特征值所对应的k个特征向量组成的矩阵，然后将U的每一行视为一个点，并使用k-means对这些点进行聚类。最后，按照U中相应行的聚类结果将V中的原始点进行划分。

2)、自动簇数目估计

原始的Ncut方法需要将簇数作为输入参数。在本发明实施例中，采用以下方法确定簇数：

假设数据本质上具有m个簇，并且

是由特征求解器生成的k个特征向量组成的矩阵，令z_i是U的第i行的向量。当k＝m时，有：

否则，如果k<m或k>m，则当v_i和v_j属于同一簇时，cos(z_i,z_j)偏离1，或者当v_i和v_j属于不同簇时，cos(z_i,z_j)偏离0。利用这种方式，可以寻找使cos(z_i,z_j)在相应条件下最大集中接近于1或0的k，作为固有簇数m的估计。在实际计算中，可以让k从2增加到n，并为每个k生成

令z_i是U的第i行的向量，计算矩阵C，其元素c_ij＝cos(z_i,z_j)。然后，将用很小的正数的阈值ε(例如ε＝0.05)对C二值化得到B＝b_ij：

令s为矩阵B所有元素的总和：

s＝∑_i,jb_ij(9)

让k从2开始递增，满足s_k≤s_k-1and s_k≤s_k+1的最小的k即为簇数m的估计。

3)、静息态网络检测

本发明实施例提供的功能近红外光谱成像静息态脑网络检测方法，使用两步聚类过程从fNIRS数据中检测RSN。对于每个受试者，首先使用上述自动确定簇数的Ncut方法，根据通道之间的RSFC对fNIRS测量区域进行个体水平的划分。其次，进一步对所有受试者的聚类结果进行组水平聚类。具有高组一致性的那些组水平簇被确定为RSN簇。

具体地，在个体水平，对于M个受试者中的第i个参与者，首先计算每对通道之间的Pearson相关系数，以构建RSFC矩阵。然后，使用Ncut将所有通道划分为K_i个簇

其中簇数K_i是自动估计的并且是受试者特异的。每个聚类

可以表示为一个地形图：

在组水平中，首先通过计算{K_i}(i＝1,2,…,M)的众数来确定组水平簇数K_group。然后，根据地形图之间的相似度，将所有个体水平的地形图

划分为K_group个簇{ξ_i}(i＝1,2,…,K_group)。然后，如果ξ_i中有多个地形图来自同一参与者，则将它们视为同一聚类的子单元并将其合并。通过计算每个ξ_i中所有个体水平地形图的平均得到群组水平一致性地形图。

为了确定ξ_i是否表征RSN，可以检查ξ_i中的各个地形图在所有受试者中是否高度一致。为此，在一个实施例中，使用二项分布B(M，0.5)对ξ_i中每一个通道上的所有个体水平地形图进行检验(M为受试者数量)，因为该分布反映了虚无假设条件下随机分配每个受试者的该通道是否属于某一个簇而得到的群组水平一致性地形图的取值分布。因为有n个通道，故可使用多重比较校正(例如Bonferroni校正、FDR校正等)对显著性水平进行校正。挑选出包含显著结果通道的那些ξ_i，ξ_i中的这些显著通道集合即构成一个RSN。

应理解的是，在上述实施例中，也可采用其他现有方法来进行个体水平的聚类和组水平的聚类，本发明对具体的聚类方法不作限制。

相应地，本发明还提供一种功能近红外光谱成像静息态脑网络检测装置，用于实现上述方法的一个方面或多个方面。例如该装置包括：第一聚类单元，其用于对于个体水平，通过自动确定簇数进行聚类，将功能近红外光谱测量通道划分为不同的簇；第二聚类单元，其用于对于个体水平的聚类结果进一步划分为不同的聚类，获得组水平的聚类结果；静息态脑网络检测单元，其用于基于组水平的聚类结果中簇的一致性，挑选出构成静息态脑网络的簇集合。该装置中的各单元可采用逻辑电路或处理器实现。

为了验证本发明的效果，在多种条件下进行了一系列模拟实验。具体地，模拟了一个196通道fNIRS测量(14行x14列)，其中添加了一个或多个RSN。每个RSN通道中的时间g(t)是合成的静息态神经活动信号h(t)和高斯白噪声n(t)的线性叠加。非RSN通道中的时间过程是纯白色高斯噪声。合成的静息态神经活动信号是一系列随机相位正弦函数的组合，频率范围为0.01Hz至0.1Hz，步长为0.001Hz，其幅度服从1/f分布。静息态信号在相同的RSN中相同，但在不同RSN之间不同。时间过程的长度为300s，采样率为10Hz。对于每种条件，均会生成100次数据以模拟个体差异。

图2示出了不同信噪比的影响。信噪比定义为

在数据中模拟了两个不同的5×5RSN团块(参见图2(A))，并操纵白高斯噪声的幅度使SNR从0.01到0.05变化，步进为0.0025。在图2中采用

作为拟合优度(GOF)指标来衡量性能，其中n为RSN的数量，C_i是第i个RSN中的通道集合，而

是本发明方法中找到的相应RSN中的显著通道集合。GOF指标的范围是[0，1]，如果检测到的RSN与实际的RSN完全匹配，则GOF将为1；如果该方法未检测到任何真实的RSN，则GOF将为0。由图2(B)可以看出，当SNR在0.02和0.04之间时，拟合优度指数迅速增加，当SNR大于或等于0.04时，本发明的方法可以完美检测所有RSN。

图3是对包含不同数量RSN的模拟数据的验证结果。设置了四个条件，每个条件分别包含一个，两个，三个和四个RSN，其中四行分别对应于不同RSN数量的四个条件。第一列显示空间配置。第二列显示了模拟100次后自动估算的簇数目的直方图。第三列和第二列显示所有组水平簇的群体一致性地形图。RSN群集大小为5×5，SNR为0.05，并在下文保持不变。由图3可以看出，本发明的方法可以完美地找到所有条件下的所有RSN簇(在所有条件下GOF＝1，组一致性检验的显著性水平为0.05，Bonferroni校正)。当实际RSN数小于或等于3时，本发明倾向于将数据划分为四个簇，包括RSN分量和噪声分量。当实际RSN数为4时，本发明将数据分为四个RSN簇和一个噪声簇。

图4是不同RSN体积模拟的结果，图中三行分别表示小(3×3)，中(5×5)和大(7×7)RSN体积的三个条件，分别对应于总测量面积的5％，13％和25％。第一列显示空间配置。第二列显示了模拟100次后自动估算的簇数直方图。第三列和第二列显示所有组水平簇的群体一致性地形图。由图4可以看出，对于所有体积大小，本发明的方法都会检测到所有RSN，小，中和大RSN量的GOF分别为0.67、1和1。

图5是不平衡的RSN体积条件的结果，示出了两个体积比为3：1的不同RSN，以测试不同RSN之间大小不平衡是否会影响性能。对于这种情况，本发明的方法将数据分为两个RSN簇和两个噪声簇。

进一步地，还进行了真实实验，以使用真实的静息fNIRS数据来验证本发明的方法。真实数据包含40位健康的右利手大学生参与者(21.7±2.5岁，男性22位，女性18位)。他们进行了10分钟的静息态扫描和5.6分钟的定位任务扫描。在休息状态下，指示参与者闭眼保持静止并放松思想。在任务态中，他们执行了一个顺序的双边按键任务，其中包含七个任务块，块的长度为20-30s。使用NIRScout连续波fNIRS系统(NIRx Medical Technologies，美国)进行fNIRS测量。以7.8125Hz的采样率测量在两个波长(785nm和830nm)处的近红外光的吸收。在实验中使用了两个4×4探测器组，每个探测器组由八个激光源和八个检测器组成，形成24个测量通道(总共48个通道)。光源与探测器距离为30毫米。将两个探测器组分别放在头上，其中心在国际10-20系统的C3和C4处，以覆盖双侧感觉运动区域。通道的皮层定位是通过使用3维数字化仪和NIRS-SPM软件形成感觉运动模板获得的。含氧血红蛋白(HbO)和脱氧血红蛋白(HbR)信号是根据修正的比尔-兰伯特定律计算获得，其中785nm和830nm的偏路径长度因子(DPF)分别为7.25和6.38。

在本发明的实验研究中，使用了基于谱聚类的RSN检测方法来检测感觉运动网络，首先，将一个0.01-0.08Hz带通滤波器应用于静止状态的HbO数据以提取自发神经活动。然后使用本发明上述方法检测RSN。此外，还使用了传统的基于种子的相关方法，并比较了从这两种方法得出的结果。根据任务激活结果(感觉运动区域，左通道5和9)确定种子通道。使用两个种子通道的平均时间进程作为种子时间进程，以计算每个通道的时间进程与皮尔逊的相关系数。将Fisher-z变换应用于Pearson的相关系数以提高其分布的高斯性。然后，通过使用单样本t检验得出组水平的RSFC t-map。实验中通过使用ROC曲线和ROC曲线下面积(AUC)比较了两种方法在检测以感觉运动模板为特征的感觉运动RSN方面的性能。

图6是真实的fNIRS数据的结果，包括本发明的基于谱聚类的方法(第一行)和现有技术的基于种子的相关方法(第二行)得到的组水平RSN地形图。第三行显示感觉运动模板。本发明的方法获得了四个成分，包括一个感觉运动RSN成分和三个噪声成分。感觉运动RSN成分的组水平一致性图如图6所示(第一行)。通过显著性检验的通道包括左半球的第5、8、9、12、13、16、17、20通道和右半球的第30、33、36、37、40、43通道(p<0.05，Bonferroni校正)。可以看出，本发明的基于谱聚类的方法得出的结果高度特异于预定义的感觉运动模板(图6，第三行)，GOF指数为0.625。与基于种子的相关方法得出的结果明显偏向种子半球(图6，第二行)相比，基于谱聚类的方法得出的结果显示双边感觉运动网络更加对称。

图7显示了ROC曲线分析结果。基于谱聚类的方法显示出比种子相关方法(AUC＝0.73，下方曲线)更高的AUC指数(AUC＝0.87，上方曲线)。

通过上述的仿真实验和真实实验可知，本发明提出的基于谱聚类的功能近红外光谱成像静息态脑网络检测方法不需要任何先验信息和人为干预，即可获得显著效果。一方面，仿真实验结果表明本发明的方法对数据噪声具有鲁棒性，能够有效地从信噪比为0.04的数据中检测出RSN(图2(B))，这一信噪比远低于实际的fNIRS数据的信噪比。此外，本发明的方法在具有不同RSN数目(图3)和大小(图4和图5)的模拟数据上也显示出良好的性能。另一方面，从实际fNIRS数据中检测到的RSN与感觉运动模板具有高度一致性(图6)。图7的ROC曲线表明本发明的方法具有较高的灵敏度和特异性(AUC＝0.87)。这些结果表明本发明所提出的方法对于fNIRS-RSN检测是有效且鲁棒的。

综上，本发明的优势体现在：首先，完全由数据驱动的，对于无法获得用于种子选择的先验信息或难以执行定位任务的情况同样适用，例如，在研究婴幼儿的RSN发育时；易于使用，本发明能够自动计算簇数并确定显著RSN地形图，而无需任何人工干预，也不需要使用者具有丰富的经验；本发明方法得出的结果不会因种子选择或手动成分选择过程而产生偏差，因此可以提供更客观的结果。

从技术上讲，本发明基于谱聚类算法。从理论上讲，本发明比传统的聚类方法(如k-means聚类)更健壮，后者在fNIRS研究中用于检测任务激活和RSFC。此外，现有技术的一些功能磁共振成像研究也采用了谱聚类方法来进行脑的子单元分割和RSN检测。与这些研究相比，本发明的方法具有几个特点。首先，采用一种自动估计算法来确定谱聚类的簇数参数。其次，使用两阶段聚类来计算作为RSN候选者的组水平聚类。第三，使用基于显著性检验的方法从组水平聚类中自动确定RSN簇。

综上所述，本发明提出了一种用于fNIRS-RSN检测的全自动，数据驱动和无模型方法。仿真实验和实际实验结果均验证了其有效性和鲁棒性。本发明的方法可以补充基于种子的相关性和ICA方法，为静息态fNIRS研究提供新工具。

需要说明的是，虽然上文按照特定顺序描述了各个步骤，但是并不意味着必须按照上述特定顺序来执行各个步骤，实际上，这些步骤中的一些可以并发执行，甚至改变顺序，只要能够实现所需要的功能即可。

本发明可以是系统、方法和/或计算机程序产品。计算机程序产品可以包括计算机可读存储介质，其上载有用于使处理器实现本发明的各个方面的计算机可读程序指令。

计算机可读存储介质可以是保持和存储由指令执行设备使用的指令的有形设备。计算机可读存储介质例如可以包括但不限于电存储设备、磁存储设备、光存储设备、电磁存储设备、半导体存储设备或者上述的任意合适的组合。计算机可读存储介质的更具体的例子(非穷举的列表)包括：便携式计算机盘、硬盘、随机存取存储器(RAM)、只读存储器(ROM)、可擦式可编程只读存储器(EPROM或闪存)、静态随机存取存储器(SRAM)、便携式压缩盘只读存储器(CD-ROM)、数字多功能盘(DVD)、记忆棒、软盘、机械编码设备、例如其上存储有指令的打孔卡或凹槽内凸起结构、以及上述的任意合适的组合。

以上已经描述了本发明的各实施例，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施例。在不偏离所说明的各实施例的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。本文中所用术语的选择，旨在最好地解释各实施例的原理、实际应用或对市场中的技术改进，或者使本技术领域的其它普通技术人员能理解本文披露的各实施例。