一种基于pca和gca的磁共振脑成像因果连接强度的检测方法

文档序号:868033阅读:599来源:国知局
专利名称:一种基于pca和gca的磁共振脑成像因果连接强度的检测方法
技术领域
本发明属于图像处理领域,具体涉及一种磁共振脑成像因果连接强度的检测方法。尤其涉及使用主成分分析PCA和格兰杰因果分析GCA进行磁共振脑成像因果连接强度的检测。
背景技术
功能磁共振成像(functionalMagnetic Resonance Imaging, fMRI)以其高时空分辨率,非侵入式等特点在神经疾病诊断治疗方面得到了广泛应用。fMRI —般指基于血氧水平依赖(blood oxygen level-dependent, BOLD)的磁共振成像,它通过测量由神经活动引起的脑血流和脑血氧等成分变化而造成的磁共振信号变化来反应脑活动。人脑是一个复杂的系统,在受到刺激条件或经历病变时各个脑区间存在的相互作用也在发生变化,反映在脑区间的功能连接强度和方向变化上。脑功能数据连接强度分析一般是根据测得的fMRI 数据寻找其活动能显著区分不同试验条件(如,刺激条件和基线条件)的脑区,再检测脑区之间的相互连接强度。合理的脑功能数据连接强度分析方法需同时考虑到人脑功能一般组织原则以及fMRI数据本身一些特性。本发明为脑区间功能连接的分析和可视化提供一种有效途径,为神经疾病的早期诊断和临床治疗提供帮助。脑功能一般遵循两个基本组织原则功能集成化和功能特异化。在大空间刻度 (规模)上,一个复杂的脑功能可能会由许多功能特异的脑区通过相互作用(集成)来完成;反过来,一个特异性脑区也会对许多不同的刺激任务进行表示或加工,通过精细空间刻度上不同的分布式脑活动来对外部不同刺激进行表示。另一方面,随着MRI技术的进步, fMRI图像空间分辨率一直在提高。传统fMRI数据在体素各维宽度约为4毫米时可获得较高的信噪比。当前,各维宽度为2毫米的体素已经可以在标准3T临床磁共振(MRI)机器上可靠地获得。随着超高场(》7T)磁共振技术的使用,fMRI图像的空间分辨率正在向亚毫米的刻度上逼近。fMRI图像空间分辨率的提高,为我们在更精细的空间刻度上研究脑功能提供了可能。目前,传统的检测脑功能因果连接强度方法大多采用动态因果模型(dynamic causal modeling, DCM),它需要预先选定相互作用的区域,并假设这些区域的任意两个之间存在影响,这种预先假设的模型在验证一些有关大脑系统之间的假设时是起作用的,但如果出现了对模型的错误指定,如漏掉了一个间接作用或新起作用的区域,就会导致错误的结论。格兰杰因果分析(Granger causality analysis, GCA)检测方法可以用于研究区域间的因果关系,而不需要先验知识,克服了 DCM的局限性,且其突出的优越性在于无需事先假设两点之间存在解剖结构的连接性。同时GCA可反映因果方向性,这为脑神经网络回路的研究提供了方便,这是普通脑功能连接强度检测方法所无法做到的。另外,以往对于同一脑激活区所包含多个体素的时间序列,采取简单平均处理得到一个平均的时间序列,这样会丢掉大量有用脑活动信号信息。本发明采用主成分分析,对多体素的时间序列进行空间上的降维,尽量消除各时间点上所有体素脑活动信号值所构成的向量之间的相关性,能够利用脑活动信号信息的主要成分,有效地提高磁共振脑功能图像因果连接强度的检测性能,更完整更准确地检测脑激活区之间的因果连接。

发明内容
为了克服已有技术的不足,本发明的目的在于设计一种提高磁共振脑成像因果连接强度的检测结果准确性及完整性的检测方法。为实现上述目的,本发明提出一种基于主成分分析PCA和格兰杰因果分析GCA的磁共振脑成像因果连接强度的检测方法,包括以下步骤步骤M 对经过数据预处理的脑功能图像进行激活区域多体素时间序列的提取, 得到激活区域多个体素的时间序列矩阵;步骤Sb 使用PCA对每个激活区多个体素的时间序列矩阵进行空间降维得到保留大部分信息量的主要成分,再平均主要成分的激活值得到一个时间序列;步骤k 建立所有激活区域的时间序列间的多变量自回归模型;步骤Sd 依据所述多变量自回归模型,计算各个时间序列间的偏相关系数;步骤% 依据偏相关系数,通过直接传递函数DTF方法计算dDTF值,得到激活脑区间的因果连接强度和方向;步骤Sf 使用替代数据法统计检验连接强度的显著性,并利用有向网络图显示显著的连接。本发明针对磁共振脑成像因果连接强度的检测,通过激活区多体素时间序列矩阵提取、PCA降维、使用基于多变量自回归模型及传递函数方法的格兰杰因果分析计算激活区之间的因果连接强度,再用替代数据的方法检验连接强度的统计显著性,提高了因果连接强度检测的准确性和鲁棒性。真实数据试验结果证明该算法比传统采用基于简单平均激活区所有体素激活值的GCA方法能更完整更准确地检测到脑激活区之间的因果连接,为脑功能数据分析方法提供了一种新途径。


图1是本发明算法计算流程示意图。图2是根据本发明算法的实施例获得的因果连接强度图。图3是使用基于简单平均体素方法获得的因果连接强度图。
具体实施例方式为使本发明的目的、技术方案和优点更加清楚明白,以下结合具体实施例,并参照附图,对本发明进一步详细说明。参照图1,本发明所涉及的一种磁共振脑成像因果连接强度的检测方法,尤其涉及使用主成分分析和格兰杰因果分析进行磁共振脑成像因果连接强度的检测的具体实施步骤如下步骤Μ,对经过数据预处理的脑功能图像进行激活区域多体素时间序列的提取, 得到激活区域多个体素的时间序列矩阵;1.数据预处理
5
由于磁共振扫描过程中各种各样的噪声的影响,个体自身存在尺度和位置上的差异,非常有必要在分析数据之前对数据做一定的预处理。在整个的实验的数据获取中,主要的噪声信息来源有(1)物理头动;(2)图像内层间扫描时间差别;(3)外在磁场的不均勻性等。脑功能图像预处理是在保留脑功能图像细节的同时,使用脑功能图像与标准模板进行仿射配准变换方式的预处理,并提高脑功能图像的信噪比。常见的预处理步骤有切片扫描时间对齐,图像序列对齐,联合配准,标准化(或称均一化),空间平滑滤波和时间平滑滤波等。2.独立成分分离一般来说,分离独立成分的方法有基于广义线性模型(General Linear Models, GLM)的方法、模式分类方法和独立成分分析方法(Incbpendent component analysis, ICA)。其中,基于ICA的分离独立成分模型如下令ν为脑功能观测量(维度),其由多个信源S经混合矩阵A组合而成(ν = AS)。现在的任务是在S与A均未知的条件下,求取解混矩阵W,使得ν通过W后所得输出y (y = W · ν)是S的最优逼近。本发明采用独立成分分析方法来分离独立成分,分离步骤如下假设分离出的各分量互相独立,由中心极限定理可知,一定条件下相互独立的随机变量之和趋向于高斯分布,所以两个相互独立的随机变量之和比任何一个参与求和的随机变量更加接近于高斯分布,这样最后求出的非高斯性最大的组合随机变量,即距离高斯分布最远的随机变量y就可以认为是源信号ν的一个独立成份,这种使距离高斯分布最远的随机变量之间相互独立的方法称为非高斯最大化方法。另外,还可用很多其它方法来定义这种接近程度,如基于信息论的最小化互信息法和最大似然估计法也是对ICA进行估计的方法。并由此发展出很多ICA算法,如最大(小)化一些相关的代价函数,基于随机梯度的自适应算法,还有现在常被采用的快速ICA算法以及hfomax算法等。在处理fMRI数据时,我们一般采用^fomax算法寻找解混矩阵,这样也就能有效地分离多个非高斯性最大的源信号。hfomax算法是一种基于信息最大化传输的单层前反馈神经网络算法,其基本思想是利用了似然函数极大法估计中的一个准则,即估计随机变量ν独立性的最大化等价于估计ν的某一非线性函数独立性的最大化。它的主要特点是引入了非线性函数f( ·),使得新向量y = f (U)的熵极大化等价于u的各分量独立性的最大化,其中u为输入ν各分量的线性组合。在该算法中选用Sigmoid函数作为f ( ·),因为Sigmoid函数具有超高斯信号 Pdf(概率密度函数)的性质,所以能有效地分离多个超高斯分布的源信号。Sigmoid 函数的形式为:f(x) = (1+e-T1。所谓信息最大化传输(熵极大化)就是使输入V与输出y的互信息最大,互信息定义为I(y, v) = H(y)-H(y/v),其中H(y)是输出信号y的边缘熵,H(y/v)为条件熵。为求互信息I (y,ν)的最大值,将等式两边对W求导 ——/(ν,ν)=——H(y),
dW ^ J dw k^j 其中,条件熵H(y/v)不依赖于解混矩阵W。
输出熵H(y)定义为
权利要求
1.一种基于主成分分析PCA和格兰杰因果分析GCA的磁共振脑成像因果连接强度的检测方法,其特征在于,包括以下步骤步骤Μ 对经过数据预处理的脑功能图像进行激活区域多体素时间序列的提取,得到激活区域多个体素的时间序列矩阵;步骤Sb:使用PCA对每个激活区多个体素的时间序列矩阵进行空间降维得到保留大部分信息量的主要成分,再平均主要成分的激活值得到一个时间序列;步骤& 建立所有激活区域的时间序列间的多变量自回归模型;步骤Sd 依据所述多变量自回归模型,计算各个时间序列间的偏相关系数;步骤% 依据偏相关系数,通过直接传递函数DTF方法计算dDTF值,得到激活脑区间的因果连接强度和方向;步骤Sf 使用替代数据法统计检验连接强度的显著性,并利用有向网络图显示显著的连接。
2.如权利要求1所述的检测方法,其特征在于,所述脑功能图像的数据预处理是在保留脑功能图像细节的同时,使用脑功能图像与标准模板进行仿射配准变换方式的预处理。
3.如权利要求1所述的检测方法,其特征在于,提取激活区域多体素时间序列的步骤为(1)利用独立成分分析ICA方法得到不同激活模式的独立成分;(2)获得所述独立成分中激活值最大点体素的坐标,以此坐标为中心作一个半径为3 倍体素棱长的球体;(3)依据预处理后的脑功能图像数据,提取包含于此球体内部的各个体素在不同时间点上激活值的时间序列。
4.如权利要求3所述的检测方法,其特征在于,所述获得独立成分中激活值最大点体素的坐标的步骤具体包括(1)在总体独立成分数据上进行单样本t检验,获得平均成分;(2)对每个平均成分中每个体素的时间序列求得方差值,作为此体素的激活值;(3)从感兴趣的成分中确定出激活值最大点体素的坐标。
5.如权利要求1所述的检测方法,其特征在于,使用PCA对每个激活区多个体素的时间序列矩阵进行空间降维得到保留大部分信息量的主要成分,再平均主成分得到一个时间序列的步骤为(1)使用总和标准化法标准化具有D个体素的时间序列矩阵Y',得到新矩阵Y";(2)计算Y"的相关系数矩阵P= (Pij),矩阵P的维数为DXD, Cov(Y:,Y;)V^cBvW'其中,P ij为相关系数矩阵元素,CovO表示两时间序列的协方差,Y" i,Y" j分别表示矩阵Y"的第i,j行元素,DO表示每行矩阵元素的方差;(3)计算相关系数矩阵P的特征值,并将特征值按从大到小的顺序排序 A1^ λ2彡…彡λ D,其相应的特征向量为Vl,v2,...,Vd;(4)计算累计贡献率达到85%的前d个特征值,其中累计贡献率按下式计算
6.如权利要求1所述的检测方法,其特征在于,所有激活区的时间序列间的多变量自回归模型为
7.如权利要求6所述的检测方法,其特征在于,步骤Sd具体包括以下步骤(1)将步骤&中得到的多变量自回归模型从时域转换到频域,得到频率转移矩阵 H(f)X(f) = A"1 (f) E (f) = H(f)E(f),其中,f表示频率,X (f)为X(t)在频域上的表示,A (f)为A (η)在频域上的表示,E (f) 为E(t)在频域上的表示,H(f) =A-^f)为频域转移矩阵;(2)计算各个时间序列间的偏相关系数eu(f)。
8.如权利要求7所述的检测方法,其特征在于,所述偏相关系数的计算包括(1)计算交叉频谱S(f) S(f) = H(f) VHh (f),其中,V表示零均值白噪声矩阵E (t)的方差,上标H表示共轭转置;(2)计算偏相关系数θ^.( )
9.如权利要求1所述的检测方法,其特征在于,步骤%中通过直接传递函数DTF方法计算dDTF值具体为
10.如权利要求1所述的检测方法,其特征在于,步骤Sf中使用替代数据的方法进行统计检验连接强度的显著性,并利用有向网络图显示显著的连接。
全文摘要
本发明公开了一种基于主成分分析PCA和格兰杰因果分析GCA的磁共振脑成像因果连接强度的检测方法,该方法包括以下步骤对经过数据预处理的脑功能图像进行激活区域多体素时间序列的提取,得到激活区域多个体素的时间序列矩阵;使用PCA对每个激活区多个时间序列矩阵空间降维得到主要成分,平均主要成分的激活值得到一时间序列;建立所有激活区域的时间序列间的多变量自回归模型;计算各个时间序列间的偏相关系数;通过直接传递函数DTF方法计算dDTF值,得到激活脑区间的因果连接强度和方向;使用替代数据法统计检验连接强度的显著性,在有向网络图上将结果显示。在真实数据集上的实验说明,本发明所述方法是一种有效的磁共振脑成像因果连接强度的检测方法。
文档编号A61B5/055GK102366323SQ20111029957
公开日2012年3月7日 申请日期2011年9月30日 优先权日2011年9月30日
发明者田捷, 白丽君, 钟崇光 申请人:中国科学院自动化研究所
网友询问留言 已有0条留言
  • 还没有人留言评论。精彩留言会获得点赞!
1