生理結構先驗引導下的虹膜精確分割算法

年炳坤，丁建睿，史夢蝶，黃子晨

(哈爾濱工業(yè)大學(威海)計算機科學與技術學院,山東 威海 264209)

隨著生物識別技術的發(fā)展，人們逐漸將某些具有區(qū)分能力的生物特征運用到身份識別中。其中虹膜是所有生物識別系統(tǒng)中應用最為廣泛的生物特征。虹膜是位于黑色瞳孔和白色鞏膜之間的圓環(huán)狀部分[1]，擁有高于人臉和指紋識別幾個數(shù)量級的精確性，被廣泛應用在身份識別的各個領域，如：美國-墨西哥邊境控制，2004年雅典奧運會身份鑒定，阿聯(lián)酋驅逐人員的管理等[2]。

20世紀90年代Daugman[3]在假設瞳孔與虹膜為同心圓的前提下提出了虹膜分割算法，取得了一定的效果。而后針對各種不同的情況，許多學者針對某一種特殊情況(瞳孔的中心偏離[4-5]、遮擋[6-7]、鏡面反射[6])，提出了特定的分割方法。為了進一步提高分割準確度，Daugman[8]又提出了基于傅里葉變換并假設瞳孔與虹膜為非同心的Snake主動輪廓模型，隨后Sutra等[9]徹底拋棄掉了同心以及輪廓外形為圓的假設，提出基于Viterbi的主動輪廓探測算法，可以較為精確地探測虹膜與瞳孔，但是算法對非虹膜區(qū)域(睫毛、眼瞼)的干擾較為敏感。Llano等[10]利用主成分分析(principal component analysis,PCA)進行圖像融合，隨后提出了基于拉普拉斯金字塔的融合方式(modified Laplacian pyramid,MLP)[11]，進一步提高了虹膜系統(tǒng)的性能。但Daugman的算法將瞳孔與虹膜看作同心圓，這就會導致?lián)p失大量的有效信息；而Snake主動輪廓模型和Viterbi算法雖然拋棄了非圓、非同心的假設，但是算法本身受睫毛、眼瞼等噪聲影響較大；而PCA融合算法和MLP融合算法雖然有效地提高了系統(tǒng)的識別性能，但其本質上是對已有分割算法進行融合，讓分割結果更容易被特征提取和編碼，因此其性能受分割算法的制約。

同時國內學者劉琳等[12]對于主動輪廓算法進一步研究，提出了結合KL散度和RSF模型的主動輪廓圖像分割方法[12]，提高了圖像分割的精度，可以一定程度上提高識別精度。

本文基于CLBP(completed local binary patterns)思想[13]，提出一種形狀敏感的成分分析算法，可以有效檢測虹膜區(qū)域，同時減少其他成分(如：睫毛與眼瞼區(qū)域)的干擾，從而避免分割過程中虹膜的信息丟失，提高整體的分割精度，進而提高虹膜識別系統(tǒng)的性能。

1 本文方法

本文提出的算法由粗分割與精細分割兩部分組成，整個流程見圖1，在分割基礎上利用Daugman[13]提出的方法進行匹配。

圖1 本文方法流程圖Fig.1 Flow chart of proposed method

首先利用瞳孔灰度值低于其他區(qū)域的生理特性，對圖片進行二值化。直接使用二值圖像進行后續(xù)處理存在以下問題：1)二值化圖像存在噪聲點，并不是所有像素點都屬于瞳孔，見圖2(a)；2)由于拍攝條件和瞳孔本身具有反光的生理特性，屬于瞳孔的像素點在局部區(qū)域會出現(xiàn)異常，致使二值化圖像中瞳孔可能出現(xiàn)斷裂現(xiàn)象，不利于后續(xù)處理，見圖2(b)。

由于噪聲點大多孤立，或者小部分局部聚集；同時瞳孔局部灰度值異常面積不會太大。為解決上述噪聲問題，本文對二值圖像進一步進行形態(tài)學操作，去除噪聲的同時，可以在一定程度上對瞳孔斷裂區(qū)域進行連接。形態(tài)學腐蝕、膨脹、閉操作分別見式(1)～(3)[14]：

A?B={z|(B)z∩Ac=?}

(1)

(2)

A·B=(A⊕B)?B

(3)

式中A為輸入的二值化圖像，B為結構元。由于本文所處理的形狀多為圓形，所以選擇3*3的結構元。處理結果見圖2(c)。

圖2 二值化圖像示例Fig.2 Examples of binary images

1.1 瞳孔粗分割

通過上述操作，可以認定圖像中的像素值為0的點絕大多數(shù)屬于瞳孔。先假定瞳孔為圓，首先利用式(4)計算圖像中黑色點集的幾何中心，如下所示：

(4)

式中：I(xk,yk)為坐標點(xk,yk)的像素灰度值，N為圖像中灰度值為0的像素個數(shù)，(x,y)為黑色點集的幾何中心，為瞳孔內一點。

對具有最大面積的連通區(qū)域補全后，求得其最小外接矩形，該矩形的內接圓邊界即為瞳孔的粗邊界，該邊界用于補全晶狀體、眼鏡等反光介質導致的局部灰度值異常所帶來的瞳孔斷裂問題。

1.2 瞳孔精確分割

利用連通域分析，找到最大連通域，并標記出其邊界，這是瞳孔的精確邊界，記作A(I)。進一步采取如下規(guī)則將瞳孔粗邊界R(I)和瞳孔細邊界A(I)進行融合，以解決由于瞳孔自身局部灰度值異常可能導致的問題。

1)如果粗邊界完全包含在精確邊界之內，則采用精確邊界；

2)如果粗邊界部分包含在精確邊界之內，則最終邊界為粗邊界和精確邊界的并集。

其數(shù)學描述如式(5)所示：

F(I)=

(5)

式中：I代表輸入的原始虹膜圖像，R(I)為瞳孔粗邊界，A(I)為瞳孔精邊界，F(xiàn)(I)代表融合后的掩碼圖，融合后的邊界即為瞳孔的最終外邊界,結果見圖3(a)。

1.3 虹膜外圓粗分割

在假設虹膜外邊界為圓形的分割算法中，大多采用Daugman的球形積分算子，在理想情況下，即沒有眼瞼、睫毛遮擋，虹膜完整呈現(xiàn)在眼球當中時，該方法的準確度較高，如圖3(b)-(1)所示。

但是大多數(shù)情況下，虹膜圖像大多存在遮擋問題，此時，球形積分算子達不到預期效果，如圖3(b)-(2)所示。

經(jīng)過分析，影響精度最大的因素是上眼瞼的存在。另外，球形積分算子本質上是探測梯度變化，當檢測到眼皮和眼球分界線時，邊界線灰度值變化遠比虹膜和鞏膜分界線的要大，所以影響較為顯著，從而導致算子檢測錯誤，如圖3(b)-(3)所示。

圖3 過程化分割結果示例Fig.3 Examples of procedural segmentation

針對以上問題，本文對其進行改進。根據(jù)眼球的生理特性可知：眼瞼、睫毛的遮擋只會發(fā)生在瞳孔的上下區(qū)域，而左右區(qū)域一般不會被遮擋，只存在鞏膜和虹膜的分界，而這有利于我們的探測，利用此先驗知識為Daugman的球形積分算子增加約束，得到新的弧形積分算子，如式(6)所示：

(6)

式中：r為所求虹膜的半徑，(x0,y0)為瞳孔幾何中心的坐標，I(xk,yk)為圖像中坐標為(xk,yk)的像素點的灰度值。

利用式(6)可以得到粗分割的虹膜外邊界,同時可以有效解決不同程度的遮擋問題，結果見圖3(c)。

1.4 虹膜外圓精分割

根據(jù)CLBP中的思想[13]，一幅圖像可以初步處理成dp，將原圖像劃分成多個3*3的區(qū)域集合,dp的定義如式(7)所示：

dp=gp-gc

(7)

其中gc為3*3區(qū)域的中心點，gp為其相鄰的點，而dp可以進一步被分解為兩部分，sp和mp，定義如式(8)所示：

(8)

文獻[13]中同時證明了mp中含有圖像絕大部分信息，而sp中含有較少的信息。本文認為：mp中的信息為結構信息，而sp中的信息為關系信息。為充分利用兩種信息，提出以下規(guī)則對其進行融合，起到類似濾波的效果，但比傳統(tǒng)意義上的濾波更加具有導向性，能夠在探測指定形狀的同時，減小其他噪聲成分的干擾，同時可以提高精度。

式(1)中,Ai為不同濃度樣品與DPPH溶液的混合液的吸光度；Aj為不同濃度樣品與無水乙醇的混合液的吸光度；A0為DPPH溶液與無水乙醇的混合液的吸光度。

1)如果需要檢測的形狀為橫向長條形，則采用上下4連通域，如果為豎向長方形，則采用左右2連通域；

2)如果需要檢測的形狀為圓形，則采用8連通域；

3)如果想擴大搜索區(qū)域，則可以擴大連通域。

1.4.1 睫毛分割

睫毛不同于眼瞼，可以近似看成一個各向異性的長方體，由此可以構造一個對四周敏感的融合式。本文構造了一個8連通域區(qū)域，中心區(qū)域的坐標點為(0,0)，集合G={(1,0),(0,1),(-1,0),(0,-1),(1,1),(1,-1),(-1,1)，(-1,-1)}。為了減少計算量，本文定義先前虹膜外圓粗分割時得到區(qū)域為ROI區(qū)域，在此區(qū)域內進行運算。對ROI區(qū)域中所有c∈C和ε∈G使用式(9)、(10)

(9)

(10)

得到濾波后的結果作為睫毛的掩碼圖層，見圖3(d)。

圖3(d)-(1)證明本文提出的算法可以檢測出較為稀疏的睫毛，而圖3(d)-(2)和圖3(d)-(3)則證明算法對于密集的睫毛同樣具有很好的檢測效果。在ROI區(qū)域內，無論是上眼瞼的睫毛，還是下眼瞼的睫毛，都被完整的檢測。同時，結果也符合先前的構造規(guī)則。

采用較高運算量的長條狀8連通域進行檢測，檢測結果基本都是長條狀物體，而上下眼瞼邊緣處陰影由于拍攝原因，可以近似認為是黑色的不連續(xù)陰影，算法也可以準確檢測，這進一步證明了式(9)和式(10)對于給定形狀的敏感性和穩(wěn)定性，能夠較大程度上去掉其他噪聲點干擾。

通過虹膜外圓粗分割，可以得到虹膜ROI區(qū)域，但是這個ROI區(qū)域包含了大量不需要的信息，例如：眼球以外的眼瞼區(qū)域。探測出整個眼皮區(qū)域的難度很大，可以探測眼球和眼瞼的邊界區(qū)域，然后和虹膜和鞏膜邊界一起構成一個閉合區(qū)域，這個區(qū)域就是去除眼瞼后的掩碼圖層。

從灰度值來看，這個邊界區(qū)域不一定是連續(xù)的，可能存在眼淚等反光物質，在邊界線的某一點處，存在灰度值異常，從而斷裂開來；同時眼瞼的探測還需解決如何在睫毛的干擾下，有效探測出準確邊緣。

根據(jù)眼瞼的生理學特點可知，人的眼瞼邊界可以近似看成一個橫向長條狀，同時為了和睫毛區(qū)分，提高探測精度，本文選擇上下4連通域，只探測橫向的相關區(qū)域，來屏蔽豎向睫毛的干擾。構造一個4連通域，連通域中心點的坐標為(0,0)，G={(0,1),(0,2),(0,-1),(0,-2)}。對所有c∈C在ROI區(qū)域中的像素點使用式(9)和(10)。

通過連通域分析后，結果見圖3(e)。可以看出，新構造的算式，能夠準確地探測出眼瞼的邊界區(qū)域，如圖3(e)-(1)結果所示，邊界上存在局部反光導致的灰度值異常區(qū)域可以被準確避開，最終形成一個閉合區(qū)域；圖3(e)-(2)中，在如此多睫毛的干擾下，仍能準確地檢測出邊界；圖3(e)-(3)中，算法可以有效避開睫毛密集的地方，雖然損失了睫毛間隔之間小部分屬于虹膜的區(qū)域，但成功避免了大量睫毛帶來的誤差。

最后將多個掩碼圖進行合并，可以得到最終的虹膜分割結果，對其進行特征編碼并識別。

2 實驗設計

本實驗采用開放數(shù)據(jù)集CASIA-IrisV3-Interval和CASIA-IrisV1[15-16]。虹膜區(qū)域分割后，采用OSIRIS 4.1版[17]的匹配模塊進行識別，它是在Biecure項目框架下的一個開源虹膜識別系統(tǒng)。

兩種數(shù)據(jù)集中照片格式均為BMP，同時數(shù)據(jù)集中的圖像質量并不相同。本文所提出的算法采用的是傳統(tǒng)分割技術，不涉及到機器學習內容，不存在過擬合等現(xiàn)象。

CASIA-IrisV3-Interval數(shù)據(jù)集中包含了來自247個不同的個體，共計2 639張圖片；CASIA-IrisV1數(shù)據(jù)集中包含了來自104個個體，共計709張圖片。數(shù)據(jù)集中的圖像示例見圖4。

圖4 虹膜數(shù)據(jù)集示例Fig.4 Examples of iris datasets

考慮到虹膜識別本質上是一個分類問題，本文采取衡量分類性能的指標，PGAR(正確接受率)在PFAR(錯誤接受率)為10-4時的值,PEER(等誤差率)以及PFRR(拒識率)的值，如公式(11)～(13)所示。對于一個二分類模型的單個樣本一共有4種預測結果：RTP(true positive)、RFP(false positive)、RFN(false negative)、RTN(true negative)。

(11)

(12)

(13)

PEER值為PFRR與PFAR相等時的值。同時采用ROC曲線下的面積(AAUC)來衡量算法的分類性能。

3 實驗結果

圖5選取了對比實驗中涉及到分割的算法：CHT[18](以同心圓為假設并部分改進的算法)、Viterbi以及本文所提出算法的內外虹膜分割結果樣例。然后將分割結果基于Daugman提出的匹配方式輸入到OSIRIS Version 4.1[15]系統(tǒng)中，得到任意兩幅圖片的相似程度，然后設定閾值，計算出本閾值下的PGAR、PEER值，并求出對應的AAUC值，實驗結果見表1。

從表1中可以看出，CHT分割算法的性能最低，原因是這個算法基于圓假設，包含了一部分無用區(qū)域或者少包含了一部分有用區(qū)域，如圖5(b)-(2)，去掉了眼瞼遮擋區(qū)域的同時，也去掉了很大一部分有用的虹膜區(qū)域，如圖5(c)-(2)中所示。Viterbi算法很大程度上提高了準確性,對于數(shù)據(jù)集CASIA-IrisV3-Interval，Viterbi算法將PGAR在PFAR=10-4的值從95.32%提升到了97.04%。這表明Viterbi算法具有更好的分割性能，提升的原因在于Viterbi算法采用了主動探測輪廓，能夠較為精確地定位出瞳孔、虹膜的外邊界區(qū)域，而不是寬泛的規(guī)定和去掉某些區(qū)域，從而有效區(qū)域損失和無用區(qū)域的包含遠比CHT算法要小的多，但Viterbi算法也存在缺點，如圖5(d)-(3)中顯示，算法對于眼瞼邊界部分檢測不夠準確，當邊界不太明顯且存在局部灰度值異常的時候，存在較大范圍誤差，同時如圖5(c)-(3)圖像中所示，由于假設分割區(qū)域灰度值成正態(tài)分布，在給定2.35個σ(標準差)外都認為是雜質，并加以去除，所以導致丟失了部分虹膜信息，影響了分類性能。

表1 在數(shù)據(jù)集CASIA-IrisV3-Interval上的結果Tab.1 Results on database CASIA-IrisV3-Interval

圖5 分割結果示例Fig.5 Examples of segmentation results

基于PCA主成分分析和MLP的融合算法，融合了CHT和Viterbi算法的分割結果，一定程度上彌補了缺陷，但是受CHT影響較大，并不能很好地提高算法性能，甚至有可能會下降。PCA主成分分析的原理是利用數(shù)學上的PCA降維手段，將兩個對象的數(shù)據(jù)，轉換到K個特征向量構成的新空間中，而MLP融合算法采用拉普拉斯金字塔進行融合，本質上合成新的特征向量，這就導致融合結果強烈依賴于融合對象，也就是分割算法的性能，如果兩個分割算法能夠互相彌補自身不足，那么融合表現(xiàn)就好，反之如果缺陷都較大，那么融合算法一定程度上也會放大缺陷。如在表2中，PCA融合后的PEER值反而上升。而本文提出的算法對于同種結構的圖片，有較強的魯棒性，可以根據(jù)不同目的構造不同的算子，將睫毛、眼瞼等虹膜特有噪聲檢測出來，在精確分割的同時，減少虹膜信息的損失。

表2 在數(shù)據(jù)集CASIA-IrisV1上的結果Tab.2 Results on database CASIA-IrisV1

從表2中可以看出，表現(xiàn)最好的MLP融合算法和本文的算法結果相差不大，PGAR、PEER值基本一致，而AAUC值也只是略高一點。通過分析CASIA-IrisV1數(shù)據(jù)集的特點，得出原因在于數(shù)據(jù)集CASIA-IrisV1拍攝條件較好，數(shù)據(jù)集中的虹膜圖像，無論是眼瞼、睫毛等遮擋均比CASIA-IrisV3-Interval中要少，圖片整體質量較高。所以Viterbi和CHT算法的缺陷(對于睫毛、眼瞼的探測精度不足)被彌補。PCA和MLP融合算法融合了Viterbi和CHT的優(yōu)點，一定程度上特征被增強，使得融合后的圖像更利于分辨。導致本文提出的算法在新性能上并沒有顯著的提高。同時這也說明了本文提出的算法，在分割性能上具有優(yōu)越性，在沒有融合的情況下達到了相似的性能，這也意味著如果將PCA與MLP融合算法的融合對象CHT換成本文算法，那么精度會進一步提高。

算法在保證精度的同時，時間復雜度較低。算法中的大部分運算可以表示為簡單的加減運算。這大大提高了算法的性能。表3給出了本文算法與CHT和Viterbi算法的運行時間比較。這些算法在帶有Intel(R)Core(TM)i7-7700hq、16 GB內存、GTX1060顯卡的計算機上運行。在Windows 10上使用C++編程語言和OpenCV4.1.0。

表3 算法平均運行時間Tab.3 Average running time of different algorithms s

4 結 論

本文基于CLBP思想，提出了一種成分分析算法，并且在開源數(shù)據(jù)集CASIA-IrisV3-Interval與CASIA-IrisV1上進行了實驗，與CHT、Viterbi、PCA和MLP融合算法進行了比較，從實驗數(shù)據(jù)上看，本文所提出的算法，在兩個數(shù)據(jù)集上都具有較高的AAUC值和較低的PEER值,說明本文算法具有更優(yōu)秀的分割性能。

同時，在實驗過程中發(fā)現(xiàn)本文所提算法對于低質量的圖片，如有模糊、散焦等問題的圖片分割效果略差，并且算法較強依賴于進行二值化的閾值，這導致算法對于不同攝像設備所拍攝的不同結構的虹膜圖像效果不顯著。如何根據(jù)圖像質量自動化確定合理的閾值，是今后算法改進的主要方向，可以借助深度學習方法，自主決定閾值，來提高算法的魯棒性。