魏雨汐，伍岳慶，陶　攀,2，姚　宇

(1.中國(guó)科學(xué)院 成都計(jì)算機(jī)應(yīng)用研究所，成都 610041；　2.中國(guó)科學(xué)院大學(xué)，北京 100049)(*通信作者電子郵箱weiyuxi272@163.com)

0　引言

在醫(yī)學(xué)超聲圖中，經(jīng)食道的超聲心動(dòng)檢測(cè)是常見(jiàn)的超聲探測(cè)方式之一，也是人體左心室功能評(píng)價(jià)的重要參考[1]。醫(yī)學(xué)圖像分割是醫(yī)學(xué)圖像處理與分析領(lǐng)域復(fù)雜而關(guān)鍵的步驟，其目的是分割出醫(yī)學(xué)圖像中具有某些特殊含義的部分，并提取相關(guān)特征，為臨床診療和病理學(xué)研究提供可靠的依據(jù)，輔助醫(yī)生作出更為準(zhǔn)確的診斷[2]。

近年來(lái)，研究者們采用各類圖像分割算法應(yīng)用于醫(yī)學(xué)圖像處理領(lǐng)域，如：Pieciak[3]基于主動(dòng)輪廓模型(Active Contour Model, ACM)方法[4]對(duì)核磁共振成像(Magnetic Resonance Imaging, MRI)的左心室進(jìn)行分割；Hajiaghayi等[5]利用改進(jìn)的ACM方法將心臟MRI數(shù)據(jù)重建為3D模型，提出了一種三維分割的算法；O’Brien等[6]基于主動(dòng)形狀模型(Active Shape Model, ASM)[7]在3D圖像中對(duì)左心室進(jìn)行多輪廓分割；劉金清等[8]基于主動(dòng)表觀模型(Active Appearance Model, AAM)[9-10]方法對(duì)膀胱MRI醫(yī)學(xué)圖像進(jìn)行分割。

由于超聲心動(dòng)圖的復(fù)雜性和多變性，其分割過(guò)程需利用圖像中特征點(diǎn)位置的迭代變換，最終達(dá)到對(duì)目標(biāo)區(qū)域進(jìn)行輪廓提取的目的，其實(shí)質(zhì)上可視為特征點(diǎn)定位的問(wèn)題。而特征點(diǎn)的定位過(guò)程可以理解為非線性最小二乘近似問(wèn)題[11]。傳統(tǒng)的基于統(tǒng)計(jì)分布的ASM方法或AAM方法在建立形狀模型或表觀模型時(shí)使用的牛頓法或擬牛頓法在算法中計(jì)算代價(jià)較高，為此文獻(xiàn)[12]提出了一種監(jiān)督下降方法(Supervised Descent Method， SDM)來(lái)替代搜索過(guò)程中雅克比矩陣和海森矩陣的復(fù)雜計(jì)算，用于解決計(jì)算機(jī)視覺(jué)中的非線性最小二乘問(wèn)題，在圖像的特征點(diǎn)定位上取得了很好的效果。

SDM的目標(biāo)函數(shù)定義為最優(yōu)化兩個(gè)特征向量的相似程度，傳統(tǒng)的方法一般采用兩個(gè)特征向量的歐氏距離來(lái)作為特征點(diǎn)之間的相似性判定度量。而文獻(xiàn)[13]中已經(jīng)通過(guò)實(shí)驗(yàn)證明了基于巴氏系數(shù)(Bhattacharyya coefficient)的直方圖相似性判定相比傳統(tǒng)的歐氏距離的準(zhǔn)確性優(yōu)勢(shì)。同時(shí)，在超聲圖像處理領(lǐng)域，文獻(xiàn)[14]提出了一種基于巴氏系數(shù)的主動(dòng)輪廓模型圖像分割方法，提升了原方法的魯棒性。本文采用一種基于巴氏系數(shù)改進(jìn)的尺度不變特征變換(Bhattacharyya-Scale Invariant Feature Transform, B-SIFT)特征的SDM，并結(jié)合多尺度圖像金字塔建立了心內(nèi)膜與心外膜的聯(lián)合形狀模型，對(duì)左心室區(qū)域進(jìn)行準(zhǔn)確分割提取，取得了良好效果。

1　基于監(jiān)督下降的左心室超聲圖像分割方法

在特征點(diǎn)定位過(guò)程中，傳統(tǒng)的基于表觀模型的AAM算法的優(yōu)化目標(biāo)可定義為一個(gè)非線性最小二乘問(wèn)題：

(1)

其中：x=[x1,y1,x2,y2,…,xN,yN]T∈R2N×1為圖像d(x)中N個(gè)特征點(diǎn)組成的形狀向量；p為變換參數(shù)；f(x,p)表示對(duì)x進(jìn)行相應(yīng)的仿射變換[16]和非剛性變換[17]；d(f(x,p))表示變換后的圖像；Ua和ca為表觀模型和表觀參數(shù)。

若用傳統(tǒng)的牛頓法解決上述最優(yōu)化問(wèn)題，在每一步的迭代中都要計(jì)算相應(yīng)的海森矩陣與雅克比矩陣，那么將花費(fèi)巨大的計(jì)算代價(jià)[18]。而SDM則無(wú)需進(jìn)行上述復(fù)雜的計(jì)算，也無(wú)需建立表觀模型，僅需利用提取到的局部特征從訓(xùn)練數(shù)據(jù)學(xué)習(xí)一組梯度下降方向，并應(yīng)用到測(cè)試集圖像上。其核心過(guò)程主要涉及以下步驟：

3)計(jì)算對(duì)應(yīng)特征點(diǎn)新的位置，更新位置迭代方向及其對(duì)應(yīng)的偏差項(xiàng)。

4)重復(fù)2)、3)直至收斂，最終得到方向序列{Rk}和偏差項(xiàng)序列{bk}。

5)在測(cè)試過(guò)程中，將訓(xùn)練得到的{Rk}和{bk}作用于測(cè)試集圖像的初始位置，得到最終的特征點(diǎn)定位結(jié)果。

1.1　監(jiān)督下降方法(SDM)

SDM的優(yōu)化目標(biāo)和AAM有所不同，定義為特征函數(shù)之間的最優(yōu)化問(wèn)題，其目標(biāo)為最優(yōu)化式(2)所示的非線性最小二乘問(wèn)題：

f(x0+Δx)=‖h(d(x0+Δx))-h(d(x*))‖2

(2)

其中:d∈Rm×1表示一個(gè)具有m個(gè)像素點(diǎn)的圖像；d(x)∈Rp×1表示圖d中的p個(gè)特征點(diǎn)組成的形狀向量；h為特征提取函數(shù)，h(d(x))表示在圖像上對(duì)p個(gè)特征點(diǎn)的局部特征進(jìn)行提?。幌蛄縳*表示p個(gè)特征點(diǎn)的最佳位置，且在訓(xùn)練階段是已知的；x0表示初始形狀位置，通常由訓(xùn)練樣本的平均形狀表示。左心室特征點(diǎn)優(yōu)化的目標(biāo)即為當(dāng)前的特征點(diǎn)位置找到一個(gè)趨近于最優(yōu)值的位置偏移Δx。

為求解式(2)中的Δx，將其進(jìn)行二階泰勒展開(kāi)，得：

(3)

令φk=h(d(xk))，在式(3)對(duì)Δx求導(dǎo)并令其為0得到Δx的第一個(gè)更新：

Δx1=-H(x0)-1Jf(x0)=

-2H(x0)-1Jh(x0)T(φ0-φ*)

(4)

其中：H(x0)為函數(shù)f在x0處的海森矩陣；Jf(x0)為函數(shù)f在x0處的雅克比矩陣。在SDM中，可通過(guò)監(jiān)督學(xué)習(xí)找到表示迭代方向的R0和表示偏差項(xiàng)的b0替代式(4)中的-2H(x0)-1Jh(x0)T與2H(x0)-1Jh(x0)Tφ*，由此便可構(gòu)造出Δx1的求解式如式(5)所示：

Δx1=R0φ0+b0

(5)

所以最初最優(yōu)化問(wèn)題轉(zhuǎn)化為求解R0和b0，其求解方式由式(6)給出，上標(biāo)i表示訓(xùn)練樣本中的第i幅圖像：

(6)

而對(duì)于非線性的最小二乘問(wèn)題，需要多次迭代才能得到最優(yōu)結(jié)果。對(duì)于每一次迭代都會(huì)產(chǎn)生一個(gè)Rk和bk，便組成一個(gè)方向序列{Rk}和一個(gè)偏差項(xiàng)序列{bk}。一般地，對(duì)于第k次迭代，結(jié)合前一步得到的Rk-1和bk-1，有：

xk=xk-1+Rk-1φk-1+bk-1

(7)

(8)

在測(cè)試階段，給定待測(cè)試圖像的初始位置x0，將訓(xùn)練得到的{Rk}和{bk}作用于x0上，以得到最終的真實(shí)位置估計(jì)。圖1簡(jiǎn)要描述了SDM從初始位置到最佳位置的大致過(guò)程。

圖1　SDM收斂過(guò)程Fig. 1　Convergence process of SDM

由式(2)不難看出，SDM的核心在于最優(yōu)化兩個(gè)特征向量的歐氏距離，所以特征提取函數(shù)的選取對(duì)算法迭代過(guò)程中{Rk}和{bk}的生成有直接影響。

1.2　基于SDM的改進(jìn)算法

由于圖像特征提取函數(shù)的選取會(huì)影響SDM的最終分割效果，因此對(duì)左心室超聲圖像準(zhǔn)確分割的實(shí)現(xiàn)有著重要意義。由于受光照程度不同、心內(nèi)組織密度各異等客觀影響，傳統(tǒng)SDM采用的SIFT特征提取函數(shù)在左心室輪廓定位上取得的效果欠佳。為弱化上述客觀影響，本文采用B-SIFT特征，其效果比傳統(tǒng)的SIFT特征要好；同時(shí)，在模型中引入了多尺度圖像金字塔模型來(lái)進(jìn)一步提升最終左心室輪廓分割效果。

1.2.1基于SIFT特征的優(yōu)化特征B-SIFT

SIFT特征是一種局部特征檢測(cè)方法，廣泛應(yīng)用于多類計(jì)算機(jī)視覺(jué)領(lǐng)域。其優(yōu)點(diǎn)在于：具有尺度、旋轉(zhuǎn)不變性，對(duì)噪聲、圖像亮度變換有很好的抵抗力，所蘊(yùn)含的局部特征信息量豐富等。

對(duì)于SIFT特征向量的提取，在初始形狀向量x0中的每一個(gè)特征點(diǎn)周圍建立16×16的像素框，計(jì)算其中每一個(gè)像素的梯度；然后每4×4的像素上計(jì)算8個(gè)方向的梯度方向直方圖如圖2(a)所示；繪制每個(gè)梯度方向的累加值，如圖2(b)所示。于是便得到4×4×8=128維的SIFT特征向量(又稱描述子)。

圖2　SIFT特征建立Fig. 2　Construction of SIFT features

本文采用基于巴氏系數(shù)的B-SIFT特征替代原SIFT特征以輔助特征點(diǎn)定位，取得了較好的效果。其構(gòu)造方法如下：

對(duì)于兩個(gè)n維單位特征向量x、y，定義x與y的巴氏系數(shù)為：

(9)

(10)

(11)

1.2.2多尺度圖像金字塔模型建立

在左心室超聲圖的探測(cè)過(guò)程中，受左心室的形狀、大小、心尖朝向等客觀因素的影響，在對(duì)訓(xùn)練集或測(cè)試集中的圖像統(tǒng)一歸一化后，可能會(huì)造成局部紋理信息丟失，影響最終結(jié)果的準(zhǔn)確性。于是，本文采用了基于多尺度圖像金字塔模型的方式，將訓(xùn)練集中的圖像歸一至多個(gè)分辨率上，在低分辨率圖像上進(jìn)行步長(zhǎng)較大的迭代，在高分辨率圖像上進(jìn)行步長(zhǎng)較小的迭代，使特征點(diǎn)更為精確地貼近真實(shí)值。分辨率從上至下依次遞增，完成每一層的特征點(diǎn)匹配后，轉(zhuǎn)至下一層。本文采用的四層金字塔結(jié)構(gòu)如圖3所示。

圖3　多分辨率金字塔模型Fig. 3　Multi-resolution pyramid model

2　實(shí)驗(yàn)及其結(jié)果分析

本文實(shí)驗(yàn)數(shù)據(jù)來(lái)自于華西醫(yī)院專家標(biāo)注的心動(dòng)周期不同階段的106張超聲心動(dòng)圖，來(lái)自于25位不同的患者，提取其中81張作為訓(xùn)練集，25張作為測(cè)試集。在特征點(diǎn)標(biāo)注方面，選取左心室心內(nèi)膜和心外膜作為參考輪廓，分別標(biāo)注17個(gè)特征點(diǎn)，其中0～16號(hào)特征點(diǎn)為心外膜，17到33號(hào)特征點(diǎn)為心內(nèi)膜，如圖4所示。

圖4　手工標(biāo)注的特征點(diǎn)Fig. 4　Mannually labeled feature points

實(shí)驗(yàn)將本文采用的方法對(duì)測(cè)試集的經(jīng)食道的超聲心動(dòng)圖像進(jìn)行了特征點(diǎn)定位操作，分割出心內(nèi)膜與心外膜的輪廓，并與傳統(tǒng)的SDM進(jìn)行了對(duì)比，描述如下。

2.1　B-SIFT與SIFT實(shí)驗(yàn)效果對(duì)比

在同時(shí)基于4層多尺度金字塔模型的基礎(chǔ)上，從上至下縮放尺度依次設(shè)為0.25、 0.5、 1.0和2.0，每一層迭代次數(shù)為2，特征提取函數(shù)分別設(shè)為B-SIFT和SIFT，實(shí)驗(yàn)結(jié)果如圖5所示，圖5(a)為采用SIFT特征的特征點(diǎn)定位輪廓，圖5(b)為采用B-SIFT特征的特征點(diǎn)定位輪廓。其中黑色特征點(diǎn)為真實(shí)值位置，白色特征點(diǎn)為實(shí)驗(yàn)結(jié)果位置。觀察兩者的心室輪廓的高曲率區(qū)域，發(fā)現(xiàn)基于B-SIFT特征的實(shí)驗(yàn)結(jié)果比SIFT更加接近真實(shí)值。

圖5　利用SIFT與B-SIFT特征進(jìn)行特征點(diǎn)定位的結(jié)果對(duì)比Fig. 5　Comparison of feature point location results using B-SIFT and SIFT features

2.2　構(gòu)建多尺度圖像金字塔后改進(jìn)效果

本實(shí)驗(yàn)構(gòu)建了4層多分辨率金字塔迭代模型，在采用相同特征提取函數(shù)的情況下，其相對(duì)于沒(méi)有金字塔結(jié)構(gòu)的SDM模型在精準(zhǔn)度上有較大提升。選取測(cè)試集中的一張圖像作演示：圖6(a)為初始化特征點(diǎn)定位x0；圖6(b)為非金字塔結(jié)構(gòu)的監(jiān)督下降方法SDM的分割效果，其只在單一尺度(scale=1.0)上進(jìn)行迭代；圖6(c)是基于4層金字塔模型的特征點(diǎn)定位效果。其中黑色特征點(diǎn)為真實(shí)值位置，白色特征點(diǎn)為實(shí)驗(yàn)迭代結(jié)束位置?？梢钥闯觯夯诙喑叨葓D像金字塔模型的特征點(diǎn)定位更接近真實(shí)輪廓。

綜上所述，結(jié)合多尺度金字塔模型和B-SIFT特征提取函數(shù)的左心室輪廓分割方法在精準(zhǔn)度上要優(yōu)于傳統(tǒng)的SDM特征點(diǎn)標(biāo)注方法。本文方法與傳統(tǒng)SDM的累計(jì)誤差分布圖如圖7所示，其中橫坐標(biāo)為測(cè)試集序列號(hào)，縱坐標(biāo)為其誤差累計(jì)。從圖7可以看出，基于B-SIFT特征提取函數(shù)的方法效果優(yōu)于傳統(tǒng)SDM。

圖6　多分辨率金字塔優(yōu)化效果Fig. 6　Optimization result of multi-resolution pyramid

圖7　傳統(tǒng)SDM和本文方法的累計(jì)誤差分布對(duì)比Fig. 7　Cumulative error distribution comparison of the traditional SDM and the proposed method

將傳統(tǒng)SDM與本文方法作用于測(cè)試集中，得到的平均誤差結(jié)果如表1所示?？梢钥闯觯贐-SIFT特征及多尺度圖像金字塔的分割方法將平均誤差從0.055 73降低至0.029 30，準(zhǔn)確度提升了47.43%。

表1　傳統(tǒng)SDM與本文方法的平均誤差對(duì)比Tab. 1　Average error comparison of the traditional SDM and the proposed method

3　結(jié)語(yǔ)

本文基于SDM對(duì)經(jīng)食道超聲心動(dòng)圖中左心室輪廓進(jìn)行了分割，針對(duì)傳統(tǒng)方法在單一尺度上定位的精確性問(wèn)題，建立了4層金字塔特征點(diǎn)定位搜索模型，并結(jié)合基于巴氏系數(shù)改進(jìn)的尺度不變特征變換(B-SIFT)特征提取函數(shù)，進(jìn)一步提高了算法的準(zhǔn)確性,使最終的特征點(diǎn)定位輸出更加接近真實(shí)輪廓。在下一步研究中，計(jì)劃建立左心室特征點(diǎn)標(biāo)注圖像庫(kù)，增加標(biāo)注圖像的樣本數(shù)量，利用局部特征描述子匹配，將算法應(yīng)用于動(dòng)態(tài)超聲視頻追蹤上。

參考文獻(xiàn):

[1]朱文玲.超聲心動(dòng)圖評(píng)價(jià)左心室功能[J].中國(guó)心血管雜志,2008,13(4):241-243. (ZHU W L. Evaluation of left ventricular function by echocardiography[J]. Chinese Journal of Cardiovascular Medicine, 2008, 13(4):241-243.)

[2]黃文博,燕楊,王云吉.醫(yī)學(xué)圖像分割方法綜述[J].長(zhǎng)春師范大學(xué)學(xué)報(bào)(自然科學(xué)版),2013,32(2):22-25. (WANG W B, YAN Y, WANG Y J. Survey of medical image segmentation methods [J]. Journal of Changchun Normal University (Natural Science), 2013, 32(2): 22-25.)

[3]PIECIAK T. Segmentation of the left ventricle using active contour method with gradient vector flow forces in short-axis MRI [M]// Information Technologies in Biomedicine, LNCS 7339. Berlin: Springer, 2012: 24-35.

[4]DAVATZIKOS C, PRINCE J L. An active contour model for mapping the cortex [J]. IEEE Transactions on Medical Imaging, 1995, 14(1): 65-80.

[5]HAJIAGHAYI M, GROVES E M, JAFARKHANI H, et al. A 3-D active contour method for automated segmentation of the left ventricle from magnetic resonance images [J]. IEEE Transactions on Biomedical Engineering, 2016, 64(1): 134-144.

[6]O’BRIEN S, GHITA O, WHELAN P F. Segmenting the left ventricle in 3D using a coupled ASM and a learned non-rigid spatial model [C/OL]// 3D Segmentation in the Clinic: A Grand Challenge III [Workshop], MICCAI 2009: Proceedings of the 12th International Conference on Medical Image Computing and Computer Assisted Intervention [2017- 03- 16]. http://doras.dcu.ie/18626/1/whelan_2009_52.pdf.

[7]COOTES T F, TAYLOR C J, COOPER D H, et al. Active shape models—their training and application [J]. Computer Vision and Image Understanding, 1995, 61(1): 38-59.

[8]劉金清,吳慶祥.一種AAM模型的MRI醫(yī)學(xué)圖像分割算法[J].電子測(cè)量與儀器學(xué)報(bào),2013,27(3):236-240. (LIU J Q, WU Q X. Algorithm of AAM model for MRI medical image segmentation [J]. Journal of Electronic Measurement and Instrument, 2013, 27(3): 236-240.)

[9]COOTES T F, EDWARDS G J, TAYLOR C J, et al. Active appearance models [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2001, 23(6): 681-685.

[10]MATTHEWS I, BAKER S. Active appearance models revisited [J]. International Journal of Computer Vision, 2004, 60(2): 135-164.

[11]XIONG X. Supervised descent method [C]// Proceedings of the 2015 IEEE Computer Vision and Pattern Recognition. Piscataway, NJ: IEEE, 2015: 2664-2673.

[12]XIONG X, DE LA TORRE F. Supervised descent method and its applications to face alignment [C]// CVPR 2013: Proceedings of the 2013 IEEE Conference on Computer Vision and Pattern Recognition. Washington, DC: IEEE Computer Society, 2013: 532-539.

[13]ARANDJELOVIC R. Three things everyone should know to improve object retrieval [C]// CVPR 2012: Proceedings of the 2012 IEEE Conference on Computer Vision and Pattern Recognition. Washington, DC: IEEE Computer Society, 2012: 2911-2918.

[14]YUAN J. Active contour driven by region-scalable fitting and local Bhattacharyya distance energies for ultrasound image segmentation [J]. IET Image Processing, 2013, 6(8): 1075-1083.

[15]LOWE D G. Object recognition from local scale-invariant features [C]// ICCV ’99： Proceedings of the Seventh IEEE International Conference on Computer Vision. Washington, DC: IEEE Computer Society, 2002, 2: 1150.

[16]KENDALL R A, DUNNING H T, Jr., HARRISON R J. Electron affinities of the first-row atoms revisited. Systematic basis sets and wave functions [J]. Journal of Chemical Physics, 1992, 96: 6796-6806.

[17]譚志國(guó),劉滕沖,孫即祥,等.一種非剛性醫(yī)學(xué)圖像點(diǎn)匹配方法[J].武漢大學(xué)學(xué)報(bào)(理學(xué)版),2007,53(6):759-764. (TAN Z G, LIU T C, SUN J X, et al. A non-rigid medical image point matching method [J]. Journal of Wuhan University (Natural Science Edition), 2007, 53(6): 759-764.)

[18]鄧梁.基于ASM與AAM的人臉特征定位與匹配算法研究[D].長(zhǎng)沙:中南大學(xué),2009:15-20. (DENG L. Research on face feature positioning and matching algorithm based on ASM and AAM [D]. Changsha: Central South University, 2009: 15-20.)