一种功能核磁共振时间序列匹配方法组成比例

本发明专利技术公开了一种功能核磁共振时间序列匹配方法，包括以下步骤：计算系数矩阵β的初始值；计算观测信号时间序列Y与设计时间序列X间的范式距离Fdist；随机获取新的记录点；计算新的范式距离Fdist；输出β，完成Y与X的β匹配。本发明专利技术在fMRI时间序列匹配中将时域信号经快速离散傅立叶模变换后得到频域序列，从而完全消除相位信息，以达到消除要匹配的fMRI时间序列间相位差的目的。本发明专利技术与当前相位校正的方法相比，更为简单，消耗的计算量少。本发明专利技术对离散傅里叶变换进行加权约束，以降低高频部分的影响，优先考虑了低频有效信号，更为明确地确定“最感兴趣频率”的位置，从而降低激活检验结果出现“伪阳性”体素的可能性。

【技术实现步骤摘要】

本专利技术涉及一种数据处理技术，特别是一种功能核磁共振时间序列同频不同相的匹配方法。
技术介绍
功能磁共振(fMRI)技术，是无损地研究人脑行为活动的重要影像技术，具有广泛而重要的应用前景。fMRI信号是一种难以直观解读的观测型信号，因此针对fMRI信号的合理分析具有重要意义，这既是热点但也是难点问题。fMRI技术的通过大脑的血氧水平依赖效应呈现大脑的活动情况，通过该效应，可以实时地检测大脑活动的脑血流变化情况。fMRI数据有以下几个特点：1)延迟特性。一个脑功能刺激，fMRI信号往往并不能做出及时的反馈，而是在几秒后达到峰值。2)混淆特性。在以往的信号分析研究中，fMRI信号通常被看作一种由不同刺激经过线性叠加而形成的观测信号。3)采样点多。fMRI的扫描精度比较高，若以2秒每帧的速度进行扫描，可以获得64*64*32分辨率的大脑图像。因此对维数比较敏感的算法，在应用到fMRI数据中往往会遇到困难。3)信噪比低。大脑是一个复杂的整体，在进行某些特定认知任务的时候，大脑还肩负着人体的调节工作，比如呼吸、心跳等一系列维持生命所必须的生理活动。而在进行有意识的活动的同时，大脑也会出现许多伴随着诸多因长期的生活训练而出现的无意识的活动。fMRI实验数据按照实验刺激的不同主要分为静息态数据和功能实验数据。静息态数据无需被扫描者做出任何活动，主要用于观察人处于无意识状态下的大脑活动。而功能数据则是在进行脑科学实验时，实验设计者给予受试者外界刺激，诱导受试者有意识地进行思维活动所产生的fMRI数据，主要用于观察当人处于有意识的行为时，大脑的活动状态。目前用于分析fMRI行为实验数据的常用软件主要有统计参数映射模型(SPM)、fMRI组独立成分工具箱(GIFT)、静息态数据处理方法(Dparsf)等，这些工具已经用于进行人脑行为认知的研究当中，最为广泛的是SPM(Statistical Parametric Maps)，它可以提取特定刺激因素引起大脑活动的范围，即脑激活区。SPM主要利用刺激时间序列这一先验条件，提出一种线性叠加模型，即广义线性模型(GLM)，来分析不同脑区的fMRI时间序列信号与刺激时间序列的线性叠加关系，从而判定该区域与特定刺激的相关性。SPM框架包括多个组成部分，包括时间校正、头动校正、标准化预处理等数据预处理方法，并集成GLM、pair&no pair T检验、F检验和贝叶斯等多种数据分析方法。Rest、Dparsf、GIFT、MICA等fMRI数据分析工具包直接集成或加载SPM的部分功能，或要求数据由SPM进行预处理后再进行数据分析。SPM虽然被广泛应用，但也有部分fMRI信号并没有得到很好地转换为脑区激活信息，因为其依赖的广义线性模型(SPM-GLM)对刺激时间序列的要求比较苛刻，这大大限制了fMRI实验设计的灵活性，而且该方法将fMRI时间序列视为一个向量，通过评价拟合参考序列与待测时间序列信号的范数距离作为目标函数，即待测fMRI信号与模型参考时间序列的匹配程度，来评价优化的目的。但是作为时间序列，相邻的时间点间联系是很紧密的，而fMRI信号作为时间序列，由刺激引发的信号发生延迟的个体差异较大，即使同一个体对不同的刺激反应延迟也存在差异，因此当脑神经元活动的时间与刺激呈现的时间存在较大不确定性的差异时，即大脑活动的时间序列本身难以准确描述的时候，传统的SPM-GLM方法很难得到理想的结果。
技术实现思路
为解决现有技术存在的上述问题，本专利技术要提出一种简便的、可以克服延迟的功能核磁共振时间序列匹配方法。为了实现上述目的，本专利技术的技术方案如下：一种功能核磁共振时间序列匹配方法，包括以下步骤：A、计算系数矩阵β的初始值将广义线性模型GLM表述为如下形式：Y＝Xβ+ε..............(1)其中，Y为观测信号时间序列，X为设计时间序列，β为系数矩阵，ε为误差，则系数矩阵β的初始值的计算公式如下：β＝(XTX)-1XTy………(2)式中y是Y的一个分量。B、计算观测信号时间序列Y与设计时间序列X间的范式距离Fdist:Fdist(Y,Xβ)＝norm(Y-Xβ)*norm(|fft(Y)|-|fft(Xβ)|*w).....(3)其中，norm()表示二范数，fft()表示快速傅立叶变换，|·|表示取绝对值，w为权重向量，其j分量w(j)由下列公式计算：w(j)＝1/(e-ηj-eθj)...........(4)w＝w/max(w)............(5)其中η为正向感兴趣频率常数，θ为负向感兴趣频率常数，它们共同通过公式(5)来确定“最感兴趣频率”的位置。C、设置初始值定义初始记录器bold＝bmin＝β，设置ε＝0.1，Δ＝1，设置迭代计算器t＝0，sj＝Fdist(y,x*β)，这里j＝1。D、设置随机方向正确次数记录器设置随机方向正确次数记录器k＝0。E、随机获取新的记录点更新步长Δ＝Δ/2，在当前记录器bold邻域内，从Nm(bold,Δ)分布中随机获取新的记录点bnew，使得bnew满足(bmin-bold)*(bnew-bold)+ε>0。F、计算新的范式距离Fdistsnew＝Fdist(y,x*bnew)............(6)其中x是X的一个分量。G、如果snew<sj，则k增加1位。H、如果k>v*p，返回步骤D。这里v是采样数，p是0～0.5之间的随机数。I、更新sj+1＝snew；βj+1＝bmin；j增加1位。J、转步骤D，迭代直到j到达迭代最大次数，所述的最大次数为90-100次。K、输出β，完成Y与X的β匹配。与现有技术相比，本专利技术具有以下有益效果：1、频域模变换：如公式(3)所示，本专利技术在fMRI时间序列匹配中将时域信号经快速离散傅立叶模变换后得到频域序列，从而完全消除相位信息，以达到消除要匹配的fMRI时间序列间相位差的目的。本专利技术与当前相位校正的方法相比，更为简单，消耗的计算量少。2、加权约束：考虑离散傅里叶变换，高频部分的离散傅里叶变换相当于高频噪声，本专利技术对离散傅里叶变换进行加权约束，如公式(4)-(5)所示，以降低高频部分的影响，优先考虑了低频有效信号，更为明确地确定“最感兴趣频率”的位置，从而降低激活检验结果出现“伪阳性”体素的可能性。3、变步长启发式寻优：在fM本文档来自技高网...

【技术保护点】
一种功能核磁共振时间序列匹配方法，其特征在于：包括以下步骤：A、计算系数矩阵β的初始值将广义线性模型GLM表述为如下形式：Y＝Xβ+ε……………(1)其中，Y为观测信号时间序列，X为设计时间序列，β为系数矩阵，ε为误差，则系数矩阵β的初始值的计算公式如下：β＝(XTX)‑1XTy………(2)式中y是Y的一个分量；B、计算观测信号时间序列Y与设计时间序列X间的范式距离Fdist:Fdist(Y,Xβ)＝norm(Y‑Xβ)*norm(|fft(Y)|‑|fft(Xβ)|*w)…..(3)其中，norm()表示二范数，fft()表示快速傅立叶变换，|·|表示取绝对值，w为权重向量，其j分量w(j)由下列公式计算：w(j)＝1/(e‑ηj‑eθj)………..(4)w＝w/max(w)…………(5)其中η为正向感兴趣频率常数，θ为负向感兴趣频率常数，它们共同通过公式(5)来确定“最感兴趣频率”的位置；C、设置初始值定义初始记录器bold＝bmin＝β，设置ε＝0.1，Δ＝1，设置迭代计算器t＝0，sj＝Fdist(y,x*β)，这里j＝1；D、设置随机方向正确次数记录器设置随机方向正确次数记录器k＝0；E、随机获取新的记录点更新步长Δ＝Δ/2，在当前记录器bold邻域内，从Nm(bold,Δ)分布中随机获取新的记录点bnew，使得bnew满足(bmin‑bold)*(bnew‑bold)+ε>0；F、计算新的范式距离Fdistsnew＝Fdist(y,x*bnew)…………(6)其中x是X的一个分量；G、如果snew<sj，则k增加1位；H、如果k>v*p，返回步骤D；这里v是采样数，p是0～0.5之间的随机数；I、更新sj+1＝snew；βj+1＝bmin；j增加1位；J、转步骤D，迭代直到j到达迭代最大次数，所述的最大次数为90‑100次；K、输出β，完成Y与X的β匹配。...

【技术特征摘要】
1.一种功能核磁共振时间序列匹配方法，其特征在于：包括以下步骤：
A、计算系数矩阵β的初始值
将广义线性模型GLM表述为如下形式：
Y＝Xβ+ε……………(1)
其中，Y为观测信号时间序列，X为设计时间序列，β为系数矩阵，ε为误
差，则系数矩阵β的初始值的计算公式如下：
β＝(XTX)-1XTy………(2)
式中y是Y的一个分量；
B、计算观测信号时间序列Y与设计时间序列X间的范式距离Fdist:
Fdist(Y,Xβ)＝norm(Y-Xβ)*norm(|fft(Y)|-|fft(Xβ)|*w)…..(3)
其中，norm()表示二范数，fft()表示快速傅立叶变换，|·|表示取绝对值，w
为权重向量，其j分量w(j)由下列公式计算：
w(j)＝1/(e-ηj-eθj)………..(4)
w＝w/max(w)…………(5)
其中η为正向感兴趣频率常数，θ为负向感兴趣频率常数，它们共同通过
公式(5)来确定“最感兴趣频率”的位置；<...

【专利技术属性】
技术研发人员：刘洪波，陈亮，张维石，冯士刚，
申请(专利权)人：大连海事大学，
类型：发明
国别省市：辽宁;21