基于多尺度特征提取的t 检验方法分析fMRI数据

支联合 周文刚 谭素敏

1(周口师范学院植物遗传与分子育种重点实验室, 河南 周口 466001)2(周口师范学院计算机科学学院, 河南 周口 466001)3(周口师范学院校医院, 河南 周口 466001)

基于多尺度特征提取的t检验方法分析fMRI数据

支联合1*周文刚2谭素敏3

1(周口师范学院植物遗传与分子育种重点实验室, 河南 周口 466001)2(周口师范学院计算机科学学院, 河南 周口 466001)3(周口师范学院校医院, 河南 周口 466001)

采用基于多尺度特征提取的相关分析方法来分析fMRI数据性能优良，但需要已知实验响应信号的精确波形，因此适用范围受到限制。研究基于多尺度特征提取的t检验新方法分析fMRI数据，以利用前者去噪彻底和后者不需要已知实验响应信号的特性，并用视觉fMRI组块型和事件相关型试验进行性能验证，且与基于多尺度特征提取的相关分析和目前权威的SPM8方法进行对比。分析结果显示，3种方法检测的激活区位置相同，但激活体素个数分别是2945、3811和1454。结果表明，由于综合了多尺度特征提取与t检验各自的优点，新方法的适用范围较基于多尺度特征提取的相关分析更广，分析性能较SPM8更好。

小波变换；特征提取；功能磁共振成像；t-检验；数据分析

引言

设计性能优良的数据分析方法，是充分发挥fMRI技术探索人脑功能[1]的重要条件。通常， 根据是否需要已知实验响应信号(paradigm responsive signal, PRS)这一标准，fMRI数据分析方法可分为模型驱动和数据驱动两类，其中前者需要PRS。现有方法中，基于一般线性模型(generalized linear model, GLM)的统计参数图(statistical parametric mapping, SPM)[2]是公认的权威方法之一，其分析结果常作为“金标准”，目前最成熟的版本是SPM8。

小波变换是现代信号分析技术之一，在fMRI数据干扰去除和激活检测等生物医学信号处理中有着广泛的应用[3-4]。比如，利用阈值法去除fMRI数据里的高频噪声[5]，基于多尺度分析特性去除fMRI数据里的低频漂移[6]，在小波域检测激活[7]，等等。

上述方法[2, 5-7]虽各有优点，但都不能同时去除高频噪声和低频漂移。为此，笔者设计了多尺度特征提取方法(multiscale feature extraction, MFE)，同时去除这两种性质不同的干扰成分[8]。但是，方法中用来检测激活的相关分析需要已知PRS的精确波形,许多fMRI试验无法满足这一条件[9]，应用上受到限制。鉴于此，本课题利用视觉fMRI组块型和事件相关型试验研究基于MFE的t检验方法，以拓展MFE的适用范围。为简单，文中将新方法称为wttest，原有利用相关分析的MFE方法称为wcor。

1 材料和方法

1.1 视觉fMRI试验数据

本研究采用文献[8]的数据。两名被试分别参加了视觉组块型和事件相关型试验,均签署知情同意书。组块型任务刺激为棋盘格闪光(持续10 s，频率8 Hz)，静息刺激为空白屏幕。分别在试验第82、122、162 和202 s给予任务刺激，其余时间静息。1.5 T Siemens 扫描机共获取全脑120帧BOLD/EPI 图像，每帧图像大小为64×64×20体素，体素大小3.4 mm×3.4 mm×6 mm(层间距1.2 mm)。扫描重复时间2 s，事件相关型试验的任务刺激持续时间2 s，其余与组块型相同。

在分析前，上述数据用SPM8进行了格式转换、运动校正、标准化到MNI解剖空间和高斯空间平滑(半高宽9 mm×9 mm×12 mm)等预处理。

1.2 方法

1.2.1 基于MFE的t检验方法

基于MFE的t检验方法，即wttest，与文献[8]的wcor相仿，包括3个步骤。

步骤1：构造特征提取矩阵。先由Matlab软件小波变换深度函数，根据体素数据长度，计算变换深度参数，将此参数加1作为小波变换深度J。再构造行数和列数都等于体素数据长度的单位矩阵，对矩阵每一列分别进行J级小波分解、剔除最高频(高频噪声)及最低频尺度(低频漂移)、小波重构共3项运算，得到特征提取矩阵。

步骤2：特征提取，即去除体素数据里的干扰。将fMRI数据排成矩阵形式，其中每个体素时间序列数据构成矩阵的一列，与特征提取矩阵相乘。

步骤3：t检验。去除干扰后，将体素时间序列数据中分别对应任务状态和静息状态的数据点进行双样本t检验，并将得到的t值用SPM8软件自带函数转化为Z值。

显然，新方法不涉及PRS的具体波形。

1.2.2 视觉试验数据处理

除wttest外，笔者还用wcor和SPM8分析以进行对比。为理解上的方便，按SPM8、wcor和wttest的顺序阐述处理流程。

SPM8采用该软件提供的标准流程，先用实验刺激函数与双伽马函数模拟的血液动力学函数卷积获得PRS，基于离散余弦变换用缺省的1/128 Hz的截止频率去除低频漂移。再用GLM建模，估计超参数及参数，用t检验逐体素地计算任务与静息两状态的t值。最后用该软件自带函数，将t值转化为Z值。

关于wcor的详细阐述见文献[8]。首先，对比PRS的频谱和PRS分解到各个小波尺度的频谱，确定PRS存在于哪些平稳小波尺度，即特征尺度；文献[8]基于sym10小波确定为d1234，这样就剔除了低频漂移尺度a4。其次，对单位矩阵的每一列提取d1234里的成分，构建特征提取矩阵。再次，用特征提取矩阵与排成矩阵形式的体素数据相乘，去除低频干扰成分。最后，将PRS与特征提取矩阵相乘，再与去除了低频干扰的体素数据进行相关分析，通过Fisher’s Z变换将相关系数转化为Z值。

因与wcor对比，故wttest也采用sym8小波和平稳小波变换，然后依据本文1.2.1节方法分析数据。其中，构造特征提取矩阵时计算得到J=4，但没有剔除最高频尺度d1，这样对wcor和SPM更公平，因为也都没有去除这种噪声。t检验时，虽然组块型试验任务刺激在第41～45、61～65、81～85、101～105次扫描时呈现，但考虑到血液动力学响应的延迟，而且PRS是刺激函数与血液动力学函数卷积的结果，因此任务样本数据点选取为42～49、62～69、82～89、102～109。同样，事件相关型数据的任务样本点选取为42～45、62～65、82～85和102～105。

为判断激活，需要选择一定的显著性水平作为阈值。由于组块型刺激的激活强度和范围较事件相关型为大，故前者选定阈值为P< 0.001，后者为P< 0.05。因多重统计检验的原因，这些阈值都由SPM8基于高斯随机场理论修正。为了直观，按放射学约定的方式，随机选取视觉皮层所在的轴向坐标z=-8 mm、-4 mm和0 mm上的激活体素显示在MNI解剖图上。在本研究中，运算都以Matlab7.8软件为平台完成。

图1 3种方法分析被试1组块型数据的激活图Fig.1 The activation maps generated by the three methods for the block-designed data of subject 1

2 结果

图1和图2显示了3种方法分析2个被试组块型数据的激活图。由图可知，3种方法检测的激活区位置一样，都在距状沟、舌回等视觉区，因此特异度相同。但wcor检测的最大，wttest次之，SPM8的最小。进一步分析可知，体素个数分别为2 535、2 029和1 152，所以两种MFE的方法较SPM8更灵敏。

图3和图4显示了3种方法分析2个被试事件相关型数据的激活图。同样的结论也可以从这两个图中得出，但更明显，因为wcor和wttest检测的激活体素个数分别是1 276和916，远多于SPM8的302。

图4 3种方法分析被试2事件相关型数据的激活图Fig.4 The activation maps generated by the three methods for the event-related data of subject 2

3 讨论和结论

本研究表明，MFE显著提高了后续统计方法分析fMRI数据的性能，因为基于这一方法的相关分析和t检验远优于SPM8。尽管3种方法都去除了低频漂移，但SPM8依据的是离散余弦变换，变换的基函数是三角函数，其定义域为无穷大，没有局域性质，只适合分析稳态信号。相反，MFE去除干扰依据的是小波变换，变换基函数是小波，具有局域特性，能够精确捕捉非稳态信号的瞬变信息，非常适合处理像fMRI这样的非稳态信号数据[3]。实际上，本研究中MFE还能非常容易地去除fMRI数据里的高频噪声，只需在构造特征提取矩阵的同时剔除掉最高频尺度d1，所以如果这样的话，去除干扰的效果会更好。更重要的是，本研究给出的去除干扰的方法与文献[8]的不同，选择特征尺度不依赖PRS，既具有数据驱动性质，也是一个较为普适的方法，所以对于fMRI数据去除干扰这一重要课题具有一定的参考价值。

本研究还表明，在MFE去除干扰框架下，相关分析比t检验的分析性能要好，因为其灵敏度更高。究其原因，t检验需要准确界定任务样本和对照样本的数据点。但对于fMRI试验，虽然实验刺激函数是方波，血液动力学函数却不是，因此两者卷积后得到的PRS也不是，这样就很难界定两个样本的数据点，从而影响t检验的分析结果。相反，相关分析依赖于PRS的精确波形，通过血液动力学函数与实验刺激函数卷积在数学上很容易实现，因为实验刺激函数是方波，任务刺激时间点和静息时间点能根据实验设计精确界定。尽管如此，由于t检验不需要PRS，所以能够处理无法模拟PRS波形的fMRI试验数据，而相关分析和SPM8对此则无能为力。而且，正如所验证的那样，由于基于MFE良好的去干扰性能,t检验方法检测激活体素的性能优于SPM8。所以，MFE去噪框架下的相关分析和t检验互有优缺点，需要在实际中针对具体问题，选择最合适的方法。

此外，还需要指出的是，虽然本方法分析fMRI数据不需要PRS，因此具有数据驱动性质，但仍然需要用已知实验刺激函数，以确定任务样本数据点和对照样本数据点，因此不是完全意义上的数据驱动方法。还有，MFE会出现混频现象[10]，而这一现象对t检验分析结果的影响有待研究。最后，本研究通过t检验检测激活体素，而正如上一段讨论的那样，怎样界定最优的任务样本数据点和对照样本数据点也是需要解决的关键问题。所以，本研究还有进一步改进的空间。

总之，本研究利用MFE去噪彻底和t检验方法检测激活体素不需要PRS精确波形的优点，将两者组合起来分析fMRI数据，设计了一种具有数据驱动性质的新方法。新方法的适用范围较基于MFE的相关分析方法更广，分析性能较权威方法SPM8更好，而且去噪方法对于fMRI数据分析具有普适性。

[1] 尧德中, 罗程, 雷旭, 等. 脑成像与脑连接 [J]. 中国生物医学工程学报, 2011, 30(1): 6-10.

[2] Friston KJ, Holmes AP, Worsley KJ,etal. Statistical parametric maps in functional imaging: a general linear approach [J]. Hum Brain Mapping, 1995, 2(4):189-210.

[3] Bullmore E, Fadili J, Maxim V,etal. Wavelets and functional magnetic resonance imaging of the human brain [J]. NeuroImage, 2004, 23(1): S234- S249.

[4] 王永轩, 邱天爽, 刘蓉. 基于小波分析方法的脑电诱发电位单导少次提取 [J]. 中国生物医学工程学报, 2011, 30(1): 34-39.

[5] Wink AM, Roerdink JBTM. Denoising functional MR images: a comparison of wavelet denoising and Gaussian smoothing [J]. IEEE Trans Med Imaging, 2004, 23 (3): 374- 387.

[6] Meyer FG. Wavelet-based estimation of a semiparametric generalized linear model of fMRI time-series [J]. IEEE Trans Med Imaging, 2003, 22(3):315-322.

[7] Van De Ville D, Seghier ML, Lazeyras F,etal. WSPM: wavelet-based statistical parametric mapping [J]. NeuroImage, 2007, 37(4): 1205-1217.

[8] 支联合, 谭素敏, 支羽光. 基于快速平稳小波变换的特征提取方法分析fMRI数据 [J]. 中国生物医学工程学报, 2012, 31(4): 620-624.

[9] 白保东, 刘健, 郭红宇. 基于独立成分分析的时间自相关方法在功能磁共振激活区自适应提取中的应用 [J]. 中国生物医学工程学报, 2014, 33(2): 194-201.

[10] 杨建国. 小波分析及其工程应用 [M]. 北京: 机械工业出版社, 2005: 107-128.

t-Test Analysis of FMRI Data Based on Multiscale Feature Extraction

Zhi Lianhe1*ZHou Wengang2Tan Sumin3

1(KeyLaboratoryofPlantGeneticsandMolecularBreeding,ZhoukouNormalUniversity,Zhoukou466001,Henan,China)2(CollegeofComputerScienceandTechnology,ZhoukouNormalUniversity,Zhoukou466001,Henan,China)3(TheHospitalAffiliatedtoZhoukouNormalUniversity,Zhoukou466001,Henan,China)

wavelet transform; feature extraction; fMRI;t-test; data analysis

10.3969/j.issn.0258-8021. 2015. 04.015

2014-10-10， 录用日期:2015-05-04

国家自然科学基金(31272168)； 河南省科技厅基础与前沿项目(132300410276)

R445.2

D

0258-8021(2015) 04-0492-04

*通信作者(Corresponding author)， E-mail: zhilianhe2008@163.com