基于圖像分割的視網膜血管圖像配準研究

劉志東, 林江莉, 陳 科*

(1. 四川大學 材料科學與工程學院, 四川 成都 610065; 2. 四川城市職業學院 汽車與信息工程學院, 四川 成都 610010)

基于圖像分割的視網膜血管圖像配準研究

劉志東1,2, 林江莉1, 陳 科1*

(1. 四川大學 材料科學與工程學院, 四川 成都 610065; 2. 四川城市職業學院 汽車與信息工程學院, 四川 成都 610010)

利用不同波長的視網膜圖像可以對視網膜血管血氧飽和度進行計算,但需進行配準處理.提出一種基于視網膜圖像血管分割的互信息圖像配準方法.為了充分利用血管的輪廓信息和灰度信息,提高配準精度,首先對配準圖像進行圖像分割,提取視網膜圖像中的血管輪廓信息；然后對分割后圖像中的血管進行相似度計算,并采用Powell優化算法中的黃金分割法一維搜索算法來提升運算速度；最后根據計算的相似度值來完成不同波長圖像的配準.研究中算法配準獲得變換參數(角度、X方向、Y方向)的誤差的均值分別為2.00%、2.53%和2.52%,誤差的方差分別為0.57、2.09和0.34,均優于直接互信息配準方法.實驗表明：算法可以自動、有效地對不同波長的視網膜血管圖像進行配準,并具有良好的可重復性和穩定性.

視網膜圖像; 圖像分割; 配準; 互信息

人眼的很多疾病,如青光眼、視網膜病變、血管阻塞等,都會引起視網膜血管血氧飽和度值的變化,因此對視網膜血管中血氧飽和度值的無損測量是目前研究的熱點問題.A. Harris等[1]根據組織對光譜的吸收特性的研究,提出了基于光譜法的無損視網膜血氧測量方法,該方法根據含氧血紅蛋白和脫氧血紅蛋白對2個波長的光吸收差異,同時獲取2個不同波長的視網膜圖像,但由于2幅圖像是非共光路引入,因此圖像位置會出現偏差,在計算血氧值之前,必須對2幅圖像進行配準.

到目前為止,圖像配準的方法很多,按照配準過程基本上可以分為兩類：基于圖像灰度信息和基于圖像特征[2-3].基于圖像特征的配準方法有一定的不足之處,在特征提取過程中很多情況需要人工干預,同時圖像特征的提取本身就是一個比較困難的問題,而且特征提取的精確度也無法保證,會直接影響配準的效果[4].基于圖像灰度信息的方法主要利用的是圖像中像素的灰度信息,根據概率密度計算公式對2幅圖像之間的相似度進行計算,通過優化搜索算法獲得配準圖像間的變換參數.

基于灰度信息配準的方法,目前主要有聯合熵法、相關熵法、最大互信息法、條件熵法等.與基于圖像特征配準方法相比,基于灰度信息的配準方法是一種自動的方法,無需人工干預,依賴于圖像本身信息進行配準,同時可以適用于不同模態圖像的配準[5-7],避免了主觀因素的影響.如洪明堅等[8]結合圖像二維信息與互信息對圖像進行配準；翟利志等[9]利用鄰域信息與互信息相結合處理紅外與可見光圖像配準的問題.越來越多的研究將灰度信息與圖像目標結構信息相結合,充分利用圖像中存在的紋理、形狀、結構等信息,提高配準的準確性是目前研究的熱點.

本文針對2幅不同波長視網膜圖像進行配準,充分利用圖像中血管輪廓信息,采用基于圖像分割的互信息配準方法進行配準研究.為了充分利用血管的輪廓信息和減少圖像配準中無關信息對配準結果的影響,首先對配準圖像進行圖像分割,提取視網膜圖像中的血管信息；對分割后圖像中的血管進行相似度計算,并采用Powell優化算法中的黃金分割法一維搜索算法來提升運算速度；根據計算的相似度值來完成不同波長圖像的配準.本文的算法流程如圖1所示.

1 圖像分割

1.1 圖像預處理 根據視網膜血管不同組織結構對光譜的吸收特性不同,血管內部和血管壁在圖像中呈現不同的灰度級別,尤其在一些大尺寸的血管中,存在中央反射的問題,為血管分割帶來了一定的困難.本文采用數學形態學中的開運算來對血管進行填充,將血管內部和血管壁的灰度值調整到一個級別上去,可以有效去除中央反射問題.同時,開運算還會對圖像起到濾波作用,將圖像中存在的一些椒鹽噪聲濾除.視網膜血管在圖像中的像素值低于周圍眼部其他組織,為了在血管分割過程中提高視網膜血管的可見性,將待處理圖像進行圖像反轉[10].在灰度圖像中像素灰度值為0～255,而所謂反轉就是指將原灰度值為g的像素值轉變為255-g,即原來像素值為0,反轉后就變為255,原來為255,反轉后就變為0.通過反轉,可以將視網膜血管的亮度進行很大程度的提升,為后續處理奠定基礎.

1.2 血管增強 視網膜圖像中的視盤和其他相關組織,在血管分割過程中都會對分割結果產生一定的影響.,為了避免這些影響,本文采用top-hat變換對血管對比度進行增強[11].

top-hat變換定義為

(1)

其中,IA為反轉后圖像,O(IA,B)為圖像IA數學形態學開運算后結果,運算的結構元素為B,Itop-hat為top-hat變換后圖像.top-hat變換的實現原理為反轉后圖像的血管的灰度值要大于視盤及周邊組織的灰度值.通過選擇一定尺度的結構元素B對圖像進行開運算,圖像中的血管信息將被周邊組織的信息所覆蓋,得到血管被模糊掉的圖像,保留了非血管圖像信息,再在原圖像的基礎上減去非血管圖像信息,使得圖像血管信息與周邊其組織信息的對比度進一步提高,從而增強了血管信息.在變換中,結構元素B的尺寸選擇就很重要,要保證結構元素B在圖像中做開運算時,能夠覆蓋直徑最大的血管,這樣才能將圖像中的血管信息全部去除,提高變換效果.通過測量本文研究中采用的圖像血管最大直徑一般在14像素左右,根據結構元素直徑要大于血管直徑的原則,本文選擇的結構元素直徑為16像素.血管增強效果如圖2(b)所示.

1.3 多窗口檢測 單窗口檢測首先需要定義一個寬和高都為H像素的檢測窗口,檢測窗口在被檢測圖像上滑動并遍歷整幅圖像,在遍歷過程中通過獲取最有意義檢測線來獲取圖像血管的特征.單窗口檢測的過程如下：

1) 定義檢測窗口,寬和高同為H像素;

(2)

其中,P(i,j)為檢測窗口覆蓋下圖像像素灰度值,NH為檢測窗口覆蓋的像素個數總數;

(3)

其中,P(i,j)為檢測窗口檢測線上的像素灰度值;

(4)

為了使血管特征函數更好地反映出血管特征,檢測窗口的大小選擇就比較重要,檢測窗口的中心處于血管中心線上時,窗口內的血管像素和背景像素數據近似,因此窗口H一般為血管寬度的2倍,而檢測線長度一般與血管寬度相當.但是該方法某些情況下效果不好,如血管緊鄰時,容易將二者合并為一條血管.

針對單窗口檢測對相鄰血管處理效果不好等存在的問題,本文采用多窗口檢測的方法,通過改變窗口的寬度,控制檢測線的長度,構建多窗口檢測器.特定檢測窗口長度h下的血管特征響應函數[12]為

(5)

其中,H為最大檢測窗口的寬度,h的范圍為1≤h≤H,表示檢測窗口的寬度.不同寬度的檢測窗口在血管特征提取過程中具有不同的作用,長寬度的檢測窗口可以有效的對大尺寸的血管進行特征提取,而短寬度的檢測窗口則可以有效的處理緊鄰的血管,避免了單窗口檢測的相鄰血管合并的問題.采用線性組合的方式,對不同寬度窗口檢測的優點進行合并,獲得多窗口檢測的最終響應函數

(6)

其中,Nh為不同寬度的檢測窗口的數量.

1.4 后處理 圖2展示了圖像血管分割過程.

通過多窗口檢測操作,將圖像血管特征進行提取,這里所謂提取并不是已經對圖像中的血管進行了分割,而是通過檢測操作,獲取了圖像每一個像素點屬于血管或背景的可能性,即檢測后的圖像像素值集中在2個區域,即圖像的直方圖集中呈現為2個分明的峰值,有了這個基礎,就可以很容易利用閾值分割方法對圖像進行分割處理.但是分割結果中存在很多的小的分割區域,這些小的區域在后續配準中存在一定的影響,必須將其去除,如圖2(c)所示.本文根據血管的連通性,將連通區域的面積小于100的連通區域去除,如圖2(d)所示.

2 圖像配準

在血管分割的基礎上,獲取分割后的血管圖像,利用血管輪廓信息與灰度信息,采用互信息配準方法對不同波長下的視網膜圖像進行配準[13].互信息用來描述2個系統間的信息相關性,用熵進行定量的數學描述.

圖像A的熵定義為

(7)

圖像B的熵定義為

(8)

圖像A、B的聯合熵定義為

3）可以通過搜索引擎，在百度、搜狐上搜索一些熟知的英文新聞網站和英文學習網站，了解最新國內外大事和與四、六級考試相關的資訊。

(9)

其中,PA(i)、PB(j)分別表示圖像A、B的灰度概率密度分布,PAB(i,j)表示圖像A、B的聯合灰度概率分布.

圖像A、B的最大互信息為

(10)

根據統計學分析,如果圖像A和圖像B是完全獨立的,則

(11)

同時I(A,B)=0；如果圖像A和圖像B完全依賴,則PAB(i,j)=PA(i)=PB(j),即H(A)=H(B)=H(A,B),則I(A,B)最大.對于圖像A和圖像B,當二者的空間位置一樣時,2幅圖像間的聯合概率密度分布最集中.從2幅圖像之間關聯性角度分析,則指2幅圖像間的相關性最大.而從熵的角度來看,即最終計算出來的2幅圖像之間的聯合熵的值最小,表明2幅圖像之間的互信息最大,這個時候表明2幅圖像已經配準了.

本研究中互信息配準流程如圖3所示.

圖像配準目的就是尋找合適的變換參數,所以其過程就是尋找最優變換參數的過程.所謂最優變換參數就是使得2幅圖像之間的互信息最大.但是,互信息最大值的計算和搜索是比較慢的,目前應用于互信息計算的優化算法很多,Powell優化算法就是其中一種[14].本文采用Powell優化算法中的一維黃金分割搜索算法.一維黃金分割搜索算法的思路就是,先利用黃金分割法確定一個較小的包含極小點的不確定區間,然后利用拋物線法獲得一個極小點,若此極小點落在此不確定區間,則利用該極小點繼續進行二次插值；否則放棄該點,改用黃金分割法繼續搜索,經過多次迭代后便可求得互信息的最大值[15-16].

3 實驗結果

表1中共列舉了5組實驗數據.將原圖與旋轉后的圖像直接采用直接互信息配準方法和本文中的配準方法進行配準,并計算每一組的配準誤差.將配準誤差進行均值和方差計算,其中采用直接互信息配準方法的角度、X方向和Y方向參數誤差的均值分別為9.68%、9.14%和10.40%,誤差的方差分別為2.56、4.30和15.03.采用本文算法的誤差均值分別為2.00%、2.53%和2.52%,誤差的方差分別為0.57、2.09和0.34.通過對比可以看出,本文算法無論準確性和算法的魯棒性都優于直接互信息配準方法.

表 1 方式-圖像配準參數表

數據實際變換參數角度/(°)XY直接互信息參數角度/(°)XY本文算法參數角度/(°)XY直接方法誤差/%角度/(°)XY本文方法誤差/%角度/(°)XY15．0010．0013．004．5811．0314．235．079．8513．288．4010．309．461．401．502．1029．00-12．00-8．009．72-13．43-9．899．11-11．69-8．158．0011．9111．131．202．581．8836．00-14．007．005．38-13．018．165．83-13．786．8910．307．0716．572．831．571．574-7．0011．0013．00-6．329．9811．87-7．1311．5513．369．719．278．691．865．002．775-8．00-12．00-16．00-8．96-11．14-16．99-8．22-11．76-16．4712．007．176．192．752．002．94

第二種方式采用設備采集的570nm圖像和600nm圖像作為配準數據,但是該組數據的旋轉和平移的變換參數無法未知,無法采用方式一中的方法進行配準結果的定量判斷.基于此,本文采用直觀觀察的方式進行效果驗證,即通過將570nm圖像分割得到的血管輪廓按照配準獲得的參數加載到600nm圖像上,如圖5中的輪廓線,該輪廓線就是570nm圖像中分割獲得的血管輪廓,通過加載到600nm圖像中,觀察600nm圖像中的血管輪廓與輪廓線的重合程度,來判斷配準結果的準確性.

如圖5顯示了直接互信息配準算法與本文算法在配準后570 nm圖像中的血管輪廓線疊加到600 nm圖像下的疊加效果.

通過圖5中2幅血管輪廓的效果疊加圖來看,圖5(a)和(b)中分別畫出了3組不同的區域,通過同種顏色框區域內容的對比,可以看出圖5(b)中各個方框中的輪廓線和血管輪廓線的重合度要優于圖5(a)中的各區域.通過直觀觀察的方式,可以證明本文配準方法優于直接互信息的配準方法.

本文方法充分考慮到不同波長視網膜圖像背景的復雜性和灰度分布的不均勻性,充分利用圖像中血管的輪廓信息和灰度信息,采用多窗口檢測的方法對視網膜圖像血管進行分割處理,再利用Powell優化算法的黃金分割搜索方法完成圖像的配準.通過實驗驗證,本文方法相對于直接互信息的配準方法,可以很好地提高配準精確度,證明本文方法在不同波長視網膜圖像的配準方法是可行的.但是本文還存在一定的問題,如配準過程中,分割的耗時相對要長一點,因此,在后續的研究中,分割方法的運行效率方面是研究的重點之一.

[1] HARRIS A, DINN R B, KAGEMANN L. A review of methods for human retinal oximetry[J]. Ophthalmic Surgery Lasers & Imaging,2003,34(2):152-164.

[2] 李新勝,劉宇. 基于邊緣特征的虛實圖像精確配準[J]. 四川大學學報(工程科學版),2014,46(5):96-103.

[3] 何敬,李永樹,李歆,等. 基于點特征和邊緣特征的無人機影像配準方法[J]. 西南交通大學學報(自然科學版),2012,46(6):955-961.

[4] 王紅玉,馮筠,崔磊,等. 應用顯著紋理特征的醫學圖像配準[J]. 光學精密工程,2015,23(9):2656-2665.

[5] PLUIM J, MAINTZ J, VIERGEVER M A. Image registration by maximization mutual information and gradient information[J]. IEEE Trans on Medical Imaging,2001,19(8):809-819.

[6] DEEPA M, SARAVANAN T. Automatic image registration using 2D-discrete wavelet transform[J]. Indian J Science & Technology,2016,9(5):1-3.

[7] ARAKI T, IKEDA N, DEY N, et al. A comparative approach of four different image registration techniques for quantitative assessment of coronary artery calcium lesions using intravascular ultrasound[J]. Computer Methods & Programs in Biomedicine,2015,118(2):158-172.

[8] 洪明堅,呂建斌,楊丹,等. 一種新的基于互信息的圖像配準方法[J]. 重慶大學學報(自然科學版),2009,32(6):697-700.

[9] 翟利志,王敬東,李鵬. 基于鄰域信息的紅外與可見光圖像互信息配準[J]. 計算機技術與發展,2008,18(10):151-154.

[10] 華梓錚,華澤璽. 結合NSCT高低頻特征的圖像邊緣檢測算法[J]. 四川師范大學學報(自然科學版),2014,37(4):612-617.

[11] 高向軍. 基于多尺度線性檢測的視網膜血管分割[J]. 科學技術與工程,2013,13(23):6820-6824.

[12] NGNYEN U T V, BHUIYAN A, PARK L A F. An effective retinal blood vessel segmentation method using multi-scale line detection[J]. Pattern Recognition,2013,13(3):703-715.

[13] YAMADA M, SIGAL L, RAPTIS M, et al. Cross-domain matching with squared-loss mutual information[J]. IEEE Transactions on Pattern Analysis & Machine Intelligence,2015,37(9):1764-76.

[14] 潘婷婷,陸麗婷,顧綺芳. QPSO算法和Powell法結合的多分辨率醫學圖像配準[J]. 計算機應用與軟件,2014,31(7):237-240.

[15] 喬玉龍,趙源萌,張存林,等. 基于混合Powell法的太赫茲/可見光雙波段圖像配準[J]. 光學學報,2015,35(2):151-157.

[16] SHAMES R , SADEGHI P, KENNEDY R. Parallel computation of mutual information on the GPU with application to real-time registration of 3D medical images[J]. Computer Methods & Programs in Biomedicine,2010,99(2):133-46.

(編輯 李德華)

Registration of Retinal Vessel Blood Image Based on Image Segmentation

LIU Zhidong1,2, LIN Jiangli1, CHEN Ke1

( 1.CollegeofMaterialsScienceandEngineering,SichuanUniversity,Chengdu610065,Sichuan; 2.DepartmentofAutomobileandInformationEngineering,UrbanVocationalCollegeofSichuan,Chengdu610010,Sichuan)

In order to calculate the blood oxygen saturation of retinal images, different wavelength images should be registered. This paper presents an image registration method based on the segmentation of blood vessel image and mutual information. In the study, in order to reduce the impact of the information on the registration result, the registration of image segmentation, extraction of retinal vessels information in the image; calculating the vascular similarity in the segmented image, and using the Powell optimization algorithm 0.618 one-dimensional search algorithm to improve the speed of operation; the different wavelength of the image registration based on the calculated similarity value. In the study, the error average of parameters (angle,Xdirection,Ydirection) calculated from the registration algorithm is 2.00%, 2.53% and 2.52%, and the variance of the error is 0.57, 2.09 and 0.34, were better than the direct mutual information registration method. Experiments show that the algorithm can automatically and effectively register retinal images with different wavelengths, and has good repeatability and stability.

retinal image; image segmentation; registration; mutual information

2017-02-21

國家自然科學基金(81301286)和四川省科技支撐項目(2014GZ0005)

TP391.4

A

1001-8395(2017)04-0554-07

10.3969/j.issn.1001-8395.2017.04.020

*通信作者簡介：陳 科(1982—),男,博士,主要從事醫學圖像處理的研究,E-mail:chenke@scu.edu.cn