基于改進(jìn)能量泛函模型的噪聲圖像分割算法

韓 明，吳朔媚，王敬濤，孟軍英

石家莊學(xué)院 計算機科學(xué)與工程學(xué)院，石家莊 050035

1 引言

圖像分割技術(shù)是計算機視覺領(lǐng)域一項重要的圖像分析技術(shù)，圖像分割技術(shù)的優(yōu)劣直接影響后續(xù)的圖像處理效果[1]。但是由于噪聲、前景與背景對比度低、邊界特征模糊、灰度不均勻等現(xiàn)象的出現(xiàn)，使得圖像分割變得越來越復(fù)雜，因此在圖像分割中如何抑制噪聲，提高分割精度成為主要研究內(nèi)容[2]。

為了克服圖像灰度不均勻以及噪聲對圖像分割帶來的影響，國內(nèi)外學(xué)者提出了很多解決方案[3-4]，在一定程度上實現(xiàn)了對噪聲的抑制，在一定程度上解決了灰度不均勻造成的影響，但是有些方法的抗噪性較差。

隨著基于動態(tài)曲線演化的水平集方法在圖像分割中的成功應(yīng)用，及其在圖像分割中的良好性能，基于水平集的圖像分割越來越具有獨特的優(yōu)勢，它將圖像分割的曲線演化問題轉(zhuǎn)化為了偏微分方程的求解，避免了對演化過程中的跟蹤。文獻(xiàn)[5]提出了局部二值擬合算法（Local Binary Fitting，LBF），該算法利用圖像的局部信息進(jìn)行分割，對于灰度不均勻的圖像分割效果較好，但是對于噪聲圖像，尤其是噪聲較強時處理效果較差，并且采用局部信息，容易導(dǎo)致水平集演化過程中陷入局部極小值。一些學(xué)者在此基礎(chǔ)上提出了全局局部融合圖像擬合（LGIF）模型[6-7]，融合了圖像的局部信息和全局信息，對圖像分割的準(zhǔn)確度而言具有較大的提升，但是模型對噪聲處理能力較差。因此針對噪聲圖像的分割方法相繼提出，例如文獻(xiàn)[8]提出局部相似性系數(shù)（Region-Based Mode via Local Similarity Factor）模型實現(xiàn)了對椒鹽噪聲圖像分割，但是當(dāng)噪聲較強時分割精度低或者出現(xiàn)了失敗的現(xiàn)象。文獻(xiàn)[9]提出了結(jié)合局部灰度差異的噪聲圖像分割算法，該算法利用不同噪聲條件下像素的灰度差異不同，通過增加噪聲修復(fù)函數(shù)提高圖像分割的抗噪性，但是算法復(fù)雜度較高，并且容易出現(xiàn)分割的“孤島”現(xiàn)象。

在圖像預(yù)處理中存在噪聲去除和目標(biāo)圖像邊緣保持之間的矛盾，各向異性擴散方程的應(yīng)用成為解決這一矛盾的合適方法[10]。該方法在圖像去除噪聲和保持目標(biāo)邊緣細(xì)節(jié)特征方面具有一定的優(yōu)勢，但是傳統(tǒng)的各項異性擴散方程擴散函數(shù)在梯度較大的區(qū)域容易隨著時間的變化而產(chǎn)生明顯的“階梯現(xiàn)象”，為了減少“階梯現(xiàn)象”的發(fā)生，并且在去噪的同時更好地保持圖像的邊緣信息，本文提出了基于非凸泛函的能量泛函模型實現(xiàn)圖像去噪。

在圖像去噪的同時為了提高模型的分割速度以及分割精度，本文從能量泛函的最小化入手，對能量泛函模型進(jìn)行改進(jìn)提高分割速度以及分割精度。綜合考慮圖像去噪和分割速度和分割精度等方面的因素，本文提出了基于改進(jìn)的能量泛函模型的噪聲圖像分割算法。該算法在圖像預(yù)處理階段采用非凸泛函模型構(gòu)建自適應(yīng)各向異性擴散方程取代傳統(tǒng)的高斯低通濾波去噪，在去噪的同時保持圖像中邊緣信息的完整性，然后將各向異性擴散方程得到的光滑圖像與水平集模型相結(jié)合，得到改進(jìn)的能量泛函模型，通過對新的能量泛函模型求解最小值，得到水平集演化的前景目標(biāo)的最佳位置。通過對比實驗證明本文算法在分割精度和分割速度上具有較高的性能。

2 自適應(yīng)各項異性擴散方程

2.1 P-M模型

Perona和Malik提出的P-M模型[11]是經(jīng)典的基于各向異性擴散方程的圖像去噪模型，該模型的擴散系數(shù)隨圖像梯度的改變而改變，旨在去除噪聲的同時盡可能多地保留邊界的細(xì)節(jié)信息。P-M模型為：

其中，u是隨時間變化的光滑圖像，div為散度算子，?為梯度算子，為擴散系數(shù)并且為非增函數(shù)，與圖像梯度成反比。

由兩個擴散方程可得：當(dāng)s→∞時，c(s)→0；當(dāng)s=0時，c(s)=1，即在梯度較大的地方擴散較小，在圖像梯度較小的地方，擴散較大。這些性質(zhì)使得這兩種擴散系數(shù)在去除圖像噪聲的同時一定程度上保護圖像邊緣特征[12]。但是P-M模型的平滑效果較差，在梯度較大的區(qū)域容易隨著時間的變化而產(chǎn)生明顯的“階梯現(xiàn)象”，這是由于方程的病態(tài)導(dǎo)致在離散問題中的不穩(wěn)定性。因此在應(yīng)用過程中為了使得圖像的處理結(jié)果更加光滑，需要減弱“階梯現(xiàn)象”的影響。

上述模型的改進(jìn)均是對偏微分方程的求解，都涉及到對初始水平集的選擇與初始化的問題，針對不同的圖像需要選擇合適的水平集，因此為了避免水平集的初始化，本文直接在離散的灰度水平集上計算能量泛函，而不是求解偏微分方程，從而大大提高分割效率。以下證明不存在全局最小值，通過求解能量泛函的極小值得到對應(yīng)的歐拉方程。

2.2 改進(jìn)的能量泛函模型

要得到圖像的光滑模型可以將方程（1）看作是一個求能量泛函極小化的能量耗散過程[13]，能量泛函的形式如下：

圖1 高斯低通濾波與不同灰度水平集進(jìn)行去噪對比

由于噪聲圖像的灰度水平集缺乏正則性，而光滑圖像的灰度水平集則是正則化的，圖1給出了不同圖像利用高斯低通濾波進(jìn)行去噪處理之后的結(jié)果，以及采用不同灰度水平集進(jìn)行去噪之后的對比結(jié)果。從左至右依次為高斯低通濾波，灰度值為60，80，100，120和160的水平集去噪結(jié)果。由圖可以看出經(jīng)過高斯低通濾波之后的圖像，圖像邊緣變得比較模糊，但是灰度水平集函數(shù)處理之后仍舊保留了邊緣信息，灰度水平集選擇灰度值為100時非常接近圖像的邊緣。

本文中采用灰度水平集函數(shù)杜圖像進(jìn)行處理，當(dāng)灰度值選擇適當(dāng)時圖像的邊緣分割準(zhǔn)確，并且保持了邊緣信息，該方法相較于高斯低通濾波在于能夠通過水平集模型使得圖像像素灰度與區(qū)域均值距離函數(shù)相互制約，保持噪聲圖像的邊緣光滑，在曲線演化過程中逼近目標(biāo)輪廓。

由于圖像中的噪聲類型不同，因此圖像的去噪方法也不盡相同，但是通過對比凸泛函模型與非凸泛函模型在去噪中的特點，發(fā)現(xiàn)非凸泛函模型能夠更好地保持圖像的邊界信息[14]。

因此本文提出了基于非凸泛函的能量泛函模型實現(xiàn)圖像的去噪處理。模型如下：

為了保證模型為非凸泛函，因此0＜α＜1，從而保證目標(biāo)邊緣信息被保留，f為噪聲圖像，λ＞0為調(diào)整參數(shù)，起到平衡平滑和保真的作用。

以下證明式（3）不存在全局最小值。

假設(shè)E(u)存在全局最小值點，則有E(u)=0，可得：

由式（4）可見，如果式（4）成立，則需要滿足 u=f和u′=0，而這兩者在圖像空間Ω上幾乎處處都成立，但是又由于 f不可能為常值函數(shù)，因此在E()u的解空間中不存在全局最小值。

通過求解式（3）能量泛函模型的極小值，即歐拉方程可得：

在圖像u的等灰度線的切線方向T=(-uy,ux)和法線方向N=(ux,uy)上的擴散，分別用uT和uN標(biāo)記如下：

由式（5）可得，當(dāng)0＜α＜1時，在圖像u的法方向的擴散系數(shù)為：

由式（8）可得對于式（5）的擴散模型而言總是在法方向N=(ux,uy)上發(fā)生倒向擴散，從而更好地保護圖像中目標(biāo)邊緣信息。

圖像的邊緣信息是區(qū)分圖像不同區(qū)域的最直觀的準(zhǔn)則，作為圖像分割的最重要特征信息，由于非凸泛函可以減小Heaviside函數(shù)在目標(biāo)輪廓處的模糊程度，一方面非凸函數(shù)使得正則項可以近似表示為的加權(quán)l(xiāng)0范數(shù)，從而使得Hε(?)的支撐更加靠近Heaviside函數(shù)的支撐，另一方面采用非凸泛函使得水平集演化的停止函數(shù)在光滑區(qū)域更加接近于1，在圖像的邊緣區(qū)域更加接近于0，因此通過非凸泛函模型使得演化過程中更好地保持圖像的邊界信息。

通過上述證明過程，通過基于非凸泛函的能量泛函模型獲得的擴散模型能夠得到光滑的圖像邊界，并且保持了圖像中目標(biāo)邊緣信息的完整性。

3 CV模型能量泛函實現(xiàn)圖像分割

3.1 CV模型表示與分析

C-V水平集模型[15-16]是基于圖像區(qū)域信息的集合活動輪廓模型，通過最小化能量泛函的方式來演化曲線，設(shè)圖像I(x,y)被閉合曲線C劃分為目標(biāo)和背景兩個同質(zhì)區(qū)域，兩個同質(zhì)區(qū)域的平均灰度值分別為Co和Cb，圖像的擬合能量泛函定義如公式（9）所示：

式中，μ,ν≥0,λ＞0為各個能量項的權(quán)重系數(shù)，均為常數(shù)。L(C)為閉合曲線的長度，A(C)為曲線C的內(nèi)部區(qū)域的面積。 δ(?)是一維的Dirac函數(shù)，H(?)表示Heaviside函數(shù)，?(x,y,t)為水平集函數(shù)。當(dāng)閉合曲線C達(dá)到兩個區(qū)域的邊界時，Ecv(C,Co,Cb)取得最小值，由于CV水平集模型是水平集思想與Mumford-Shah模型的結(jié)合，即C-V水平集模型是Mumford-Shah模型的特例，又因為Mumford-Shah模型的泛函存在極小值，所以C-V水平集模型的能量泛函Ecv(C,Co,Cb)一定存在極小值。即最優(yōu)化式（9）可得最終的目標(biāo)圖像的分割曲線C的位置。

CV模型的示意圖如圖2所示，圖2（a）表示演化曲線包含目標(biāo)輪廓圖；圖2（b）表示目標(biāo)輪廓曲線包含演化曲線；圖2（c）表示目標(biāo)輪廓曲線與演化曲線相交；圖2（d）表示演化曲線與輪廓曲線相重合，只有在圖2（d）所示的情況下能量泛函Ecv(C,Co,Cb)取極小值。

圖2 CV模型示意圖

為了計算式（9）中的數(shù)值解進(jìn)行能量泛函的逼近計算，Dirac函數(shù)δ(x)采用正則化函數(shù)δε(x)，其定義如下：

在實際應(yīng)用中利用數(shù)值計算的方法計算H(x)函數(shù)比較困難，因此在本文中選擇非緊支撐、光滑并嚴(yán)格單調(diào)的正則化形式：

與本文的式（11）進(jìn)行對比實驗的是式（12）的CV模型中的Heaviside函數(shù)：

對比實驗結(jié)果如圖3所示。

圖3 Hε(x)函數(shù)對比

實驗中ε的取值均為0.1。由圖可見當(dāng)ε的取值相等時，本文的非緊支撐、光滑并嚴(yán)格單調(diào)的正則化Hε(x)函數(shù)相對于式（12）中傳統(tǒng)CV模型中選擇的Hε(x)函數(shù)而言能夠更好地逼近Heaviside函數(shù)，因此可見本文的正則化函數(shù)對在分割精度上更高。

Dirac函數(shù)的正則化函數(shù)δε(x)和Heaviside函數(shù)的正則化函數(shù)Hε(x)對于不同的ε取值如圖4所示。

圖4 不同ε值的δε(x)函數(shù)

由圖4可見如果ε取值太?。ㄈ纾害?0.5）則δε(x)的值趨向于0，由圖可見其作用范圍很窄，從而導(dǎo)致能量函數(shù)陷入局部極小化，如果初始輪廓距離目標(biāo)較遠(yuǎn)，容易出現(xiàn)目標(biāo)分割失敗的現(xiàn)象。如果ε取值太大（如：ε=2），由圖可見雖然δε(x)可以得到全局最小化，最終目標(biāo)輪廓可能會不精確[17]。

即當(dāng) ε→0 時，δε(?)→δ(?)，Hε(?)→H(?)，此時能量泛函Ecv(C,Co,Cb)的逼近形式為：

通過求解能量泛函對應(yīng)的Euler-Lagrange方程實現(xiàn)對式（13）的極小化，利用有限差分法得到如下的水平集演化方程，將?表示為梯度流下降的形式為：

通過上述演化方程可得到水平集函數(shù)?(x ,y,t)的演化過程，同時可得演化曲線的演化過程，因此經(jīng)過若干次迭代演化之后初始圖像通過演化曲線被分割，即如圖2（d）所示的分為了曲線內(nèi)部的部分和曲線外部的Ωb=部分。

3.2 模型的灰度水平集算法實現(xiàn)

利用第2章中提出的各項異性擴散模型得到的光滑圖像u對式（13）的能量泛函進(jìn)行改進(jìn)得到新的模型，并且求解能量泛函的最小值，得到水平集演化的前景目標(biāo)的最佳位置。

其中：

由式（17）和式（18）可見改進(jìn)之后的新模型是只依賴于K的函數(shù)。

根據(jù)變分法和Euler-Lagrange方程，求解式（18）對Euler-Lagrange方程表示的水平集演化方程為：

求解式（18）對應(yīng)的能量泛函最小值問題即演化為求解式（19）中目標(biāo)零邊界條件的Euler-Lagrange方程的解：

即

則有，如果L(C)為非常值函數(shù)，由式（24）可見，存在極小值點K∈(Kmin,Kmax)使得 E′cv(K)=0。

本文的改進(jìn)能量泛函模型主要做了以下幾方面的改進(jìn)：

（1）采用基于非凸泛函的能量泛函模型實現(xiàn)圖像的去噪處理，在圖像去噪的同時，利用非凸泛函的圖像邊緣保持特性，保持圖像的邊緣信息，使得水平集演化時更加準(zhǔn)確。

（2）利用式（11）的正則化函數(shù)代替式（12）的正則化函數(shù)使得目標(biāo)函數(shù)更好地逼近Heaviside函數(shù)，從而提高圖像分割的準(zhǔn)確性。

（3）利用通過各項異性擴散模型得到的光滑圖像u代替水平集函數(shù)?，實現(xiàn)對噪聲圖像的非同質(zhì)區(qū)域的邊緣抑制外，能夠更好地進(jìn)行全局優(yōu)化。

（4）在算法的時間復(fù)雜度上本文算法主要在于通過各項異性擴散模型得到的光滑圖像u與水平集函數(shù)?之間的差異，本文算法在能量泛函的迭代求解過程中由于采用了各項異性擴散方程實現(xiàn)光滑圖像的獲取，在每一步的迭代中每個時間步長總是保持減小，直到這種狀態(tài)保持到模型達(dá)到穩(wěn)定為止，而傳統(tǒng)的方法中時間步長保持不變，因此相對而言本文算法的運算復(fù)雜度相對較小，通過第4章實驗中的運行時間的計算同時證明了本文算法運算的復(fù)雜度較小。

3.3 本文算法具體步驟

在上述論述的基礎(chǔ)上，本文基于各向異性擴散的灰度水平集圖像分割算法步驟如下：

步驟1輸入原始圖像I(x,y)利用式（3）的非凸泛函能量泛函模型實現(xiàn)圖像的去噪處理，得到去噪之后的光滑圖像u。

步驟2對步驟1中得到的光滑圖像u給出初始化曲線C，初始參數(shù)設(shè)置如下：

離散網(wǎng)格的間隔為h=1，時間步長為Δt=0.1，正則化參數(shù)ε=1，輸入的原始圖像u中的待分割目標(biāo)較小，長度懲罰項L(C)的權(quán)值μ取較小的值，當(dāng)待分割目標(biāo)較大時則 μ取較大的值，其取值形式可表示為μ=O×2552，其中O∈[ ]0,1。

步驟3根據(jù)式（19）的水平集函數(shù)的演化方程，求出演化的水平集函數(shù)u-K。

步驟4根據(jù)步驟3中的水平集函數(shù)u-K ，提取出零水平集。

步驟5求解式（18）對應(yīng)的能量泛函最小值，判斷演化是否停止，如果停止，則為水平集演化的目標(biāo)最佳邊緣位置，給出分割結(jié)果，算法結(jié)束，否則，轉(zhuǎn)到步驟3繼續(xù)。

4 實驗分析

為了驗證本文算法的有效性，本文選取了不同的視頻圖像序列進(jìn)行去噪和分割實驗，與本文進(jìn)行對比實驗的算法是考慮局部圖像信息的LBF算法和文獻(xiàn)[9]的結(jié)合局部灰度差異的噪聲圖像分割算法。

4.1 實驗環(huán)境和參數(shù)設(shè)置

本文的實驗環(huán)境如下：實驗采用Matlab R2009a作為仿真環(huán)境，實驗計算機配置為：Win10操作系統(tǒng)，Intel Core i7-6700 CPU 3.40 GHz，64 GB內(nèi)存。無特別說明，實驗中參數(shù)設(shè)置如下：離散網(wǎng)格的間隔為h=1，時間步長為Δt=0.1，正則化參數(shù)ε=1，λ=1.0，標(biāo)準(zhǔn)方差隨處理圖像的不同而不同。對算法性能的定量評價，采用迭代次數(shù)和運行時間評價算法的計算效率，采用True Positive Rate（TPR），F(xiàn)alse Positive Rate（FPR）和Similarity Index（SI）[18]進(jìn)行分割精度的度量。

4.2 實驗處理與分析

4.2.1 算法噪聲處理能力實驗

首先測試本文算法對噪聲的處理能力。如圖5所示為三幅具有高強度噪聲的圖像，運用本文算法進(jìn)行噪聲處理和分割，由圖4可見本文算法在高強度噪聲下的目標(biāo)分割取得了較好的效果。

圖5 高強度噪聲下本文算法的分割結(jié)果

4.2.2 不同算法對比實驗

為了實現(xiàn)對算法性能的定量評價，本文采用真陽性比率（True Positive Rate，TPR），假陽性比率（False Positive Rate，F(xiàn)PR），相似性系數(shù)（Similarity Index，SI）三個度量進(jìn)行評估。三個度量的定義如下：

其中，ST表示圖像中待分割圖像的真實前景區(qū)域像素集，SA表示通過模型分割算法所獲取的前景區(qū)域像素集，理想情況下，SI和TPR接近于1，F(xiàn)NR接近于0，表示算法性能越好，TPR值越大，F(xiàn)N R越小說明更多的真實目標(biāo)區(qū)域被分割區(qū)域所覆蓋，SI值越大說明目標(biāo)區(qū)域與分割區(qū)域相似度越高。

另外對于本文算法的分割精確性的度量，通過測量分割之后的目標(biāo)輪廓與真實目標(biāo)輪廓之間的距離進(jìn)行比較，即均方根誤差（Root Mean Square Error，RMSE）進(jìn)行比較[19]。

（1）增加不同噪聲的算法比較

如圖6所示為在本文算法和文獻(xiàn)[9]算法對原始圖像上增加了不同程度的噪聲之后的圖像的分割結(jié)果，其中后兩幅圖像分別增加了標(biāo)準(zhǔn)差為0.01和0.02的高斯噪聲。在實驗中兩個算法的μ=0.06×2552。

圖6 不同噪聲的圖像分割結(jié)果

由圖可以看出對于原始圖像和噪聲污染程度較低的圖像本文算法和文獻(xiàn)[9]算法均能夠?qū)崿F(xiàn)正確的圖像分割，并且SI和TPR的值都接近于1，F(xiàn)PR接近于0，但是對于增加標(biāo)準(zhǔn)差為0.02的噪聲時，由圖可以看出當(dāng)噪聲增強時分割結(jié)果的準(zhǔn)確性不同，本文算法比文獻(xiàn)[9]算法的分割結(jié)果更加精確。此時文獻(xiàn)[9]算法和本文算法SI分別為0.935和0.991，TPR分別為0.982和0.995，F(xiàn)PR分別為0.008和0，通過定性分析可見本文算法抗噪性更強。

圖中標(biāo)準(zhǔn)差0.01和0.02情況下的均方根誤差比較如表1所示。

表1 不同噪聲圖像分割的均方根誤差比較

（2）邊緣模糊圖像對比實驗

該實驗采用腫瘤圖像作為分割圖像，旨在比較不同的算法對前景目標(biāo)灰度和背景相似的圖像的分割結(jié)果的比較，在實驗中μ=0.02×2552。該圖像噪聲強，目標(biāo)圖像的邊緣特征模糊，由圖7可見，LBF算法能夠有效地進(jìn)行分割，但是對于第三列的背景與前景目標(biāo)邊界不明顯的圖像來說分割失敗，對于文獻(xiàn)[9]算法而言由分割結(jié)果可見能夠準(zhǔn)確地定位目標(biāo)的邊界，但是出現(xiàn)了部分“孤島”現(xiàn)象，本文算法對于所實驗的圖像均能夠達(dá)到較好的分割效果。由表2可以看出本文的算法分割精度相對較高，算法分割性能較好，在相同的迭代次數(shù)下用時較短。

圖7 不同算法對腫瘤圖像分割結(jié)果

（3）MR圖像分割實驗

該實驗采用真實的醫(yī)學(xué)超聲圖像進(jìn)行分割實驗，實驗?zāi)康氖欠指畛鰣D像中的白質(zhì)邊界，三幅均帶有較大程度噪聲的醫(yī)學(xué)超聲圖像，并且邊界模糊，在實驗中μ=0.02×2552。由圖可見在噪聲干擾下，LBF算法和文獻(xiàn)[9]算法均沒有能夠成功地提取出白質(zhì)輪廓，文獻(xiàn)[9]算法效果比LBF算法效果較好，但是仍然存在一些不正確的邊界定位，而本文模型則能夠準(zhǔn)確地演化到白質(zhì)邊界，效果較好。三個圖像分割的均方差誤差和重疊區(qū)域誤差比較如表3所示，由表數(shù)據(jù)可見，在圖像中LBF算法效果較差，重疊區(qū)域誤差較大，運行時間相對較長，本文算法相對于其他兩種算法對于噪聲圖像具有更好的分割性能。

表2 圖7中三幅圖像分割結(jié)果的均方誤差、重疊區(qū)域誤差、運行時間比較

圖8 不同算法對MR圖像分割結(jié)果

表3 圖像分割結(jié)果的均方誤差、重疊區(qū)域誤差、運行時間比較

5 結(jié)束語

為了提高水平集對圖像分割的抗噪性，本文在分析了各向異性擴散方程和CV水平集模型的基礎(chǔ)上，提出了一種基于各項異性擴散方程和灰度水平集的噪聲圖像分割算法。與LBF和文獻(xiàn)[9]中的模型相比，本文的算法在強噪聲環(huán)境下具有更好的魯棒性，并且能夠?qū)崿F(xiàn)對強噪聲圖像的高精度分割。本文模型通過引入非凸泛函到能量泛函模型，實現(xiàn)對噪聲圖像平滑處理的同時保持前景目標(biāo)邊緣信息。然后將平滑之后的圖像函數(shù)與水平集能量泛函模型相結(jié)合實現(xiàn)對強噪聲圖像中目標(biāo)邊界的水平集準(zhǔn)確演化，從而保證在強噪聲環(huán)境下圖像分割的準(zhǔn)確性。最后本文通過多組對比實驗進(jìn)行分析和驗證，證明了本文提出的算法對噪聲圖像的分割具有優(yōu)越性。