专利名称：一种基于移动平均时间序列模型的脑功能激活区检测方法
技术领域：
本发明属于生物信息、图像识别技术领域，涉及一种基于时间序列状态分析的核 M^Ml^it (functional Magnetic Resonance Imaging, fMRI) g.M^biiii^。
背景技术：
功能磁共振成像(fMRI)技术是在传统磁共振成像(Magnetic Resonance Imaging,MRI)技术的基础上发展而来。fMRI技术基于血氧水平依赖(Blood Oxygenation Level-Dependent, BOLD)，通过测量由脑血流或脑血氧变化而引起的磁共振信号变化来反映脑活动。由于其具有无创性、可重复性等优点，并可提供较高的图像空间分辨率，该技术已被广泛应用于神经科学、临床医学等领域的研究。作为一种非侵入性检测手段，fMRI可为一些重要脑功能区(如运动、语言等)定位。脑功能定位研究是神经科学研究发展极为迅速的领域之一，为人类的认知与思维活动及疾病的康复治疗与预测带来了新的前景。目前， 脑功能定位研究已广泛应用于神经外科、认知功能相关研究、康复治疗与预测、社会心理等方面。随着基础研究与临床方面的应用深入，合理有效的分析fMRI数据并提取与外部刺激相关的脑激活区对于脑功能活动区的精确定位愈发重要。现有fMRI脑功能激活区检测方法大多基于广义线性模型(General Linear Model, GLM)方法。该模型本质上是一种模型驱动(Data Driven)的一元回归分析方法， 一般先指定一个先验模型(即一个包含实验设计、人为增加的心跳、呼吸、头动等因素的模型)，将该模型分别与大脑中的每一个体素的时间序列进行匹配，确定每一个体素的“激活” 状态。因此，该方法的应用前提是fMRI时间序列中“激活”体素的信号强度时间序列变化需符合实验设计中的指定的先验模型。然而，越来越多的影像学研究及临床应用表明，神经激活信号的时间变化模式并无准确先验。因此，基于GLM模型的先验假设过于单一，并不能够很好地表示fMRI数据。脑活动的分布模式、血液动力学模型等不能由GLM模型恰如其分地模拟。
发明内容
为了避免现有方法的不足，本发明设计一种基于移动平均时间序列模型的脑功能激活区检测方法，以克服现有方法依赖先验模型的缺点，达到有效提取脑激活区目的。为实现上述目的，本发明提出的一种基于移动平均时间序列模型的磁共振成像脑激活区检测方法，其特征在于，该方法包括以下步骤步骤1，对采集到的原始功能磁共振成像(fMRI)图像进行预处理；步骤2，对于每个被测试者逐体素计算每个被测试者的个体体素信号强度时间序列模型；步骤3，对于组群，计算每个个体在组群中的权重值，并根据全脑个体体素信号强度时间序列模型，逐体素构建组群体素信号强度时间序列模型；步骤4，根据全脑组群体素信号强度时间序列模型，应用区域增长方法实现脑激活区的提取。本发明针对现有脑影像fMRI数据分析方法过于依赖先验模型这一问题，采用移动平均方法构建大脑中每个体素的时间序列状态变化模型，该时间序列动态变化模型不依赖与先验知识；通过区域增长方法度量各体素时间序列动态变化的活动相似性，可自适应地提取脑功能激活区。将该方法应用于fMRI影像数据集的结果亦表明，该方法是一种有效的脑功能激活区提取方法。
图1是本发明方法流程示意图。图2是本发明中个体体素信号强度时间序列动态变化模型。图3是本发明中组群体素信号强度时间序列动态变化模型。图4是本发明中提取的脑功能区激活图。
具体实施例方式为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明进一步详细说明。参照图1，本发明提出的一种基于移动平均时间序列模型的脑功能激活区检测方法，具体包括以下步骤步骤1、对采集到的原始fMRI图像进行预处理。首先对输入的三维fMRI图像进行时间点校正和头动校正；再将校正后的图像进行标准化处理采用SPM(Statistical Parametric Mapping, SPM，统计参数映射)软件自带的EPI (Echo Planar Image，平面回波成像)模板，将校正后的图像归一化到蒙特利尔神经科学研究所(Montreal Neurological Institute, MNI)标准脑模板上，并把三维fMRI 图像中的每个体素重新采样至3X3X3mm3 ；然后对标准化后的图像使用高斯核函数进行平滑，其中，高斯核函数的全宽半高值(FWHM)为6X6X6mm3;最后对平滑后的图像去除线性漂移，将扫描基线的漂移进行去除。步骤2、对于每个被测试者(individual subject)，逐体素计算每个被测试者的个体(individual subject)体素信号强度时间序列模型。步骤2进一步包括以下几个步骤步骤2. 1，构建个体体素信号强度时间序列预测统计量；对于每个体素，假设第i个个体体素的fMRI神经响应的时间序列分布状态X满足以下两种数据分布形式基线状态& Ν( θ μ Σ D，以及激活状态& Ν( θ ” Σ D，两种状态分布的均值不同，分别对应于无刺激状态和刺激状态。其中，θ ^，Q1为均值，Ei为协方差。那么，实测的信号强度时间序列统计模型为X1t =S1t= \,...,m, t = l，...，n其中，t为时间序列取值的时间点，x丨为第i个个体在时刻t测得的信号强度，勻为第i个个体在时刻t所对应的噪音项，f为第i个个体对应于整个时间序列的噪音向量，
t服从均值为0，协方差为Σ i的正态分布，m为所有被测试者人数，η为时间点个数，4表示第i个个体在t时刻的信号强度期望值，^满足,；={ θ°, = 1,…,T ，其中，τ为从基线状态到激活状态的状态改变所对应的时间点，在τ时刻，fMRI 时间序列的状态发生改变，如由外界刺激条件引起大脑皮层神经活动的增强(或减弱)，进而使时间序列中的数据从均值为θ ^转变到均值为Q1的分布状态。利用移动平均法构建个体体素信号强度时间序列预测统计量ζ丨，ζ丨满足Z1t = Ax1t + (1 — λ)ζ\_λ， /\<t<n,其中，λ为权重系数，0彡λ彡1，λ决定了该模型检测fMRI时间序列信号突然发生变化的能力，即时效性。λ越大，时效性越强，反之越弱。这里选择λ =0.2。步骤2. 2，基于个体体素信号强度时间序列预测统计量，统计推断个体体素信号强度时间序列变化状态，完成个体体素信号强度时间序列模型的建立。通过对数据分布形式加以判断，可以判断特定时刻时间序列信号的状态是否发生变化。为此做出如下假设零假设(Htl)为只存在一种均值状态，均值μ为常数，即信号状态没有发生改变；备择假设(H1)为存在两种状态(基线状态与激活状态)，均值μ对应不同值，即H0 μ = θ ，t = 1，. · ·，η ；H1 μ = θ ，t = 1，· · ·，τ ；μ = θ j, t = τ +1, . . . , η ；为了检验个体体素信号强度的状态是否发生变化，定义检验统计量T为
Γ π Τι Z't —θ0Λ = ι .，
^jVar (ζ)其中，Far^)为统计量ζ丨在t时刻的方差。检验统计量T服从自由度为d的t分布t ( Σ，df)，控制阈限(Control Limits, 即置信度区间)定义为代±广 ^Kar(X)，其中，t*为假阳性率(一般取值为95% )。一旦统计量ζ丨超过该阈限值，则表明t 时刻该体素处于“激活”状态。由于检验统计量T服从自由度为d的t分布t( Σ, df)，因此fMRI时间序列中每个时间点t都对应着一个概率P值。为遍历整个时间序列，需对P值进行校正。 令？=识,7；,...,7；丨为对应于所有时间点的检验统计量向量，f服从多变量t分布，即 ？ 《Σ,O。优选地，可采用蒙特-卡罗(Monte Carlo integration)方法在时间维度上对 P值进行FWER (Fami ly-wise Error Rrate)校正，以估计t*值。具体来说，可由t ( Σ , df) 分布随机生成2000个长度为η的特征向量，用Tk表示其中一个特征向量。对于Tk，可求得 Tk各分量的绝对值最大，表示为『腿=Jaxfrj。因此，对应于随机生成的2000个特征向量，可求得2000个Tmax值。将Tmax按照由
6大到小的顺序排序=Tmaxil > Tfflax,2Tmax,2_，由于蒙特-卡罗方法的校正过程本质上是
多重检验的过程，一般通过控制假阳性率(FDRJalse Discovery Rate)来决定ρ值。实践过程中FDR的值一般选取0. 05。因此t*取为Tmax,95，t*即为校正后的ρ值。对于每个个体i的任一体素，用f表示第i个个体对应于整个时间序列的信号强度预测统计量，^模拟了信号强度时间序列动态变化，通过与阈限值
代土广·^/沿r(z丨)比较，可以判断某时刻该体素是否处于“激活”状态。至此，通过计算个体体素信号强度时间序列预测统计量，并在此基础上估计信号强度时间序列变化状态，完成了个体体素信号强度时间序列模型的建立。上述描述是针对单一体素举例，将上述步骤应用于全脑体素，可以得到全脑个体体素信号强度时间序列模型。步骤3、对于组群(group subjects)，计算每个个体在组群中的权重值，并根据全脑个体体素信号强度时间序列模型，逐体素构建组群(group subjects)体素信号强度时间序列模型。步骤3进一步包括以下几个步骤步骤3. 1，构建组群体素信号强度时间序列预测统计量；对于组群，考虑到不同被测试者之间的个体差异，实测的信号时间序列统计模型可写为x't =s't+ ε\, i = 1，···，/ηs't =sfop+η η t = l,---,n其中，χ丨为第i个个体在时刻t测得的信号强度 ’ε;为第i个个体在时刻t所对应的噪音项，^ =( ，...，。,S1-N(OX1)表示第i个个体在时刻t的信号强度期望值；
nt为时刻t组群内不同被测试者个体差异引起的噪音项^~#(0，25；)，这里[B为不同个体间的噪音协方差矩阵，m为个体个数，η为时间点个数，&胃满足Sfop=I，将该模型写为矩阵形式，可得 ―1 ―1 ―1χ = S +εm由上面2个公式，该模型又可写为一1 一一一1 一 -^1 X _ S pop + Tj~\~ ￡ — S pop H- ζ j其中，用^表示第土个个体的总噪音项，^ ^^，么+乙^，令^^^+乙^^为第 i个个体的体素信号强度时间序列预测统计量。则第i个个体信号强度时间序列预测统计量-的噪音协方差矩阵f可写为V； =AV1At =Σ*+Σ；
其中，Λ为
权利要求
1.一种基于移动平均时间序列模型的磁共振成像脑激活区检测方法，其特征在于，该方法包括以下步骤步骤1，对采集到的原始功能磁共振成像(fMRI)图像进行预处理；步骤2，对于每个被测试者逐体素计算每个被测试者的个体体素信号强度时间序列模型；步骤3，对于组群，计算每个个体在组群中的权重值，并根据全脑个体体素信号强度时间序列模型，逐体素构建组群体素信号强度时间序列模型；步骤4，根据全脑组群体素信号强度时间序列模型，应用区域增长方法实现脑激活区的提取。
2.如权利要求1所述的方法，其特征在于，所述步骤1中对采集到的原始fMRI图像进行预处理进一步包括对fMRI图像进行时间点校正、头动校正、标准化处理、重新采样、平滑以及去除线性漂移。
3.如权利要求2所述的方法，其特征在于，对fMRI图像进行标准化处理时采用统计参数映射软件自带的平面回波成像(EPI)模板对fMRI图像进行标准化处理。
4.如权利要求1所述的方法，其特征在于，所述步骤2进一步包括以下步骤步骤2. 1，构建个体体素信号强度时间序列预测统计量；步骤2. 2，基于个体体素信号强度时间序列预测统计量，统计推断个体体素信号强度时间序列变化状态，完成个体体素信号强度时间序列模型的建立。
5.如权利要求4所述的方法，其特征在于，所述步骤2.1进一步使用移动平均算法来构建所述个体体素信号强度时间序列预测统计量。
6.如权利要求4所述的方法，其特征在于，所述步骤2.2进一步采用蒙特-卡罗方法在时间维度上对概率P进行FWER校正来统计推断个体体素信号强度时间序列变化状态，其中，概率P为fMRI时间序列中时间点t所对应的概率值。
7.如权利要求1所述的方法，其特征在于，所述步骤3进一步包括以下步骤步骤3. 1，构建组群体素信号强度时间序列预测统计量；步骤3. 2，计算每个个体在组群中的权重值D ；步骤3. 3，基于组群体素信号强度时间序列预测统计量，统计推断组群体素信号强度时间序列变化状态；步骤3. 4，基于全脑个体体素信号强度时间序列模型、每个个体在组群中的权重值D以及组群体素信号强度时间序列变化状态，逐体素构建组群体素信号强度时间序列模型。
8.如权利要求7所述的方法，其特征在于，所述步骤3.2进一步使用期望最大化算法来计算每个个体在组群中的权重值D。
9.如权利要求1所述的方法，其特征在于，所述步骤4进一步包括以下步骤步骤4. 1，度量不同体素间的活动相似性；步骤4. 2，使用体素Vtl作为种子点，使用区域增长算法选择同种子点活动最相似的邻域体素加入局部一致脑区；迭代进行区域增长算法，直到局部一致脑区达到指定局部区域大小后停止；步骤4. 3，将上述步骤4. 1和步骤4. 2应用到全脑体素，根据脑激活区内体素信号强度时间序列动态变化基本一致这一原理，实现脑激活区的提取。
10.如权利要求9所述的方法，其特征在于，所述步骤4. 1进一步使用组群体素信号强度时间序列模型的皮尔逊相关系数作为区域一致性准则，度量不同体素间的活动相似性。
全文摘要
本发明公开了一种基于移动平均时间序列模型的磁共振成像脑激活区检测方法，该方法包括以下步骤使用移动平均算法，对于每个个体，逐体素计算个体(individual subject)体素信号强度时间序列模型；在此基础上使用EM-algorithm计算每个个体在组群中的权重值，逐体素构建组群(group subjects)体素信号强度时间序列动态变化模型；最后通过区域增长方法，度量组群体素间信号强度时间序列模型的活动相似性，确定脑功能激活区。在真实fMRI数据集上的实验说明，本发明所述方法是一种有效的脑功能激活区提取方法。
文档编号A61B5/055GK102508184SQ20111032939
公开日2012年6月20日 申请日期2011年10月26日 优先权日2011年10月26日
发明者田捷, 白丽君, 薛婷 申请人:中国科学院自动化研究所