【作 者】吳曉唯，周傳清＊，任秋實(shí)

上海交通大學(xué)生物醫(yī)學(xué)工程系，上海，200240

腎臟是人體的重要器官，它負(fù)責(zé)生成尿液，排除代謝產(chǎn)物；維持體液平衡及體內(nèi)酸堿平衡；還具有內(nèi)分泌功能[1]?，F(xiàn)在越來(lái)越多的人受到腎臟疾病的困擾，而臨床醫(yī)生需要各類腎臟組織參數(shù)來(lái)評(píng)估疾病的嚴(yán)重程度。除了作為參數(shù)之一的腎臟體積之外，其內(nèi)部皮質(zhì)與髓質(zhì)等組織的區(qū)分也顯得非常重要[2]，所以需要合適的圖像分割算法來(lái)提取這些目標(biāo)組織。

除了被認(rèn)為是“黃金標(biāo)準(zhǔn)”但需要更多人力與時(shí)間投入的“手動(dòng)勾勒”方法之外[3]，Kass M.提出的“動(dòng)態(tài)輪廓”[4]作為一種基于模型的方法被廣泛地應(yīng)用。它是一個(gè)由內(nèi)部連續(xù)性、光滑性作用力和外部的基于圖像特征的作用力（或者叫限制性能量成分）來(lái)引導(dǎo)的最小化能量的樣條。但是，傳統(tǒng)動(dòng)態(tài)輪廓模型的限制之一，就是初始輪廓的位置必須盡量接近最終輪廓的位置，否則輪廓在變化過(guò)程中可能會(huì)被某些離散的噪聲點(diǎn)給吸引住。針對(duì)這一點(diǎn)，Lapeer R.J.等人[5]在運(yùn)用動(dòng)態(tài)輪廓模型提取腹部器官的過(guò)程中，采用基于標(biāo)記的分水嶺算法，得到接近最終目標(biāo)邊緣位置的初始輪廓。而Cohen L.D.等人[6]在動(dòng)態(tài)輪廓模型中添加了一個(gè)“氣球力”，使初始輪廓可以手動(dòng)設(shè)置在遠(yuǎn)離最終輪廓的位置。

本文基于腎臟的核磁共振圖像，提出了一種結(jié)合閾值分割與動(dòng)態(tài)輪廓模型的方法，用來(lái)提取腎臟內(nèi)部的皮質(zhì)與髓質(zhì)區(qū)域。首先，閾值分割得到一個(gè)初始的掩碼圖像；其次，模仿氣球力的簡(jiǎn)化動(dòng)態(tài)輪廓模型得到所需組織的初始輪廓；最后，初始輪廓在傳統(tǒng)動(dòng)態(tài)輪廓模型的作用下被真正的由梯度圖像定義的組織邊緣給吸引住。

1 算法與實(shí)現(xiàn)

圖1 算法流程圖Fig.1 The fl ow of the analysis process

圖1顯示了算法的整個(gè)過(guò)程。首先，通過(guò)預(yù)處理從原始腎臟圖像中得到對(duì)比度增強(qiáng)的圖像，繼而就能通過(guò)梯度運(yùn)算得到梯度圖像，同時(shí)通過(guò)閾值分割得到掩碼圖像。接著，手動(dòng)設(shè)置的起始輪廓在簡(jiǎn)化的動(dòng)態(tài)輪廓模型的氣球力作用下變形，得到所需的組織初始輪廓，同時(shí)也得到了內(nèi)部作用力。最后，初始輪廓在傳統(tǒng)的動(dòng)態(tài)輪廓模型的內(nèi)外部作用力下，變形得到最終的皮質(zhì)與髓質(zhì)的區(qū)域。

1.1 圖像預(yù)處理與閾值分割

圖像預(yù)處理的目的是為了提高原始圖像不同區(qū)域間的對(duì)比度，改善之后的梯度圖像和閾值分割圖像的質(zhì)量。文中采用基于圖像形態(tài)學(xué)操作[7]的tophat濾波和bottomhat濾波結(jié)果相減的方法。

對(duì)閾值分割[7]來(lái)說(shuō)，它是非模型方法中最簡(jiǎn)易最快的自動(dòng)分割算法。在一個(gè)最小閾值設(shè)定后，圖像中任何灰度值高于這個(gè)閾值的像素都被認(rèn)為是感興趣區(qū)域的一部分。閾值分割在文中的目的是為了得到一個(gè)腎臟二值掩碼圖像，基于它，可以得到所需的皮質(zhì)與髓質(zhì)部分的初始輪廓。

1.2 傳統(tǒng)的離散動(dòng)態(tài)輪廓模型

動(dòng)態(tài)輪廓模型是一個(gè)連續(xù)和光滑的可變形樣條，它受到由內(nèi)，外部作用力組成的能量最小化函數(shù)的作用[4]。在實(shí)際應(yīng)用中，采用公式(1)中的離散形式：

“貪婪算法”[8]被用來(lái)作為動(dòng)態(tài)輪廓上離散點(diǎn)v1,v2,......vN運(yùn)動(dòng)的依據(jù)，見(jiàn)圖2所示。每個(gè)點(diǎn)的運(yùn)動(dòng)被限制在以它為中心的3×3矩陣中，移動(dòng)到的那個(gè)像素點(diǎn)具有公式（1）中的三個(gè)能量和的最小值。

圖2 vi以為中心的3×3鄰域矩陣。點(diǎn)vi將會(huì)移動(dòng)到具有最小能量和的像素位置。Fig.2 The 3×3 neighborhood in the discrete active contour model.Pointwill be moved to the pixel vi with the minimal energy value in this 3×3 matrix.

內(nèi)部作用力參考Lobregt S.等人提出的模型[9]，對(duì)于點(diǎn)來(lái)說(shuō)，

i一個(gè)離散點(diǎn)指向當(dāng)前離散點(diǎn)的差向量。同時(shí)，

定義為。在當(dāng)前點(diǎn)的兩個(gè)內(nèi)部作用力的3×3領(lǐng)域矩陣創(chuàng)建后，需要除以最大值對(duì)他們進(jìn)行歸一化。

圖3 點(diǎn)v的內(nèi)部作用力。曲率定義為。是從前一點(diǎn)v指ii-1向當(dāng)前點(diǎn)vi的差向量Fig.3 Internal forces of the point vi. The curvature s de fi ned as 1. is the difference vector pointing from vi-1 to vi

外部作用力采用傳統(tǒng)的梯度圖像，即EgradiertimgelI是對(duì)比度增強(qiáng)后的圖像。在3×3領(lǐng)域矩陣創(chuàng)建后，同樣需要進(jìn)行歸一化。

在輪廓上所有離散點(diǎn)每變動(dòng)一次位置后 ，需要用“重采樣”的機(jī)制[9]來(lái)控制相鄰點(diǎn)間的距離。如果距離小于用戶定義的最小距離lmin，那么兩個(gè)點(diǎn)將由它們的中間點(diǎn)來(lái)替代；如果距離超過(guò)了用戶定義的最大距離lmax，這兩個(gè)點(diǎn)的正中位置將插入一個(gè)新的點(diǎn)。

1.3 模仿氣球力的簡(jiǎn)化動(dòng)態(tài)輪廓模型

為了從腎臟的二值掩碼圖像中得到所需組織的初始輪廓，文中模仿“氣球力”[6]的作用，對(duì)公式(1)的傳統(tǒng)離散動(dòng)態(tài)輪廓模型進(jìn)行了簡(jiǎn)化。不同于Cohen L.D.[6]提出的在模型中添加一個(gè)氣球力的方式，我們將表示內(nèi)部作用力的“曲率”在作用上轉(zhuǎn)化成使輪廓不斷膨脹變形的“氣球力”，同時(shí)去掉對(duì)“氣球膨脹”沒(méi)有作用的Econtinous。簡(jiǎn)化公式表示如下：

首先，切向量和法向量[9]的定義方式見(jiàn)圖4。結(jié)合圖3，是向量和的和，而是順時(shí)針轉(zhuǎn)90度的結(jié)果。

圖4 切向量和法向量r的定義。是與d的和的單位向量，它順ii時(shí)針轉(zhuǎn)90度就成為ri。Fig.4 The tangential component whichis the unit vector of the sum of and he radial componentiss therotated 90 degrees clockwise.

外部作用力在公式（4）中不是原始圖像的梯度圖像，而是二值化掩碼圖像的梯度圖像。

2 結(jié)果

2.1 二值掩碼圖像

在matlab的編程環(huán)境中，首先對(duì)原始MRI腎臟圖像用imtophat-imbothat的形態(tài)學(xué)方式來(lái)提高圖像的對(duì)比度(結(jié)構(gòu)元素采用半徑為40的“disk”)。之后用graythresh自動(dòng)檢測(cè)閾值的方式進(jìn)行圖像的閾值分割。最后對(duì)分割結(jié)果，先用參數(shù)為3的“square”結(jié)構(gòu)元素通過(guò)imreconstruct去除圖像背景上的一些噪聲點(diǎn)，接著用參數(shù)為4的“square”結(jié)構(gòu)元素對(duì)圖像進(jìn)行imclose的閉操作得到最后的腎臟掩碼圖像。圖5分別顯示原始圖像、對(duì)比度增強(qiáng)后的圖像和閾值分割處理后的掩碼圖像。

圖5 (A) 原始MRI腎臟圖像； (B) 對(duì)比度增強(qiáng)圖像；(C)閾值分割掩碼圖像Fig.5 (A)Original MRI kidney image (B)Contrast enhanced image(C)Thresholding mask

2.2 皮質(zhì)與髓質(zhì)區(qū)域的初始輪廓

根據(jù)2.3節(jié)中對(duì)原始動(dòng)態(tài)輪廓模型的簡(jiǎn)化，公式(4)中參數(shù)設(shè)置為：B2=0.5, R2=5。在髓質(zhì)和皮質(zhì)區(qū)域設(shè)置的起始三角型輪廓，會(huì)在氣球力的作用下迅速變形，?？康浇M織的邊緣處(見(jiàn)圖6)，作為下一步驟的初始輪廓。

圖6 手動(dòng)設(shè)置的起始輪廓與輪廓“變形過(guò)程”。(A) 5個(gè)三角型分別是髓質(zhì)與皮質(zhì)的手動(dòng)設(shè)置的起始輪廓；(B)、(C) 是髓質(zhì)中兩個(gè)部分的輪廓變化過(guò)程；(D) 是皮質(zhì)部分的輪廓變化過(guò)程Fig.6 Manually set starting contours and the deforming process. (A)Five triangles are the starting contours of the cortex and medulla (B)and(C) show the deforming process of the contour in two regions of the medullas (D)displays one deforming phase of the contour in the cortex region.

2.3 最終分割的結(jié)果

在得到皮質(zhì)與髓質(zhì)的初始輪廓之后，根據(jù)公式(1)中傳統(tǒng)的離散動(dòng)態(tài)輪廓模型(設(shè)定)，就能獲取目標(biāo)組織的光滑的邊緣輪廓，位置是由預(yù)處理之后的梯度圖像來(lái)決定的，見(jiàn)圖7所示。

圖7 最終輪廓 (A)髓質(zhì)區(qū)域的分割結(jié)果；(B)皮質(zhì)區(qū)域的分割結(jié)果Fig.7 Final contours. (A) Segmentation of medullary regions(B) Extraction of cortical regions.

3 結(jié)論

本文實(shí)現(xiàn)了一種結(jié)合非模型與基于模型的分割算法。非模型的閾值分割先得到的是具有目標(biāo)組織形態(tài)的掩碼圖像。通過(guò)仿照氣球膨脹作用而簡(jiǎn)化的動(dòng)態(tài)輪廓模型，可以手動(dòng)設(shè)置最簡(jiǎn)單的起始輪廓（三角型），并在很短時(shí)間內(nèi)得到所需組織的初始輪廓。它很接近最終輪廓位置，因而在傳統(tǒng)的動(dòng)態(tài)輪廓模型作用下就能準(zhǔn)確地被組織邊緣吸引住，從而完成相關(guān)區(qū)域的分割與提取。未來(lái)的工作中，除了會(huì)嘗試更多的非模型與基于模型算法的組合方式，如分水嶺和動(dòng)態(tài)輪廓之外，也會(huì)對(duì)更多的醫(yī)學(xué)圖像數(shù)據(jù)進(jìn)行測(cè)試。

[1] 百度百科-腎臟 [OL]. http://baike.baidu.com/view/66040.htm. 最近更新2010. 12.25.

[2] de Priester JA, den Boer JA, Giele EL, et al. MR renography:an algorithm for calculation and correction of cortical volume averaging in medullary renographs [J]. J Magn Reson Imaging,2000, 12(3): 453-459.

[3] Chow TW, Takeshita S, Honjo K, et al. Comparison of manual and semi-automatic delineation of regions of interest for radioligand PET imaging analysis [J]. BMC Nuclear Medicine, 2007, 7:2.

[4] Kass M, WitKin A, Terzopoulos D. Snake: Active contour model [J].International Journal of Computer Vision, 1988: 321-331.

[5] Lapeer R.J., Tan A.C., Aldridge R. Active watersheds: combining 3D watershed segmentation and active contours to extract abdominal organs from MR images [A]. MICCAI [C] 2002, 2488:596-603.

[6] Cohen L D On active contour models and balloons [J]. Cvgip:Image Understanding, 1991, 53(2): 211-218.

[7] Gonzalez R C, Woods R E 著, 阮秋綺, 等, 譯. 數(shù)字圖像處理(第二版)[M]. 北京: 電子工業(yè)出版社. 2007.

[8] Williams DJ, Shah M. A fast algorithm for active contours and curvature estimation [J]. Computer Vision, Graphics, and Image Processing: Image Understanding. 1992, 55(1):14-26.

[9] Lobregt S, Viergever M. A discrete dynamic contour model [J].IEEE Transactions on medical imaging. 1995, 14(1): 12-24.