CN103824295A - 一种肺部ct图像中粘连血管型肺结节的分割方法 - Google Patents

Abstract

本发明涉及一种肺部CT图像中粘连血管型肺结节的分割方法，包括：（Ⅰ）提取肺实质；（Ⅱ）预分割肺结节：通过对肺实质轮廓自动分割获得其中的肺结节的初始轮廓，再通过内核图割法获取初始肺结节内轮廓；（Ⅲ）提取肺结节标记点：采用射线法在肺结节图像上提取图割的目标点和背景点；（Ⅳ）图割法分割肺结节：采用基于最大流算法的图割法对肺结节作进一步的分割，获取最终的肺结节轮廓。采用本发明的方法可以将肺结节实质和模糊外围进行区分，获取精确的肺结节边界和特征，有利于医生准确确定肺结节实质。

Description

一种肺部CT图像中粘连血管型肺结节的分割方法

技术领域

本发明属于图像分割技术领域，具体涉及一种肺部CT图像中粘连血管型肺结节的分割方法。

背景技术

对肺结节的精确分割是基于CT图像的早期肺癌计算机辅助诊断的关键步骤，能否从CT图像中精确地分割出肺结节，最终会影响到技术机辅助诊断系统的性能。

在CT图像中，由于结节等病变部位的密度与心脏及胸廓等软组织的密度极为接近，故传统的基于阈值的分割方法，虽然能有效地提取肺实质，却常常会丢失如结节等病变部位的细节信息。现有的成像技术还无法确切地区分出结节和胸廓这样密度相似的软组织，要解决这个问题只能依赖图像分割算法。针对这一问题，现有的算法多是先利用阈值分割、区域增长、边缘检测等方法得到CT图像中肺实质的大致区域，然后对分割结果进行修正，以弥补因各种病变造成的肺实质的缺失部分。然而，这种分割不够精确，不能满足临床诊断的需要。

由Boykov等人提出的Graph cuts（图割法）是目前应用最为广泛的图像分割方法之一，其基于图的最大流算法实现了图像的快速分割，尤其它采用的是用户提取目标和背景点，通过最小化Gibbs能量函数E来对图像中未标记点进行标记，获得比较精准的图像分割结果。这种分割利用图像像素间的距离和灰度相似性对像素进行分割，但是单纯地考虑两个像素之间的欧式距离对有些像素间有很大灰度差的像素不准确；此外，由于图割方法需要确定标记点，而肺结节图像比较特殊，在肺实质区域医生首先要确定感兴趣区域（ROI）的肺结节，然后医生手动标记点，这种获取标记点的方式对Graph cuts的分割结果有很大影响。

发明内容

针对现有技术中存在的缺陷，本发明的目的是提供一种肺部CT图像中粘连血管型肺结节的分割方法，采用该分割方法可以将肺结节实质和模糊外围进行区分，获取精确的肺结节边界和特征，有利于医生准确确定肺结节实质。

为达到以上目的，本发明采用的技术方案是：一种肺部CT图像中粘连血管型肺结节的分割方法，包括以下步骤：

（Ⅰ）提取肺实质：在输入的肺部图像上提取包含肺结节在内的肺实质轮廓，得到肺实质边缘曲线；

（Ⅱ）预分割肺结节：对步骤（Ⅰ）获得的肺实质进行自动分割，获得其中的肺结节的初始轮廓，再采用基于混合高斯模型的内核图割法对得到的肺结节进行初始分割，获取初始肺结节内轮廓；

（Ⅲ）提取肺结节标记点：在步骤（Ⅱ）得到的肺结节内部取一中心点，通过该中心点向初始肺结节内轮廓发射一组等角度的射线，取射线与初始肺结节内轮廓相交的第一点作为初始标记点，并基于该组初始标记点的平均距离进行初始标记点的插值操作，然后沿每一射线上的初始标记点分别向中心点提取目标点，向内轮廓外围提取背景点；

（Ⅳ）图割法分割肺结节：利用步骤（Ⅲ）提取的目标点和背景点，采用基于最大流算法的图割法对肺结节作进一步的分割，获取最终的肺结节轮廓。

进一步，步骤（Ⅰ）中，采用区域生长法并结合边缘检测法提取肺实质。

进一步，步骤（Ⅰ）中得到肺实质边缘曲线后，根据曲率和曲线延伸的惯性规律修正肺实质边缘曲线。

再进一步，根据曲率和曲线延伸的惯性规律修正肺实质边缘曲线，按照以下方式进行：

（1）根据式①求得肺结节边界曲线和与其相交的肺实质边界曲线的曲率K，并提取曲率方向改变的点；

$K=\left|\frac{I_y^2{I}_{\mathrm{xx}}-2I_x{I}_y{I}_{\mathrm{xy}}+I_x^2{I}_{\mathrm{yy}}}{\sqrt{I_x^2+I_y^2}}\right|$ ①；

（2）根据式②定义的曲率半径ρ判断曲线各个部分的凹凸性；

$\mathrm{\ρ}=-M\frac{{[{\left(I_x\right)}^2+{\left(I_y\right)}^2]}^{\frac{3}{2}}}{(I_x{I}_{\mathrm{yy}}-I_y{I}_{\mathrm{xx}})}$ ②；

式①、②中,I_x，I_y，I_xx，I_yy分别为曲线的坐标x,y对参变量的一阶、二阶导数，M是经验系数；

（3）求得肺实质边界曲线和肺结节边界曲线相交时的拐点，将曲率小的曲线段去掉；

（4）利用曲率方向改变的点，结合曲线延伸的惯性规律，曲线拟合法修正曲线边缘，得到包含肺结节在内的肺实质轮廓。

进一步，步骤（Ⅱ）中采用边缘检测法自动分割肺实质以获取肺结节的初始轮廓。

进一步，步骤（Ⅲ）中，目标点和背景点到初始标记点的距离分别为中心点到初始标记点距离的1/3。

进一步，步骤（Ⅳ）的实现过程如下：

（1）以步骤（Ⅲ）中标记的一组背景点所构成的曲线作为当前轮廓线；

（2）构造s-t网络：使当前轮廓线对应s-t网络的源点s，使步骤（Ⅲ）中标记的一组目标点所形成的内边界对应s-t网络的汇点t，其余每一个像素对应s-t网络中的一个结点，相邻结点间以8邻域方式带权值连接；

（3）最大流最小代价切割：计算式③中Gibbs能量函数E的区域项和边界项，用最大流算法对s-t网络进行最小代价切割，得到一条新当前轮廓线；

$E\left(L\right)=\mathrm{\λ}\underset{i\∈V}{\mathrm{\Σ}}R(i,L_i)+\underset{\{i,j\}\∈N}{\mathrm{\Σ}}B(i,j)\·\mathrm{\δ}(L_i\≠L_j)$

③；

$\mathrm{\δ}(L_i\≠L_j)=\left\{\begin{array}{cc}1;& L_i\≠L_j\\ 0;& L_i=L_j\end{array}\right.$

式③中：V为图像中像素点集合；N表示图像中所有邻接像素对{i,j}；L_i是一个二进制标记，L_i=1时，像素点i标记为前景标记，L_i=0时，像素点i标记为背景标记；而λ表示区域项R(i,L_i)和边界项B(i,j)之间的权重；

区域项R(i,L_i)定义如式④所示：

R(i,L_i)=-lnP(I_i|L_i)             ④；

式④中：I_i表示像素点i的灰度值；P(I_i|L_i)表示灰度值I_i被标记为L_i的概率；

边界项B(i,j)定义如式⑤所示：

$B(i,j)=\exp (-\frac{{(I_i-I_j)}^2}{2{\mathrm{\σ}}^2})\·\frac{1}{D_g(i,j)}\·|L_i-L_j|$ ⑤；

式⑤中：D_g(i,j)表示相邻像素点i和j间的测地距离；

（4）使用步骤（3）的新当前轮廓线代替当前轮廓线，重复步骤（2）-（3），使新当前轮廓线迭代收缩，直至达到图像中确切的肺结节轮廓线；

（5）输出图像中被分割的肺结节图像。

再进一步，测地距离D_g(i,j)定义如式⑥所示：

$D_g(i,j)=\underset{P_{i.j}}{\min }{\∫}_0^{1}D\left(P_{i,j}\left(t\right)\right)\left|\right|{\stackrel{\·}{P}}_{i,j}\left(t\right)\left|\right|\mathrm{dt}$ ⑥；

式⑥中：P_i,j(t)是连接像素i,j的路径；灰度函数D(x)为距离增量，定义如式⑦所示：

$D\left(x\right)=e^{\frac{E\left(x\right)}{v}},E\left(x\right)=\frac{\left|\right|\▿I\left|\right|}{G_{\mathrm{\σ}}*\left|\right|\▿I\left|\right|+\mathrm{\γ}}$ ⑦；

式⑦中：v是尺度参数；E(x)为边缘测度；为图像的梯度幅值；G_σ为标准偏差，为σ的高斯函数，其中σ表示方差；γ是确保弱边界不被灰度计算影像的一个参数。

本发明提供的一种肺部CT图像中粘连血管型肺结节的分割方法，首先采用经典的区域增长提取肺实质，利用曲率信息对与胸膜壁相连接的肺结节进行修正，获得精确的肺实质感兴趣区域，然后采用混合高斯函数内核图割方法对肺结节图像进行初始自动分割，确定肺结节初始内轮廓，然后通过肺结节中心点的位置向肺结节初始内轮廓等角度地发射出一组射线，使该组射线与初始肺结节内轮廓相交，获得初始标记点，并通过插点均匀提取肺结节的标记点，并沿射线上的初始标记点分别向中心点提取目标点，向内轮廓外围提取背景点，然后再通过Graph cuts法并结合测地距离对肺结节进行精细分割。本发明在提取肺结节区域的时候对肺结节的内部结构进行精确的分析，最终可以将肺结节实质和模糊外围进行区分，有利用医生准确确定肺结节实质，获取精确的肺结节边界和特征，为医生提供有效的诊断依据。本发明通过对两组肺结节真实数据的实验结果表明，该方法可以快速有效地解决困难结节的分割问题。

附图说明

图1是本发明提供的一种肺部CT图像中粘连血管型肺结节的分割方法的流程图；

图2是本发明步骤（Ⅰ）中根据曲率和惯性规律填补肺结节表面凹陷的方法示意图；

图3是根据本发明步骤（Ⅱ）得到的初始肺结节内轮廓的效果示意图；

图4是根据本发明步骤（Ⅲ）提取肺结节标记点的过程示意图，为清楚显示的目的，图4（c）、（d）中去掉了射线；

图5示出了本发明步骤（Ⅳ）中测地距离的改进效果；

图6-10是根据本发明的方法分别分割各种不同肺结节的过程示意图。

具体实施方式

下面结合附图和具体实施方式对本发明作进一步描述。

如图1所示，本发明所提供的一种肺部CT图像中粘连血管型肺结节的分割方法，包括以下步骤：

（Ⅰ）提取肺实质：在输入的肺部图像上提取包含肺结节在内的肺实质轮廓，得到肺实质边缘曲线；

（Ⅱ）预分割肺结节：对步骤（Ⅰ）获得的肺实质进行自动分割，获得其中的肺结节的初始轮廓，再采用基于混合高斯模型的内核图割法对得到的肺结节进行初始分割，获取初始肺结节内轮廓；

（Ⅲ）提取肺结节标记点：在步骤（Ⅱ）得到的肺结节内部取一中心点，通过该中心点向初始肺结节内轮廓发射一组等角度的射线，取射线与初始肺结节内轮廓相交的第一点作为初始标记点，并基于该组初始标记点的平均距离进行初始标记点的插值操作，然后沿每一射线上的初始标记点分别向中心点提取目标点，向内轮廓外围提取背景点；

（Ⅳ）图割法分割肺结节：利用步骤（Ⅲ）提取的目标点和背景点，采用基于最大流算法的图割法对肺结节作进一步的分割，获取最终的肺结节轮廓。

本发明的步骤（Ⅰ）中，可以采用区域生长法并结合边缘检测法提取肺实质时，具体可按照以下方式进行：

（1）利用大律法（OTSU）求得CT图像最优阈值作为终止条件，并利用四角置点区域生长法得到背景区域，进而连通域标记合并小面积区域去除背景。

具体流程如下：

a.阈值以上文得到的最优阈值作为迭代终止条件；

b.分别从图像的四个角点（作为种子点）向中心区域扫描，若扫描的像素值小于阈值即认为是背景区域并将像素值置为0（黑色），不断扫描直至所有背景像素点合并完毕；

c.生长过程：从种子点开始遍历其邻域像素，若像素灰度值小于阈值，则标记为背景，若大于阈值且其四邻域有小于阈值的像素，则标记为边界；

d.将标记为背景的像素置为白色。

（2）以步骤（1）的阈值分割，并将分割的图像取反，然后进行闭运算使肺部气管与肺实质连成一体，得到粗略的肺实质模板。

（3）采用区域生长法生成肺实质掩膜。

具体步骤为：在分割后肺实质模板的每一个区域内，选取图像的中心点为种子点，生长准则为像素的相似性，在种子点的四邻域范围内按照顺时针方向进行搜索，与种子点像素一致的四邻域点被纳入生长区域内，同时作为新的种子点，直到没有满足条件的四邻域点被纳入生长区域内的时候，生长结束。

依据上述方法生长出肺实质掩膜是基于原来CT图像上去除了周围器官、支气管、血管等不期望得到的图像，得到分割后的肺实质图像。

（4）采用边缘检测法提取肺实质边缘点集，并将提取的边缘点集连接为线，得到肺实质轮廓曲线。

此处采用的边缘检测法可使用canny边缘算子，算法步骤可参考期刊IEEE Transactions on Pattern Analysis and Machine Intelligence于1986年第6期第679-698页报道的文献名称为“A computational approachto edge detection”中所记载的内容。具体的算法步骤可简述为：

step1:用高斯滤波器平滑图像；

step2:用一阶偏导的有限差分来计算梯度的幅值和方向；

step3:对梯度幅值进行非极大值抑制；

step4:用双阈值算法检测和连接边缘。

进而，本发明步骤（Ⅰ）中获取肺实质边缘曲线后，根据曲率和惯性规律填补肺实质边缘曲线上的凹陷部分，以修正肺实质边缘曲线，可以按照以下方式进行：

（1）根据式①求得肺结节边界曲线和与其相交的肺实质边界曲线的曲率K，并提取曲率方向改变的点（例如，提取4-6个点）；

$K=\left|\frac{I_y^2{I}_{\mathrm{xx}}-2I_x{I}_y{I}_{\mathrm{xy}}+I_x^2{I}_{\mathrm{yy}}}{\sqrt{I_x^2+I_y^2}}\right|$ ①；

（2）根据式②定义的曲率半径ρ判断曲线各个部分的凹凸性；

$\mathrm{\ρ}=-M\frac{{[{\left(I_x\right)}^2+{\left(I_y\right)}^2]}^{\frac{3}{2}}}{(I_x{I}_{\mathrm{yy}}-I_y{I}_{\mathrm{xx}})}$ ②；

式①、②中,I_x，I_y，I_xx，I_yy分别为曲线的坐标x,y对参变量的一阶、二阶导数，M是经验系数；

（3）求得肺实质边界曲线和肺结节边界曲线相交时的拐点，将曲率小的曲线段去掉；

（4）利用曲率方向改变的点，结合曲线延伸的惯性规律，曲线拟合法修正曲线边缘，得到包含肺结节在内的肺实质轮廓。

以下以图2（a）为例具体说明上述对肺实质边缘曲线的修正过程。

图2（a）为肺部图像的某一帧，图2（b）通过两次区域增长的方法获得肺实质部分，从图像中可以看出在图像的左侧肺实质底部有一个肺结节，由于跟胸膜壁灰度相似，在进行肺实质提出时将其分出了肺实质区域，图2（c）为2（b）的边缘图像。我们将感兴趣区域放到被提出的肺结节边缘上，如图2（d）所示，图2（e）上的点本文中通过曲率获得的曲率改变方向的点，通过这几个点的确定，将原来的边缘图像进行修正，通过惯性规律获得大范围图像的曲率来修补肺实质边缘，图2（f）为修补后的肺实质边缘图像，可以看出将肺结节充分的包含在肺实质内部，完成了肺实质的提取。

进一步地，在本发明的步骤（Ⅱ）中，可以采用边缘检测法（算法步骤可参考期刊IEEE Transactions on Pattern Analysis and MachineIntelligence在1986年第6期第679-698页报道的文献名称为“Acomputational approach to edge detection”中所记载的内容。）进行自动分割，以获取肺实质中肺结节的初始轮廓。具体的实现过程如下：

step1:用高斯滤波器平滑图像；

step2:用一阶偏导的有限差分来计算梯度的幅值和方向；

step3:对梯度幅值进行非极大值抑制；

step4:用双阈值算法检测和连接边缘。

进而，本发明的步骤（Ⅱ）中获取肺结节初始轮廓后，采用基于混合高斯模型的内核图割法（算法步骤可参考期刊IEEE Transactions on ImageProcessing在2011年第20卷第2期第545-557页报道的文献名称为“Salah M B,Mitiche A,Ayed I B.Multiregion image segmentation byparametric kernel graph cuts”中所记载的内容）对得到的肺结节进行自动的初始分割，以获取初始肺结节内轮廓时。

图3给出了5组分割比较困难的肺结节，如图3（a）-图3（e）所示，大部分肺结节为与血管相连或者与胸膜壁相连，肺结节周围的边界比较模糊，肺结节的特征多为毛刺或不规则形状等特征而造成分割结果不理想。本发明首先采用自动分割，获得肺结节初始轮廓，然后采用基于混合高斯模型的内核图割分割法对肺结节进行初始分割，获得初始肺结节内轮廓。后续通过该初始肺结节内轮廓来获得肺结节的初始标记点。

本发明在步骤（Ⅲ）中提供了一种图割法的标记点（种子点）的提取方法，这是本发明的主要创新点。以下以图4为例，具体说明本发明步骤（Ⅲ）中提取肺结节标记点的具体实现方法（本文简称“射线法”）：

如图4所示，图4（a）为采用基于混合高斯模型的内核图割分割法分割的初始肺结节内轮廓，从中可以看出，血管部分被图割分割到肺结节区域；然后在肺结节区域确定一个中心点C，通过中心点C向初始肺结节内轮廓发射一组射线，如图4（b）所示；射线与初始肺结节内轮廓线相交，由于本发明选取射线第一次与初始肺结节内轮廓相交的点，这样就可以把部分血管的点（即射线与图中曲线相交的其他点）去掉，从而获得肺结节的初始标记点，如图4（c）所示；由于在提取初始标记点时，有些标记点因含有一些凸起或者血管区域，导致产生的标记点之间的距离远远大于其他相邻标记点距离，此时可参考该组初始标记点的平均距离对距离大的标记点之间进行插点操作，如图4（d）所示，图4（d）中的点I即为新添加的标记点；然后，对每一射线，以该射线上的初始标记点为基准沿射线向中心点C区域平移一定距离作为图割的目标点（如图4（e）中位于肺结节内部的标记点），并向肺结节外围平移一定距离作为图割的背景点（如图4（e）中肺结节外部的标记点）。

在上述提取目标点和背景点时，优选情况下，所提取的目标点和背景点到初始标记点的距离分别为中心点到初始标记点距离的1/3。

本发明在步骤（Ⅲ）中通过标记点的精确确定，就可以通过下面步骤（Ⅳ）的图割法对肺结节图像进行更精细的分割。

本发明步骤（Ⅳ）中采用基于最大流算法的图割法分割肺结节时，可参考期刊IEEE Transactions on Pattern Analysis and MachineIntelligence在2001年第23卷第11期第1222-1239页报道的文献名称为“Fast approximate energy minimization via graph cuts”所记载的方法步骤。

在本发明的一种实施方式中，步骤（Ⅳ）中采用基于最大流算法的Graph cuts法分割肺结节的具体实现过程如下：

（1）以步骤（Ⅲ）中标记的一组背景点所构成的曲线作为当前轮廓线；

（2）构造s-t网络：使当前轮廓线对应s-t网络的源点s，使步骤（Ⅲ）中标记的一组目标点所形成的内边界对应s-t网络的汇点t，其余每一个像素对应s-t网络中的一个结点，相邻结点间以8邻域方式带权值连接；

（3）最大流最小代价切割：计算式③中Gibbs能量函数E的区域项和边界项，用最大流算法对s-t网络进行最小代价切割，得到一条新当前轮廓线；

$E\left(L\right)=\mathrm{\λ}\underset{i\∈V}{\mathrm{\Σ}}R(i,L_i)+\underset{\{i,j\}\∈N}{\mathrm{\Σ}}B(i,j)\·\mathrm{\δ}(L_i\≠L_j)$

$\mathrm{\δ}(L_i\≠L_j)=\left\{\begin{array}{cc}1;& L_i\≠L_j\\ 0;& L_i=L_j\end{array}\right.$ ③；

式③中：V为图像中像素点集合；N表示图像中所有邻接像素对{i,j}；L_i是一个二进制标记，L_i=1时，像素点i标记为前景标记，L_i=0时，像素点i标记为背景标记；而λ表示区域项R(i,L_i)和边界项B(i,j)之间的权重；

区域项R(i,L_i)定义如式④所示：

R(i,L_i)=-ln P(I_i|L_i)        ④；

式④中：I_i表示像素点i的灰度值；P(I_i|L_i)表示灰度值I_i被标记为L_i的概率；

边界项B(i,j)定义如式⑤所示：

$B(i,j)=\exp (-\frac{{(I_i-I_j)}^2}{2{\mathrm{\σ}}^2})\·\frac{1}{D_g(i,j)}\·|L_i-L_j|$ ⑤；

式⑤中：D_g(i,j)表示相邻像素点i和j间的测地距离；

（4）使用步骤（3）的新当前轮廓线代替当前轮廓线，重复步骤（2）-（3），使新当前轮廓线迭代收缩，直至达到图像中确切的肺结节轮廓线；

（5）输出图像中被分割的肺结节图像。

本发明在上述步骤（3）的式⑤中用测地距离D_g(i,j)代替了原来的欧氏距离来得到边界项B(i,j)，并定义如式⑥所示：

$D_g(i,j)=\underset{P_{i.j}}{\min }{\∫}_0^{1}D\left(P_{i,j}\left(t\right)\right)\left|\right|{\stackrel{\·}{P}}_{i,j}\left(t\right)\left|\right|\mathrm{dt}$ ⑥；

式⑥中：P_i,j(t)是连接像素i,j的路径；灰度函数D(x)为距离增量，定义如式⑦所示：

$D\left(x\right)=e^{\frac{E\left(x\right)}{v}},E\left(x\right)=\frac{\left|\right|\▿I\left|\right|}{G_{\mathrm{\σ}}*\left|\right|\▿I\left|\right|+\mathrm{\γ}}$ ⑦；

式⑦中：v是尺度参数；E(x)为边缘测度；为图像的梯度幅值；G_σ为标准偏差，为σ的高斯函数，其中σ表示方差；γ是确保弱边界不被灰度计算影像的一个参数。

对于灰度值比较接近的像素点，本发明通过对测地距离D_g(i,j)的约束，可以更有效的计算标记点到未知像素的距离。

本发明步骤（Ⅳ）中用测地距离代替欧氏距离的改进效果如图5所示。图5中肺结节的两个毛刺部分（位于图5（a）中正方形方框的左下角）的端点距离（即图5（b）中P1、P2之间的距离），若用欧式距离表示，应为线段A的长度，若用测地距离表示，则为曲线B的长度。从图5上看，肺结节的两个毛刺部分的端点距离用测地距离表示时，更能准确反映毛刺部分边缘特征。

应用例

采用本发明提供的方法对一组实验数据（图像）进行分割。

实验数据由两部分构成：

第一部分:由天津医科大学总医院放射科提供26个病例，数据采用16排或者64排薄层螺旋CT扫描，扫描参数层厚10mm或者5mm，扫描间隔10mm或者5mm，120kv，280mA，机架旋转一周时间0.35-0.8s，DFOV36.0cm，螺距1.375:1，扫描矩阵512*512、重建层厚1.25mm，间距1.25mm，基本实现图像体素各相同性以确保容积分割的准确性。

第二部分：来自于美国国家癌症研究会颁布的肺部图像数据库联盟（Lung Image Database Consortum，LIDC），从该数据库中提取68个病例的138张全肺CT扫描图像。

此处从以上病例中选取了5组比较有代表性的全肺CT扫描图像，采用本发明提供的方法进行分割。以下结合图6-10进行说明。

图6-10中分别给出了分割连续帧肺结节的分割结果。

图6-10中的第一行图像分别是肺结节原图像，图像中肺结节的特点是，肺结节部分与血管连接区域较大，血管灰度与肺结节灰度一致，即使利用医生交互分割，也需要在肺结节周围点击很多的标记点才能完成精确分割。

图6-10中的第二行图像分别是肺结节通过初始的内核图割法分割所得到的结果，可以看出，图6和图10这两帧肺结节与血管连接紧密，初始分割时把血管部分体积分割到肺结节的内核区域。

图6-10中的第三行图像分别给出了通过射线获取的目标点和背景点。

图6-10中的第四行图像分别给出了通过上文描述的插点操作所得到的目标点和背景点，可以看出通过改进提取的标记点为理想的标记点。

图6-10中的第五行图像分别给出了通过基于最大流算法的图割法并结合测地距离所获取的肺结节轮廓。

从图6-10显示的分割结果可以看出，通过本发明提供的上述方法获取的标记点对肺结节进行分割，能有效去除血管等其他组织，获得精确的肺结节。

上述实施例只是对本发明的举例说明，本发明也可以以其它的特定方式或其它的特定形式实施，而不偏离本发明的要旨或本质特征。因此，描述的实施方式从任何方面来看均应视为说明性而非限定性的。本发明的范围应由附加的权利要求说明，任何与权利要求的意图和范围等效的变化也应包含在本发明的范围内。

Claims (8)

1.一种肺部CT图像中粘连血管型肺结节的分割方法，包括以下步骤：

（Ⅰ）提取肺实质：在输入的肺部图像上提取包含肺结节在内的肺实质轮廓，得到肺实质边缘曲线；

（Ⅱ）预分割肺结节：对步骤（Ⅰ）获得的肺实质进行自动分割，获得其中的肺结节的初始轮廓，再采用基于混合高斯模型的内核图割法对得到的肺结节进行初始分割，获取初始肺结节内轮廓；

（Ⅲ）提取肺结节标记点：在步骤（Ⅱ）得到的肺结节内部取一中心点，通过该中心点向初始肺结节内轮廓发射一组等角度的射线，取射线与初始肺结节内轮廓相交的第一点作为初始标记点，并基于该组初始标记点的平均距离进行初始标记点的插值操作，然后沿每一射线上的初始标记点分别向中心点提取目标点，向内轮廓外围提取背景点；

（Ⅳ）图割法分割肺结节：利用步骤（Ⅲ）提取的目标点和背景点，采用基于最大流算法的图割法对肺结节作进一步的分割，获取肺结节的最终轮廓。

2.根据权利要求1所述的一种肺部CT图像中粘连血管型肺结节的分割方法，其特征在于，步骤（Ⅰ）中，采用区域生长法并结合边缘检测法提取肺实质。

3.根据权利要求1或2所述的一种肺部CT图像中粘连血管型肺结节的分割方法，其特征在于，步骤（Ⅰ）中得到肺实质边缘曲线后，根据曲率和曲线延伸的惯性规律修正肺实质边缘曲线。

4.根据权利要求3所述的一种肺部CT图像中粘连血管型肺结节的分割方法，其特征在于，根据曲率和曲线延伸的惯性规律修正肺实质边缘曲线，按照以下方式进行：

（1）根据式①求得肺结节边界曲线和与其相交的肺实质边界曲线的曲率K，并提取曲率方向改变的点；

$K=\left|\frac{I_y^2{I}_{\mathrm{xx}}-2I_x{I}_y{I}_{\mathrm{xy}}+I_x^2{I}_{\mathrm{yy}}}{\sqrt{I_x^2+I_y^2}}\right|$ ①；

（2）根据式②定义的曲率半径ρ判断曲线各个部分的凹凸性；

$\mathrm{\ρ}=-M\frac{{[{\left(I_x\right)}^2+{\left(I_y\right)}^2]}^{\frac{3}{2}}}{(I_x{I}_{\mathrm{yy}}-I_y{I}_{\mathrm{xx}})}$ ②；

式①、②中,I_x，I_y，I_xx，I_yy分别为曲线的坐标x,y对参变量的一阶、二阶导数，M是经验系数；

（3）求得肺实质边界曲线和肺结节边界曲线相交时的拐点，将曲率小的曲线段去掉；

（4）利用曲率方向改变的点，结合曲线延伸的惯性规律，曲线拟合法修正曲线边缘，得到包含肺结节在内的肺实质轮廓。

5.根据权利要求1所述的一种肺部CT图像中粘连血管型肺结节的分割方法，其特征在于，步骤（Ⅱ）中采用边缘检测法自动分割肺实质以获取肺结节的初始轮廓。

6.根据权利要求1所述的一种肺部CT图像中粘连血管型肺结节的分割方法，其特征在于，步骤（Ⅲ）中，目标点和背景点到初始标记点的距离分别为中心点到初始标记点距离的1/3。

7.根据权利要求1所述的一种肺部CT图像中粘连血管型肺结节的分割方法，其特征在于，步骤（Ⅳ）的实现过程如下：

（1）以步骤（Ⅲ）中标记的一组背景点所构成的曲线作为当前轮廓线；

（2）构造s-t网络：使当前轮廓线对应s-t网络的源点s，使步骤（Ⅲ）中标记的一组目标点所形成的内边界对应s-t网络的汇点t，其余每一个像素对应s-t网络中的一个结点，相邻结点间以8邻域方式带权值连接；

（3）最大流最小代价切割：计算式③中Gibbs能量函数E的区域项和边界项，用最大流算法对s-t网络进行最小代价切割，得到一条新当前轮廓线；

$E\left(L\right)=\mathrm{\λ}\underset{i\∈V}{\mathrm{\Σ}}R(i,L_i)+\underset{\{i,j\}\∈N}{\mathrm{\Σ}}B(i,j)\·\mathrm{\δ}(L_i\≠L_j)$

$\mathrm{\δ}(L_i\≠L_j)=\left\{\begin{array}{cc}1;& L_i\≠L_j\\ 0;& L_i=L_j\end{array}\right.$ ③；

式③中：V为图像中像素点集合；N表示图像中所有邻接像素对{i,j}；L_i是一个二进制标记，L_i=1时，像素点i标记为前景标记，L_i=0时，像素点i标记为背景标记；而λ表示区域项R(i,L_i)和边界项B(i,j)之间的权重；

区域项R(i,L_i)定义如式④所示：

R(i,L_i)=-ln P(I_i|L_i)            ④；

式④中：I_i表示像素点i的灰度值；P(I_i|L_i)表示灰度值I_i被标记为L_i的概率；

边界项B(i,j)定义如式⑤所示：

$B(i,j)=\exp (-\frac{{(I_i-I_j)}^2}{2{\mathrm{\σ}}^2})\·\frac{1}{D_g(i,j)}\·|L_i-L_j|$ ⑤；

式⑤中：D_g(i,j)表示相邻像素点i和j间的测地距离；

（4）使用步骤（3）的新当前轮廓线代替当前轮廓线，重复步骤（2）-（3），使新当前轮廓线迭代收缩，直至达到图像中确切的肺结节轮廓线；

（5）输出图像中被分割的肺结节图像。

8.根据权利要求7所述的一种肺部CT图像中粘连血管型肺结节的分割方法，其特征在于，

测地距离D_g(i,j)定义如式⑥所示：

$D_g(i,j)=\underset{P_{i.j}}{\min }{\∫}_0^{1}D\left(P_{i,j}\left(t\right)\right)\left|\right|{\stackrel{\·}{P}}_{i,j}\left(t\right)\left|\right|\mathrm{dt}$ ⑥；

式⑥中：P_i,j(t)是连接像素i,j的路径；灰度函数D(x)为距离增量，定义如式⑦所示：

$D\left(x\right)=e^{\frac{E\left(x\right)}{v}},E\left(x\right)=\frac{\left|\right|\▿I\left|\right|}{G_{\mathrm{\σ}}*\left|\right|\▿I\left|\right|+\mathrm{\γ}}$ ⑦；

式⑦中：v是尺度参数；E(x)为边缘测度；

为图像的梯度幅值；G_σ为标准偏差，为σ的高斯函数，其中σ表示方差；γ是确保弱边界不被灰度计算影像的一个参数。

Images