基于融合相位特征的視網(wǎng)膜血管分割算法

李媛媛，蔡軼珩，高旭蓉

(北京工業(yè)大學(xué) 信息學(xué)部，北京 100124)(*通信作者電子郵箱caiyiheng@bjut.com)

0 引言

眼底圖像具有獲取簡便、成像清晰等優(yōu)勢，可以幫助人們盡早地檢測眼部疾病，及時治療。視網(wǎng)膜血管作為眼底圖像的重要組成部分，是唯一可以通過非入侵的方式觀察到的微小血管，而且一些全身性的疾病，如糖尿病、高血壓等，往往會引起視網(wǎng)膜血管發(fā)生病變，所以通過對眼底視網(wǎng)膜血管圖像的分析來進(jìn)行一些疾病的早期檢測和輔助治療是一種非常好的醫(yī)療手段。

鑒于視網(wǎng)膜血管分割的重要性，國內(nèi)外學(xué)者提出了許多有效的視網(wǎng)膜血管分割算法，主要分成非監(jiān)督(例如，血管追蹤[1]、匹配濾波[2]、基于模型[3]等)和監(jiān)督[4-9]兩類。其中，監(jiān)督的方法是將眼科專家手工標(biāo)記的視網(wǎng)膜血管圖像作為標(biāo)準(zhǔn)樣本，并基于給定的血管特征通過訓(xùn)練學(xué)習(xí)優(yōu)化分類器參數(shù)得到最終的分類器，然后將此分類器用于眼底視網(wǎng)膜血管的分割。由于監(jiān)督的方法依據(jù)已有的先驗標(biāo)記信息來設(shè)計，所以基于監(jiān)督方法的血管分割效果普遍優(yōu)于基于非監(jiān)督方法的血管分割效果；但是基于監(jiān)督的分割方法需要根據(jù)提取的圖像特征及標(biāo)準(zhǔn)標(biāo)記信息構(gòu)造大量訓(xùn)練樣本，在訓(xùn)練學(xué)習(xí)的過程中逐步優(yōu)化分類器參數(shù)，所以其訓(xùn)練時間較長，分割效率不高，因此，在基于監(jiān)督的分割方法中，提取優(yōu)質(zhì)的圖像特征用于訓(xùn)練和分割，不僅可以提高分割的準(zhǔn)確率，而且可以顯著減少訓(xùn)練樣本數(shù)和訓(xùn)練時長，提高算法效率。

Soares等[4]提取眼底視網(wǎng)膜圖像的Gabor變換和像素灰度特征，利用貝葉斯分類器實現(xiàn)血管檢測。其中，不同尺度和不同角度的Gabor變換可以對不同寬度不同走向的血管進(jìn)行檢測，而且Gabor變換可以濾除噪聲且對光照變化不敏感；但是，該方法只考慮了圖像Gabor特征和灰度特征，對細(xì)小血管的檢測不夠準(zhǔn)確。吳奎等[5]在Gabor特征的基礎(chǔ)上加入組合線性算子提高了分割的準(zhǔn)確率。吳奎[6]提取了眼底視網(wǎng)膜圖像的組合位移濾波響應(yīng)(Combination Of Shifted FIlter REsponses， COSFIRE)特征，采用貝葉斯模型來進(jìn)行血管分割。其中，對稱和非對稱COSFIRE濾波器的使用可以檢測出連續(xù)血管及血管末端，得到更優(yōu)質(zhì)的血管特征。slani等[7]提取眼底視網(wǎng)膜圖像的Hessian矩陣、高帽變換、Gabor變換、條帶選擇組合位移濾波響應(yīng)(Bar-selective Combination Of Shifted FIlter REsponses, B-COSFIRE)特征以及像素灰度特征用于分類。該方法中作者評估了每個特征在混合特征向量中的重要性，其中B-COSFIRE濾波和Hessian矩陣的預(yù)測重要性指標(biāo)較高。Zhu等[8]提取眼底視網(wǎng)膜圖像的形態(tài)學(xué)、相位一致性、Hessian矩陣、向量散場以及局部特征用于分類。其中，向量散場特征具有較高的特異性，形態(tài)學(xué)特征有較高的敏感性。朱承璋等[9]提取眼底視網(wǎng)圖像的不變矩、灰度共生矩陣、高斯拉普拉斯(Laplace Of Gaussian， LOG)算子結(jié)合高斯二階導(dǎo)、相位一致性及Hessian矩陣特征用于分類。其中，LOG算子結(jié)合高斯二階導(dǎo)特征有較高的敏感性，Hessian特征有較高的特異性，相位一致性特征有較高的準(zhǔn)確性；但由于相位一致性特征主要用于血管邊緣檢測，提取出的相位一致性特征中血管中心部分未能檢測出來，所以其敏感性低。

上述基于監(jiān)督的眼底視網(wǎng)膜血管分割算法大多提取較多的圖像特征，對每個像素構(gòu)造高維的特征向量用于分類器訓(xùn)練與分割，耗時長。此外，現(xiàn)有文獻(xiàn)中大多直接使用眼底視網(wǎng)膜圖像的相位一致性特征。相位一致性特征雖對血管邊緣信息有很好的檢測效果，但其對血管中心部位檢測不足，不適用于基于像素分類的視網(wǎng)膜血管分割算法。通過對常用的眼底視網(wǎng)膜圖像特征的分析和研究，本文基于預(yù)測重要性指標(biāo)較高的B-COSFIRE特征、Hessian矩陣特征以及對方向性檢測敏感的Gabor變換特征，針對相位一致性特征對血管中心部位檢測不足的問題，提出了基于融合相位特征的視網(wǎng)膜圖像血管分割算法，該算法克服相位一致性特征對血管中心檢測的不足，且提取較少的特征，實現(xiàn)高效的視網(wǎng)膜血管的分割。

1 基于融合相位特征的視網(wǎng)膜血管分割算法

本文算法是由預(yù)處理、特征提取、分類器分類以及后處理4個模塊組成，具體算法流程如圖1所示。預(yù)處理后，提取視網(wǎng)膜圖像的Hessian矩陣特征、Gabor特征、B-COSFIRE特征以及相位特征，對圖像中每個像素構(gòu)造一個4D的特征向量，采用支持向量機(jī)(Support Vector Machine, SVM)分類器進(jìn)行訓(xùn)練與分類，實現(xiàn)了高效的視網(wǎng)膜血管分割。

圖1 本文算法流程

1.1 預(yù)處理

本文通過預(yù)處理過程增強眼底視網(wǎng)膜圖像中血管與背景的對比度，便于后續(xù)的血管提取。在原始的彩色視網(wǎng)膜圖像中：如圖2(a)，血管與背景部分對比度不高，且存在圖像噪聲，不適合進(jìn)行后續(xù)處理。而各類研究表明，綠色通道分量圖像的對比度高于其他兩個通道分量，且完整地保存了視網(wǎng)膜圖像中的血管結(jié)構(gòu)，故本文選取綠色通道分量圖像進(jìn)行后續(xù)圖像處理；圖2(b)為提取的綠色通道分量的圖像。為進(jìn)一步增強視網(wǎng)膜圖像的對比度，抑制噪聲，本文采用對圖像增強效果較好的受限制對比度自適應(yīng)直方圖均衡化(Contrast Limited Adaptive Histogram Equalization, CLAHE)[6]作圖像預(yù)處理；圖2(c)為對綠色通道分量進(jìn)行CLAHE的結(jié)果。預(yù)處理后，眼底視網(wǎng)膜圖像中血管與背景的對比度明顯增強，更適合后續(xù)的處理。

圖2 預(yù)處理過程結(jié)果

1.2 特征提取

1.2.1 Hessian矩陣

根據(jù)眼底視網(wǎng)膜圖像中目標(biāo)血管呈線性管狀結(jié)構(gòu)的特點，本文采用對線性結(jié)構(gòu)敏感的Hessian矩陣來提取像特征。Hessian矩陣是一個由多元函數(shù)二階偏導(dǎo)數(shù)構(gòu)成的方陣，它的特征值和特征向量可以很好地描述視網(wǎng)膜血管類似樹杈形狀的拓?fù)浣Y(jié)構(gòu)，并且通過濾波函數(shù)將眼底視網(wǎng)膜圖像中的噪聲去除。Frangi等[10]最先將Hessian矩陣用于眼底視網(wǎng)膜血管增強，并提出線性濾波函數(shù)。

由于視網(wǎng)膜血管的直徑存在變化，單一尺度的Hessian矩陣的血管檢測效果不佳，故將Hessian矩陣的差分運算與高斯核函數(shù)結(jié)合，通過改變高斯函數(shù)的標(biāo)準(zhǔn)差，來獲得不同尺度下的Hessian矩陣。根據(jù)高斯函數(shù)的卷積性質(zhì)，空間尺度導(dǎo)數(shù)Fab可由輸入圖像I與高斯函數(shù)的二階倒數(shù)卷積得到：

(1)

通過分析Hessian矩陣的各特征值與其對應(yīng)的線性目標(biāo)關(guān)系，F(xiàn)rangi等[10]提出在σ尺度下的線性濾波公式：

(2)

1.2.2 Gabor變換

在視網(wǎng)膜血管圖像特征中，Gabor變換具有良好的方向和尺度選擇特性，可以檢測不同方向與不同尺寸的血管，而且Gabor變換對光照變化不敏感，可以克服眼底視網(wǎng)膜圖像光照不均的缺點，因此本文選擇Gabor變換進(jìn)行圖像特征的提取。Gabor變換即窗函數(shù)為高斯函數(shù)的加窗傅里葉變換，有時可歸為小波變換的一種。本文中Gabor變換的定義同文獻(xiàn)[4]，即分析小波為Gabor小波的連續(xù)小波變換。

2D Gabor變換為：

(3)

(4)

在眼底視網(wǎng)膜圖像血管檢測中，Gabor變換參數(shù)設(shè)置大多同文獻(xiàn)[1]，ε=4，K0∈[0,3]，旋轉(zhuǎn)角度θ以10°為步長在10°～170°內(nèi)旋轉(zhuǎn)，并取不同方向變換的最大值為尺度a下的輸出響應(yīng)。膨脹尺度a的選擇應(yīng)與血管相匹配，包括可能的血管寬度。經(jīng)過實驗，本文取尺度a∈[2,5]，并將各尺度下的Gabor變換結(jié)果進(jìn)行疊加作為最后輸出。

1.2.3 B-COSFIRE

B-COSFIRE濾波模型由Azzopardi等[11]提出，這種模型模擬復(fù)雜感受野細(xì)胞的工作機(jī)制，對血管類的帶狀結(jié)構(gòu)有很好的檢測效果。在眼底視網(wǎng)膜血管檢測中，可以通過對稱和非對稱的B-COSFIRE濾波器來檢測出連續(xù)血管和血管末端，提取更完整的視網(wǎng)膜血管特征。B-COSFIRE濾波響應(yīng)就是由一組高斯差分(Different of Gaussian, DoG)濾波響應(yīng)的乘積加權(quán)幾何平均得到。

濾波公式為：

(5)

B-COSFIRE濾波器通過自動配置過程來實現(xiàn)其對管狀結(jié)構(gòu)的選擇性，圖3為B-COSFIRE濾波器自動配置模型。在濾波器自動配置過程中，給定一個標(biāo)準(zhǔn)差為σ的DoG濾波器，給定中心點，如圖3(a)、(b)中點1，該濾波器會沿著同心圓進(jìn)行濾波響應(yīng)，其中響應(yīng)最大的位置就是表征強度變化大的點，即感興趣的點，如圖3(a)中的點2、3、4、5，圖3(b)中的點2、3。這些感興趣點的數(shù)目取決于所考慮的同心圓的數(shù)目和指定的原型模型，并用(σi,ρi,φi)來表示這些感興趣的點，S={σi,ρi,φi|i=1,2,…,n}為這些點的集合，其中：σi為最大的DoG濾波器標(biāo)準(zhǔn)差，(ρi,φi)為這些點相對于中心點的極坐標(biāo)，n為DoG濾波器個數(shù)。

圖3 B-COSFIRE濾波器配置模型

為了提高各個點位置的容錯性，對DoG濾波響應(yīng)進(jìn)行模糊和移位操作。模糊操作即計算每一個位置極坐標(biāo)鄰域的DoG濾波器的最大加權(quán)閾值響應(yīng)，權(quán)重為DoG濾波響應(yīng)與高斯函數(shù)Gσ′(x′,y′)的系數(shù)乘積，σ′=σ0′+αρi，σ0、α為常數(shù)，ρi為與中心濾波器之間的距離。移位操作即將模糊后的DoG濾波響應(yīng)向φi的相反方向移動距離ρi，模糊移位后的DoG濾波響應(yīng)為：

yi′)Gσ′(x′,y′)}

(6)

其中：Δxi=-ρicosφi；Δyi=-ρisinφi；-3σ′≤x′，y′≤3σ′。

B-COSFIRE濾波響應(yīng)公式為：

(7)

1.2.4 相位特征

1)相位一致性特征。

對于一個由方波或三角波組成的傅里葉級數(shù)，在方波的階躍或三角波的波峰波谷處，其各次諧波分量的相位相同。故，相位一致性認(rèn)為圖像信號相位一致性最大處就是圖像特征[12]。在眼底視網(wǎng)圖像血管檢測中，相位一致性特征不受光照強度和圖像對比度的影響，且對于血管的邊緣有很好的檢測效果。

Kovesi[13]提出了二維相位一致性特征的計算公式：

(8)

其中：o表示方向的索引；n為尺度信息；To為噪聲補償；Ano為加權(quán)函數(shù)；ε為常數(shù)，避免分母為零。具體推導(dǎo)過程及參數(shù)設(shè)置見文獻(xiàn)[13]。

2)融合相位特征。

如圖4(a)所示，相位一致性特征對血管邊緣包括細(xì)微的分支血管都有很好的檢測效果，但是無法檢測出血管中心部分。故單一的相位一致性特征不適用于基于像素分類的眼底視網(wǎng)膜血管分割算法。為此，本文提出了一種新的融合相位特征，即將分別提取的相位一致性特征與Hessian矩陣特征進(jìn)行小波融合，融合后的相位特征不僅保留了良好的邊緣信息，其血管中心部分也被填充為白色，更適用于基于像素分類的視網(wǎng)膜血管的檢測。

在小波融合過程中，源圖像被分解成高頻部分和低頻部分，其中，高頻部分包括源圖像的細(xì)節(jié)和紋理信息，低頻部分包括源圖像的概貌和平均信息，因此，小波融合更容易提取圖像的結(jié)構(gòu)和細(xì)節(jié)信息。本文小波融合的具體步驟如下：

步驟1 對兩幅源圖像，如圖4(a)相位一致性特征、圖4(b)Hessian矩陣特征，分別進(jìn)行二層小波分解，得到各自的高頻子圖像系數(shù)與低頻子圖像系數(shù)。

步驟2 將兩幅源圖像的高頻系數(shù)和低頻系數(shù)按照加權(quán)平均的融合規(guī)則進(jìn)行融合，得到融合后的高頻系數(shù)與低頻系數(shù)。

步驟3 對融合后得到的系數(shù)進(jìn)行小波重構(gòu)，得到融合后圖像，如圖4(c)融合相位特征。融合相位特征即保留了相位一致性特征的良好血管邊緣信息，且血管中間部分被填充為白色。

圖4 相位特征融合過程

本文中提取的4D特征如圖5所示。對視網(wǎng)膜圖像中的每一個像素都構(gòu)造了4維的特征向量，用于后續(xù)的分類。

1.3 SVM分類

本文使用公共的用于血管提取的數(shù)字視網(wǎng)膜圖像(Digital Retinal Images for Vessel Extraction, DRIVE)數(shù)據(jù)庫中的訓(xùn)練與測試集合，通過LIBSVM工具箱進(jìn)行SVM訓(xùn)練和測試。

圖5 本文方法的各特征圖

訓(xùn)練過程 對DRIVE數(shù)據(jù)庫的訓(xùn)練集中的每幅視網(wǎng)膜圖像的每個像素構(gòu)造4維的特征向量，并分別提取1 500個血管點和非血管點，歸一化后，作為訓(xùn)練樣本。通過LIBSVM工具箱對SVM進(jìn)行訓(xùn)練，得到SVM模型。訓(xùn)練過程中選擇徑向基函數(shù)(Radial Basis Function, RBF)核函數(shù)作為核函數(shù)。為節(jié)省訓(xùn)練時間，本文中懲罰系數(shù)C選擇常用默認(rèn)值1，gamma選擇0.01，訓(xùn)練得到SVM模型。

測試過程 從DRIVE數(shù)據(jù)庫的測試集中任意選擇視網(wǎng)膜圖像并對其像素點提取特征向量，歸一化后，作為測試樣本，輸入訓(xùn)練后得到的SVM模型，輸出二值化的視網(wǎng)膜血管結(jié)構(gòu)圖。

1.4 后處理

采用連通域度量的方法，將輸出圖像中連通域面積小于25像素的孤立噪聲點去除。

2 實驗結(jié)果與分析

2.1 評價指標(biāo)

為了對算法的分割結(jié)果進(jìn)行評價，將本文的實驗結(jié)果與DRIVE數(shù)據(jù)庫中專家手工分割的結(jié)果進(jìn)行比較，使用準(zhǔn)確率(Accuracy, Acc)和受試者特性工作曲線面積(Area Under ROC Curve, AUC)對算法性能進(jìn)行定量分析。

準(zhǔn)確率是評價對視網(wǎng)膜圖像中血管點和非血管點的正確分割的概率，計算公式為：

(9)

其中：TP(True Positive)為真正類，TN(True Negative)為真負(fù)類，F(xiàn)P(False Positive)為假正類，F(xiàn)N(False Negative)為假負(fù)類。

受試者工作曲線(Receiver Operating Characteristic, ROC)是反映敏感性和特異性連續(xù)變量的綜合指標(biāo)，ROC曲線以真陽性率為縱坐標(biāo)，假陽性率為橫坐標(biāo)。受試者特性工作曲線面積(AUC)可由梯形規(guī)則計算得到，其值越接近1，表明算法的性能就越好。

2.2 實驗結(jié)果

本文以DRIVE數(shù)據(jù)庫中第1位專家手工分割結(jié)果為評價標(biāo)準(zhǔn)，實驗結(jié)果如表1所示，本文算法進(jìn)行血管分割所得平均Acc為0.957 4，平均AUC為0.970 2，結(jié)果較好。圖6為本文方法在DRIVE數(shù)據(jù)庫所有圖像中分割準(zhǔn)確率最高的最好分割結(jié)果，圖中大量血管均被分割出來，可用于臨床參考。圖7為本文方法在DRIVE數(shù)據(jù)庫所有圖像中分割準(zhǔn)確率最低的最壞分割結(jié)果，圖中血管主干部分被分割出來，在發(fā)生病變的區(qū)域，病變部分被錯誤檢測成血管。雖然圖7中出現(xiàn)成塊的錯誤血管區(qū)域，但這依然可以提醒醫(yī)生此處異常，為醫(yī)生進(jìn)行疾病檢測節(jié)省時間。

表1 本文算法分割結(jié)果

圖6 本文方法在DRIVE數(shù)據(jù)庫中分割準(zhǔn)確率最高的最好分割結(jié)果

圖7 本文方法在DRIVE數(shù)據(jù)庫中分割準(zhǔn)確率最低的最壞分割結(jié)果

為了分析融合相位特征用于像素分類的血管檢測的效果，本文進(jìn)行了不同特征組合的實驗，結(jié)果如表2所示。組合1提取了眼底視網(wǎng)膜圖像的Hessian矩陣、B-COSFIRE濾波和Gabor變換特征，構(gòu)造3維特征向量用于分類器分割，所得Acc和AUC分別為0.957 3、0.965 8，該特征組合分割效果好。組合2在上述3維特征組合中加入相位一致性特征，構(gòu)造4維的特征向量用于分類，但血管分割的Acc反而降低。因為相位一致性特征雖對血管邊緣檢測效果好，但其對血管中心部位檢測不足，直接使用相位一致性特征反而降低了分割的準(zhǔn)確率。組合3則在上述3維特征組合中加入融合相位特征，即本文方法，實驗得到的血管分割的Acc和AUC均有所提升。結(jié)果表明，在組合特征中，融合相位特征比相位一致性特征更適用于基于像素的視網(wǎng)膜血管分割。

表2 不同特征組合分割結(jié)果比較

此外，為進(jìn)一步驗證融合相位特征的有效性，本文將提取的單一的融合相位特征與相位一致性特征，分別通過SVM進(jìn)行分類，結(jié)果如表3所示，單一的相位一致性特征用于像素分類提取出血管的Acc和AUC分別為0.919 1、0.935 9，單一融合相位特征用于像素分類提取出血管的Acc和AUC分別為0.947 8和0.957 8，均高于相位一致性分割結(jié)果。實驗結(jié)果表明，融合的相位特征比相位一致性特征更適用于基于像素分割的視網(wǎng)膜血管分割算法。

表3 單一特征分割結(jié)果比較

為了分析本文算法的有效性，將本文算法與近期基于特征提取且算法準(zhǔn)確率較優(yōu)的文獻(xiàn)算法[6-9,11]，以及近期較熱且具有相似復(fù)雜度的基于神經(jīng)網(wǎng)絡(luò)的文獻(xiàn)算法[14-15]進(jìn)行比較。表4為本文算法與文獻(xiàn)[6-9,11,14-15]算法在DRIVE數(shù)據(jù)庫上分割準(zhǔn)確度的比較，圖8為本文算法與文獻(xiàn)[6,9]算法在DRIVE數(shù)據(jù)庫上分割結(jié)果圖的比較，圖9為本文算法與文獻(xiàn)[7,14]算法在DRIVE數(shù)據(jù)庫上分割結(jié)果圖的比較。如表4所示，與文獻(xiàn)[6-9,11,14-15]算法相比，本文算法分割所得到的平均Acc與AUC總體較高。文獻(xiàn)[9]算法的Acc較高，但其AUC低于本文算法，且如圖8所示，相較于文獻(xiàn)[9]分割得到的血管圖像，本文算法分割得到更多的細(xì)小血管且血管連續(xù)性較好。文獻(xiàn)[14]算法的AUC較高，但其Acc略低，且如圖9所示，相較于文獻(xiàn)[14]，本文分割雖然部分血管末端連續(xù)性不夠，但本文算法分割出更多的細(xì)小血管。故從總體指標(biāo)來看，本文算法結(jié)果較優(yōu)。此外，相較于文獻(xiàn)[8]算法構(gòu)造的39維特征向量、文獻(xiàn)[9]算法構(gòu)造的23維特征向量，文獻(xiàn)[14]算法采用的完全卷積神經(jīng)網(wǎng)絡(luò)，本文算法僅提取了4維特征通過SVM分類，降低了算法的復(fù)雜度，提高了分割效率。

表4 不同算法分割結(jié)果比較

3 結(jié)語

針對相位一致性特征對血管中心檢測不足的問題，本文提出基于融合相位特征的視網(wǎng)膜血管分割算法。不同于文獻(xiàn)[8-9]算法，本文算法沒有直接使用相位一致性特征，而是將相位一致性特征同Hessian矩陣特征進(jìn)行小波融合，得到一種新的融合相位特征。該特征既保留了相位一致性特征中良好的血管邊緣信息，又克服了相位一致性特征中血管中心部位與背景相似的缺點。由于基于像素分類的視網(wǎng)膜血管分割算法需要對每一個像素的特征進(jìn)行識別分類，因此，本文中的融合相位特征比相位一致性特征更適用于基于像素分類的血管分割算法。此外，本文對每一個像素提取Hessian矩陣、Gabor小波、B-COSFIRE、相位特征，構(gòu)造4維特征向量用于SVM分類，提取較少的特征與樣本，實現(xiàn)快速、高效的眼底視網(wǎng)膜血管的分割；但是，與現(xiàn)有大多的眼底視網(wǎng)膜血管分割算法一樣，本文對于眼底視網(wǎng)膜中毛細(xì)血管的分割不足，損失了一些細(xì)小血管的信息。在接下來的工作中，如何對細(xì)微小的血管進(jìn)行更高準(zhǔn)確度的分割仍需要繼續(xù)研究。

圖8 文獻(xiàn)[6,9]算法與本文算法結(jié)果比較

圖9 文獻(xiàn)[7,14]算法與本文算法結(jié)果比較