基于Level Set方法的圖像分割

郭元卡

(西安電子科技大學電子工程學院，陜西西安 710071)

圖像分割是由圖像處理到圖像分析的關鍵步驟，分割的準確性直接影響后續任務的有效性，因此具有重要意義。變分水平集方法［1］是直接通過極小化關于水平集函數自身的能量函數而得到演化的偏微分方程(PDE)，故具有參數化Snake模型的優點，又具有水平集方法的優點。相對于傳統的由純粹PDE驅動的水平集圖像分割方法，基于變分水平集的圖像分割方法可以在能量函數中自然地融入附加約束信息，因而可以產生魯棒性更強的結果。

1 Level Set方法

水平集的基本思想是將曲線看成高一維空間中某一函數φ的零水平集，同時曲線的演化也擴充到高一維的空間中。將水平集函數按照它所滿足的演化方程進行演化或迭代，由于水平集函數不斷進行演化，所以對應的零水平集也在不斷變化，當水平集演化趨于平穩時，演化停止，得到曲線形狀。

考慮零水平集x(t)所對應的水平集函數φ，則有

對方程兩邊求關于時間的偏導數，有

假設F為外法向方向的速度，那么

基本方程式是水平集函數及相應的水平集在法向力F推動下的演化方程。

2 算法原理

傳統的水平集方法在每次迭代后，都需要對水平集函數進行重新初始化，以此修正迭代造成的水平集位置偏離。雖然這樣保證了演化計算的穩定性，但花費時間長，而且重新初始化較復雜，嚴重阻礙了水平集方法的廣泛應用。

Gomes J與 Faugeras O［2］指出，重新初始化水平集函數是水平集理論與實踐的矛盾之處。理論上，重新初始化沒有必要;實踐中，為保證水平集在演化過程中不會偏離符號距離函數(SDF)太遠，時間步長必須選得足夠小，這樣做又增加了迭代次數，增加了演化的時間。下面介紹一種無需重新初始化的變分水平集算法［3］。

(1)內部能量。設φ為水平集函數，SDF滿足一個非常好的性質，即反之，滿足1的任意函數φ都可表示成一個SDF加一個常數［4］。所以引出如下積分

作為函數φ與Ω中SDF接近程度的一個衡量標準。

根據上述定義的函數P(φ)，提出以下變分公式

式中，μ＞0為常數，它控制著“懲罰”水平集φ偏離SDF的力度，P(φ)稱為函數φ的內部能量;基于圖像數據的能量εm(φ)用于驅動φ的零水平集曲線向圖像的邊緣演化，稱為外部能量。

他是最小化函數ε的梯度流［5］。

函數φ根據式(7)表示的梯度流進行演化，在φ的演化過程中，零水平集曲線通過外部能量εm(φ)進行演化;與此同時，由于內部能量“懲罰”的作用，演化函數將自動保持為一個近似的SDF。因此，在提出的公式中重新初始化步驟被完全消除了。

(2)外部能量。

設I為原始圖像，g為邊緣檢測函數，定義如下

式中，Gσ是標準差為σ的高斯核。當輪廓線位于圖像的邊緣處時，圖像的高斯梯度較大，邊緣檢測函數g較小，曲線的演化速度幾乎為零，輪廓線就會停止在圖像的邊緣處。

定義函數φ(x，y)的外部能量如下

式中，λ＞0和v均為常數。Lg(φ)和Ag(φ)定義如下

式中，δ為單變量Dirac函數，H為Heaviside函數。

假設φ的零水平集用可微的參數化曲線C(p)，p∈［0，1］表示，那么Lg(φ)表示的是在正形投影公式中φ的零水平集曲線的長度［6］。

能量函數Ag(φ)用來加速曲線的演化，可以看出，當邊緣檢測函數g為常數1時，Ag(φ)即為Ωφ=區域的加權面積。一般來說，當初始輪廓線位于目標外部時，v為正值，這樣可以促使輪廓線快速收縮;反之，當初始輪廓線位于目標內部時，v為負值，可以加快輪廓線的膨脹速度。

于是，得到完整的能量函數

式中，Δ為Laplacian算子;▽φ是φ的梯度;div為散度。這一梯度流即是零水平集函數的演化方程。第2、3項由外部能量而來，用于驅動φ的零水平集曲線向目標邊緣運動。第1項中的

3 算法實現

(1)偏微分方程的求解。該算法采用有限差分的方法來求解式(13)。在實際應用中，Dirac函數采用正則化的 δε(x)進行計算，δε(x)定義如下

文中所有實驗選取ε=1.5。

(2)時間步長的選擇。在實現水平集方法時，可以選擇的時間步長τ的范圍明顯大于在傳統水平集方法中使用的時間步長范圍，在這個算法中，τ可以設在0.1～100之間。實驗中，為保持水平集演化的穩定性，時間步長τ和系數μ必須滿足

采用較大的時間步長可以加速演化，但是，如果時間步長的選擇過大，在邊緣定位時可能出錯。一般地，對于大多數圖像選擇τ≤10.0。

(3)水平集函數的初始化。在傳統的水平集方法中，必須將水平集函數φ初始化為符號距離函數φ0，如果初始水平集函數明顯不同于SDF，那么重新初始化步驟也不能將這個函數重新初始化為一個SDF。在算法中，不僅重新初始化步驟得以完全消除，而且水平集函數也不再需要初始化為SDF。

這里提出以下函數作為初始函數φ0。令Ω0是圖像區域Ω的一個子集，?Ω0是Ω0的所有邊界點的集合，它們能夠被一些簡單的形態學操作有效識別。于是，初始函數φ0定義如下

式中，ρ＞0是一個常數。ρ一般選取ρ≥2ε，ε為正則化Dirac函數δε(x)定義中的寬度。

與傳統由輪廓計算得來的SDF不同，上述初始水平集函數是從圖像區域Ω中任意一個子區域Ω0中計算而來，水平集函數的這種基于區域的初始化方法不僅計算高效，而且在一些情況下應用較為靈活。例如，如果感興趣區域可以通過一些途徑粗略和自動獲取，那么可以使用這些粗略獲得的區域作為區域Ω0來構造初始水平集函數φ0。接下來，初始水平集函數將根據演化方程穩定演化，同時它的零水平曲線可以收斂于感興趣區域的精確邊緣。

雖然這種初始函數φ0明顯偏離了SDF，但由于能量函數中的懲罰函數P(φ)的調節作用，使水平集函數φ在其零水平集附近能夠近似等于SDF，這對于邊緣的定位沒有太大影響。

(4)水平集函數收斂的判斷。水平集函數收斂的判斷條件為

算法的實現步驟為:

(1)初始化:設置初始輪廓，可以設置為圓形、矩形或是一條非閉合曲線，根據式(18)計算水平集函數的初始值φ0。

(2)利用迭代公式(17)，計算φn+1。

(3)檢查本次迭代是否收斂，如果收斂，則停止計算;否則，轉至步驟(2)，繼續計算。

(4)更新輪廓線，獲得分割圖像。

4 實驗結果與分析

用多幅圖像對上述算法進行了實驗，獲得了較好的分割效果。所有實驗都采用式(18)所定義的初始化方法。

圖1中，參數λ=5.0，μ=0.04，v=1.5，時間步長τ=5給出了一幅像素大小為83×65的顯微鏡下觀測到的兩個細胞圖像的輪廓演化過程。可以看到，兩個細胞的某些邊緣非常模糊，演化曲線自動分裂，將兩個細胞分割出來，實驗結果表明這種方法在處理弱目標邊緣時有較強的魯棒性。實驗中，選擇了一個矩形區域作為區域Ω0來計算初始水平集函數φ0，如圖1(a)中所示。

圖1 細胞圖像(83×65)分割結果

圖2中，參數λ=5.0，μ=0.04，v=3.0，時間步長τ=5給出了一幅像素大小為84×84的兩個目標物圖像的輪廓演化過程。這個時間步長明顯大于在傳統水平集方法中使用的時間步長。在傳統的水平集方法中，初始輪廓線一般選取一條簡單閉合曲線。而圖2(a)中使用了一條非閉合的直線作為初始輪廓線，演化曲線將兩個目標物都提取了出來，這是傳統水平集方法做不到的。

圖2 兩個目標物圖像(84×84)分割結果

圖3與圖2參數相同給出了一幅像素大小為128×128的腦部CT圖像的去骨分割過程。

圖3 腦部CT圖像(128×128)的去骨分割結果

這種方法具有幾個顯著優點:它使得真正影響水平集計算量的重新初始化步驟得以完全消除;靈活的水平集函數初始化方法，使得初始輪廓線的選擇更加自由，并且計算也簡便了很多。但該模型中的初始輪廓線的選取依然受到很大限制，初始輪廓線需包圍所有待分割的物體，而且不能跨越同質區域。

5 結束語

本文提出的無需重新初始化的變分水平集模型，獲得了滿意的分割效果。但由于這個模型僅利用了圖像的邊緣梯度信息，而沒有考慮圖像的全局區域信息，使得模型對邊緣模糊圖像、噪聲圖像、紋理圖像的分割效果不理想。此外，該模型的迭代次數過多，初始輪廓線的選取依然受到很大限制。

［1］ ZHAO H K，CHAN T F，MERRIMAN B，et al.A variational level set approach to multiphase motion［J］.Journal Of Computational Physics，1996，167:179 －195.

［2］ GOMES J，FAUGERAS O.Reconciling distance functions and level sets［J］.J.Visiual Communic and Image Representation，2000，11(2):209－223.

［3］ LI Chunming，XU Chengyang，GUI Changfeng，et al.Level set evolution without reinitialization:a new variational formulation［C］.San Diego:IEEE International Conference on Computer Vision and Pattern Recognition(CVPR)，2005:430－436.

［4］ ARNOLD V I.Geometrical methods in the theory of ordinary differential equations［M］.New York:Springer－Verlag，1983.

［5］ EVANS L C.Partial differential equations［M］.Providence R I:American Mathematical Society，1998.

［6］ VEMURI B，CHEN Y.Joint image registration and segmentation［M］.New York:Geometric Level Set Methods in Imaging，2003.

［7］ OSHER S，FEDKIW R.Level set methods and dynamic implicit surfaces［M］.New York:Spring－Verlag，2002.