CN107392900B - 一种肺结节图像的多尺度增强方法 - Google Patents

Abstract

本发明公开了一种用多尺度增强滤波器和3D特征的肺结节检测方法。首先，构建结节和血管的理想模型，并构建增强滤波器分别对肺部3D图像中的结节图像和血管图像进行增强；接着提取疑似肺结节图像，并提取其SNOAH特征；最后采用不同核函数的SVM分类器对特性进行分类。实验结果表明，本发明中提出的肺结节特征描述子‑表面法线方向角度直方图是有效的，可以很好地区分结节图像和血管等其它的图像，能够有效地降低检测结果的假阳性。

Description

一种肺结节图像的多尺度增强方法

技术领域

本发明涉及肺部结节图像的检测，具体涉及一种肺结节图像的多尺度增强方法。

背景技术

在肺结节的计算机CAD系统中，由于CT图像中血管的密度、CT值等都与结节相似，往往会存在较高的假阳性，影响肺结节的检测准确率。在本发明中，提出了一种基于多尺度增强滤波器和3D形状特征的肺结节自动检测方法。

发明内容

本发明所要解决的技术问题是针对现有技术的不足提供一种用多尺度增强滤波器和3D特征的肺结节检测方法。

本发明的技术方案如下：

一种肺结节图像的多尺度增强方法，包括以下步骤：

A，构建结节和血管理想模型；

构建三维结节、血管和曲面的理想模型的表达式如下：

海森矩阵是一个由多元函数的二阶偏导数构成的方阵，描述了函数的局部曲率；在理想球体模型中，其每一个体素V(x,y,z)对应的Hessian矩阵H的表达式为：

通过计算发现，H是三阶对称矩阵，且其六个混合偏导f_xy，f_xz，f_yx，f_yz，f_zx，f_zy的值都为0，f_xx，f_yy，f_zz计算公式如下：

另外，对每个Hessian矩阵都按照以下公式进行矩阵分解，其中λ₁，λ₂，λ₃是矩阵的特征值，且|λ₁|≥|λ₂|≥|λ₃|，对应的特征向量分别为e₁，e₂和e₃；

由于上述公式中的Hessian矩阵为对角矩阵，所以对角线上的元素就是这个矩阵的特征值，计算公式如下：

Hessian矩阵的特征值对应的特征向量两两正交，而且特征向量的方向对应着三维椭球各轴的主方向，特征值的大小对应着各轴的长度，特征值和特征向量共同反映了物体的形状和大小；在球体内部，任何一点对应的海森矩阵的特征值都与该点的坐标值相关，且当该点离球心距离越近，其特征值大小的绝对值越大，当点为球心时，其3个特征值的绝对值都达到最大且相等，如下公式所示：

同样的，某点离圆柱中心线距离越近，其特征值大小的绝对值越大，当点在中心线上时，其2个特征值的绝对值都达到最大且相等，第3个特征值为0；

对于不同的理想模型，其特征值的对应关系为：球状结构，λ_b1≈λ_b2≈λ_b3＜0；管状结构，λ_l1≈λ_l2＜0,λ_l3≈0；面状结构λ_p1＜0,λ_p2≈λ_p3≈0，其中，λ_b1,λ_b2,λ_b3、λ_l1,λ_l2,λ_l3和λ_p1,λ_p2,λ_p3分别表示不同理想模型下的特征值。

B,构建多尺度结节图像和血管图像增强滤波器，采用基于高斯函数的多尺度滤波方法，利用高斯函数对图像进行卷积运算，去除图像中的噪声同时平滑图像，再对图像进行多尺度增强，增强不同尺度大小的结节图像；

步骤B中构建了结节图像和血管图像增强滤波器，其增强因子的表达式分别如下；

其中，λ_b1,λ_b2,λ_b3和λ_l1,λ_l2,λ_l3分别表示球状和管状理想模型下的特征值。

所述的肺结节图像的多尺度增强方法，候选结节的直径范围为3mm-30mm，尺度的个数N为5，当N确定后，每个尺度大小的计算公式如下：

实验结果表明，本发明方法能够有效地检测出肺结节，并降低假阳性，能达到的0.9681敏感度和0.9477特异度，肺结节的检测准确率为0.9501，平均每个病例的假阳性数目0.96。

附图说明

图1是本发明进行肺结节检测的方法流程图。

图2是本发明中理想球体模型(a)和圆柱体模型(b)。

图3是本发明特征值模型。(a)球形模型(b)圆柱模型(c)特征向量坐标系。

图4是本发明中结节增强函数Eball输出曲线。

图5是本发明中血管增强函数Eline输出曲线。

图6是不同尺度的肺结节图像增强结果对比图。

图7是本发明18-邻域3D区域增长法示意图。

图8是本发明中结节(a)和血管(b)梯度方向示意图。

图9是本发明中仰角θ和方位角φ示意图。

图10结节SNOAH特征。(a)结节(b)仰角(c)方位角。

图11血管SNOAH特征。(a)血管(b)仰角(c)方位角。

图12胸膜SNOAH特征。(a)胸膜(b)仰角(c)方位角。

图13胸膜牵拉SNOAH特征。(a)胸膜牵拉(b)仰角(c)方位角。

图14血管交叉SNOAH特征。(a)血管交叉(b)仰角(c)方位角；

具体实施方式

以下结合具体实施例，对本发明进行详细说明。

参照图1，本发明方法的实现流程如下：

(1)结节和血管模型

由于结节和血管在形态上表现出球状和管状特性，在三维图像中可以将它们分别近似看作一个球体和圆柱体。为此，首先构造三个理想模型，在三维空间中他们分别代表点(球体)，线(柱状)，面，表达式如下所示，其中理想球体模型和圆柱体模型示意图如图2所示。

海森矩阵是一个由多元函数的二阶偏导数构成的方阵，描述了函数的局部曲率。在理想球体模型中，其每一个体素V(x,y,z)对应的Hessian矩阵H的表达式为。

通过计算很容易发现，H是三阶对称矩阵，且其六个混合偏导f_xy，f_xz，f_yx，f_yz，f_zx，f_zy的值都为0，f_xx，f_yy，f_zz计算公式如下。

另外，对每个Hessian矩阵都可以按照以下公式进行矩阵分解，其中λ₁，λ₂，λ₃是矩阵的特征值，且|λ₁|>＝|λ₂|>＝|λ₃|，其对应的特征向量分别为e₁，e₂和e₃。

由于上述公式中的Hessian矩阵为对角矩阵，所以对角线上的元素就是这个矩阵的特征值，计算公式如下：

它的特征值对应的特征向量两两正交，而且特征向量的方向对应着三维椭球各轴的主方向，特征值的大小对应着各轴的长度，它们共同反映了物体的形状和大小。球形模型和圆柱模型的特征值模型如图3所示；由上述公式和图3可知，在球体内部，任何一点对应的海森矩阵的特征值都与该点的坐标值相关，且当该点离球心距离越近，其特征值大小的绝对值越大，当点为球心时，其3个特征值的绝对值都达到最大且相等，如下公式所示：

同样的，某点离圆柱中心线距离越近，其特征值大小的绝对值越大，当点在中心线上时，其2个特征值的绝对值都达到最大且相等，第3个特征值为0。

而且对于不同的理想模型，其特征值的对应关系总结如表1所示。

表1不同理想模型对应的海森矩阵特征值关系

(2)增强滤波器

由于结节和血管在形态上表现出球状和管状特性，根据不同理想模型对应的海森矩阵特征值关系，我们可以设计增强滤波器用来增强图像中的结节和血管。

根据三种不同的理想模型对应的海森矩阵特征值的特点，为了提高区分结节和血管的能力，我们充分考虑特征值之间的关系，将特征值之间的比值关系作为增强函数的输入，提出了新的结节增强函数E_ball和血管增强函数E_line。

首先，我们将增强函数的输入值和输出值都限定在[0,1]之间，我们定义如下函数。

对上述公式进行分析可知，函数的输入n/m和输出

都[0,1]之间，且单调递增，当m和n的值越接近时，函数的输出也就越大，当二者相等时，输出最大为1。为了充分考虑特征值两两之间的关系，我们定义新的结节增强函数为。

通过化简得到最终的结节增强函数E_ball如下。

由于在血管模型中，特征值λ₃的值趋近于0且值很小，同样的，为了充分考虑到特征值之间的关系，我们定义新的血管增强函数E_line如下：

通过化简得到最终的血管增强函数E_line为

以λ₂/λ₁和λ₃/λ₂为坐标输入，绘制了结节增强函数E_ball输出曲线，如图4所示。由曲线可知，λ₂/λ₁和λ₃/λ₂对E_ball具有相同的贡献值，而且E_ball的值域为[0,1]，且单调递增。体素点的E_ball值能够表达该体素点属于结节的可能性，E_ball值越大，则该体素点属于结节的可能性也越大。某体素点的E_line值能够表达该体素点属于血管的可能性，E_line值越小则该体素点属于血管的可能性也越小。

类似的，我们以λ₂/λ₁和λ₃为坐标输入，绘制了血管增强函数E_line输出曲线，如图5所示。由曲线可知，E_line值在[0,1]单调递增。因此，某体素点的E_line值能够表达该体素点属于血管的可能性，E_line值越大，则该体素点属于血管的可能性也越大；E_line值越小则该体素点属于血管的可能性也越小。因此可以通过E_line增强函数的大小来对血管图像进行增强。

(3)高斯函数多尺度计算

在肺部CT图像中，结节的尺寸大小具有不确定性，而且存在许多图像噪声，而体素的二阶偏导计算过程对图像噪声具有很强的敏感性，因此，如果将增强滤波器直接应用在图像上，将不会产生好的结果。为了有效检测出不同大小的结节，在本发明采用了一种基于高斯函数的多尺度滤波方法。利用高斯函数对图像进行卷积运算，去除图像中的噪声同时平滑图像，再对图像进行多尺度增强，可以增强不同尺度大小的结节图像。

对于标准差为σ的高斯函数，在范围为[-4σ,4σ]内包含了高斯函数95％以上的面积区域。因此在设计增强函数的尺度σ时，当需要增强的结节直径为d时，函数的标准差σ则应为d/4。另外，在结节增强之前，应先用标准差为σ的高斯函数对图像进行卷积运算，这样不仅可以消除图像噪声还可以提高滤波器的增强效果。

假设待增强的肺结节的直径范围为[d₀,d₁]，可将高斯滤波器的尺度范围设定在[d₀/4,d₁/4]之间。然后选取N个不同的σ值，分别对图像进行卷积和增强运算。N越大，结节的检测精度越高，计算越复杂；N越小，结节的检测精度较低，计算越简单。因此，N的选择对结节的检测结果精度影响很大。在本发明中，候选结节的直径范围为[3mm-30mm]，尺度的个数N为5，当N确定后，每个尺度大小的计算如下公式：

综上所述，使用多尺度滤波器进疑似结节检测的具体的步骤如下：

1)确定待增强物体大小的范围[d₀,d₁]和过滤器的尺度个数N，并利用下述公式计算每一个尺度σ_N的大小。

2)对每一个尺度σ_N，进行步骤(3)-(6)；

3)对图像进行高斯滤波平滑图像，并对每个体素进行步骤(4)-(6)；

4)构建Hessian矩阵并计算三个特征值λ₁，λ₂，λ₃，计算每个体素的结节增强函数E_ball和血管增强函数E_line。

5)分别使用本发明提出的结节增强滤波器和血管增强滤波器对体素进行增强；

6)结束对体素和尺度的循环。

7)对每个体素点，最终的输出是不同尺度计算得到的最大值，得到结节增强图像Nmax和血管增强图像Vmax。

8)将Nmax图像中减去与Vmax相交的图像，得到疑似肺结节的肺部3D图像；

9)输出最终的结节检测结果。

对每个体素点，最终的输出是不同尺度计算得到的最大值，得到结节增强图像和血管增强图像。增强的结果如图6所示。

(4)疑似肺结节提取

每个体素的最终输出结果是所有不同尺度下的增强函数的最大值，体素的结节增强函数值越大，其隶属于肺结节的可能性就越大，因此，候选结节的坐标位置可以通过分析体素点的结节增强函数值得到。另外，对每个体素点，我们采用一个5维的特征向量[E,g,x,y,z]^T来表示，其中E代表体素点的结节增强函数值E_ball，g代表体素点的灰度值，x和y代表点的坐标，z代表CT序列号。两个像素的特征向量之间的欧式距离用来衡量像素之间的相似性。

在本发明中，为了提取所有的候选结节，提出了一种基于多种子点的18-邻域3D区域增长法，用于提取每一个结节图像，利用种子点进行18-邻域区域增长的过程示意图如图7所示。

本发明中提取肺结节具体的步骤如下：

1)设置结节增强函数E_ball的阈值T_E，获取疑似结节的体素种子点集合{Seed₁,Seed₂,…,Seed_i}；

2)在每个种子点的18邻域内进行扩展；

3)对18领域内的每个体素点，计算其与当前种子点的欧氏距离D；

4)如果D<T_d，将该体素包含到种子区域内，并将该体素作为新的种子点；

5)重复步骤(2)-(4)直到种子区域大小不在发生变化为止；

6)输出检测的所有的疑似结节3D图像。

首先设置一个阈值T_E，当体素点的结节增强函数值E_ball大于阈值T_E时，选取这些体素点作为结节初始种子点。然后，对所有的种子点，计算其18邻域内的体素点与种子点特征向量的欧氏距离D，若D小于阈值T_d，则将其扩入种子区域中，并以该体素点作为新的种子点进行循环，直到种子区域的大小不再发生改变为止。

(5)特征提取

图8所示为结节和血管各像素的梯度方向，结节像素的梯度方向均匀分布，没有明显的主方向；血管像素的梯度方向基本一致。像素法线方向和梯度方向垂直，则结节的法线方向分布也呈均匀的放射状，没有主方向；血管的法线方向和血管中血液流向一致。因此，利用像素的法线方向分布可以区分结节和血管。

本发明提出的肺结节特征描述子-表面法线方向直方图SNOAH(Surface normalorientation angles histogram)能够反映肺结节表面法线方向角度的概率分布情况。

在本发明中，表面法线方向角度直方图(SNOAH)特征提取的步骤如下：首先，输入候选结节图像，通过Hessian矩阵分解求解特征值和特征向量，然后计算每个体素的表面法线方向求解其方向角度，最后进行统计分析得到其体素表面法线方向角度分布直方图。

在上述步骤(1)中，对Hessian矩阵进行了分解，对分解后的表达式可进行一步变形，我们可以得到下列公式。

H＝(λ₁-λ₂)e₁e₁ ^T+(λ₂-λ₃)(e₁e₁ ^T+e₂e₂ ^T)+λ₃(e₁e₁ ^T+e₂e₂ ^T+e₃e₃ ^T)

其中，(λ₁-λ₂)e₁e₁ ^T，(λ₂-λ₃)(e₁e₁ ^T+e₂e₂ ^T)，λ₃(e₁e₁ ^T+e₂e₂ ^T+e₃e₃ ^T)分别表示表面粗糙度surface-ness，曲度curve-ness和球点度point-ness，而且它们系数λ₁-λ₂、λ₂-λ₃、λ₃分别表示对应的权值。

surface-ness：方向是e₁，权值为λ₁-λ₂；

curve-ness：方向是e₃，权值为λ₂-λ₃；

point-ness：无方向，权值为λ₃；

从定义中我们发现，表面粗糙度的法向量方向与特征值λ₁对应的特征向量e₁的方向一致，因此，我们可以通过计算特征向量e₁的方向来确定平面法向量的方向。

对候选结节图像表面的每个体素，进行海森矩阵分解，计算其特征值和特征向量。在空间坐标系中，空间内任意一点都可以由半径r、仰角θ和方位角

的组合

来表示，如图9所示。由于surface-ness(λ₁-λ₂)e₁e₁ ^T的方向是法向量e₁的方向，其权值大小为λ₁-λ₂，所以仰角θ和方位角

可以根据以下公式进行计算得到。

其中，e₁ ^(x),e₁ ^(y),e₁ ^(z)分别表示表面法线方向e₁在x、y和z方向上的分量，仰角的范围为[0，180]度，方位角的范围为[0,360]度，然后通过对仰角θ和方位角

进行统计，得到二维统计直方图。

在进行统计表面法线方向俯仰角和方位角时，将俯仰角θ的范围平均分为m份，每份覆盖的角度范围为(180/m)度。同理，也将方位角

的范围平均分为m份，每份覆盖的角度范围为(360/m)度。也就是说，最终提取得到的特征向量的维度是2m。

分别对结节、血管和胸膜提取其SNOAH特征的结果如图10、图11和图12所示。从图中可以看出，结节的表面法向量的方位角没有明显的主方向，在0-360°之间均匀分布；血管的表面法向量的方位角集中在180°和360°上，在仰角上也有明显的主方向；而胸膜的表面法向量的方位角和仰角都仅仅分布在几个特殊角度上。因此，本发明提出的特征描述子-SNOAH能够进行物体形状的区分。此外，我们还对几种比较特殊的图像进行了角度直方图的统计。胸膜牵拉图像如图13所示，血管相交图像如图14所示。

(6)分类

SVM考虑寻找一个满足分类要求的超平面，并且使训练集中的点距离分类面尽可能的远，也就是寻找一个分类面使它两侧的空白区域(margin)最大。

在本发明中，首先构建N个疑似肺结节图像的样本集X＝{(x_i,y_i)}，其中特征向量x_i∈Rⁿ，i＝1,2，…N,y_i∈{1,-1}。SVM训练样本的目的就是要最小化误差函数，也就是解决下面公式的优化问题。

Subject to

其中f(x_i)是一个决策函数，κ是一个函数映射，能够将特征向量x_i映射到高维特征空间，SVM的目的是寻找一个满足分类要求的超平面(w，b)，并且使样本点距离超平面尽可能的远。C大于0，是一个惩罚因子。

对于解决线性不可分问题，SVM通过引入核映射，把样本空间映射到一个高维空间中，这样就使得原来样本空间中线性不可分的问题变得高维样本空间中线性可分，从而很好地解决了非线性分类问题。其中核函数恒定义为K(x_i,x)≡κ(x_i)^Tκ(x)，在本发明中，用到了我们使用RBF核K_r(x_i,x)，多项式核K_p(x_i,x)和混合核K_mix(x_i,x)分别作为核函数，它们的公式分别如下。

K_p(x_i,x)＝||x_i·x+1||^p

K_mix(x_i,x)＝βK_p(x_i,x)+(1-β)K_r(x_i,x)

其中β为混合权重因子。

通过实验表明，采用k-折交叉验证分类结果时，当k＝6时，各个分类器的最佳结果对比如下表2所示。

表2k＝6时的结节检测结果对比

由表2可知，基于混合核函数的SVM分类器SVM_M具有最佳的分类结果，能达到96.81％的敏感度，假阳性率为0.96，效果较好。就分类敏感度而言，SVM_P分类器和SVM_R分类器的分类准确率也比较好，能达到93.62％和91.49％，这也从另外一个角度也可以说明我们提取的特征是比较有效的。因此，整体而言，本发明的方法能够有效的降低检测结果的假阳性，同时也能保持较高的检测敏感度和准确率。

应当理解的是，对本领域普通技术人员来说，可以根据上述说明加以改进或变换，而所有这些改进和变换都应属于本发明所附权利要求的保护范围。

Claims (2)

1.一种肺结节图像的多尺度增强方法，其特征在于，包括以下步骤：

A，构建结节和血管理想模型；

构建三维结节、血管和曲面的理想模型的表达式如下：

海森矩阵是一个由多元函数的二阶偏导数构成的方阵，描述了函数的局部曲率；在理想球体模型中，其每一个体素V(x,y,z)对应的Hessian矩阵H的表达式为：

通过计算发现，H是三阶对称矩阵，且其六个混合偏导f_xy，f_xz，f_yx，f_yz，f_zx，f_zy的值都为0，f_xx，f_yy，f_zz计算公式如下：

另外，对每个Hessian矩阵都按照以下公式进行矩阵分解，其中λ₁，λ₂，λ₃是矩阵的特征值，且|λ₁|≥|λ₂|≥|λ₃|，对应的特征向量分别为e₁，e₂和e₃；

由于上述公式中的Hessian矩阵为对角矩阵，所以对角线上的元素就是这个矩阵的特征值，计算公式如下：

Hessian矩阵的特征值对应的特征向量两两正交，而且特征向量的方向对应着三维椭球各轴的主方向，特征值的大小对应着各轴的长度，特征值和特征向量共同反映了物体的形状和大小；在球体内部，任何一点对应的海森矩阵的特征值都与该点的坐标值相关，且当该点离球心距离越近，其特征值大小的绝对值越大，当点为球心时，其3个特征值的绝对值都达到最大且相等，如下公式所示：

同样的，某点离圆柱中心线距离越近，其特征值大小的绝对值越大，当点在中心线上时，其2个特征值的绝对值都达到最大且相等，第3个特征值为0；

对于不同的理想模型，其特征值的对应关系为：球状结构，λ_b1≈λ_b2≈λ_b3＜0；管状结构，λ_l1≈λ_l2＜0,λ_l3≈0；面状结构λ_p1＜0,λ_p2≈λ_p3≈0，其中，λ_b1,λ_b2,λ_b3、λ_l1,λ_l2,λ_l3和λ_p1,λ_p2,λ_p3分别表示不同理想模型下的特征值；

B,构建多尺度结节图像和血管图像增强滤波器，采用基于高斯函数的多尺度滤波方法，利用高斯函数对图像进行卷积运算，去除图像中的噪声同时平滑图像，再对图像进行多尺度增强，增强不同尺度大小的结节图像；

步骤B中构建了结节图像和血管图像增强滤波器，其增强因子的表达式分别如下；

其中，λ_b1,λ_b2,λ_b3和λ_l1,λ_l2,λ_l3分别表示球状和管状理想模型下的特征值。

2.根据权利要求1所述的肺结节图像的多尺度增强方法，其特征在于，候选结节的直径范围为3mm-30mm，尺度的个数N为5，当N确定后，每个尺度大小的计算公式如下：

Images