楊 剛 徐學(xué)軍

(長沙理工大學(xué)物理與電子科學(xué)學(xué)院 湖南 長沙 410000)

0 引 言

無論是靜脈穿刺,還是中醫(yī)針灸,臨床上,醫(yī)護(hù)人員主要通過目測或經(jīng)驗(yàn)定位血管,但對部分血管過細(xì)、脂肪較厚或血管塌陷的特殊患者(例如兒童、體胖者、老人等),血管難以辨別,嚴(yán)重影響靜脈穿刺成功率[1]。而在中醫(yī)針灸中,醫(yī)生需要將針灸針一根一根施于患者的相應(yīng)穴位中,在施針的時(shí)候,醫(yī)生用雙眼很難精確地去確定每個(gè)穴位的具體位置,同時(shí)還得避免扎錯(cuò)地方,若是扎到了血管或其他要害部位,還會(huì)危害患者的健康[2-3]。為了使靜脈能夠被準(zhǔn)確識別,人們通過光照顯影、紅外顯像、靜脈投影等輔助方式,提高靜脈精準(zhǔn)定位。而近年來根據(jù)相關(guān)研究,主要利用近紅外成像技術(shù),來區(qū)分靜脈和其他組織,有較好的實(shí)施效果和適用性[1],其靜脈圖像提取系統(tǒng)如圖1所示。

圖1 近紅外靜脈圖像提取系統(tǒng)

然而,基于近紅外成像技術(shù)采集的靜脈圖像,由于受到光照的不均勻性和個(gè)體差異性的影響,導(dǎo)致近紅外圖像出現(xiàn)對比度低、邊緣模糊和細(xì)節(jié)信息丟失等缺陷,造成了圖像分割后出現(xiàn)很多虛假的靜脈信息。因此,如何增強(qiáng)靜脈圖像的對比度以及將靜脈紋路與皮膚背景分割開來成為了研究的熱點(diǎn)。魯周迅等[4]提出了基于子直方圖均衡化算法,結(jié)合NiBlack分割方法對圖像進(jìn)行靜脈提取。林劍等[5]提出了基于Hessian矩陣的靜脈圖像分割方法,通過Hessian矩陣的跡初次濾除了非靜脈區(qū)的像素點(diǎn)。王定漢等[6]提出基于靜脈灰度值特征的圖像分割與對比度增強(qiáng)算法。Malutan等[7]提出了將LLBP描述子與形態(tài)學(xué)算子相結(jié)合用于提取靜脈特征,有效地提取出了靜脈紋理特征。以上算法都取得了一定的效果,但是也存在一些不足。因此,探索新的靜脈圖像的對比度增強(qiáng)和分割方法是非常必要的。

由于人體靜脈特征網(wǎng)絡(luò)具有復(fù)雜性和方向高度隨機(jī)性,本文算法在基于方向?yàn)V波器的基礎(chǔ)上,首先提出多方向Gabor濾波算法進(jìn)行靜脈圖像增強(qiáng),然后將每個(gè)方向所得的靜脈圖像進(jìn)行重組得到增強(qiáng)后的圖像。最后利用馬爾可夫隨機(jī)場分割方法對靜脈網(wǎng)絡(luò)進(jìn)行分割,該靜脈特征提取流程如圖2所示。實(shí)驗(yàn)結(jié)果表明,該方法能夠得到可靠有效的靜脈網(wǎng)絡(luò)特征。

1 靜脈圖像增強(qiáng)算法

靜脈沒有虹膜那么豐富的紋理,也不像指紋具有較多的細(xì)節(jié)點(diǎn)特征,其主要體現(xiàn)在具有復(fù)雜變化關(guān)系的血管網(wǎng)絡(luò)結(jié)構(gòu)中[8]。人體靜脈特征通常要通過靜脈圖像區(qū)域分割來提取,而靜脈圖像分割的好壞會(huì)直接影響后續(xù)的靜脈特征提取質(zhì)量。因此,本文提出一種基于多方向的Gabor濾波器的靜脈圖像增強(qiáng)算法。

Gabor濾波是一種眾所周知的特征提取方法,已在機(jī)器視覺領(lǐng)域得到廣泛研究和應(yīng)用[9]。Gabor 濾波器對于方向和角度的表達(dá)類似于人類視覺系統(tǒng)對于方向和角度的表達(dá)[1],而在本文中引入的多方向Gabor濾波器目的是能夠檢測不同方向的靜脈血管。在空間域中,一般二維(2D)Gabor濾波器被定義為包含復(fù)正弦信號的高斯函數(shù)[10]。其表達(dá)式如下:

(1)

(2)

式中:G(x,y)是像素的坐標(biāo);σx、σy分別為沿著x和y方向的高斯包絡(luò)的標(biāo)準(zhǔn)偏差,xθ、yθ表示二維Gabor濾波器旋轉(zhuǎn)θ角后更新后的x、y坐標(biāo);*表示卷積運(yùn)算符;f0代表二維Gabor濾波器的空間中心頻率。

Gabor濾波函數(shù)又分為實(shí)部和虛部,經(jīng)研究發(fā)現(xiàn),Gabor濾波器的實(shí)部也叫做偶對稱 Gabor濾波器,特別適合圖像中脊線的檢測,且增強(qiáng)了Gabor濾波器處理時(shí)間的有效性[1]。偶對稱Gabor濾波器如式(3)所示。

(3)

基于前人的經(jīng)驗(yàn)可知, 當(dāng)濾波器的方向角θ個(gè)數(shù)超過8個(gè)時(shí), 對靜脈的增強(qiáng)效果很差,故在選擇方向個(gè)數(shù)的時(shí)候, 通常會(huì)選擇4個(gè)或6個(gè)方向,一般來說, 利用Gabor濾波器的多個(gè)方向可以挖掘出更多的靜脈信息。然而, 選擇方向過多時(shí)會(huì)產(chǎn)生更多的冗余信息,這會(huì)大大降低Gabor 濾波器的性能,選擇方向較少時(shí)又不能挖掘出更多的靜脈信息[8]。所以本文選擇6個(gè)方向,其中i(i=1,2,3,4,5,6)表示6個(gè)方向的索引值,對應(yīng)的6個(gè)方向θi為π/6、2π/6、3π/6、4π/6、5π/6。fi表示根據(jù)第i個(gè)方向的偶數(shù)對稱Gabor濾波器的空間中心頻率,其值的選擇決定了輸出濾波圖像的質(zhì)量。圖3顯示了6個(gè)方向的偶對稱Gabor濾波器的空間響應(yīng)。

圖3 6個(gè)方向的偶對稱Gabor濾波器的空間響應(yīng)

設(shè)I(x,y)為一幅近紅外成像設(shè)備獲得的靜脈圖像,將每一方向?qū)?yīng)的偶對稱Gabor濾波函數(shù)與原圖像進(jìn)行卷積運(yùn)算,獲得濾波后的靜脈圖像Oi(x,y)如下:

(4)

由于近紅外靜脈圖像中靜脈線比皮膚區(qū)域暗,因此,選擇6幅圖像中每個(gè)位置像素點(diǎn)灰度值最低的點(diǎn)構(gòu)成濾波后的靜脈圖像,如圖4所示。

圖4 濾波后的靜脈圖像

2 MRF圖像分割

眾所周知,靜脈圖像分割提供了具有分離的靜脈區(qū)域的圖案分割圖像,這對于特征提取很重要[11]。對于靜脈圖像而言,分割使得靜脈紋路和皮膚背景區(qū)分開來,以便于特征提取等后續(xù)處理。傳統(tǒng)的圖像分割方法主要分為基于像素的方法、基于區(qū)域的方法、基于邊界的方法、基于模型的方法和多種組合算法[12]。而本文在使用Gabor濾波器靜脈圖像增強(qiáng)后,采用MRF隨機(jī)場用于靜脈圖像分割。實(shí)驗(yàn)表明,該算法在靜脈圖像分割中取得很好的分割效果。

基于馬爾可夫隨機(jī)場理論是指圖像中的像素點(diǎn)與其相鄰像素相關(guān),即相鄰像素的強(qiáng)度值變化不大。任何數(shù)字圖像都是由離散的像素集合構(gòu)成的,這些像素可以被建模為一個(gè)隨機(jī)場。圖像中的每個(gè)像素都是一個(gè)站點(diǎn),每個(gè)站點(diǎn)都有一個(gè)標(biāo)簽X,它通常表示一個(gè)像素的強(qiáng)度值[13]。

由上文可知靜脈紋路暗,皮膚背景亮,可分別記為標(biāo)簽“0”和“1”,這樣將大小為m×n的整幅圖像分為了兩類,表示為:

W={w1,w2,…,wm×n}

式中:W為所得點(diǎn)集;w(w=0,1)表示圖像各點(diǎn)像素值為s所屬類。

為了使得分割圖像更精準(zhǔn),我們需要將像素值更好地分為所屬類。P(w|s)表示在已知像素信息s的條件下屬于一個(gè)標(biāo)簽的概率,也叫后驗(yàn)概率。因此,該過程就是通過最大化后驗(yàn)概率來尋找像素點(diǎn)s所屬最優(yōu)標(biāo)簽w。由于后驗(yàn)分布是一個(gè)條件分布,我們可以用Bayes規(guī)則來估計(jì)它,根據(jù)這個(gè)規(guī)則可得后驗(yàn)概率表達(dá)式[14]:

(5)

式中:P(w)是先驗(yàn)概率,在此認(rèn)為其符合MRF模型;P(s|w)稱為似然函數(shù);P(s)是一個(gè)常量。

根據(jù)Hammersley-Clifford[15]定理,用 Gibbs 隨機(jī)場的概率密度函數(shù)描述馬爾可夫隨機(jī)場,則先驗(yàn)概率P(w)表示為:

(6)

(7)

本文采用雙點(diǎn)勢團(tuán)的 8 鄰域系統(tǒng),像素點(diǎn)t為像素點(diǎn)鄰域s內(nèi)的點(diǎn)。其中β為耦合系數(shù),控制著勢團(tuán)內(nèi)兩點(diǎn)相互作用的強(qiáng)弱,一般取值在[0.5,1.0]區(qū)間之間,耦合系數(shù)β越小,對圖像的分割越細(xì)膩,隨著β的增大,對細(xì)碎紋理的檢測能力逐漸降低[16]。

而P(s|w)為似然函數(shù),表示在已知分類標(biāo)簽w情況下,點(diǎn)像素灰度值為s的概率。因此,我們需要將圖像隨機(jī)初始化分類為w,并將屬于w類的點(diǎn)找出來。一般認(rèn)為每一類里面的所有點(diǎn)分布是具有不同均值和標(biāo)準(zhǔn)差的高斯函數(shù)[16],對每一類像素點(diǎn)建立高斯密度函數(shù)。對于靜脈圖像而言,我們得到兩類高斯密度函數(shù)。

那么已知像素點(diǎn)灰度值則能計(jì)算出該像素點(diǎn)分別屬于“0”“1”類的概率值P(s|w0)、P(s|w1)。由式(5)可知,當(dāng)P(s|w)×P(w)達(dá)到最大值時(shí),圖像分割效果最理想。具體分割算法描述如下:

(1) 獲取Gabor濾波增強(qiáng)后的靜脈圖像。

(2) 使用k-means聚類算法初始化圖像標(biāo)簽“0”和“1”。

該礦體為浸染(星點(diǎn))狀輝鉬礦，礦體呈灰白色—淺肉紅色—褐紅色，地表出露長度為142 m，厚度1.72～2.86 m，產(chǎn)狀310°∠69°。礦體規(guī)模不大，沿走向露頭不好。礦石金屬礦物主要有黃鐵礦、褐鐵礦、輝鉬礦、鉬華，蝕變主要為硅化。脈石礦物主要為二長花崗巖及碎裂石英巖。圍巖為二長花崗巖，礦體與圍巖界線不明顯，沿礦體兩側(cè)圍巖產(chǎn)生鉀化、硅化蝕變。礦石鉬品位0.038%～0.137%。

(3) 統(tǒng)計(jì)每個(gè)標(biāo)號場8領(lǐng)域的相同標(biāo)簽個(gè)數(shù),獲得圖像標(biāo)號場的先驗(yàn)概率。

(4) 計(jì)算似然函數(shù)P(s|w)。

(5) 求取最大化后驗(yàn)概率,并更新標(biāo)號場。

(6) 判斷是否達(dá)到最大迭代次數(shù),若沒有,重復(fù)步驟(2)-步驟(4),否則,算法結(jié)束。

本文在 Windows 10、64位操作系統(tǒng)下,通過 MATLAB 2014a實(shí)現(xiàn)上述算法。相關(guān)函數(shù)來自MATLAB自帶函數(shù)庫,部分主要代碼如下:

Class_num=2;

%設(shè)置分類數(shù)

Max_iter=20;

%最大迭代次數(shù)

Label=kmeans(img(:),Class_num);

for i=1:Class_num

Label_i=i*ones(size(label));

Temp=～(label_i-label_u)+～(label_i-

Label_d)……;

p_c=(i,:)=temp(:)/8;

%計(jì)算最大后驗(yàn)概率

end

for j=i:Class_num

MU=repmat(mu(j),size(img,1)*

size(img,2),1);

p_sc(j,:)=1/sqrt(2*pi*sigma(j))*

exp(-(img(:)-MU).^2/2/sigma(j));

end

%計(jì)算似然概率P(s|w)

[M,label]=max(log(p_c)+log(p_sc));

%尋找最大后驗(yàn)概率

分割后的靜脈圖像如圖5所示。

圖5 分割后的靜脈圖像

3 實(shí)驗(yàn)結(jié)果分析

為了驗(yàn)證本文方法的有效性,實(shí)驗(yàn)采用Intel(R)Core(TM)i5-5200U CPU @2.20 GHz 2.20 GHz,內(nèi)存為8 GB的PC機(jī),在MATLAB 2014 a環(huán)境下完成的。

本實(shí)驗(yàn)選擇了3幅實(shí)驗(yàn)室自行設(shè)計(jì)的近紅外靜脈成像設(shè)備采集到的噪聲未經(jīng)處理的靜脈紋路高度隨機(jī)性的分辨率為640×480像素的靜脈圖像進(jìn)行測驗(yàn)。在進(jìn)行Gabor濾波時(shí),為避免圖像欠增強(qiáng)和產(chǎn)生更多的冗余信息,因此將這3幅低復(fù)雜度背景的靜脈圖像選用六方向分別為30°、60°、90°、120°、150°、180°,x、y方向的高斯包絡(luò)偏差σx、σy均為2的Gabor濾波器進(jìn)行靜脈圖像增強(qiáng),如圖6所示,經(jīng)實(shí)驗(yàn)所得靜脈圖像得到了顯著增強(qiáng),克服了GHE等圖像增強(qiáng)算法出現(xiàn)的局部過增強(qiáng)和欠增強(qiáng)等問題,同時(shí)也很好地抑制了處理過程中產(chǎn)生的噪聲、邊緣偽影等信息,為后續(xù)進(jìn)行靜脈網(wǎng)絡(luò)特征提取提供了好的靜脈環(huán)境。

圖6 實(shí)驗(yàn)所得靜脈圖像

其次,得到增強(qiáng)后的靜脈圖像用馬爾可夫隨機(jī)場進(jìn)行分割,耦合系數(shù)β越小,對靜脈圖像的紋理與背景分離越明顯,隨著β的增大,許多靜脈細(xì)節(jié)信息丟失。在進(jìn)行中值濾波去毛刺處理時(shí),經(jīng)多次實(shí)驗(yàn)對比,選擇窗口為5×5的濾波模板時(shí)去毛刺效果最好。將本文方法與其他三種常用靜脈分割方法進(jìn)行比較,結(jié)果如圖7所示。Otsu分割算法對靜脈圖像做分割時(shí),分割后的圖像仍然保留了過多的噪聲。使用Niblack分割算法得到的二值圖像中對于皮膚背景信息較復(fù)雜的靜脈圖像而言,容易引進(jìn)較多的偽靜脈信息。K-means聚類算法對于較高質(zhì)量的靜脈圖像能夠提取到很好的靜脈信息,且殘留的噪聲較少,對于復(fù)雜背景的靜脈圖像分割依舊有噪聲存在。

圖7 各種算法對比

本文提出的靜脈圖像分割方法對低背景復(fù)雜度的靜脈圖像做分割時(shí),能夠提取到更好地靜脈信息,沒有噪聲、偽靜脈信息的干擾。而對于背景信息相對復(fù)雜的靜脈圖像而言,能夠更好地抑制噪聲信息的產(chǎn)生。采用主觀評價(jià)法通過人眼來比較不同的分割結(jié)果圖像,其穩(wěn)定性較差,難以定量描述,在應(yīng)用中受到很大限制。為了更好地評價(jià)本文提出的靜脈特征提取模型性能,因此選擇無監(jiān)督評價(jià)指標(biāo)中的區(qū)域內(nèi)一致性指標(biāo)對圖像分割性能進(jìn)行客觀評價(jià)。而區(qū)域一致性指標(biāo)又分為最大對比度評價(jià)和方差評價(jià)一致性,本文使用最大對比度評價(jià)(Zeb)[18],結(jié)果如表1所示?？梢钥闯鍪褂帽疚乃惴ㄋ肸eb值相比幾種傳統(tǒng)分割算法而言,其值更小,所得靜脈分割區(qū)域內(nèi)一致性更好,所得分割圖像效果更好。因此本文方法提取出來的靜脈較為完整,提取的靜脈結(jié)構(gòu)非常適合原始圖像的血管拓?fù)浣Y(jié)構(gòu)。

表1 不同分割算法的最大對比度評價(jià)指標(biāo)

4 結(jié) 語

本文為了解決利用近紅外成像設(shè)備采集到的靜脈所存在的問題,提出一種基于多方向Gabor濾波器的靜脈圖像增強(qiáng)方法,系統(tǒng)有效地解決了靜脈信息的噪聲大、對比度低、圖像清晰度不夠的實(shí)際困難,主要采取了兩種技術(shù)手段。(1) 根據(jù)靜脈的紋理特征,采用了六維度Gabor濾波器進(jìn)行濾波,獲得了滿意的增強(qiáng)靜脈圖像;(2) 為了實(shí)現(xiàn)靜脈特征提取和針灸穴位精準(zhǔn)定位,用馬爾可夫隨機(jī)場對增強(qiáng)后的靜脈進(jìn)行二值分割,并采用中濾波取毛刺處理。通過范例實(shí)驗(yàn)結(jié)果表明,本文算法能夠有效精準(zhǔn)地提取皮下近紅外靜脈特征和針灸穴位,可靠地增強(qiáng)靜脈圖案質(zhì)量,為中醫(yī)針灸、西醫(yī)穿刺、生物特征識別等領(lǐng)域提供了新思路。