胡 敏，梁 挺，張光華，黃宏程

(1.重慶郵電大學(xué) 通信與信息工程學(xué)院，重慶 400065；2.太原學(xué)院 計算機科學(xué)與工程系，山西 太原 030000)

0 引 言

甲狀腺相關(guān)性眼病(thyroid associated ophthalmopathy，TAO)是一種與自身免疫有關(guān)的眼部疾病，是一種常見的眼部疾病，占Graves病的40%～50%，免疫性甲亢的2～5%。關(guān)于眼瞼退縮等定義見參考文獻[1-3]。目前已知的眼瞼與虹膜外圓的位置信息可分為全局和局部兩種。

目前國內(nèi)外有關(guān)虹膜位置的相關(guān)文獻多以某一特定的虹膜采集設(shè)備為基礎(chǔ)，收集到的圖像如下：①僅包括人類的眼睛；②對虹膜及瞳仁進行了擴大；③增加了對眼、眉等的干擾作用；④虹膜、瞳孔和鞏膜有顯著的區(qū)別。傳統(tǒng)的、通用的人臉圖片都是由便攜式設(shè)備獲取的，例如手機和普通相機，而不是由儀器獲取，已有的研究手段無法對其進行有效應(yīng)用，因此，開展基于面部的虹膜定位的技術(shù)十分重要。人類的虹膜和瞳孔都位于眼眶的正中央，上眼瞼的形狀和位置都有相同的特點，而瞳孔、虹膜、眼瞼等都比正常的瞳孔、虹膜、眼瞼等情況更加復(fù)雜。首先，要解決TAO病人影像上的虹膜位置問題，然后將此算法應(yīng)用到便攜式上。

全采集采用便攜式設(shè)備，通常為全貌、局部采集使用專業(yè)的虹膜儀器，僅含有眼睛的部分信息。虹膜的定位一般是確定虹膜的內(nèi)外邊緣和上、下4個分界線[4-6]。基于全局采集的虹膜定位算法包括：Daugman等[5]采用積分微分算子，通過圓型模板對虹膜的中心及半徑進行搜索，得到了較好的定位效果。Wilde等[7]提出了利用霍夫變換進行虹膜內(nèi)、外邊界的方法，但是要對整個圖像進行檢索，不僅費時費力，而且耗費較大。目前，基于局部提取的虹膜定位方法主要有：田子林等[8]所提出的基于最小二乘法的虹膜定位方法；王琪等[9,10]所提出的基于SDM和SIFT的虹膜定位算法，通過對上瞼和虹膜外緣的比例變換不變特征進行定位，從而使其定位精度得到了提高，但其計算量大、耗時長。采用特殊的攝影器材進行局部采集，可以放大眼睛的特點，但由于設(shè)備的成本較高，不宜大規(guī)模購買；而且，近距離獲得的眼球信息不能達到眼瞼退縮的判別準(zhǔn)則。

為了解決上述問題，本文在利用便攜式攝像機采集到的完整的面部影像基礎(chǔ)上，采用區(qū)域霍夫變換進行虹膜定位，并與EUGOGO眼瞼退縮準(zhǔn)則相結(jié)合。實驗結(jié)果表明，該算法能夠很好地實現(xiàn)對圖像的輔助診斷。

1 區(qū)域Hough變換的虹膜定位方法

1.1 特征點提取

樣本圖像由便攜式設(shè)備拍攝，以第一眼位圖為標(biāo)準(zhǔn)(當(dāng)兩眼在真正水平面注視6 m以上的目標(biāo)時，兩眼注視線同時向前方相并平行)，背景顏色選擇白色，圖像包含了從頭發(fā)到雙肩之間的整個面部區(qū)域。面部對齊是根據(jù)圖像中的區(qū)域和尺寸來確定面部形狀[11]。臉部特征點是面部對齊的重要指標(biāo)，被看作是學(xué)習(xí)一個回歸函數(shù)的過程[12,13]，將原始圖像作為輸入，利用多個簡化的小波變換算法對所抽取的特征進行分類，目標(biāo)的輸入取決于上一步的輸出，反復(fù)地進行殘差擬合直至獲得最優(yōu)的估計形態(tài)。初始的特征點定位共包含68個特征點，并沒有覆蓋額頭區(qū)域，如圖1所示。

圖1 68特征點效果

在此基礎(chǔ)上，本文利用拓展后的81特征點法，如圖2所示，該方法在前額頭上新增13個特征點，彌補了68個特征點的不足之處。所有的特征點依次由1-81標(biāo)記，圖2中1至17的標(biāo)號與69至81的標(biāo)號，30個特征點共同組成面部外輪廓，眼瞼部位一次標(biāo)號為37至48。

圖2 81特征點效果

1.2 Bezier曲線擬合

本文采用Bezier曲線擬合面部輪廓特征點，截取面部圖像，排除頭發(fā)、背景等干擾因素。Bezier曲線是法國工程師Bezier[14]在1962年為了設(shè)計汽車車身形狀提出的，應(yīng)用于二維圖形程序，在計算機圖形學(xué)中是相當(dāng)重要的參數(shù)曲線。與最小二乘法等擬合方法相比，其優(yōu)勢在于Bezier曲線實現(xiàn)簡便、需要的條件更加簡單且擬合后的曲線更光滑。

貝塞爾曲線的定義嚴(yán)格依賴于確定該段曲線控制點的個數(shù)，控制點也就是曲線上的點，N次多項式的曲線由N+1個頂點確定。貝塞爾曲線上各點的參數(shù)方程為

(1)

其中，Pi為第i個頂點的坐標(biāo)值，基函數(shù)Bi,N(t) 為Bernstein多項式，該多項式公式為

(2)

由上面兩式可得Bezier曲線，根據(jù)N+1個頂點 (A0,A1,…,An)， 得到其貝塞爾曲線公式為

(3)

該公式的原始數(shù)學(xué)表達式為伯恩斯坦多項式，即在 [A0,An] 區(qū)間上所有的連續(xù)函數(shù)都可以用多項式來逼近，并屬于一致收斂。曲線的起點是A0，終點是An,在所有控制點都是同線的情況，曲線表示為一條從A0到An的直線。該曲線可以在任何控制點上進行分割，每條被分割的曲線都屬于Bezier曲線，其兩端的點分別為分割處的控制點。

在1.1節(jié)中，圖2通過標(biāo)號1-81標(biāo)注了所有的特征點，面部最外一圈的特征點構(gòu)成了擬合區(qū)域，包含了標(biāo)號1-17，69-81，一共30個特征點。將30個特征點用A0-A29表示，即用于人臉特征點擬合的Bezier公式為

B(t)=P0(1-t)30t0+30P1(1-t)29t+…+
30P29(1-t)t29,t∈[0,1]

(4)

通過式(4)的擬合效果如圖3所示，圖3閉合線條區(qū)域為30個特征點擬合后的效果，該區(qū)域切割后的效果如圖4所示。對比圖3和圖4可知經(jīng)過Bezier曲線擬合后頭發(fā)、雙肩以及左右耳朵等明顯的干擾區(qū)域被處理。

圖3 Bezier曲線擬合效果

圖4 處理后效果

1.3 三次B樣條插值

本文應(yīng)用B樣條插值方法[15]獲得了眼瞼輪廓曲線，并充分利用了它的局部支撐性和插值足夠光滑的特性優(yōu)勢，從而獲取了更加精準(zhǔn)的上下眼瞼輪廓曲線。B樣條曲線包含一系列的控制點，并將其劃分為節(jié)點和系數(shù)，每一個系數(shù)對應(yīng)一個控制點，其核心思想是采用分段低階多項式將高次多項式進行連續(xù)拼接。由于Bezier曲線不能進行局部修改，而且在拼接成長線時，會產(chǎn)生一些波形，對結(jié)果產(chǎn)生影響，因此，我們選擇了B樣條進行輪廓擬合；而上下眼皮分別用4個點來進行定位，它們符合三次B樣條曲線的對應(yīng)條件，也就是三次B樣條曲線的插值。

在對上下眼瞼進行劃分時，左右兩個眼部被均分為上下兩個部分，包括4條眼瞼線，分別稱之為左上眼瞼、左下眼瞼、右上眼瞼、右下眼瞼。在1.1節(jié)圖2中，可將4條眼瞼分別對應(yīng)4個標(biāo)號，即三次B樣條插值所需的4個控制點。左上眼瞼對應(yīng)的4個特征點標(biāo)號為：(37、38、39、40)。左下眼瞼對應(yīng)的4個特征點標(biāo)號為：(37、42、41、40)。右上和右下分別對應(yīng)的標(biāo)號為：(43、44、45、46),(43、48、47、46)。

B樣條曲線的總方程為

(5)

其中，Pi是控制曲線的特征點，F(xiàn)i,k(t) 則是K階B樣條基函數(shù)。

三次B樣條曲線方程中基函數(shù)為

(6)

(7)

(8)

(9)

(10)

由式(6)到式(10)可得三次B樣條曲線方程

P(t)=P0*F0,3(t)+P1*F1,3(t)+
P2*F2,3(t)+P3*F3,3(t)

(11)

下面是三次B樣條曲線的擬合法，如圖5所示，圖5中封閉線條區(qū)域表示上下眼瞼經(jīng)過式(11)插值后的效果區(qū)域。結(jié)果表明，采用三次B樣條插值算法進行上下眼瞼邊緣定位有很好的效果，同時有很高的準(zhǔn)確率。

圖5 眼瞼擬合

1.4 Hough圓檢測

Hough變換是Paul Hough在1962年首次提出的一種線條檢查算法，在1972年，經(jīng)過Richard Duda和Peter Hart推廣使用，經(jīng)典的Hough變換通常被用來檢測圖像中的直線，經(jīng)過眾多學(xué)者的不斷研究和發(fā)展，Hough變換被擴展到可以進行任意形狀物體的檢測，現(xiàn)應(yīng)用最為廣泛的形狀是圓和橢圓的檢測[16,17]。人眼的虹膜區(qū)域通常位于鞏膜和瞳孔之間，其外邊界和鞏膜相切，實際為眼白和棕色線。虹膜在形狀上無限接近于圓形，并且具有不會輕易受到其它因素影響而改變形狀的特性。針對虹膜這種特性，本文利用Hough變換圓檢測方法進行虹膜外圓的定位。

Hough變換圓檢測原理是對于一個半徑為r，圓心坐標(biāo)位置為 (a,b) 的圓，將其表示為

(x-a)2+(y-b)2=r2

(12)

圓檢測方法步驟如下：

(1)建立一個由各個象素單位組成的累積空間。開始時，將每一個單元格都設(shè)定為0。

(2)針對各圖象的邊界點座標(biāo)位置(x,y)，可以是一個圓形的中央的網(wǎng)孔數(shù)值，如方程式(12)所示。方程式中的單位用a代表。

(3)在上述的步驟中，從各個可能發(fā)現(xiàn)的數(shù)值a中，得到符合方程的全部可能的數(shù)值b。

(4)在累積空間內(nèi)搜尋局部極大值。這些單元格表示算法檢測到的所有圓圈。

以上所述的圓檢法所涉及到的參數(shù)空間是三維的，若要在三維空間上進行證據(jù)累計，則需要耗費大量的時間和空間復(fù)雜性。因此本文對該Hough圓檢測方法進行改進，主要考慮優(yōu)化空間開銷和時間開銷。為此，對這種Hough圓法進行了改進，重點是在優(yōu)化的空間和時間上。本文提出了一種基于局部霍夫圓法的新算法，該算法可分為兩個階段：

(1)選擇目標(biāo)區(qū)域?；贐ezier曲線方法，擬合圖像，然后根據(jù)人臉五官的分布，獲取目標(biāo)眼部區(qū)域信息。

(2)確定搜索區(qū)域范圍。首先，利用改進后的Canny算子對圖像進行邊緣信息提取，然后對搜索區(qū)域進行去噪處理，這樣更好地突出虹膜外圓的線條信息，并規(guī)定了經(jīng)典Hough圓檢測方法中的搜索半徑r指定范圍值，而非盲目進行區(qū)域搜索，避免了三維的參數(shù)空間，只需要在目標(biāo)區(qū)域內(nèi)對檢測圓心(a,b)進行證據(jù)累計，從而縮短了時間和空間，增加了虹膜外圓定位的準(zhǔn)確率。

目標(biāo)區(qū)域選擇通過縮小Hough圓圈來探測和降低空間消耗。該算法采用了一種改進的Canny邊界提取算法來突出圖像的邊緣，同時根據(jù)虹膜的直徑來限定搜索半徑r的范圍，從而降低了搜索的耗時。

1.4.1 選擇目標(biāo)區(qū)域

“三庭五眼”的概念最早在明末王繹的《寫像秘訣》中記錄[18]。“大三庭”通常是把頭部上下分為三等分，從前額發(fā)際線至眉骨，從眉骨至鼻底，從鼻底至下頦，各占臉長的三分之一。“五眼”從左側(cè)發(fā)際至右側(cè)發(fā)際，為5只眼形。兩只眼睛之間有一只眼睛的間距，兩眼外側(cè)至側(cè)發(fā)際各為一只眼睛的間距，各占比例的1/5。如圖6所示。

圖6 “三庭五眼”

本文截取眼部區(qū)域位于圖6中的“中庭”處。直接取“中庭”區(qū)域為眼部區(qū)域會包含部分眉毛信息，為了避免其干擾，本文通過大量的圖片測量和驗證，確定目標(biāo)區(qū)域為眉線以下六分之一以及鼻底線以上二分之一。同時對目標(biāo)區(qū)域進行歸一化處理，統(tǒng)一像素值為278×59。以4位TAO患者圖像為例，截取目標(biāo)區(qū)域圖，分別如圖7(a)、圖7(b)、圖7(c)、圖7(d)所示。

圖7 TAO患者目標(biāo)區(qū)域

1.4.2 確定搜索區(qū)域范圍

搜索區(qū)域處理主要分為兩個部分：

(1)區(qū)域圈定：Canny算子能檢測真正的弱邊緣，可使用兩種不同的閾值分別檢測強邊緣和弱邊緣，不容易受噪聲的干擾，具有準(zhǔn)確的邊緣定位功能，所以本文采用Canny算子對得到的目標(biāo)區(qū)域進行邊緣檢測。經(jīng)過Canny算子處理后的圖像存在無用的小面積區(qū)域，本文將結(jié)合連通域處理方法進行圖像篩選，過濾無用區(qū)域。改進后的Canny算子工作流程如圖8所示。

圖8 改進后Canny算子工作流程

圖9 原始Canny算子處理后效果

圖10 改進后Canny算子處理效果

圖9(a)、圖9(b)、圖9(c)、圖9(d)這4幅圖像經(jīng)過原始Canny算子處理后眼瞼與虹膜基本輪廓較為明顯，每位患者的左右眼虹膜區(qū)域內(nèi)存在少數(shù)小圓圈(小圓圈是拍攝時人眼部位的反光聚焦點，表現(xiàn)為白色亮斑，經(jīng)Canny算子處理后表現(xiàn)為小圓圈)，眼瞼處存在斷斷續(xù)續(xù)的不連接線段。由于TAO患者眼部區(qū)域的特殊性和每個人眼部區(qū)域的差異性，部分圖像經(jīng)過Canny算子檢測后虹膜輪廓不能呈現(xiàn)出完整的輪廓信息。圖10(a)、圖10(b)、圖10(c)、圖10(d)是經(jīng)過改進后Canny算子處理結(jié)果圖。對比圖9和圖10，小圓圈和不連接線段均被去除，保留了大部分虹膜外圓的邊緣信息。

(2)區(qū)域搜索：由于人類的虹膜直徑比較相近，平均虹膜直徑12 mm，而中國人的平均虹膜直徑大約為11.4 mm。本文對隨機的500張TAO患者圖像和100張非TAO患者圖像進行虹膜直徑的估計，98.6%的TAO患者眼直徑在10 mm～13 mm范圍內(nèi)，1.4%的TAO患者眼直徑在該范圍外，所有的非TAO患者圖像眼徑在10.5 mm到11.6 mm之間。直徑分布見表1。

表1 直徑分布

根據(jù)直徑分布表信息，1000張TAO患者圖像中只有14張圖像的眼直徑范圍在10 mm～13 mm外，即將搜索半徑r的值縮小在10 mm～13 mm內(nèi)。同樣以4位TAO患者圖像為例，原始Hough圓檢測實驗效果如圖11所示，圖中4位患者的定位結(jié)果存在兩個主要問題，第一，檢測到左右眼虹膜外圓但外圓的大小和位置相比于真實的虹膜外圓位置和大小都出現(xiàn)了一定偏差；第二，部分圖像只能檢測左眼或者右眼的虹膜外圓輪廓。本文的區(qū)域Hough圓檢測實驗效果如圖12所示，圖中本文方法對每張圖都能準(zhǔn)確找到兩只眼的虹膜外圓，保證虹膜外圓大小和位置符合其實際值，另外對虹膜外圓邊緣信息不足的情況，補全了擬合圓的輪廓信息，完成虹膜外圓的定位。

圖11 經(jīng)典Hough圓檢測定位結(jié)果

圖12 本文Hough圓檢測結(jié)果

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

本文實驗的數(shù)據(jù)集由兩部分組成，分為私有數(shù)據(jù)集和公有數(shù)據(jù)集。私有數(shù)據(jù)集主要用來完成診斷眼瞼退縮實驗，由眼科醫(yī)院提供的2500張圖像組成，其包括1822張已知存在眼瞼退縮的TAO患者圖與678張確認(rèn)為非TAO患者的圖像。圖片尺寸統(tǒng)一為1915×1080，圖片為jpg格式。公有數(shù)據(jù)集主要用來驗證本文所提出的方法同樣適用于專業(yè)設(shè)備拍攝的局部圖像，在定位精度和時間上有一定的提升，數(shù)據(jù)集從中科院CASIA虹膜庫中隨機抽取的100幅圖像組成，分辨率為640×480，圖片為png格式。

2.1 精度分析

為了進一步驗證所提出的方法的可行性，進行了兩組對比實驗，實驗基于中科院CASIA虹膜庫隨機抽取的2500幅圖像。第一組是通過對虹膜的定位精度、虹膜外圈的定位時間和最小二乘法與霍夫變換的虹膜定位算法相比較，后者采用形態(tài)學(xué)、自適應(yīng)閾值算法、Canny邊界提取、最小二乘法等技術(shù)，對虹膜外圓進行了定位。第二組是將本文方法與文獻[9]所提出方法進行對比，后者根據(jù)SIFT技術(shù)和SDM技術(shù)，結(jié)合文獻[9]中所提的虹膜定位技術(shù)，采用徑向?qū)ΨQ變換、微積分算子、SIFT特性以及SDM等方法，實現(xiàn)了對目標(biāo)的精確定位。并給出了式(13)中的位置精度的計算公式

(13)

式中：O表示實際虹膜面積，N表示測得的虹膜面積，C表示重復(fù)的虹膜面積。根據(jù)實際虹膜與測得的虹膜面積覆蓋率表示其定位精度。實驗結(jié)果見表2。

表2 CASIA虹膜庫圖像定位結(jié)果對比

表2在CASIA虹膜庫數(shù)據(jù)集的基礎(chǔ)上一共進行了5組實驗，根據(jù)歐洲眼病專家組(European group of Graves’ orbitopathy，EUGOGO)[19]進行精度和方差的對比，Daugman提出的基于微積分的虹膜定位算法、Wilde提的Hough變換與Canny相結(jié)合的定位算法、文獻[8]中的最小二乘法與Hough變換相結(jié)合的虹膜定位算法、文獻[10]提出的基于SDM和SIFT虹膜定位算法以及本文提出的基于區(qū)域Hough變換的虹膜定位方法。實驗數(shù)據(jù)集從CASIA數(shù)據(jù)庫中隨機抽取2500張，從表格可以看出，定位精度、方差和外圓定位耗時三方面進行對比，本文所提方法在精度、穩(wěn)定性以及時間上都優(yōu)于其它算法。其中Daugman和Wilde所提的定位方法由于在早期就已經(jīng)提出，所以其精度較其它算法都較低。文獻[8]中二值化過程比較粗糙，不夠精確，其利用灰度直方圖的波峰波谷取得閾值，很容易導(dǎo)致實驗結(jié)果出現(xiàn)一定偏差，但是在區(qū)域Hough變換法中就通過對目標(biāo)區(qū)域和搜索半徑的范圍進行確定，在基于CASIA虹膜庫數(shù)據(jù)集的實驗中表明其具有更高精度和更低耗時的優(yōu)勢。文獻[10]中方法主要思想是以徑向?qū)ΨQ和微積分算子為基礎(chǔ)確定虹膜中心和半徑，但是對圖像質(zhì)量要求較高，初步的計算結(jié)果將直接影響最終的虹膜定位結(jié)果。綜上所述，本文所提區(qū)域Hough變換法進行虹膜定位在較短的時間消耗情況下具有較高的精度。

2.2 眼瞼退縮診斷結(jié)果分析

為了驗證本文方法在實際應(yīng)用場景中能有效輔助醫(yī)生診斷患者是否存在眼瞼退縮，共進行了5組對比實驗，同時計算了實際誤診率。在實驗中對2500張圖像進行眼瞼退縮診斷，其中1822張確認(rèn)為TAO患者，678張確認(rèn)為非TAO患者圖像。最終實驗定位效果如圖13所示(圖片已經(jīng)獲得可公開授權(quán))。

圖13 虹膜外圓與眼瞼定位

眼瞼退縮診斷結(jié)果見表3。

表3 眼瞼退縮檢測實驗結(jié)果

如表3所示，進行了5組實驗，在2500張數(shù)據(jù)集上通過準(zhǔn)確率、誤診率和平均耗時3個指標(biāo)來評估算法。其中Daugman和Wilde這兩種方法相比其它算法存在準(zhǔn)確率低且誤診率高的問題，這是因為這兩種算法對圖像質(zhì)量要求較高，但是患者圖像面部的特征不清晰并且不規(guī)則，所以在檢測時存在較大誤差，耗時長。文獻[8]和文獻[10]準(zhǔn)確率與誤差值較為相近，耗時也相近，這是由于兩種算法均是在原始方法上進行部分特征改進，所以在準(zhǔn)確率和誤診率以及耗時上比Daugman和Wilde方法更好，但是由于患者圖像的特殊性，其準(zhǔn)確率等都存在一定上限，而本文的算法是專門針對固定的患者圖像提出的，所以在準(zhǔn)確率等指標(biāo)上比其它方法更具有優(yōu)勢與針對性。本文算法準(zhǔn)確率為93.5%，且誤診率低，達到了實際應(yīng)用的基本要求。綜上所述，本文所提出的區(qū)域Hough法在準(zhǔn)確率、誤診率以及耗時三方面均優(yōu)于其它幾種算法。

3 結(jié)束語

針對傳統(tǒng)的人工眼瞼退縮診斷準(zhǔn)確率低以及耗費時間長的缺點，本文提出了一種基于區(qū)域Hough變換新的虹膜定位方法。與傳統(tǒng)的Hough方法相比，該算法采用特征點、Bezier曲線和三次B樣條擬合出眼瞼輪廓和目標(biāo)檢測區(qū)域，從而縮小了搜索范圍，并減少大量的時間開銷，通過改進的Canny邊緣檢測方法和改進的Hough圓檢測方法提升了鞏膜外圓定位的精度，再結(jié)合EUGOGO眼瞼退縮標(biāo)準(zhǔn)進行診斷。實驗結(jié)果表明本文改進的方法對眼瞼退縮的診斷有較大的提升，滿足實際計算機輔助診斷的需求，又能在一定程度上有效提高在利用便攜式設(shè)備局部獲取圖像下的精度，縮短定位時間[22]。