CN103996196B - 一种基于多变量的dti图像分析方法 - Google Patents

Abstract

本发明公开了一种基于多变量的DTI图像分析方法，用于疾病大脑机制研究中对脑部白质病变区域确认和提取，从而为临床治疗提供影像学证据，具体步骤为：对数据进行预处理，包括去涡流、头动矫正、去脑壳、弥散拟合、构建白质骨架；对预处理后的数据进行特征提取,通过置换检验并以年龄为协变量得到病人组和正常人组的显著差异区域，并分别计算显著差异区域中特定变量的平均值得到特征值；用留一法进行交叉验证，判断是否满足停止条件，若不满足则计算每个特征的平均权值，并去除平均权值最小的特征，直到满足停止条件为止；最后得到的脑部区就是脑部病变区域。本发明利用DTI成像模式，通过多变量的研究方法，为临床治疗找到了病变区域提供了影像学证据。

Description

一种基于多变量的DTI图像分析方法

技术领域

本发明涉及一种图像分析方法，具体涉及一种基于多变量的DTI图像分析方法。

背景技术

扩散张量成像(DTI)是一种非侵入性的能提供活体内水分子扩散运动的成像技术，它能检测传统MRI所不能观察到的组织的微观变化，是MR成像技术的重大突破。基于脑影像信息的模式分类是目前脑影像研究中的热点课题。利用图像分类方法，计算DTI图像具有某种属性的可能性大小，或者自动的判别图像的类别属性，是计算机辅助分析的一个重要应用。

文献“Alexander AL,Lee JE,et al.Diffusion tensor imaging of the corpuscallosum in Autism.Neuroimage.2007；34(1):61–73.”使用基于感兴趣区域的方法研究胼胝体在孤独症中的作用，但是这种方法需要有关某种病理或者病变区域的先验知识，所以没有很好地推广性。

文献“Ridgway GR,Henley SM,et al.Ten simple rules for reporting voxel-based morphometry studies.Neuroimage.2008；40(4):1429–1435.”使用基于体素的双样本t检验统计方法研究病人和正常人之间的组间差异。这种方法是在假设变量符合正态分布的情况下仅仅考虑单个变量的影响没有同时考虑多个变量的影响，而且无法判断单个人是否是病人或者是正常人。

文献“Madhura Ingalhalikar,et al.Diffusion based Abnormality Markersof Pathology:Towards Learned Diagnostic Prediction of ASD.Neuroimage.2011；57(3):918–927”将地图集中的各向异性(FA)和平均扩散率(MD)的值作为特征，然后将所得特征加入支持向量机(SVM)中，通过留一法不断地选择可以使支持向量机得到最好的正确率和推广性的特征。然而，这种基于地图集的方法无法提取地图集下子区域的相关变量作为特征，这样就无法找到子区域中病变的区域。同时，这种方法没有考虑年龄因素对白质的影响。

发明内容

针对现有技术的不足，本发明旨在提供一种推广性更强的基于多变量的DTI图像分析法，本发明首先经过预处理和弥散张量计算得到各向异性值(FA)、平均扩散率值(MD)、径向弥散率值(RD)，针对计算得到的三个量进行置换检验后得到病人与正常人间的组间差异，然后将差异区域的各向异性值、平均扩散率值、径向弥散率值分别平均后加入线性支持向量机中用留一法训练支持向量机，直至找到使分类器准确率高、推广性好的特征量，最终得到的这些特征量所在的区域就很有可能是病变区域。

为了实现上述目的，本发明采用如下技术方案：

一种基于多变量的DTI图像分析方法，包括图像数据的采集，其中，所述图像数据分为正常人组和病人组两组，通过3T的磁共振设备采集，再采用弥散加权序列得到弥散图像数据，所述方法还包括如下步骤：

步骤1，将所述弥散图像数据进行预处理，具体如下：

步骤1.1，用基于Linux的FSL软件对所述弥散图像数据进行去涡流处理和头动矫正；

步骤1.2，对经所述步骤1.1处理得出的图像进行去脑壳处理；

步骤1.3，对经所述步骤1.2处理得出的图像进行弥散张量拟合得出各向异性图像、平均扩散率图像、径向弥散率图像；

步骤1.4，将所述各向异性图像通过非线性配准方法配准到标准空间；

步骤1.5，将所有配准到标准空间的各向异性图像进行平均得到平均各向异性图像；

步骤1.6，将所述平均各向异性图像进行骨架化，得到白质骨架；

步骤1.7，将步骤1.3得出的每个个体的各向异性图像、平均扩散率图像、径向弥散率图像分别投射到所得所述白质骨架上，得到每个个体的各向异性骨架图像、平均扩散率骨架图像、径向弥散率骨架图像；

需要说明的是，所述病人组中的图像数据采集自神经退化类的病人。

需要进一步说明的是，所述正常人组和病人组的量的最优比例为1:1。

步骤2，提取步骤1得出的所述各个个体的各向异性骨架图像、平均扩散率骨架图像、径向弥散率骨架图像的特征：

步骤2.1，将所述各个个体的各向异性骨架图像、平均扩散率骨架图像、径向弥散率骨架图像分别进行逐体素置换检验并以年龄作为协变量得到检验后的图像；

步骤2.2，将所述检验后的图像设置阈值，得到病人组和正常人组在白质骨架上各向异性值、平均扩散率值、径向弥散率值存在显著差异的区域；

步骤3，对步骤2所得的所述存在显著差异的区域进行多变量分析：

步骤3.1，分别对所述存在显著差异的区域中的各向异性值、平均扩散率值、径向弥散率值进行平均，得到所述存在显著差异的区域中的平均各向异性值、平均平均扩散率值、平均径向弥散率值；

步骤3.2，基于MATLAB软件将步骤3.1中所述平均的各向异性值、平均扩散率值、径向弥散率值作为特征输入到线性支持向量机中，通过留一法对线性支持向量机进行训练，最终得到特征所在区域，从而得到和病变有关的区域。

需要说明的是，所述步骤3.2中，通过留一法对线性支持向量机进行训练的具体实施步骤如下：

步骤3.2.1，用n表示数据样本中的个体总数，而各个体都有m个特征量，而各个体的类属性都是已知的，即病人或正常人；将得到的数据样本分为两组，一组是测试集，包含一个个体，一组是训练集，包括除测试集中所含个体外的所有人，共n-1个个体；

步骤3.2.2，用所述训练集训练所述线性支持向量机，得出所述训练后支持向量机：根据以下公式计算得到权值向量w，w为一个m维列向量，其中的每个元素对应一个特征量；

y_i(w^Tx_i+b)-1+ξ_i≥0

s.t.ξ_i≥0；

其中，γ是惩罚参数，用来实现算法复杂度和错分样本数的折中；ξ_i测量错分程度；y_i每个人的类属性；x_i为每个个体的特征向量；b为常数；

需要说明的是，b是通过拉格朗日乘子法将有约束的最小化问题转化为无约束的最小化化问题，求解该无约束优化得到，而ξ_i大于1时，则第i个人被错分，通过拉格朗日乘子法将有约束的最小化问题转化为无约束的最小化化问题，求解该无约束优化得到。

步骤3.2.3，用已知类属性的所述测试集来评估所述训练后的支持向量机的性能：用所述训练后支持向量机来判断所述测试集的类属性，所述训练后支持向量机会给出属性标签1或-1，其中1为病人，-1为正常人，通过所述训练后支持向量机得出的判断结果与所述测试集的实际类属性相比较，两者若一致，则所述训练后支持向量机分类正确，否则，则分类错误；

步骤3.2.4，重新将n个个体分为测试集和训练集，所述测试集包含一个个体，且该个体与上一轮的测试集中的个体不相同，余下的所有个体作为训练集，然后按照步骤3.2.2的方法训练所述线性支持向量机，得出训练后支持向量机，然后再按照所述步骤3.2.3的方法评估得出的训练后支持向量机的性能；重复步骤3.2.4n-1次后停止；

步骤3.2.5，将每个特征的n次权值求平均权值，并依据平均权值将特征由大到小进行排序，去除排序最低的特征量；

步骤3.2.6，重复步骤3.2.1至步骤3.2.4，接着执行步骤3.2.7；

步骤3.2.7，依据步骤3.2.6中重复步骤3.2.1至步骤3.2.4后所得出的该轮n次测试中的分类正确率和上一轮n次测试的分类正确率比较的结果判断是否停止：如果该轮n次测试的分类正确率大于或等于上一轮n次测试的分类正确率，则返回执行步骤3.2.5至步骤3.2.6，否则停止。

需要进一步说明的是，上一轮n次测试得到特征排序和一个分类正确率，去除排序最低的特征后，用余下的特征量进行下一轮n次测试，下一轮n次测试又得到特征排序和一个分类正确率，如果该轮n次测试的分类正确率大于或等于上一轮n次测试的分类正确率，则去除该轮中排序低的特征量，然后继续下一轮n次测试，如果该轮n次测试的分类正确率小于上一轮n次测试的分类正确率则终止，将上一轮n次测试中剩余的特征量作为最后选择的特征量。

本发明的有益效果在于：

1、使用置换检验并以年龄因素为协变量进行初步的特征提取，克服了双样本t检验对于变量要服从正态分布的假设，考虑了年龄对白质的影响；

2、在体素水平进行置换检验，克服了基于感兴趣区域(ROI)的方法需要有先验知识的缺点，同时克服了基于地图集方法无法观察地图集中子区域的缺点；

3、使用线性支持向量机对初选的特征进行二次筛选，去除不是由于病变引起的差异而是由于图像预处理或者噪声引起的差异，同时考虑了多个不同区域不同变量的相互作用，克服了t检验只考虑单个变量的缺点。

附图说明

图1为本发明的实施步骤示意图；

图2为图1中预处理的子流程图；

图3为图1中特征提取的子流程图；

图4为图1中留一法交叉验证的子流程图。

具体实施方式

以下将结合附图对本发明作进一步的描述，需要说明的是，本实施例以本技术方案为前提，给出详细的实施步骤，但并不限于本实施例。

如图1所示，所述一种基于多变量的DTI图像分析方法包括：对磁共振扫描后采集并采用弥散加权序列得到的弥散图像数据进行数据预处理，采集的图像数据分为病人组和正常人组，所述病人组中的病人是指神经退化类的病人，与所述病人相对应的为正常人；所述正常人组和病人组的量的比例为1:1；再采用置换检验提取特征，最后用留一法进行交叉验证得到可能的病变区域。具体步骤如下：

步骤1，数据预处理分析阶段(基于FSL软件)，如图2所示：

步骤1.1，对所述弥散图像数据进行涡流矫正，涡流矫正就是矫正在磁共振扫描过程中由于人身体切割磁场产生的电流以增强图像质量。

由于在磁共振扫描过程中，被试的呼吸、血流脉动等生理因素造成的头部运动在所难免，所以对涡流矫正过的数据进行头动矫正，也就是将一个序列的每一帧图像都和这个序列的第一帧图像进行配准，配准到同一坐标系下，以矫正头动。对齐后再对头动矫正后的数据进行人工检查，如果平移和旋转分别超过1毫米和1度就排除掉、不予分析。

步骤1.2，对头动矫正后的图像进行去脑壳，在数据分析中只关注大脑，所以要将脑壳部分去除。

步骤1.3，对去脑壳后的图像进行弥散张量拟合，从而得到弥散张量的3个特征值和3个特征向量。然后对得到的特征值进行算术运算得到各向异性图像、平均扩散率图像、径向弥散率图像。

步骤1.4，将得到各向异性图像配准到蒙特利尔神经病学研究所(MNI)空间，实验中存在多个被试，被试和被试之间的脑部大小形状存在一定的差异，为了后续的统计分析，须进行大脑形状的归一化，将被试的大脑配准归一化到标准的大脑模板上。

步骤1.5，将MNI空间上的所有个体的各向异性图像求平均得到平均的各向异性图像；

步骤1.6，由于只是分析白质，所以将平均的各向异性图像骨架化以去除灰质，这样就得到了白质骨架了。

步骤1.7，将每个个体的各向异性图像、平均扩散率图像、径向弥散率图像投射到白质骨架得到骨架化的各向异性图像、平均扩散率图像、径向弥散率图像。

步骤2，特征提取的步骤如下，如图3所示：

步骤2.1，将骨架化的各向异性图像、平均扩散率图像、径向弥散率图像分别进行逐体素置换检验并以年龄因素为协变量，并进行多重比较的校正，得到校正后的p值图；

步骤2.2，将p值取0.05，得到病人组和正常人组的骨架化的各向异性图像、平均扩散率图像、径向弥散率图像的存在显著差异的区域。

步骤3，对步骤2所得的所述存在显著差异的区域进行多变量分析：

步骤3.1，分别对所述存在显著差异的区域中的各向异性值、平均扩散率值、径向弥散率值进行平均，得到所述存在显著差异的区域中平均的各向异性值、平均扩散率值、径向弥散率值；

步骤3.2，基于MATLAB软件将步骤3.1中所述平均的各向异性值、平均扩散率值、径向弥散率值作为特征输入到线性支持向量机中，通过留一法对线性支持向量机进行训练，最终得到特征所在区域，从而得到和病变有关的区域。

如图4所示，为了验证分类器的推广性和准确性，采用留一法进行交叉验证。留一法交叉验证(基于MATLAB软件)的具体步骤如下：

步骤3.2.1，用n表示样本中的个体总数，而每个个体都有m个特征量，每个个体的类属性都是已知的，即病人或正常人，将得到的特征量分为两组，一组是测试集，只有一个个体，一组是训练集，包括除测试集中所含个体外的所有个体，也就是n-1个个体；

步骤3.2.2，用所述训练集训练所述线性支持向量机，得出训练后支持向量机：根据以下公式计算得到权值向量w，w为一个m维列向量，其中的每个元素对应一个特征量；

y_i(w^Tx_i+b)-1+ξ_i≥0

s.t.ξ_i≥0；

其中，γ是惩罚参数，用来实现算法复杂度和错分样本数的折中；ξ_i测量错分程度；y_i每个人的类属性；x_i为每个个体的特征向量；b为常数；

步骤3.2.3，用已知类属性的所述测试集来评估所述训练后支持向量机的性能：用所述训练后支持向量机来判断所述测试集的类属性，所述训练后支持向量机会给出属性标签1或-1，其中1为病人，-1为正常人，通过所述训练后支持向量机得出的判断结果与所述测试集的实际类属性相比较，两者若一致，则所述训练后支持向量机分类正确，否则，则分类错误；

步骤3.2.4，重新将n个个体分为测试集和训练集，所述测试集包含一个个体，且该个体与上一轮的测试集中的个体不相同，余下的所有个体作为训练集，然后按照步骤3.2.2的方法训练所述线性支持向量机，得出训练后支持向量机，然后再按照所述步骤3.2.3的方法评估得出的训练后支持向量机的性能；重复步骤3.2.4n-1次后停止；

步骤3.2.5，将每个特征的n次权值求平均权值，并依据平均权值将特征由大到小进行排序，去除排序最低的特征量；

步骤3.2.6，重复步骤3.2.1至步骤3.2.4，接着执行步骤3.2.7；

步骤3.2.7，依据步骤3.2.6中重复步骤3.2.1至步骤3.2.4后所得出的该轮n次测试中的分类正确率和上一轮n次测试的分类正确率比较的结果判断是否停止：如果该轮n次测试的分类正确率大于或等于上一轮n次测试的分类正确率，则返回执行步骤3.2.5至步骤3.2.6，否则停止。

对于本领域的技术人员来说，可以根据以上的技术方案和构思，作出各种相应的改变和变形，而所有的这些改变和变形，都应该包括在本发明权利要求的保护范围之内。

Claims (2)

1.一种基于多变量的DTI图像分析方法，包括图像数据的采集，其中，所述图像数据分为正常人组和病人组两组，通过3T的磁共振设备采集，再采用弥散加权序列得到弥散图像数据，其特征在于，所述方法还包括如下步骤：

步骤1，将所述弥散图像数据进行预处理：

步骤1.1，用基于Linux的FSL软件对所述弥散图像数据进行去涡流处理和头动矫正；

步骤1.2，对经所述步骤1.1处理得出的图像进行去脑壳处理；

步骤1.3，对经所述步骤1.2处理得出的图像进行弥散张量拟合，得出各向异性图像、平均扩散率图像、径向弥散率图像；

步骤1.4，将所述各向异性图像通过非线性配准方法配准到标准空间；

步骤1.5，将所有配准到标准空间的所述各向异性图像进行平均得到平均各向异性图像；

步骤1.6，将所述平均各向异性图像进行骨架化，得到白质骨架；

步骤1.7，将步骤1.3中得出的所述各向异性图像、平均扩散率图像、图像按个体分别投射到所述白质骨架上，得到每个个体的各向异性骨架图像、平均扩散率骨架图像和径向弥散率骨架图像；

步骤2，提取步骤1得出的每个个体的所述各向异性骨架图像、平均扩散率骨架图像和径向弥散率骨架图像的特征：

步骤2.1，将各个体的各向异性骨架图像、平均扩散率骨架图像、径向弥散率骨架图像分别进行逐体素置换检验并以年龄因素作为协变量得到检验后的图像；

步骤2.2，将所述检验后的图像设置阈值，得到所述病人组和所述正常人组在白质骨架上各向异性值、平均扩散率值、径向弥散率值存在显著差异的区域；

步骤3，对步骤2所得的存在显著差异的区域进行多变量分析：

步骤3.1，分别对所述存在显著差异的区域中的各向异性值、平均扩散率值、径向弥散率值进行平均，得到所述存在显著差异的区域中的平均各向异性值、平均平均扩散率值、平均径向弥散率值；

步骤3.2，基于MATLAB软件将步骤3.1中所述的平均各向异性值、平均平均扩散率值、平均径向弥散率值作为特征输入到线性支持向量机中，通过留一法对线性支持向量机进行训练，最终得出特征所在的区域，从而得到与病变有关的区域；

通过留一法对线性支持向量机进行训练的具体实施步骤如下：

步骤3.2.1，用n表示数据样本中的个体总数，而各个体都有m个特征量，而各个体的类属性都是已知的，即病人或正常人；将得到的数据样本分为两组，一组是测试集，包含一个个体，一组是训练集，包括除测试集中所含个体外的所有人，共n-1个个体；

步骤3.2.2，用所述训练集训练所述线性支持向量机，得出训练后支持向量机：根据以下公式计算得到权值向量w，w为一个m维列向量，其中的每个元素对应一个特征量；

$\begin{array}{cc}min& \frac{1}{2}{w}^{T}w+\mathrm{\γ}\underset{i=1}{\overset{n}{\Σ}}{\mathrm{\ξ}}_i\end{array}$

y_i(w^Tx_i+b)-1+ξ_i≥0

s.t. ξ_i≥0；

其中，γ是惩罚参数，用来实现算法复杂度和错分样本数的折中；ξ_i测量错分程度；y_i为每个人的类属性；x_i为每个个体的特征向量；b为常数；

步骤3.2.3，用已知类属性的所述测试集来评估所述训练后支持向量机的性能：用所述训练后支持向量机来判断所述测试集的类属性，所述训练后支持向量机会给出属性标签1或-1，其中1为病人，-1为正常人，通过所述训练后支持向量机得出的判断结果与所述测试集的实际类属性相比较，两者若一致，则所述训练后支持向量机分类正确，否则，则分类错误；

步骤3.2.4，重新将n个个体分为测试集和训练集，所述测试集包含一个个体，且该个体与前一次测试的测试集中的个体不相同，余下的所有个体作为训练集，然后按照步骤3.2.2的方法训练所述线性支持向量机，得出训练后支持向量机，然后再按照所述步骤3.2.3的方法评估得出的所述训练后支持向量机的性能；重复步骤3.2.4n-1次后停止；

步骤3.2.5，将每个特征的n次权值求平均权值，并依据平均权值将特征由大到小进行排序，去除排序最低的特征量；

步骤3.2.6，重复步骤3.2.1至步骤3.2.4，接着执行步骤3.2.7；

步骤3.2.7，依据步骤3.2.6中重复步骤3.2.1至步骤3.2.4后所得出的该轮n次测试中的分类正确率和上一轮n次测试的分类正确率比较的结果判断是否停止：如果该轮n次测试的分类正确率大于或等于上一轮n次测试的分类正确率，则返回执行步骤3.2.5至步骤3.2.6，否则停止；

所述病人组中的图像数据采集自神经退化类的病人。

2.根据权利要求1所述的一种基于多变量的DTI图像分析方法，其特征在于，所述正常人组和病人组的量的最优比例为1:1。