基于水平集的图像分割方法

一、引言

借鉴一些流体中的重要思想, 1988年,Osher和Sethian首次提出了水平集算法[1],这是一种有效解决曲线演化问题的数值方法,并且计算稳定,适宜任意维数空间。随后,Osher等人对水平集算法做出扩展和总结[2,3], Giga也做了相关的理论扩展[4]。近年来这种算法已被广泛地应用在图像处理领域[5]中 ,尤其在图像分割中已取得了很大的进展。事实上,用水平集来解决图像分割问题的实质就是与活动轮廓模型结合,用水平集方法来求解这些模型得到的偏微分方程PDE(partial differential equation) ,属于边缘检测的分割方法。

二、水平集方法

水平集方法将 n维曲面的演化问题转化为 n +1维空间的水平集函数曲面演化的隐含方式来求解,主要包括 3个要素:超曲面的数据表示,控制曲面演化的一系列 PDE以及相应的数值解法。Osher和Fedkiw对水平集方法做了归纳和总结[6]。为简化问题,这里选择n=2来讨论。水平集的主要思想是将2维平面闭合曲线 C(t)表示为水平集函数φ的零水平集φ(t=0),即将界面嵌入到一个曲面中,将2维曲线的演化转化为 3维曲面演化。具体情况如图 1所示。

图1 水平集方法原理示意图

首先,定义符号距离函数(signed distance function, SDF):

(1)

其中， 表示点 到初始闭合曲线 的最短距离，其符号取决于点 在曲线的内部还是外部，通常在外部取正号，内部取负。

水平集函数曲面的演化遵循下面的Hamilton-Jacobi方程[1]：

(2)

其中， 表示曲线上各点的演化速度，方向沿着曲线的法线方向，通常与图像梯度和曲线曲率有关。我们想要分析和计算的就是在速度 的作用下，曲面随后的演化情况。从式(2)可以看出，只要速度 变化平滑，则 始终保持是一个光滑函数，零水平集始终与运动曲面相对应，曲面的拓扑变化可以很容易的描述。对于不同的分割模型，速度 表达式也不相同，从而出现了多种基于水平集的图像分割方法。另外，可以看出式(2)对任意维数的曲面演化都适用。曲线的几何特性也可以很容易的从水平集函数 得出。计算过程中，可以利用有限差分法结合离散网格在数值上近似求解，水平集函数的梯度也可以利用空间导数得到很好的近似。

通常情况下，速度 只定义在零水平层，而扩展实在整个曲面进行，因此需要将 扩展到整个函数曲面。最初的扩展方法就是全局扩展算法，但这种方法的运算区域是整个平面，计算量达，效率低。针对这一问题，Malladi等人提出了窄带算法及自适应窄带算法，该方法将计算区域由整个图像平面缩减到曲线周围的一个窄带内，计算量大大减少，提高了效率。

这种基于水平集的数值方法以一种隐含的方式表述曲面的运动变化。所有曲面运动相关的信息和运动曲面本身的情况都体现在水平集函数之中。

三、基于梯度信息的活动轮廓模型的水平集图像分割方法

3.1 活动轮廓模型

活动轮廓模型即Snake模型，其基本思想是使初始曲线在一系列的外部约束力和图像内在能量的相互作用下进行演化，直至它满足一定的收敛条件停止在图像的边缘，实现图像的分割。Snake模型的原理就是要计算能量最小，其具体形式为：

(3)

其中，参数 为非负常数， 为初始曲线。式(3)前两项的作用是用来控制轮廓的平滑性，输入内部能量；最后一项是控制分割物体周围的轮廓的，属于外部能量。可以看出，要是能量达到最小，就要把演化曲线 置于梯度 最大的地方，即目标物体的边缘，从而实现图像分割。

随后，基于平均曲率的变化情况，Caselles等人提出了水平集的几何活动轮廓模型，其水平集的曲线演化方程为：

(4)

式中， 是大于零的常数。 是边缘检测子，是一个正的递减函数，它依赖于图像梯度信息，为：

(5)

式中， 是图像 与方差为 的高斯滤波器的卷积，显然在梯度较大处 趋近于0.

同时，Malladi等人提出了两种基于梯度的水平集活动轮廓模型。他们将速度函数 分为两部分： ， 是基于符号控制曲线扩张和收缩的常量， 是依赖于曲线的几何特性(如曲率K)的变量。第一种活动轮廓模型的曲线演化方程为：

(6)

其中， 和 分别表示图像梯度的最大值和最小值。

第二种活动轮廓模型的曲线演化方程为

(7)

可以看出，第二种模型与基于几何的活动轮廓模型式(4)很相似。随后，Caselles等人在基于能量最小化的Snakes模型和基于曲线进化理论的几何活动轮廓模型的基础上提出了测地线活动轮廓模型，对分割具有裂口、缝隙的图像有很好的效果。此模型是建立在几何活动轮廓模型的基础上进行改进的，其曲线演化方程为：

(8)

对于依赖于梯度有缝隙、裂口的图像，采用以上的水平集方法进行图像分割具有一定的优势。

3.2 数值算法

基于PDE的水平集求解的基本步骤是：（1）初始化距离函数及设置各种参数（如时间间隔等）。（2）利用上述曲线演化模型通过迭代方式计算水平集函数直到收敛为止。

该算法的第一步是初始化水平集函数。Adalsteinsson和Sethian提出基于窄带的快速推进法在迭代过程进行水平集函数的重新初始化，大大提高了曲线演化速度。快速推进法（fast marching method）是一种求解边界值问题——Eikonal方程的数值解法。该算法结合了Hamilton-Jacobi方程粘性解法中的上推算法、窄带算法及最小数据堆的数据结构。随后，Zhao提出了快速扫描法（fast sweeping method）来求解Eikonal方程以及Hamilton-Jacobi方程。快速扫描方法将计算复杂度由快速推进发的 降到 。而且因为快速扫描算法不需要堆栈，所以比快速推进法更易实现。

2005年，ChunmingLi提出一种不需要重新初始化的水平集算法，将算法实现复杂度大大降低，同时提高了曲线演化速度，进而提高了算法的效率[8]。

四、基于Mumford-Shad模型的水平集图像分割方法

对于比较模糊或者噪声较大的图像，采用以上的基于梯度信息的图像分割模型可能出现不正确的分割结果。针对这种情况，Chan-Vese提出一种基于Mumford-Shah模型的水平集图像分割方法[7]。此方法不依赖梯度信息，利用计算能量函数最小实现图像分割，因此该方法具有全局最优。

4.1 Chan-Vese图像分割模型

Mumford-Shah模型是一种基于能量最小化的用于图像分割或降噪模型。Mumford-Shah模型的基本形式如下：

（9）

式中， 是非负常数， 表示图像区域，曲线C表示区域的边界， 表示初始图像， 表示逼近原始图像 一个分段光滑图像。要解决的主要问题是最小化能量函数 ，Chan-Vese提出用水平集方法来解决Mumford-Shah问题。此外，他们将Mumford-Shah模型又分为分段光滑和分段常数两种情况。

首先介绍分段光滑（piecewise smooth）Chan-Vese模型。假设C是一条封闭曲线（活动轮廓），它将图像区域 分为曲线的内部和外部，分别表示为 和 。利用水平集方法求解该模型，用水平集函数 来代替未知曲线 ，如果当前点在曲线C外，则 ，反之， ，在曲线上，则 。那么式（9）将变为：

（10）

将 最小化，可得到相应的欧拉-拉格朗日方程：

（11）

上式中， 是Heaviside函数， 是Dirac函数，定义为：

（12）

该模型进行图像分割时，具有平滑和降噪作用。

Marquina-Osher将 替换成 ，此法不依赖于水平集函数的选择消除了Dirac函数对检测远离活动轮廓线边缘的抑制，问题就变成形态学问题。

分段常数（piecewise constant）Chan-Vese模型是Mumford-Shah模型的一个特例。现在以两相（two-phase）情况为例讨论其水平集表示形式。这里只涉及一个水平集函数 ，其能量函数为：

（13）

通过能量最小化，得到相应的欧拉-拉格朗日方程为：

（14）

其中， 和 分别表示区域内部和区域外部的图像灰度平均值， 是正常数。基于Mumford-Shah模型的水平集分割方法初始曲线可选在图像区域的任意位置，而且不依赖图像梯度大小，尤其对含有噪声的图像具有很好的分割效果。

4.2 数值解法

从上述基于Mumford-Shah模型和水平集方法的曲线演化方程可以看出，基于水平集的图像分割方法都需要求解欧拉-拉格朗日方程。为解决欧拉-拉格朗日方程计算量打的问题，人们在Chan-Vese模型的基础上提出一些快速算法。Gibou-Fedkiw定义了 来求解欧拉-拉格朗日方程，其中当 时，相应的像素点值就会改变符号，从而实现水平集曲线的演化。

Chan-Vese模型算法有计算量大、降噪能力差和初始曲线很难给定等缺陷。Gao-Buiti提出一种基于模型的分级算法，该算法将曲线演化过程和降噪过程分开进行，先曲线演化，接着根据水平集曲线分割的图像区域分别降噪处理，这样既提高了运算速度，也可利用其他方法进行降噪增加了算法灵活性。

五、基于水平集的多相位图像分割方法

以上讨论都是基于水平集的两相位图像的分割情况，研究者还对水平集的多相位图像的分割做了大量工作。

Zhao-Chan等提出用n个水平集函数表示n个相位的水平集方法来分割多相位的图像，也就是说每一个分割部分 ，由一个水平集函数 来表示[9]。这种方法要求每个子区域之间互不重叠，其分割是将图像中的每个像素按所在的类归到不同的子区域中。对于包含很多相位的图像，该算法效率往往很低。

Vese-Chan提出用n个水平集函数表示 个相位，该算法食用水平集函数的符号来作为相位的二进制代码[10]。假设图像有4个相位 ， ，那么就需要两个水平集函数 和 来表示：

可以估计多个相位的个数，多余相位会随着收敛的过程而消失。该方法计算效率比Zhao的有所提高，但存在相位错误辨识和存储量大等缺陷。

针对计算存储量大的问题，Lie-Lysaker提出用一个水平集函数的不同水平层来表示不同的相位。定义基本函数：

其中，

问题就变为求解以下能量函数最小化：

式中, 。 如果满足 ,那么存在唯一的 使得每一个 有 。Lie-Lysaker采用了扩展的拉格朗日方法解决以上最小化问题。但此方法在表示曲线单位法矢及曲率方面有一定的困难，而且还用高阶多项式来作为基本函数，正确的分割往往依赖于较好的初始曲线。

六、结论

水平集的图像分割算法已广泛地应用与各个领域，基于梯度的水平集图像分割方法适用与梯度变化很大的图像；基于Mumford-Shah模型的水平集图像分割方法适用与模糊或者噪声较大的图像。多相位图像的所有相位可用多个水平集函数来表示，也可用一个水平集函数的不同水平层来表示，从而达到分割目的。水平集方法是一种简单、精确、灵活的数值方法，主要优势在于处理外形复杂、拓扑结构变化的图像。目前的研究主要基于水平集方法不断的降低计算复杂度，提高算法速度和分割准确度。