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

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

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

0　引言

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

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

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

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

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

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

(1)

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

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

3)計算對應(yīng)特征點新的位置，更新位置迭代方向及其對應(yīng)的偏差項。

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

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

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

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

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

(2)

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

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

(3)

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

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

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

(4)

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

Δx1=R0φ0+b0

(5)

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

(6)

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

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

(7)

(8)

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

圖1　SDM收斂過程Fig. 1　Convergence process of SDM

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

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

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

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

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

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

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

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

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

(9)

(10)

(11)

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

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

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

2　實驗及其結(jié)果分析

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

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

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

2.1　B-SIFT與SIFT實驗效果對比

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

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

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

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

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

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

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

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

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

3　結(jié)語

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

參考文獻(xiàn):

[1]朱文玲.超聲心動圖評價左心室功能[J].中國心血管雜志,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].長春師范大學(xué)學(xué)報(自然科學(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].電子測量與儀器學(xué)報,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]譚志國,劉滕沖,孫即祥,等.一種非剛性醫(yī)學(xué)圖像點匹配方法[J].武漢大學(xué)學(xué)報(理學(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].長沙:中南大學(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.)