融合Harris角點檢測算法的肺實質分割方法

孫 紅，李 晶

(上海理工大學，上海 200093) (上?，F代光學系統重點實驗室，上海 200093)

1 引 言

由于大氣污染和吸煙，肺癌是世界上對人類健康和生命威脅最大的惡性腫瘤之一，發病率和死亡率目前均占所有惡性腫瘤的首位.降低肺癌死亡率、提高病人生存質量，最關鍵的是早期發現、早期診斷.早期肺癌往往表現為肺部小結節，肺部結節成為腫瘤是個漸進的過程，因此我們需要及時發現并鑒別診斷出肺結節.目前主要通過電子計算機斷層掃描(CT掃描)來檢測肺結節.

肺結節是指小于等于75px的肺部圓形病變.然而，由于肺內器官組織的復雜性，胸膜表面的肺結節呈現不同的特征，有些甚至僅通過肉眼無法輕易識別[1].CT掃描使得早期肺癌的發現率大大提高，而薄層CT技術有助于檢測更小的肺結節.但是，在使用薄層CT技術時，CT成像數量明顯增加，影像醫師的讀片量增加，工作量也增大，這樣很可能出現對肺結節的漏檢和錯誤判斷的情況.

本文著重針對胸膜表面粘連的肺結節，運用凸包與凸缺陷算法與harris角點檢測的方法對肺實質的輪廓進行修正，從而對這些的肺結節做出分割.本文采用的圖像樣本來源于公開數據庫The Cancer Imaging Archive (TCIA)[2]中的SPIE-AAPM Lung CT Challenge樣本集.

目前對肺結節自動檢測的研究，主要分為：肺部圖像的預處理方法、結節的特征提取、運用監督或非監督學習的方法對可疑結節進行檢測.在對肺實質分割之前，通常運用最大類間方差法進行自動閾值分割得到易于分析的二值化的圖像.對于與胸膜粘連的結節，Armato[3]提出了一種滾球算法.在對肺實質分割之前，對肺實質的輪廓進行補正，用來避免這些位于胸膜表面的結節的丟失.滾球法的主要的問題是需要對滾球模型設定半徑，不同的半徑會對分割結果產生很大的影響，容易產生過分割或邊界補正不完全.Pei等[4]采用基于凸包操作和射線追蹤的邊界校正算法進行邊界校正，此方法同樣可能會造成過分割的情況.Shen等[5]提出了一種無參數的肺結節分割算法.采用雙向鏈編碼方法，結合支持向量機(SVM)分類器，選擇性地平滑肺邊界，同時盡量減少相鄰區域的過度分割.HUANG[1]提出了一種基于邊界逼近的肺實質分割方法.針對自然彎曲的肺實質內輪廓與外輪廓，通過判斷自然凹陷最大深度與凹陷寬度的關系，區分自然凹陷與病灶點.由于肺部結構的復雜性以及邊界的不確定性.以上兩種算法[1，5]的缺點在于對于一些比較特殊的將連通的肺實質斷開的胸膜結節，該算法將無法對邊界進行補正.本文通過harris角點檢測的方法來彌補這一缺陷.

2 方 法

本文提出的方法主要涵蓋4個步驟：

1)運用最大類間方差法將原CT圖像轉化為二值圖像；

2)初步提取肺實質部分的輪廓；

3)運用凸包與凸缺陷以及角點檢測方法對邊界進行校正，得到完整的模板；

4)遍歷整個圖像根據得到的校正后的模板分割出肺實質內部的所有結節候選點.

2.1 預處理

本文運用最大類間方差法(又稱Otsu算法或大津法)對圖像進行預處理以得到二值化的圖像以便進一步分析.如圖1(b)所示，在對一幅圖像進行二值化處理時，希望得到一個最佳閾值將原圖像分成前景和背景兩個圖像.在Otsu算法中，通過使得類間方差最大來衡量前后背景的差異，即當類間方差取得最大值時，得到的閾值就是最佳閾值.因為方差是灰度分布均勻性的一種度量，那么前景與背景間的類間方差越大就說明灰度差別越大，也就越能區分前景和背景.

圖1 Otsu算法處理得到的二值化圖像Fig.1 Processed by the otsu algorithm

為了分割出圖1(a)中肺實質內部的肺結節可疑點，需要得到一個二值模板，這個模板對應于肺實質的部分的灰度值為255，而其余部分為0.用這個模板與原圖像相乘即可得到需要分割的可疑結節點.如圖2所示，先提取圖1(b)中肺實質的輪廓，再對輪廓內部做區域填充，即可得到二值模板.

2.2 凸包與凸缺陷

在2.1節中，初步得到的模板與原圖像相乘之后，圖(1)b中圓圈位置所標注的凹陷，即胸膜表面的結節無法被分割.因為與胸膜相粘連的結節灰度值與胸膜接近，在基于灰度的分割方法中，會被誤判為肺實質外的區域.為了修正這個問題，本文引入了輪廓的凸包與凸缺陷.

那么凸包是什么呢？給定平面上有一個(有限)點集，這個點集的凸包就是包含點集中所有點的最小面積的凸多邊形.本文通過這種方法來分割出圖1(b)中圓圈標注的凹陷.凸包最常用的凸包算法是Graham掃描法[6]和Jarvis步進法[7].本文運用的是Graham掃描法.將該方法運用在圖2中的輪廓圖像即可得到修正以后的模板圖像圖3(a)，將模板圖像與原圖像做乘運算后得到所有結節候選點的結果，如圖3(b)所示.處于胸膜表面的結節以及內部的所有組織都被分割出來，以便進行下一步分析.

圖2 輪廓提取與區域填充Fig.2 Contour extraction and area filling

圖3 修正以后的結果.Fig.3 Result of correction

2.3 harris角點檢測方法

針對文獻[1，5]中提出的問題，由于肺實質內部結構復雜，同一名患者的CT圖像在不同層的肺實質輪廓都不相同.

圖4 方法的失效Fig.4 Method failed

對于一些特殊位置的結節，如圖4所示.粘連于胸膜表面的結節使得左肺邊緣出現“斷裂”.

本文運用harris角點檢測算法對該問題進行修正，得到期望的理想輪廓線與模板，從而使位于邊緣的病灶不產生丟失.對于圖4所產生的問題，根據大量的樣本顯示，這類情況多出現在肺實質區域狹長的層，因病灶的緣故使得肺實質之間產生“斷裂”.

通過對圖像進行行掃描可尋找到該斷裂位置，并以該位置為中心設定感興趣區域.在該斷裂區域尋找到肺實質區域邊緣的角點并將這些點用一個最小多邊形將其連接起來，這樣填充以后即可得到可以分割出病灶的模板.

Harris角點檢測是由Harris在文獻[8]中提出的，其主要思想是用一個模板窗口在圖像上滑動，在垂直和水平方向上都產生較大變化的就是角點.首先定義自相關函數：

(1)

令D(x，y)=Iu+x，v+y-Iu，v，將D(x，y)在(0，0)按泰勒級數展開得：

(2)

有D(0，0)=Iu+0，v+0-Iu，v，高階無窮小忽略不計：

(3)

所以(1)可以化為：

(4)

(5)

(5)式本質是一個橢圓函數.橢圓的扁率和尺寸是由矩陣H的特征值決定的，橢圓函數的特征值與圖像中的角點，直線和平滑區域的關系可分為：

1)直線.一個特征值較大，另外一個特征值較小.

2)平滑區域.兩個特征值都很小，且近似相等.自相關函數值在個方向上都很小.

3)角點.兩個特征值都很大，且近似相等，自相關函數在所有方向上都很大.

為了判定圖像中角點，根據hessian矩陣中行列式的值、矩陣的跡，將它們代入角點計算公式：

(6)

當R大于一個特定閾值時，該特征點(u，v)被認為是角點.由于通常情況下一個連通域分裂為兩個區域，為了補正它們之間的縫隙，尋找曲率低的平滑的點，將這些點連接所圍成的區域不足以使缺失區域補全，因此考慮尋找斷裂區域周圍的角點作為候選點進行補正.

3 實驗結果與分析

3.1 harris角點檢測方法

如圖5所示，左肺的病灶出現在斷裂的位置，運用凸缺陷的方法無法得到包含病灶點的模板.運用上文中的harris角點檢測算法，首先通過行掃描的方法對每一行像素點的個數進行積分，若某一行積分值圖變為0，則認為該位置存在斷裂的區域，并將上一行中首個不為0的像素所在的列數作為ROI區域的中心點，并以該點為中心設置一個感興趣40×40像素的區域.即僅對該斷裂區域檢測肺實質部分的角點，將滿足條件的候選點按逆時針排列后將它們兩兩連結成一個多邊形，如圖5所示.

圖5(b)是原圖中感興趣區域補正的結果，大圖是結果全體.對該補正后的空洞區域做填充即得到整個連通的模板，用此模板與原圖像相乘可得病灶候選點.本文對The Cancer Imaging Archive (TCIA)的SPIE-AAPM Lung CT Challenge樣本集中的263張CT圖像(來自不同的6例病例)，對于樣本集中存在如文中所說凹陷的圖像均可做出正確的補正.圖6是實驗結果，一例在左肺產生斷裂，另一例是在右肺邊緣產生病灶的樣本圖像，本結果來自于windows 8環境下、CPU主頻2.5GHz、內存3.9G、Visual Studio 2013、opencv 2.4.9.

圖5 原圖像與補正后的結果Fig.5 Original image and result of the harris algorithm

圖6 原圖像與修正結果對比Fig.6 Original image and result of correction

3.2 滾球算法結果對比與分析

文獻[3]中首次提出運用一個固定半徑的圓繞著輪廓滾動一周.在這過程中，平滑區域應當與圓相切即只有一個交點.在凹陷處，即表面存在結節處應當存在一個以上的交點.將這兩個交點用直線相連即可完成補正的方法.本結果來自于Image J圖像處理軟件自帶的rolling ball算法.

圖7為本文修正結果與文獻[3]的滾球法的效果對比，圖示說明中r表示滾球法中設定的圓半徑參數.對照圖2(b)修正前的圖像，圖(a)-(d)均正確地將結節處在的邊緣進行了修復.圖(b)的結果與本文方法相近.而圖(c)和圖(d)當增大參數半徑時就會產生過度分割，即將非實質區域誤認為是結節區域進行了不必要的修正.又如圖(f)-(h)，在不同r參數下有不同的分割結果.由此可見，文獻[3]中的分割方法非常依賴參數r的設定，該參數r并不具有自適應性.

圖7 本文方法與滾球法對比Fig.7 Comparison between rolling ball algorithm and the method in this paper

3.3 DRLSE(Distance Regularized Level Set Evolution)模型的水平集方法結果與分析

DRLSE模型是由LI在文獻[9]中提出的，運用該方法所產生的結果來自于 windows 8 環境下、CPU主頻2.5GHz、內存3.9G、Matlab 2014b.

如圖8所示，該方法相比于傳統的水平集方法，提出了距離正則化項用來修正水平集函數與符號距離函數之間的差，從而保證最終演化結果精確地在邊緣處停止并舍棄了重新初始化水平集函數的步驟.本文根據文獻[9]通過水平集方法尋找到包含病灶的肺實質輪廓，從而得到分割結節候選點的二值模板.而由結果可見，該方法主要有二點問題，第一，水平集方法需要設定初始的輪廓位置，而CT影像在不同層的位置和大小都有所不同，所以曲線演化的初始位置不易選擇.第二，由圖8(c)可見，在曲線演化至接近目標邊緣時，病灶位置的曲線的拓撲結構發生了變化，并未將該病灶包含在輪廓內部.

圖8 水平集方法Fig.8 Level set method

3.4 邊界逼近法結果與分析

將邊界逼近法運用于肺實質分割是由HUANG在文獻[1]中提出的，該方法中選取的樣本大小均為512像素×512像素，實驗環境為Windows 8.1、主頻2.5 GHz、內存4GB，開發平臺是VS2013、OPENCV2.4.9.

如圖9(b)所示，該算法對于內輪廓上自然凹陷區域分割較好，同時也較好保留外輪廓存在病灶的區域.但是圖9(d)所示結果，未檢測出右肺的病灶點，粘連于胸膜表面的結節使得右肺邊緣出現“斷裂”.該算法主要的問題是：對如圖9所示樣本，因病灶的緣故使得肺實質之間產生“斷裂”，將肺實質割裂成兩塊獨立的連通域，從而并未將該病灶包含在輪廓內部.

圖9 邊界逼近法Fig. 9 Border approximation

4 討論和結論

本文主要針對CT圖像中位于胸膜表面的結節的檢測提出了解決方法，運用graham掃描算法對肺部圖像的輪廓進行修補，從而得到完整的二值模板來分割位于邊緣的病灶.并通過rolling ball、水平集方法、邊界逼近法與本文角點檢測方法做了對比.由于rolling ball與水平集方法都需要手動設定固定的參數或初始化位置，而同一位病人的CT圖像在不同層上的肺實質圖像差異很大，因此難以設定初始值和初始輪廓，從而不易將邊緣的病灶檢測出來.同時邊界逼近法也不能很好的檢測出粘連于胸膜表面的結節，可能會出現“斷裂”現象.Harris角點檢測的方法對文獻[1]和文獻[5]提出的問題具有針對性，如實驗結果可見，該區域無需手動設定參數，當尋找到該ROI區域后自動進行補正.因此，本文方法針對將連通的肺實質斷開的胸膜結節檢測具有有效性.本研究的下一步工作是運用機器學習的方法根據形態、大小、密度等特征對分割出的胸膜表面的結節進行綜合分析，從而對其進行分類，以完善計算機輔助診斷系統.