CN104867143A - 基于局部引导核拟合能量模型的水平集图像分割方法 - Google Patents

Abstract

一种基于局部引导核拟合能量模型的水平集图像分割方法，主要包括水平集函数的定义、分割模型能量泛函的构建、能量泛函的形式简化、水平集函数的演化和水平集函数的光滑处理。本发明通过将基于局部引导核拟合能量泛函提取图像局部信息，并在每次迭代过程中进行高斯滤波，以避免周期性重新初始化水平集函数。本发明所提出的分割方法不仅提高了灰度不均匀场景中弱边界目标的分割精度，而且完全避免了周期性的重新初始化问题，同时还降低了算法的计算复杂度。

Description

基于局部引导核拟合能量模型的水平集图像分割方法

技术领域

本发明涉及的是一种图像处理技术领域图像分割的方法，具体是一种局部引导核拟合(Local Kernel-induced Fitting，LKF)能量模型水平集图像分割方法。

技术背景

基于变分水平集(Variational Level Set Method)的图像分割方法，凭借其自由拓扑性及多信息共融性，被广泛应用于图像处理领域，例如图像分割、目标提取、目标跟踪。在成像过程中，由于受成像设备固有缺陷和照度不均等因素的影响，实际中的图像广泛存在灰度分布不均匀性。为了能够分割灰度不同质图像，Vese和Chan提出了分段光滑(Piecewise Smooth，PS)模型来解决分段常量(Piecewise Constant，PC)模型不能分割灰度不均匀图像的问题；李纯明等提出了局部二值拟合(Local Binary Fitting，LBF)图像分割模型。但是，该模型在每次迭代过程，需要计算水平集函数的偏微分方程以及周期性重新初始化过程会增大计算量，且对初始化十分敏感，从而阻碍该模型的实际应用价值。近年来，基于改进的水平集方法由于充分考虑综合利用了区域与边界信息，已成为分割灰度分布均匀图像的研究热点。

经对现有技术文献的检索发现，袁克虹等提出基于水平集函数的灰度不均匀图像分割方法(专利号：CN102354396A)；王爽等提出了基于邻域概率密度函数特征的水平集图像分割方法(专利号：CN101571951)；曹宗杰提出了基于概率论模型的水平集方法(专利号：CN101221239)。这些方法均是在变分理论基础上进行改进，提出不同模型的图像分割方法，但是这些模型分割复杂灰度不均匀图像效果不佳且非常耗时。

发明内容

本发明的目的是为了更准确有效地分割灰度不均匀图像，而提供了一种局部引导核拟合能量模型的水平集图像分割方法。本发明提出基于局部引导核拟合能量模型的水平集图像分割方法，通过引入局部核引导拟合能量提取图像局部信息，并在每次迭代过程中进行高斯滤波以避免周期性重新初始化水平集函数，因而能有效分割灰度不均匀图像。

本发明基于局部引导核拟合能量模型的水平集图像分割方法，实现该方法的步骤主要包括：水平集函数的定义，分割模型能量泛函的构建，能量泛函的形式简化，水平集函数的演化和水平集函数的光滑处理。具体步骤如下：

步骤1：水平集函数的定义。水平集函数φ(x,t)表示在t时刻点x到轮廓曲线的最短距离，即符号距离函数(Signed Distance Function，SDF)。定义为符号距离函数：

$\mathrm{\φ}(x,t=0)=\left\{\begin{array}{cc}\mathrm{\ρ}& x\∈R_1\\ 0& x\∈\∂C\\ -\mathrm{\ρ}& x\∈R_2\end{array}\right.$

参数ρ是正常量，R₁和R₂分别表示零水平集函数φ(x,t)＝0内部区域和外部区域，C表示曲线的边界。假设H(□)表示Heavide函数，在实际应采用规范化的Heaviside函数H_ε(x)和Dirac函数δ(x)(正规化参数ε)，其表达式如下：

$H_{\mathrm{\ϵ}}\left(x\right)=\frac{1}{2}[1+\frac{2}{\mathrm{\π}}\arctan \left(\frac{x}{\mathrm{\ϵ}}\right)],\mathrm{\δ}\left(x\right)=H^{\′}\left(x\right)=\frac{d}{\mathrm{dx}}H\left(x\right)$

步骤2：分割模型能量泛函的构建。拟合能量泛函表达式如下：

$E^{\mathrm{LKF}}\left(\mathrm{\φ}\right)=\frac{1}{4}{\∫}_{\mathrm{\Ω}}{\left|\right|\mathrm{\θ}\left(I\left(x\right)\right)-\mathrm{\θ}\left(I^{\mathrm{LKF}}\left(x\right)\right)\left|\right|}^2$

θ(·)是从观察空间到高维特征空间的一个映射，I(x)和I^LKF(x)分别表示原始图像灰度与局部引导核拟合图像灰度，其中I^LKF＝u₁H(φ)+u₂(1-H(φ))，u₁和u₂分别表示曲线C内部R_C(s)和外部区域的拟合函数，u₁和u₂其表达式如下：

$u_1=\frac{K_{\mathrm{\α}}*\left[I\left(x\right)H\left(\mathrm{\φ}\right)\right]}{K_{\mathrm{\α}}*H_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right)},u_2=\frac{K_{\mathrm{\α}}*\left[I\left(x\right)(1-H_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right))\right]}{K_{\mathrm{\α}}*(1-H_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right))}$

其中K_α是的高斯核函数，方差α是尺度参数，α大小的选取与图像的特征有关，如图1所示。

步骤3：能量泛函的形式简化。假设非线性高维数据空间核函数K(y,z)可以用非线性映射函数θ(·)来表示，其表达关系为K(y,z)＝θ(y)^T·θ(z)。能量泛函中的数据项核函数用原始数据空间的非欧氏距离向量来替换，使二维的数据空间数据转化为一维的数据空间，表达式如下：

J_K＝‖θ(I(x))-θ(I^LKF(x))‖²

＝[θ(I(x))-θ(I^LKF(x))]^T·[θ(I(x))-θ(I^LKF(x))]

＝θ(I(x))^T·θ(I(x))^T-θ(I(x))^T·θ(I^LKF(x))-θ(I^LKF(x))^T·θ(I(x))

+θ(I^LKF(x))^T·θ(I^LKF(x))

＝K(I(x),I(x))+K(I^LKF(x),I^LKF(x))-2K(I(x),I^LKF(x))

本发明选取径向基函数(Radial Basis Functio，RBF)中的高斯核函数，RBF核函数的表达式K(y,z)＝exp(-(y-z)²/σ²)，σ为方差。其能量泛函简化形式如下：

$\begin{array}{c}{E}^{\mathrm{LKF}}\left(\mathrm{\φ}\right)=\frac{1}{4}{\∫}_{\mathrm{\Ω}}{\left|\right|\mathrm{\θ}\left(I\left(x\right)\right)-\mathrm{\θ}\left(I^{\mathrm{LKF}}\left(x\right)\right)\left|\right|}^2\mathrm{dx}\\ =\frac{1}{4}{\∫}_{\mathrm{\Ω}}[K(I\left(x\right),I\left(x\right))+K(I^{\mathrm{LKF}}\left(x\right),I^{\mathrm{LKF}}\left(X\right))-2K(I\left(x\right),I^{\mathrm{LKF}}\left(x\right))]\mathrm{dx}\\ =\frac{1}{4}{\∫}_{\mathrm{\Ω}}[2-2K(I\left(x\right),I^{\mathrm{LKF}}\left(x\right))]\mathrm{dx}={\∫}_{\mathrm{\Ω}}[1-K(I\left(x\right),I^{\mathrm{LKF}}\left(x\right))]\mathrm{dx}\\ =\frac{1}{2}{\∫}_{\mathrm{\Ω}}[1-\exp (-\frac{{(I\left(x\right)-I^{\mathrm{LKF}}\left(x\right))}^2}{{\mathrm{\σ}}^2})]\mathrm{dx}\end{array}$

步骤4：水平集函数的演化。采用梯度下降流方法，求取演化曲线的水平集演化方程式。假设关于水平集函数φ的变量η＝δφ，水平集函数变化量φ'＝φ+εη，先固定拟合函数u₁和u₂，通过对参数φ进行微分运算求取最小化能量泛函E^LKF，当自变量ε→0，可得：

$\begin{array}{c}\frac{\mathrm{\δ}{E}^{\mathrm{LKF}}\left(\mathrm{\φ}\right)}{\mathrm{\δ\φ}}=\underset{\mathrm{\ϵ}\→0}{\lim }\frac{1}{2}{\∫}_{\mathrm{\Ω}}[1-\exp (-\frac{{(I\left(x\right)-I^{\mathrm{LKF}}\left(x\right))}^2}{{\mathrm{\σ}}^2})]\mathrm{dx}\\ =\underset{\mathrm{\ϵ}\→0}{\lim }\frac{1}{2}{\∫}_{\mathrm{\Ω}}[1-\exp (-\frac{{(I\left(x\right)-u_1{H}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right)-u_2(1-H_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right)))}^2}{{\mathrm{\σ}}^2})]\mathrm{dx}\\ =\underset{\mathrm{\ϵ}\→0}{\lim }(-{\∫}_{\mathrm{\Ω}}\frac{(u_1-u_2[I\left(x\right)-u_1{H}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right)-u_2(1-H_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right))]}{{\mathrm{\σ}}^2})\\ \exp (-\frac{{[I\left(x\right)-u_1{H}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right)+u_2(1-H_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right))]}^2}{{\mathrm{\σ}}^2}{\mathrm{\δ}}_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right)\mathrm{\ηdx})\\ =-\underset{\mathrm{\ϵ}\→0}{\lim }{\∫}_{\mathrm{\Ω}}\frac{(u_1-u_2)[I\left(x\right)-I^{\mathrm{LKF}}\left(x\right)]}{{\mathrm{\σ}}^2}\\ \exp (-\frac{{(I\left(x\right)-I^{\mathrm{LKF}}\left(x\right))}^2}{{\mathrm{\σ}}^2}{\mathrm{\δ}}_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right)\mathrm{\ηdx}\end{array}$

通过嵌入水平集函数φ(s,t):[0,1]→Ω，根据Euler-Lagrange方程求解关于水平集函数的能量最小值，通过解决如下的偏微分方程式：

$\frac{\∂\mathrm{\φ}}{\∂t}=-\frac{\mathrm{\δ}{E}^{\mathrm{LKF}}\left(\mathrm{\φ}\right)}{\mathrm{\δ\φ}}$

用梯度下降法可以得到曲线演化方程式为：

$\frac{\∂\mathrm{\φ}}{\∂t}=\frac{(u_1-u_2)[I\left(x\right)-I^{\mathrm{LKF}}\left(x\right)]}{{\mathrm{\σ}}^2}\exp (-\frac{{(I\left(x\right)-I^{\mathrm{LKF}}\left(x\right))}^2}{{\mathrm{\σ}}^2}){\mathrm{\δ}}_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right)$

步骤5：水平集函数的光滑处理。水平集函数的符号距离函数作为该曲线的初始化轮廓。因此，我们在水平集函数每次迭代后对符号距离函数进行高斯滤波操作，其表达式如下：

φ_i ⁿ⁺¹＝G_ξ*φⁿ⁺¹

其中G_ξ是高斯核函数，协方差为ξ，协方差ξ应该满足ξ∈[0.45,1]。协方差ξ应满足 表示时间步长。

本发明的优点:本发明通过将局部引导核拟合能量提取图像局部信息，并用高斯函数对水平集函数进行光滑处理，不仅提高了灰度不均匀场景中弱边界目标的分割精度，而且完全避免了周期性的重新初始化问题，降低了算法的计算复杂度，从而提高了分割的精度和分割效率。

附图说明

图1表示本发明实施例中基于局部引导核拟合能量模型的水平集图像分割方法流程图。

图2表示基于高斯核权重的窗口函数。

图3为LKF模型分割医学图像结果效果。

其中：图(a_i)表示分割图像的初始轮廓图像；图(b_i)分割方法在迭代过程中，演化曲线中间结果；图(c_i)曲线最终演化的位置，其中i＝1，2。

图4为单分辨率和多分辨率多区域水平集方法的分割结果。

其中：图(a)和图(e)显示了两组图像的初始轮廓；图(b)和图(f)显示局部二值拟合模型分割该组图像；图(c)和图(g)显示LKF模型的分割结果；图(d)和图(h)显示Chan-Vese模型分割的结果。

具体实施方式

本发明具体实施步骤包括如下：

(1)输入分割图像，设置初始化参数：给定的尺度参数α，时间步长Δt，Heavide函数的正规化参数ε，符号距离函数常量ρ，协方差为ξ；

(2)初始化演化曲线的水平集函数φ，定义为符号距离函数φ(x,t)＝0。

(3)根据描述步骤4中所描述的曲线演化方程式进行计算；

(4)将演化后的水平集函数φ按照步骤5中的高斯滤波方程式进行计算，即φⁿ＝G_ξ*φⁿ；

(5)判断描述的水平集演化曲线是否得到满足终止，如果是，则输出各分割区域的图像和分割结果。否则，将高斯滤波后的水平集函数φⁿ⁺¹＝φⁿ作为下一次迭代的初始水平集函数，转到第三步。

图3显示了LKF模型分割灰度不均匀医学图像结果，试验中尺度参数均为α＝5。图(a₁)半径R＝15像素的圆；图(a₂)为长和宽分别为20和40像素的矩形。该组图像经过迭代后，初始轮廓曲线迅速扩张并逐渐包围图像图标区域，图(b_i)显示了迭代的中间过程。图(c_i)显示了演化曲线最终停止位置,其中i＝1,2。

图4选取了灰度不均匀的MR图像，比较了LKF模型与Chan-Vese模型和LBF模型的分割效果。实验中，Chan-Vese模型的最佳的长度控制参数μ＝0.01×255²，LKF模型中尺度参数均设置为α＝10。在两组图像试验中，每组图像的初始化轮廓相同，图(a)和图(e)显示了两组图像的初始轮廓。图(b)和图(f)显示LBF模型分割该组图像。图(c)和图(g)显示本文提出的LKF模型的分割结果。从分割的结果看，LKF分割结果几乎接近LBF模型，但LKF模型具有更快的收敛速度以及更有效的计算效率，LBF模型和LKF模型迭代次数以及运行时间如表1所示。

表1LBF模型与LKF模型迭代次数与运行时间

Claims (3)

1.一种基于局部引导核拟合能量模型的水平集图像分割方法，其分割模型的能量泛函表达式如下：

$E^{\mathrm{LKF}}\left(\mathrm{\φ}\right)=\frac{1}{4}{\∫}_{\mathrm{\Ω}}{\left|\right|\mathrm{\θ}\left(I\left(x\right)\right)-\mathrm{\θ}\left(I^{\mathrm{LKF}}\left(x\right)\right)\left|\right|}^2$

θ(·)表示从观察空间到高维特征空间的一个映射，I(x)和I^LKF(x)分别表示原始图像灰度与局部引导核拟合图像灰度，其中I^LKF＝u₁H(φ)+u₂(1-H(φ))，u₁和u₂分别表示曲线C内部R_C(s)和外部区域的拟合函数，u₁和u₂其表达式如下：

$u_1=\frac{K_{\mathrm{\α}}*\left[I\left(x\right)H\left(\mathrm{\φ}\right)\right]}{K_{\mathrm{\α}}*H_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right)},u_2=\frac{K_{\mathrm{\α}}*\left[I\left(x\right)(1-H_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right))\right]}{K_{\mathrm{\α}}*(1-H_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right))}$

其中K_α是的高斯核函数，方差α是尺度参数；

其分割方法的具体步骤如下：

步骤1：输入分割图像，定义水平集函数：用符号距离函数来表示水平集函数φ(x,t)，即表示在t时刻点x到轮廓曲线的最短距离；

步骤2：分割模型能量泛函的构建：通过考虑图像的局部拟合信息，来建立图像分割模型的能量泛函；

步骤3：能量泛函的形式简化：通过将能量泛函中的数据项核函数用原始数据空间的非欧氏距离向量来替换，使二维的数据空间数据转化为一维的数据空间，以便能量泛函能用梯度下降法进行求解；

步骤4：水平集函数的演化：采用梯度下降流方法，求取演化曲线的水平集演化方程式；

步骤5：水平集函数的光滑处理：在水平集函数每次迭代后对符号距离函数进行高斯滤波操作，使演化曲线能变光滑。

2.根据权利要求1所述的基于局部引导核拟合能量模型的水平集图像分割方法，其特征在于:水平集函数的演化，采用梯度下降流方法，求取演化曲线的水平集演化方程式；假设关于水平集函数φ的变量η＝δφ，水平集函数变化量φ'＝φ+εη，先固定拟合函数u₁和u₂，通过对参数φ进行微分运算求取最小化能量泛函E^LKF，当自变量ε→0，可得：

$\begin{array}{c}\frac{\mathrm{\δ}{E}^{\mathrm{LKF}}\left(\mathrm{\φ}\right)}{\mathrm{\δ\φ}}=\underset{\mathrm{\ϵ}\→0}{\lim }\frac{1}{2}{\∫}_{\mathrm{\Ω}}[1-\exp \left(\frac{{(I\left(x\right)-I^{\mathrm{LKF}}\left(x\right))}^2}{{\mathrm{\σ}}^2}\right)]\mathrm{dx}\\ =\underset{\mathrm{\ϵ}\→0}{\lim }\frac{1}{2}{\∫}_{\mathrm{\Ω}}[2-\exp (-\frac{{(I\left(x\right)-u_1{H}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right)-u_2(1-H_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right)))}^2}{{\mathrm{\σ}}^2}]\mathrm{dx}\\ =\underset{\mathrm{\ϵ}\→0}{\lim }(-{\∫}_{\mathrm{\Ω}}\frac{(u_1-u_2)[I\left(x\right)-u_1{H}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right)-u_2(1-H_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right))]}{{\mathrm{\σ}}^2}\\ \exp (-\frac{{[I\left(x\right)-u_1{H}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right)+u_2(1-H_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^{\′}\right))]}^2}{{\mathrm{\σ}}^2}{\mathrm{\δ}}_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right)\mathrm{\ηdx})\\ =-\underset{\mathrm{\ϵ}\→0}{\lim }{\∫}_{\mathrm{\Ω}}\frac{(u_1-u_2)[I\left(x\right)-I^{\mathrm{LKF}}\left(x\right)]}{{\mathrm{\σ}}^2}\\ \exp (-\frac{{(I\left(x\right)-I^{\mathrm{LKF}}\left(x\right))}^2}{{\mathrm{\σ}}^2}{\mathrm{\δ}}_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right)\mathrm{\ηdx}\end{array}$

通过嵌入水平集函数φ(s,t):[0,1]→Ω，根据Euler-Lagrange方程求解关于水平集函数的能量最小值，通过解决如下的偏微分方程式：

$\frac{\∂\mathrm{\φ}}{\∂t}=-\frac{\mathrm{\δ}{E}^{\mathrm{LKF}}\left(\mathrm{\φ}\right)}{\mathrm{\δ\φ}}$

用梯度下降法可以得到曲线演化方程式为：

$\frac{\∂\mathrm{\φ}}{\∂t}=\frac{(u_1-u_2)[I\left(x\right)-I^{\mathrm{LKF}}\left(x\right)]}{{\mathrm{\σ}}^2}\exp (-\frac{{(I\left(x\right)-I^{\mathrm{LKF}}\left(x\right))}^2}{{\mathrm{\σ}}^2}){\mathrm{\δ}}_{\mathrm{\ϵ}}\left(\mathrm{\φ}\right).$

3.根据权利要求1所述的基于局部引导核拟合能量模型的水平集图像分割方法，其特征在于：水平集函数的光滑处理，其水平集光滑表达式如下：

${\mathrm{\φ}}_i^{n+1}=G_{\mathrm{\ξ}}*{\mathrm{\φ}}^{n+1}$

其中G_ξ是高斯核函数，协方差为ξ，协方差ξ应该满足ξ∈[0.45,1]；协方差ξ应满足 $\mathrm{\ξ}>\sqrt{\mathrm{\Δt}},$ 表示时间步长。