復(fù)雜背景下基于AD-GAC模型和最大熵閾值法的葉片病斑分割

趙輝 芮修業(yè) 岳有軍 王紅君

摘要：旨在研究復(fù)雜背景下葉片病斑的分割。由于復(fù)雜背景會帶來巨大的噪聲，產(chǎn)生過多的邊緣和灰度值不均勻的區(qū)域，很容易導(dǎo)致過分割的現(xiàn)象，因此在復(fù)雜背景下，很難通過1次分割就完成對葉片病斑的分割。為了解決復(fù)雜背景下過分割的現(xiàn)象，提出兩步分割的策略。第1步先用筆者提出的各向異性擴散測地線活動輪廓模型（anisotropic diffusion geodesic active contour model，簡稱AD-GAC模型）進(jìn)行預(yù)分割，在此過程中構(gòu)造新的邊緣檢測函數(shù)（edge stop function，簡稱ESF）;第2步通過最大熵閾值法完成最終的分割。隨后，提取并計算預(yù)分割部分各像素灰度值的最大熵，以得到病斑部分與葉片部分的灰度值閾值，通過閾值來完成最后1步的分割。通過MATLAB仿真，可以證明該算法可以有效地將病斑從復(fù)雜背景下的葉片上分割出來。研究結(jié)果后續(xù)的病斑識別作了鋪墊。

關(guān)鍵詞：各向異性擴散;測地線活動輪廓;復(fù)雜背景;最大熵閾值法;病斑分割

中圖分類號： S126文獻(xiàn)標(biāo)志碼： A

文章編號：1002-1302（2019）18-0136-05

收稿日期：2018-06-14

基金項目：天津市科技計劃（編號：15ZXZNGX00290）;天津市農(nóng)業(yè)科技成果轉(zhuǎn)化與推廣項目（編號：201203060、201303080）。

作者簡介：趙?輝（1963—），男，天津人，博士，主要研究方向為復(fù)雜系統(tǒng)智能控制理論與應(yīng)用、農(nóng)業(yè)信息與精準(zhǔn)農(nóng)業(yè)智能監(jiān)測理論與技術(shù)。E-mail：zhaohui3379@126.com。

通信作者：芮修業(yè)，碩士研究生，主要研究方向為數(shù)字圖像處理與模式識別。E-mail：292318849@qq.com。

我國是一個農(nóng)業(yè)大國。在設(shè)施農(nóng)業(yè)中，如果我們能準(zhǔn)確地確定作物疾病的時間線，并及時作出處理，可以極大地提高國家農(nóng)業(yè)的經(jīng)濟收益。對于經(jīng)濟作物的病害診斷而言，一個巨大挑戰(zhàn)就是作物復(fù)雜的自然生長環(huán)境帶來的影響[1]。圖像分割技術(shù)對于后期關(guān)于病斑特征的提取和識別任務(wù)是至關(guān)重要的，并且分割效果的好壞將直接影響病斑的識別效果[2]。

圖像分割技術(shù)目前已經(jīng)在很多方面被提出，但是面對各種復(fù)雜的分割條件，還沒有一種方法是普遍適用的，尤其是對于復(fù)雜背景下的葉片病斑分割[3]。大多數(shù)主流分割算法，無論是基于區(qū)域的分割算法[4]還是基于邊緣的分割算法[5]，對于復(fù)雜背景下的葉片病斑分割效果都不盡如人意。

作為策略的第一步，筆者選擇活動輪廓模型進(jìn)行預(yù)分割。活動輪廓模型已經(jīng)被證明可以成功地運用于核磁共振成像（MRI）圖像的病斑分割[6]。然而，這種方法尚不能被成熟地運用于復(fù)雜背景下的植物葉片病斑分割?；顒虞喞难芯块_始于Snake模型的提出[7]。此模型通過在圖像中的最小化能量函數(shù)來完成對圖像中目標(biāo)的分割，最小化能量函數(shù)的方法是水平集方法（level set method，簡稱LSM）。在這個過程中，我們可以得到模型演化的偏微分方程。其中1種LSM是基于邊緣信息的[8]，GAC（測地線活動輪廓）模型[9]就是其中1個基于邊緣信息的模型。GAC模型最終可以得到1個閉合輪廓，這是1個可以利用的特征。GAC模型利用圖像的梯度信息去構(gòu)造能量函數(shù)，開始設(shè)置的初始輪廓隨著演化可以收斂在目標(biāo)物體的邊緣。但是GAC模型對噪聲很敏感，復(fù)雜背景下的高噪聲會影響曲線的演化，使得活動輪廓不能靠近葉片和病斑的邊緣。如果不能成功地獲得1個合適的預(yù)分割結(jié)果，那也不能得到合適的灰度值閾值，這也會導(dǎo)致最終的分割失敗。

為了克服GAC模型在復(fù)雜背景下不能完成對葉片病斑預(yù)分割的局限性。本研究提出的各向異性擴散測地線活動輪廓模型（anisotropic diffusion geodesic active contour model，簡稱AD-GAC模型）比傳統(tǒng)GAC模型的邊緣檢測函數(shù)（edge stop function，簡稱ESF）增加了各向異性擴散濾波器。各向異性擴散濾波器能夠在去除噪聲的同時增強邊緣信息[10]，本研究將原來的外部能量函數(shù)轉(zhuǎn)變成了可以自動去除噪聲、增強圖像梯度信息的新的外部能量函數(shù)。因此，用于AD-GAC模型內(nèi)部運算的是強化過的梯度信息，各向異性擴散作為外部能量函數(shù)的一部分，是一種技術(shù)的融合而不僅僅是簡單的堆積。此外，筆者用懲罰項作為內(nèi)部能量函數(shù)，可以避免重新初始化帶來的時間上的巨大消耗。最后，通過將預(yù)分割部分像素的灰度值放入與原圖維度相同的全零矩陣來提取預(yù)分割結(jié)果，然后，將提取部分的像素灰度值作為第2部分的輸入數(shù)據(jù)，通過計算最大熵得到閾值，完成最終的病斑分割。

1?相關(guān)工作基礎(chǔ)

1.1?各向異性擴散

在圖像處理和計算機領(lǐng)域，各項異型擴散的核心是Perona-Malik方程，所以各向異性擴散模型也稱為P-M模型。這是一種在保留圖像重要部分的同時消除圖像噪聲的技術(shù)。特別是對于圖像中的邊緣、線條或其他重要部分。從原理上看，各向異性擴散是在考慮目標(biāo)邊界之后平均邊界一側(cè)的像素值。各向異性擴散可以表示為如下公式[10]：

It=div[c（‖I‖）·I]I（t=0）=I0。（1）

式中：c（‖I‖）是擴散因子。

P-M模型的擴散方程可以表示為如下公式：

c（‖I‖）=exp[-（‖I‖/k）2]。（2）

式中：k是擴散門限，根據(jù)本試驗情況，將k值設(shè)為15。本研究將利用各向異性擴散后的圖像梯度信息來構(gòu)造新的ESF。

1.2?水平集方法

LSM的基本思想是將n維曲面演化為在（n+1）維空間中求水平集曲面的隱函數(shù)解。用1條平面閉合曲線隱含地表示水平集函數(shù)的零水平集：

C（t）={（x，y）|（x，y，t）=0}。（3）

水平集演化方程[11]如下：

t=F||（x，y，0）=0（x，y）。（4）

式中：F控制著水平集函數(shù)演化的速度;為梯度算子。

2?AD-GAC模型

在復(fù)雜背景下，傳統(tǒng)GAC模型的ESF在邊緣檢測方面存在對噪聲敏感的缺點。復(fù)雜背景下的一些陰影、水珠、交錯的遮擋的葉片都會產(chǎn)生很多無法定義的邊緣。由于梯度信息的變化，輪廓將不會收斂在預(yù)期的目標(biāo)邊緣。

AD-GAC的核心在于各向異性擴散濾波器，各向異性擴散濾波器的離散表達(dá)式可以寫成如下形式：

It+1p=Itp+λ1|ηp|∑q∈ηpc（Itp，q）Itp，q。（5）

式中：控制系數(shù)λ1是1個控制著總體擴散強度的常量，經(jīng)過多次試驗，在本研究中將其設(shè)置為0.2;|ηp|表示鄰域空間的大小，一般選擇4-型鄰域空間。以Ip作為邊緣檢測因子的一部分：

Is=Ip×Gσ。（6）

式中：Ip表示各向異性擴散后的圖像強度;Is表示對Ip進(jìn)行高斯濾波后的結(jié)果;Gσ是帶有標(biāo)準(zhǔn)差的高斯核函數(shù)。x軸、y軸方向上的梯度可以按如下公式進(jìn)行計算：

Ix，Iy=Ist。（7）

因此，新的ESF可以被定義為下式：

gd=11+f， f=I2x+I2y。（8）

利用gd，1個新的外部能量函數(shù)就可以自然地被定義為如下形式：

ζgd，λ，v（）=λ2∫Ωgd·δ（）||dxdy+v∫Ωgd·H（-）dxdy。（9）

式中：δ是單變量狄拉克函數(shù);H是海氏函數(shù);λ2是加權(quán)長度項的系數(shù);v是加權(quán)區(qū)域項的系數(shù)，如果初始輪廓在目標(biāo)的外面，v選擇正值，反之，v選擇負(fù)值。

經(jīng)過很多輪迭代之后，水平集函數(shù)會失去平滑性和距離特性[12]，為了避免這個系統(tǒng)缺陷，本研究利用懲罰項作為內(nèi)部能量函數(shù)[13]。用于最小化能量函數(shù)的最陡梯度下降流公式定義為如下形式：

t=μ[Δ-div||]+vgdδ（）+λ2δ（）divgd||。（10）

式中：μΔ-div||是懲罰項，也就是整個能量函數(shù)的內(nèi)部能量函數(shù)，其中，μ是懲罰項系數(shù)，當(dāng)時間步長T滿足T·μ<0.25時，水平集函數(shù)的演化可以保持穩(wěn)定，時間步長T可以根據(jù)具體情況進(jìn)行選擇。

經(jīng)過大量對比試驗，筆者發(fā)現(xiàn)，當(dāng)T=5，μ=0.04，迭代次數(shù)為1 000次時，AD-GAC模型的能量函數(shù)收斂效果最好。

3?分割的第2步

在分割的第2步，筆者將完成對于復(fù)雜背景下葉片病斑分割的最后一部分，利用預(yù)分割結(jié)果，本研究可以獲得葉片和病斑的灰度值閾值，從而區(qū)分葉片部分和病斑部分。

3.1?移除背景

為了移除背景，筆者首先構(gòu)造1個與原圖維度相同的全零矩陣，同時，獲得預(yù)分割結(jié)果部分的每個像素的坐標(biāo)。根據(jù)坐標(biāo)，將對應(yīng)的各向異性擴散后的圖像的像素灰度值按照坐標(biāo)對應(yīng)地放進(jìn)全零矩陣中。

整個移除背景的過程如下所示：

（1）假設(shè)各向異性擴散后的圖像（維度與原圖相同）是1個維度為5×5的矩陣，假設(shè)預(yù)分割結(jié)果如圖1所示。

（2）假設(shè)預(yù)分割部分像素的坐標(biāo)為（1，2），（1，3），（2，2），（2，3），（2，4），（3，2），（3，3），（3，4），（4，2），（4，3），（4，4），（5，2）。

（3）構(gòu)造1個維度為5×5的全零矩陣（圖2）。

（4）將上述坐標(biāo)對應(yīng)各向異性擴散的像素灰度值放進(jìn)全零矩陣中（圖3）。

算法描述——移除背景：

u：預(yù)分割結(jié)果的二值圖像

U：原圖像的單通道圖

p：各向異性擴散結(jié)果

步驟：

（1） 找到預(yù)分割區(qū)域所有像素的坐標(biāo)：

[row，col，v]=find（u>0）

（2） 將坐標(biāo)存在a1中：

a1=[row，col]

（3） 構(gòu)造1個與原圖維度相同的全零矩陣b1：

[x，y]=size（U）

b1=zeros（x，y）

（4） 將p中對應(yīng)坐標(biāo)的灰度值放進(jìn)b1中：

[t，t1]=size（a1）

for q=1 to t do

m1=a1（q，1）

n1=a1（q，1）

b1（m1，n1）=p（m1，n1）

end for

3.2?最大熵閾值

在去除背景后，只剩下主要病斑和一些帶有病斑的葉片。所以現(xiàn)在矩陣b1中含有3種類型的元素：（1）零元素（背景區(qū)域）;（2）葉片部分的灰度值;（3）病斑部分的灰度值。

將矩陣b1中的數(shù)據(jù)轉(zhuǎn)化為一維數(shù)據(jù)，去除零元素，作為最大熵閾值法[14]的輸入。

在去除背景后的圖中，含有2種不同的灰度值區(qū)域，那必然存在1個灰度值t，作為閾值，可以得到最佳分割結(jié)果。假設(shè)葉片區(qū)域的灰度值范圍是a～t，該區(qū)域的概率分布為∑ti=ah（i），病斑區(qū)域的灰度值范圍是（t+1）～b，該區(qū)域的概率分布為∑bi=t+1h（i）。概率分布應(yīng)該滿足：

∑tah（i）+∑bt+1h（i）=1。（11）

葉片區(qū)域像素灰度值的熵為

Hlesion（t）=-∑ti=ah（i）∑tj=ah（j）log2h（i）∑tj=ah（j）。（12）

病斑區(qū)域像素灰度值的熵為

Hleaf（t）=-∑bi=t+1h（i）∑bj=t+1h（j）log2h（i）∑bj=t+1h（j）。（13）

根據(jù)信息論，當(dāng)閾值能夠最佳區(qū)分2個區(qū)域的時候，熵應(yīng)該最大，所以最佳閾值可以通過最大化2種類型像素的熵來選擇：

T=ArgMaxt=a，…，bHlesion（t）+Hleaf（t）。（14）

根據(jù)最佳閾值，將b1矩陣中灰度值小于T的值設(shè)置為0，剩下的就是病斑區(qū)域。至此，本研究就完成了復(fù)雜背景下對于葉片病斑的分割。

4?結(jié)果與分析

為了評估筆者提出的算法框架的效果，本研究選擇幾幅復(fù)雜背景下的葉片病斑圖進(jìn)行測試，數(shù)據(jù)圖片來自本研究所在的試驗農(nóng)場，圖片的像素值分別為395×395和640×427。

根據(jù)本試驗內(nèi)容，參數(shù)的具體選擇如下：各向異性擴散光滑系數(shù)λ1=0.1，時間步長T=5，懲罰項系數(shù)μ=0.04，加權(quán)長度項系數(shù)λ2=5，加權(quán)區(qū)域項系數(shù)v=3，AD-GAC模型迭代次數(shù)為1 000次。程序運行的硬件環(huán)境如下：Windows7旗艦版（32-bit），內(nèi)核為i5-4590 3.30-GHz處理器，運行內(nèi)存為4 GB。

圖4展示了AD-GAC模型相對于GAC模型在復(fù)雜背景下對于葉片病斑預(yù)分割的優(yōu)越性，可以看出，由于原圖的復(fù)雜背景帶來過多的邊緣信息，導(dǎo)致GAC模型無法很好地收斂

到病斑以及帶有病斑的葉片周圍，過多的預(yù)分割區(qū)域意味著進(jìn)入下一步分割的像素灰度值有更大的范圍，會直接影響最大熵閾值法求得的最佳閾值的準(zhǔn)確度。

完整的分割過程如圖5所示。從圖5-b可以看出，利用gd構(gòu)建的外部能量函數(shù)能使水平集輪廓較為準(zhǔn)確地收斂在病斑以及帶有病斑的葉片周圍。如圖5-c所示，為了避免復(fù)雜背景帶來的干擾，將預(yù)分割結(jié)果從復(fù)雜背景中提取出來。從圖5-d可以看出，最終的分割結(jié)果有效避免了過分割現(xiàn)象，并且能夠得到準(zhǔn)確的病斑區(qū)域。

在另一個試驗中，病斑的灰度值小于葉片的灰度值，完整的分割過程如圖6所示。

為了更好地展現(xiàn)最終的分割效果，本研究將最終分割結(jié)果的背景設(shè)置為白色。用過漏檢率（R1）、過檢率（R2）來評估試驗結(jié)果。作為試驗結(jié)果的對比，本研究用基于高斯分布改進(jìn)C-V模型對同樣的圖片進(jìn)行分割[15]。對于基于高斯分布的C-V模型，筆者仔細(xì)選擇參數(shù)，讓其分割效果達(dá)到最好。

真實葉片病斑區(qū)域面積Sa由人工計算得到，算法分割出的病斑區(qū)域面積Sq通過size（data）得到，Si為Sa和Sq的交集區(qū)域。可以通過以下公式計算R1和R2：

R1=|Sa-Si|Sa;（15）

R2=|Sq-Si|Sa。（16）

由表1可知，用本研究提出的算法得到的R1和R2值遠(yuǎn)小于基于高斯分布改進(jìn)C-V模型得到的R1和R2值，這說明本研究算法得到的分割結(jié)果準(zhǔn)確度更高。

5?結(jié)論

為了解決復(fù)雜背景對于葉片病斑分割帶來的過分割問題，本研究提出了兩步分割的策略，通過加入各向異性擴散濾波器改造外部能量函數(shù)，得到了魯棒性更強的AD-GAC模型，從而克服了GAC模型對噪聲敏感的缺點。然后，通過提取預(yù)分割結(jié)果來去除復(fù)雜背景，將提取出來的像素灰度值作為一維數(shù)據(jù)，通過計算最大熵得到最佳閾值，將最終的病斑區(qū)域分割出來。試驗結(jié)果證明，本研究提出的算法能有效地完成復(fù)雜背景下對葉片病斑的分割。

參考文獻(xiàn)：

[1]Mai X C，Meng M Q H，Automatic lesion segmentation from rice leaf blast field images based on random forest[C]. RCAR，2016：255-259.

[2]Hu Q X，Tian J，He D J. Wheat leaf lesion color image segmentation with improved multichannel selection based on the ChanVese model[J]. Computers and Electronics in Agriculture，2017，135（C）：260-268.

[3]Zaitoun N M，Aqel M J. Survey on image segmentation techniques[J]. Procedia Computer Science，2015，65：797-806.

[4]Pratondo A，Ong S H，Chui C K，Region growing for medical image segmentation using a modified multiple-seed approach on a multi-core CPU computer[C]//The 15th International Conference on Biomedical Engineering，2013：112-115.

[5]Vincent L，Soille P，Watersheds in digital spaces：an efficient algorithm based on immersion simulations[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence，1991，13（6）：583-598.

[6]Zhou S P，Wang J J，Zhang S，et al. Active contour model based on local and global intensity information for medical image segmentation[J]. Neurocomputing，2016，186：107-118.

[7]Kass M，Witkin A，Terzopoulos D. Snakes：active contour models[J]. International Journal of Computer Vision，1988，1（4）：321-331.

[8]Zeng L，Chen J，Xu Y F，et al. Level set method for image segmentation based on local variance and improved intensity inhomogeneity model[J]. IET Image Processing，2016，10（12）：1007.

[9]Csaelles V，Kimmel R，Sapiro G. Geodesic active contours[C]//Proceedings of IEEE International Conference on Computer Vision，1995：694-699.

[10]Perona P，Malik J. Scale-space and edge detection using anisotropic diffusion[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence，1990，12（7）：629-639.

[11]Caselles V，Catt F，Coll T，et al. A geometric model for active contours in image processing[J]. Numerische Mathematik，1993，66（1）：1-31.

[12]Osher S，F(xiàn)edkiw R. Level set methods and dynamic implicit surfaces[M]. New York：Springer-Verlag，2002：158-174.

[13]Li C M，Xu C Y，Gui C F，et al. Level set evolution without re-initialization：a new variational formulation[C]//IEEE Conf Comput Vis Pattern Recognit，2005：430-436.

[14]Long J W，Zhang J X，Xiang N，et al. An iterative maximum entropy thresholding algorithm[C]//2016 International Conference on Cyberworlds，2016：171-174.

[15]Tian J，Hu Q X，Ma X Y. Color image segmentation of plant lesion using improved C-V model based on Gaussian distribution[J]. Transactions of the Chinese Society of Agricultural Engineering，2013，29（16）：166-173.