自适应的基于先验形状的图像分割方法

摘要

一种图像处理技术领域的基于先验形状的图像分割方法，采用整数符号函数克服由于噪声的干扰对图像分割的影响，并针对其需要手动调节先验形状模型和传统活动轮廓模型的权值系数，提出约束变分模型使得该权值系数可以自适应的收敛到稳定值，同时以识别为基础的形状模板选择用以在分割时候确定采用哪个先验形状的形状模板，避免现有技术中得不到基于先验形状模型的分割结果的问题。

说明书

技术领域

本发明涉及的是一种图像处理的方法，具体是一种自适应的基于先验形状的图像分割方法。

背景技术

活动轮廓模型和水平集方法已经广泛的应用在图像处理和机器视觉当中。然而，当目标的边缘和背景间的差异不是很大的时候，演化曲线将会出现泄漏，导致分割失败。而且，如果图像中存在噪声、杂波或者遮挡的时候，普通的分割模型不能提取出有意义的目标。一个有效的办法是将目标的一些先验知识融入活动轮廓模型的框架中。低层的先验知识，例如：灰度级、颜色、纹理和运动等信息往往不足以表示目标特征，在实际应用中通常不能保证这些特征的不变性。近些年来，高层的知识，尤其是形状被很多学者所重视，并将其融入到分割模型中。这类模型往往被称为基于先验形状的图像分割方法。该方法有两个难点：首先如何将先验形状表示、描述成可以被分割模型利用的形状模型，其次如何将该形状模型融入到传统的分割模型中。

很多学者对这一问题进行了研究。目前对该问题的研究大致分为两类：统计先验模型和基于形状模板匹配的模型。

在统计先验模型中，M.E.Levento等在CVPR(计算机视觉与模式识别国际会议)(2000，1：1316-1323)上发表了“Statistical shape influence in geodesic active contours”(“统计先验模型对测地线活动轮廓模型的影响”)一文，该文率先提出了通过估计形状的概率密度，然后将形状的概率密度和活动轮廓模型融入到贝叶斯框架中的方法。在该文中，主成分分析和核主成分分析被用来提取形状特征，用于估计形状的概率密度函数。另外，M.Rousson和D.Cremers在MICCAI(医学影像计算与计算机辅助介入国际会议)(2005：757-764)上发表了“Eficientkernel density estimation of shape and intensity priors for level set segmentation”(“引入对形状核密度估计的水平集分割方法”)一文，该文对比了高斯分布、一致分布和核密度估计，证明了核密度估计对形状的概率密度函数估计的有效性。

在基于形状模板匹配的模型中，N.Paragios等在ECCV(欧洲计算机视觉国际会议)(2002：775-789)上发表了“Matching distance functions：a shape-to-area variational approach forglobal-to-local registration”(“匹配距离函数-用于全局到局部配准的结合形状和区域的配准方法”)一文，该文提出采用单一形状模板作为形状模型，通过配准演化曲线和形状模板，从而使得演化曲线朝着形状模板演化，通过手动调整该配准项和传统活动轮廓模型间的权值系数即可实现基于形状模板的图像分割。

上诉基于先验形状的分割模型虽然在当背景中存在噪声、杂波和遮挡的时候也能保证一定的分割精度，但是仍然存在以下问题没有解决：将先验形状模型和传统活动轮廓模型相结合时，需要手动调节这两项间的权值系数，并且没有任何准则来指导如何调整该系数。如果该系数调节不当，往往得不到基于先验形状模型的分割结果；大多数的形状模型采用符号距离函数来描述图像，然而由于符号距离函数所在空间是非线性的，通过线性加权得到的统计形状模型和参数形状模板所描述的图像已经不再满足符号距离函数。

发明内容

本发明针对现有技术存在的上述两点不足，提供一种自适应的基于先验形状的图像分割方法，采用整数符号函数克服由于噪声的干扰对图像分割的影响，并针对其需要手动调节先验形状模型和传统活动轮廓模型的权值系数，提出约束变分模型使得该权值系数可以自适应的收敛到稳定值，同时以识别为基础的形状模板选择用以在分割时候确定采用哪个先验形状的形状模板，避免现有技术中得不到基于先验形状模型的分割结果的问题。

本发明是通过以下技术方案实现的，本发明包括以下步骤：

第一步、将目标图像的形状表示为形状模板作为先验形状模板库，初始化演化曲线。

所述的整数符号函数的形状模板是指：以双边的链表结构C表示形状轮廓，该双边的链表结构的内边缘为L_in，外边缘为L_out，在C内部且不在L_in中的点为C_in，在C外部且不在L_out中的点为C_out，则有作为形状模板的整数符号函数：其中：(x，y)为目标图像中任一一点的坐标，a和b为整数值，本发明中a＝1，b＝3。

所述形状模板是指：将先验形状的轮廓作为L_out，根据整数符号函数的定义将先验形状表示为整数符号函数。

所述的演化曲线是指：在图像平面中的闭合曲线，初始化的时候一般需要手动。在本发明中同样采用整数符号函数描述在图像平面中的闭合曲线，其中将演化曲线作为L_out。

第二步、将演化曲线和先验形状模板库中的形状模板进行逐一配准，并根据配准结果更新形状模板库：

所述的配准是指：使得演化曲线和形状模板配准。最小化演化曲线和形状模板的误差平方和，优化尺度、旋转和平移使得误差平方和最小，具体步骤为：

2.1)演化曲线和形状模板的平方误差和表示为：E_S(φ，f)＝∫∫_Ω(φ-ψ(f))²dxdy，其中：ψ和φ分别为形状模板和演化曲线的整数符号函数，具体为：其中(u，v)为形状模板平面上的点，其中(x，y)为图像平面上的点；E_S为形状能量，f为刚体变换函数，具体为包括平移项T_x，T_y、旋转项θ和缩放项s，(T_x，T_y表示在x方向和y方向的平移)；

2.2)作为尺度、旋转和平移的梯度下降方程具体为：其中：p∈{T_x，T_y，θ，s}。

2.3)根据平移项T_x，T_y、旋转项θ和缩放项s更新形状模板。

第三步、利用形状相似性度量从所述先验形状模板库中选择形状模板：采用部分豪斯道夫距离方法度量演化曲线与更新后的形状模板的相似性程度，并取出相似性程度最大的形状模板作为识别结果；

所述的相似性程度是指：H_LK(A，B)＝max(h_L(A，B)，h_K(B，A))，其中：为改进后的有向豪斯道夫距离，代表距离集中的第K个值，同样的意义用于L，点集A＝{a₁，...，a_m}为演化曲线上的所有点，m表示A中点的数目，点集B＝{b₁，...，b_n}为形状模板轮廓上的点，n表示B中点的数目。定义p₁＝K/m和p₂＝L/n，则p₁和p₂的范围在[0，1]，且p₁和p₂的选择取决于目标被遮挡的程度。通过选定p₁和p₂就可以确定L和K。

第四步、利用约束变分模型求得识别结果中的形状和灰度之间的权值系数，并结合形状模板和图像灰度信息分割目标图像得到分割结果。

所述的约束变分模型为：

$\min {\int \int }_{\Omega }{\lambda }_i{(I-c_i)}^2(-\frac{\phi -b}{2b})\mathrm{dxdy}+{\int \int }_{\Omega }{\lambda }_o{(I-c_o)}^2\left(\frac{\phi +b}{2b}\right)\mathrm{dxdy},$

s.t.∫∫_Ω(φ-ψ(f))²dxdy＝α∫∫_Ω(2b)²δ(ψ(f)-a)dxdy，

上述约束变分模型的解析解满足欧拉方程(EE)和LMR，具体如下：

$\frac{{\lambda }_i}{2b}{(I-c_i)}^2-\frac{{\lambda }_o}{2b}{(I-c_o)}^2-2\xi (\phi -\psi \left(f\right))=0,$

∫∫_Ω(φ-ψ(f))²dxdy＝α∫∫_Ω(2b)²δ(ψ(f)-a)dxdy，

其中：ζ为拉格朗日乘子常数，即为将C-V能量函数和形状能量E_S结合起来的权值系数，α是常量，α＝2.5或0.2。

所述的拉格朗日乘子常数通过交替计算以下迭代方程，直到φ和ζ达到稳定：

$\frac{\partial \phi }{\partial t}=\frac{{\lambda }_i}{2b}{(I-c_i)}^2-\frac{{\lambda }_o}{2b}{(I-c_o)}^2-2\xi (\phi -\psi \left(f\right)),$

$\frac{\partial \xi }{\partial t}={\int \int }_{\Omega }{(\phi -\psi \left(f\right))}^2\mathrm{dxdy}-\alpha {\int \int }_{\Omega }{\left(2b\right)}^2\delta (\psi \left(f\right)-a)\mathrm{dxdy},$

其中：为曲线的演化方程，为拉格朗日乘子的迭代方程。

本发明技术效果在于：解决了传统的基于先验形状的图像分割方法需要手动调节权值系数的问题(该权值系数用于平衡形状知识和图像本身的灰度级)；并且，本发明可以自适应的选择最相似的形状模板用来指导图像分割；再者，当图像背景存在噪声、杂波，甚至待分割目标被部分遮挡时，本发明仍然能够分割出想要的目标。

附图说明

图1不同λ取值对分割结果的影响。

图2为待分割海星和先验形状模板库。

图3为海星分割检测结果。

图4为海星分割检测的拉格朗日乘子和部分豪斯道夫距离。

图5为左心室先验形状模板库。

图6为左心室分割检测结果。

图7为左心室分割检测的部分豪斯道夫距离和拉格朗日乘子。

图8为行人分割检测结果。

图9为行人分割检测的部分豪斯道夫距离和拉格朗日乘子。

图10为不同视点观察下的兵马俑分割检测先验形状模板库。

图11为不同视点观察下的兵马俑分割检测结果。

具体实施方式

下面对本发明的实施例作详细说明，本实施例在以本发明技术方案为前提下进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。

本实施例包括以下步骤：

第一步、初始化演化曲线(分割结果)，并将目标图像的形状表示为整数符号函数的形状模板作为先验形状模板库：

将形状表示为整数符号距离函数：假设轮廓为C，将其表示为双边的链表结构：内边缘L_in和外边缘L_out。定义在C内部且不在L_in中的点为C_in；在C外部且不在L_out中的点为C_out。

$\phi (x,y)=\left(\begin{array}{c}-b,(x,y)\in C_{\mathrm{in}}\\ -a,(x,y)\in L_{\mathrm{in}}\\ a,(x,y)\in L_{\mathrm{out}}\\ b,(x,y)\in C_{\mathrm{out}}\end{array}\right),$

其中(x，y)为当前点的坐标，a和b为整数值，本实施例中令a＝1，b＝3。

第二步、将图像分割结果和预置的形状库中的形状模板逐一配准，并根据配准结果更新形状模板库：

假设分别用整数符号函数ψ和φ分别代表某个形状模板和演化曲线。通过最优化其误差平方和(最优化旋转、平移和缩放项)，可将其配准。ψ和φ的误差平方和定义如下：

E_S(φ，f)＝∫∫_Ω(φ-ψ(f))²dxdy

其中，E_S为形状能量，f为刚体变换函数，包括平移项T_x，T_y、旋转项θ和缩放项s。

配准可以通过梯度下降法来实现，并且，刚体变换参数的最大梯度下降方向如下：

$\frac{\partial p}{\partial t}=2{\int \int }_{\Omega }(\psi \left(f\right)-\phi )▿\psi \frac{\partial f}{\partial p}\mathrm{dxdy}$

其中，p∈{T_x，T_y，θ，s}。

根据所求得的平移项T_x，T_y、旋转项θ和缩放项s，更新该形状模板。

第三步、利用形状相似性度量从所述先验形状模板库中选择形状模板：部分豪斯道夫距离在比较两个形状的部分相似性时具有很好的效果。当目标物体存在部分遮挡时，使用部分豪斯道夫距离同样可以对物体进行识别。在前面所述的有向豪斯道夫距离的定义中将全局最大值改为第K个值，即得到改进后的有向豪斯道夫距离：

$h_K(B,A)=K_{a\in A}^{\mathrm{th}}\underset{b\in B}{\min }\left|\right|a-b\left|\right|$

其中，代表距离集中的第K个值。

于是，部分豪斯道夫距离的定义如下：

演化曲线上的所有点构成点集A，形状模板轮廓上的点构成点集B。部分豪斯道夫距离的定义式中含有两个参数p₁＝K/m和p₂＝L/n，p₁和p₂的范围是[0，1]。p₁和p₂的选择取决于目标被遮挡的程度。先验形状模板库里的每一个形状模板都将与演化曲线分别计算部分豪斯道夫距离，与演化曲线最相似的形状模板，即部分豪斯道夫距离值最小的形状模板，将被选择用于指导曲线的演化。

本实施例中采用形状识别方法从形状模板库中选择形状模板，避免了由于符号距离函数的非线性性带来的形状表示的不准确。

在第一阶段中，令p₁＝p₂＝0.8；在第二阶段中，令p₁＝p₂＝0.95；

第四步、利用提出的约束变分模型求得形状和灰度间的权值系数，并结合所述形状模板和图像灰度信息分割所述目标图像得到新的分割结果：本实施例提出的约束变分模型如下：

$\min {\int \int }_{\Omega }{\lambda }_i{(I-c_i)}^2(-\frac{\phi -b}{2b})\mathrm{dxdy}+{\int \int }_{\Omega }{\lambda }_o{(I-c_o)}^2\left(\frac{\phi +b}{2b}\right)\mathrm{dxdy}$

s.t.∫∫_Ω(φ-ψ(f))²dxdy＝α∫∫_Ω(2b)²δ(ψ(f)-a)dxdy

其解析解须满足欧拉方程(EE)和LMR，具体如下：

$\frac{{\lambda }_i}{2b}{(I-c_i)}^2-\frac{{\lambda }_o}{2b}{(I-c_o)}^2-2\xi (\phi -\psi \left(f\right))=0$

ζ＝constant

∫∫_Ω(φ-ψ(f))²dxdy＝α∫∫_Ω(2b)²δ(ψ(f)-a)dxdy

其中，ζ为拉格朗日乘子(LM)，α是一个常量。

这里的拉格朗日乘子ζ即为将C-V能量函数和形状能量E_S结合起来的权值系数。但是，由于难以得到φ和ζ的解析解，本方法给出演化曲线和拉格朗日乘子的迭代计算方法，具体如下：

$\frac{\partial \phi }{\partial t}=\frac{{\lambda }_i}{2b}{(I-c_i)}^2-\frac{{\lambda }_o}{2b}{(I-c_o)}^2-2\xi (\phi -\psi \left(f\right))$

$\frac{\partial \xi }{\partial t}={\int \int }_{\Omega }{(\phi -\psi \left(f\right))}^2\mathrm{dxdy}-\alpha {\int \int }_{\Omega }{\left(2b\right)}^2\delta (\psi \left(f\right)-a)\mathrm{dxdy}$

其中，为曲线的演化方程，为拉格朗日乘子的迭代方程。通过交替地计算这两个方程直到φ和ζ达到稳定，就可以得到最终的演化曲线。

本实施例中，ζ的值是自适应变化的，能够根据图像和形状模板识别的情况自动收敛到稳定状态，这就避免了该权值系数的手动调节。

在第一阶段中，令α＝2.5；在第二阶段中，令α＝0.2；

第五步、如果拉格朗日乘子ζ和水平集函数φ都稳定，即可得到分割结果。

下面以海星分割检测、心脏超声图像中左心室分割检测、行人分割检测和不同视点观察下的兵马俑分割检测来演示本实施例的技术效果。并且，在心脏超声图像中左心室分割检测中，将本实施例的分割效果与传统C-V模型的分割效果进行对比，以显示本实施例的显著优点。

海星分割检测的目的是将图2.(a)中海星分割出来，使用的是如图2.(b)，(c)，(d)所示的三个不同形状模板的先验形状库，其中图2.(b)和图2.(c)分别是两个不同的先验海星形状，而图2.(d)是手形状。海星分割结果如图3所示，其中图3.(b)是第一阶段分割的结果(相应迭代次数为75)，图3.(c)-(f)是分割和配准结果。海星与图像背景及海星自己的影子分离开来，并且其被沙滩遮挡的触手亦被重建出来。图3.(d)和图3.(e)显示：海星形状的形状模板能够与演化曲线相配，而人手形状的形状模板则不能。图4.(a)显示在本实施例方法的两个阶段自适应地计算得到拉格朗日乘子，并且最终收敛于稳定状态。图4.(b)显示部分豪斯道夫距离的值。在图4.(b)中可以看出图2.(c)中的形状模板与演化曲线的相似性最大。即使目标物体存在部分遮挡、或者受杂波、噪声的影响，利用豪斯道夫距离都能够从先验形状集中选择出合适的形状模板。

心脏超声图像中左心室分割检测的目的是从心脏超声图像中分割出左心室。检测所用到的左心室先验形状模板库如图5。图6是采用本实施例方法的分割结果，其中图6.(e)显示了仅利用C-V模型，没有结合先验形状分析的检测结果，图中演化曲线在弱边界处出现了泄漏问题；并且，由于C-V模型以分段常数的形式表示M-S模型，很显然C-V模型不能应对例如图6.(e)中心室目标灰度分布不均匀的情况。而采用本方法的分割方法能够在灰度分布不均匀的区域重建丢失的形状，如图6.(d)所示。图6中(f)-(j)显示了配准结果。配准在第一阶段完成，将其结果投影到第二阶段的分割结果上。在图6.(j)中，图5.(e)的形状模板与演化曲线失配，这是图6.(c)中的边界泄漏造成的影响，从图7.(a)中可看出其相应的部分豪斯道夫距离(标注为e的曲线)相对比较大。图7显示了该检测的部分豪斯道夫距离和拉格朗日乘子。图7.(a)中，标注为“d”的曲线在经过若干次迭代之后部分豪斯道夫距离相对较小，因此它对应的形状模板(图5.(d))被选择出来，用于指导图6中的曲线演化。

行人分割检测的目的是分割出图像中的行人。本实施例从一个人行走的连续图像中选出11张图像作为先验形状库，并将这11张图像重叠在图8.(b)中的初始演化曲线上。虽然目标物身着大衣，采用本实施例的分割方法仍然可以分割出行人本身，如图8.(d)。图9显示了行人分割检测的部分豪斯道夫距离和拉格朗日乘子。图9.(a)中标注为(e)-(o)的曲线的稳定状态对应的是图8.(e)-(o)中的部分豪斯道夫距离。在经过329次迭代之后，拉格朗日乘子的值也将达到稳定，如图9.(b)。需要特别指出的是：在本实施例中，第一阶段α＝2.5，p₁＝p₂＝0.8；考虑到第一阶段的配准结果比较粗略，并且一些部分豪斯道夫距离值比较接近，所以在检测第二阶段没有一次性将参数设置为固定值，而是让这些参数逐渐由第一阶段的值变化为α＝0.3，p₁＝p₂＝0.95，这样做可以让识别的过程和形状模板的选择更加精确。

不同视点观察下的兵马俑分割检测选用如图10所示的25张形状模板构成兵马俑的先验形状模板库，这些形状模板是通过对3D兵马俑多视点观察得到的2D图像。在图11中，人为地在原始待分割图像上添加了高斯噪声，并且对兵马俑添加了部分遮挡。在第一阶段令部分豪斯道夫距离值为0.6，这是因为这种遮挡将对采用部分豪斯道夫距离的形状模板选择带来较为严重的干扰。对于传统的活动轮廓模型，此类噪声和遮挡将会严重影响模型的分割效果，但是使用本实施例的分割方法则可以重建出目标物体遗失的部分。

上述实施例结果证明：只要提供足够多的先验形状模板，即使是对不同视点观察得到的图像，使用本方法的方法仍然可以识别和分割出目标物体。