杏仁核是人类情绪记忆网络中最为核心的脑区，它在情绪的产生和表达、社会判断和面孔识别等认知活动中扮演了重要的角色。它位于内侧颞叶，是一个双侧结构，每个脑半球包含一个杏仁核。杏仁核接受视觉、听觉、其他感觉和痛觉的输入，并将不同的刺激和这些刺激造成的结果联系起来，然后发送信息至其他情绪相关脑区^[1]。杏仁核结构和功能的异常会引起孤独症、抑郁症、帕金森症等精神疾病^[2,3]，它已经成为神经病理学研究的热点，但目前国内外对杏仁核内部结构和功能的细化研究相对较少，如何快速精确自动化对杏仁核脑区进行划分已经显得尤为重要。人脑是中枢神经系统最高级的部分，它是由140亿个神经细胞构成的具有高度组织性的器官，且具有相似功能的细胞团体其活动具有显著的一致性，同样，在杏仁核内部，具有相同功能的细胞核团其功能也会体现出某种程度上的一致性。杏仁核中具有相同功能的细胞核团，其静息态功能也具有高度的相似性。基于这样的理论基础，本文通过静息态功能成像，根据特定脑区内各体素全脑的功能连接模式的相似性，以谱聚类对杏仁核脑区进行分割，验证提出的分割方法的有效性和稳定性，为杏仁核脑区分割方法提供新的策略。

1　材料与方法

1.1　研究对象与采集参数

       采用千人功能连接组脑计划(1000 Functional Connectomes Project)和国际影像数据联盟(International Neuroimaging Data-Sharing Initiative)提供的静息态数据，选取其中内森克莱研究中心(Nathan Kline Institute，NKI)研究组的数据，该数据包含大量的志愿者，年龄跨度大(4~85岁)，静息态数据采集过程中志愿者被要求闭眼休息但不能睡着，本文选取其中年龄在19~29岁的30名志愿者[18名男性和12名女性，年龄(23.2±6.4)岁]的静息态数据来进行脑区分割。志愿者在试验前对其交代了试验完成过程，签署了试验协议，同意进行本试验，并通过了西安交通大学伦理委员会审查。

       磁共振数据使用西门子3.0 T核磁共振设备采取。结构像采集参数如下：扫描矩阵为256×256×176，体素大小为1 mm×1 mm×1 mm。功能像扫描参数如下：TR=2 s，体素大小为3 mm×3 mm×4 mm，扫描矩阵为64×64，持续时间465 s，采集顺序为倒序采集，采集帧数(volume)为200。考虑到功能磁共振成像信号达到稳定和志愿者对环境的适应需要一定的时间，前5个volume在数据发布的时候已经去除，层数设置为34层以覆盖全脑。

1.2　研究方法

1.2.1　杏仁核脑区分割方法

       杏仁核是人脑边缘系统的核心脑区，对于情绪处理具有极其重要的作用。但目前对该脑区内部功能和结构方面的研究很少，杏仁核内部分区的研究对于了解该脑区神经机制以及抑郁症等精神疾病的预防和治疗有着重要意义。基于静息态功能磁共振成像^[4]对杏仁核进行谱聚类分割的算法流程如图1所示。首先采用自动解剖标记(automated anatomical labeling，AAL)模板将皮层和皮层下的灰质划分成90个脑区。左右杏仁核分别为41和42号脑区，将杏仁核感兴趣区域(region of interest，ROI)以及单个志愿者的静息态数据都配准到MNI空间。对于志愿者的静息态磁共振图像数据进行测试，表明一致性和重复性良好。通过计算杏仁核中每个体素的全脑连接模式，进而得到每两个体素之间的全脑连接模式相似度矩阵，以此矩阵进行子区划分。具体做法是首先提取杏仁核每个体素的时间序列信号，对其他88个脑区的静息态信号采用脑区内所有体素信号平均的方式提取。然后通过计算每个杏仁核内部的体素信号与其他88个脑区的信号之间的相关系数，得到杏仁核体素的全脑连接矩阵^[5,6,7,8]M，其中元素(i，j)表示第i个杏仁核体素到第j个脑区的功能连接强度。该矩阵的每一行描述了一个杏仁核体素在静息态下和全脑其他脑区之间的连接模式，值越接近1，相关性越大。最后计算每两行连接模式之间的相关系数就可以度量两个体素在全脑功能连接模式上的相似度，记为相似度矩阵N，其中每个元素表示杏仁核第i个体素和第j个体素在全脑功能连接模式上的相似度。在此相似度矩阵的基础上，本文就可以采用聚类算法对杏仁核的体素进行谱聚类研究^[9]。

点击查看大图

图1  应用谱聚类对杏仁核进行分割的算法流程
Fig. 1  Algorithm flow of applying spectral clustering algorithm for segmentation of the amygdala.

1.2.2　算法实现

       整个分割算法编程基于python 2.7语言实现，它是一种面向对象、解释型计算机程序设计语言，由Guido van Rossum于1989年发明，Python具有丰富和强大的库来满足各种开发需求，本算法科学计算部分主要使用numpy (Numeric Python)模块实现，对磁共振数据的读取和写入则依赖NiBabel模块，NiBabel提供了对常见医学影像数据格式读取的接口。而聚类算法则使用SciKit-Learn模块，SciKit-Learn是一个基于NumPy，SciPy，Matplotlib的开源机器学习工具包。

1.2.3　组分析方法

       本文采用最大概率图谱法(maximum probabilistic map，MPM)对单个志愿者杏仁核脑区的分割结果进行组分析。MPM计算是基于MATLAB中Anatomy 1.7工具箱实现的，其过程是先将所有分割结果配准到MNI空间，然后以将每个像素分配给最可能属于的分割子区域为原则来计算最大概率图谱，如果某些特殊像素属于两个不同子分割区域的概率相等，则这个像素会分给与之相邻的26个像素以求出平均概率更高的那个子区域。

2　结果

       根据Amunts等^[10]在2005年基于细胞构筑的研究，Julich图谱计算并构建了这些杏仁核亚区的概率模板。该图谱将杏仁核分成了3个子区域，分别为基底外侧核(laterobasal，LB)、浅表核(superficial，SF)和中央内侧核(centromedial，CM)分别对应图2A中蓝色、绿色和红色脑区。LB包括副层状核、底内侧核、基底外侧核和外侧核，CM由内侧核和中央核组成，SF包括杏仁核梨状皮质移行区、杏仁前区、杏仁核-海马区以及腹侧和后侧皮质核。图2A中Julich图谱显示出杏仁核子区域的空间位置关系，图2B为本文基于杏仁核内部体素与其他各脑区功能连接进行谱聚类对杏仁核划分的结果，二者进行对比可发现本文的聚类结果和图谱具有高度一致性，表明这一方法的可行

       卡尔•皮尔逊基于弗朗西斯•高尔顿的算法提出了皮尔逊线性相关系数^[11]，通常该相关系数也被称作"皮尔森相关系数"(pearson's linear correlation)。相关系数ρ考察两个变量的相关程度，取值范围在-1~1之间，其中1表示变量完全正相关，0表示无关，-1表示完全负相关。

       式中p(X,Y)：变量X和Y的相关系数；cov(X,Y)：变量X和Y的协方差；：变量X的标准差；：变量Y的标准差；μx：变量X的均值；μy：变量Y的均值。

       本文将大小m×n×o的待测三维图像A和参考图像B看作为随机变量，则两个图像之间的相关系数为：

       式中A表示待评价三维图像，B表示参考图像，Ai,j,k表示待测图像中的像素点，Bi,j,k表示参考图像中的像素点，表示图像A的均值，表示图像B的均值，m表示图像行数，n表示图像的列数，o表示图像高度。图像A和B的相关系数越小(越接近0)，说明图像A和B的相似度越小；相关系数越大(越接近1)，说明图像A和B的相似度越大。

       本文的聚类组分析结果分别与Julich图谱对应的杏仁核子区域的皮尔森相关系数如图3所示，结果越接近1表示二者空间分布越相似。由图3可知右侧CM亚区分割结果与图谱最接近(CC=0.78)，右侧SF亚区与图谱相似度较低(CC=0.65)。

       本文将每个志愿者杏仁核亚区的分割结果与该数据集的组分析结果进行空间相关，以各子区域与组分析结果的空间相似度的变化剧烈程度来反映本文分割算法的稳定性，图4描述了左侧和右侧杏仁核子区单个志愿者与组分析结果的重叠情况，通过计算右侧CM空间相关系数的方差为0.0303(最大)，表明该亚区分割的结果稳定性比较差，而右侧LB脑区空间相关系数的方差为0.0017(最小)，说明这个脑区聚类算法的分割结果可重复性(最高)，稳定性也最好(其他脑区的方差分别为var (CML)=0.0041，var (LBL)=0.0185，var (SFR)=0.0098，var (SFL)=0.0194)。

       在讨论了分割算法的精确性和稳定性的基础上，叠加不同强度的噪声到杏仁核体素的时间序列中再进行脑区划分来考察本文算法的健壮性，过程如图5所示。首先将杏仁核左右脑区各体素的时间序列提取出来，将信噪比分别为10~90 db的高斯白噪声叠加到各时间序列中，然后基于谱聚类进行杏仁核子区域分割，整个实验仿真过程在python 2.7平台实现。

       叠加噪声后杏仁核6个亚区的聚类分割如图6所示，图6反映了不同噪声下基于ROI体素时间序列与AAL模板其余88个脑区平均序列计算的功能连接矩阵进行分割的结果与没有加入噪声的情况下各子区空间相似性的变化情况。由图6可知不同信噪比的图像数据均取得了比较满意的分割结果，在信噪比大于70 db的情况下，各子区都保持了近乎与没有加入噪声的结果100%的相似性，在10 db的低信噪比下各区皮尔森相关系数也能保持在0.7左右。

点击查看大图

图2  杏仁核生理解剖图谱与本文谱聚类分割结果对比。蓝色、绿色和红色脑区分别为LB、SF和CM子脑区
Fig. 2  Comparison of the physiological anatomy of the amygdala and the results of this paper. Blue, green and read brain regions denote LB, SF and CM subregions, respectively.

点击查看大图

点击查看大图

点击查看大图

图3  谱聚类分割杏仁核结果与Julich图谱相似度
Fig. 3  Similarity between spectral clustering segmentation of amygdala and Julich patterns.

点击查看大图

图4  单个志愿者与左右侧杏仁核组分析结果空间相关系数分布图
Fig. 4  Spatial correlation coefficient distribution of single and left and right amygdala group.

点击查看大图

图5  健壮性分析：叠加噪声到杏仁核体素的时间序列中的谱聚类脑区划分。蓝色、绿色和红色脑区分别为LB、SF和CM子区
Fig. 5  Robustness analysis: spectral clustering time series of brain regions when adding noise into amygdala. Blue, green and read brain regions denote LB, SF and CM subregions, respectively.

点击查看大图

图6  加入噪声对基于功能连接和基于血氧水平依赖信号分割结果的影响
Fig. 6  Effects of adding noise on the functional connectivity and blood oxygen level dependent-based signal segmentation results.

3　总结

       本文提出了一种基于静息态功能成像的杏仁核谱聚类分割方法，即一种基于谱聚类算法，根据杏仁核内部体素功能的相似性对该脑区进行自动化高效率划分的方法，先对静息态磁共振数据预处理、再对杏仁核脑区提取、然后对杏仁核内部体素全脑功能连接计算以及最后对功能连接矩阵进行谱聚类分割。本文提出的自动分割算法与杏仁核临床解剖结果取得很大程度的一致性，并且在稳定性和抗噪声干扰方面也取得较满意的结果，相对于传统手动分割方法，具有更加简单、方便、高效、可重复性高等多方面的优点。