楊名宇

（中國科學院 長春光學精密機械與物理研究所，中國科學院航空光學成像與測量重點實驗室，吉林 長春 130033）

1 引 言

自Kass等人提出主動輪廓模型［1］以來，基于曲線演化的形變模型已被廣泛地應(yīng)用于圖像分割、跟蹤等領(lǐng)域。其中，在醫(yī)學圖像分割領(lǐng)域的成功應(yīng)用［2－6］，使得該模型備受關(guān)注。該模型的優(yōu)點在于將圖像分割問題轉(zhuǎn)化為求解泛函的極小值問題，并能有效地融合利用圖像自身的低層視覺屬性與對待分割物體已有的先驗知識，因此可較好地獲得分割區(qū)域的完整表示，在一定程度上克服了傳統(tǒng)非模型分割方法由于自身的局限性使得分割區(qū)域的邊界可能不完整，以及缺乏結(jié)合先驗知識能力等缺陷。Malladi、Caselles等人后又提出幾何主動輪廓模型［7－8］，該模型基于曲線演化理論和水平集方法，通過水平集函數(shù)的零水平集間接的表達物體輪廓，該方式雖不如參數(shù)主動輪廓［7］模型直觀，但易于處理曲線的拓撲形變。

2001年Chan和Vese提出了基于水平集的Chan－Vese模型，該模型一旦初始化后便可自動檢測待分割物體的內(nèi)部輪廓，并且曲線的初始化與位置、梯度無關(guān)，因此該模型既可在很大程度上減少噪聲的干擾，又適用于待分割物體邊緣模糊甚至不連續(xù)的情況。但上述模型在實際應(yīng)用中都有共同的缺點，即沒有明確的算法終止條件。目前常用的做法是預(yù)先設(shè)定一個迭代次數(shù)，但這種做法有其缺點：如果迭代次數(shù)設(shè)置過小，則不能得到圖像的完整分割；反之如果迭代次數(shù)設(shè)置過大的話，則不能正確地反映出算法性能，造成計算資源的浪費。因此如何設(shè)置一個恰當?shù)耐Ｖ箺l件對于此類基于能量泛函最小值的分割顯得至關(guān)重要。

本文提出了一種改進的Chan－Vese模型，并給出了該模型相應(yīng)算法的終止條件。實驗表明改進后的模型在新的算法終止條件下，分割結(jié)果正確，并且可以大幅減少原始Chan－Vese模型的迭代次數(shù)，且速度更快。

2 曲線演化理論與水平集方法

曲線演化是解決靜止或運動圖像中分割和目標檢測問題的一種有效方法［8－9］。該方法利用閉合曲線或曲面形變的特定規(guī)律，定義度量閉合曲線或曲面形變的能量函數(shù)，最小化此能量函數(shù)可使閉合曲線逐漸逼近圖像中指定的目標邊緣。該方法將圖像分割、運動檢測轉(zhuǎn)化為求解泛函的極值問題，因而有大量成熟的數(shù)學工具可以直接使用，使得該方法在處理紋理、結(jié)構(gòu)復(fù)雜的情況下有較好的效果。曲線演化在具體求解中常以偏微分方程的形式來表達，其表達方式可分為顯式和隱式兩類。顯式表達直接通過曲線的各種參數(shù)變化來描述曲線的演化過程；隱式表達則通過隱式函數(shù)來描述曲線的演化過程。

水平集方法（Level Set Method）最早由Osher和Sethian［9－10］提出，該方法的基本原理是將運動界面作為零水平集嵌入到高一維的水平集函數(shù)中，由閉超曲面的演化方程可以得到水平集函數(shù)的演化方程，而嵌入的閉超曲面總是其零水平集，所以最終只要確定零水平集即可確定移動界面的演化結(jié)果。經(jīng)過多年的發(fā)展，水平集方法已經(jīng)在圖像分割、拓撲優(yōu)化、閉合界面演化等領(lǐng)域得到廣泛的應(yīng)用［8－10］。

在水平集方法中，二維曲線的演化問題轉(zhuǎn)化為三維空間中水平集函數(shù)曲面演化的零水平集求解問題［10］，這里的水平集指的是由水平集函數(shù)值相等的點組成的集合。設(shè)φ（t）＝φ（x，y，t）表示t時刻的二維閉合曲線，引入水平集函數(shù)后，φ（t）被隱式表示為三維空間上一簇具有相同函數(shù)值的曲線，如圖1所示。

圖1 引入水平集后的二維曲線在三維空間表示Fig.1 2Dcurve displayed in 3Dspace with level set

由圖1可知，當t＝0時三維曲線C（x，y，t＝0）即為二維曲線c（x，y），此時的封閉曲線為零水平集。不斷更新水平集函數(shù)即可實現(xiàn)隱含在水平集函數(shù)中閉合曲線的演化，而最終水平集函數(shù)的零水平集即為曲線的演化結(jié)果。因此，水平集方法其實質(zhì)就是將n維空間的問題映射到n＋1維空間來求解，從而避免了因拓撲結(jié)構(gòu)變化而需做的處理。

在曲線演化過程中，零水平集上的點始終滿足

對式（1）兩端對t求導得：

φ為水平集函數(shù)，Vn為水平集在法線方向上的速度，曲線的內(nèi)外區(qū)域通常用φ≤0和φ＞0來表示。水平集函數(shù)φ通常初始化為符號距離函數(shù)，即

其中：設(shè)Ω表示圖像區(qū)域，Ω0為閉合曲線內(nèi)部區(qū)域，?Ω0為閉合曲線的邊界。d（x，y）表示點x和點y的歐氏距離。

3 Chan－Vese模型

Mumford－Shah模型（以下簡稱 M－S模型）是20世紀80年代提出并被廣泛使用研究的一種圖像處理模型［11］，該模型通過求解一個廣義能量泛函的最小值將圖像分割、噪聲去除和圖像重建三個問題巧妙的融合在一起。Mumford和Shah認為，可以通過最小化下述泛函得到圖像的一個分段光滑的近似：

其中：Ω是圖像定義域，u0為實際觀察到的圖像，C是圖像u0邊緣的連續(xù)封閉曲線，u是u0的分段光滑近似。第一項是長度規(guī)整項，第二項是能量項，第三項是面積平滑項。

通過求解式（3）的最小值，可以一次性地實現(xiàn)圖像的分割、去噪和重建；而且，由于該模型不依賴于圖像的梯度信息，因此對具有弱邊界、甚至不連續(xù)邊界的圖像也能有較好的分割效果。

但由于M－S模型在具體求解中存在較大難度，因此眾多學者分別提出了簡化的M－S模型。Chan和Vese提出將原圖像分成兩個分段光滑的區(qū)域，并用每個區(qū)域的均值來表示該區(qū)域［12］，同時省略了M－S模型中的面積項，模型如下：

結(jié)合水平集方法，設(shè)Ф是根據(jù)初始輪廓線C構(gòu)造的水平集函數(shù)，即｛C｜Φ（x，y）＝0｝，式（4）可寫成：

式（5）關(guān)于Φ的梯度下降方程為：

由于Chan－Vese模型采用了非緊支且光滑的Hε來近似H，故式（6）與下式具有相同的靜態(tài)解：

而式（7）又是式（8）對應(yīng)的梯度下降方程：

由于式（8）是關(guān)于Φ的一次齊次函數(shù)，一般來說不存在最小值，即隨著演化的持續(xù)，水平集函數(shù)Φ將趨于＋∞或－∞。因此，Chan和Vese提出了將水平集函數(shù)Φ約束在［0，1］以求取式（6）的最小值［13－14］的方法。

4 改進的Chan－Vese模型

與上述出發(fā)點不同，新模型不是對Φ的取值范圍加以約束，而是通過對原Chan－Vese模型進行改動，從而得到算法的終止條件。從Chan－Vese模型式（5）中可以看出，由于Heaviside函數(shù)的特性，只有當水平集函數(shù)Φ改變符號時，才會對最小化有影響。因此，當演化曲線附著在物體的外層輪廓時，水平集函數(shù)Φ在下一步?jīng)]有改變符號，此時算法就會終止。因此，本文提出如下改進的Chan－Vese模型：

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

其中：｜·｜2表示2范數(shù)，此處定義矩陣設(shè)時間步長為Δt，可得如下迭代方程：

A＝（αij）的2范數(shù)為

對于式（9）的求解，首先固定Φ，可得c1，c2如下：

再由梯度下降法可得關(guān)于φ的偏微分方程，形式如下：

其中：Hε為正則化的 Heaviside函數(shù)，δε為正則化的Dirac函數(shù)，即Heaviside函數(shù)之導數(shù)，表示分別如下：

使用有限差分的隱式策略求解式（14），（15）。

5 實驗結(jié)果對比

下面的實驗對比Chan－Vese算法和本文方法對同一幅圖像的分割效果。實驗環(huán)境：CPU為Intel T5500，內(nèi)存1G，軟件平臺為 Matlab 7.0，操作系統(tǒng)為 Windwos XP SP3。在以下實驗中，λ1和λ2均取作0.5，ε取作0.01，α取作3。

圖2 Chan－Vese模型與本文模型分割結(jié)果對比Fig.2 Comparison between Chan－Vese model and proposed model

圖2是對一幅合成圖像的分割實驗，其中圖2（a）是原始圖像及初始輪廓。這里，初始輪廓Φ取作以圖像中心為分界，左邊為－1，右邊為1，圖像尺寸為200×200；圖2（b）為Chan－Vese模型的分割結(jié)果，迭代次數(shù)為84次；圖2（c）為新模型的分割的效果。由圖中可以看出，該圖像兩種方法的分割效果類似，但新模型在其停止條件下，僅用了18次迭代便自動停止。

圖3 Chan－Vese模型與本文模型分割結(jié)果對比Fig.3 Comparison between Chan－Vese model and proposed model

圖3是一幅人體腳骨的圖像，其中圖3（a）是原始圖像及初始輪廓。這里，初始輪廓Φ簡單的取作以圖像中心為分界，左邊為－1，右邊為1，圖像尺寸為150×125；圖3（b）為Chan－Vese模型的分割結(jié)果，迭代次數(shù)為184次；圖3（c）為新模型分割的效果?？梢钥闯?，Chan－Vese模型在局部出現(xiàn)了一些錯誤，相比之下，新模型的分割結(jié)果更為準確，且在新的停止條件下，僅用了24次迭代便自動停止。表1給出了Chan－Vese模型與本文方法的迭代次數(shù)與運行時間對比，可以看出，對于圖2和圖3來說，新模型的收斂速度比Chan－Vese模型快3～6倍。

表1 Chan－Vese方法與新方法迭代次數(shù)和運行時間對比Tab.1 Comparison of iteration number and running timebetween Chan－Vese method and proposed method

6 結(jié) 論

針對基于水平集方法的圖像分割沒有統(tǒng)一、明確的算法終止條件，本文提出一種改進的Chan－Vese模型，并在此基礎(chǔ)上明確給出了相應(yīng)算法的終止條件。在新的算法終止條件下，無需手工設(shè)定迭代次數(shù)，擴大了算法的實用性。實驗結(jié)果表明，新模型在新的終止條件下，分割結(jié)果正確，與傳統(tǒng)Chan－Vese模型相比，新模型的手鏈速度比Chan－Vese模型快3～6倍，并且具有較好魯棒性。進一步的研究將集中在水平集函數(shù)初始化方面。

［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］ 王衛(wèi)星，蘇培垠.基于顏色、梯度矢量流活動輪廓及支持向量機實現(xiàn)白細胞的提取和分類［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］ 田豐，夏雪，王鶴.真三維顯示在醫(yī)學教育與仿真中的應(yīng)用［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多分辨率分析分割血管內(nèi)超聲圖像［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，F(xiàn)luid 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.