一种基于fMRI探测动态脑自发活动的数据方法
【专利摘要】该发明公开了一种基于fMRI探测动态脑自发活动的数据方法,属于生物医学图像模式识别技术领域，具体涉及基于功能磁共振时间序列的动态脑自发活动探测方法。本发明根据脑自发活动的动态变化特性，利用短时间间隔内脑信号是平稳(伪平稳)属性，提出一种组合滑动窗(sliding?window)及信号振幅(ALFF)方法，并结合给变异程度(Variance)值，给出了动态脑自发活动定量化计算方法。本发明的方法具有鲁棒性高，稳定性强等特点，对于磁共振数据模式识别等领域脑自发活动的动态性探测提供了新的有效技术。
【专利说明】
一种基于f MR I探测动态脑自发活动的数据方法
技术领域
[0001]本方法属于生物医学图像模式识别技术领域，具体涉及基于功能磁共振时间序列的动态脑自发活动探测方法。
【背景技术】
[0002]大脑处于静息态是指大脑无特定外在任务刺激下的“静息态”脑活动，这种“静息态”活动仍然持续地进行血液循环和消耗氧。研究发现即使是在静息状态，大脑脑区间的血氧依赖功能磁共振(BOLD-fMRI)信号的振荡保持着高度同步和相关。大脑静息状态下固有的自发神经活动具有特有的时空组织模式。
[0003]目前，对静息状态下的大脑自发脑活动的模式探测主要采用低频振荡方法、fMRI信号强度变异及功率谱的方法。这些方法都是在时域计算信号的变异程度或者在频域计算特殊频段的功率谱或者幅度。这几种方法在本质上是一致的。利用这些方法，我们可以探测发现大脑自发神经活动是大脑在最原始的默认状态下，自主产生的最基本的神经活动，主要分布于大脑的网络。这一指标稳定可靠，可以作为良好的生物学指标、从多个方面对大脑的活动、疾病、塑形、发育及进化进行研究，促进加深对大脑机理的理解。
[0004]但是，上述算法计算是均利用了fMRI整体时间序列(一般5?8分钟)，忽略了大脑本身的动态性(Dynamic)。近来，研究表明，即使是在静息态采集或者执行任务这一特定时间段内，大脑的功能连接也不是稳态，而是一种非稳态，其连接模式会随着时间变化。目前，还没有稳定可靠的方法用以探测探测动态脑自发活动。所以，本发明采用一种滑动窗信号振幅(dynamic amplitude of low-frequency fluctuat1ns ,dALFF)，代替常规脑自发脑活动方法采用所有时间信息，而忽略这种中尺度时间段内的动态变化，用以揭示脑动态脑自发活动。
【发明内容】
[0005]针对目前尚无方法实现动态脑自发活动探测，本发明提供了一种针对功能磁共振数据的动态脑自发活动的数据方法，利用组合滑动窗(si iding-window)及信号振幅(ALFF)计算时间过程的脑自发活动，最后采用变异程度(Variance)刻画动态脑自发活动。
[0006]为实现上述目的，本发明的拘束方案为一种基于fMRI探测动态脑自发活动的数据方法，该方法具体步骤如下:
[0007]步骤1:获取原始功能磁共振(fMRI)数据矩阵，并将原始功能磁共振数据矩阵转化为V*T矩阵，其中V为体素数目，1~为时间点长度；
[0008]步骤2:针对每个体素，可得其时间序列(1*T的向量)。并给定滑动窗长度L(1彡L彡T)，滑动步长n(l<n彡L)，将每个体素时间序列1*T的向量划分为[(T-L)/n] + l段的1*L向量，其中[]代表向下取整函数；
[0009]步骤3:利用傅里叶变化，计算第I个体素中各滑动窗j(j= l，2，…[(T-L)/n] + l)内时间序列为I *L的向量的频谱幅度；
[0010]步骤4:对滑动窗内时间序列进行滤波，选择符合fMRI生理信号的频段；
[0011 ]步骤5:利用傅里叶反变换将频域信号变成时域信号；
[0012]步骤6:计算其信号幅度的平均值，作为该滑动窗内的脑自发活动强度；
[0013]步骤7:采用步骤3?步骤6的方法计算其余各体素中各滑动窗内的脑自发活动强度；
[0014]步骤8:根据步骤6和步骤7获得的脑自发活动强度，利用变异程度(Variance)计算每个体素的动态脑自发活动变异程度；
[0015]步骤9:循环每个体素i(i= l，2，…，V)，进而得到全脑的每个体素的动态脑自发活动变异程度脑图谱。
[0016]进一步的，所述步骤2中的滑动窗长度为40-80。
[0017]进一步的，所述步骤2中滑动窗选择为汉明窗或变长滑动窗。
[0018]进一步的，所述步骤4中滤除滑动窗内频率在0.01-0.08Hz之外的频段。
[0019]本发明根据脑自发活动的动态变化特性，利用短时间间隔内脑信号是平稳(伪平稳)属性，提出一种组合滑动窗(sliding-window)及信号振幅(ALFF)方法，并结合给变异程度(Variance)值，给出了动态脑自发活动定量化计算方法。本发明的方法具有鲁棒性高，稳定性强等特点，对于磁共振数据模式识别等领域脑自发活动的动态性探测提供了新的有效技术。
【附图说明】
[0020]图1为特征区域(脑区)位置图。
[0021]图2为特征区域的时间信号。
[0022]图3为在每个滑动窗内的频谱及信号幅度。
[0023]图4为整个时间序列计算得到的脑自发活动模式。
[0024]图5为利用滑动窗计算得到的动态脑自发活动变异程度模式。
[0025]图6为特征区域的动态脑自发活动变异程度和滑动窗长度关系。
[0026]图7为滑动窗窗宽范围优化。
[0027]图8为不同随机序列的动态动态脑自发活动计算。
[0028]图9为不同滑动步长的动态动态脑自发活动计算。
【具体实施方式】
[0029]下面结合附图和实施例对本发明的【具体实施方式】做进一步详细地描述，以下实施例用于说明本发明，但不用来限制本发明的范围。
[0030]一种基于fMRI探测动态脑自发活动的方法，具体实施步骤如下:
[0031]步骤A:仿真数据构造。
[0032]生成一个32X 32体素的fMRI图像，其中每一个体素包含了一个240个时间点的时间序列。在A、B、C、D区域构造各自特定特征的时间序列(左上为A区域，右上为B区域，左下为C区域，右下为D区域)。特征区域位置标识见图1，特征信号显示于图2。为模拟功能磁共振信号，我们将时间采样间隔定为2s，采样频率为0.5Hz。fMRI图像局部区域内相似性较高，因此在A、B、C、D区域内各模拟了9个信号，这9个信号具有98%的相似性。区域A、B的信号设计为规律组块间隔设计，组块长度为40个时间点，共6个组块。区域C、D的信号设计为随机短组块设计，组块长度为20。仿真信号为正弦信号叠加。区域A的信号设计方案为组块同频率段(0.01-0.08Hz)，不同幅度(Al = I，A2 = 0.5)，高幅度、低幅度间隔出现。区域B的信号设计方案为组块同幅度(A=I)，不同频率段(0.01-0.08Hz、0.1-0.2Hz)，低频段、高频段间隔出现。区域C的信号设计方案为同频率段(0.01-0.08Hz)，不同幅度(Al = I，A2 = 0.5)，高幅度、低幅度随机出现。区域D的信号设计方案为组块同幅度(A=I)，不同频率段(0.01-0.08Hz、0.1-0.2Hz)，低频段、高频段随机出现(图2，图3)。生成的信号均叠加了 10%的白噪声，共生成9X4X240(36 X 240) (S卩 V = 36，T = 240)的有效信号。
[0033 ] 步骤B:全时间段(240时间点)信号振幅(ALFF)计算
[0034]对于每个体素，利用傅里叶变换，计算全时间段240时间点的信号频谱，带通滤波截取低频段(0.01-0.08Hz)，再利用傅里叶逆变换，返回时域空间，计算滤波后时域信号平均幅度(ALFF)，并在体素位置标记(图4)。
[0035]步骤C:动态信号振幅(dALFF)计算
[0036]对于给定体素，给定给滑动窗长度(30时间点)(即L= 30)，滑动步长为I时间点，逐段截取体素信号。利用傅里叶变换，计算窗内信号的频谱，带通滤波截取低频段(0.01-
0.08Hz)，并利用傅里叶逆变换，回到时域空间，并计算滤波后时域信号平均幅度(ALFF)(图
3)。滑动窗口，并遍历整段信号，记录滑动序列动态变化ALFF序列。图3共分为四个部分，每个部分均再分为三个子图，原始信号图，组块(N时间点)频域分布图，以及dALFF序列图。原始信号子图中纵向线为组块分割标记，框线为滑动时间窗，三条水平线为信号平均值线及加减一倍方差线，组块频域分布子图中纵线为组块分割标记，框线为0.01-0.08Hz范围。dALFF序列子图中框线为滑动窗，点为滑动窗中点，代表了滑动窗内平均ALFF值。区域A在图3左上部分，其原始信号频率范围固定(0.01-0.08Hz)，组块大小为40时间点，信号幅度呈现为高-低-高-低-高-低分布，dALFF序列子图也出现了对应的峰值和谷值变化。区域B在图3的右上部分，其原始信号信号幅度固定，组块大小为40时间点，信号频率呈现为低(0.01-0.08Hz)-高(0.1 -0.2Hz) _低-高_低-高分布，dALFF序列子图也出现了对应的峰值和谷值变化。区域C在图3的左下部分，其原始信号频率范围固定(0.01-0.08Hz)，组块大小为20时间点，信号幅度随机(低-低-低-高-高-低-低-高-高-高-低-高)，dALFF序列子图也出现了较明显对应的峰值和谷值变化。区域D在图3的左下部分，其原始信号幅度固定，组块大小为20时间点，信号频率随机(低-低-高-高-高-低-低-高-高-低-高-低)，dALFF序列子图也出现了较明显对应的峰值和谷值变化。
[0037]步骤D:动态ALFF变异程度计算
[0038]对步骤C得到的滑动序列动态变化ALFF序列(dALFF)，逐体素计算dALFF的变异系数，并将变异系数在体素位置进行标记(图5)。
[0039]步骤E:不同滑动窗长度dALFF变异性变化趋势
[0040]采用改变滑动窗长度(10:5:50)的方式，重复步骤C、D，并计算不同滑动窗情况下，信号dALFF的变异性(图6)。滑动窗长度越大，变异性变小，检测能力受到限制，但滑动窗长度很小的时候起稳定性受到考验，结合BOLD-fMRI的自身特征，20-40时间点(40-80s)的窗宽适合探测特异性变化(图7)。
[0041]步骤F:不同构造的随机序列(区域C、D)的dALFF序列计算[0042 ]采用改变随机序列区块长度(I O、12、15、24、30时间点)的方式，用固定30时间点的滑动窗，重复步骤C(图8)。
[0043]步骤G:不同滑动步长dALFF序列计算
[0044]采用固定窗宽30时间点，不同步长(1:5:41)的滑动窗模式，重复步骤C(图9)。
[0045]综上所述，本发明提出的方法在不同信噪比，不同窗长度，不同滑动步长选择时，仍然能够有效探测到脑自发活动模式的动态变化特性，有效克服了现有方法存在的局限性。
【主权项】
1.一种基于fMRI探测动态脑自发活动的数据方法，该方法具体步骤如下: 步骤I:获取原始功能磁共振(fMRI)数据矩阵，并将原始功能磁共振数据矩阵转化为V*T矩阵，其中V为体素数目，1~为时间点长度； 步骤2:针对每个体素，可得其时间序列(1*T的向量)。并给定滑动窗长度L(KLST),滑动步长n(l<n<L)，将每个体素时间序列1*Τ的向量划分为[(T-L)/n]+l段的1*L向量，其中[]代表向下取整函数； 步骤3:利用傅里叶变化，计算第I个体素中各滑动窗j(j = l，2，…[(T-L)/n]+l)内时间序列为1*L的向量的频谱幅度； 步骤4:对滑动窗内时间序列进行滤波，选择符合fMRI生理信号的频段； 步骤5:利用傅里叶反变换将频域信号变成时域信号； 步骤6:计算其信号幅度的平均值，作为该滑动窗内的脑自发活动强度； 步骤7:采用步骤3?步骤6的方法计算其余各体素中各滑动窗内的脑自发活动强度；步骤8:根据步骤6和步骤7获得的脑自发活动强度，利用变异程度(Variance)计算每个体素的动态脑自发活动变异程度； 步骤9:循环每个体素i (i = I，2，…，V)，进而得到全脑的每个体素的动态脑自发活动变异程度脑图谱。2.如权利要求1所述的一种基于fMRI探测动态脑自发活动的数据方法，其特征在于:所述步骤2中的滑动窗长度为40-80。3.如权利要求1所述的一种基于fMRI探测动态脑自发活动的数据方法，其特征在于:所述步骤2中滑动窗选择为汉明窗或变长滑动窗。4.如权利要求1所述的一种基于fMRI探测动态脑自发活动的数据方法，其特征在于:所述步骤4中滤除滑动窗内频率在0.01-0.08Hz之外的频段。
【文档编号】A61B5/055GK106037741SQ201610524807
【公开日】2016年10月26日
【申请日】2016年7月4日
【发明人】廖伟, 段旭君, 陈华富
【申请人】电子科技大学