一種改進的變分水平集的圖像分割算法

張桂梅， 周飛飛， 儲 珺

（南昌航空大學，江西省圖像處理與模式識別重點實驗室，江西 南昌 330063）

一種改進的變分水平集的圖像分割算法

張桂梅， 周飛飛， 儲 珺

（南昌航空大學，江西省圖像處理與模式識別重點實驗室，江西 南昌 330063）

為了解決灰度不均勻現(xiàn)象對醫(yī)學圖像的干擾問題，提出了基于局部極性信息的活動輪廓模型。通過引入局部圖像信息，該模型能有效地分割灰度不均勻圖像。在規(guī)則化項中增加的能量懲罰項，使得水平集函數(shù)在演化過程中保持為近似的符號距離函數(shù)。該算法將圖像分割問題歸結為曲線能量泛函的最小化，首先建立包含局部灰度信息（極性信息）和改進的符號距離函數(shù)的曲線演化能量泛函；然后采用變分水平集方法求解能量函數(shù)的最小值，得到最終的分割結果。真實醫(yī)學圖像和人工合成圖像的實驗結果表明，此方法對灰度不均勻的醫(yī)學圖像有較高的分割精確度，在圖像分割速度上有較大提高。由于利用了局部灰度信息，可以有效地分割灰度不均勻的醫(yī)學圖像，而改進后的變分水平集可以完全避免重新初始化，使得圖像分割效率大大提高了。

變分水平集；圖像分割；符號距離函數(shù)；極性信息

圖像分割是圖像處理領域中基礎的工作，也是進行后續(xù)圖像處理和分析的重要步驟，隨著計算機對圖像處理領域的研究發(fā)展，圖像分割在工業(yè)自動化、遙感衛(wèi)星圖像處理、生物醫(yī)學圖像處理、視頻圖像處理，以及軍事、農業(yè)工程等方面已經(jīng)得到廣泛地應用。

目前基于PDF的圖像分割方法有兩種：基于邊緣的圖像分割方法和基于區(qū)域的分割方法。基于邊緣的分割方法[1]主要是依賴于圖像的梯度信息，與圖像像素在鄰域的灰度變化有關。此方法在圖像梯度變化不明顯或者梯度無意義的弱邊緣圖像中的分割效果并不夠理想。基于區(qū)域的分割方法[2]主要依賴圖像的區(qū)域信息，如圖像的顏色、紋理和平均灰度值等。此方法在無先驗知識的情況下，可以分割較復雜的圖像，但是計算復雜度較高。能量函數(shù)最小化進行圖像分割的方法[3]被提出來了，其主要思路是首先將圖像分割問題看作能量函數(shù)最小化問題，當能量函數(shù)最小時，實現(xiàn)最優(yōu)圖像分割。然后，對給定的能量函數(shù)定義一個圖，讓圖割集的代價正好等于給定的能量函數(shù)。最后，通過求圖的最小割（即圖割），實現(xiàn)給定能量函數(shù)的最小化，也就實現(xiàn)了最優(yōu)圖像分割。Osher和Sethian[4]在1988年提出水平集的圖像分割方法，該方法實現(xiàn)對一個隨時間變化的輪廓曲線求解其能量泛函的偏微分方程的過程。針對這類由曲線的能量泛函最小化所導出的曲線演化問題，Cai等[5]提出了一種新的水平集方法，稱之為變分水平集方法。在數(shù)值實現(xiàn)時變分水平集方法比水平集方法的穩(wěn)定性更高。為了完全避免重新初始化，Zhan等[6]提出了一種改進的變分水平集方法，即在能量泛函中加入一項懲罰項，該方法比變分水平集的分割效率好。Chan和Vese[7]在 Mumford-Shah模型的基礎之上提出了C-V模型，該模型假設圖像的前景和背景灰度是均勻的，通過待分割圖像的統(tǒng)計區(qū)域特征來驅動輪廓曲線的演化。該模型對噪聲具有較好的魯棒性，但是該模型難以從復雜背景中提取目標邊緣，對灰度不均勻的圖像處理效果不夠理想，需要定期初始化曲線輪廓，計算復雜度高。文獻[8]提出了基于最小類內差和最大類間差的圖像分割算法；文獻[9]提出了基于流體靜力學原理活動輪廓模型，但是該模型在灰度較復雜的圖像中分割時間較長；Li等[10]提出了符號距離函數(shù)的約束項，該約束項作為內部能量來約束水平集函數(shù)逼近符號距離函數(shù)，可以解決重新初始化輪廓的計算復雜度問題，提高了計算效率。但是此方法存在局部不穩(wěn)定，從而影響演化的結果。此后，Li等[11]又在此基礎上提出了一種改進的符號距離函數(shù)約束項，稱之為Double-Well型的符號距離函數(shù)，其可以解決局部不穩(wěn)定性問題。DRLSE模型是Li等[11]提出的基于梯度信息和改進的符號距離函數(shù)的活動輪廓模型，DRLSE模型可以有效地解決演化曲線初始化問題，避免了定期初始化的復雜過程，提高了計算的效率。但是該模型利用圖像梯度信息將初始曲線牽引到目標邊界，無法提取以非梯度形式定義的邊界。不能較好地分割灰度不均勻和含有噪聲的圖像。上述方法的曲線演化能量泛函都是基于梯度信息來控制演化的結果，但是梯度信息對紋理和噪聲都比較敏感，從而將部分紋理細節(jié)或噪聲信號誤判為目標邊界。

在圖像分割中采用水平集演化的方式使得零水平集自動收斂于目標物體的邊界。在采用該方法進行分割之前，通常都先對圖像進行濾波，但是濾波方法在對物體進行平滑操作的同時，也會造成物體邊界的弱化，部分消除了物體內部的灰度不一致性，但是不能完全消除不均勻灰度分布，從而影響了分割的效果，而本文方法可以克服這一缺點。

針對以上問題，本文主要研究灰度不均勻的醫(yī)學圖像和人工合成的圖像。利用像素的極性信息[12]來替代上述方法中的梯度信息，由于局部極性信息是衡量像素點附近區(qū)域的梯度矢量方向，增大了局部像素點對輪廓的影響，減少了輪廓距離較遠的像素點對其影響，有效地分辨出像素點是否位于目標的邊界上。而改進后的符號距離函數(shù)可以避免重新初始化水平集的步驟，極性信息作為局部信息對灰度不均勻的圖像處理效果較好。本文將這兩者引入活動輪廓模型中，通過變分法最小化能量函數(shù)得到圖像分割結果。應用真實醫(yī)學圖像和人工合成的圖像進行實驗，結果表明該分割方法有較高的精確性和較低的計算復雜度。

1 相關理論

1.1 水平集方法

活動輪廓模型的基本思想是將圖像分割問題歸結為最小化一個封閉曲線的能量函數(shù)，通過變分水平集的方法求解能量函數(shù)的最小值，此時得到的演化曲線就是分割圖像。Sethian[13]在曲線演化過程做了較多地研究工作，Osher和Sethian[4]在1988年提出了用來描述曲線演化過程的水平集方法，其將二維的演化曲線隱含地看作為更高維的曲面水平集。設演化的曲線集為 C (p,t)，初始的輪廓曲線為C0(p)。C(p, t)表示高維函數(shù) z =φ (x,y,t)的零水平集，即 z =φ (x,y,t)= 0 ，z函數(shù)定義為 z =φ (x ,y,t)。其中，z為點 p (x,y)到曲線上的最短代數(shù)距離，一般函數(shù)符號為負數(shù)時表明該點在演化曲線的內部，函數(shù)符號為正數(shù)時表明該點在演化曲線的外部。而所有在曲線上距離為零的點則組成了演化曲線，即零水平集。用水平集函數(shù)φ表示曲線演化的方程為：

其中，F(xiàn)為演化曲線沿著法線方法的速率。

1.2 改進的符號距離函數(shù)

活動輪廓模型一般要求水平集函數(shù)在演化的過程中保持為符號距離函數(shù)，即|? φ |=1。如果在演化過程中水平集函數(shù)不保持為符號距離函數(shù)，將導致數(shù)值迭代過程不穩(wěn)定，從而得到錯誤的分割。為阻止這種退化，一般在若干次迭代后重新初始化水平集函數(shù)，以將其重新置為符號距離函數(shù)。通過如下的偏微分方程來實現(xiàn)水平集的重新初始化：

其中，φ為水平集函數(shù)，sign(φ)為符號函數(shù)，當曲線演化穩(wěn)定后，即將初始化后的水平集函數(shù)作為新的初始水平集函數(shù)，開始新的演化。但究竟需要經(jīng)過多少次迭代進行重新初始化，目前還沒有準確的結論，故采用此方法在重新初始化上需要花費大量的計算時間，從而導致分割效率低。為避免重新初始化水平集函數(shù)，文獻[14]提出一種可避免重新初始化的變分水平集方法，該方法引入二次函數(shù)型的符號距離函數(shù)約束如下式：

式(3)的函數(shù)表達如圖1所示，從圖1可以看出，約束項是根據(jù)當前水平集函數(shù)偏離符號距離函數(shù)的程度來計算函數(shù)的值。當前水平集函數(shù)偏離符號距離函數(shù)越多，該函數(shù)值越大。

其中，差分率的計算函數(shù)為：

圖1 式(3)的曲線圖

圖2表示符號距離函數(shù)的差分率(difference rate)從圖中可以看出，當|? φ |≥ 1，d(| |)p?φ 為負數(shù)，差分為后向差分，使得|? φ | 增加。前向、后向差分使得符號距離函數(shù)保持在|? φ |=1，但是在|? φ | 接近于0時，差分率d(| |)p?φ 為負的無限大，這使得符號距離函數(shù)在零點附近出現(xiàn)震蕩的現(xiàn)象，從而影響水平集演化的穩(wěn)定性。

圖2 式(5)的曲線圖

針對上面存在的震蕩問題，Li等[11]提出了Double-Well型的符號距離函數(shù)約束項，如下式：

式(6)的函數(shù)表達如圖3所示，類似圖1，也是表明約束項根據(jù)當前水平集函數(shù)偏離符號距離函數(shù)的程度計算函數(shù)的值。當前水平集函數(shù)偏離符號距離函數(shù)越大，該函數(shù)值越大。而且可以看出函數(shù)在0點和1點時取得最小值。同樣采用式(4)差分率計算函數(shù)，計算差分率的式子如下：

圖3 式(6)的曲線圖

從圖4可以看出改進后的符號距離函數(shù)有如下特點：

圖4 式(7)的曲線圖

1.3 極性信息

像素的極性值可以精確地分辨出像素點是分布在紋理或噪聲中還是位于目標的邊界上，以此來分離目標的兩個顯著區(qū)域。相比于梯度信息，極性信息更能準確地區(qū)別出顯著目標的邊界。極性值是一量個方局向部，的設特定征→為量主，方用向于，衡量為x附垂近直區(qū)于域的的梯單度位矢向量。設像素x的極性值為 P (x)，其表達式為：

其中，圖像亮度的梯度 ? I為圖像與高斯濾波器的一階導數(shù)在每一維的卷積。Gσ是寬度為σ的平滑高斯內核高斯窗函數(shù)。和分別代表正極和負極部分參數(shù)的調整，E+表示在窗函數(shù) Gσ( u,v)中有多少梯度矢量在 的正極部分→，E-表示為窗函數(shù)Gσ( u,v)中有多少梯度矢量在η的負極部分。極性值一般在0-1之間，所有尺度中，當極性值在接近于1表示像素點在顯著目標分開的邊界上。

2 本文的圖像分割模型

活動輪廓模型主要用于分割圖像中的目標物體。其基本思想是將圖像分割問題歸結為最小化一個封閉曲線的能量泛函，能量方程由兩部分組成：①以規(guī)范化曲線形狀為目的的項，稱為內能量(internal force)；②以靠近目標物體邊緣為目的的項，稱為外能量(external force)。最小化內能量使得曲線不斷向內部緊縮且保持平滑，外能量是保證曲線緊縮到目標物體邊緣停止。通過變分水平集的方法最小化能量函數(shù)，演化曲線收斂到目標物體邊緣，即分割圖像。

通過以上分析可知，基于改進的符號距離函數(shù)可以完全避免重新初始化輪廓曲線問題，數(shù)值計算穩(wěn)定，計算復雜度低。基于局部灰度信息的幾何活動輪廓模型能較好地處理灰度不均勻的圖像，因此引入極性信息作為模型中的局部灰度信息。本文將兩者的優(yōu)點有效結合起來，基于改進的符號距離函數(shù)作為組合模型能量函數(shù)的內能量，而局部灰度的極性信息引入到外部能量函數(shù)中。本文的能量函數(shù)式如下：

其中，E表示整個區(qū)域的能量函數(shù)， μ,λ,α為大于零的常量， H (φ)對應為Heaviside函數(shù)，δ (φ)對應為 Dirac函數(shù)，φ 對應于水平集函數(shù)。正則化后的Hε(φ)和 δε(φ)的函數(shù)表達式為：

gp為演化曲線停止的條件。其表達式如下：

其中， P(I)表示像素的極性值，通過變分法得到上式的梯度下降流，即水平集函數(shù)的演化方程：

式中右側第1項是曲線演化中保持水平集函數(shù)為符號距離函數(shù)的約束項，第2項和第3項主要控制演化曲線到達目標的分割邊界。

本文算法的主要步驟如下：

步驟1. 輸入待分割的目標圖像，在目標圖像中設置演化曲線的大小和位置；

步驟2. 將水平集函數(shù)初始化為二值函數(shù)；

步驟3. 對式(13)中的參數(shù) μ ,λ,α進行設置；

步驟4. 根據(jù)式(11)、(12)計算 Hε(φ),δε(φ),gp的值；

步驟5. 根據(jù)式(13)更新水平集函數(shù)并判斷水平集函數(shù)是否收斂，如果收斂則輸出分割圖像，否則返回步驟3。

3 實驗與結果分析

3.1 實驗

本實驗是在Intel(R)Core(TM)i5 CPU 4.00 GB內存的環(huán)境下進行，程序通過matlab2012a實現(xiàn)的。為了驗證本文算法的有效性，本文實驗的參數(shù)設置為 μ = 0 .04， λ =5， α = 1 .5；實驗圖像為灰度不均勻的醫(yī)學圖像和高斯噪聲圖像。

實驗1. 本實驗的圖像為4幅灰度不均勻的醫(yī)學圖像，為了驗證本文算法對灰度不均勻醫(yī)學圖像的分割效果，將本文模型與C-V模型[7]、DRLSE[9]模型進行比較。實驗結果如圖5所示，從該實驗結果可以看出C-V模型和DRLSE模型都不能正確地分割圖像，甚至有部分分割錯誤，圖中綠色的圓圈代表錯誤的分割結果。而本文算法能夠精確地分割出灰度不均勻的醫(yī)學圖像。

圖5 C-V模型、DRLSE模型和本文模型對灰度不均勻的醫(yī)學圖像的分割圖

實驗 2. 為了驗證本文方法的分割效率，將選用圖6的人工合成圖像進行實驗，驗證在不同程度噪聲干擾的情況下，本文算法的分割時間。第一列的圖像是均值為零，方差為0.005的高斯噪聲圖。第二列的圖像是均值為零，方差為0.010的高斯噪聲圖。從分割的結果可以看出每個模型都能較好地分割噪聲圖像，但是迭代次數(shù)和處理的 CPU時間不同，從表1可知，對第一行圖像的處理，C-V模型需要迭代200次才能分割圖像，而本文的模型只需迭代 10次就能達到較好地分割結果，大大減少了迭代的次數(shù)。在第二行的噪聲圖像中，C-V模型能較好地處理這種情況，但是因為需要定期初始化輪廓，迭代的時間會稍慢一些。而DRLSE模型相對C-V模型處理的時間快點。而本文模型迭代次數(shù)少，能較快地收斂，計算效率高。相比于C-V模型 和DRLSE模型，本文方法能高效率的分割圖像。

圖6 對含有高斯噪聲圖像的分割結果

表1 本文模型與文獻[7]、[11]模型對噪聲圖像分割的效率

實驗 3. 本實驗主要是驗證本文改進的變分水平集方法與普通水平集方法的分割效率，由于改進的變分水平集不僅完全避免了重新初始化的問題，也使得嵌入函數(shù) φ0( x,y)的工作大大的簡化。而且模型本身具有自動促使嵌入函數(shù)φ近似為距離函數(shù)的能力，故初始化的 φ0( x,y)就不一定要嚴格的讓其是C0的距離函數(shù)。圖 7是對兩組血管圖像分割的結果，從實驗中可以看出，兩種水平的方法都能分割出血管，但使用普通水平集方法需要迭代200次才能完成，CPU的處理時間為5.0 s，如圖7(b)所示，而本文用改進后的變分水平集方法只需要迭代 90次就得到完成，CPU的處理時間為0.8 s，如圖7(c)所示。所以改進后的變分水平集函數(shù)有更高地分割效率。

圖7 血管圖像的分割結果

3.2 性能評價

本文主要從分割正確率和分割時間來評價本文算法分割的精確度和有效性。為了驗證本文算法分割結果的精確度，采用Jaccard相似性比率[15]進行性能分析。Jaccard相似性比率是兩個集合 SA、SB的交集和并集的比率，Jaccard相似性比率的計算公式為：

其中，J(SA,SB)表示正確的分割率，SA表示分割圖像的像素點集合，SB表示真實圖像的像素點集合，利用式(14)對分割圖像的正確率做計算，正確率越高，說明錯誤分割率越小，分割的效果越好。表2給出了圖5中C-V模型，DRLSE模型和本文模型正確分割率的數(shù)據(jù)，從表2中可以得出，本文模型的正確分割率最高，分割效果最好。

表2 本文模型與文獻[7]、[11]模型分割結果正確率(%)

4 結 論

為了解決醫(yī)學圖像灰度分布不均勻及分割圖像時的速度問題，本文結合圖像的極性信息以及DRSLE模型，提出新的分割模型。將極性信息引入活動輪廓模型中，并且融合改進的符號距離函數(shù)，建立了曲線演化的能量泛函，通過變分法最小化能量泛函。由于極性信息對紋理圖像和噪聲圖像比較穩(wěn)定，而改進的符號距離函數(shù)約束項可以避免演化曲線的定期初始化，提高演化的穩(wěn)定性和速度，基于這兩種信息的活動輪廓模型的優(yōu)點組合成新的能量函數(shù)，通過能量最小化來使演化曲線停止于目標邊界。由于本文引入的符號距離函數(shù)可以有效避免輪廓的定期初始化，大大提高了圖像分割的速度，而極性信息能有效地分辨出像素點位于目標的邊界，從而提高了分割的精確度。用真實醫(yī)學圖像和人工合成圖像進行了實驗，結果表明，本文模型能有效地處理灰度不均勻醫(yī)學圖像，同時對不同類型的圖像分割效率也有較大地提高。本文使用的極性信息比梯度信息的分割結果更精確，但還是達不到最精確的效果，在以后的研究工作中，將考慮加入被研究圖像的相關區(qū)域和形狀信息來提高分割精確度。

[1] Xu Chenyang, Prince J L. Snakes, shapes, and gradient vector flow [J]. Image Processing, IEEE Transactions on, 1998, 7(3): 359-369.

[2] Ronfard R. Region-based strategies for active contour models [J]. International Journal of Computer Vision, 1994, 13(2): 229-251.

[3] 劉松濤, 殷福亮. 基于圖割的圖像分割方法及其新進展[J]. 自動化學報, 2012, 38(6): 911-922.

[4] Osher S, Sethian J A. Fronts propagating with curvature-dependent speed: algorithms based on Hamilton-Jacobi formulations [J]. Journal of Computational Physics, 1988, 79(1): 12-49.

[5] Cai Li, Gao Hao, Xie Wenxian. Variational level set method for left ventricle segmentation [C]//TENCON 2013-2013 IEEE Region 10 Conference (31194). Xi′an, China, 2013: 1-4.

[6] Zhan Tianming, Zhang Jun, Xiao Liang, et al. An improved variational level set method for MR image segmentation and bias field correction [J]. Magnetic Resonance Imaging, 2013, 31(3): 439-447.

[7] Chan T F, Vese L A. Active contours without edges [J]. Image Processing, IEEE Transactions on, 2001, 10(2): 266-277.

[8] 吳 昊, 汪榮貴, 方 帥, 等. 基于最小類內差和最大類間差的圖像分割算法研究[J]. 圖學學報, 2011, 32(1): 67-75.

[9] 黃文鈞, 方麗菁. 應用于灰度圖像分割的基于流體靜力學原理活動輪廓模型[J]. 圖學學報, 2014, 35(3): 350-357.

[10] Li Chunming, Xu Chenyang, Gui Changfeng, et al. Level set evolution without re-initialization: a new variational formulation [C]//Computer Vision and Pattern Recognition, CVPR 2005. IEEE Computer Society Conference on. Storrs, USA, 2005: 430-436.

[11] Li Chunming, Xu Chenyang, Gui Changfeng, et al. Distance regularized level set evolution and its application to image segmentation [J]. Image Processing, IEEE Transactions on, 2010, 19(12): 3243-3254.

[12] Ksantini R, Shariat F, Boufama B. An efficient and Fast active contour model for salient object detection [C]// Canadian Conference on Computer and Robot Vision. Kelowna, CA, 2009: 124-131.

[13] Sethian J A. Curvature and the evolution of fronts [J]. Communication of Matllematical Physics, 1985, 101(4): 487-502.

[14] 崔玉玲. 基于改進符號距離函數(shù)的變分水平集圖像分割算法[J]. 模式識別與人工智能, 2013, 26(11): 1033-1040.

[15] Li Chunming, Kao Chiuyen, Gore J C, et al. Minimization of region-scalable fitting energy for image segmentation [J]. Image Processing, IEEE Transactions on, 2008, 17(10): 1940-1949.

An Improved Variational Level Set Used in Image Segmentation Algorithm

Zhang Guimei, Zhou Feifei, Chu Jun
(Nanchang Hangkong University, Key Laboratory of Image Processing and Pattern Recognition, Nanchang Jiangxi 330063, China)

To solve the problem caused by intensity inhomogeneity in medical images, this paper proposed an active contour model based on the local polarity information. By incorporating the local information of image, the proposed model can efficiently segment the image which intensity is nonuniform. Through introducing a penalizing energy term into the regularization, the level set function can be approximated as a signed distance function all the time in the process. This algorithm regarded the image segmentation as the task to minimize the curve of energy functional. Firstly, to establish the evolution curve of energy functional including the local gray level information (polarity information) and the improved signed distance function. Then to solve the minimization of energy function by using the variational level set method to be able to get the final segmentation result. The experiment results of real medical images and artificial synthesis images have shown that our method can get a higher precision of segmentation for the medical images with uneven grayscale. In addition, the image processing speed has been improved greatly. Because of the utilization of local gray information, we can segment the medical image with uneven grayscale effectively, the improved variational level set can completely avoid the re-initialization, and greatly improve the efficiency of image segmentation.

variational level set; image segmentation; symbolic distance function; polarity information

TP 391.41

A

2095-302X(2015)05-0740-07

2015-01-20；定稿日期：2015-03-25

國家自然科學基金資助項目(61462065, 61263046, 61165011)

張桂梅(1970-)，女，江西臨川人，教授，博士。主要研究方向為圖形圖像處理與模式識別。E-mail：guimei.zh@163.com