劉曉軍 董恩清 呂成林 李曉陽(yáng) 崔 波 鄭 強(qiáng)
(山東大學(xué)(威海)機(jī)電與信息工程學(xué)院,威海 264209)
引言
對(duì)腦部基底節(jié)區(qū)及黑質(zhì)核團(tuán)中的幾個(gè)微小區(qū)域的磁共振成像(MRI)圖像進(jìn)行精確而穩(wěn)定的分割是準(zhǔn)確測(cè)定發(fā)生帕金森?。≒arkinson’s disease,PD)病變MRI影像特征,評(píng)價(jià)PD患者患病程度最為關(guān)鍵性的工作。
目前,已經(jīng)提出了大量的計(jì)算機(jī)輔助分割方法[1],其中包括可變模型法[2]、動(dòng)態(tài)輪廓線演化法[3]、基于圖譜的彈性配準(zhǔn)法[4]、基于手工標(biāo)記訓(xùn)練集的Bayesian方法[5]和信息融合[6]等。盡管在這些文獻(xiàn)中有很多自動(dòng)和半自動(dòng)的方法,但很多臨床醫(yī)生因?yàn)閾?dān)心這些方法的可靠性、對(duì)圖像數(shù)據(jù)變化的魯棒性、高的計(jì)算代價(jià)及對(duì)初始化敏感性等仍然會(huì)選擇手工分割。目前由李偉等提出的一種以解剖先驗(yàn)知識(shí)為約束基于動(dòng)態(tài)曲面模型和自適應(yīng)區(qū)域增長(zhǎng)的自動(dòng)3D分割方法[7],對(duì)黑質(zhì)結(jié)構(gòu)體進(jìn)行精確三維分割和提取具有一定的效果。目前還未見(jiàn)采用其他算法對(duì)PD發(fā)病的其他區(qū)域進(jìn)行分割研究的報(bào)道。
蟻群算法(ant colony algorithm,ACA)不僅能夠智能搜索、全局優(yōu)化,而且具有魯棒性、正反饋、分布式計(jì)算,且易與其他算法相結(jié)合等特點(diǎn)[8]。它的聚類特性與圖像的識(shí)別過(guò)程有相似之處,已廣泛地用于醫(yī)學(xué)圖像分割中[9]。
王曉年等根據(jù)主動(dòng)輪廓模型的特點(diǎn)構(gòu)建了一類新的蟻群求解算法[10],其實(shí)質(zhì)是把圖像分割問(wèn)題轉(zhuǎn)化成最優(yōu)路徑的搜索問(wèn)題。白楊等針對(duì)基本蟻群算法易出現(xiàn)早熟和停滯現(xiàn)象,提出了一種動(dòng)態(tài)自適應(yīng)蟻群算法[11],通過(guò)自適應(yīng)調(diào)整初始聚類中心和動(dòng)態(tài)更新局部信息素濃度,使其收斂性和穩(wěn)定性有一定的提高。
近些年來(lái),在國(guó)內(nèi)外開(kāi)展MRI圖像分割研究逐漸增多,但大多數(shù)學(xué)者致力于采取某些算法來(lái)對(duì)整個(gè)腦部MRI圖像進(jìn)行完全精細(xì)分割的目的是不現(xiàn)實(shí)的,僅能滿足宏觀分割的要求。況且在實(shí)際醫(yī)學(xué)檢測(cè)診斷中,常常只關(guān)注其中某個(gè)局部區(qū)域,無(wú)需關(guān)注整體的精細(xì)分割。由于目前這些宏觀的整體分割算法并不是借助于局部腦部解剖結(jié)構(gòu)來(lái)開(kāi)展有監(jiān)督的精細(xì)分割,這勢(shì)必導(dǎo)致分割結(jié)果與實(shí)際解剖結(jié)構(gòu)之間的差異增大。所以說(shuō),基于局部微小組織的精細(xì)分割將是以后醫(yī)學(xué)圖像分割研究的方向。
所以,在分析多年來(lái)國(guó)內(nèi)外對(duì)MRI圖像分割研究結(jié)果的基礎(chǔ)上,開(kāi)展針對(duì)某一微小目標(biāo)區(qū)域具有一定輪廓約束的,且輪廓的邊界可以是模糊不連續(xù)的自動(dòng)分割研究,本研究提出基于解剖結(jié)構(gòu)特征,采用蟻群聚類算法進(jìn)行有監(jiān)督精細(xì)分割算法研究。
蟻群算法的基本思想是通過(guò)模仿螞蟻在一區(qū)域內(nèi)更新信息素而找到最短路徑的方法來(lái)求解實(shí)際最優(yōu)化問(wèn)題[12]。對(duì)于數(shù)字圖像的邊緣提取,初始狀態(tài)時(shí)螞蟻隨機(jī)分布在像素點(diǎn)上,在螞蟻行經(jīng)的節(jié)點(diǎn)上釋放信息素,根據(jù)其鄰域像素點(diǎn)的信息素與梯度值,以較大的概率選擇信息素分布多、梯度值高的像素點(diǎn),并增強(qiáng)經(jīng)過(guò)像素點(diǎn)的信息素。由于圖像的邊緣點(diǎn)梯度值高,因此,螞蟻逐漸向邊緣匯聚,從而得到圖像邊緣。
2.1.1 提出的搜索空間方案
這里介紹的搜索空間方案1和方案2分別由文獻(xiàn)[12]采用的和由本研究提出的。
方案1,對(duì)于給定的模板輪廓線,先求其質(zhì)心,然后在以質(zhì)心與輪廓線上某一點(diǎn)的連線方向上,在輪廓線上的點(diǎn)的內(nèi)外分別取等間隔若干個(gè)點(diǎn)作為該輪廓線上點(diǎn)的橫向搜索范圍,這樣對(duì)于輪廓線上的每個(gè)點(diǎn)總是有一個(gè)可以搜索選擇的點(diǎn)集。
方案2,對(duì)于所給模板輪廓上的每個(gè)點(diǎn),求其法線,然后在法線方向取搜索空間。
為了減少運(yùn)算量,上面提出的搜索空間方案都是在初始輪廓的基礎(chǔ)上構(gòu)造的。無(wú)論初始輪廓是否封閉,都可以按照某種規(guī)則分別對(duì)模板輪廓線上總點(diǎn)數(shù)為M的每個(gè)點(diǎn)鄰域取L個(gè)離散點(diǎn),形成由M個(gè)集合構(gòu)成的一維搜索空間Ω,那么整個(gè)搜索空間就整理成M行L列的矩陣。其中一行向量空間為v(i)=[x(i,1),x(i,2),…,x(i,L)],i=(1,2,…,M)。
2.1.2 選取搜索空間方案的討論
對(duì)于這兩類方案,主要是針對(duì)腦部殼核區(qū)域的實(shí)驗(yàn)為例進(jìn)行的。
方案1,由于殼核的形狀比較特殊,它的質(zhì)心離邊界比較近,導(dǎo)致搜索空間的點(diǎn)集中有許多重復(fù)點(diǎn),甚至搜索空間的波動(dòng)明顯不符合殼核的大致形狀,導(dǎo)致分割結(jié)果差異較大。
方案2,既可以逼近殼核的形狀,操作起來(lái)又比較簡(jiǎn)單,效果也比較好,所以采取了方案2。
為了尋找最優(yōu)路徑,針對(duì)第k只螞蟻,我們定義如下目標(biāo)函數(shù)
式中,xk(i,lo)∈v(i)和xT(i)分別位于第k只螞蟻的路徑和初始模板輪廓線上,lo∈(1,2,…,L)隨i變化。于是提出的有監(jiān)督蟻群算法求解過(guò)程就是在所有的螞蟻中去尋找極小值的過(guò)程,圖像分割問(wèn)題就轉(zhuǎn)變?yōu)橄率剿镜膬?yōu)化問(wèn)題
式中,S為最終輪廓線,N為螞蟻的總數(shù)。
信息素場(chǎng)記憶信息素的強(qiáng)度[12]是引導(dǎo)螞蟻進(jìn)行路徑選擇的重要依據(jù),也是螞蟻之間交換信息的途徑,初始時(shí)各路徑的信息素相等,但其會(huì)隨著蟻群算法的迭代而改變,設(shè)第k只螞蟻在任意位置xk(i,l1)∈v(i),該螞蟻從v(i+1)集合中選擇任意元素,共有L種可能。那么集合v(i)中的每個(gè)元素與集合v(i+1)中任意元素的序?qū)Χ伎赡苁俏浵佭x擇的路徑,即每?jī)蓚€(gè)集合之間都有L×L種連接的可能。如果初始輪廓離散成M個(gè)控制點(diǎn)(每個(gè)點(diǎn)必須且只有唯一的屬于某集合v(i)),那么所有可能的連接權(quán)值構(gòu)成信息素矩陣τ(i,xk(i,l1),xk(i+1,l2))(信息素場(chǎng)),其中i=(1,2,…,M),k=(1,2,…,N),l1,l2,l3∈(1,2,…,L)。
為了能夠不斷地積累有用的搜索信息,使螞蟻向著圖像邊緣移動(dòng),在信息素正反饋的作用下,使螞蟻趨向于最優(yōu)路徑,這里仍然采用基于基本蟻群系統(tǒng)的信息素更新策略[13]。
當(dāng)蟻群完成一次迭代之后,相應(yīng)路徑上的信息素濃度必須進(jìn)行更新處理,模仿人類記憶的特點(diǎn),對(duì)舊的信息進(jìn)行削弱。設(shè)0<ρ<1為信息激素殘留系數(shù),則剩余信息素為ρ·τ(i,xk(i,l1),xk(i+1,l2)),這里只增加本次迭代最優(yōu)路徑對(duì)應(yīng)的信息素。設(shè)S*為當(dāng)前所有螞蟻找到的最優(yōu)路徑,則其對(duì)應(yīng)的信息素的增量為Δτ,整個(gè)過(guò)程如式(3)描述
式中,Sk′為第k′∈(1,2,…,N)只螞蟻的最優(yōu)路徑S*,Δτ為信息素增量,這里選Δτ為一常量1。
設(shè)螞蟻在xk(i,l1)∈v(i)位置,選擇xk(i+1,l2)∈v(i+1)的概率為[13]:
式中,τ(i,xk(i,l1),xk(i+1,l2))表示連接xk(i,l1)∈v(i)和xk(i+1,l2)∈v(i+1)的信息素值。式(5)定義的η(xk(i+1,l2))為xk(i+1,l2)點(diǎn)的啟發(fā)式信息,考慮當(dāng)前梯度信息的影響,α,β表示螞蟻在路徑選擇過(guò)程中對(duì)所積累的信息和其啟發(fā)信息的權(quán)重
顯然當(dāng)螞蟻選擇概率值較大的點(diǎn)時(shí),能夠保證算法快速收斂,而當(dāng)選擇小概率其他點(diǎn)時(shí),保證全局收斂。因此由式(4)得到的概率直接影響算法的收斂速度。
筆者所提出的求解有監(jiān)督的蟻群算法模型就是尋找式(2)的極小值。通過(guò)構(gòu)造搜索空間Ω,把尋優(yōu)能量函數(shù)極小值的過(guò)程轉(zhuǎn)換成在搜索空間Ω尋找式(2)定義的最優(yōu)路徑S的過(guò)程。這里同樣采用常規(guī)的蟻群算法的優(yōu)化過(guò)程。如果假設(shè)定義的理想輪廓線的目標(biāo)函數(shù)為最小值,蟻群算法就是在信息素的影響下,第k只螞蟻從第i行向量空間中選擇第l1點(diǎn)xk(i,l1)∈v(i)到第i+1行向量空間中選擇第l2點(diǎn)xk(i+1,l2)∈v(i+1),使最終行進(jìn)路徑的目標(biāo)函數(shù)最小。在一次迭代后,每只螞蟻都得到一條行進(jìn)路徑,然后增加所有螞蟻行進(jìn)路徑中的最優(yōu)路徑所對(duì)應(yīng)的信息素。這樣新的信息素場(chǎng)又會(huì)誘導(dǎo)下次螞蟻迭代的過(guò)程,在整個(gè)蟻群的協(xié)同下,直到目標(biāo)函數(shù)的不再下降為止,最后便得到目標(biāo)函數(shù)為最小的路徑。
文中處理的是DICOM國(guó)際標(biāo)準(zhǔn)格式的磁共振腦部圖像。預(yù)處理主要是對(duì)數(shù)據(jù)進(jìn)行線性歸一化生成了普通格式的圖像數(shù)據(jù);然后采用中值濾波及層疊濾波去噪處理,并進(jìn)行灰度均衡,獲得了質(zhì)量相對(duì)較好的圖像供后續(xù)分割。后期處理主要是對(duì)算法獲得的輪廓線進(jìn)行插值平滑,使其更接近于真實(shí)殼核的形狀;同時(shí)采用模板來(lái)對(duì)得到的結(jié)果輪廓進(jìn)行約束,去除偏離模板比較大的輪廓點(diǎn)。
為了分析和驗(yàn)證提出的有監(jiān)督蟻群分割算法(supervised ant colony segmentation,SACS)的實(shí)際效果,同時(shí)采用了模糊聚類分割算法(fuzzy c-mean segmentation,F(xiàn)CMS)、區(qū)域生長(zhǎng)分割算法(region growth segmentation,RGS)、GVF Snake模型分割算法(gradient vector flow snake segmentation,GVFS)及基本蟻群分割算法(basic ant colony segmentation,BACS)一同來(lái)對(duì)MRI腦部殼核進(jìn)行自動(dòng)分割,從而有利于非??陀^地分析說(shuō)明提出的新分割算法的優(yōu)越性。在實(shí)驗(yàn)中采用了McGill University的MRI腦圖像庫(kù)中的數(shù)據(jù)進(jìn)行算法對(duì)比測(cè)試[14]。在這個(gè)MRI圖像數(shù)據(jù)庫(kù)中,共有20個(gè)正常人腦部MRI掃描數(shù)據(jù)體圖像,且每個(gè)人腦部切片序列圖像能夠顯示殼核區(qū)域的切片層數(shù)為15~25層。從每個(gè)人腦圖像序列中的相同位置選取3張切片,共計(jì)60張切片進(jìn)行測(cè)試。
根據(jù)60張測(cè)試圖像及人腦解剖結(jié)構(gòu)圖譜會(huì)同專家預(yù)先設(shè)計(jì)相應(yīng)位置的殼核輪廓“金標(biāo)準(zhǔn)”模板(putamen templet,PT)表1為對(duì)這60張測(cè)試MRI圖像進(jìn)行分割處理得到的分割面積和分割邊界周長(zhǎng)分別與殼核模板的差異百分比。為表述方便,測(cè)試中的分割區(qū)域的面積和分割邊界的周長(zhǎng)都采用像素?cái)?shù)表示。由于殼核左右大小和形態(tài)的差異,所以對(duì)同一位置切片左右殼核分別預(yù)定模板。對(duì)其中的某一個(gè)分割算法的分割結(jié)果,分別統(tǒng)計(jì)同一切片左右殼核分割的面積(邊界周長(zhǎng))相對(duì)于相應(yīng)的模板的面積(邊界周長(zhǎng))的絕對(duì)差異值與該相應(yīng)模板的面積(邊界周長(zhǎng))的百分比。用同樣的方法得到所有測(cè)試切片左右殼核區(qū)域分割的面積和邊界周長(zhǎng)差異百分比并取平均。這兩個(gè)指標(biāo)分別從分割的平面區(qū)域面積和分割邊界的光滑度來(lái)說(shuō)明各種算法的差異性。從這兩個(gè)指標(biāo)說(shuō)明提出的算法非常接近真實(shí)的“金標(biāo)準(zhǔn)”,誤差都在10%之內(nèi)。
圖1為選定該圖像庫(kù)中的一張含有殼核區(qū)域的整個(gè)腦部切片圖。圖2為上面提到的幾種分割算法對(duì)圖1中所示區(qū)域內(nèi)的殼核進(jìn)行分割效果對(duì)比圖。從這些對(duì)比結(jié)果來(lái)看,RGS、FCMS和BACS與實(shí)際的分割輪廓相差太遠(yuǎn),不能滿足實(shí)際的要求。提出的SACS的結(jié)果與實(shí)際分割結(jié)果最為接近,由于GVFS結(jié)果在分割邊界處過(guò)于平滑,不能夠反映真實(shí)邊界的細(xì)節(jié)部分,也不如提出的BACS的分割效果,但比其他幾個(gè)方法要好。
對(duì)同一數(shù)據(jù)體連續(xù)切片進(jìn)行分析,在上下相鄰層之間的形態(tài)相似和單一層內(nèi)的連續(xù)性,只有提出的SACS的結(jié)果非常接近真實(shí)的分割邊界,而其他分割算法的分割結(jié)果與真實(shí)的分割邊界具有非常大的差異性,完全無(wú)法與提出的分割結(jié)果相比。
表1 分割面積和分割邊界周長(zhǎng)分別與殼核模板的差異率Tab.1 The difference rate of the segmentation area and edge perimeter with PT
圖1 含有殼核區(qū)域的整個(gè)腦部切片圖Fig.1 Whole brain MRI slice with putamen region
圖2 幾種算法對(duì)殼核區(qū)域分割對(duì)比圖Fig.2 Comparison chart of putamen segmentation using a few algorithms
由于采用了搜索輪廓與模板的方差作為目標(biāo)函數(shù),所以提出的有監(jiān)督的蟻群算法模型的求解過(guò)程就等價(jià)于把圖像分割過(guò)程變成求與目標(biāo)函數(shù)極小值的過(guò)程。也可以理解為采用蟻群算法來(lái)優(yōu)化初始模板的目標(biāo)函數(shù),形成基于初始模板的蟻群算法,從而把圖像分割問(wèn)題轉(zhuǎn)變成尋找搜索空間的最優(yōu)路徑問(wèn)題。
基于解剖結(jié)構(gòu)特征對(duì)腦部基底節(jié)區(qū)內(nèi)的殼核這樣小規(guī)模、復(fù)雜解剖結(jié)構(gòu)的區(qū)域進(jìn)行精細(xì)分割是一項(xiàng)具有挑戰(zhàn)性的研究課題,目前在國(guó)內(nèi)外未見(jiàn)報(bào)道。
[1]Caviness VSJ,Lange NT,Makris N,et al.MRI-based brain volumetrics:emergence of a developmental brain science[J].Brain and Development,1999,21(5):289-295.
[2]Tohka J,Wallius E,Hirvonen J,et al.Automatic extraction of caudate and putamen in[11C]Raclopride PET using deformable surface models and normalized cuts[J].IEEE Transactions on Nuclear Science,2006,53(1):200-227.
[3]Yushkevich P,Piven J,Hazlett H,et al.User-guided 3-D active contour segmentation of anatomical structures:Significantly improved efficiency and reliability[J].NeuroImage,2006,31:1116-1128.
[4]Kelemen A,SzekelyG,Gerig G.Elastic model-based segmentation of 3-D neuroradiological data sets[J].IEEE Transactions on Medical Imaging,1999,18(10):828-839.
[5]Fischl B,Salat DH,Busa E,et al.Whole brain segmentation:Automated labeling of neuroanatomical structures in the human brain[J].Neuron,2002,33(3):341-355.
[6]Barra V,Boire JY,Automatic segmentation of subcortical brain structures in MR images using information fusion[J].IEEE Transactions on Medical Imaging,2001,20(7):549-558.
[7]李偉,陳武凡.人腦黒質(zhì)神經(jīng)核團(tuán)的精確三維自動(dòng)分割[J].計(jì)算機(jī)工程與應(yīng)用,2008,44(25):206-209.
[8]Dorigo M,Blum C.Ant colony optimization theory:A survey[J].Theoretical Computer Science,2005,344:243-278.
[9]白楊,孫躍,胡銀萍.蟻群算法在磁共振圖像分割中的應(yīng)用[J].中國(guó)醫(yī)學(xué)影像技術(shù),2007,23(9):1402-1404.
[10]王曉年,馮遠(yuǎn)靜,馮祖仁.一種基于主動(dòng)輪廓模型的蟻群圖像分割算法[J].控制理論與應(yīng)用,2006,23(4):515-522.
[11]白楊,孫躍,王君,等.基于動(dòng)態(tài)自適應(yīng)蟻群算法的MRI圖像分割[J].計(jì)算機(jī)科學(xué),2008,35(2):226-229.
[12]Nezamabadi-Pour H,Saryazdi S,and Rashedi E.Edge detection using ant algorithms[J].Soft Computing,2006,10(7):623-628.
[13]曹會(huì)治,王晨.結(jié)合蟻群算法的Snake模型的醫(yī)學(xué)圖像分割方法[J].北京生物醫(yī)學(xué)工程,2007,6(26):245-248.
[14]The McConnell Brain Imaging Centre(BIC)of the Montreal Neurological Institute,McGill University,BrainWeb:Simulated Brain Database[DB].http://mouldy.bic.mni.mcgill.ca/brainweb/anatomic_normal_20.html,2006-06-12/2010-05-12.