楊名宇
(中國科學院 長春光學精密機械與物理研究所,中國科學院航空光學成像與測量重點實驗室,吉林 長春 130033)
自Kass等人提出主動輪廓模型[1]以來,基于曲線演化的形變模型已被廣泛地應用于圖像分割、跟蹤等領域。其中,在醫學圖像分割領域的成功應用[2-6],使得該模型備受關注。該模型的優點在于將圖像分割問題轉化為求解泛函的極小值問題,并能有效地融合利用圖像自身的低層視覺屬性與對待分割物體已有的先驗知識,因此可較好地獲得分割區域的完整表示,在一定程度上克服了傳統非模型分割方法由于自身的局限性使得分割區域的邊界可能不完整,以及缺乏結合先驗知識能力等缺陷。Malladi、Caselles等人后又提出幾何主動輪廓模型[7-8],該模型基于曲線演化理論和水平集方法,通過水平集函數的零水平集間接的表達物體輪廓,該方式雖不如參數主動輪廓[7]模型直觀,但易于處理曲線的拓撲形變。
2001年Chan和Vese提出了基于水平集的Chan-Vese模型,該模型一旦初始化后便可自動檢測待分割物體的內部輪廓,并且曲線的初始化與位置、梯度無關,因此該模型既可在很大程度上減少噪聲的干擾,又適用于待分割物體邊緣模糊甚至不連續的情況。但上述模型在實際應用中都有共同的缺點,即沒有明確的算法終止條件。目前常用的做法是預先設定一個迭代次數,但這種做法有其缺點:如果迭代次數設置過小,則不能得到圖像的完整分割;反之如果迭代次數設置過大的話,則不能正確地反映出算法性能,造成計算資源的浪費。因此如何設置一個恰當的停止條件對于此類基于能量泛函最小值的分割顯得至關重要。
本文提出了一種改進的Chan-Vese模型,并給出了該模型相應算法的終止條件。實驗表明改進后的模型在新的算法終止條件下,分割結果正確,并且可以大幅減少原始Chan-Vese模型的迭代次數,且速度更快。
曲線演化是解決靜止或運動圖像中分割和目標檢測問題的一種有效方法[8-9]。該方法利用閉合曲線或曲面形變的特定規律,定義度量閉合曲線或曲面形變的能量函數,最小化此能量函數可使閉合曲線逐漸逼近圖像中指定的目標邊緣。該方法將圖像分割、運動檢測轉化為求解泛函的極值問題,因而有大量成熟的數學工具可以直接使用,使得該方法在處理紋理、結構復雜的情況下有較好的效果。曲線演化在具體求解中常以偏微分方程的形式來表達,其表達方式可分為顯式和隱式兩類。顯式表達直接通過曲線的各種參數變化來描述曲線的演化過程;隱式表達則通過隱式函數來描述曲線的演化過程。
水平集方法(Level Set Method)最早由Osher和Sethian[9-10]提出,該方法的基本原理是將運動界面作為零水平集嵌入到高一維的水平集函數中,由閉超曲面的演化方程可以得到水平集函數的演化方程,而嵌入的閉超曲面總是其零水平集,所以最終只要確定零水平集即可確定移動界面的演化結果。經過多年的發展,水平集方法已經在圖像分割、拓撲優化、閉合界面演化等領域得到廣泛的應用[8-10]。
在水平集方法中,二維曲線的演化問題轉化為三維空間中水平集函數曲面演化的零水平集求解問題[10],這里的水平集指的是由水平集函數值相等的點組成的集合。設φ(t)=φ(x,y,t)表示t時刻的二維閉合曲線,引入水平集函數后,φ(t)被隱式表示為三維空間上一簇具有相同函數值的曲線,如圖1所示。

圖1 引入水平集后的二維曲線在三維空間表示Fig.1 2Dcurve displayed in 3Dspace with level set
由圖1可知,當t=0時三維曲線C(x,y,t=0)即為二維曲線c(x,y),此時的封閉曲線為零水平集。不斷更新水平集函數即可實現隱含在水平集函數中閉合曲線的演化,而最終水平集函數的零水平集即為曲線的演化結果。因此,水平集方法其實質就是將n維空間的問題映射到n+1維空間來求解,從而避免了因拓撲結構變化而需做的處理。
在曲線演化過程中,零水平集上的點始終滿足

對式(1)兩端對t求導得:

φ為水平集函數,Vn為水平集在法線方向上的速度,曲線的內外區域通常用φ≤0和φ>0來表示。水平集函數φ通常初始化為符號距離函數,即

其中:設Ω表示圖像區域,Ω0為閉合曲線內部區域,?Ω0為閉合曲線的邊界。d(x,y)表示點x和點y的歐氏距離。
Mumford-Shah模型(以下簡稱 M-S模型)是20世紀80年代提出并被廣泛使用研究的一種圖像處理模型[11],該模型通過求解一個廣義能量泛函的最小值將圖像分割、噪聲去除和圖像重建三個問題巧妙的融合在一起。Mumford和Shah認為,可以通過最小化下述泛函得到圖像的一個分段光滑的近似:

其中:Ω是圖像定義域,u0為實際觀察到的圖像,C是圖像u0邊緣的連續封閉曲線,u是u0的分段光滑近似。第一項是長度規整項,第二項是能量項,第三項是面積平滑項。
通過求解式(3)的最小值,可以一次性地實現圖像的分割、去噪和重建;而且,由于該模型不依賴于圖像的梯度信息,因此對具有弱邊界、甚至不連續邊界的圖像也能有較好的分割效果。
但由于M-S模型在具體求解中存在較大難度,因此眾多學者分別提出了簡化的M-S模型。Chan和Vese提出將原圖像分成兩個分段光滑的區域,并用每個區域的均值來表示該區域[12],同時省略了M-S模型中的面積項,模型如下:

結合水平集方法,設Ф是根據初始輪廓線C構造的水平集函數,即{C|Φ(x,y)=0},式(4)可寫成:

式(5)關于Φ的梯度下降方程為:

由于Chan-Vese模型采用了非緊支且光滑的Hε來近似H,故式(6)與下式具有相同的靜態解:

而式(7)又是式(8)對應的梯度下降方程:

由于式(8)是關于Φ的一次齊次函數,一般來說不存在最小值,即隨著演化的持續,水平集函數Φ將趨于+∞或-∞。因此,Chan和Vese提出了將水平集函數Φ約束在[0,1]以求取式(6)的最小值[13-14]的方法。
與上述出發點不同,新模型不是對Φ的取值范圍加以約束,而是通過對原Chan-Vese模型進行改動,從而得到算法的終止條件。從Chan-Vese模型式(5)中可以看出,由于Heaviside函數的特性,只有當水平集函數Φ改變符號時,才會對最小化有影響。因此,當演化曲線附著在物體的外層輪廓時,水平集函數Φ在下一步沒有改變符號,此時算法就會終止。因此,本文提出如下改進的Chan-Vese模型:

通過在原Chan-Vese模型中加入(ε+φ)與(ε-φ),可以明確地反映出在每一次迭代中能量項的變化,同時,還可以加快模型的收斂速度;其中ε為一小正數,作用在于控制Φ的大小。Hε(α+Φ)和Hε(α-Φ)為 Heaviside函數的一個變形,即在原始Chan-Vese模型的Hε(Φ)和Hε(-Φ)分別加上正常數α。由于Heaviside函數的二值性,所以加入正常數α可以使Hε(α+Φ)中的Φ當Φ>-α時,隨著曲線的演化逐漸趨近于-α;同理,可以使Hε(α-Φ)中的Φ當Φ<α時,隨著曲線的演化逐漸趨近于α。故改進后的能量項將水平集函數Φ的取值范圍約束在[-α,α],隨著曲線的演化,Φ將趨向于-α或α,即|Φ|將趨近于α。這樣,算法的終止條件(Stop Criterion)可以通過下式來確定:

其中:|·|2表示2范數,此處定義矩陣設時間步長為Δt,可得如下迭代方程:

A=(αij)的2范數為

對于式(9)的求解,首先固定Φ,可得c1,c2如下:

再由梯度下降法可得關于φ的偏微分方程,形式如下:

其中:Hε為正則化的 Heaviside函數,δε為正則化的Dirac函數,即Heaviside函數之導數,表示分別如下:

使用有限差分的隱式策略求解式(14),(15)。
下面的實驗對比Chan-Vese算法和本文方法對同一幅圖像的分割效果。實驗環境:CPU為Intel T5500,內存1G,軟件平臺為 Matlab 7.0,操作系統為 Windwos XP SP3。在以下實驗中,λ1和λ2均取作0.5,ε取作0.01,α取作3。

圖2 Chan-Vese模型與本文模型分割結果對比Fig.2 Comparison between Chan-Vese model and proposed model
圖2是對一幅合成圖像的分割實驗,其中圖2(a)是原始圖像及初始輪廓。這里,初始輪廓Φ取作以圖像中心為分界,左邊為-1,右邊為1,圖像尺寸為200×200;圖2(b)為Chan-Vese模型的分割結果,迭代次數為84次;圖2(c)為新模型的分割的效果。由圖中可以看出,該圖像兩種方法的分割效果類似,但新模型在其停止條件下,僅用了18次迭代便自動停止。

圖3 Chan-Vese模型與本文模型分割結果對比Fig.3 Comparison between Chan-Vese model and proposed model
圖3是一幅人體腳骨的圖像,其中圖3(a)是原始圖像及初始輪廓。這里,初始輪廓Φ簡單的取作以圖像中心為分界,左邊為-1,右邊為1,圖像尺寸為150×125;圖3(b)為Chan-Vese模型的分割結果,迭代次數為184次;圖3(c)為新模型分割的效果。可以看出,Chan-Vese模型在局部出現了一些錯誤,相比之下,新模型的分割結果更為準確,且在新的停止條件下,僅用了24次迭代便自動停止。表1給出了Chan-Vese模型與本文方法的迭代次數與運行時間對比,可以看出,對于圖2和圖3來說,新模型的收斂速度比Chan-Vese模型快3~6倍。

表1 Chan-Vese方法與新方法迭代次數和運行時間對比Tab.1 Comparison of iteration number and running timebetween Chan-Vese method and proposed method
針對基于水平集方法的圖像分割沒有統一、明確的算法終止條件,本文提出一種改進的Chan-Vese模型,并在此基礎上明確給出了相應算法的終止條件。在新的算法終止條件下,無需手工設定迭代次數,擴大了算法的實用性。實驗結果表明,新模型在新的終止條件下,分割結果正確,與傳統Chan-Vese模型相比,新模型的手鏈速度比Chan-Vese模型快3~6倍,并且具有較好魯棒性。進一步的研究將集中在水平集函數初始化方面。
[1] Kass M,Witkin A,Terzopoulos D.Snakes:Active contour models[J].International Journal of Computer Vision,1988,1(4):321-332.
[2] 王醒策,張美霞,武仲科,等.基于全局LBF水平集模型的腦血管層次粗分割[J].光學 精密工程,2013,21(12):3283-3297.Wang X C,Zhang M X,Wu Z K.Level coarse brain vessel segmentation based on global LBF model[J].Opt.Precision Eng.,2013,21(12):3283-3297.(in Chinese)
[3] 王衛星,蘇培垠.基于顏色、梯度矢量流活動輪廓及支持向量機實現白細胞的提取和分類[J].光學 精密工程,2012,20(12):2781-2790.Wang W X,Su P Y.Blood cell image segmentation on color and GVF snake for Leukocyte classification on SVM[J].Opt.Precision Eng.,2012,20(12):2781-2790.(in Chinese)
[4] Chunming L,Chiu-Yen K,Gore C,et al.Minimization of region-scalable fitting energy for image segmentation[J].IEEE Trans.Image Processing,2008,17(10):1940-1949.
[5] 田豐,夏雪,王鶴.真三維顯示在醫學教育與仿真中的應用[J].液晶與顯示,2012,27(4):535-538.Tian F,Xia X,Wang H.Applications of volumetric three-dimensional display in medical simulation and education[J].Chinese Journal of Liquid Crystals and Displays,2012,27(4):535-538.(in Chinese)
[6] 張麟,汪源源,王威琪,等.活動輪廓模型和Contourlet多分辨率分析分割血管內超聲圖像[J].光學 精密工程,2008,16(11):2303-2311.Zhang L,Wang Y Y,Wang W Q,et al.Intravascular ultrasound image segmentation based on active contour model and contourlet multiresolution analysis[J].Opt.Precision Eng.,2008,16(11):2303-2311.(in Chinese)
[7] Malladi R,Sethian J A,Vernuri B C.Shape modeling with front propagation:A level set approach [J].IEEE Trans.Pattern Analysis and Machine Intelligence,1995,17(2):158-175.
[8] Caselles V,Kimmel R,Sapiro G.Geodeisic active contours[J].International Journal of Computer Vision,1997,22(1):61-79.
[9] Sethian J A.Level Set Methods and Fast Marching Methods:Evolving Interfaces in Computational Geometry,Fluid Mechanics,Computer Vision,and Materials Science [M].London:Cambridge University Press,1999.
[10] Osher S,Sethian J A.Fronts propagating with curvature dependent speed:algorithms based on Hamilton-Jacobi formulation[J].Journal of Computer Physics,1988,79(1):12-49.
[11] Mumford D,Shah J.Optimal approximation by piecewise smooth functions and associated variational problems[J].Comm.Pure Appl.Math.,1989,42(5):577-685.
[12] Chan F,Vese L.Active contours without edges[J].IEEE Trans.Image Processing,2001,10(2):266-177.
[13] Chan F,Esedoglu S,Nikolova M.Algorithms for finding global minimizers of image segmentation and denoising models[J].SIAM J.Appl.Math.,2006,66(5):1632-1648.
[14] Bresson X,Esedoglu S,Vandergheynst P,et al.Global minimizer of the active contour/snake model[J].J.Math.Image Vis.,2007,28(1):151-167.