王 南,許道云
貴州大學 計算機科學與技術學院,貴陽 550025
左心室MR 圖像分割是運用圖像分割的方法對左心室內外膜進行準確分割,其結果可以為醫生定量地分析左心室收縮功能,提供參數指標。目前,關于左心室的分割方法主要分為圖像驅動模型、基于圖譜的方法[1]、基于深度學習的方法[2-5]以及主動輪廓模型等。圖像驅動模型包括基于閾值、基于區域信息以及像素相似性質分類等方法,該類方法通常具有分割速度快、實現較為簡單的特點,但其不能與左心室的先驗信息有效結合,對左心室邊界弱邊界、對比度不高的圖像分割能力有限。基于圖譜的方法在分割給定的圖像時,需要與大型訓練集匹配分割圖像,該方法的分割結果受配準圖像數據集的大小所影響[6]。基于深度學習的方法,在處理左心室圖像時,具有魯棒提取圖像特征以及實現左心室全自動分割的優點,但完成該方法的神經網絡訓練,需要大量圖像數據,并且訓練時間較長。主動輪廓模型相較于其他方法,其優點是在高噪聲的情況下,也能得到較好的分割效果,并且容易與其他幾種方法進行結合,分割速度較快。
主動輪廓模型分為參數輪廓模型和幾何輪廓模型。參數輪廓模型是Kass 等提出的傳統Snake 模型及其改進模型[7-9],其主要是將輪廓曲線向檢測物邊緣演化的過程與求解能量泛函最小值的過程相關聯。該方法在給定適當的初始輪廓后,可以自動收斂到能量最小狀態,具有對噪聲和對比度不敏感的優點,但其對初始輪廓的位置具有很強的依賴性,在圖像的弱邊界處容易陷入局部最優。
幾何活動輪廓模型是基于曲線演化理論與水平集的方法,將平面曲線隱含的表達為三維曲面函數,通過曲面的演化來隱含地求解曲線演化[10]。C-V模型[11]是幾何輪廓模型中一個具有代表性的基于區域信息進行圖像分割的算法,其具有易于處理拓撲結構變化、可分割較為復雜的對象、容易實現低維到高維的擴展、對初始輪廓不敏感等優勢,但同時具有一定局限性,該模型對強度不均勻圖像會產生不良分割結果,需要反復初始水平集,造成計算成本高等限制,眾多學者在該算法的基礎上進行改進與創新[12-17]。在文獻[13]中,盧小鵬等人將目標區域與背景區域的灰度平均值之間的差值添加到C-V模型面積項中,使得C-V模型對亮度不均勻圖像的分割能力有所提升。在文獻[14]中,Wang等人提出一種將全局和局部統計信息結合的水平集新方法,克服了C-V模型不能成功分割強度不均勻圖像的缺點。Murtaza等人在文獻[14]基礎上,使用變異系數代替方差,提出新的圖像變分模型PSM[15],在對圖像分割實驗中取得更好的分割結果。Li 等人在文獻[16]中提出一種對原圖像進行擬合的方法。文獻[17]中,盧振泰等人提出了一種基于熵和局部鄰域信息的高斯約束C-V 模型對左心室進行分割,該模型具有算法執行速度快、避免重新水平集初始化、提高分割精度。文獻[18]中,Li 等人提出距離正則化水平集演化(DRLSE)方法在能量函數中增加了一個正則項,以避免在水平集演化過程中出現形狀不規則和不穩定性等情況,并且在該模型不需要重新初始化水平集,但該方法對初始輪廓敏感且在對比度不高的圖片中演化過程容易陷入局部最小值。在運用DLRSE算法對左心室MR圖像內膜進行分割時,會因為受到乳突肌與小梁組織的干擾而造成曲線陷入局部最優的結果,為了避免這種情況的出現,文獻[19]中,根據左心室內膜的先驗知識,劉宇提出一種保持凸性的水平集演化模型,該方法能保證曲線演化過程中具有良好的凸性,最終演化到左心室內膜邊緣,但該算法存在對初始輪廓敏感,對圖像弱邊緣分割不精確的缺點。
綜上所述,在圖像分割中DRLSE 算法具有對初始輪廓敏感,演化過程僅依賴圖像梯度信息的缺陷。為降低DRLSE 算法對初始輪廓的敏感程度,提升分割弱邊緣、強度不均勻圖像的能力,提高算法對左心室MR 圖像內外的分割精度,提出一種適用于弱邊緣信息左心室MR圖像的DRLSE分割算法,實現對左心室內膜與外膜分割。在本文的內膜算法中,首先,使用Li提出的圖像擬合方法,擬合得到新圖像,計算出原圖像與新圖像的圖像差,替換原PSM模型中的圖像差,形成PSM新的局部項;其次,將文獻[13]中C-V 面積項與上一步得到的PSM 新的局部項添加進 DRLSE 中;再次,將 CPLSE 模型的保持曲線凸性的曲率因子添加到DRLSE 模型中;最終實現對左心室MR 圖像的內膜分割。在對左心室外膜進行分割時,首先,運用內膜分割結果,結合左心室心肌各方向厚度基本一致的形狀先驗,構造外膜形狀約束項,形成左心室外膜形狀約束項;其次,將外膜形狀約束項與內膜分割算法中的PSM新的局部項、文獻[13]中的C-V 模型面積項嵌入到DRLSE 算法中;最終實現對左心室MR圖像外膜的分割。最后通過實驗,對本文提出的模型進行測試,對比實驗結果進行分析。
在C-V 模型中,將Ω?R2作為圖像范圍,u:Ω→R作為給出的灰度圖像,定義Ωin是在圖像u上被閉合曲線C所劃分的內部區域,Ωout是在圖像u上被閉合曲線C所劃分的外部區域。C-V模型的能量泛函定義如下:

其中,μ≥0 ,v≥0 ,λ1>0 ,λ2>0 是各個項的系數,Length(C)是曲線長度項有平滑曲線的作用,Area(inside(C))是Ωin面積項,一般設其系數v=0,后兩項為保真項是用于驅動曲線C向真實邊界演化。
水平集定義表示如下:

E(C)表示為水平集形式:

其中,c1、c2分別為Ωin、Ωout的灰度平均值,其計算公式如下:

δ是Dirac函數,H是Heaviside函數。

當曲線演化到目標邊緣的時,泛函式(3)將達到最小,其最小化可以通過以下的梯度下降流來求解:

C-V模型可以很好地對圖像進行分割,但對灰度不均勻、均勻強度重疊的圖像分割效果很差。2010 年Wang等人提出LCV模型,該模型對圖像灰度不均勻圖像有很好的分割效果,但在圖像均勻強度一致部分所產生的分割結果較差。
PSM 模型是 2011 年由 Murtaza 等人在 LCV 模型基礎上提出,該模型具體如下:

其中,u*=gk?u-u,gk是高斯核,u*:Ω*→R作為給出的新灰度圖像,是曲線C將u*內外兩部分區域,m1、m2分別表示區域灰度的平均值。水平集下的PSM公式如下:

其中:

PSM模型的最小化求解過程可以通過以下方程:

PSM 模型在灰度非均勻圖像的分割結果要優于LCV、C-V模型。
DRLSE 模型作為一種基于區域的水平集算法,是由Li等人于2010年提出的模型。
定義水平集函數:

DRLSE模型的能量泛函定義為:

其中,μ≥0 ,λ >0 ,α∈R,Rg(φ)是距離正則項,Lg(φ)是加權長度項可以驅使邊緣向零水平集向著目標邊緣方向演化,Ag(φ)是加權面積項可以改變水平集的演化速度,它們分別被定義為:

其中,P為勢函數,是邊緣指示函數,Gσ是高斯核;? 是梯度算子。在DRLSE方法中,選取:

能量泛函式(14)的最小化可以通過以下梯度下降流來求解:

DRLSE模型允許有更靈活、更有效的初始化,在曲線演化過程中顯著地減少了迭代次數與計算時間,但該方法對初始輪廓敏感,僅依靠梯度信息驅動曲線演化,模型容易陷入局部最優。
對左心室圖像分割時,DRLSE 方法存在對初始位置敏感且無法有效分割弱邊界圖像的問題,為改善以上情況,本文運用Li提出的圖像擬合方法得到新圖像,將新圖像與原圖像之間的圖像差,替換新的PSM 模型中局部項的圖像差,將新的局部統計信息與梯度信息結合改進DRLSE模型。為降低新局部項受到圖像灰度不均勻情況的干擾,將對分割灰度不均勻圖像有一定貢獻的文獻[13]中C-V的面積項加入到DRLSE模型中,得到初步改進DRLSE,其水平集下公式為:

其中,Hε1(x)是選取與DRLSE中相同的Hε(x) ,Hε2(x)選取與PSM中相同的Hε(x),u*=f-u,m1、m2是曲線劃分u*產生區域內外的灰度平均值,擬合圖像f的公式為:

其中:

gε是卷積核函數。
運用|?φ|來消除傳統C-V模型中的Dirac函數對檢測遠邊緣的抑制,并將邊停函數g加入到演化函數中。基于此,改進的DRLSE 能量泛函方程的最小化求解對應的梯度下降流方程:

在內膜分割中,內膜形狀類似于橢圓,而且DRLSE算法受到乳突肌或小梁的影響,內膜分割結果可能是凹形曲線,出現欠分割情況,將保持曲線凸性的曲率因子加入到上述改進的DRLSE 方法中,這里對文獻[19]中kmean公式作修改,定義:

最終內膜分割算法水平集最小化演化方程如下:

在分割外膜時,外膜與肌肉組織的灰度接近,在使用DRLSE 模型對圖像外膜進行分割時因受到弱邊緣、對比度低等因素的影響,造成過分割或欠分割等不好結果,觀察圖像中左心室的特性,發現左心室心肌各方向厚度基本一致[10],本文設計的外膜形狀約束項是采用已分割好的內膜與已知左心室最小厚度Rmin結合先驗信息[20],結合先驗信息形成形狀約束,以實現對演化過程中曲線的約束,具體公式如下:

其中dis是外膜分割曲線c上的點與內膜曲線的最短距離。

最終外膜分割算法水平集最小化演化方程如下:

其中,ω >0 ,γ >0 是兩項牽引力的系數。
本文方法的具體步驟如下:
(1)內膜分割
①設定模型參數,對初始水平集φ0,即φt(t=0),最大迭代次數T,0 ≤t≤T;
②根據當前φt,結合公式(4),計算更新當前區域內外平均值c1、c2;
③ 根據公式(20),擬合出圖像ft,計算更新u*、m1、m2;
④ 據公式(22)、(23)、(24),計算更新k、s(k)in、kmean;
⑤根據公式(25),更新φt+1;
⑥ 判斷t >T,若是,則停止;否則判斷φt+1與φt輪廓是否相同,若是,則停止;否則,φt=φt+1,轉回②。
內膜分割流程圖如圖1。

圖1 內膜分割算法流程圖
(2)外膜分割
①設定模型參數,對初始水平集φ0,即φt(t=0),最大迭代次數T,0 ≤t≤T;
②根據當前φt,結合公式(4),計算更新當前區域內外平均值c1,c2;
③根據公式(20),擬合出圖像ft,計算更新u*、m1、m2;
④ 根據公式(26)、(27)、(28),計算更新dmean、s(k)out;
⑤根據公式(29),更新φt+1;
⑥ 判斷t >T,若是,則停止;否則判斷φt+1與φt輪廓是否相同,若是,則停止;否則,φt=φt+1,轉回②。
外膜分割流程圖如圖2。

圖2 外膜分割算法流程圖
本文對提出的方法通過實驗進行驗證,對其結果進行定量分析,并與DRLSE、CPLSE方法進行對比。近年來,基于深度學習的方法在圖像分割領域應用成為主流方法,對左心室也有較好的分割效果,選擇U-Net 網絡加入到對比方法中。實驗進行驗證所使用的數據集來源于http://www.cse.yorku.ca/~mridataset/的Cardiac MRI dataset,該數據集是從33 個受試者獲得共7 980 張的心臟MR 圖像,并對圖像進行手工分割獲得5 011 個分段MR 圖像和10 022 個輪廓。每個實驗受測對象由20 幀和8~15個切片組。對U-Net網絡進行訓練時,從數據集中選取4 200 張帶有標簽的左心室圖像,將數據圖像分為兩部分,總數量10%作為測試數據集,90%作為訓練數據集輸入到U-Net網絡之中進行訓練,采用二分類交叉熵作為損失函數,經過訓練100 次循環,得到訓練好的U-Net網絡,其內膜與外膜訓練的損失函數值分別為0.001 與0.002。在數據集中選取不同受試者、不同幀、不同切片的左心室MR圖像,分別使用DRLSE、DRLSE-shape、CPLSE、訓練好的U-Net 網絡和本文模型對MR圖像進行分割。選擇的實驗環境為Windows 7 操作系統,CPU為AMD Athlon II X4 640,5.75 GB內存,實驗平臺是Python3.5。
從數據集中選取不同受試者選取不同幀、不同層、像素為256×256的左心室MR圖像進行實驗。分別賦予每幅圖各自不同的內膜初始化位置,運用本文左心室內膜分割算法對MR圖像進行分割,該算法涉及的參數設置為μ=0.04,λ1=2,α=-1.5,λ2=0.1,v=0.001,η=1。為每幅圖賦予各自不同的外膜初始化輪廓,運用本文左心室外膜分割算法對MR 圖像進行分割,其涉及的參數設置為μ=0.04,λ1=2,α=5,λ2=0.1,v=0.001,η=0.25,γ=5。
(1)給定一張左心室的MR圖像,分別運用DRLSE、CPLSE及本文算法,在給與三種方法相同初始化位置的情況下,對MR 圖像進行內膜分割,其分割結果如圖3,其中每張圖的紅色虛線為手工輪廓,綠色曲線為算法分割結果。

圖3 同一MR圖像不同初始位置的分割結果
(2)運用(1)中的三種方法對不同類型左心室MR圖片的內膜進行分割,結果如圖4。

圖4 正常的MR圖像分割結果
圖5~圖7 中,情形 1、2、3 依次代表圖像含有乳突肌、弱邊界、乳突肌和弱邊界的情況。紅色虛線為手工輪廓結果,綠色曲線為算法分割結果。

圖5 含情形1的MR圖像分割結果

圖6 含情形2的MR圖像分割結果

圖7 含情形3的MR圖像分割結果
(3)運用DRLSE、DRLSE-shape、U-Net 網絡以及本文算法對六個不同左心室圖像進行外膜分割,最終MR圖像內外膜分割結果如圖8。

圖8 不同左心室內外膜的分割結果
圖8(a)~(f)中的紅色曲線代表的是專家手工輪廓;白色虛線代表的是DRLSE算法對MR圖像內外膜分割結果;紫色虛線代表的是DRSLE-shape對MR圖像外膜分割結果;藍色曲線代表的是U-Net 網絡對MR 圖像內外膜的分割結果;綠色曲線代表的是本文算法對左心室內與外膜的分割結果。
(4)將專家手工輪廓、本文算法分割好的內外膜結果進行二值化后得到結果如圖9。
為了對本文提出的模型進行定量評價,采用相似系數(Dice Similarity Coefficient,DSC)、平均垂直距離(Average Perpendicular Distance,APD)。DSC是對兩個區域重疊度的測量,越接近1,兩幅圖的重疊度越高;APD 是自動分割所得輪廓曲線上的全部點與手動輪廓曲線的最短垂直距離的平均值,其值越小代表輪廓曲線重合越高;質心距離(Centroid Distance,CD)是自動分割環形質心與手工分割環形質心的歐氏距離。兩者指標定義如下:

圖9 不同左心室心肌分割結果的二值圖

其中,A代表專家已分割輪廓,Aor代表專家手工分割內外膜輪廓之間區域,B代表經過算法自動分割輪廓,Bre代表自動分割內外膜輪廓之間區域,|Ω|代表計算Ω區域內像素點的總數,例如,|A|+|B|代表兩個區域內全部像素點個數,|A?B|代表兩個區域相交部分像素點個數,xi代表B中的像素點,(xor,yor)代表Aor的質心,(xre,yre)代表Bre的質心。
表1 給出的是運用DRLSE、DRLSE-shape、CPLSE、U-Net 網絡與本文算法對左心室內外膜分割所用平均時間。

表1 各分割算法的平均時耗s
表2給出的是運用DRLSE、CPLSE和U-Net網絡與本文算法分割不同心臟左心室、不同層MR 圖像內膜,所得分割結果的DSC與APD值。
表3給出的是運用DRLSE、DRLSE-shape、U-Net網絡與本文算法分割不同心臟左心室、不同層MR圖像外膜,所得分割結果的DSC和APD值。
表4給出的是,針對不同圖像,分別采用本文所提算法對內膜與外膜的分割結果、U-Net網絡對內膜與外膜的分割結果、DRLSE對內膜與外膜的分割結果、DRSLE對內膜分割的結果與DRLSE-shape 對外膜分割的結果進行組合得到左心室心肌結果,并將其結果與專家手工分割的心肌結果進行比較,所得到的DSC和CD值。

表2 不同圖片內膜分割結果的APD和DSC值

表3 不同圖片外膜分割結果的APD和DSC值

表4 不同左心室心肌分割結果的DSC和CD值
由圖3知DRLSE、CPLSE依據梯度信息進行MR圖像分割時,存在對初始輪廓敏感的問題,本文算法同時依據梯度信息和圖像局部區域統計信息來驅動曲線演化,降低了算法對初始位置的依賴性。
由表1可以直觀看出各算法的時間復雜度,運用各種方法對左心室內外膜分割,U-Net的平均消耗時間較短,但該方法對網絡進行訓練時,不僅需要大量圖片標簽數據,還會消耗較長時間。結合表2、3 可知,相較CPLSE、DRLSE-shape 及本文方法,DRLSE 雖然平均時間最短但其分割效果最不理想。本文算法、CPLSE 及DRLSE-shape方法的運行時間相近,本文算法在運行速度方面表現并不突出,主要因為本文內膜與外膜算法在DRLSE 基礎上,兩者除了都增加新局部信息及輪廓內外區域灰度平均值的計算之外,內膜算法添加了曲率因子的計算,外膜算法添加了形狀約束項的計算,從而影響了內外膜算法的運行速度。
由圖4~圖7 與表2 知,對比 DRLSE、CPLSE 模型、U-Net網絡,本文算法加入圖像局部信息,對邊緣處有加強的作用,在處理弱邊緣、對比度低且有乳突肌影響的左心室圖像時,可以得到更高的DSC平均值,本文算法所得輪廓包圍的內部區域與專家手工分割所得輪廓包圍的內部區域更加重合。分割結果APD平均值相對于對比中的其他模型分割所得APD 平均值更低,本文算法所得的輪廓更加靠近專家手工分割輪廓,有更高的分割精度。
由圖8 與表3 知,深度學習方法對左心室MR 圖像分割時,可能會無法對圖像做出分割。在對DRLSE 算法加入本文先驗形狀約束項后,提高了DRLSE 方法的分割能力。本文外膜分割算法基于本文內膜分割效果的基礎上進行的,內膜分割結果優于DRLSE、CPLSE模型、U-Net網絡分割結果,并且在外膜分割方法中依靠了圖像差的局部信息,所以本文外膜分割效果是對比方法中最好。
由圖9 與表4 知,在對左心室整體的分割結果對比中,U-Net 網絡對左心室心肌整體分割的質心CD 平均值最低,更接近手工輪廓質心,但其DSC均值在對比方法中低于本文算法的DSC均值,且本文的質心CD值較傳統方法分割結果相較于傳統分割方法更接近手工輪廓的質心,本文算法有更好的分割效果。
針對運用DRLSE 模型分割存在弱邊界、對比度低的左心室圖像時不能得到很好分割效果,為了提高左心室內外膜的分割精度,本文提出改進的DRSLE 算法內膜分割算法和帶有形狀約束的外膜分割算法。通過區域擬合圖像與原圖像比較,提供局部區域信息,將改進的C-V 面積項與局部區域信息與保持凸性的方法添加到DRLSE算法中,以此驅動曲線進行演化。實驗表明,本文算法在對弱邊界、對比度低的左心室圖像進行內膜分割時,相對DRLSE、CPLSE算法、U-Net網絡有更高的分割精度;在對弱邊界外膜分割時,也取得相對穩定、較高的分割精度。但本文算法,需要手動給出初始輪廓,以及調節部分參數,在分割過程中仍依賴人工經驗。未來的工作將是,需求新的方法、減少變參數量、降低模型對人工依賴,以達到低成本人工成本同時達到更高的分割精度。