基于DBN與對象融合的遙感圖像變化檢測方法

竇方正,, ,,

(1.中國科學(xué)院電子學(xué)研究所,北京 100190; 2.中國科學(xué)院大學(xué),北京 100190; 3.北京跟蹤與通信技術(shù)研究所,北京 100094)

0 概述

圖像變化檢測是指通過對同一區(qū)域不同時(shí)間觀測的圖像進(jìn)行處理和分析確定目標(biāo)變化情況的過程。遙感圖像因其覆蓋地域大的特點(diǎn)被廣泛應(yīng)用于地物目標(biāo)的變化檢測。遙感圖像的變化檢測在各領(lǐng)域具有廣泛的應(yīng)用,如城市規(guī)劃、災(zāi)害預(yù)測和戰(zhàn)場情報(bào)分析等,因此,開展遙感圖像變化檢測方法的研究具有重要的應(yīng)用價(jià)值。

隨著遙感觀測技術(shù)的快速發(fā)展和成像分辨率的提高,單位地表面積的遙感圖像所承載的信息不斷增加。在高分辨率遙感影像中,目標(biāo)細(xì)節(jié)信息更加豐富,這些目標(biāo)本身的顏色和紋理變化較為復(fù)雜,且背景地物有時(shí)會呈現(xiàn)出與目標(biāo)相似的顏色紋理特征,單純利用圖像中的顏色、紋理、局部邊緣等低層特征難以取得理想的變化檢測結(jié)果。為此,科研人員開始分析和提取高級的圖像特征以改善分類器性能,以提高算法在海量高分辨率遙感圖像數(shù)據(jù)下的魯棒性。文獻(xiàn)[1]針對高斯過程分類器在變化檢測中的不足,提出了基于空間上下文相關(guān)的高斯過程變化檢測方法。文獻(xiàn)[2]提出了一種融合多種特征的變化檢測方法,該方法利用光譜直方圖、LBP等多種特征提高變化檢測精度。然而傳統(tǒng)的基于像素的方法易受到局部噪聲和復(fù)雜背景的干擾,而基于面向?qū)ο蠓治龅淖兓瘷z測能夠較好地解決此問題。文獻(xiàn)[3]利用圖像背景先驗(yàn)知識提出一種基于對象的分析方法,通過圖像變化檢測結(jié)果實(shí)現(xiàn)目標(biāo)的檢測與分析。文獻(xiàn)[4]采用面向?qū)ο蟮乃枷霕?gòu)造對象的特征向量用于SVM的變化分類。文獻(xiàn)[5]提出一種基于對數(shù)變差函數(shù)紋理增強(qiáng)的圖像變化檢測方法,采用方向?qū)?shù)變差函數(shù)算法對圖像的紋理特征進(jìn)行增強(qiáng)。文獻(xiàn)[6]在面向?qū)ο蟮幕A(chǔ)上提出一種基于直推式支持向量機(jī)的變化檢測方法,該方法能夠利用未標(biāo)注的歷史數(shù)據(jù)進(jìn)行支持向量機(jī)訓(xùn)練,從而提高檢測方法的泛化能力。然而,由于圖像中的對象通常具有大小不一的尺寸和形狀,因此對對象進(jìn)行特征設(shè)計(jì)和提取面臨很大挑戰(zhàn),而現(xiàn)有的面向?qū)ο蟮姆治龇椒ǘ鄶?shù)僅能利用一些簡單特征組合來得到不同尺寸對象的特征。

針對上述問題,本文提出一種基于深度置信網(wǎng)絡(luò)(Deep Belief Network,DBN)和面向?qū)ο蠼Y(jié)果融合的圖像變化檢測方法?？紤]到DBN能夠?qū)W習(xí)樣本數(shù)據(jù)集中低層到高層本質(zhì)特征,因此,將圖像變化檢測看作二分類問題,利用DBN模型進(jìn)行分類。為降低圖像和目標(biāo)局部噪聲的對檢測結(jié)果的影響,設(shè)計(jì)面向?qū)ο蟮姆诸惾诤戏椒?以圖像像素為單元進(jìn)行圖像特征提取和分類,將DBN分類結(jié)果作為該像素的變化檢測結(jié)果進(jìn)行面向?qū)ο蟮姆诸惤Y(jié)果融合。此外,提出一種多尺度的圖像特征學(xué)習(xí)和分類方法,在特征學(xué)習(xí)和分類階段利用圖像的金字塔變換得到多個(gè)尺度上的數(shù)據(jù),并利用DBN學(xué)習(xí)得到多尺度的目標(biāo)和上下文特征,從而增強(qiáng)本文方法對環(huán)境變化因素的魯棒性,提高檢測準(zhǔn)確率。

1 深度置信網(wǎng)絡(luò)

隨著深度學(xué)習(xí)方法的研究與發(fā)展以及該方法在自然場景圖像解譯、大數(shù)據(jù)分析[7]等一系列應(yīng)用中取得的成功,使得深度學(xué)習(xí)逐步引起遙感圖像解譯研究人員的重視,并開始研究如何將深度學(xué)習(xí)算法引入到遙感圖像的特征提取和分析中。文獻(xiàn)[8]利用深度置信網(wǎng)絡(luò)對遙感圖像中的飛機(jī)目標(biāo)進(jìn)行檢測,文獻(xiàn)[9]利用棧式自編碼器實(shí)現(xiàn)了SAR圖像中目標(biāo)的變化檢測。本文將深度學(xué)習(xí)引入遙感圖像特征提取和分類中,提出一種基于深度置信網(wǎng)絡(luò)的光學(xué)遙感圖像變化檢測方法。深度置信網(wǎng)絡(luò)(DBN)[10]是一種產(chǎn)生式的深度神經(jīng)網(wǎng)絡(luò)模型,其能利用大量的未標(biāo)注歷史圖像數(shù)據(jù)和一種逐層的非監(jiān)督訓(xùn)練方法實(shí)現(xiàn)特征自學(xué)習(xí)和提取,并通過少量標(biāo)注數(shù)據(jù)進(jìn)行模型的分類訓(xùn)練和優(yōu)化。DBN的基本結(jié)構(gòu)是一個(gè)多層的神經(jīng)網(wǎng)絡(luò),其概率圖模型如圖1所示。其中,最上面一層網(wǎng)絡(luò)是無向圖模型,即限制玻爾茲曼機(jī)模型,其余各層是自上而下的有向連接圖模型。基于這種網(wǎng)絡(luò)結(jié)構(gòu),文獻(xiàn)[10]提出一種通過逐層訓(xùn)練限制玻爾茲曼機(jī)進(jìn)行參數(shù)初始化的訓(xùn)練方法,用于深度置信網(wǎng)絡(luò)特征學(xué)習(xí)。

圖1 深度置信網(wǎng)絡(luò)圖模型

限制玻爾茲曼機(jī)的圖模型如圖2所示。它是一個(gè)包含2層節(jié)點(diǎn)的無向圖模型,即可視層與隱藏層,分別以v與h表示。不同層節(jié)點(diǎn)的連接由權(quán)值矩陣W表示。

圖2 限制玻爾茲曼機(jī)圖模型

限制玻爾茲曼機(jī)的能量函數(shù)表示如下:

E(v,h)=-hTWv-aTv-bTh

(1)

v=(v1,v1,…,vNv)T表示可視層的狀態(tài)向量,h=(h1,h1,…,hNh)T表示隱藏層的狀態(tài)向量,其中vi和hj分別表示第i個(gè)可視層節(jié)點(diǎn)與第j個(gè)隱藏層節(jié)點(diǎn)的狀態(tài),wij表示權(quán)值矩陣W中連接這兩個(gè)節(jié)點(diǎn)的權(quán)值,ai和bj分別表示可視層偏置項(xiàng)a中第i個(gè)元素與隱藏層偏置項(xiàng)b中的第j個(gè)元素。W、a和b統(tǒng)稱為限制玻爾茲曼機(jī)的網(wǎng)絡(luò)參數(shù)。利用式(1)定義的能量函數(shù),可以得到限制玻爾茲曼機(jī)中可視層節(jié)點(diǎn)與隱藏層節(jié)點(diǎn)的聯(lián)合概率分布:

(2)

(3)

在限制玻爾茲曼機(jī)的參數(shù)求解中,為了計(jì)算的方便,通常會對對數(shù)形式的似然函數(shù)進(jìn)行梯度計(jì)算,即:

(4)

其中,<·>data和<·>model分別表示模型關(guān)于數(shù)據(jù)經(jīng)驗(yàn)分布和模型真實(shí)分布的期望。為解決吉布斯采樣效率較低的問題,文獻(xiàn)[11]提出一種近似采樣算法,稱為對比散度(Contrastive Divergence,CD)算法。目前,對比散度算法已成為限制玻爾茲曼機(jī)的標(biāo)準(zhǔn)訓(xùn)練算法。

2 本文方法

本文方法流程如圖3所示,主要分為模型訓(xùn)練與變化檢測2個(gè)主要步驟。在第1個(gè)步驟中,經(jīng)過預(yù)處理的訓(xùn)練圖像對(待檢測圖像和基準(zhǔn)圖像)將用來訓(xùn)練一個(gè)深度置信網(wǎng)絡(luò),經(jīng)過限制玻爾茲曼機(jī)非監(jiān)督訓(xùn)練和監(jiān)督微調(diào)訓(xùn)練2個(gè)過程就得到了用于變化檢測的深度神經(jīng)網(wǎng)絡(luò)模型。在變化檢測步驟中,對經(jīng)過預(yù)處理的測試圖像將利用上一步驟中訓(xùn)練得到的深度神經(jīng)網(wǎng)絡(luò)模型進(jìn)行基于像素的分類和變化檢測,同時(shí)本文對預(yù)處理后的待測圖像進(jìn)行基于對象的分割,隨后進(jìn)行基于對象分割的分類結(jié)果融合,以降低局部噪聲的影響,經(jīng)過融合后的結(jié)果即最終的變化檢測結(jié)果。

圖3 基于深度置信網(wǎng)絡(luò)的遙感圖像變化檢測流程

2.1 圖像預(yù)處理

預(yù)處理的過程主要包括圖像配準(zhǔn)、相對輻射校正和濾波去噪。圖像配準(zhǔn)的目的是消除成像設(shè)備、拍攝位置、角度不同所帶來的成像差異,本文采用了基于SIFT[12]特征點(diǎn)匹配的圖像配準(zhǔn)方法。相對輻射校正的目的是消除2幅或多幅圖像由于光照、天氣等變化帶來的成像明暗差異,這些差異會對變化檢測的結(jié)果造成較大影響,本文中采用直方圖匹配的方法實(shí)現(xiàn)相對輻射校正。濾波去噪的目的是消除圖像中局部噪聲的影響,本文采用了保邊濾波中的雙邊濾波[13]實(shí)現(xiàn)噪聲去除,與傳統(tǒng)的高斯濾波等方法不同,雙邊濾波在去噪的同時(shí)能夠很好地保留圖像中的邊緣信息,這些信息對于此后基于對象的分析和結(jié)果融合非常重要。經(jīng)過圖像預(yù)處理之后,待檢測圖像和基準(zhǔn)圖像組成的圖像對將用于深度置信網(wǎng)絡(luò)模型訓(xùn)練和最終變化檢測結(jié)果的產(chǎn)生。

2.2 深度置信網(wǎng)絡(luò)模型訓(xùn)練

本文方法將經(jīng)過預(yù)處理后的訓(xùn)練圖像對作為訓(xùn)練數(shù)據(jù),來進(jìn)行深度置信網(wǎng)絡(luò)的參數(shù)求解,從而得到一個(gè)用于分類的深度神經(jīng)網(wǎng)絡(luò)模型。其中,基于深度置信網(wǎng)絡(luò)模型的分類變化檢測按像素進(jìn)行分類,即圖像中的每個(gè)像素被賦予一個(gè)標(biāo)簽,標(biāo)簽值用離散的0和1表示,標(biāo)簽為1表示圖像對中該像素位置產(chǎn)生了變化。為了充分利用遙感圖像目標(biāo)的多尺度上下文信息,本文提出了一種基于圖像金字塔變換的訓(xùn)練方法,即一個(gè)像素通過結(jié)合其周圍一定范圍內(nèi)的多尺度上下文信息來得到分類結(jié)果。如圖4所示,本文獲取以待分類像素為中心周圍一定范圍內(nèi)的金字塔變換后的圖像(黑框范圍內(nèi)圖像)作為深度置信網(wǎng)絡(luò)的輸入數(shù)據(jù)。

圖4 利用多尺度圖像金字塔訓(xùn)練深度置信網(wǎng)絡(luò)模型示意圖

深度置信網(wǎng)絡(luò)的訓(xùn)練包括逐層非監(jiān)督預(yù)訓(xùn)練和監(jiān)督微調(diào)訓(xùn)練,逐層非監(jiān)督預(yù)訓(xùn)練是特征學(xué)習(xí)和提取的過程,訓(xùn)練利用限制玻爾茲曼機(jī)和對比散度算法實(shí)現(xiàn)。預(yù)訓(xùn)練完成后,本文利用反向傳播算法[14]進(jìn)行網(wǎng)絡(luò)參數(shù)的有監(jiān)督微調(diào)訓(xùn)練,網(wǎng)絡(luò)的輸出層是一個(gè)二維向量,二維輸出向量中的最大值位置表示了檢測分類結(jié)果。本文利用事先標(biāo)注的標(biāo)簽進(jìn)行監(jiān)督訓(xùn)練,最后得到一個(gè)深度置信網(wǎng)絡(luò)模型,并將其用于變化檢測分類中。

2.3 變化檢測和面向?qū)ο蟮慕Y(jié)果融合

在變化檢測階段,本文方法首先將待檢測的圖像對進(jìn)行預(yù)處理,按2.2節(jié)所述的方法獲取多尺度的深度置信網(wǎng)絡(luò)輸入數(shù)據(jù),然后將數(shù)據(jù)輸入到訓(xùn)練好的深度置信網(wǎng)絡(luò)模型中,模型的輸出即為變化檢測結(jié)果。

在高分辨率遙感圖像中,由于光照、天氣等變化因素,在同一位置不包含變化的2幅圖像中會出現(xiàn)局部紋理、明暗的變化,這些變化會使得最后的檢測結(jié)果產(chǎn)生大量虛警。為了解決這一問題,本文在基于像素的深度置信網(wǎng)絡(luò)變化檢測基礎(chǔ)上,提出一種面向?qū)ο蟮慕Y(jié)果融合方法。對于經(jīng)過預(yù)處理的測試圖像對,首先利用訓(xùn)練好的深度置信網(wǎng)絡(luò)模型進(jìn)行基于像素的變化檢測,同時(shí)對測試圖像對利用圖像分割方法進(jìn)行對象分割,然后利用對象分割結(jié)果對深度置信網(wǎng)絡(luò)的分類結(jié)果進(jìn)行融合,得到最終的變化檢測結(jié)果,具體過程如圖5所示。

圖5 面向?qū)ο蟮淖兓瘷z測結(jié)果融合過程

常用的對象分割方法包括基于均值漂移的方法[15]、水平集分割[16]、分水嶺算法[17]、金字塔算法[18]等。本文采用基于均值漂移的分割方法獲取對象。均值漂移是一種非參數(shù)估計(jì)密度函數(shù)方法[19],它使每一個(gè)點(diǎn)通過有效的統(tǒng)計(jì)迭代“漂移”到概率密度函數(shù)的局部極大值點(diǎn)上。均值漂移算法計(jì)算簡單,且不需要任何先驗(yàn)知識,在實(shí)際應(yīng)用中通常具有較好的穩(wěn)定性和較高的效率。由于待檢測圖像和基準(zhǔn)圖的分割結(jié)果往往會有差別,因此,本利用文獻(xiàn)[3]中的方法進(jìn)行分割結(jié)果進(jìn)行合并。

進(jìn)行圖像分割后,本文方法利用分割得到的對象對深度置信網(wǎng)絡(luò)分類得到的變化檢測結(jié)果進(jìn)行融合,其原理如下:

(5)

其中,label表示對象Ω中某個(gè)像素經(jīng)過結(jié)果融合后的標(biāo)簽,pi表示對象Ω中像素i的深度置信網(wǎng)絡(luò)預(yù)測標(biāo)簽,其取值為1或0,Area表示對象Ω的像素面積分割后,T是本文方法定義的自適應(yīng)融合閾值。

T的計(jì)算公式如下:

(6)

其中,MaxArea表示圖像中最大面積對象的像素面積值。通過面向?qū)ο蟮娜诤?可以很好地解決由于對象尺寸大小不一和形狀變化而難以進(jìn)行高層特征設(shè)計(jì)和提取的難題。同時(shí),該方法利用對象分割結(jié)果,對基于像素的深度置信網(wǎng)絡(luò)分類結(jié)果進(jìn)行校正,以抑制局部噪聲對檢測結(jié)果的影響,從而提高檢測精度。

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

3.1 實(shí)驗(yàn)設(shè)置與數(shù)據(jù)集描述

為驗(yàn)證本文方法的有效性,選取20組QucikBird光學(xué)遙感圖像對組成實(shí)驗(yàn)數(shù)據(jù)集。其中,圖像分辨率為0.61 m,尺寸在1 000×1 000像素與10 000×10 000像素之間。這些圖像中包含全球部分民用機(jī)場、港口、工廠等大型設(shè)施。圖6為其中一組測試數(shù)據(jù)。

圖6 測試數(shù)據(jù)

在深度置信網(wǎng)絡(luò)訓(xùn)練中,本文利用2.2節(jié)所述方法在5組圖像對中隨機(jī)產(chǎn)生20 000組多尺度數(shù)據(jù),作為深度置信網(wǎng)絡(luò)的訓(xùn)練數(shù)據(jù)。每組訓(xùn)練數(shù)據(jù)包含待測圖和基準(zhǔn)圖的各3層金字塔圖像,每層圖像選擇16×16大小的切片尺寸。其余15組圖像對作為變化檢測方法的測試數(shù)據(jù)集。

為了評價(jià)算法的性能,本文選擇準(zhǔn)確率、召回率和F1值作為定量評價(jià)指標(biāo)。準(zhǔn)確率即檢測結(jié)果為變化的像素里實(shí)際為變化的比例,召回率是檢測結(jié)果中檢測出為變化的像素本占所有實(shí)際為變化像素的比例。F1值的計(jì)算方法為:

(7)

其中,P為準(zhǔn)確率,R為召回率。

3.2 變化檢測結(jié)果比較與分析

通過觀察F1值得變化,本文對深度置信網(wǎng)絡(luò)的參數(shù)進(jìn)行調(diào)整,以獲取最優(yōu)的變化檢測結(jié)果,這些參數(shù)主要包括隱藏層數(shù)量以及每層的節(jié)點(diǎn)數(shù)量。通過實(shí)驗(yàn)發(fā)現(xiàn),當(dāng)選擇節(jié)點(diǎn)數(shù)量分別為1 536(輸入層,16×16×6=1 536)、300、200、100、300、2(分類層)的包含4個(gè)隱藏層的深度置信網(wǎng)絡(luò)時(shí),可以取得最優(yōu)的檢測結(jié)果。因此,下面的對比實(shí)驗(yàn)中均選擇這種參數(shù)組合。

實(shí)驗(yàn)1將本文方法與利用深度置信網(wǎng)絡(luò)在不加入多尺度上下文信息情況下僅加入本文對象融合策略的方法(方法1)和像素分類的方法(方法2)進(jìn)行比較。圖7顯示了不同方法檢測結(jié)果的定性比較,其中圖7(c)表示僅利用深度置信網(wǎng)絡(luò)按像素分類得到的檢測結(jié)果,可以發(fā)現(xiàn)結(jié)果受局部噪聲影響明顯,與之相比,利用本文提出的融合方法進(jìn)行對象融合后(如圖7(d)所示),能夠顯著降低局部噪聲的影響。圖7(e)展示了本文方法的檢測結(jié)果,與圖7(d)的結(jié)果相比可以看出,通過在訓(xùn)練和檢測中引入多尺度的上下文信息,可以很好地抑制由光照、紋理等因素變化而產(chǎn)生的干擾,提高檢測準(zhǔn)確率。

圖7 不同方法的檢測結(jié)果比較

對以上方法的定量統(tǒng)計(jì)結(jié)果進(jìn)行比較,如表1所示。實(shí)驗(yàn)結(jié)果表明,本文方法可使檢測準(zhǔn)確率和召回率得到明顯提升。

表1 檢測準(zhǔn)確率、召回率和F1值的比較1 %

本文方法與文獻(xiàn)[1]方法和文獻(xiàn)[6]方法的檢測結(jié)果比較如表2所示。實(shí)驗(yàn)結(jié)果表明,與其他2種方法相比,本文方法在準(zhǔn)確率和召回率上均有提高。2種對比方法雖然較好地利用了面向?qū)ο蟮姆治龇椒▉硪种圃肼暤挠绊?但由于對象尺寸大小不一的問題,僅能對對象進(jìn)行了簡單特征的提取,因此限制了方法的進(jìn)一步提高。而本文方法在解決上述問題的基礎(chǔ)上通過利用多尺度上下文信息,顯著地提升了檢測性能。

表2 檢測準(zhǔn)確率、召回率與F1值的比較2 %

4 結(jié)束語

本文提出一種基于深度置信網(wǎng)絡(luò)和對象融合的圖像變化檢測方法。將變化檢測問題轉(zhuǎn)化為二分類問題,并在分類問題中把圖像像素作為分類單元。在特征學(xué)習(xí)和分類階段充分利用圖像目標(biāo)的上下文信息,設(shè)計(jì)了多尺度的圖像特征學(xué)習(xí)和分類方法。在此基礎(chǔ)上,本文還提出一種基于對象的分類融合方法,通過基于對象融合的方法處理深度置信網(wǎng)絡(luò)分類得到的結(jié)果,從而降低圖像和目標(biāo)局部噪聲的影響。在QucikBird影像數(shù)據(jù)集上的實(shí)驗(yàn)結(jié)果表明,本文方法具有較高的變化檢測準(zhǔn)確率。下一步將在本文工作的基礎(chǔ)上研究深度學(xué)習(xí)以及不同形式的融合方法在遙感圖像變化檢測中的應(yīng)用。

[1] 陳克明,周志鑫,盧漢清,等.基于高斯過程的高分辨率遙感圖像變化檢測[J].遙感學(xué)報(bào),2012,16(6):82-94.

[2] 李 亮,舒 寧,王 凱,等.融合多特征的遙感影像變化檢測方法[J].測繪學(xué)報(bào),2014,43(9):945-953.

[3] 陳麗勇,閆夢龍,孫 顯,等.一種基于背景先驗(yàn)的飛機(jī)目標(biāo)檢測方法[J].測繪科學(xué),2016,41(3):69-74.

[4] HUO C,ZHOU Z,LU H,et al.Fast object-level change detection for VHR images[J].IEEE Geoscience and Remote Sensing Letters,2010,7(1):118-122.

[5] 王 威,劉 洋,劉 婧,等.基于對數(shù)變差函數(shù)紋理增強(qiáng)的圖像變化檢測[J].計(jì)算機(jī)工程,2016,42(2):224-228.

[6] 陳麗勇,孫 顯,王宏琦.一種基于TSVM-MRF的變化檢測方法[J].國外電子測量技術(shù),2015,3(7):32-36.

[7] BENGIO Y.Learning deep architectures for AI[J].Foundations and Trends in Machine Learning,2009,2(1):1-127.

[8] CHEN X,XIANG S,LIU C L,et al.Aircraft detection by deep belief nets[C]//Proceedings of the 2nd Asian Conference on Pattern Recognition.Washington D.C.,USA:IEEE Press,2013:54-58.

[9] LIU J,GONG M,ZHAO J,et al.Difference representation learning using stacked restricted Boltzmann machines for change detection in SAR images[J].Soft Computing,2016,20(12):4645-4657.

[10] HINTON G E,OSINDERO S,TEH Y W.A fast learning algorithm for deep belief nets[J].Neural Computation,2006,18(7):1527-1554.

[11] HINTON G E.Training products of experts by minimizing contrastive divergence[J].Neural Computation,2002,14(8):1771-1800.

[12] LOWE D G.Object recognition from local scale-invariant features[C]//Proceedings of the 7th IEEE International Conference on Computer Vision.Washington D.C.,USA:IEEE Press,1999:1150-1157.

[13] TOMASI C,MANDUCHI R.Bilateral filtering for gray and color images[C]//Proceedings of the 6th International Conference on Computer Vision.Washington D.C.,USA:IEEE Press,1998:839-846.

[14] LeCUN Y,BOSER B,DENKER J S,et al.Back-propagation applied to handwritten zip code recognition[J].Neural Computation,1989,1(4):541-551.

[15] COMANICIU D,MEER P.Mean shift:a robust approach toward feature space analysis[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,2002,24(5):603-619.

[16] OSHER S,SETHIAN J A.Fronts propagating with curvature-dependent speed:algorithms based on Hamilton-Jacobi formulations[J].Journal of Compu-tational Physics,1988,79(1):12-49.

[17] GAUCH J M.Image segmentation and analysis via multiscale gradient watershed hierarchies[J].IEEE Transactions on Image Processing,1999,8(1):69-79.

[18] REZAEE M R,van der ZWET P M,LELIEVELDT B P,et al.A multiresolution image segmentation technique based on pyramidal segmentation and fuzzy clustering[J].IEEE Transactions on Image Processing,2000,9(7):1238-1248.

[19] KEINOSUKE F,HOSTETLER L.The estimation of the gradient of a density function,with applications in pattern recognition[J].IEEE Transactions on Information Theory,1975,21(1):32-40.