CN104881680A - 一种基于二维特征和三维特征融合的阿尔茨海默病及轻度认知功能障碍识别方法 - Google Patents

Abstract

本发明提出了一种基于三维形态学特征和二维纹理特征的阿尔茨海默病及轻度认知功能障碍计算机辅助识别方法，具体包括：1)对医学图像进行预处理，包括预分割和配准等过程；2)对医学图像进行二维纹理特征提取，包括灰度共生矩阵的二次统计量和Gabor小波变换的多尺度多方向特征值；3)对医学图像进行三维形态学特征提取，即提取感兴趣区域的体积特征；4)对三维形态学特征和二维纹理特征进行特征融合；5)构建支持向量机分类器实现对阿尔茨海默病和轻度认知功能障碍的识别。本发明将三维形态学特征和二维纹理特征结合，能够更加全面准确地表达医学图像的内容，所以该方法能够提高对阿尔茨海默病和轻度认知功能障碍的识别，从而提供更有效的临床辅助诊断。

Description

一种基于二维特征和三维特征融合的阿尔茨海默病及轻度认知功能障碍识别方法

技术领域

本发明属于计算机辅助诊断领域，更为具体地讲，涉及一种基于脑部医学图像的阿尔茨海默病及轻度认知功能障碍计算机辅助识别的方法。

背景技术

阿尔茨海默病(Alzheimer’s Disease,AD)是一种致死性神经退行性疾病，严重损害了老年人的身体健康，我国八十岁以上老人发病率在30％，轻度认知功能障碍(Mild Cognitive Impairment,MCI)则被认为是介于正常老化与AD之间的一种临床状态。而随着影像医学和计算机领域相关技术的发展，如今诊断疾病很大程度上依赖于医学图像的识别来辅助医生进行诊断，所以通过实现基于医学图像的计算机辅助识别分类来辅助医生诊断AD和MCI具有重大的作用和意义。其中，图像识别的前期工作和关键环节就是图像特征提取与选择。

传统的特征提取则主要分为底层特征提取和语义特征提取。底层特征主要有灰度，纹理和形状三大类特征。语义特征则属于高层特征，该类特征需要研究者基于对图像内容的先验知识以及实际应用需求从语义的角度对其信息进行描述刻画。而与普通图像相比，医学图像本身具有高精度，信息量巨大，且与临床密切相关等专业性极强的特点。所以近年来，国内外学者对医学图像的特征提取展开了大量研究。

经过对现有文献分析可知，纹理特征非常适用于病变的纵向检测，目前被广泛应用于医学图像领域，G.Castellano等人在文献(Castellano G,Bonilha L,Li L M,et al.Texture analysis of medical images[J].Clinical radiology,2004,59(12):1061-1069.)中阐明了医学图像中常用的纹理分析技术，可分为四种形式：结构型，基于模型，统计分析和基于变换，这四种形式各有优缺点。但是，随着技术的发展采用单一特征已经难以准确地表达医学图像的特征，从而不利于后期的识别分类。

本发明方案基于上述不足，提出了一种基于二维特征和三维特征融合的阿尔茨海默病及轻度认知功能障碍计算机辅助识别方法，不但保留了纹理特征用于疾病识别的优势，而且利用三维特征弥补了二维特征的不足，能有效提高疾病的识别率。

发明内容

本发明的目的在于克服现有辅助识别阿尔茨海默病和轻度认知功能障碍方法的不足，提出了一种基于二维纹理特征和三维形态学特征相结合的识别方法，突破了传统方法中单个维度或者单一特征类型在反映医学图像特征上的局限性，有助于全面、准确地表达医学图像特性，有效地提高疾病识别的准确性，从而更有效的辅助医生的临床诊断。

为了实现上述目的，本发明提出了一种基于二维纹理特征和三维形态学特征相结合的识别方法，包括：

对多个样本个体对应的结构磁共振图像进行预处理，采用基于体素的脑形态学方法(Voxel-based morphometry，VBM)和通过取幂李代数的微分同胚解剖配准(diffeomorphic anatomical registration through exponentiated lie algebra，DARTEL)方法，包括将原始图像配准到标准空间，去除非脑组织和预分割为脑白质、脑白质以及脑脊液三部分等过程；

选取每个样本的二维图像横断面序列图像，针对每幅图像，构建不同方向和尺度的多个Gabor滤波器，然后提取系数幅度序列的均值及其标准方差作为特征值，以及构建多个灰度共生矩阵，提取能量，对比度，相关度，以及熵共4维特征量；

利用二维图像序列及其二维序列的位置信息，按其空间位置关系进行曲面重建，使用双样本T检验对其进行统计分析，比较样本之间的脑区体积和密度的显著差异区；

将显著差异区作为感兴趣区域(region of interest，ROI)，提取每个样本的所有感兴趣区域的体积参数作为三维形态学特征；

以上提取的特征信息均采用维特征向量表示，然后利用主成分分析法对上述提取的所有特征向量进行特征级融合，作为每个样本的多维特征向量；

将样本分为训练集和测试集，选择径向基核函数的支持向量机(SupportVector Machine，SVM)，用训练集的特征向量对SVM分类器进行训练，训练完成后，将测试集的特征向量输入SVM分类器，最后评价测试结果。

本发明的特点在于基于三维形态学特征和二维纹理特征融合的阿尔茨海默病及轻度认知功能障碍识别方法，一方面，灰度共生矩阵思想发展成熟且易于实现，在图像的局部模式上有较好的鉴别能力，而Gabor小波能较好地描述人类的视觉神经元的感受野，具有时域和频域的联合最佳分辨率，另一方面，提取ROI的形态学特征有助于表现医学图像特有的生理结构信息，将上述特征相融合可以得到一种具有更强鉴别力的特征向量，能够有效地提高阿尔茨海默病及轻度认知功能障碍的识别率。

附图说明

图1是本发明提出的基于三维形态学特征和二维纹理特征的阿尔茨海默病及轻度认知功能障碍识别方法的总体流程图；

图2是基于VBM-DARTEL的预处理流程图；

图3是SVM的示意图。

具体实施方式

下面结合附图对本发明的具体实施方式进行描述，以便本领域的技术人员更好地理解本发明。需要注意的是，以下描述仅仅用于解释本发明，并不用于限定本发明。

本发明提出的基于三维形态学特征和二维纹理特征的阿尔茨海默病及轻度认知功能障碍识别方法的总体流程图如图1所示，具体如下：

步骤(1)是对结构磁共振图像进行预处理，包括患者图像和正常人图像。采用了基于体素的脑形态学方法和通过取幂李代数的微分同胚解剖配准方法。VBM方法是一种以体素为单位的形态测量学方法，通过计算单位体积内包含某种体素的密度变化来显示脑组织的形态学改变，DARTEL算法则是一种VBM分析的算法，通过对流场取幂的方法获得相应的变形场，从而保证映射的结果是微分同胚的。

基于VBM-DARTEL的预处理流程如图2所示：

步骤(1.1)首先是在原始空间上根据大脑的先验概率分布模板对图像进行预分割，包括去除非脑组织即颅骨剥离，以及分割出灰质，白质和脑脊液，这一过程保证了后期统计分析的数据不会受到白质和灰质的相互影响。

步骤(1.2)利用DARTEL方法将正常人的灰质图像迭代配准到对应的平均图像上，从而得到多个模板图像，选取最为精确的模板图像作为标准模板，再通过DARTEL配准方法把所有灰质图像配准到该标准模板上，得到不同图像对应的变形场。

步骤(1.3)基于步骤(1.2)得到的变形场将所有灰质图像配准到MNI空间进行标准化，然后对图像进行调制从而保证体素的总数量不变，该过程中包括仿射变换和平滑处理，平滑处理一般采用高斯平滑核。

步骤(2)是对图像进行特征提取，首先选取每个样本的二维图像横断面序列图像，然后分别对每幅图像提取二维纹理特征和三维形态学特征，具体如下：

步骤(2.1)对图像进行二维纹理特征提取，具体如下：

步骤(2.1.1)采用灰度共生矩阵统计分析图像中灰度的空间相关特性，提取其二次统计量作为纹理特征。灰度共生矩阵的定义为从灰度级i的点离开某个固定位置关系d＝(Dx,Dy)达到灰度为j的概率，用P_d(i,j)(i,j＝0,1,2,…，L-1)表示，其中L是图像的灰度级，i和j表示像素的灰度，d表示两个像素间的空间位置关系，包括了距离和方向，一般方向取0°,45°,90°,135°四个方向。通过图像构建灰度共生矩阵之后将各个元素P_d(i,j)除以各元素之和S得到归一化值由此得到归一化共生矩阵。考虑到便于计算以及能够提供较高的分类精度，一般提取如下四个特征值作为纹理特征：

1)能量：

$f_1=\underset{i=0}{\overset{L-1}{\Σ}}\underset{j=0}{\overset{L-1}{\Σ}}{P}_d^2(i,j)$

能量反映了图像灰度分布均匀程度和纹理粗细程度。

2)对比度：

$f_2=\underset{n=0}{\overset{L-1}{\Σ}}{n}^2\left\{\underset{i=0}{\overset{L-1}{\Σ}}\underset{j=0}{\overset{L-1}{\Σ}}{\hat{P}}_d(i,j)\right\}$

对比度反映了图像的清晰度和纹理沟纹深浅的程度。

3)相关：

$f_3=\frac{\underset{i=0}{\overset{L-1}{\Σ}}\underset{j=0}{\overset{L-1}{\Σ}}ij{\hat{P}}_d(i,j)-{\mathrm{\μ}}_1{\mathrm{\μ}}_2}{{\mathrm{\σ}}_1^2{\mathrm{\σ}}_2^2}$

其中

${\mathrm{\μ}}_1=\underset{i=0}{\overset{L-1}{\Σ}}i\underset{j=0}{\overset{L-1}{\Σ}}{\hat{P}}_d(i,j)$

${\mathrm{\μ}}_2=\underset{i=0}{\overset{L-1}{\Σ}}j\underset{j=0}{\overset{L-1}{\Σ}}{\hat{P}}_d(i,j)$

${\mathrm{\σ}}_1^2=\underset{i=0}{\overset{L-1}{\Σ}}{(i-{\mathrm{\μ}}_1)}^2\underset{j=0}{\overset{L-1}{\Σ}}{\hat{P}}_d(i,j)$

${\mathrm{\σ}}_2^2=\underset{i=0}{\overset{L-1}{\Σ}}{(i-{\mathrm{\μ}}_2)}^2\underset{j=0}{\overset{L-1}{\Σ}}{\hat{P}}_d(i,j)$

相关用于衡量灰度共生矩阵的元素在行或列方向上的相似程度。

4)熵：

$f_4=\underset{i=0}{\overset{L-1}{\Σ}}\underset{j=0}{\overset{L-1}{\Σ}}{P}_d^2(i,j){\mathrm{logP}}_d^2(i,j)$

熵反映图像中纹理的复杂程度或非均匀度。

步骤(2.1.2)是基于Gabor多通道滤波技术，提取多尺度多方向的纹理特征。Gabor函数是目前唯一能够达到时频测不准关系下届的函数，能最好地兼顾信号在时域和频域中的分辨能力。二维Gabor函数及其傅里叶变换如公式所示：

$g(x,y)=\left(\frac{1}{2{\mathrm{\π\σ}}_x{\mathrm{\σ}}_y}\right)\exp \{-\frac{1}{2}\[\frac{x^2}{{\mathrm{\σ}}_x^2}+\frac{y^2}{{\mathrm{\σ}}_y^2}\]+2\mathrm{\π}jwx\}$

但是由于Gabor函数本身不具有小波函数的正交特性，所以需要将g(x,y)作为母小波对其进行伸缩和变换处理，得到Gabor小波如公式所示：

g_mn(x,y)＝a^-mg(x,y)

其中

x＝a^-m(xcosθ+ysinθ)

y＝a^-m(ycosθ-xsinθ)

a^-m即是尺度因子，m和n分别是指定小波的尺度和方向，即是说，m＝0,1,2…,s-1，n＝0,1,2…,t-1,s和t分别是总的尺度数和总的方向数。即是说，针对W×H的医学图像H(x,y)的离散Gabor小波变换是：

$G_{mn}(x,y)=\underset{i}{\Σ}\underset{j}{\Σ}H(x-i,y-j)g_{mn}^{*}(i,j)$

其中g*_mn(i,j)是gmn(i.j)的共轭复数，i和j即是窗口大小(即滤波器模板大小)的变量，此时设计的滤波器需要确保Gabor小波滤波器族的响应在频谱半高宽能相互接触但不重叠。然后由于系数幅度值能够表征不同方向和尺度上的图像能量，所以，计算出不同方向和尺度上变换后所得的系数幅度序列为：

$E(m,n)=\underset{x}{\Σ}\underset{y}{\Σ}\left|G_{mn}(x,y)\right|$

因此，纹理特征分别是系数幅度序列的均值μmn及其标准方差σmn：

μ_mn＝E(m,n)/WH

${\mathrm{\σ}}_{mn}=\sqrt{\underset{x}{\Σ}\underset{y}{\Σ}{\left(\right|G_{mn}(x,y)|-{\mathrm{\μ}}_{mn})}^2}/WH$

步骤(2.2)对图像进行三维形态学特征提取，具体如下：

步骤(2.2.1)是选取感兴趣区域。利用二维图像序列及其二维序列的位置信息，按其空间位置关系进行曲面重建，然后将不同组图像中的每个体素作为参数，进行双样本T检验，找出具有显著差异的区域，则是具有统计学意义的簇，将这些簇作为感兴趣区域。其中，为了使簇具有显著差异，在统计过程中可设置p≤0.05，并且进行FDR校正。

步骤(2.2.2)是提取每个样本的感兴趣区域的体积特征。将在步骤(2.2.1)中得到的感兴趣区域制作成多个mask，然后用每个mask点乘每个样本的灰质图像，从而计算出每个样本的感兴趣区域的体积。所以，每个样本的感兴趣区域的体积就是一组特征向量。

步骤(3)是对步骤(2)提取所得的特征进行选择融合，以便消除特征间的冗余信息，且在保留了所需识别信息的基础上降低特征维数，减少计算代价。主成分分析法是一种把多个实测变量转换为少数几个不相关的综合指标(即主成分)的多元统计分析方法，其基本思想是基于实测变量之间必定存在一定相关关系，因此，可以用综合各变量信息的互不重叠的综合指标来代表这些实测变量。

每幅图像通过步骤(2)提取了n维特征，一共有k幅图像，那么这个原始特征集就是：

$X=\left[\begin{array}{cccc}{x}_{11}& x_{12}& ...& x_{1k}\\ x_{21}& x_{22}& ...& x_{2k}\\ \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& & \begin{array}{c}.\\ .\\ .\end{array}\\ x_{n1}& x_{n1}& ...& x_{nk}\end{array}\right]=\[x_1{x}_2{x}_3...x_k\]$

对所有特征值采用标准差变换和极差转换进行标准化处理之后，使用x₁，x₂，…x_k线性表示综合变量y₁，y₂，…y_k：

$\left\{\begin{array}{c}{y}_1=u_{11}{x}_1+u_{12}{x}_2+u_{13}{x}_3+...+u_{1k}{x}_k\\ y_2=u_{21}{x}_1+u_{22}{x}_2+u_{23}{x}_3+...+u_{2k}{x}_k\\ ...\\ y_k=u_{k1}{x}_1+u_{k2}{x}_2+u_{k3}{x}_3+...+u_{kk}{x}_k\end{array}\right.$

其中，需满足如下条件：(2)y_i和y_j互不相关；(3)y₁是x₁，x₂，…x_k所有线性组合中的方差最大值，y₂则是与前者不相关的x₁，x₂，…x_k所有线性组合中的方差最大值，以此类推，y_k是与之前所有综合变量不相关的x₁，x₂，…x_k所有线性组合中的最大方差值。

综上所述，得到k个新的综合变量(即主成分)，然后根据矩阵代数，把特性向量u₁，u₂，…u_k所对应的非零特征根设为λ₁，λ₂，…λ_k，则可以定义前p个主成分y₁，y₂，…y_k的累计贡献率为：

$\underset{i=0}{\overset{p}{\Σ}}{\mathrm{\λ}}_i{\left(\underset{i=0}{\overset{k}{\Σ}}{\mathrm{\λ}}_i\right)}^{-1}$

当p个主成分的累计贡献率超过85％时，即选取这p个主成分作为特征融合的结果。

步骤(4)基于SVM方法对样本图像构建分类器，实现对阿尔茨海默病，轻度认知功能障碍和正常人之间的自动识别。SVM的基本思想是将输入数据通过内积函数费线性映射到一个高维空间，然后在这个空间建立一个最大间隔超平面同时求解出最优的线性分类面，其示意图如图3所示，每个中间节点对应一个支持向量，输出是中间节点的线性组合。

构建SVM分类器时，将步骤(3)融合所得的特征向量作为输入，采用径向基核函数以及leave-one-out交叉验证法设置参数。由于是三分类问题，所以采用一对一的投票策略，即总共构建3个二分类器，训练集分别对应于(AD，正常人)，(MCI，正常人)，(AD，MCI)，测试时将测试样本依次送入这三个分类器，最后得到一组结果，取分类结果投票数最多的类别作为最终结果，若出现同样多的投票，一般取第一个最大值对应的类别。最后以平均预测正确率作为评价结果标准。

以上对本发明所提出的计算机辅助识别阿尔茨海默病以及轻度认知功能障碍方法进行了详细的介绍，但是具体实施方式的描述仅用于解释本发明的方法及其核心思想，以便于本技术领的技术人员理解本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。

Claims (3)

1.一种基于二维纹理特征和三维形态学特征结合的阿尔茨海默病和轻度认知功能障碍的识别方法，包括以下步骤：

步骤一：对多个样本个体对应的结构磁共振图像进行预处理，采用基于体素的脑形态学方法和通过取幂李代数的微分同胚解剖配准方法，包括配准到标准空间，去除非脑组织和预分割为脑白质、脑白质以及脑脊液三部分等过程；

步骤二：选取每个样本的二维图像横断面序列图像，针对每幅图像，构建不同方向和尺度的多个Gabor滤波器，然后提取系数幅度序列的均值及其标准方差作为特征值，以及构建多个灰度共生矩阵，提取能量，对比度，相关度，以及熵作为特征值；

利用二维图像序列及其二维序列的位置信息，按其空间位置关系进行曲面重建，使用双样本T检验对其进行统计分析，将样本之间的脑区体积和密度的显著差异区作为感兴趣区域，提取每个样本的所有感兴趣区域的体积参数作为三维形态学特征；

步骤三：以上提取的特征信息均采用特征向量表示，然后利用主成分分析法对上述提取的所有特征向量进行特征级融合，作为每个样本的多维特征向量；

步骤四：将样本分为训练集和测试集，选择径向基核函数的支持向量机(SVM)，用训练集的特征向量对SVM分类器进行训练，训练完成后，将测试集的特征向量输入SVM分类器，最后评价测试结果。

2.根据权利要求1所述的基于二维纹理特征和三维形态学特征结合的阿尔茨海默病和轻度认知功能障碍识别方法，其特征在于，所述步骤二和三中将二维纹理特征和三维形态学特征相结合，具体包括如下步骤：

步骤(1.1.1)选取每个样本的二维图像横断面序列图像，为每幅图像构建多个灰度共生矩阵，提取如下四个特征值作为二维纹理特征:

能量：

$f_1=\underset{i=0}{\overset{L-1}{\mathrm{\Σ}}}\underset{j=0}{\overset{L-1}{\mathrm{\Σ}}}{P}_d^2(i,j)$

对比度：

$f_2=\underset{n=0}{\overset{L-1}{\mathrm{\Σ}}}{n}^2\left\{\underset{i=0}{\overset{L-1}{\mathrm{\Σ}\stackrel{}{}}}\underset{j=0}{\overset{L-1}{\mathrm{\Σ}}}{\hat{P}}_d(i,j)\right\}$

相关：

$f_3=\frac{\underset{i=0}{\overset{L-1}{\mathrm{\Σ}}}\underset{j=0}{\overset{L-1}{\mathrm{\Σ}}}\mathrm{ij}{\hat{P}}_d(i,j)-{\mathrm{\μ}}_1{\mathrm{\μ}}_2}{{\mathrm{\σ}}_1^2{\mathrm{\σ}}_2^2}$

熵：

$f_4=\underset{i=0}{\overset{L-1}{\mathrm{\Σ}}}\underset{j=0}{\overset{L-1}{\mathrm{\Σ}}}{P}_d^2(i,j)\log P_d^2(i,j)$

步骤(1.1.2)为每幅图像构建不同方向和尺度的多个Gabor滤波器，计算不同尺度和方向上的系数幅度序列：

$E(m,n)=\underset{x}{\mathrm{\Σ}}\underset{y}{\mathrm{\Σ}}\left|G_{\mathrm{mn}}(x,y)\right|$

然后计算如下两个特征值作为二维纹理特征：

系数幅度序列的均值：

m_mn＝E(m,n)/WH

系数幅度序列的标准差：

${\mathrm{\σ}}_{\mathrm{mn}}=\sqrt{\underset{x}{\mathrm{\Σ}}\underset{y}{\mathrm{\Σ}}{\left(\right|G_{\mathrm{mn}}(x,y)|-{\mathrm{\μ}}_{\mathrm{mn}})}^2}/\mathrm{WH}$

步骤(1.2.1)利用二维图像序列及其二维序列的位置信息，按其空间位置关系进行曲面重建，然后将不同组的图像中的每个体素作为参数，进行双样本T检验，找出具有显著差异的区域作为感兴趣区域，计算出每个样本的感兴趣区域的体积作为三维形态学特征；

步骤(2)对上述得到的特征进行标准化处理，然后采用主成分分析法对其进行融合，选取其中m个最具判别性的综合指标，即主成分作为每个样本的多维特征向量。

3.如权利要求1所述的一种基于二维纹理特征和三维形态学特征结合的阿尔茨海默病和轻度认知功能障碍识别方法，其特征在于，所述步骤四中的基于SVM的阿尔茨海默病及轻度认知功能障碍识别的具体步骤如下：

构建三个支持向量机分类器，采用径向基核函数，将经过特征融合得到的特征向量作为输入，采用一对一投票策略，对阿尔茨海默病，轻度认知功能障碍和正常人进行分类识别，最后以平均预测正确率作为评价结果标准。