CN107240114B - 一种基于点距离函数形状约束的半自动医学图像分割方法 - Google Patents

Abstract

本发明公开一种基于点距离函数形状约束的半自动医学图像分割方法，包括以下步骤：选取待处理医学图像，并定义图像中平面上N点距离函数；将点距离函数的形状约束融入到变分框架，得到基于点距离函数形状约束的活动轮廓模型；求解活动轮廓模型；求解梯度流方程，得到图像的分割。本发明能够灵活地描述如圆、类圆、椭圆、超椭圆、曲边三角形、曲边四边形和心形等形状，能够有效地分割带有边界缺失的医学图像，并且不需要建立目标形状数据库或者优化形状参数。

Description

一种基于点距离函数形状约束的半自动医学图像分割方法

技术领域

本发明属于医学图像处理技术领域，具体涉及一种基于点距离函数形状约束的半自动医学图像分割方法。

背景技术

医学图像中，常常存在各种原因引起的图像数据缺失的情况，表现在图像上为器官、组织和病灶的边界缺失。此时，仅利用图像的灰度信息很难分割出感兴趣的组织或者器官轮廓。针对目标图像数据丢失的图像分割，现有一种常有的分割方法是结合图像的先验形状信息。

现有图像分割中先验形状的表示方法主要有两种：

一种是通过建立目标形状的数据库，利用数据分析(如聚类、机器学习等方法)得到目标先验形状的表示方法，此种方法的缺陷是建立目标先验形状数据库需要耗费大量的人力物力，并且对于大部分试验而言，这样的数据库往往很难得到；

另一种图像先验形状的表示方法是利用待分割目标所具有的共有的形状特征，将目标的先验形状用参数化的曲线或曲面来描述，该方法需要估计先验形状的参数。

Cootes等人提出了一种主动形状模型(ASM)。主动形状模型的基本思想是以目标物体的轮廓作为训练样本，用点分布模型(point-distribution model)进行描述，构造样本的先验模型。点分布模型主要体现了训练集的平均形状和变形方式。在先验模型建立后，ASM方法通过灰度匹配搜索过程寻找图像中的目标轮廓，在迭代搜索的同时调整先验模型的参数，使模型与目标物体边界吻合。Leventon等[5]是第一个提出将形状先验信息融合到变分框架中来实现图像分割的。随后，Chen等人[6]将心脏的先验形状融入到GAC模型中得到：

上式中第一项为几何测地活动轮廓线的水平集形式，第二项为待分割目标的先验形状约束项。其中C^*为一条曲线，代表待分割目标的先验形状。ν,R,T分别表示伸缩、旋转和平移参数。d(C^*,(x,y))表示点(x,y)与曲线C^*的距离。形状约束项约束分割曲线的形状，使其通过伸缩、旋转或平移能够逼近给定的先验形状C^*，并且在能量泛函的约束下与几何测地活动轮廓线分割的结果达到一个平衡，得到分割结果。Saroul等人提出基于Rayleigh分布和以广义超椭圆为先验形状的超声前列腺组织分割模型；但是超椭圆隐函数中有许多形状参数要随着水平集的演化而优化，使算法的计算时间和计算量大大增加。

[1]X.Bresson,S.Esedoglu,P.Vandergheynst,J.-P.Thiran and S.Osher,Fastglobal minimization of the active contour/snake model[J],Journal ofMathematical Imaging and Vision,2007,28:151-167.

[2]J.Huang,X.Yang and Y.Chen.A fast algorithm for global minimizationof maximum likelihood based on ultrasound image segmentation[J].InverseProblems and Imaging,2011,5(3):645-657.

[3]J.Huang,X.Yang,Y.Chen,L.Tang.Ultrasound kidney segmentation with aglobal prior shape[J].Journal of Visual Communication and ImageRepresentation,2013,24(7):937-943.

[4]L.Gong,S.D.Pathak,D.R.Haynor,P.S.Cho,Y.Kim.Parametricshapemodeling using Deformable superellipses for prostate segmentation[J].IEEE Transactions on Medical Imaging,2004,23(3),340-349.

[5]M.E.Leventon,O.Faugeras,W.E.L.Grimson,et al.Level set basedsegmentation with intensity and curvature priors[C],IEEE Workshop onMathematical Methods in Biomedical Image Analysis,2000:4-11.

[6]Y.Chen,H.Tagare,S.R.Thiruvenkadam,Using prior shapes in geometricactive contours in a variational framework,International journal of computervision,50(3),(2002),315-328.

[7]T.F.Cootes,C.J.Taylor,D.H.Cooper,J.Graham.Active shape models-their trainin and application[J].Computer Vision and Image Understanding,1995,61(1),38-59.

以上文献或多或少公布了部分类似方法，但是现有的这些先验形状或是通过对收集的目标物体形状的分析而得到的，或者是需要在求解的过程中优化大量的形状参数；然而收集的目标形状是通过有经验的放射学专家人工手绘得到的，通常不易获得此类先验形状的数据库，而且人工手绘非常耗时，且针对同一器官或组织，不同的放射学专家也可能得到不一样的分割结果。而且上述的形状模型可以描述的形状多为单一形状，训练、优化一个模型后往往只能分割一种单一形状的目标。

发明内容

发明目的：本发明的目的在于解决现有技术中存在的不足，提供一种基于点距离函数形状约束的半自动医学图像分割方法，利用点距离函数来隐式的表达形状信息，能够灵活地描述如圆、类圆、椭圆、超椭圆、曲边三角形、曲边四边形和心形等形状。

技术方案：本发明一种基于点距离函数形状约束的半自动医学图像分割方法，包括以下步骤：

(1)选取待处理医学图像，并定义图像中平面上N点距离函数；

(2)将点距离函数的形状约束融入到变分框架，得到基于点距离函数形状约束的活动轮廓模型；

(3)快速算法求解步骤(2)中建立的基于点距离函数形状约束的活动轮廓模型；

(4)利用主对偶算法求解步骤(3)中建立的梯度流方程，得到图像的分割。

进一步的，所述步骤(1)中的N点距离函数定义如下：

假设{P₁,P₂,,P_N}为R²中给定的点，{D₁,D₂,,D_N}是每个点P_b＝(x_b,y_b)对应的欧氏距离：

则N点距离函数定义为：

为平面的中点P和给定点{P₁,P₂,,P_N}之间距离的线性组合，代表平面中点P和给定点{P₁,P₂,,P_N}之间的某种距离，α_b为实常数且b＝1,,N，对任意正常数T记 表示一个区域，的边界与某些形状吻合；如果所有的α_b都为正数，则的边界为凸的，如果α_b中某些为负，的边界可能与某些凹边界吻合。例如：针对不同器官、组织或病灶的医学图像自身特点，此处，点距离函数包括N＝1,2,3,4和α_b＝±1。

进一步的，所述步骤(2)中将步骤(1)定义的点距离函数融入全局极小化的变分框架中，进而得到基于点距离函数约束的活动轮廓模型：

其中λ₁,λ₂,μ为权重参数,H_ε(φ)是正则化的Heaviside函数,即

D是步骤(1)中定义的点距离函数的归一化，即D∈(0,1)。d₁和d₂是数据项，不同类型的图像可以选择不同的d₁和d₂，由于医学图像中经常存在图像灰度不均的情况，需要更加鲁棒的纹理描述算子来刻画图像的纹理；因此此处取d₁和d₂如下：

d₁＝1-Q(t(y),t_in(y)),d₂＝1-Q(t(y),t_out(y)),

Q是一个高斯核函数，t＝[m(x),s(x)]是一个简单的纹理描述算子，用以描述图像的纹理信息。Q具有以下形式：

m(x),s(x)分别代表图像灰度的均值和方差。t_in和t_out表示轮廓内外的图像灰度均值和方差组成的纹理描述算子，具体如下：

Ω₁代表图像中轮廓内的区域，I(y)表示待分割图像上像素点y的灰度值。

进一步地，能量泛函(#)关于φ的变分，得到控制水平集函数演化的Euler-Lagrange方程：

其中δ_ε(φ(x))表示H_ε(φ(x))的导数，即

进一步的，所述步骤(3)的具体方法为：

去掉δ_ε(φ)项，得到与公式(*)稳态解相吻合的简化的梯度流方程：

基于上述梯度流，提出新的能量泛函：

对某些α∈(0,1)，通过对最小化子φ进行阈值化得到待分割的目标区域：

Ω₁＝{x,φ(x)≥α},Ω₂＝{x,φ(x)<α}；

D,d₁和d₂分别是能量泛函中最主要的部分，d₁和d₂主要关注的重点是图像中的纹理信息，具有相似纹理特征的像素点将会被归为一类。另一方面，能量泛函中的另外两项D·d₁和(1-D)·d₂，通过将D与d₁相乘，对于给定的常数T>0，算法给D＝T上的像素点以相同的权重。位于点距离函数的越小的T水平集上的像素点会越靠近用户选定的点。因此在迭代的过程中，这些点被当成目标的可能性越大。此时，能量泛函中的D·d₁会变得很小。相似的，通过乘以1-D，远离用户选定点的像素点则更有可能被认为是背景。此时，能量泛函中的(1-D)·d₂也会很小。水平集的演化过程会在能量极小的时候停止。由于在能量极小化的过程中不需要优化形状参数，所以此处基于点距离函数的形状先验在实际操作中更为灵活。

此处，基于点距离函数的先验形状表现为指定位于D＝T上的像素点的权重，相同的T上的像素点的权重相同。分割目标的过程中，水平集的演化不仅由先验形状驱动，还由图像的灰度和纹理信息驱动，所以此处基于点距离函数的形状先验的活动轮廓分割模型在分割给定目标的时候能够得到较好的分割结果。

进一步的，所述步骤(4)的具体方法为：

采取主对偶算法来求解步骤(3)得到的梯度流方程(**)，(**)中的完全变分项可以写成如下对偶形式：

其中，A＝{ω:|ω(x,y)|_∞≤μ,ω|＝0}，因此对于给定的t_in和t_out，得到以下的极小化问题：

其中

R(x,d₁,d₂,D)＝λ₁d₁D-λ₂d₂(1-D)。

给定任意的第k迭代的中间解(φ^k,ω^k)，主对偶算法的具体求解过程为：

1)给定φ＝φ^k,考虑极大化问题：

上升方向为：按照如下式子更新ω：

这里

2)给定ω＝ω^k+1,考虑极小化问题：

下降方向为-div(ω^k+1)+R(x,d₁,d₂,D)，按照下式更新φ：

φ^k+1＝max{min{φ^k+1,1},0},

φ^k+1＝φ^k+β^k(-div(ω^k+1)+R(x,d₁,d₂,D))

其中τ^k,β^k为第k次迭代的步长。

有益效果：本发明中的点距离形状先验模型能够灵活地描述圆、椭圆、超椭圆、曲边四边形和心形等常见图像，而且本发明中的点距离形状模型不仅能够描述凸的形状还可以描述凹的形状。另外，通过将此形状约束融合到水平集框架中，本发明能够分割多种器官组织和病灶如：肾脏、前列腺以及甲状腺结节等。

与传统的方法相比，本发明具有以下优点：

(1)轮廓的演化由图像灰度驱动的同时还有形状的约束。

(2)本发明不需要优化形状参数或训练大量的目标形状。

(3)本发明中的模型可分割的形状更加的灵活，不同的形状只需要选择不同的初始点来解决。

(4)本发明的图像灰度的均值和方差被融合到变分模型中用来描述图像的纹理信息。

附图说明

图1为本发明的整体流程图；

图2为本发明与现有技术在甲状腺结节超声图像的分割结果的比较示意图；

图3为本发明与医生手动和现有技术在肾脏超声图像的分割结果的比较示意图；

图4为本发明与现有技术在前列腺MR图像分割结果的比较示意图；

其中，图4(a)为原始的CT图像，图4(b)为实施例中使用三点距离形状约束模型的分割结果示意图，图4(c)为采用文献[1]中算法的分割结示意图。

具体实施方式

下面对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。

如图1所示，发明的一种基于点距离函数形状约束的半自动图像分割方法，包括以下步骤：

(1)选取待处理医学图像，并定义图像中平面上N点距离函数；

(2)将点距离函数的形状约束融入到变分框架，得到基于点距离函数形状约束的活动轮廓模型；

(3)快速算法求解步骤(2)中建立的基于点距离函数形状约束的活动轮廓模型；

(4)利用主对偶算法求解步骤(3)中建立的梯度流方程，得到图像的分割。

实施例1：

本实施例将基于一个点的距离函数形状约束的半自动图像分割方法应用于甲状腺结节超声图像的分割，如图2所示，第一列是原始超声图像，第二列是本实施例中对应的分割结果图，最后一列是采用文献[1]中方法的分割结果图。由于文献[1]中的能量泛函没有形状约束项，所以得到的都是过分割的结果。

实施例2：

本实施例将基于两个点的距离函数形状约束的半自动图像分割方法应用于肾脏超声图像的分割。

如图3所示，第一行图像为原始的肾脏超声图像，第二行为医生手动分割的结果图，第三行为实施例中使用两点距离形状约束模型的对应分割结果图。正常人类肾脏的形状类似于豌豆。本发明的基于两个点的点距离函数生成的形状先验为一系列的椭圆。第四行是文献[1]中算法的分割结果示意图。最后一行是文献[3]中算法的分割结果示意图。文献[3]在变分框架中引入了一个参数化的超椭圆形状先验来分割超声肾脏图像，该方法中有许多参数需要优化如：旋转、尺度、伸缩等，算法的时间复杂度较高。但本发明的方法，仅需要人工选取两个点来生成形状约束。

实施例3：

本实施例将基于三个点的距离函数形状约束的半自动图像分割方法应用于肝脏CT图像中肿瘤的分割。

如图4所示，图4(a)为原始的CT图像，图4(b)为实施例中使用三点距离形状约束模型的分割结果图，图4(c)为采用文献[1]中算法的分割结果图。肝脏CT图像中，肿瘤内部存在钙化使得灰度不均匀，并且与周边正常肝组织对比度低。肝脏周围还有肋骨、软组织等，使得图像灰度的异质性高。因此，对于不加形状约束的算法(如[1])即使我们将初始轮廓设置在肿瘤边界附近，仍然难以得到肿瘤的正确边界。但是由于实施例中的算法的能量泛函中结合了形状约束，能够准确找到肿瘤的边界。

综上所述，本发明能够有效地分割带有边界缺失的医学图像，并且不需要建立目标形状数据库或者优化形状参数，本发明在多种器官组织和病灶图像的分割上做了验证，本发明大大优于已有分割方法。

Claims (1)

1.一种基于点距离函数形状约束的半自动医学图像分割方法，其特征在于：包括以下步骤：

(1)选取待处理医学图像，并定义图像中平面上N点距离函数；

(2)将点距离函数的形状约束融入到变分框架，得到基于点距离函数形状约束的活动轮廓模型；

(3)快速算法求解步骤(2)中建立的基于点距离函数形状约束的活动轮廓模型；

(4)利用主对偶算法求解步骤(3)中建立的梯度流方程，得到图像的分割；

所述步骤(1)中的N点距离函数定义如下：

假设{P₁,P₂,…,P_N}为R²中给定的N个点，{D₁,D₂,…,D_N}是每个点P_b＝(x_b,y_b)对应的欧氏距离：

则N点距离函数定义为：

为平面上点P和给定点{P₁,P₂,…,P_N}之间距离的线性组合，代表平面上点P和给定点{P₁,P₂,…,P_N}之间的某种距离，α_b为实常数，b＝1,…,N,对任意正常数T，记 表示一个区域，的边界与某些形状吻合；如果所有的α_b都为正数，则的边界为凸的，如果α_b中某些为负，的边界可能与某些凹边界吻合；

其中，步骤(2)中将步骤(1)定义的点距离函数融入全局极小化的变分框架中的具体方法为：

其中λ₁,λ₂,μ为权重参数,H_ε(φ)是正则化的Heaviside函数,即

D是点距离函数的归一化，即D∈(0,1)，d₁和d₂是数据项，不同类型的图像可以选择不同的d₁和d₂，由于医学图像中经常存在图像灰度不均的情况，需要更加鲁棒的纹理描述算子来刻画图像的纹理，因此此处取d₁和d₂如下：

d₁＝1-Q(t(y),t_in(y)),d₂＝1-Q(t(y),t_out(y)),

Q是一个高斯核函数，t＝[m(x),s(x)]是一个简单的纹理描述算子，用以描述图像的纹理信息，Q具有以下形式：

m(x),s(x)分别代表图像灰度的均值和方差；t_in和t_out表示轮廓内外的图像灰度均值和方差组成的纹理描述算子，具体如下：

Ω₁代表图像中轮廓内的区域，I(y)表示待分割图像上像素点y的灰度值；

能量泛函(#)关于φ的变分，得到控制水平集函数演化的Euler-Lagrange方程：

其中δ_ε(φ(x))表示H_ε(φ(x))的导数，即

上述步骤(3)中使用快速算法求解活动轮廓模型的具体方法为：

去掉δ_ε(φ)项，得到与公式(*)稳态解相吻合的简化的梯度流方程：

基于上述梯度流，提出新的能量泛函：

对某些α∈(0,1)，通过对最小化子φ进行阈值化得到待分割的目标区域：

Ω₁＝{x,φ(x)≥α}

Ω₂＝{x,φ(x)<α}

对于给定的常数T>0，给D＝T上的像素点以相同的权重；

上述步骤(4)的具体方法为：

采取主对偶算法来求解步骤(3)得到的梯度流方程(**)，(**)中的完全变分项可以写成如下对偶形式：

其中，因此对于给定的t_in和t_out，得到以下的极小化问题：

其中

R(x,d₁,d₂,D)＝λ₁d₁D-λ₂d₂(1-D)；

给定任意的第k迭代的中间解(φ^k,ω^k)，主对偶算法的具体求解过程为：

1)给定φ＝φ^k,考虑极大化问题：

上升方向为：按照如下式子更新ω：

这里

2)给定ω＝ω^k+1,考虑极小化问题：

下降方向为-div(ω^k+1)+R(x,d₁,d₂,D)，按照下式更新φ：

其中τ^k,β^k为第k次迭代的步长。