劉維平,楊新
前列腺癌是我國男性中發(fā)病率呈增長趨勢的惡性腫瘤,早期診斷對于及時(shí)治療有重要的意義。醫(yī)院確診前列腺癌的主要手段是經(jīng)直腸的超聲活組織穿刺取樣檢驗(yàn),而結(jié)合多種成像模式的穿刺活檢能夠提高前列腺癌的檢出率。核磁共振成像對前列腺及其內(nèi)部結(jié)構(gòu)的成像效果遠(yuǎn)好于超聲圖像,便于醫(yī)生觀察以確定癌癥位置。若將核磁共振成像中疑似癌癥的位置映射到超聲圖像中作有針對性的穿刺,就需要對超聲和核磁共振成像的前列腺分別分割,然后通過配準(zhǔn)分割所得的輪廓使得超聲和核磁共振成像前列腺的切片對應(yīng),方便與醫(yī)生觀察和比較。因此,超聲前列腺分割對前列腺癌癥的診斷有著重要的意義。除此之外,其在前列腺特異抗原計(jì)算、超聲聚焦、超聲引導(dǎo)的前列腺短程放療等應(yīng)用領(lǐng)域,同樣具有極大的潛力。
近年來,基于先驗(yàn)知識的分割方法在超聲前列腺分割中受到了廣泛的關(guān)注。該方法的特點(diǎn)在于在傳統(tǒng)的分割方法中引入了先驗(yàn)知識,對于待分割目標(biāo)在復(fù)雜背景下產(chǎn)生的分割結(jié)果的不確定性,可依靠目標(biāo)的先驗(yàn)知識得到解決,即利用先驗(yàn)知識對目標(biāo)進(jìn)行恢復(fù)和重構(gòu)的過程。
Shen[1]提出處理局部形變的形狀模型,采用旋轉(zhuǎn)不變特征作為超聲圖像屬性驅(qū)動前列腺形狀模型的形變。為了避免一些優(yōu)化過程中出現(xiàn)的局部能量最小,Shen采用多分辨率的策略,在低分辨率中粗略分割,在高分辨率中精確分割。
Xu[2]提出利用超聲前列腺的先驗(yàn)直方圖和先驗(yàn)形狀模型實(shí)現(xiàn)分割。先驗(yàn)直方圖和形狀都是通過機(jī)器學(xué)習(xí)的方法從訓(xùn)練集中學(xué)習(xí)出來,通過一個(gè)手動調(diào)節(jié)的參數(shù)平衡這兩項(xiàng)對曲線演化的影響。
Yan[3]采用穩(wěn)定、正確的特征提取的結(jié)果來引導(dǎo)形狀模型的形變。在特征提取的時(shí)候?qū)μ卣鬟M(jìn)行初步的篩選,那些不好的或者不穩(wěn)定的特征提取結(jié)果是不對形狀模型的形變產(chǎn)生任何影響的。
本文采用基于先驗(yàn)形狀的活動輪廓模型對超聲前列腺圖像分割。將前列腺的形狀通過機(jī)器學(xué)習(xí)的方法求得形狀模型,并將其作為先驗(yàn)知識引入活動輪廓模型,使得曲線滿足先驗(yàn)形狀的約束。
活動輪廓模型是一種圖像分割方法,其主要思想是在圖像域定義一個(gè)曲線或曲面,并且在由曲線或曲面自身相關(guān)的內(nèi)力,以及由圖像數(shù)據(jù)定義的外力的作用下移動。通常情況下,外力的定義是基于邊緣特征或者區(qū)域特征的,它們分別采用了待分割圖像的局部和全局信息。為了兼顧全局和局部特征對演化曲線的影響,近些年來,一些學(xué)者提出對背景和目標(biāo)做局部統(tǒng)計(jì)的方法。最為代表性的思路是基于局部擬合方法的能量函數(shù),用以解決灰度分布不均勻的圖像分割問題。
Wang[4]提出基于局部高斯擬合的活動輪廓模型,對于圖像上的某點(diǎn),定義該點(diǎn)處的能量函數(shù)為
其中,Kσ(x)為截?cái)嗟母咚购撕瘮?shù),σ為方差,pi,x(I(y )),i= 1,2,為像素的灰度在區(qū)域Ω1,Ω2中的概率密度分布。點(diǎn)x處能量函數(shù)被限制在以該點(diǎn)為中心的某個(gè)窗口中,截?cái)喔咚购说姆讲顩Q定了窗口ρ的大小,并且,圖像中的點(diǎn)y對能量函數(shù)的貢獻(xiàn)隨著y與x距離的增大而減小,該截?cái)喔咚购撕瘮?shù)如下所示:
其中,c是歸一化的常量,通常情況下,令ρ=2σ。
在Wang提出的方法中,假設(shè)概率密度分布 pi,x(I(y ))滿足高斯分布,即
其中, ui( x)和 σi(x)2分別為點(diǎn)x附近第i個(gè)局部區(qū)域灰度的均值和方差。
整個(gè)目標(biāo)的邊界通過最優(yōu)化圖像域Ω中所有局部高斯擬合能量函數(shù)得到,即
其中,φ為水平集函數(shù),H(·)為 Heaviside函數(shù),公式(5)是公式(4)的水平集表示。對該能量函數(shù)求極值并結(jié)合梯度下降方法可得到演化方程和某點(diǎn)處演化曲線內(nèi)部和外部灰度值的局部擬合方程。該方法能夠有效地克服圖像亮度不一致的問題,并且在紋理圖像分割中有著較好的效果。
影響超聲前列腺圖像分割的主要原因有如下幾個(gè)方面:1)斑點(diǎn)噪聲是超聲圖像中常見的干擾,導(dǎo)致分割困難。2)在超聲前列腺區(qū)域中,灰度或者紋理分布是不均勻的,有些部位灰度級較低,而另外一些部位灰度值較高;在超聲前列腺外部,有些地方的灰度分布和前列腺區(qū)域內(nèi)部的分布一樣,故難以采用全局的紋理或者灰度特征區(qū)別前列腺的內(nèi)外部分。3)超聲偽影也會對前列腺分割產(chǎn)生嚴(yán)重的影響。因此,超聲圖像的分割迄今仍然是一個(gè)難題。
在超聲前列腺圖像分割中,單純采用基于邊緣或者區(qū)域的傳統(tǒng)活動輪廓模型的方法很難得到較好的分割結(jié)果,因此引入前列腺的先驗(yàn)性形狀作為傳統(tǒng)活動輪廓模型的約束。對于先驗(yàn)形狀約束的圖像分割方法,首先要將待分割目標(biāo)形狀的訓(xùn)練集配準(zhǔn),然后將配準(zhǔn)后的形狀采用機(jī)器學(xué)習(xí)的方法提取形狀特征并且對形狀建模,最后將形狀模型嵌入到傳統(tǒng)活動輪廓模型當(dāng)中。
主成分分析方法(Principle Component Analysis, PCA)常常用于提取形狀模型的形狀特征,而形狀形變能力主要由形狀模型形變參數(shù)的概率密度估計(jì)表示。當(dāng)目標(biāo)的形狀變化不是很大的情況下,可以假設(shè)形變參數(shù)的概率密度函數(shù)滿足高斯分布,而當(dāng)形變很大的情況下,就需要采用非參數(shù)的方法估計(jì)形變參數(shù)[5]。但是超聲前列腺形狀的形變不是嚴(yán)重,因此針對這個(gè)具體問題,采用高斯分布的形狀模型。
形狀配準(zhǔn)是機(jī)器視覺中的一個(gè)非常復(fù)雜的問題,其基本思想可以描述如下:給定目標(biāo)形狀、當(dāng)前要配準(zhǔn)的形狀和相似性度量,求解某變換使得變換后兩個(gè)形狀的相似性度量最小。Paragios[6]采用符號距離函數(shù)描述形狀,通過最優(yōu)化兩形狀的平方誤差和從而實(shí)現(xiàn)形狀的剛性配準(zhǔn)。本章采用這種方法將形狀訓(xùn)練集中的圖像配準(zhǔn)。
假設(shè)目標(biāo)形狀為ψT(X , Y ),當(dāng)前形狀為ψC(U , V ),目的是找到兩者間的剛性變換
其中,F(xiàn)(·)為前向映射;s,θ,Tx,Ty分別為從坐標(biāo)系(U, V)前向變換到坐標(biāo)系(X, Y)的尺度參數(shù),旋轉(zhuǎn)參數(shù),x軸方向的平移參數(shù),y軸方向的平移參數(shù),使得兩個(gè)形狀的平方誤差和最小,即最優(yōu)化如下能量函數(shù)
其中,Φ為坐標(biāo)系(X, Y)空間;B(·)為后向映射,即為 F(·)的逆映射
采用梯度下降的方法對公式(12)求最優(yōu),得到配準(zhǔn)參數(shù)的迭代方程如下
通過迭代公式(13)~(16)直至穩(wěn)定,可將兩個(gè)以符號距離函數(shù)表示的形狀配準(zhǔn)。在訓(xùn)練形狀模型的時(shí)候,給定形狀的訓(xùn)練集,選取其中任意一個(gè)形狀的作為目標(biāo)形狀,將其余的形狀與該形狀采用上述方法配準(zhǔn)即可得到配準(zhǔn)后的訓(xùn)練集。
其中,U的列向量描述了形狀模型變化的正交系,Λ是相應(yīng)奇異值的對角矩陣。新形狀可以用k個(gè)主分量和k維列向量系數(shù)γ計(jì)算如下:
其中,Uk是d×k維矩陣,包含了k個(gè)最大特征值對應(yīng)的U的列向量。
一個(gè)確定k的方法是選擇滿足公式
最小的k。假設(shè)參數(shù)向量γj滿足高斯分布,并且滿足-3≤γj≤3,γj是 γ的第j個(gè)分量,λj是(1 n) MMT的第 j個(gè)特征值。通過調(diào)整γ使形狀模型產(chǎn)生形變。
定義基于先驗(yàn)形狀模型的活動輪廓的能量函數(shù)如下:
其中,ETotal為總的能量函數(shù),ELGDF為局部高斯擬合的活動輪廓模型的能量函數(shù),EShape為形狀模型的能量函數(shù)。分割的目的是通過最優(yōu)化以上的總的能量函數(shù)來驅(qū)動水平集曲線演化。λ為權(quán)值系數(shù),用以平衡兩個(gè)能量項(xiàng)對曲線演化的影響。
在基于先驗(yàn)形狀的活動輪廓模型中,需要最優(yōu)化形狀模型與演化曲線的平方誤差和,使得形狀模型與演化曲線最匹配。因此,整個(gè)系統(tǒng)中有兩套參數(shù)需要被優(yōu)化:配準(zhǔn)參數(shù)和形變參數(shù)。構(gòu)造的能量函數(shù)如下所示:
其中,B(·)為后向映射,包括尺度s,旋轉(zhuǎn)θ,x和y方向上的平移參數(shù)分別為Tx和Ty。配準(zhǔn)參數(shù)的迭代方程如下
其中F(·)為前向映射。通過梯度下降方法得到形變參數(shù)的迭代方程:
其中,jγ為形狀模型的第j個(gè)形變參數(shù),j=1,…,k。
對能量函數(shù)求最優(yōu)解并且利用梯度先將方法可以得到如下公式
其中,δ(φ)為Dirac函數(shù)。公式 (23) 為曲線演化方程,公式 (24)~ (27) 為計(jì)算局部高斯擬合方法的參數(shù)計(jì)算公式。
其中,α=0.1為設(shè)定的參數(shù)。
給出先驗(yàn)形狀約束的活動輪廓模型的具體計(jì)算步驟:
1.給定演化方程的終止條件
2.先驗(yàn)形狀約束的前列腺分割算法步驟:
3.初始化活動輪廓和形狀模型;
4.將形狀模型與演化曲線配準(zhǔn),公式(17)~(20);
5.優(yōu)化形變參數(shù),公式(22);
6.迭代曲線演化方程,公式(23)~(27);
7.判斷曲線是否滿足終止條件(28),如果不滿足則轉(zhuǎn)到第4步;如果滿足則結(jié)束算法。
超聲數(shù)據(jù)庫有12個(gè)病例,每個(gè)病例有6個(gè)2維超聲前列腺圖像的切片,每個(gè)切片的大小為 700×400。這些超聲圖像數(shù)據(jù)來自于新華醫(yī)院超聲科,并且由醫(yī)生手動分割完畢。采用留一法,即每個(gè)病歷留下一個(gè)作為分割結(jié)果的測試集,其余的圖像的手動分割結(jié)果作為形狀模型的訓(xùn)練集。因此,在這個(gè)模型當(dāng)中,采用60幅前列腺的形狀作為形狀模型的訓(xùn)練集,12幅圖像作為算法的測試。
將配準(zhǔn)后的前列腺形狀作為訓(xùn)練集來訓(xùn)練前列腺的形狀模型。通過分別調(diào)整公式(13) 中的jγ,j=1,…,5,得到前列腺形狀模型的形變結(jié)果,如圖1所示:
圖1 前列腺形狀模型的形變
圖1可以看出,這種形狀模型能夠較好地模擬前列腺形變。當(dāng)然,有學(xué)者提出更為復(fù)雜的形狀模型,用以提取形變較大的形狀,例如基于核密度估計(jì)的形狀模型。但是在二維前列腺分割中,基于PCA和高斯分布的形狀模型足以應(yīng)付前列腺的形變,因此,針對特定問題采用了該種形狀模型。
在超聲前列腺的分割中,通過人為在前列腺圖像中選取上下左右四點(diǎn)擬合成一個(gè)閉合區(qū)域作為活動輪廓模型演化曲線的初始化,如圖2(a)所示。先將λ賦予較小的值,再將該值在每次迭代曲線演化方程的時(shí)候乘以一個(gè)大于 1的參數(shù)。這樣使得形狀約束逐漸作用于曲線演化,避免了能量函數(shù)的局部最小。圖2(b)~(d)為曲線演化和形狀模型配準(zhǔn)、形變的過程,可以看到曲線最終收斂到形狀模型附近,得到了較好的分割效果,如圖2(a) ~(d)所示:
圖2 前列腺分割過程
從另外一些超聲前列腺的分割結(jié)果中可以看出,盡管超聲圖像中的噪聲、偽影、灰度分布不均勻等問題很嚴(yán)重,提出的方法仍然可以較好地提取前列腺的邊緣,如圖3所示:
圖3 其它分割結(jié)果
為了對分割算法的準(zhǔn)確程度有一個(gè)客觀的評測,醫(yī)學(xué)圖像中,普遍采用的是和臨床醫(yī)生的手工分割的標(biāo)準(zhǔn)輪廓進(jìn)行比較。采用 Dice相似性系數(shù) (Dice Similarity Coefficient,DSC) 作為衡量標(biāo)準(zhǔn)。Dice相似性系數(shù)定義如下
其中,As和Ag分別是算法分割結(jié)果和醫(yī)生手動分割結(jié)果的面積。算法的分割結(jié)果與醫(yī)生手動分割結(jié)果的Dice相似性系數(shù),如表1所示:
表1 算法分割結(jié)果和醫(yī)生手動分割結(jié)果的Dice相似性系數(shù)
提出一種半監(jiān)督的、基于先驗(yàn)形狀約束和高斯擬合活動輪廓的分割模型,主要應(yīng)用于超聲前列腺圖像分割。通過將分割結(jié)果與醫(yī)生手動分割結(jié)果的對比,驗(yàn)證了所提出方法的有效性。
[1]Dinggang Shen, Yiqiang Zhan, Christos Davatzikos.Segmentation of prostate boundaries from ultrsound images using statistical shape model [J].IEEE Transactions on Medical Imaging, 2003, 22 (4): 539-551.
[2]Robert S.Xu, Oleg V.Michailovich, Igor Solovey, etc.A probability tracking approach to segmentation of ultrasound prostate images using weak shape priors [C].In Proceeding of SPIE Medical Imaging, 2010.
[3]Pingkun Yan, Sheng Xu, Baris Turbey, etc.Discrete deformable model guided by partial active shape model for TRUS image segmentation [J].IEEE Transactions on Biomedical Engineering, 2010, 57 (5): 1158-1166.
[4]Li Wang, Lei He, Arabinda Mishra, etc.Active contours driven by local Gaussian distribution fitting energy [J].Signal Processing, 2009, 89 (12): 2435-2447.
[5]Weiping Liu, Yanfeng Shang, Xin Yang, etc.A shape prior constraint for implicit active contours, Pattern Recognition Letters [J].2011, 32(15): 1937-1947
[6]Nikos Paragios, Mikael Rousson, Visvanathan Ramesh.Non-rigid registration using distance functions [J].Computer Vision and Image Understanding, 2003, 89(2-3): 142-165.