劉昊,王瑜,王怡寧,徐橙
1.北京工商大學人工智能學院,北京100048;2.中國醫學科學院北京協和醫院,北京100730
作為世界人口主要的死亡原因之一,心血管疾病早期體征不明顯且易導致患者猝死。在所有心血管疾病中,多數是冠狀動脈問題,主要由左或右動脈狹窄引起[1]。因此,冠狀動脈疾病的早期診斷、評估風險和治療規劃至關重要。
近年來,冠狀動脈計算機斷層血管成像檢查(Coronary Computed Tomography Angiography,CCTA)由于其無創性、多視角等特點,已成為心臟病診斷和治療的重要工具。但是由于CCTA 數據量大,復雜性高,因此即使專家進行圖像診斷也通常需要花費大量時間。此外,由于CCTA 圖像的噪聲影響,冠脈結構復雜且被其他組織環繞,迅速準確的冠狀動脈分割難以實現[2]。冠狀動脈序列圖像分割方法主要包括區域生長算法[3]、水平集[4]、圖割方法[5]、分水嶺算法[6]和基于深度學習的方法[7]等。其中,分水嶺算法在醫學圖像分割中應用廣泛,具有運算速度快、計算簡單、邊緣定位精準等優點,但易受到微弱邊緣的影響,很容易因噪聲和對比度低等原因造成過分割和易丟失重要輪廓等問題。為了克服分水嶺算法的缺陷,可以結合有效預處理方法先行去除噪聲,或對分割后的圖像進行區域融合。
為此,本文提出一種新穎的冠狀動脈分割方法,基本思路是,利用分水嶺算法將冠脈圖像分割為子區域,然后交互式分割出種子區域,接著根據形態學和同一組織在斷層圖像的連續性,以分水嶺子區域為單位進行區域生長,最終實現冠狀動脈的三維分割。本算法可解決傳統分水嶺算法造成的過分割及低對比度圖像易丟失重要輪廓的問題,且運算速度快、分割精準,同時利用形態學思想解決了少數相鄰層冠狀動脈區域不重疊的現象。
本文提出的方法主要包括3 個步驟:(1)數據預處理;(2)利用分水嶺算法將CCTA 圖像分割為子區域;(3)利用以子區域為單位的區域生長方法,得到最后的冠狀動脈分割結果。
有效的圖像預處理是冠狀動脈準確分割的基礎,針對CCTA 圖像結構的復雜性,本文使用閾值分割方法[8]對原始圖像進行預處理,以灰度直方圖中冠狀動脈最高灰度值為閾值,低于閾值的區域灰度值保持不變,其余區域灰度值設為0。具體公式如下所示:

其中,f(x,y,z)表示原圖像(x,y,z)點的灰度值,T表示冠狀動脈最高灰度值。
本文根據經驗和調試,設置閾值為1 500(數據格式為DICOM)。閾值分割方法能很大程度剔除冠狀動脈周圍高密度組織的干擾。此外,針對CCTA 圖像包含噪聲的問題,使用空間域高斯平滑方法[9]對圖像進行平滑處理,為后續分水嶺分割提供高質量數據。
分水嶺算法[10]是一種流行且簡單的圖像分割方法,其本質源于地理。假設景觀被雨水淹沒,雨水自然沿著最陡峭的路徑流下,并最終落入許多盆地區域。分水嶺是這些區域的分界線,是景觀的自然分割輪廓。將圖像看作是景觀,從而利用分水嶺算法將圖像分割成許多細小子區域。但是,傳統分水嶺算法對圖像中存在的噪聲、物體表面細微的灰度變換較為敏感,通常會出現過分割現象,導致分割結果不理想。為此本文根據CCTA 圖像中細節豐富的特點,提出一種改進的分水嶺算法。首先對于預處理后圖像利用Sobel 算子[11]得到梯度幅值圖像,然后增強圖像中極值的對比度,同時保持圖像相對平滑區域,最后利用形態學修正梯度圖像,并提取圖像的標記修正梯度圖像進行分水嶺變換,得到分割結果。
1.2.1 梯度圖像計算預處理后的冠狀動脈CCTA 圖像,雖然大量非冠狀動脈組織被剔除,且圖像噪聲得到抑制,但是如果僅僅通過傳統分水嶺算法的梯度計算方法得到圖像的邊緣,容易出現信息丟失現象,影響算法后續的分割效果。因此,本文采用sobel 算子計算梯度圖像,計算公式如下:

其中,f(x,y)表示圖像(x,y)點的灰度值,Gx表示橫向梯度檢測的灰度值,Gy表示縱向梯度檢測的灰度值。該算子較傳統梯度計算方法有明顯優勢,能夠找出圖像中的細小邊緣,并且可以生成噪聲水平較低的圖像,更適用于CCTA圖像。
1.2.2 前景與背景標記雖然利用sobel 算子代替傳統的梯度計算方法一定程度上解決了邊緣信息丟失以及圖像噪聲現象,但是圖像中仍然存在與目標無關的極小值,導致過分分割的現象。為了解決這個問題,本文使用標記控制的分水嶺分割。如果可以提前識別并標記前景對象和背景位置,則分水嶺算法的分割效果會更好。本文使用形態學運算和最大類間方差算法[12]來提取圖像的標記,具體步驟如下:(1)計算圖像的局部極大值;(2)利用公式f=(fΘB)⊕B對局部極大值圖像進行關操作,并作為前景標記;(3)利用公式f=f⊕B對圖像進行膨脹操作;(4)利用最大類間方差算法將圖像二值化;(5)計算二值化圖像距離變換,并利用分水嶺變換求取其相鄰區域間的分界線,以此作為背景標記。
傳統區域生長算法[13]以像素單位作為生長的種子點,因此此類算法不僅計算時間長,且對噪聲和灰度不均較為敏感,適合一些分割灰度均勻的平滑區域,無法準確高效地提取醫學影像中的血管組織。針對上述問題,本文以分水嶺子區域為單元,交互式分割種子區域,根據斷層圖像上下層間的相似性,以及冠狀動脈的連通性,利用相鄰層圖像中對應區域的面積重疊關系以及密度相似性,自上而下進行區域生長,最終實現冠狀動脈的三維分割。
1.3.1 種子區域的選擇種子區域的選擇直接影響最終的分割結果,尤其是在斷層圖像分割中。CCTA 圖像對于同一組織在斷層圖像的相鄰層之間普遍具有相似性。因此,可以根據冠狀動脈斷層間的相似性來選取種子區域。在本文研究的CCTA 圖像中,除了斷層圖像中部分冠狀動脈目標區域像素在前一張斷層圖像的目標區域中沒有重疊的像素點外,其余所有冠狀動脈目標區域均能夠在前一張斷層圖像中周圍找到重疊的區域像素。因此,選擇無法通過上層投影而分割的目標區域作為種子區域,進行序列圖像的三維區域生長分割,通過對圖像的種子區域集合設置向下投影,從而實現序列圖像的三維分割,如此就不會錯過部分斷層圖像中的冠狀動脈區域。
1.3.2 區域生長規則的選擇在斷層圖像中,相鄰層之間具有良好的連續性,從垂直于斷層平面的方向來看,相鄰層圖像中同一組織普遍存在面積重疊和密度相似關系。但是,在個別情況下相鄰層冠狀動脈斷面相互不重疊。盡管在這種情況下,二者并不重疊,但其相隔距離并不遠[14]。針對這種不重疊現象,本文運用公式f=f⊕B進行形態學膨脹,對預分割的兩相鄰層中的前一層目標區域進行膨脹處理后再利用面積重疊關系合并下層目標區域,其中B為結構元素。本文選用圓形結構元素,半徑為3 個像素,如圖1所示。其中,圓形虛線是冠狀動脈截面進行膨脹后的邊界,圓形虛線所包圍的區域與相鄰下層中的冠狀動脈重疊,可以看出形態學膨脹可以克服同一組織不重疊的問題。

圖1 冠狀動脈斷層膨脹處理Fig.1 Coronary artery cross-section dilation
本研究選擇面積重疊關系和密度相似性作為區域生長的生長準則,具體如下式所示:

其中,RAnew為新的分水嶺子區域,RAbefore‐dilate為上一層已分割目標區域進行膨脹后的區域,RAnew?RAbefore‐dilate為RAnew和RAbefore‐dilate相交的像素點個數;RMnew為新的分水嶺子區域的灰度均值,RMbefore為上一層已分割目標區域的灰度均值。本文通過上述生長準則,將符合條件的分水嶺子區域納入冠狀動脈區域,并計算該層目標區域,進行膨脹處理后作為新的種子區域,進而搜索下一層符合條件的分水嶺子區域,直到沒有新的分水嶺子區域納入冠狀動脈區域為止。
提出的利用交互方法分割出冠脈的步驟如下:(1)將種子區域所在層作為當前層,即層1。(2)對層1目標區域進行膨脹,在層2 中搜索與層1 中膨脹后目標區域滿足面積重疊和密度相似關系的分水嶺子區域,記錄層2 中滿足條件的子區域并進行區域合并,其結果作為層2 中目標區域的分割結果。(3)將層2作為當前層,并重復步驟(2)合并下一層中的冠狀動脈區域,直到在下一層中沒有滿足條件的子區域時,停止搜索。依照上述步驟進行循環生長,最終實現冠狀動脈各分支的三維分割。
實驗數據集來自中國醫學科學院北京協和醫院的CCTA 圖像,包含231 幅切片,并由專家標注數據,圖像為DICOM 格式,大小為512×512 像素。為了驗證提出方法的有效性,本文精心設計了一系列實驗,對比了現有算法和提出方法的結果,并與專家標注結果進行了比較。實驗是在Intel Core (TM)i5‐7400CPU 3.00 GHz 的PC 機上基于MATLAB R2019a軟件平臺實現的。經過多次優化,實驗具體的參數設定如下,冠狀動脈灰度值上限是T= 1 510,形態學運算結構元素是半徑為3的圓形結構,面積重疊閾值Thδ1= 0.48,相似性密度范圍Thδ2= 380,其中選擇初始層種子區域的平均密度和膨脹后的區域面積作為實驗的初始模式點進行計算,種子區域選擇情況如圖2所示。

圖2 交互式分割種子區域的選擇Fig.2 Selection of seed region in interactive segmentation
本文對比了傳統分水嶺方法[10]、三維區域生長法[13]和本文提出方法的冠狀動脈分割結果,3 種方法將在相同預處理和種子點選取的情況下進行實驗,具體結果如圖3和表1所示。使用有經驗的放射科醫生手動分割結果作為標準,分別利用Jaccard 指數和Dice 相似系數得分進行客觀評測,兩種評測標準的范圍都為[0,1],其中1 表示重疊較高,0 表示無重疊,具體公式如下:

表1 冠狀動脈分割的評估Tab.1 Evaluation of coronary artery segmentation results

圖3 分割結果及醫生標注Fig.3 Segmentation results obtained by different methods and manual segmentation results

其中,Y是醫生標注結果,Yp是算法分割結果。
從圖3中可以看出,傳統分水嶺算法由于過分割現象嚴重,分水嶺子區域過于粉碎,很難選取合適的面積重疊閾值,導致錯分割或欠分割現象嚴重。而三維區域生長算法的結果同樣極易發生欠分割現象,這是因為傳統區域生長算法只能從體素灰度特征上設定相似性規則,而對于冠狀動脈和周圍組織的模糊邊界,則很難選取合適的閾值進行有效分割。而本文的方法以分水嶺子區域為單位,利用面積重疊關系和密度相似性進行區域生長,可以有效抑制圖像中其他組織對于冠狀動脈的影響,雖然在細節上和手動分割結果有些誤差,但是冠狀動脈的主枝區域是十分準確的。
從表1可以看出,在兩個定量評價指標中,本文算法的Jaccard 指數為0.936 8,Dice 相似系數得分為0.967 3。兩者的指標都接近1,高于傳統分水嶺算法和三維區域生長算法的定量評價指標,說明本文的分割結果更為精確,即本文算法可以更有效地分割冠狀動脈。
本文結合分水嶺算法與區域生長算法,提出了一種新穎的冠狀動脈三維分割方法,首先利用閾值法和高斯濾波對冠狀動脈CCTA 圖像進行預處理,然后利用改進的分水嶺算法將斷層圖像分割為若干數個子區域,并利用交互式方法分割出種子區域,最后利用形態學方法和斷層圖像相鄰關系,實現冠狀動脈的三維分割。實驗結果表明,本文方法能有效地分割出冠狀動脈,但在血管細小分支部分有待提高,后續工作將進一步研究有效提高血管細小分支部分精準分割的方法。