乳腺腫瘤超聲圖像的反應(yīng)擴(kuò)散水平集分割

陳 帥 張 麒

(上海大學(xué)通信與信息工程學(xué)院，上海 200072)

乳腺腫瘤超聲圖像的反應(yīng)擴(kuò)散水平集分割

陳 帥 張 麒

(上海大學(xué)通信與信息工程學(xué)院，上海 200072)

針對乳腺腫瘤超聲圖像分割，提出一種改進(jìn)的反應(yīng)擴(kuò)散(RD)水平集分割算法。先使用Gabor各向異性擴(kuò)散模型進(jìn)行濾波，由此構(gòu)造邊界停止函數(shù)；再將該函數(shù)融入RD水平集演化方程，以控制曲線的演化得到乳腺腫瘤的邊界。采用該方法和傳統(tǒng)RD方法對77例病人的111幅乳腺超聲圖像進(jìn)行分割實(shí)驗(yàn)，分割準(zhǔn)確率分別為98.5%和98.0%，真陽性率分別為88.2%和82.7%，與金標(biāo)準(zhǔn)之間的均方根誤差分別為3.6和4.6像素。結(jié)果表明，該改進(jìn)算法可獲得更加準(zhǔn)確的乳腺腫瘤分割結(jié)果。

邊界停止函數(shù) 反應(yīng)擴(kuò)散 水平集演化 圖像分割 乳腺超聲圖像

0 引言

乳腺癌是女性面臨的最常見的惡性腫瘤?；诔晥D像的計(jì)算機(jī)輔助診斷(computer aided diagnosis, CAD)系統(tǒng)是提高乳腺腫瘤診斷客觀性的重要手段。腫瘤圖像分割是CAD的前提條件和關(guān)鍵步驟，國內(nèi)外對此提出多種方法，如灰度閾值結(jié)合動態(tài)規(guī)劃法[1-2]、一致性直方圖法[3]、活動輪廓模型[4-6]。但由于超聲圖像中斑點(diǎn)噪聲嚴(yán)重、對比度低和乳腺腫瘤形狀復(fù)雜多變等因素[5]，現(xiàn)有方法難以取得理想的分割結(jié)果。水平集方法屬于幾何活動輪廓模型，具有適應(yīng)拓?fù)浣Y(jié)構(gòu)變化、穩(wěn)定性好等優(yōu)點(diǎn)，廣泛用于醫(yī)學(xué)圖像分割[7]，成為研究熱點(diǎn)。本文結(jié)合近年提出的反應(yīng)擴(kuò)散(reaction diffusion,RD)水平集[13]和Gabor各向異性擴(kuò)散模型[16]，提出改進(jìn)的水平集算法，以提高輪廓演化的效率和細(xì)節(jié)分割的準(zhǔn)確性，實(shí)現(xiàn)乳腺腫瘤的精確快速分割。

1 算法設(shè)計(jì)與實(shí)現(xiàn)

1.1 RD水平集原理

原始的水平集方法有兩種形式[13]，即基于偏微分方程的水平集方法和基于變分的水平集方法，方程分別如式(1)和式(2)所示。

(1)

(2)

式中：φ為水平集函數(shù)(levelsetfunction,LSF)，φ0為初始的水平集函數(shù)；F為邊界力函數(shù)，控制著目標(biāo)輪廓曲線的演化；為梯度算子；δ(·)為Dirac函數(shù)。這兩種水平集方法在演化過程中，在過平滑或者過陡峭邊緣容易出現(xiàn)嚴(yán)重的數(shù)值誤差，因此需要重新初始化為一個有符號的距離函數(shù)。

反應(yīng)擴(kuò)散方程最初用于檢測動物皮毛的化學(xué)機(jī)制[15]，它包括了反應(yīng)和擴(kuò)散兩個過程，可用于描述紋理分析、自然圖像建模與相變建模等領(lǐng)域的動態(tài)過程[17]。Zhang等[13]受相變理論啟發(fā)，提出在傳統(tǒng)的水平集演化方程中增加一個擴(kuò)散項(xiàng)構(gòu)造成RD方程，以使水平集演化保持穩(wěn)定且消除重新初始化的過程。該方法不僅適用于變分水平集方法，而且適用于基于偏微分方程的水平集方法。

RD演化方程如下：

(3)

式中：ε為一個小的正常數(shù)；Δ為拉普拉斯算子；φ為LSF，φ0為初始的水平集函數(shù)；L(φ)為演化方程，在變分水平集方法中，L(φ)=-Fδ(φ)；在基于偏微分方程的水平集方法中，L(φ)= -F|φ|。

式(3)包含了兩個動態(tài)過程：εΔφ為擴(kuò)散項(xiàng)，在每個分割區(qū)域內(nèi)逐漸把LSF正則化為一個分段常數(shù)；-ε-1L(φ)為反應(yīng)項(xiàng)，使RD方程趨于穩(wěn)定，直至L(φ)=0，得到分割區(qū)域。在傳統(tǒng)的LSM中，由于無此擴(kuò)散項(xiàng)，因此，不得不通過增加一個額外的過程(如重新初始化)來正則化LSF。

Zhang等[13]在基于相變的VanderWaals-Cahn-Hilliard理論上分析了當(dāng)ε0+時(shí)，基于變分水平集方法的RD方程平衡解，構(gòu)成與偏微分水平集方法的統(tǒng)一框架，從而，式(3)可看作傳統(tǒng)水平集方法的一般形式。采用了兩步分裂法(two-stepsplittingmethod,TSSM)來迭代地求解RD方程。

1.2 引入GAD的RD改進(jìn)算法

在大噪聲情況下，傳統(tǒng)的RD算法存在分割細(xì)節(jié)不準(zhǔn)確、演化效率低的問題。本文圍繞傳統(tǒng)RD方程存在的問題，在RD方程基礎(chǔ)上增加一個邊界停止函數(shù)項(xiàng)，構(gòu)建新的演化方程：

(4)

式中：g為邊界停止函數(shù)；α為控制邊界停止函數(shù)權(quán)重的系數(shù)。

在變分水平集模型下L(φ)=-Fδ(φ)，邊界力函數(shù)F=μdiv(φ/|φ |)-λ1(I-c1)2+λ2(I-c2)2，其中，μ、λ1、λ2為系數(shù)常量；I為目標(biāo)分割圖像；c1、c2分別為曲線內(nèi)、外部常量值。在進(jìn)行目標(biāo)分割演化時(shí)。通過曲線內(nèi)外力值大小的改變，控制曲線向內(nèi)或者向外演化。當(dāng)F=0時(shí)，表示內(nèi)外力達(dá)到平衡，此時(shí)得到最終目標(biāo)輪廓。增加邊界停止函數(shù)后，邊界力函數(shù)F =μdiv(φ/|φ |)-λ1(I-c1)2+λ2(I-c2)2+αg；通過停止函數(shù)g降低圖像噪聲對力的影響，同時(shí)改變系數(shù)α，增加弱邊緣力對演化方程的影響，從而提高分割的準(zhǔn)確性。

在弱噪聲情況下，邊界停止函數(shù)可表示為[6]：

(5)

式中：Gσ為方差為σ的高斯函數(shù)；I為原始圖像；*表示卷積運(yùn)算。

由式(5)可以看出，當(dāng)圖像梯度較大時(shí)，g的取值趨于0，對邊界力函數(shù)影響相對較弱，曲線演化；反之，當(dāng)梯度較小時(shí)，g的取值趨于1，邊界力函數(shù)影響相對增強(qiáng)，控制曲線停止演化，此時(shí)得到目標(biāo)的輪廓。

由于超聲圖像中存在很強(qiáng)的乘性斑點(diǎn)噪聲，為克服其對演化曲線的影響，通常需要使用方差很大的高斯濾波器；但濾波器在濾除噪聲的同時(shí)也使邊緣變模糊，容易導(dǎo)致在弱邊界處產(chǎn)生“邊界泄漏”(即邊界不連續(xù)、斷裂)現(xiàn)象。Gabor各向異性擴(kuò)散(gabor-basedanisotropicdiffusion,GAD)模型是一種對斑點(diǎn)噪聲較為有效的濾波算法，它不僅能夠在不同程度斑點(diǎn)噪聲下保持較好的降噪魯棒性，而且具有優(yōu)越的邊緣保持特性[16]。因此，本文在RD方程中引入GAD模型，以提高RD水平集分割的準(zhǔn)確性和效率。引入GAD的邊界停止函數(shù)表示為：

(6)

式中：GAD(I)為經(jīng)過GAD濾波后的圖像。

1.3 算法實(shí)現(xiàn)

文獻(xiàn)[18]利用TTSM算法實(shí)現(xiàn)反應(yīng)擴(kuò)散，首先由反應(yīng)方程生成一個0和1的二值函數(shù);然后將擴(kuò)散方程應(yīng)用到二值函數(shù)中，得到曲率運(yùn)動;最后，得到TTSM算法實(shí)現(xiàn)RD方程。該方程可以表述為：

① 解反應(yīng)項(xiàng)φt=-ε-1L(φ0)，φ(x,t=0)=φn，在某時(shí)間Tr達(dá)到中間解,表示為φ(n+1)/2=φ(x,Tr)。

② 解擴(kuò)散項(xiàng)φt=εΔφ，φ(x,t=0)=φ(n+1)/2，在某時(shí)間Td，最后水平集表示為φn+1=φ(x,Td)。

在以上2個步驟中，通過選擇小Tr和Td，φ(n+1)/2和φn+1分別離散地近似為φ(n+1)/2=φn+Δt1[-ε-1L′×(φn)]，φn+1=φn+1/2+ Δt2(εΔφn+1/2)，其中時(shí)間步長Δt1和Δt2分別表示時(shí)間Tr和Td。將參數(shù)ε并入到時(shí)間步長Δt1和Δt2中，如Δt1Δt1(-ε-1)和Δt2Δt2ε。因此，類似文獻(xiàn)[18]，只要選擇兩個能夠保持?jǐn)?shù)值穩(wěn)定的時(shí)間步長Δt1和Δt2即可。

用n表示迭代次數(shù)，φn表示第n次迭代的LSF，基于水平集演化的RD算法實(shí)現(xiàn)過程概括如下。

① 初始化：φn=φ0，n=0；

② 計(jì)算φ(n+1)/2：

φ(n+1)/2=φn-Δt1L′(φn)

(7)

式中：L′(φn) =L(φn)+αgδ(φn);Δt1為反應(yīng)項(xiàng)的時(shí)間步長，影響著演化的穩(wěn)定性。

③ 計(jì)算φn+1：

φn+1=φn+Δt2Δφn

(8)

式中：φn=φn+1/2；Δt2為擴(kuò)散項(xiàng)的時(shí)間步長，取值范圍為0～0.25。

④ 如φn+1滿足穩(wěn)定條件，停止演化；反之，n=n+1，返回步驟②。

2 實(shí)驗(yàn)結(jié)果與討論

為驗(yàn)證算法的有效性，采集來自77例患者(年齡17～79歲)的111幅(79幅良性；32幅惡性)乳腺腫瘤B型超聲圖像進(jìn)行分割實(shí)驗(yàn)。超聲儀型號為法國SuperSonic Aixplorer，超聲探頭頻率范圍4～15 MHz。

為比較分割算法的性能，將本文提出的改進(jìn)RD算法與Zhang提出的傳統(tǒng)RD 算法及手工分割方法進(jìn)行實(shí)驗(yàn)對比。由于本文算法是在RD算法的基礎(chǔ)上提出的，為了客觀地對比算法的特性，將兩者共有的參數(shù)設(shè)置為相同的值，并選擇相同的初始輪廓(手工描記4～6個點(diǎn)連成多邊形)，避免初始輪廓對分割結(jié)果的影響。

圖1和圖2分別為兩種算法對良性乳腺腫瘤超聲圖像和惡性乳腺腫瘤超聲圖像的分割結(jié)果對比。

圖1 良性乳腺腫瘤超聲圖像分割結(jié)果

圖2 惡性乳腺腫瘤超聲圖像分割結(jié)果

圖1是一個形狀較為規(guī)則的良性乳腺腫瘤的圖像。乳腺腫瘤內(nèi)部灰度不均勻，1點(diǎn)鐘方向呈現(xiàn)出弱邊緣。由分割結(jié)果可看出，在相同參數(shù)情況下，本文算法在弱邊緣和大噪聲下體現(xiàn)出較好的分割特性。圖2是一個形狀相對復(fù)雜、邊緣凹凸不平的惡性乳腺腫瘤的圖像。在相同參數(shù)情況下，本文算法相對傳統(tǒng)RD方法能更有效提取9點(diǎn)鐘及2點(diǎn)鐘方向凸出的邊緣，與手工分割結(jié)果更加相符。

本文采用真陽性率(true positive rate,TPR)、假陽性率(false positive rate,FPR)、準(zhǔn)確率(accuracy,Acc)、特異性(specificity,Speci)、均方根誤差(root mean square error,RMSE)共5個定量指標(biāo)，并以醫(yī)生手工分割作為金標(biāo)準(zhǔn)，來度量算法的分割質(zhì)量。計(jì)算公式如下：

(9)

TPR代表實(shí)際分割腫瘤區(qū)域的準(zhǔn)確率，F(xiàn)PR是非腫瘤區(qū)域被算法分作腫瘤區(qū)域的概率，Acc表示綜合分隔腫瘤區(qū)域的準(zhǔn)確率，Speci表示非腫瘤區(qū)域被算法正確判別的概率。TPR、Acc、Speci度量了分割算法得到區(qū)域與金標(biāo)準(zhǔn)區(qū)域之間的值，取值在0～1之間，值越大表示分割越準(zhǔn)確。RMSE則反映了分割輪廓坐標(biāo)點(diǎn)與金標(biāo)準(zhǔn)輪廓坐標(biāo)點(diǎn)之間的平均距離，值越小表示分割越準(zhǔn)確。圖3給出GT、TP、TN、FP、FN等參數(shù)的物理意義。

對111幅乳腺腫瘤圖像分別使用兩種算法進(jìn)行分割，得到的量化指標(biāo)平均值和標(biāo)準(zhǔn)差如表1所示。從表1可知，RD和本文算法的Acc和TPR分別為98.0%1.0%和98.5%0.7%，82.7%7.7%和88.2%5.7%，RMSE分別為4.6341.373pixels和3.5970.813pixels。本文算法在Acc、TPR和RMSE三個指標(biāo)上均優(yōu)于傳統(tǒng)RD算法。值得注意的是，本文算法的FPR和Speci指標(biāo)略遜于RD算法，但這兩個指標(biāo)僅在分割灰度均勻的目標(biāo)時(shí)具有較高參考價(jià)值，在分割不均勻的乳腺腫瘤時(shí)并不能很準(zhǔn)確地表征分割性能。實(shí)際乳腺腫瘤病灶內(nèi)部灰度值的分布往往是很不均勻的，當(dāng)初始輪廓定義在病灶內(nèi)部、演化時(shí)被病灶內(nèi)偽邊緣所吸引而未到達(dá)真實(shí)輪廓，Speci和FPR也能達(dá)到比較好的值；比如在分割結(jié)果完全位于金標(biāo)準(zhǔn)內(nèi)部時(shí)，Speci=1，F(xiàn)PR=0。因此，本文將Acc、TPR和RMSE視為主要指標(biāo)，F(xiàn)PR和Speci視為次要指標(biāo)。綜合考慮所有指標(biāo)，可知本文算法對乳腺腫瘤超聲圖像的分割精度優(yōu)于傳統(tǒng)RD算法。

圖3 GT、TP、TN、FP、FN等參數(shù)的物理意義

表1 RD算法和本文算法對乳腺超聲圖像的分割定量指標(biāo)

Tab.1 The quantitative indicators of breast ultrasound image segmentation(RD vs. improved RD)

指標(biāo)RD算法均值標(biāo)準(zhǔn)差本文算法均值標(biāo)準(zhǔn)差TPR0．8270．0220．8820．033FPR0．0150．0100．0280．007Acc0．9800．0100．9850．032Speci0．9971．3730．9920．813RMSE/pixel4．6343．597

3 結(jié)束語

本文提出了一種改進(jìn)的RD水平集乳腺腫瘤超聲圖像分割方法。與傳統(tǒng)的RD水平集方法不同的是，本文的方法利用具有斑點(diǎn)噪聲強(qiáng)抑制力的GAD算法進(jìn)行濾波處理，得到邊界停止函數(shù)，并將該函數(shù)作為RD能量函數(shù)中的一個增加項(xiàng)，使得改進(jìn)的RD算法具有更好的演化效率以及捕捉弱邊緣的能力。此外，本文的改進(jìn)RD模型是通過水平集方法實(shí)現(xiàn)的，繼承了水平集良好的拓?fù)渥赃m應(yīng)性優(yōu)勢，可以較好地分割出形狀復(fù)雜多變的乳腺腫瘤。實(shí)驗(yàn)結(jié)果表明，本文算法在乳腺腫瘤超聲圖像分割中取得良好效果。尤其在應(yīng)用于乳腺超聲圖像中目標(biāo)區(qū)域灰度不均勻、存在凹凸不平邊界及弱邊界情況下，與傳統(tǒng)RD方法相比，本文算法顯示出優(yōu)越性能。下一步將考慮如何自動提取初始輪廓以進(jìn)一步提高分割的客觀性與自動化水平。

[1] 汪源源，沈嘉琳，王涌，等.基于灰度閾值分割和動態(tài)規(guī)劃的超聲圖像乳腺腫瘤邊緣提取[J].航天醫(yī)學(xué)與醫(yī)學(xué)工程，2005,18(4):281-286.

[2] 汪源源，沈嘉琳，王涌，等.基于形態(tài)特征判別超聲圖像中乳腺腫瘤的良惡性[J].光學(xué)精密工程，2006,14(2):333-340.

[3] 黃劍華，于瀛星，張英濤，等.基于一致性直方圖的超聲乳腺圖片分割方法[J].哈爾濱工業(yè)大學(xué)學(xué)報(bào)，2008,40(7):1103-1106.

[4] Wang W M,Zhu L,Qin J,et al.Multi-scale geodesic active contours for ultrasound image segmentation using speckle reducing anisotropic diffusion [J].Optics and Lasers in Engineering,2014,54:105-116.

[5] Liu B,Chen H D,Huang J H,el al.Probability density difference based active contour for ultrasound image segmentation [J].Pattern Recognition,2010,43(6):2028-2042.

[6] Chan T F,Vese L A.Active contours without edges [J].IEEE Transactions on Image Processing,2001,10(2):266-277.

[7] 張虎重，康志偉.融合模糊聚類的變分水平集圖像分割模型[J].電子測量與儀器學(xué)報(bào)，2011,25(4):325-330.

[8] 高粱,劉曉云，陳武凡，等.基于相位和GGVF的水平集乳腺超聲圖像分割[J].儀器儀表學(xué)報(bào)，2012,33(4)：870-877.

[9] Li L,Kruusmaa M.Ultrasound image segmentation with localized level set based curve evolution and Rayleigh distribution [A],Tuduce R.IEEE International Conference on Systems,Signal and Image Processing[C],Bucharest:IEEE,2013:139-142.

[10]Li C M,Xu C Y,Gui C F,et al.Level set evolution without re-initialization:a new variational formulation [A],Schmid C.IEEE Computer Social on Computer Vision and Pattern Recognition[C],San Diego:IEEE,2005:430-436.

[11]Li C M,Xu C Y,Konwar K M,et al.Fast distance preserving level set evolution for medical image segmentation [A].Processing of the 9th International Conference on Control,Automation,Robotics and Vision[C],Singapore:IEEE,2006:1-7.

[12]Li C M,Xu C Y,Gui C F,et al.Distance regularized level set evolution and its application to image segmentation [J].IEEE Transactions on Image Processing,2010,19(12):3243-3254.

[13]Zhang K H,Zhang L,Song H H,et al.Re-initialization free level set evolution via reaction diffusion [J].IEEE Transactions on Image Processing,2013,22(1):258-271.

[14]虞紅偉，厲力華，徐偉棟，等.一種基于水平集的多尺度乳腺腫塊分割方法[J].儀器儀表學(xué)報(bào)，2010,31(6):1418-1423.

[15]Turing A M.The chemical basis of morphogenesis,philosophical transactions of the royal society of London [J].Series B,Biological Sciences,1952,237(641):37-72.

[16]Zhang Q,Han H,Ji C H,et al.Gabor-based anisotropic diffusion for speckle noise reduction in medical ultrasonography [J].Journal of the Optical Social of America A,2014,31(6):1273-1283.

[17]Turk G.Generating textures on arbitrary surfaces using reaction-diffusion[J].Computer Graphics,1991,25(4):289-298.

[18]Merriman B,Bence,Osher S.Motion of multiple junctions:a level set approach [J].Journal of Computational Physics,1994,112(2):334-363.

增強(qiáng)安全意識，提升防爆技能——自儀院防爆技術(shù)培訓(xùn)火熱招生

隨著國民經(jīng)濟(jì)的迅速發(fā)展，石油、化工、鋼鐵、煤炭、制藥等企業(yè)的生產(chǎn)規(guī)模不斷擴(kuò)大，防止爆炸性事故的發(fā)生已得到社會普遍關(guān)注，昆山工廠粉塵爆炸等多起事故以生命和鮮血一次次地敲響了安全生產(chǎn)的警鐘。長期實(shí)踐表明，防爆生產(chǎn)裝備相關(guān)的管理、設(shè)計(jì)、安裝、使用、維護(hù)人員全面掌握防爆電氣理論和操作技能對于保障企業(yè)安全生產(chǎn)的必要性已毋容置疑。

為增強(qiáng)企業(yè)安全意識，提升防爆作業(yè)人員技術(shù)水平，上海工業(yè)自動化儀表研究院繼續(xù)教育培訓(xùn)中心將于2015年11月10日在上海舉辦防爆技術(shù)培訓(xùn)，現(xiàn)已火熱招生，誠邀您的參與。

1.特種作業(yè)人員電工(防爆電氣)安全技術(shù)培訓(xùn)班

(頒發(fā)國家安監(jiān)總局監(jiān)制《特種作業(yè)(電工：防爆電氣)操作證》，全國通行)

2.人保部“專業(yè)技術(shù)人才知識更新工程”防爆技術(shù)研修班

(頒發(fā)國家人保部“專業(yè)技術(shù)人才知識更新工程”繼續(xù)教育證書，全國通行)

招生詳情可見www.sipai.com，或詢培訓(xùn)中心姜老師、陳老師(021-64368180-324，training@sipai.cn/tccedu@sipai.cn)。

Segmentation of Breast Tumor Ultrasound Images Using Reaction Diffusion Level Set

For segmentation of breast tumor ultrasonic images, the improved reaction diffusion (RD) level set segmentation algorithm is proposed. Firstly, the Gabor anisotropic diffusion model is used for filtering, to construct edge stopping function; then the function is integrated into RD level set evolution equation, for controlling the evolution of curves to obtain the boundary of breast tumor. Segmentations of 111 breast ultrasonic images from 77 patients are conducted with the proposed method and the traditional method, the segmentation accuracies are 98.5% and 98.0% respectively; the true positive rate are 88.2% and 82.7%; and the mean square errors are 3.6 and 4.6 pixels. The results show that the proposed improved algorithm provides more accurate segmentation for breast tumor.

Edge-stopping function Reaction diffusion(RD) Level set evolution Image segmentation Breast ultrasound image

國家自然科學(xué)基金青年資助項(xiàng)目(編號：61401267);

上海市自然科學(xué)基金資助項(xiàng)目(編號：12ZR1444100);

上海市教委人才計(jì)劃“晨光計(jì)劃”資助項(xiàng)目(編號:11CG45)；

陳帥(1987- )，男，現(xiàn)為上海大學(xué)生物醫(yī)學(xué)工程專業(yè)在讀碩士研究生；主要從事醫(yī)學(xué)圖像處理方面的研究。

TP391

A

10.16086/j.cnki.issn1000-0380.201509008

上海市教委科研創(chuàng)新項(xiàng)目(編號：12YZ026)。

修改稿收到日期：2015-03-24。