改進分水嶺與K-means++結合的人體脊柱圖像分割方法

郎成洪，華云松，張嘉棋，鐘雪蓮，王麒翔

（上海理工大學光電信息與計算機工程學院，上海 200093）

0 引言

人體脊柱圖像的準確分割對于脊柱微創手術十分重要，圖像分割效果會直接影響其三維重建圖像效果以及微創手術前病灶信息的獲取。由于不同成像設備和環境條件的影響，使得CT圖像常常包含噪聲。為提高圖像分割效果，人們提出了許多圖像分割方法，主要分為基于閾值、邊緣、區域、聚類分析、能量泛函和人工神經網絡的分割方法。

基于閾值的分割方法是根據一個或幾個閾值對圖像中的灰度值進行分類，常用的有自適應閾值法、最大熵法、模糊閾值法、類間閾值法等。此類方法計算簡單、高效、速度快，但對噪聲比較敏感，且對灰度差異不明顯的目標分割效果不佳。如李長春等使用基于多閾值分割與B樣條插值的方法去除CT圖像中的金屬偽影。

基于邊緣的分割方法通過檢測不同區域邊緣來分割不同目標，邊緣檢測效果好且速度快，但無法得到好的區域結構，且抗噪性差。常見的邊緣檢測算子有Robert算子、Sobel算子、Prewitt算子、Canny算子等。如李慶忠等提出基于改進Canny算子的圖像弱邊緣檢測算法，該算法通過改進梯度倒數加權濾波的方法檢測平滑塊中細小的弱邊緣。

基于聚類分析的圖像分割方法根據相似性準則，將圖像像素在特征空間進行分類聚集并進行分割，然后將分割結果映射回原圖像空間。此類方法主要包括K均值和模糊C均值聚類算法，但其主要缺點是灰度值不均勻、對噪聲比較敏感，且運算開銷大。如姜文斌等提出基于自適應權重法的K-means模型對遙感圖像進行分割的方法，該方法利用二維高斯函數對遙感圖像進行濾波平滑處理，然后根據收斂度和自適應值找到最優權重作為初始聚類中心，以提高圖像分割精度。

基于能量泛函的分割方法是通過連續曲線表示目標邊緣，然后定義一個能量泛函使其自變量包括邊緣曲線，從而進行圖像分割，最具代表性的是Osher提出的水平集模型。如黃國鵬等利用偏移場校正和局部區域信息擬合來構建能量函數，并通過加性算子分裂策略提高圖像分割精度和效率。

基于人工神經網絡的分割方法是通過訓練多層感知機得到線性決策函數，然后通過決策函數對圖像像素進行分類，從而分割出目標的過程。其優點是有效解決了圖像灰度不均勻和噪聲問題，但需要大量數據集，且速度很慢。如唐思源等提出一種基于改進粒子群優化算法的BP神經網絡醫學圖像分割方法，提高了圖像分割精度。

基于區域的分割方法主要采用分裂與合并兩種方法，把具有相似性的像素連通，構成最終的分割區域。其優點是區域特征較好，且能有效克服圖像分割空間小的缺點，但容易出現過分割現象。其中分水嶺變換分割是最強大、應用最廣泛的基于區域的圖像分割方法之一。針對其過分割的缺點，人們提出一些改進的分水嶺分割方法。如周小明采用最大熵算法對梯度圖像進行校正，然后進行圖像分割，取得了良好的圖像分割效果。

隨著機器學習算法的快速發展，越來越多機器學習算法被應用于解決分水嶺算法的過分割問題。機器學習算法是一種監督算法，需要進行大量訓練和相應標記才能對測試圖像進行準確分割。但在大多數情況下，訓練數據和測試數據的數量不足，無法保證訓練圖像的質量和相似度，故不需要標記數據的無監督方法在數據短缺情況下表現得更好。無監督的聚類方法不需要訓練及標記數據，可自動將相似數據點分組到同一聚類中。模糊Cmeans聚類和K-means聚類算法是最著名的無監督聚類算法，已應用于各種領域的數據聚類。本文采用一種比K-means聚類算法更優秀的聚類算法——K-means++聚類算法對人體脊柱圖像進行去噪，并修改區域最小值。相較于K-means，該算法分割效率更高。為了得到更好的分割效果，采用在傳統分水嶺算法基礎上進行改進的標記符控制的分水嶺算法對K-means++聚類算法處理后的初始聚類圖像進行分割。

1 K-means++聚類算法

在原始人體脊柱CT圖像中往往存在灰度相似的區域，這些區域會產生許多不必要的局部極小值，如圖1所示。

Fig.1 CT image of human spine圖1 人體脊柱CT圖像

如果先用K-means對圖像進行初始聚類處理，可減少這種不必要的局部極小值，從而改善分水嶺算法的過分割現象。K-means++算法為解決K-means算法中聚類中心選取過于隨意，且結果對初始聚類中心過于敏感的問題，對初始聚類中心選取方式進行改進，具體步驟如下：

（1）從數據集中隨機選取一個樣本作為初始聚類中心

c

。（2）計算數據集中每個樣本

x

與已選取聚類中心的距離，用

D

(

x

)表示。

（3）計算每個樣本被選為下一個聚類中心的概率，與當前已有聚類中心距離較遠的樣本有更大概率被選為下一個聚類中心。具體概率計算公式如下：

其中，

x

表示樣本點，

D

(

x

)表示樣本點到聚類中心的距離。（4）重復步驟（2）、（3），直到選擇出

K

個聚類中心。（5）計算每個樣本分別到

K

個聚類中心的距離，將各樣本劃分到與其距離最近的聚類中心所屬類別中。

（6）針對每個類別計算其所有樣本特征均值，將該均值作為該類別新的聚類中心。其均值計算公式如下：

其中，

i

表示該類別中的樣本點數，

x

表示樣本點。

（7）重復步驟（5）、（6），直到聚類中心位置不再發生變化。

原始人體脊柱圖像的像素經過K-means++處理后分成了

K

類，由于在醫學圖像中，脊柱區域是聚類中心值最大且面積最大的連通區域，故可將聚類中心值最大且面積最大的連通分量

S

(

i

=1，2，…，

k

)作為脊柱圖像區域。

其中，argmax表示括號內的最大值，即為脊柱區域面積值。

2 形態學處理與標記符控制的分水嶺算法

分水嶺變換是一種基于形態學的圖像分割工具，其基本思想是把圖像看作測地學上的拓撲地貌，把每個像素點灰度值表示為該點的海拔高度，然后模擬泛洪過程。假設在圖像的局部極小值面刺穿一個小孔，并讓水以均勻的速率從小孔中涌出，隨著水深逐步增加，為防止不同集水盆地匯聚到一起，在兩個集水盆地間修建大壩，即形成分水嶺。

2.1 形態學技術去除噪聲

利用K-means++聚類算法初始分割后的人體脊柱圖像雖然過分割率下降了很多，但是仍然存在過分割現象，且脊柱區域連通分量會參雜部分背景像素，故對Kmeans++聚類算法初始分割后的脊柱圖像進行基于重建的開操作和閉操作，在保留圖像重要區域輪廓的同時去除細節和噪聲。

2.1.1 開操作與閉操作

形態學處理有兩種基本運算，分別為腐蝕和膨脹，開操作和閉操作就是結合腐蝕和膨脹兩種基本形態學運算而形成的。開操作是先腐蝕再膨脹，閉操作是先膨脹再腐蝕。

形態學閉操作如下：

其中，

A

表示圖像，

B

表示結構元。

形態學開操作如下：

其中，

A

表示圖像，

B

表示結構元。

膨脹可填充圖像中的小孔洞和圖像邊緣處的一些小凹陷，對圖像的外部進行濾波；腐蝕可去除圖像中一些小的“突刺”，對圖像內部進行濾波，并縮小圖像。開操作和閉操作是一種去噪處理，開操作可濾掉小的“突刺”，閉操作可填充孔洞缺口。

2.1.2 基于重建的開操作與閉操作

采用基于重建的開操作和閉操作清理使用Kmeans++聚類算法處理后的人體脊柱圖像，先對圖像采用基于重建的開操作，在開操作結果上利用基于重建的閉操作進行濾波重建。在簡化初始分割圖像的同時，準確定位輪廓分水線，從根源上消除過分割現象。在分割前需要選一個結構元素，用來對圖像進行形態學重建。利用該結構元素處理圖像后，圖像中的灰度躍變急劇增加，可消除梯度對邊緣方向的依賴性。

基于重建的開操作：

基于重建的閉操作：

2.2 標記符控制的分水嶺算法

在標記符控制的分水嶺算法中，首先將使用Kmeans++聚類算法初始分割后的人體脊柱圖像轉換成灰度圖像，應用形態學算子對圖像進行去噪處理，確定前景標記并與原始圖像疊加，然后確定背景標記，獲得原始圖像的梯度圖像并對其進行修改，最后應用分水嶺變換生成分割圖像。基于標記符控制的分水嶺算法流程如圖2所示。

Fig.2 Flow of watershed algorithm based on control marker圖2 基于控制標記符的分水嶺算法流程

具體步驟如下：

（1）由于脊柱骨骼區域在人體脊柱CT圖中一般是灰度值較高的連通像素塊，將形態學處理后的圖像記為

f

。（2）計算

f

的局部區域最大值，并將局部區域最大值作為前景標記，然后把前景標記的灰度圖像轉換成二值圖像，記為

f

。

（3）將前景圖像疊加到原始圖像上，某些對象中的前景標記會一直覆蓋到對象邊緣。

（4）通過對前景圖像執行先閉操作后侵蝕的形態學去噪處理，清理標記斑點的邊界，使得前景標記向內收縮。

（5）使用閾值分割方法對圖像

f

進行分割，并將分割結果作為背景標記圖像，記為

f

，以識別背景中的暗像素。（6）通過對背景標記圖

f

進行距離變換，使背景標記不會非常靠近想要分割的物體邊界，然后使用經過距離變換后的分水嶺算法產生圖像的分水嶺脊線。對圖像中所有像素(

i

，

j

)求距離變換公式如下：

（7）使用梯度幅值作為分割函數，并使用Sobel算子計算圖像梯度。梯度值總是在物體邊緣處高，而在物體內部低。

其中，

I

、

I

分別表示圖像在垂直和水平方向的梯度圖像，?I表示圖像梯度。

（8）梯度圖像被修改，使得在前景和背景標記處僅存在區域最小值。經過修改后的梯度圖像表示為：

其中，

IMMIN

(*)表示形態學極小值標定操作，?I表示圖像梯度，

f

為前景標記的二值化圖。（9）將修改的梯度圖像?I′作為分水嶺算法的輸入，計算得到理想的物體分割結果，其過程可表示為：

其中，

WST

(*)表示為分水嶺變換操作，?I表示圖像梯度。

3 實驗數據及結果分析

3.1 實驗數據

本文實驗數據采用合作醫院的雙源64排CT機對患者進行掃描獲取的人體脊柱CT圖像，并通過DICOM Converter軟件將人體脊柱CT圖像的DICOM格式數據轉換成JPG格式的原始實驗圖像，二維圖像分辨率大小為512×512。本文實驗數據集中共有186張不同部位的人體脊柱圖像，部分數據集如圖3所示。

Fig.3 Experimental data set of partial human spine images圖3 部分人體脊柱圖像實驗數據集

3.2 實驗結果分析

本文實驗在PyCharm Community Edition 2020.3.5x64和Matlab R2018a軟件平臺上進行，首先直接應用分水嶺變換分割方法對人體脊柱圖像進行分割，并查看結果。圖4顯示了直接對人體脊柱原始圖像進行分水嶺分割的結果。如圖所示，直接應用分水嶺變換方法的分割結果顯示出很高的過分割率，并導致錯誤的結果。為提高分割效果，采用K-means++聚類方法進行初始聚類。選用不同聚類數

K

對原始圖像進行聚類，不同

K

值的K-means++聚類效果如圖5所示。

Fig.4 Effect drawing of direct watershed圖4 直接分水嶺分隔效果

Fig.5 Effect of K-means++clustering with different K values圖5 不同K值的K-means++聚類效果

由先驗知識可知，在醫學圖像中，骨骼區域一般是完整且灰度值較高的連通像素塊。通過分析以上不同

K

值的初始聚類圖像可知，當

K=

3時，能得到原始圖像的最佳聚類效果，在分割出脊柱骨的同時，盡量減少局部極小值點。對使用

K=3

的K-means++聚類算法處理后的初始聚類圖像進行分水嶺變換，可發現其錯誤的局部極小值明顯較少，過分割率明顯降低，如圖6所示。

為進一步減少分水嶺變換的過分割率，對K-means++算法聚類后的初始分割圖像（見圖5（b））應用5×5的圓形結構體進行腐蝕，然后對腐蝕圖像進行形態學重建，再使用相同的結構體對重建后的脊柱圖像進行膨脹處理，最后對膨脹后的圖像進行形態學重建處理，使得K-means++處理后的人體脊柱圖像變得平滑，如圖7所示。

Fig.6 Effect of direct watershed transformation after K-means++processing圖6 K-means++處理后的直接分水嶺變換效果

Fig.7 Result after morphological denoising圖7 形態學去噪處理后結果

人體脊柱醫學圖像中灰度值較大的區域一般為脊柱區域，故提取經過形態學處理后的結果圖像（見圖6）的區域最大值，得到前景標記。標記的前景圖如圖8所示。

Fig.8 Tagged foreground圖8 標記的前景圖

上述經過形態學處理后的結果圖像（見圖7）中比較亮的是人體脊柱骨區域，比較暗的是背景，故采用閾值分割法對形態學處理后的結果圖像進行分割，將分割結果標記為背景圖像，如圖9所示。

Fig.9 Background image of tag圖9 標記的背景圖

為了不讓背景標記太接近目標物體，對背景標記圖（圖9）進行距離變換，以達到細化目標物體的目的。由于梯度總是在邊緣處高，而在物體內部低，利用Sobel算子計算原始圖像梯度。修改梯度圖像，使得梯度圖像的區域最小值只出現在標記的前景和背景處，如圖10所示。

Fig.10 Gradient graph calculated by Sobel operator圖10 Sobel算子計算后的梯度圖

最后將修改后的梯度圖像（見圖10（b））利用分水嶺變換算法進行分割，分割結果如圖11所示。

Fig.11 Effect drawing of watershed transformation for modified gradient image圖11 修改后的梯度圖像進行分水嶺變換效果

3.3 比較實驗分析

對實驗結果采用

Dice

、

Jaccard

和

Precision

3種圖像分割評價指標進行分析，3種指標取值范圍均在0～1之間。該值越大，表示分割效果越好，具體定義如下：

其中，

SEG

表示算法對人體脊柱圖像自動分割的結果圖，

GT

表示采用MITK軟件手動分割的結果圖。MTIK軟件手動分割結果如圖12所示。

Fig.12 MITK manual segmentation result圖12 MITK手動分割結果

對原始人體脊柱CT圖像分別采用傳統分水嶺算法與K-means++結合的分割算法（K-means++and Watershed，KM++-W）、基于控制標記符分水嶺與K-means結合的分割算法（K-means and Marker-Controlled Watershed，KMMCW）、基于神經網絡的U-NET算法及本文提出的控制標記符分水嶺算法與K-means++結合的分割算法（Kmeans++and Marker-Controlled Watershed，KM++-MCW）進行分割實驗，然后對分割結果圖進行二值化處理。4種算法分割效果如圖13所示。

Fig.13 Segmentation rendering of four algorithms圖13 4種算法分割效果

將分割結果分別代入式（12）-式（14）進行計算，得到4種分割算法的

Dice、Jaccard

和

Precision

評價指標，并對其進行比較，如表1所示。

Table 1 Three evaluation coefficient values of four segmentation algorithms表1 4種分割算法的3種評價系數值

由表1可以看出，KM++-MCW算法的3種圖像分割評價指標數值均比另外3種算法高，從而證明了本文提出的人體脊柱圖像分割算法具有較高的準確度及良好的圖像分割效果。

對100幅不同部位的人體脊柱CT切片圖像分別采用3種算法進行實驗，并計算其

Dice

、

Jaccard

和

Precision

指標，指標數值分布如圖14所示。從圖14可以看出，4種算法的

Dice

系數均值分別為0.072 1（KM++-W）、0.726（KM-MCW）、0.740（U-NET）、0.892（KM++-MCW），

Jaccard

系數均值分別為0.052 5（KM++-W）、0.593（KM-MCW）、0.696（U-NET）、0.823（KM++-MCW），

Precision

系數均值分別為0.062 8（KM++-W）、0.736（KM-MCW）、0.751（U-NET）、0.854（KM++-MCW）。從4種算法的各評價系數均值可以看出，本文算法的分割效果較另外3種算法有明顯提高。其中，U-NET是深度學習算法，需要大量數據進行訓練并打上標簽，在數據量較少且無標簽的情況下，無監督學習的圖像分割算法更具有優勢。本文提出的人體脊柱圖像分割算法相比U-NET算法，3種評價指標平均值分別高出了約15.2%、12.7%、10.3%，對人體脊柱圖像的分割具有良好的精確度，能有效分割出人體脊柱圖像中的脊柱骨骼區域。

Fig.14 Index value distribution圖14 指標數值分布

4 結語

本文提出一種結合K-means++聚類算法與改進分水嶺算法的人體脊柱圖像分割方法，如果直接對人體脊柱原始圖像應用分水嶺變換方法進行分割會導致很高的過分割率。為防止這種情況，首先利用K-means++聚類算法對原始圖像進行初始分割，然后應用形態學技術進行去噪處理，標記出前景圖像和背景圖像，計算原始圖像的梯度圖像并對其進行修改，使得前景和背景標記處僅存在區域最小值，最后對修改后的梯度圖像應用分水嶺變換得到分割結果。實驗結果證明，該方法能很好地降低分水嶺變換的過分割率，并能有效將原始圖像中的人體脊柱骨骼區域分割出來。但是本文只針對人體脊柱圖像進行了實驗，對于如何將該算法運用于人體其它組織和器官的分割尚未進行研究，這也是未來研究的重點方向。