梁義濤, 孟亞敏, 朱玲艷， 張 猛， 李永剛

(1.河南工業(yè)大學(xué)糧食信息處理與控制教育部重點(diǎn)實(shí)驗(yàn)室， 鄭州 450001; 2.河南工業(yè)大學(xué)信息科學(xué)與工程學(xué)院， 鄭州 450001)

圖像分割是圖像處理和分析中的關(guān)鍵一環(huán)，而閾值分割是圖像分割中使用較多的基本和關(guān)鍵技術(shù)之一[1]。在閾值分割技術(shù)中，Otsu法[2]因其良好的實(shí)用化模型特點(diǎn)已被廣泛應(yīng)用在分類識別[3-6]、計(jì)算機(jī)視覺[7-10]等多個(gè)領(lǐng)域。1993年，針對經(jīng)典Otsu法抗噪能力較差的問題，劉建莊等[11]提出了基于灰度圖像的二維Otsu法。此方法考慮了像素本身及其鄰域結(jié)構(gòu)信息兩方面影響因素，有很好的抗噪聲性。但該方法假定了邊界區(qū)域信息對于圖像分割的影響可以忽略。此假設(shè)使得經(jīng)典的二維Otsu法對于邊緣信息較為豐富的圖像分割效果不佳。因此，后續(xù)研究者給出了很多改進(jìn)算法：范九倫等[12]提出了曲線閾值型Otsu法及其快速算法，經(jīng)典的二維Otsu法可視為該方法的一種特例。但由于曲線閾值的解析式階次需要提前人為設(shè)定，該文獻(xiàn)快速算法依然是基于直線閾值近似地給出。如何針對不同圖像自適應(yīng)獲得最優(yōu)的曲線閾值解析表達(dá)依然是一個(gè)值得進(jìn)一步探討的問題。在直線閾值的求解方面，文獻(xiàn)[13-15]結(jié)合對二維直方圖的先驗(yàn)分析，給出了二維直方圖區(qū)域斜分方法及改進(jìn)，同時(shí)還將線閾值推廣到法線與灰度級軸成任意θ角的直線。但該方法在實(shí)際應(yīng)用時(shí)均固定選取θ=45°情況時(shí)直線閾值實(shí)現(xiàn)分割。且文獻(xiàn)[14]指出需根據(jù)實(shí)際圖像的信噪比及結(jié)果需求人為介入直線閾值θ角的值，使得該方法的自適應(yīng)能力降低。受文獻(xiàn)[16]的曲線閾值理論的啟發(fā)，梁義濤等[17]提出了灰度圖像二維Otsu折線閾值分割法。該方法本質(zhì)上是曲線閾值思想在實(shí)際應(yīng)用中的一種近似。仿真結(jié)果表明該方法對邊界信息豐富的圖像也能保持良好的分割效果。但是該方法需要事先針對不同的圖像設(shè)定二維直方圖迭代處理次數(shù)的經(jīng)驗(yàn)值，降低了實(shí)際應(yīng)用中自適應(yīng)程度，且未對噪聲污染圖像的適用情況進(jìn)行分析。

現(xiàn)在折線閾值法迭代思想的基礎(chǔ)上進(jìn)行改進(jìn)，提出二維Otsu擬合線閾值分割法。改進(jìn)和拓展主要涉及兩方面：一是明確給出迭代停止的量化條件，提高算法的自適應(yīng)能力；二是針對折線閾值法分類賦值過程煩瑣、效率低的局限性，提出使用線閾值完成圖像分割的思路。最后選取邊界信息較為豐富的圖像及含噪圖像進(jìn)行仿真實(shí)驗(yàn)，給出實(shí)驗(yàn)結(jié)果和分析，并從分割結(jié)果、抗噪性能和運(yùn)行時(shí)間與相關(guān)算法進(jìn)行比較與分析。

1 二維Otsu曲線閾值分割法

對于一幅M×N的數(shù)字圖像，用f(x,y)表示圖像上坐標(biāo)為(x,y)的像素點(diǎn)的灰度，g(x,y)表示圖像上坐標(biāo)為(x,y)的像素點(diǎn)的3×3鄰域平均灰度，g(x,y)的定義為

(1)

式(1)中：?」為取整；m為像素點(diǎn)左右位置；n為像素點(diǎn)上下位置。從g(x,y)的定義可以看出，如果圖像的灰度級為L，那么相應(yīng)的像素鄰域平均灰度的灰度級也為L，f(x,y)和g(x,y)組成的二元組記為(i,j)。在此基礎(chǔ)上定義圖像的二維直方圖，該二維直方圖定義在一個(gè)L×L大小的正方形區(qū)域，其橫坐標(biāo)表示圖像像元的灰度值,縱坐標(biāo)表示像元的鄰域平均灰度值。直方圖任意一點(diǎn)的值定義為pij，它表示二元組(i,j)發(fā)生的頻率。pij由式(2)確定，即

(2)

圖1中，對角線上的兩個(gè)區(qū)域Ⅰ和Ⅲ分別對應(yīng)于背景和目標(biāo)，遠(yuǎn)離對角線的區(qū)域Ⅱ和Ⅳ對應(yīng)于邊緣和噪聲。在二維Otsu法中認(rèn)為區(qū)域Ⅱ和Ⅳ的概率值近似為0，這種假設(shè)致使二維Otsu法忽略了邊界區(qū)域的信息，造成了它在某些場合是不適用的。

圖1 二維直方圖平面投影圖Fig.1 Plane projection of two-dimensional histogram

文獻(xiàn)[12]提出曲線閾值型Otsu方法。在該方法中，如果(s,t)是選取的閾值點(diǎn)，作過(s,t)的曲線r(i,j)將二維區(qū)域分成C0(s,t)和C1(s,t)兩塊，分別表示背景和目標(biāo)，二維Otsu曲線閾值直方圖如圖2所示。

圖2 二維Otsu曲線閾值直方圖Fig.2 Threshold histogram of two-dimensional Otsu curve

背景和目標(biāo)出現(xiàn)的概率分別為P0和P1,公式為

(3)

(4)

背景和目標(biāo)對應(yīng)的均值矢量為U0和U1，其表達(dá)式分別為

(5)

(6)

二維直方圖上總的均值矢量為UT

(7)

容易證明

P0+P1=1

(8)

UT=P0U0+P1U1

(9)

定義類間的離差矩陣SB為

SB=P0[(U0-UT)(U0-UT)′]+P1[(U1-UT)(U1-UT)′]=

(10)

使用SB的跡作為類間的離散度測度,有

trSB=P0[(U00-UT0)2+(U01-UT1)2]+

P1[(U10-UT0)2+(U11-UT1)2]=

(11)

最佳閾值(s*,t*)由式(12)確定

(12)

由圖2可以得出曲線r(i,j)取為(0,t)→(s,t)→(s,0)構(gòu)成的折線，曲線Otsu選取r(i,j)為過(s,t)且垂直于二維直方圖的定義域?qū)蔷€的直線。通過i+j=s*+t*的直線，對原始圖形進(jìn)行分割，歸類方式為

(13)

式(13)稱為基于直線進(jìn)行閾值選取的方法為二維直線閾值型Otsu法。

2 二維Otsu擬合線閾值分割法

2.1 算法思想介紹

在給出曲線Otsu法原理性描述之后，文獻(xiàn)[12]指出，如何求解線閾值是一個(gè)關(guān)鍵問題。受此啟發(fā)，基于逼近的思想和前期折線閾值法的研究，在二維Otsu折線閾值算法和曲線擬合方法的基礎(chǔ)上，提出了另外一種易于工程實(shí)現(xiàn)的近似算法——二維Otsu擬合線閾值法。該方法用原二維Otsu法對邊界信息或噪聲所屬區(qū)域的像素點(diǎn)迭代處理，獲得多個(gè)閾值點(diǎn)；再將多個(gè)閾值點(diǎn)擬合成線閾值；以此線閾值作為分割標(biāo)準(zhǔn)實(shí)現(xiàn)分割。

假設(shè)由二維Otsu法獲取最佳閾值點(diǎn)為(s*,t*)，t表示擬合的多項(xiàng)式，也用f(s)表示，a0表示常數(shù)項(xiàng)，n表示擬合曲線的階次，an擬合多項(xiàng)式的系數(shù)，則過點(diǎn)(s*,t*)的最佳線閾值解析式可表示為

(14)

此時(shí)，問題轉(zhuǎn)化為閾值曲線上點(diǎn)的分布區(qū)域和階次的選取問題。由于閾值點(diǎn)(s*,t*)的存在，Ⅰ和Ⅲ區(qū)域的像素點(diǎn)是已經(jīng)明確分類的像素點(diǎn)，而Ⅱ區(qū)域和Ⅳ區(qū)域包含的像素點(diǎn)為未分類的像素點(diǎn)。進(jìn)而可知，閾值曲線上的各點(diǎn)必然分布在Ⅱ區(qū)域或Ⅳ區(qū)域，此時(shí)假定最終獲取的最佳線閾值如圖3所示。

圖3 最佳線閾值直方圖Fig.3 Best line threshold histogram

任取曲線上一點(diǎn)，以該點(diǎn)為中心可將Ⅱ區(qū)域或Ⅳ區(qū)域細(xì)分為4個(gè)子區(qū)域。如此遍歷曲線上點(diǎn)細(xì)分二維直方圖，則可實(shí)現(xiàn)所有像素點(diǎn)分類。

2.2 算法流程介紹

輸入：一幅彩色圖像。

輸出：分割好的黑白圖像。

(1)灰度化并計(jì)算其鄰域平均灰度值，組成圖像二維直方圖。

(2)利用原始二維Otsu法，得到閾值點(diǎn)。

(3)判斷Ⅱ、Ⅳ區(qū)域像素占整幅圖像像素的比率是否小于0.02，如果比率和小于0.02，停止處理；如果比率和大于或等于0.02，則Ⅱ、Ⅳ區(qū)域成為新的待處理區(qū)域，條件不滿足則不斷迭代處理。

(4)判斷閾值個(gè)數(shù)是否大于1，如果閾值個(gè)數(shù)大于1，用最小二乘法將這些點(diǎn)擬合成階次為n的曲線，稱為線閾值,以此線閾值作為分割標(biāo)準(zhǔn)，如果閾值個(gè)數(shù)等于1，以此閾值作為分割標(biāo)準(zhǔn)。

(5)根據(jù)分割標(biāo)準(zhǔn)，得到二值化圖像，實(shí)現(xiàn)分割。算法流程如圖4所示。

圖4 算法流程圖Fig.4 Algorithm flow chart

2.3 迭代處理

本文方法在文獻(xiàn)[17]折線閾值法的迭代思想基礎(chǔ)上提出迭代停止條件為該類像素點(diǎn)概率滿足小概率事件的概率標(biāo)準(zhǔn)要求(稱為小概率設(shè)定值ε)時(shí)停止迭代。

具體的過程說明如下：用二維Otsu法處理，獲取初次閾值點(diǎn)(s*,t*)，依照點(diǎn)(s*,t*)將坐標(biāo)系對應(yīng)區(qū)域分為四個(gè)子區(qū)域：Ⅰ和Ⅲ區(qū)域?yàn)槟繕?biāo)和背景區(qū)域，Ⅱ和Ⅳ區(qū)域?yàn)檫吘壭畔⒒蛟肼?。統(tǒng)計(jì)屬于Ⅱ區(qū)域和Ⅳ區(qū)域灰度像素的概率和，判斷Ⅱ和Ⅳ區(qū)域的概率和是否小于ε。如果概率和大于或等于ε，則利用二維Otsu法再分別處理區(qū)域Ⅱ和Ⅳ區(qū)域，再次獲取對應(yīng)閾值點(diǎn)，如圖5所示。如果概率和小于ε，則停止迭代。

圖5 迭代處理圖Fig. 5 Iterative processing diagram

利用二維Otsu進(jìn)行迭代對Ⅱ、Ⅳ區(qū)域的像素逐級分類處理，最終達(dá)到使得Ⅱ、Ⅳ區(qū)域像素?cái)?shù)達(dá)到可忽略的標(biāo)準(zhǔn)(概率上可以認(rèn)為Ⅱ、Ⅳ區(qū)域未分類像素?cái)?shù)占整幅圖像總像素的比率滿足小概率事件標(biāo)準(zhǔn))。

從不同設(shè)定值ε的運(yùn)行時(shí)間進(jìn)行分析比較研究。為方便比較起見，其中，擬合曲線選擇過點(diǎn)(s*,t*)的兩段直線[以點(diǎn)(s*,t*)為界，分別在Ⅱ區(qū)域和Ⅳ區(qū)域各擬合一條直線],噪聲密度標(biāo)準(zhǔn)差設(shè)定為20。

通過對6幅圖像(圖6)進(jìn)行實(shí)驗(yàn)，統(tǒng)計(jì)了6幅圖像的Ⅱ、Ⅳ區(qū)域像素占整幅圖像像素的百分比，如表1所示，分別用n(Ⅱ)/N和n(Ⅳ)/N表示；Ⅱ區(qū)域和Ⅳ區(qū)域的像素占整幅圖像像素的百分比，用n(Ⅱ、Ⅳ)/N表示。其中，n表示該區(qū)域的像素?cái)?shù)，N表示整幅圖像的像素?cái)?shù)。

圖6 輸入圖像Fig.6 Input image

表1 圖像直方圖未分類區(qū)域像素比例

從表1統(tǒng)計(jì)的數(shù)據(jù)可看出，Ⅱ、Ⅳ區(qū)域的像素比例遠(yuǎn)大于可以忽略的條件，所以對其進(jìn)行處理是必要的。而且如果對邊界處理要求較高，Ⅱ、Ⅳ區(qū)域的處理結(jié)果也很重要。

在概率論中一般將概率在0.01～0.05的事件稱為小概率事件，所以表2中ε選取0.01、0.02、0.05做實(shí)驗(yàn)，統(tǒng)計(jì)了3幅圖像的閾值點(diǎn)和運(yùn)行時(shí)間，其中運(yùn)行時(shí)間是取3次實(shí)驗(yàn)的平均值。

表2 不同小概率設(shè)定值ε的影響比較

從表2中閾值點(diǎn)個(gè)數(shù)的數(shù)據(jù)可看出，ε為0.01時(shí)需要處理像素最多，ε為0.05時(shí)需要處理像素最少，隨著ε變大，運(yùn)行的時(shí)間變短。比如對于圖6的電路，ε選取0.01和0.02得到的閾值點(diǎn)個(gè)數(shù)相同，但對于ε選0.02來說，ε選0.01會處理像素更多，所以ε選0.01運(yùn)行時(shí)間長，ε選取0.05時(shí)迭代次數(shù)減少，運(yùn)行時(shí)間短。計(jì)算圖6的電路的運(yùn)行時(shí)間，ε取0.05比ε取0.01運(yùn)行速度提高了3%。

由圖7的實(shí)驗(yàn)結(jié)果可得，針對不同要求的實(shí)際應(yīng)用，可以選擇不同的ε值。選擇設(shè)定II、IV區(qū)域未分類像素占整幅圖像的百分比小于0.02達(dá)到可忽略不計(jì)的標(biāo)準(zhǔn)，以ε=0.02為例做實(shí)驗(yàn)。

圖7 不同小概率設(shè)定值的實(shí)驗(yàn)結(jié)果Fig.7 Experimental results of different small probability settings

2.4 擬合線閾值

實(shí)際中，擬合曲線的階次n不可能無限大，且擬合多項(xiàng)式階數(shù)n>5時(shí)，其法方程的系數(shù)矩陣是病態(tài)的。因此，考慮過擬合、實(shí)現(xiàn)效率等問題，這里假定n>3時(shí)，曲線閾值錯(cuò)分類像素點(diǎn)數(shù)可忽略。下面為分別以階數(shù)n=1、2、3進(jìn)行分割的比較研究實(shí)驗(yàn)。仿真結(jié)果如圖8所示。其中，小概率設(shè)定值為0.02，噪聲密度標(biāo)準(zhǔn)差設(shè)定為20。

圖8 不同階數(shù)的實(shí)驗(yàn)結(jié)果Fig.8 Experimental results of different orders

通過迭代處理得到的多個(gè)閾值點(diǎn)，用最小二乘法擬合這些閾值點(diǎn)得到線閾值，作為分割標(biāo)準(zhǔn)，從圖8可以看出，擬合曲線階數(shù)n=1處理效果很好，因此設(shè)定擬合曲線的階次為1。

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

實(shí)驗(yàn)硬件配置為Intel(R) Core(TM) i5-3210MCPU @2.50 GHz 4GB RAM的PC機(jī)，軟件使用MATLAB R2017b進(jìn)行實(shí)驗(yàn)。實(shí)驗(yàn)的算法包括原始二維Otsu法，改進(jìn)的灰度圖像二維Otsu折線閾值分割法，改進(jìn)的二維直方圖區(qū)域斜分算法及本文方法。通過實(shí)驗(yàn)，驗(yàn)證了本文方法的可行性，通過方法的對比，說明此方法提高了分割質(zhì)量和普適性。

3.1 分割性能實(shí)驗(yàn)

為驗(yàn)證本文方法的可行性，以下給出圖像的實(shí)驗(yàn)結(jié)果，實(shí)驗(yàn)圖像依次為：電路、電路板、細(xì)胞、西瓜、血管和馬。實(shí)驗(yàn)方法依次為原始二維Otsu法、二維Otsu折線法、二維Otsu斜分法、本文方法。

3.2 抗噪聲性能實(shí)驗(yàn)

為驗(yàn)證本文方法的抗噪性，用lean圖像做實(shí)驗(yàn)，在原圖加入不同程度的高斯噪聲，圖10為加入不同高斯噪聲不同算法實(shí)驗(yàn)結(jié)果。圖10(a)～圖10(d)、圖10(e)～圖10(h)、圖10(i)～圖10(l)分別加入標(biāo)準(zhǔn)差為10、20、30的高斯噪聲，使用的分割方法分別為原始二維Otsu法、二維Otsu折線法、二維Otsu斜分法、本文方法。

3.3 實(shí)驗(yàn)分析

圖9(e)～圖9(h)和圖9(u)～圖9(%)是存在較多噪聲點(diǎn)的一個(gè)典型情況，二維Otsu折線法和二維Otsu斜分法的處理效果比原始二維Otsu法的要好很多但是這兩種方法的自適應(yīng)能力較差，本文方法能自適應(yīng)處理并且處理效果較好。

由圖9(i)～圖9(l)可以看出，原始二維Otsu法和二維Otsu折線法處理的結(jié)果沒有分割出右下的細(xì)胞，二維Otsu斜分法和本文方法可以將弱對比度的細(xì)胞也分割出來，圖9(a)～圖9(d)的右上區(qū)域，因?yàn)楣庹赵蚨SOtsu法、二維Otsu折線法和二維Otsu斜分法處理結(jié)果不好，區(qū)域的線路沒有分割出來，本文方法將這部分區(qū)域的線路清晰分離出來。由圖9(a)～圖9(d)、圖9(i)～圖9(l)可知，對于對比度較弱或者有弱光影響下，本文方法處理效果好且自適應(yīng)性強(qiáng)。由圖9(q)～圖9(t)和圖9(m)～圖9(p)可以直觀地看出，本文方法對邊緣處理效果更好。

圖9 不同算法的實(shí)驗(yàn)結(jié)果Fig.9 Experimental results of different algorithms

圖10為本文方法的抗噪聲實(shí)驗(yàn)結(jié)果，可以看出本文方法比其他方法的抗噪性更好，并且噪聲越大越能突出本文方法的抗噪性越強(qiáng)。

圖10 加入不同高斯噪聲不同算法實(shí)驗(yàn)結(jié)果Fig.10 Experimental results of different algorithms with different Gaussian noise

與二維Otsu折線閾值法比較，二維Otsu折線閾值法需要自己確定迭代次數(shù)，降低了智能自動程度，本文方法的改進(jìn)主要是確定了迭代分割的次數(shù)，保持了原始二維Otsu法的自動性，處理時(shí)間也比二維Otsu折線閾值法快。

表3是本文方法得到的閾值點(diǎn)及線閾值，ε設(shè)定為0.02，n設(shè)定為1，噪聲標(biāo)準(zhǔn)是20，x表示像素的灰度，y表示像素的鄰域平均灰度。

表3 實(shí)驗(yàn)圖像的分割線閾值

3.4 分割質(zhì)量評價(jià)

主觀評價(jià)圖像分割質(zhì)量的方法主觀性強(qiáng)，因此對圖像分割質(zhì)量的客觀評價(jià)是必要的。

采用最大相關(guān)準(zhǔn)則(MCC)作為評價(jià)函數(shù)[17]。最大相關(guān)準(zhǔn)則是Yen等[18]基于最大熵準(zhǔn)則提的，相關(guān)數(shù)總量與分割質(zhì)量成正相關(guān)的關(guān)系。陳修橋等[19]將一維最大相關(guān)準(zhǔn)則推廣至二維。表4為四種算法分割結(jié)果的定量評價(jià)，本文算法在主客觀上一致，有較好的視覺效果和相關(guān)數(shù)指標(biāo)。

表4 不同算法分割結(jié)果的定量評價(jià)

峰值信噪比(PSNR)[20]是衡量算法抗噪性能的量化指標(biāo)，PSNR越大，表示算法的抗噪性能越好。從表5可知，本文算法和二維Otsu折線法隨著噪聲程度變大，PSNR變大；另外兩種算法隨著噪聲程度變大，PSNR變小。證明了本文算法的抗噪聲能力，并噪聲程度越大，本文算法效果越好。

表5 不同算法對噪聲圖像分割的PSNR數(shù)據(jù)對比

4 結(jié)論

(1)對于二維Otsu法及改進(jìn)的方法存在的問題：不合理忽略像素和普適性低，提出二維Otsu擬合線閾值圖像分割法，通過原始二維Otsu法對二維直方圖分區(qū)，若邊緣信息或噪聲區(qū)域占整幅圖像像素的比率達(dá)到在概率上可以忽略的小概率事件，則合適忽略，否則用原始二維Otsu法繼續(xù)迭代處理邊緣信息或噪聲區(qū)域，解決了經(jīng)典二維Otsu法在實(shí)際應(yīng)用中可能存在的假設(shè)前提不合理的問題；若邊緣信息或噪聲區(qū)域占整幅圖像像素的比率本就小于小概率事件概率，則直接用閾值進(jìn)行分割，解決了普適性低的問題。

(2)擬合閾值是對文獻(xiàn)[12]提出的曲線閾值的一種應(yīng)用近似，且解決了最優(yōu)閾值為曲線時(shí)得到解析式困難的問題。如果最終采用了一次解析式擬合閾值時(shí)，針對不同的圖像，其自動生成的線閾值斜率是不同的。這就解決了文獻(xiàn)[14]中線閾值的斜率需要人為設(shè)定的問題。

(3)在迭代處理環(huán)節(jié)，本文算法基于折線閾值法的思路進(jìn)行改進(jìn)，加入量化的迭代停止條件，提高了算法的自適應(yīng)能力。在分類賦值方式上，引入曲線擬合方法，相對折線閾值法而言，簡化實(shí)現(xiàn)過程，提高了算法執(zhí)行效率。

(4)實(shí)驗(yàn)結(jié)果表明，本文方法對邊緣信息豐富并且要求分割質(zhì)量高的圖像處理效果更好。本文方法保持了原始Otsu法的自適應(yīng)性，提升了算法的普適性。