鄒耀斌 雷幫軍 臧兆祥 王俊英 胡澤海 董方敏

圖像閾值分割由于具有簡(jiǎn)潔性、直觀性和有效性等特點(diǎn)而被持續(xù)關(guān)注[1?5].閾值分割的基本原理是將圖像中每個(gè)像素的灰度值與選擇的閾值進(jìn)行比較,然后確定像素是屬于前景還是背景.選擇閾值由此成為閾值分割的關(guān)鍵[6],而讓所選閾值盡可能地接近最優(yōu)閾值則是閾值分割的核心目標(biāo)[7].

在閾值的自動(dòng)選擇方面,一類(lèi)思想是以概率統(tǒng)計(jì)學(xué)和/或信息論為基礎(chǔ),利用灰度直方圖提供的原始信息和/或派生信息計(jì)算閾值.根據(jù)灰度直方圖的維度以及構(gòu)建灰度直方圖過(guò)程中所涉及到的圖像空間,可以將灰度直方圖分為全局一維直方圖,局部一維直方圖,全局聯(lián)合直方圖三大類(lèi).

全局一維灰度直方圖通過(guò)對(duì)整幅圖像的灰度值分布進(jìn)行統(tǒng)計(jì)而建立.利用全局一維灰度直方圖的典型方法包括Otsu方法[8]、最小誤差法[9]、最大香農(nóng)熵法[10]、最小交叉熵法[11]以及基于非香農(nóng)熵的方法[12?15].這些方法雖然適合分割灰度直方圖為雙峰模式且峰谷特征明顯的圖像,但是容易誤分割灰度直方圖為單峰模式的圖像.一些新近提出的方法[16?21]嘗試在全局一維灰度直方圖范圍內(nèi)對(duì)相關(guān)方法進(jìn)行改進(jìn),但是它們面臨著一個(gè)共同的難題:對(duì)于全局一維灰度直方圖相同但是圖像內(nèi)容不同的兩幅圖像,它們無(wú)法計(jì)算出不同的閾值.

局部一維灰度直方圖通過(guò)對(duì)圖像特定區(qū)域的灰度值分布進(jìn)行統(tǒng)計(jì)而建立.典型的特定區(qū)域是前景和背景之間的過(guò)渡區(qū)[22?25],其一維灰度直方圖具有近似對(duì)稱分布的傾向[7].這種傾向性雖然降低了局部一維灰度直方圖退化為單峰模式的風(fēng)險(xiǎn),但同時(shí)也弱化了直方圖的峰谷特征.另外,為了獲取過(guò)渡區(qū)的局部一維灰度直方圖信息,需要先提取出合理的過(guò)渡區(qū),這會(huì)引發(fā)如何自動(dòng)評(píng)估未知特征參數(shù)的新難題.

全局聯(lián)合灰度直方圖通過(guò)對(duì)原始圖像和其派生圖像的灰度值進(jìn)行聯(lián)合分布統(tǒng)計(jì)而建立.派生圖像可以是原始圖像的均值濾波或者中值濾波結(jié)果[26?33],也可以是原始圖像的閾值化分割結(jié)果[34].限于計(jì)算復(fù)雜度,目前全局聯(lián)合灰度直方圖的維度主要是二維[26?30]和三維[31?33].在利用全局聯(lián)合灰度直方圖進(jìn)行閾值分割方面,由盧振泰等[34]提出的結(jié)合k均值和互信息量(Combiningk-means and mutual information,CKMI)的方法頗具特色.然而,CKMI方法傾向于選擇這樣的二值圖像:值為0的像素個(gè)數(shù)和值為1的像素個(gè)數(shù)的差盡量小(具體分析見(jiàn)第1節(jié)).這種傾向性使得CKMI方法很容易誤分割前景和背景大小不平衡的圖像.

上述方法大多限于分割具有近似雙峰模式或者近似單峰模式直方圖的圖像,并且相對(duì)而言,它們一般更適合于處理前景和背景的灰度分布可以用正態(tài)分布近似逼近的情形.然而,灰度直方圖所呈現(xiàn)的模式是復(fù)雜多變的.除了雙峰和單峰模式外,還有無(wú)峰和多峰模式.另外,構(gòu)成每種直方圖模式的基本分布又有多種可能,除了正態(tài)分布外,還有常見(jiàn)的極值分布、瑞利分布、貝塔分布和均勻分布等.當(dāng)前景或背景的灰度分布呈現(xiàn)為極值分布、瑞利分布、貝塔分布、均勻分布或者它們的混合分布時(shí),如何自動(dòng)選擇出合理的分割閾值是個(gè)新難題.

為了在統(tǒng)一框架內(nèi)處理上述不同灰度分布情形下的閾值選擇問(wèn)題,提出了一種歸一化互信息量最大化導(dǎo)向(Maximizing normalized mutual information,MNMI)的自動(dòng)閾值選擇方法.MNMI方法比較不同閾值分割圖像所對(duì)應(yīng)的輪廓圖像和一幅多尺度梯度乘圖像的歸一化互信息量,并以最大的歸一化互信息量所對(duì)應(yīng)的閾值作為最終分割閾值.實(shí)驗(yàn)結(jié)果表明,當(dāng)前景或背景的灰度分布可以用正態(tài)分布、極值分布、瑞利分布、貝塔分布、均勻分布或者它們的混合分布逼近時(shí),MNMI方法計(jì)算出的閾值比4個(gè)新近開(kāi)發(fā)的方法所得閾值在總體上更優(yōu).

剩余內(nèi)容安排如下:第1節(jié)先簡(jiǎn)介CKMI方法,然后剖析該方法的主要缺陷;第2節(jié)提出MNMI方法,并從理論上闡述了該方法所具有的優(yōu)勢(shì);第3節(jié)討論MNMI方法和5個(gè)被比較方法的實(shí)驗(yàn)結(jié)果;第4節(jié)給出結(jié)論.

1 CKMI方法及其缺陷

1.1 計(jì)算閾值的目標(biāo)函數(shù)

CKMI方法主要包括兩個(gè)步驟:第一步用k均值算法計(jì)算出一個(gè)初始閾值tinitial,第二步是在[tinitial?δ,tinitial+δ]范圍內(nèi),以最大化原始圖像和二值圖像之間的互信息量為目標(biāo),搜索一個(gè)更優(yōu)的最終閾值tfinal.盧振泰等[34]建議δ=10.

記原始圖像為X,使用一個(gè)閾值t∈[tinitial?δ,tinitial+δ]閾值化X,并記閾值化所得二值圖像為Yt.CKMI方法的第二步使用如下目標(biāo)函數(shù)計(jì)算最終閾值tfinal:

其中,I(X;Yt)表示計(jì)算X和Yt之間的互信息量.對(duì)于兩幅具有相同大小的圖像X和Yt,可以用式(2)計(jì)算它們之間的互信息量[34]:

其中,H(X)和H(Yt)分別是X和Yt的香農(nóng)熵,而H(X,Yt)是X和Yt的聯(lián)合香農(nóng)熵.

1.2 目標(biāo)函數(shù)的化簡(jiǎn)

下面將證明CKMI方法的目標(biāo)函數(shù),即式(1),可以被化簡(jiǎn)為:

證明.由于二值圖像Yt的像素值為0或1,則圖像X和Yt的聯(lián)合灰度分布PX,Yt(x,y)與圖像X的灰度分布PX(x)有如下特殊關(guān)系:1)當(dāng)x

1.3 計(jì)算閾值時(shí)的缺陷

在CKMI方法中,第一步采用了k均值算法.如果k均值算法估計(jì)的初始閾值tinitial嚴(yán)重偏離最優(yōu)閾值,以至于一個(gè)小的δ值無(wú)法讓區(qū)間[tinitial?δ,tinitial+δ]包含最優(yōu)閾值.此時(shí)CKMI方法的第二步將無(wú)法在區(qū)間[tinitial?δ,tinitial+δ]內(nèi)搜索到最優(yōu)閾值.

退一步講,假設(shè)擴(kuò)大δ取值使得區(qū)間[tinitial?δ,tinitial+δ]包含最優(yōu)閾值,那么CKMI方法將依賴式(1)計(jì)算最終閾值tfinal.由式(1)到式(3)的化簡(jiǎn)過(guò)程,可知這等效于使用式(3)計(jì)算最終閾值tfinal.根據(jù)香農(nóng)熵的定義,式(3)中的H(Yt)可以展開(kāi)為:

其中,p表示Yt中像素值為0的概率,而1?p則是Yt中像素值為1的概率.H(p)和p的關(guān)系是:1)0

情形1.如果區(qū)間[tinitial?δ,tinitial+δ]內(nèi)的任意閾值t閾值化原始圖像X,所得二值圖像Yt對(duì)應(yīng)的p值都小于0.5,則由于0

情形2.如果區(qū)間[tinitial?δ,tinitial+δ]內(nèi)的任意閾值t閾值化原始圖像X,所得二值圖像Yt對(duì)應(yīng)的p值都大于0.5,則由于0.5

情形3.如果區(qū)間[tinitial?δ,tinitial+δ]內(nèi)的任意閾值t閾值化原始圖像X,所得二值圖像Yt對(duì)應(yīng)的p值部分小于0.5,部分大于0.5,則式(3)將使得最終分割閾值tfinal具有這樣的特性:在tfinal閾值化原始圖像所得二值圖像Ytfinal中,值為0的像素個(gè)數(shù)和值為1的像素個(gè)數(shù)的差盡量小.

對(duì)于情形1,如果最優(yōu)閾值小于tinitial,那么tinitial+δ會(huì)誤分割更多像素.對(duì)于情形2,如果最優(yōu)閾值大于tinitial,那么tinitial?δ也會(huì)誤分割更多像素.退一步講,假設(shè)情形1(或情形2)所得閾值tinitial+δ(或tinitial?δ)比tinitial更優(yōu)或者就是最優(yōu)閾值,則會(huì)引發(fā)一個(gè)新問(wèn)題:對(duì)不同的圖像如何自動(dòng)計(jì)算出不同的δ.沒(méi)有證據(jù)表明這個(gè)問(wèn)題會(huì)比自動(dòng)計(jì)算閾值本身更容易.盧振泰等[34]建議δ=10,除了減少計(jì)算代價(jià)的解釋外,沒(méi)有關(guān)于δ=10的其他深入分析.此外,在上述3種情形下,CKMI方法將不管前景像素和背景像素的實(shí)際比例,只是去力求:在最終分割結(jié)果圖像中,值為0和值為1的像素個(gè)數(shù)盡可能相等.這種傾向性使得CKMI方法只適合處理前景和背景大小較為均衡的圖像.一旦前景和背景的大小失衡,CKMI方法將產(chǎn)生明顯的誤分割,這點(diǎn)也被后面的實(shí)驗(yàn)所證實(shí).

2 提出的MNMI方法

2.1 計(jì)算閾值的目標(biāo)函數(shù)

MNMI方法采用如下目標(biāo)函數(shù)計(jì)算最終閾值tfinal:

式(6)拓展了t的取值范圍,由CKMI方法中的[tinitial?δ,tinitial+δ]拓展為[tmin,tmax],其中,tmin和tmax分別表示原始圖像的最小灰度值和最大灰度值.這種拓展使得MNMI方法:1)不需要k均值算法去估計(jì)初始閾值tinitial;2)不需要考慮額外的參數(shù)δ;3)避免了CKMI方法中區(qū)間[tinitial?δ,tinitial+δ]不包含最優(yōu)閾值的可能性.

式(6)中的MX表示對(duì)原始圖像X進(jìn)行多尺度梯度乘變換后所得圖像(圖1(b)和(e)分別給出了X和MX的直觀示意).MX在像素位置(r,w)處的值可以按如下式(7)計(jì)算[7]:

圖1 MNMI方法所涉及的關(guān)鍵概念和步驟的圖示((a)原始圖像X 的灰度直方圖;(b)原始圖像X;(c)使用t閾值化原始圖像X所得的二值圖像Yt;(d)對(duì)原始圖像X進(jìn)行多尺度梯度乘變換后所得圖像MX的灰度直方圖;(e)圖像MX;(f)從二值圖像Yt中提取輪廓后所得輪廓圖像CYt.注意,為了能更清楚地顯示灰度區(qū)間[0,255]內(nèi)灰度值出現(xiàn)的頻率,(d)中灰度直方圖在頻率為80處進(jìn)行了截?cái)?)Fig.1 Graphic illustration of crucial concepts and steps in MNMI method((a)Gray level histogram of original image X;(b)original image X;(c)binary image Ytobtained by thresholding original image X with a gray level t;(d)gray level histogram of image MXproduced by applying a multiscale gradient multiplication transformation to original image X;(e)image MX;(f)contour image CYtextracted from binary image Yt.The gray level histogram in(d)is truncated at the frequency 80 for more clearly showing the frequency of gray level in the range[0,255].)

式(6)中的CYt表示二值圖像Yt對(duì)應(yīng)的輪廓圖像(圖1(c)和(f)分別給出了Yt和CYt的直觀示意).CYt可如下獲得:首先讓CYt等于Yt,然后對(duì)Yt中每個(gè)值為1的像素進(jìn)行判斷,如果該像素的四鄰域均為1,則設(shè)置該像素在圖像CYt中對(duì)應(yīng)位置的像素值為0.

式(6)中的Inormal(MX;CYt)表示計(jì)算圖像MX和圖像CYt的歸一化互信息量.由于MX和CYt具有相同大小,可以應(yīng)用式(8)計(jì)算兩者之間的歸一化互信息量[35]:

式中,H(MX)和H(CYt)分別是MX和CYt的香農(nóng)熵,而H(MX,CYt)是MX和CYt的聯(lián)合香農(nóng)熵.

2.2 MNMI方法所具優(yōu)勢(shì)的理論分析

其次,MNMI方法將原始圖像X變換成圖像MX,實(shí)際上是實(shí)施了直方圖規(guī)范化:將圖像X的直方圖規(guī)范化為圖像MX的單峰長(zhǎng)拖尾直方圖.如圖1(d)所示,圖像MX的灰度直方圖的峰起始于灰度值0處,尾結(jié)束于灰度值255處.這種直方圖規(guī)范化有利于MNMI方法以統(tǒng)一的框架處理具有不同直方圖模式的原始圖像X.圖像MX的單峰長(zhǎng)拖尾直方圖的形成源于邊緣信號(hào)和噪聲對(duì)FDoG濾波器不同的響應(yīng)特點(diǎn).隨著濾波尺度的增加,邊緣信號(hào)的響應(yīng)值在總體上得到更好的保持,而噪聲的響應(yīng)值則更快地衰減.當(dāng)不同濾波尺度下的響應(yīng)值相乘時(shí),絕大部分噪聲的總響應(yīng)值被往0處抑制,強(qiáng)弱不同的邊緣信號(hào)的總響應(yīng)值有大有小,分布在[0,255]區(qū)間.這些綜合因素使得圖像MX的灰度直方圖具有單峰長(zhǎng)拖尾特征.

2.3 算法描述

為了更清晰地說(shuō)明如何利用式(6)～(8)計(jì)算最終閾值tfinal,這里進(jìn)一步給出了實(shí)現(xiàn)MNMI方法的8個(gè)步驟.

算法1.算法名稱:MNMI

輸入.原始灰度圖像X

輸出.閾值tfinal及閾值分割結(jié)果圖像Ytfinal

步驟1.對(duì)輸入的原始圖像X,分別計(jì)算出n個(gè)濾波尺度 (即σ1,σ2,···,σn) 下的 FDoG 濾波結(jié)果圖像,并將這n個(gè)濾波結(jié)果圖像相乘產(chǎn)生出規(guī)范圖像MX(注:FDoG濾波系數(shù)可以預(yù)先生成并按規(guī)律存儲(chǔ)在一個(gè)MAT文件中,做FDoG濾波時(shí)直接訪問(wèn)預(yù)先載入的MAT文件即可快速獲得相應(yīng)的濾波系數(shù)).

步驟2.統(tǒng)計(jì)規(guī)范圖像MX的離散灰度分布,并利用該分布計(jì)算出MX的香農(nóng)熵H(MX)(注:香農(nóng)熵H(MX)只需計(jì)算一次).

步驟 3.遍歷原始圖像X的灰度區(qū)間[tmin,tmax],對(duì)于該區(qū)間內(nèi)的每個(gè)灰度值t,如果它在原始圖像X中出現(xiàn)過(guò),則按步驟4到步驟7順序執(zhí)行.執(zhí)行時(shí),用一個(gè)變量Ifinal記錄當(dāng)前最大的歸一化互信息量,初始時(shí)令I(lǐng)final=0;用一個(gè)變量tfinal記錄當(dāng)前最大的歸一化互信息量對(duì)應(yīng)的灰度值,初始時(shí)令tfinal=tmin.

步驟4.對(duì)于灰度值t,使用如下規(guī)則產(chǎn)生一幅對(duì)應(yīng)的二值圖像Yt:對(duì)原始圖像X的每個(gè)像素進(jìn)行判斷,如果X在某個(gè)像素位置(r,w)處的灰度值大于或等于t,則將Yt(r,w)設(shè)置為1,否則將Yt(r,w)設(shè)置為0.

步驟 8.輸出最終的分割閾值tfinal,同時(shí)按步驟4的規(guī)則產(chǎn)生一幅相應(yīng)的二值圖像Ytfinal并輸出.

3 實(shí)驗(yàn)結(jié)果與討論

3.1 實(shí)驗(yàn)環(huán)境、定量評(píng)估指標(biāo)、參與比較的方法

實(shí)驗(yàn)環(huán)境配置如下: Intel Core i3-2350M 2.3GHz CPU,4GB DDR2內(nèi)存,Windows 7操作系統(tǒng),Matlab 2009a編程平臺(tái).測(cè)試圖像集包含具有不同灰度分布模式的9幅合成圖像和59幅真實(shí)世界圖像,每幅測(cè)試圖像的分割參考圖像借助于Adobe Photoshop CS軟件手工產(chǎn)生(請(qǐng)?jiān)L問(wèn)https://pan.baidu.com/s/1bo9QNH9下載測(cè)試圖像、分割參考圖像和分割參考圖像的制作說(shuō)明).

采用誤分類(lèi)率 (Misclassification error,ME)[3,7]來(lái)定量評(píng)估分割方法的分割精度.ME是定量評(píng)估分割精度的一個(gè)常用指標(biāo),反映了分割結(jié)果圖像中背景像素被誤分為前景像素,以及前景像素被誤分為背景像素的情況,其計(jì)算公式為:

其中,Fg和Bg分別表示分割參考圖像中的前景和背景,而Ft和Bt分別表示使用t閾值化原始圖像所得二值圖像中的前景和背景.符號(hào)T表示取交集運(yùn)算,符號(hào)|·|表示計(jì)算元素的個(gè)數(shù).當(dāng)閾值化所得二值圖像和分割參考圖像相同時(shí),ME等于0%;當(dāng)閾值化所得二值圖像和分割參考圖像完全相反時(shí),ME等于100%.

提出的MNMI方法和交互式閾值(Interactive thresholding,IT)方法[36]、廣義熵閾值(Generalized entropy thresholding,GET)方法[15]、迭代三類(lèi)閾值(Iterative triclass thresholding,ITT)方法[19]、過(guò)渡區(qū)閾值(Transition region thresholding,TRT)方法[7]以及CKMI方法[34]進(jìn)行了比較.值得注意的是,IT方法人工交互式選擇一個(gè)閾值,該閾值對(duì)應(yīng)的二值圖像具有最小的ME值.因此,IT方法可以作為其他比較方法在分割精度方面的參考.

3.2 在合成圖像上的比較實(shí)驗(yàn)

9幅代表性合成圖像被用于檢驗(yàn)6個(gè)分割方法對(duì)不同直方圖模式的適應(yīng)能力.在這9幅圖像中,4幅圖像的灰度直方圖呈現(xiàn)出雙峰特征(見(jiàn)圖2),3幅圖像的灰度直方圖呈現(xiàn)出單峰特征(見(jiàn)圖3),1幅圖像的灰度直方圖呈現(xiàn)出無(wú)峰特征(見(jiàn)圖4),1幅圖像的灰度直方圖呈現(xiàn)出復(fù)雜的多峰特征(見(jiàn)圖5).

在灰度直方圖呈現(xiàn)出雙峰特征的4幅圖像中,各有1幅圖像的前景或背景的灰度分布呈現(xiàn)為正態(tài)分布(見(jiàn)圖2(a))、極值分布(見(jiàn)圖2(b))、瑞利分布(見(jiàn)圖2(c)),以及貝塔分布(見(jiàn)圖2(d)).表1顯示:1)對(duì)雙峰正態(tài)分布情況,GET和CKMI方法的分割結(jié)果比其他4種方法差.特別是GET方法,其計(jì)算出的閾值和IT方法選擇的閾值相差31個(gè)灰度級(jí).2)對(duì)雙峰極值分布、瑞利分布和貝塔分布三種情況,GET、ITT、TRT和CKMI方法的分割精度明顯下降,它們各自計(jì)算出的閾值和IT方法選擇的閾值的平均絕對(duì)偏差分別達(dá)到39,69,40和73個(gè)灰度級(jí),反映在ME值方面,GET、ITT、TRT和CKMI方法在這三種情況下的平均ME值分別達(dá)到了16.22%,33.48%,13.03%和34.07%.3)對(duì)雙峰正態(tài)分布、極值分布和貝塔分布三種情況,MNMI方法和IT方法的分割結(jié)果完全一致.在雙峰瑞利分布情況下,MNMI方法和IT方法所得閾值僅相差1個(gè)灰度級(jí).

對(duì)灰度直方圖呈現(xiàn)為單峰特征的情況,所合成的3幅圖像的前景或背景的灰度分布分別呈現(xiàn)為正態(tài)分布(見(jiàn)圖3(a))、極值分布(見(jiàn)圖3(b))以及瑞利分布(見(jiàn)圖3(c)).由于前景在整幅圖像中所占比例僅為1%,因此圖像的整體直方圖呈現(xiàn)為單峰特征.表2顯示:1)GET和TRT方法僅在單峰正態(tài)分布情況下獲得較好的分割結(jié)果,在單峰極值分布和瑞利分布情況下誤分割都較嚴(yán)重.2)在單峰的三種情況下,ITT和CKMI方法的ME值都超過(guò)35%,誤分割非常嚴(yán)重.3)對(duì)單峰的三種情況,MNMI方法計(jì)算出的閾值始終和IT方法選擇的閾值保持一致,因此它們的ME值也保持一致.

圖2 灰度直方圖呈現(xiàn)出雙峰特征的4幅合成圖像.每幅子圖的左邊顯示了合成圖像,右邊顯示了對(duì)應(yīng)的灰度直方圖,黑色曲線顯示了提出的MNMI方法計(jì)算閾值的目標(biāo)函數(shù)曲線,虛線及旁邊的數(shù)字標(biāo)示了MNMI方法計(jì)算出的閾值(下同)((a)正態(tài)分布情形;(b)極值分布情形;(c)瑞利分布情形;(d)貝塔分布情形)Fig.2 4 synthetic images with bimodal gray level histogram.In each sub- figure,a synthetic image is shown on the left;on the right,the gray level histogram is shown,the objective function of MNMI method is illustrated with a black curve,and the threshold obtained by MNMI method is indicated with a dashed line and a number(the same below)((a)Normal distribution,(b)Extreme value distribution,(c)Rayleigh distribution,(d)Beta distribution)

圖3 灰度直方圖呈現(xiàn)出單峰特征的3幅合成圖像((a)正態(tài)分布情形;(b)極值分布情形;(c)瑞利分布情形)Fig.3 3 synthetic images with unimodal gray level histogram((a)Normal distribution,(b)extreme value distribution,(c)Rayleigh distribution)

圖4 灰度直方圖呈現(xiàn)出無(wú)峰特征的合成圖像Fig.4 A synthetic image with a uniform distribution of gray level

表1 6個(gè)閾值分割方法在灰度直方圖呈現(xiàn)出雙峰特征的4幅合成圖像上的閾值t和ME值Table 1 Threshold values t and ME values of 6 thresholding methods on 4 synthetic images with bimodal gray level histogram

圖5 灰度直方圖呈現(xiàn)出多峰特征的合成圖像(灰度直方圖的灰度區(qū)間[0,50]由瑞利分布和均勻分布組合而成,區(qū)間[51,100]為均勻分布,區(qū)間[101,150]由極值分布和均勻分布組合而成,區(qū)間[151,200]由貝塔分布和正態(tài)分布組合而成,區(qū)間[201,255]為正態(tài)分布)Fig.5 A synthetic image with multimodal gray level histogram(The gray level histogram in the range[0,50]is combined by a Rayleigh distribution and a uniform distribution,[51,100]by a uniform distribution,[101,150]by an extreme value distribution and a uniform distribution,[151,200]by a beta distribution and a normal distribution,and[201,255]by a normal distribution.)

表2 6個(gè)閾值分割方法在灰度直方圖呈現(xiàn)出單峰特征的3幅合成圖像上的閾值t和ME值Table 2 Threshold values t and ME values of 6 thresholding methods on 3 synthetic images with unimodal gray level histogram

表3 6個(gè)閾值分割方法在灰度直方圖分別呈現(xiàn)出無(wú)峰特征和多峰特征的合成圖像上的閾值t和ME值Table 3 Threshold values t and ME values of 6 thresholding methods on synthetic images with uniform gray level histogram and multimodal one

灰度直方圖呈現(xiàn)為無(wú)峰特征的測(cè)試圖像是通過(guò)添加均勻分布的噪聲而合成的(見(jiàn)圖4).因?yàn)闊o(wú)法利用雙峰直方圖的峰谷特征或者單峰直方圖的凹陷特征,閾值化此類(lèi)合成圖像比閾值化前面兩類(lèi)合成圖像更有挑戰(zhàn)性.表3第2列數(shù)據(jù)顯示:1)GET、ITT、TRT和CKMI方法計(jì)算出的閾值離最優(yōu)閾值230分別差50,104,56和102個(gè)灰度級(jí),結(jié)果這4種方法均存在較為嚴(yán)重的誤分割.2)MNMI方法依然能精確地計(jì)算出最優(yōu)閾值230,獲得和IT方法一樣的分割結(jié)果.

鑒于真實(shí)世界圖像的灰度直方圖經(jīng)常呈現(xiàn)為多峰特征,我們也合成了一幅具有多峰特征的測(cè)試圖像(見(jiàn)圖5).注意這幅圖像的前景只形成了一個(gè)峰,其余峰都是由背景所形成.因此,分割這幅圖像的難度并不亞于分割具有無(wú)峰直方圖特征的圖像.表3第3列數(shù)據(jù)顯示:1)MNMI方法和IT方法計(jì)算出相等的閾值204,兩者輸出一樣的分割結(jié)果.2)但是GET、ITT、TRT和CKMI方法的誤分割嚴(yán)重,它們對(duì)應(yīng)的ME值分別高達(dá)17.75%,55.60%,27.89%和55.90%,它們計(jì)算出的閾值和最優(yōu)閾值204分別差24,94,50和96個(gè)灰度級(jí).

3.3 在真實(shí)世界圖像上的比較實(shí)驗(yàn)

59幅真實(shí)世界圖像被用于進(jìn)一步檢驗(yàn)6個(gè)分割方法對(duì)不同直方圖模式的適應(yīng)能力.這些測(cè)試圖像的灰度直方圖可用正態(tài)分布、極值分布、瑞利分布、貝塔分布和均勻分布中的若干分布的混合進(jìn)行逼近.在這59幅圖像中,編號(hào)1到19的圖像的灰度直方圖呈現(xiàn)出雙峰特征,編號(hào)20到36的圖像的灰度直方圖呈現(xiàn)出單峰特征,編號(hào)37到40的圖像的灰度直方圖呈現(xiàn)出無(wú)峰特征,編號(hào)41到59的圖像的灰度直方圖呈現(xiàn)出多峰特征(大于等于3峰).

圖6顯示了6種方法在59幅測(cè)試圖像上分割結(jié)果的量化比較.從該圖可以觀察到:1)對(duì)具有雙峰特征、單峰特征、無(wú)峰特征或者多峰特征的圖像,IT方法和MNMI方法對(duì)應(yīng)的ME均值都小于0.6%;而其他方法,除了TRT在單峰特征情況外,每種方法在不同灰度直方圖特征情況下對(duì)應(yīng)的ME均值都大于4.3%;2)在灰度直方圖呈現(xiàn)為無(wú)峰特征情況下,GET、ITT、TRT和CKMI 4種方法的誤分割都很?chē)?yán)重,各自的ME均值都超過(guò)了15%;3)ITT和CKMI方法在單峰和多峰特征情況下,誤分割情況也很?chē)?yán)重,各自的ME均值都超過(guò)了16%,并且相對(duì)而言,CKMI方法的誤分割更嚴(yán)重;4)僅就GET、ITT、TRT和CKMI 4種方法自身比較而言,GET方法和ITT方法分別在多峰和雙峰直方圖特征情況下略有優(yōu)勢(shì),而TRT方法則更適合分割具有單峰直方圖特征的圖像.

IT方法以最小化ME值為標(biāo)準(zhǔn)人工交互式選擇閾值,因此它選擇的閾值是最小化ME值意義下的最優(yōu)閾值.這樣,通過(guò)比較其他方法所得閾值和IT方法所選閾值之間的差異,可以衡量其他方法在閾值計(jì)算方面的偏離程度.從圖7可以觀察到,MNMI、GET、ITT、TRT、CKMI 5種方法和IT方法的閾值差分別分布在區(qū)間[?9,9],[?60,76],[?93,85],[?53,78]和 [?105,85]內(nèi).這里,每個(gè)區(qū)間的上下限反映了相應(yīng)方法在閾值計(jì)算方面最大可能的偏離度(針對(duì)此59幅測(cè)試圖像而言).進(jìn)一步,如果將大于0的閾值差稱為正向閾值差,將小于0的閾值差稱為反向閾值差,那么可以計(jì)算出每種方法對(duì)應(yīng)的正向閾值差的均值和反向閾值差的均值,這兩個(gè)均值可以用于衡量每種方法在閾值計(jì)算方面的平均偏離度.MNMI、GET、ITT、TRT和CKMI 5種方法對(duì)應(yīng)的平均偏離度區(qū)間分別為 [?4.38,3.32],[?20.91,25.21],[?41.17,20.93],[?14.17,22.24]和[?48.1,30.26].這組數(shù)據(jù)直接反映了:1)GET、ITT、TRT和CKMI 4種方法計(jì)算出的閾值和最優(yōu)閾值總體上相差甚遠(yuǎn);2)總體而言,MNMI方法在不同灰度直方圖模式下具有更精確的閾值計(jì)算能力,或者說(shuō),它對(duì)不同灰度直方圖模式的圖像具有更強(qiáng)的適應(yīng)能力.

3.4 計(jì)算效率的比較實(shí)驗(yàn)

通過(guò)比較不同方法在相同硬件和軟件條件下的CPU耗時(shí),可以直觀地反映不同方法在計(jì)算效率方面的差異.在相同的硬件和軟件條件下,同一個(gè)程序在不同的時(shí)間點(diǎn)運(yùn)行,其CPU耗時(shí)通常會(huì)有輕微的浮動(dòng).為了降低這種浮動(dòng)所產(chǎn)生的負(fù)面效應(yīng),讓每個(gè)參與比較的分割方法在同一幅測(cè)試圖像上連續(xù)運(yùn)行20次,并以這20次運(yùn)行所耗時(shí)間的均值作為該分割方法在該測(cè)試圖像上的CPU耗時(shí).以此為基礎(chǔ),可以進(jìn)一步計(jì)算出每個(gè)分割方法在9幅合成圖像和59幅真實(shí)世界圖像上CPU耗時(shí)的均值和標(biāo)準(zhǔn)偏差(見(jiàn)表4).

表4 5個(gè)自動(dòng)閾值分割方法的計(jì)算效率比較Table 4 Quantitative comparisons of 5 automatic thresholding methods in computational efficiency

圖6 6個(gè)閾值方法在59幅真實(shí)世界圖像上分割精度的量化比較(在每幅子圖中,各條水平虛線標(biāo)示了對(duì)應(yīng)情形下ME值的平均值)Fig.6 Quantification comparisons of segmentation accuracy for 6 thresholding methods on 59 real-world images(In each sub-figure,each horizontal dashed lines indicate the corresponding average ME,respectively)

如表4所示,總體而言ITT和GET兩種方法的計(jì)算效率相對(duì)更高,而CKMI、TRT和MNMI 3種方法的計(jì)算效率相對(duì)更低.ITT方法只需預(yù)先統(tǒng)計(jì)原始圖像X的灰度分布信息,然后進(jìn)行相對(duì)簡(jiǎn)單的均值和方差的運(yùn)行,因此它的CPU耗時(shí)最少.GET方法也需要預(yù)先統(tǒng)計(jì)原始圖像X的灰度分布信息,不過(guò)由于它計(jì)算廣義熵涉及到對(duì)數(shù)運(yùn)算,因此GET方法的計(jì)算代價(jià)比ITT方法略高.CKMI方法需要一個(gè)k均值算法估計(jì)初始閾值,然后涉及到對(duì)原始圖像X和二值圖像Yt的香農(nóng)熵和聯(lián)合香農(nóng)熵的計(jì)算,因此它的計(jì)算量比GET方法多一些,相應(yīng)地也需要花費(fèi)更多的CPU處理時(shí)間.TRT方法的主要計(jì)算代價(jià)發(fā)生在多尺度梯度乘變換和平穩(wěn)過(guò)渡區(qū)計(jì)算上,而MNMI方法除了涉及多尺度梯度乘變換的計(jì)算外,還涉及二值圖像輪廓提取操作,以及圖像MX和CYt的香農(nóng)熵和聯(lián)合香農(nóng)熵的計(jì)算,因此MNMI方法的計(jì)算代價(jià)最大,CPU耗時(shí)最多.

為了提高M(jìn)NMI方法的計(jì)算效率,可以將計(jì)算和尋找最大歸一化互信息量的目標(biāo)分解為兩個(gè)過(guò)程來(lái)完成:跳躍式處理和逐級(jí)式處理.在跳躍式處理過(guò)程中,對(duì)原始圖像X的灰度區(qū)間[tmin,tmax],將算法1中步驟3的逐個(gè)灰度處理,改為以一定步長(zhǎng)ρ跳躍式處理,而計(jì)算歸一化互信息量的其他步驟不變.跳躍式處理過(guò)程執(zhí)行完畢,輸出一個(gè)灰度值ttemp,它對(duì)應(yīng)于跳躍式處理過(guò)程中產(chǎn)生的最大歸一化互信息量.然后,在逐級(jí)式處理過(guò)程中,對(duì)灰度區(qū)間[ttemp?ρ+1,ttemp+ρ?1]逐個(gè)處理,即將算法1中步驟3的灰度區(qū)間[tmin,tmax]替換為[ttemp?ρ+1,ttemp+ρ?1],而計(jì)算歸一化互信息量的其他步驟不變.逐級(jí)式處理過(guò)程執(zhí)行完畢,輸出最終的分割閾值tfinal,它對(duì)應(yīng)于灰度區(qū)間[ttemp?ρ+1,ttemp+ρ?1]內(nèi)的最大歸一化互信息量.如表5所示,不管是合成圖像還是真實(shí)世界圖像,1)當(dāng)步長(zhǎng)ρ從1變化到6時(shí),平均CPU耗時(shí)都逐漸減少,但是平均ME值基本保持不變;2)當(dāng)步長(zhǎng)ρ從7變化到15時(shí),平均CPU耗時(shí)基本維持穩(wěn)定,但是平均ME值出現(xiàn)波動(dòng).這樣,當(dāng)步長(zhǎng)ρ取4,5或6時(shí),MNMI方法可以在保持分割精度的同時(shí),在計(jì)算效率方面接近于CKMI方法和TRT方法.

圖7 對(duì)59幅真實(shí)世界圖像,MNMI,GET,ITT,TRT,CKMI方法和IT方法獲得的閾值之差(在每幅子圖中,每個(gè)黑點(diǎn)旁的數(shù)字標(biāo)示了相應(yīng)的閾值差)Fig.7 The differences between MNMI,GET,ITT,TRT,CKMI and IT methods in segmentation thresholds for 59 real-world images(In each sub- figure,the number beside each black point labels the specific difference of segmentation threshold)

表5 不同步長(zhǎng)下MNMI方法的計(jì)算效率和誤分類(lèi)率Table 5 Computational efficiency and ME of MNMI method with different steps

4 結(jié)論

對(duì)灰度直方圖呈現(xiàn)為雙峰、單峰、無(wú)峰或者多峰特征的圖像,當(dāng)前景或背景的灰度分布可以用正態(tài)分布、極值分布、瑞利分布、貝塔分布、均勻分布或者它們的混合分布逼近時(shí),MNMI方法雖然在計(jì)算效率方面不優(yōu)于GET、ITT、TRT和CKMI方法,但是MNMI方法在自動(dòng)閾值選擇方面具有更穩(wěn)健的適應(yīng)能力.MNMI方法的穩(wěn)健性和適應(yīng)性得益于:1)MNMI方法利用多尺度梯度乘變換進(jìn)行了圖像的規(guī)范化,可將具有不同灰度直方圖模式的原始圖像轉(zhuǎn)化為具有單峰長(zhǎng)拖尾分布特征的規(guī)范圖像;2)MNMI方法計(jì)算規(guī)范圖像和二值輪廓圖像之間的統(tǒng)計(jì)相關(guān)性,兼顧了規(guī)范圖像和二值輪廓圖像各自提供的信息;3)MNMI方法引入基于幾何均值的歸一化互信息量,有利于應(yīng)對(duì)二值輪廓圖像中值為1的像素個(gè)數(shù)和值為0的像素個(gè)數(shù)比例失衡的情況.