周家香，朱建軍，馬慧云，梅小明

(中南大學(xué) 地球科學(xué)與信息物理學(xué)院，湖南 長(zhǎng)沙，410083)

國(guó)內(nèi)外在圖像分割方面已取得了很多成果，現(xiàn)有自然圖像分割方法可概括為以下幾個(gè)方面：基于全局信息或局部信息；基于區(qū)域或邊緣；針對(duì)灰度影像或紋理影像[1]。由于遙感影像所記錄地物的復(fù)雜性和多樣性，數(shù)據(jù)量巨大，將這些分割方法引入到遙感影像分割中面臨很大困難[2?3]。Mean-Shift (MS)算法是一種高效的統(tǒng)計(jì)迭代算法，依靠特征空間樣本點(diǎn)的統(tǒng)計(jì)特征，不需要任何先驗(yàn)知識(shí)，近年來(lái)被廣泛地應(yīng)用于圖像分割[4?7]和跟蹤[8]等計(jì)算機(jī)視覺(jué)領(lǐng)域。傳統(tǒng)MS算法的計(jì)算復(fù)雜度是O(kdN2) (其中k為每個(gè)數(shù)據(jù)點(diǎn)的平均迭代次數(shù)，d為數(shù)據(jù)維數(shù)，N為1幅圖像的像素)[9?12]，如IKONOS和QuickBird等遙感影像，需要多次迭代才能收斂。同時(shí)，若所用的特征向量的維數(shù)多，則收斂速度比較慢，影響遙感圖像的實(shí)際應(yīng)用。近年來(lái)，一些研究者研究了MS算法的原理，提出減少迭代次數(shù)和降低每次迭代的計(jì)算量這2種提高算法效率的方法[10]。為了提高M(jìn)S算法對(duì)遙感影像分割的速度和效果，本文提出以下改進(jìn)措施：

(1) 對(duì)不同空間分辨率遙感圖像給定空間帶寬參考值，不同波段采用不同的灰度值域帶寬。基于傳統(tǒng)MS算法編寫(xiě)的EDISON軟件需要人工輸入空間帶寬和值域帶寬，且所有波段采用相同的值域帶寬。若將遙感影像分割成有意義的區(qū)域且沒(méi)有過(guò)多的細(xì)節(jié)，自適應(yīng)帶寬才能得到好的分割效果。

(2) 面對(duì)大數(shù)據(jù)量的遙感影像分割速度慢的問(wèn)題，本文主要在MS運(yùn)行效率上采取多項(xiàng)加速措施。

(3) 針對(duì) MS分割存在過(guò)分割現(xiàn)象，將檢測(cè)出的模點(diǎn)采用基于全局模點(diǎn)融合來(lái)得到穩(wěn)定的分割結(jié)果。

(4) 可用于遙感影像分割的特征很多，模點(diǎn)檢測(cè)或模點(diǎn)融合時(shí)利用全部特征勢(shì)必增加計(jì)算時(shí)間，改進(jìn)方法將可將特征分組，利用位置和顏色特征進(jìn)行模點(diǎn)檢測(cè)，利用紋理特征進(jìn)行模點(diǎn)融合。

1 改進(jìn)MS遙感影像分割算法

1.1 傳統(tǒng)高斯MS算法模點(diǎn)檢測(cè)

傳統(tǒng)高斯核MS由高斯核密度估計(jì)的梯度推導(dǎo)獲得[4]。給定d維歐拉空間R中n個(gè)采樣點(diǎn){xi，i=1，…，n}，利用高斯核函數(shù)K (x)= ( 2π)?d/2exp(?x2)，則MS迭代向量為：

式中：xs為二維位置特征；xr為多維值域特征；hs為空間帶寬；hr為值域帶寬。則傳統(tǒng)MS計(jì)算步驟：

(1) 選用高斯核函數(shù)，帶寬hs和 hr，隨機(jī)給定樣本點(diǎn)x；

(2) 根據(jù)式( 1) 計(jì)算樣本點(diǎn)的迭代向量m(x)；

開(kāi)展拓展實(shí)驗(yàn)是培養(yǎng)學(xué)生科學(xué)探究能力的有效補(bǔ)充途徑。目前高中生物教材中的很多實(shí)驗(yàn)為體驗(yàn)或驗(yàn)證性實(shí)驗(yàn)，即使是探究性實(shí)驗(yàn)，也列出了詳細(xì)的探究過(guò)程與具體內(nèi)容，學(xué)生往往在未經(jīng)充分思維后就直接參考了教材。因此，教材實(shí)驗(yàn)在培養(yǎng)學(xué)生自主探究能力上還稍顯不足，不能滿足熱愛(ài)生物學(xué)科的學(xué)生的求知需求和探究愿望，而拓展實(shí)驗(yàn)可以作為有效補(bǔ)充。為此，在教學(xué)實(shí)踐中，筆者結(jié)合校本課進(jìn)行了拓展實(shí)驗(yàn)教學(xué)，積極探索有效開(kāi)展拓展實(shí)驗(yàn)促進(jìn)學(xué)生科學(xué)探究能力提升的實(shí)施策略。

(3) 條件判斷，若|m(x)?x|＞ε(ε為閾值，如 0.001)，則移動(dòng)點(diǎn)x到m(x)，執(zhí)行步驟(2)，否則，結(jié)束循環(huán)；

(4) 將結(jié)束循環(huán)后的 m(x)賦給 x，此時(shí) m(x)稱為樣本點(diǎn)x的模點(diǎn)；

(5) 移動(dòng)到下一樣本點(diǎn)x，重復(fù)步驟(2)～(4)，直至檢測(cè)出影像內(nèi)每個(gè)數(shù)據(jù)點(diǎn)的模點(diǎn)為止；

(6) 將收斂至相同模點(diǎn)的數(shù)據(jù)點(diǎn)歸為同一區(qū)域，得到初分割結(jié)果。

這樣，整個(gè)影像被分解成m個(gè)子區(qū)域Sj(j=1，2，…，m)，每個(gè)子區(qū)域像素為nj，中心為cj，則：

1.2 全局模點(diǎn)融合

由于MS分割存在過(guò)分割問(wèn)題，經(jīng)過(guò)MS分割后，將以點(diǎn)陣表示的圖像轉(zhuǎn)化為塊狀表示。然后，采用如下過(guò)程的基于區(qū)域的MS迭代。

對(duì)每一個(gè)子區(qū)域Sj(j=1，2，…，m)，進(jìn)行以下迭代：

式中：np為區(qū)域p的像素；h為帶寬；k為迭代次數(shù)；= c 。用上式進(jìn)行迭代計(jì)算，直至j 為止(εw是閾值，如0.001)，將最終收斂值記錄為。

1.3 改進(jìn)的MS遙感影像分割算法實(shí)現(xiàn)

1.3.1 濾波

圖像噪聲會(huì)使檢測(cè)出的局部模點(diǎn)數(shù)過(guò)多，這樣使后面的模點(diǎn)融合變得困難。濾波能降低噪聲對(duì)模點(diǎn)檢測(cè)的影響。事實(shí)上，MS本身對(duì)噪聲不太敏感[11]，模點(diǎn)檢測(cè)前的濾波能很好地降低噪聲的影響，這樣能使模點(diǎn)檢測(cè)過(guò)程更簡(jiǎn)單。本文選擇3×3的中值濾波。

1.3.2 帶寬的選擇

在MS算法中最重要的參數(shù)是帶寬矩陣H，它決定迭代速度和分割效果。目前，Plug-in 和Cross-validation是廣泛使用帶寬選擇方法[1]。Plug-in規(guī)則是用估計(jì)參數(shù)來(lái)代替未知參數(shù)，漸進(jìn)積分均方誤差是應(yīng)用于Plug-in規(guī)則的通用標(biāo)準(zhǔn)。使?jié)u進(jìn)積分均方誤差最小便可估計(jì)全局值域帶寬：

式中：d為特征向量維數(shù)；j為單個(gè)波段；n為數(shù)據(jù)量；σj為特征向量的標(biāo)準(zhǔn)差。

在plug-in準(zhǔn)則下獲得的值域帶寬hr=(hR，hG，hB)能盡可能地降低估計(jì)密度函數(shù)與原始密度函數(shù)之間的變形。為了將遙感影像分割成有意義的區(qū)域而不含有過(guò)多的細(xì)節(jié)，一般來(lái)說(shuō)，過(guò)小的hr是無(wú)意義的。但是，過(guò)大的值域帶寬導(dǎo)致忽略更多的影像細(xì)節(jié)，同時(shí)降低檢測(cè)出的模點(diǎn)數(shù)；相反，小的灰度帶寬使更多的地物細(xì)節(jié)被提取出來(lái)。因此，灰度值域帶寬hr被認(rèn)為是分割的尺度。

1.3.3 快速高斯模點(diǎn)檢測(cè)

應(yīng)用傳統(tǒng)MS算法進(jìn)行聚類，可得到初始分割結(jié)果。初始分割是對(duì)影像中每一個(gè)像點(diǎn)進(jìn)行MS迭代得到其模點(diǎn)，由于遙感影像數(shù)據(jù)量大，迭代計(jì)算量非常大，本文采用以下加速策略。

(1) 以待處理點(diǎn) xs為中心，事先根據(jù)圖像的空間分辨率來(lái)選擇空間帶寬hs，以長(zhǎng)、寬均為2hs+1的正方形按下式構(gòu)造高斯空間模板；

式中：xs表示待處理點(diǎn)的二維坐標(biāo)，xsi為以xs為中心，長(zhǎng)和寬為2hs+1正方形內(nèi)的采樣點(diǎn)二維坐標(biāo)?？梢钥醋鳈?quán)值，與xs點(diǎn)近的采樣點(diǎn)權(quán)值大，遠(yuǎn)的采樣點(diǎn)權(quán)值小。則式(1)變成：

(2) 為了減少M(fèi)S迭代次數(shù)，從xn點(diǎn)開(kāi)始MS迭代，第 τ次迭代后的向量為等于或近似于另一特征向量xm，則認(rèn)為xn和xm收斂于同一個(gè)模點(diǎn)，停止迭代。

(3) 對(duì)任意一隨機(jī)樣本點(diǎn)x，沿一定的路徑收斂至一固定點(diǎn)(即模點(diǎn))，將收斂路徑一定寬度內(nèi)的點(diǎn)都認(rèn)為收斂至相同的模點(diǎn)，這些點(diǎn)不再參加迭代計(jì)算。

(4) MS迭代計(jì)算量大的1個(gè)主要原因是特征向量的維數(shù)，而遙感圖像可用的特征向量的維數(shù)很多，本文將特征向量分為2組，首次迭代時(shí)使用多光譜灰度特征，模點(diǎn)融合時(shí)使用紋理特征。

1.3.4 基于全局的模點(diǎn)融合

由于傳統(tǒng)MS分割一般存在過(guò)分割現(xiàn)象，進(jìn)一步將分割后的影像進(jìn)行全局模點(diǎn)融合，得到最終分割結(jié)果。遙感影像的特征有位置、灰度、紋理、形狀等，為了達(dá)到理想的分割結(jié)果，應(yīng)盡可能地利用遙感影像的多維特征[9]，這樣勢(shì)必會(huì)增加計(jì)算時(shí)間。為了充分利用遙感影像的多維特征又不增加計(jì)算時(shí)間，在模點(diǎn)融合時(shí)使用紋理特征。

此時(shí)，用到的紋理特征f1，f2和f3分別為：

式中：i和j為像素點(diǎn)灰度值；p(i,j)為灰度值i和j出現(xiàn)的歸一化概率；L為灰度級(jí)。則式(3)中h應(yīng)為紋理帶寬ht。對(duì)任一分割子區(qū)域Sj(j=1，2，…，m)，其紋理帶寬由下式計(jì)算：

式中：nj為第 j塊的像素，fki(k=1，2，3)為該塊內(nèi)像素i的 fk紋理特征；ckj(k=1，2，3)為第j塊的fk紋理均值。

2 實(shí)驗(yàn)結(jié)果評(píng)價(jià)與分析

2.1 分割實(shí)驗(yàn)

在改進(jìn)的分割方法中，空間帶寬hs主要影響計(jì)算速度。不同傳感器的影像有不同的空間分辨率。圖1所示為 SPOT5影像采用不同空間帶寬分割結(jié)果的對(duì)比。實(shí)驗(yàn)發(fā)現(xiàn)，空間分辨率為10.00 m的SPOT5影像選用hs=5比較合適，空間分辯率為2.44 m的Quickbird影像選用hs=9比較合適。

圖2(a)是1幅Quickbird 影像270×270的子圖像，采用前面提出的參考空間帶寬 hs=9和由式(4)計(jì)算的灰度值域帶寬進(jìn)行快速高斯MS迭代得到圖2(b)所示的分割結(jié)果(290個(gè)模點(diǎn))。圖2(c)所示是用EDSION軟件(1個(gè)基于傳統(tǒng)MS算法編寫(xiě)的軟件)分割結(jié)果。將本文的分割結(jié)果與EDISON軟件分割結(jié)果進(jìn)行對(duì)比，由圖2(c)可見(jiàn)：本文分割結(jié)果要優(yōu)于EDISON軟件分割結(jié)果。

圖1 SPOT5影像不同空間帶寬對(duì)比實(shí)驗(yàn)Fig.1 Contrast experiment of SPOT5 image with various space bandwidths

圖2 改進(jìn)分割方法與EDISON分割結(jié)果對(duì)比Fig.2 Comparative segmentation between presented approach and EDISON of Quickbird image

從圖2(b)可以看出：本文的分割結(jié)果比較理想，但與分割參照?qǐng)D(圖3(b))比較還存在過(guò)分割現(xiàn)象。對(duì)初分割結(jié)果進(jìn)行全局模點(diǎn)融合，保持原有空間帶寬不變，采用式(11)計(jì)算每塊的紋理帶寬，得到圖3(a)所示的實(shí)驗(yàn)結(jié)果。由圖3(c)可見(jiàn)：本文分割結(jié)果比較理想。

2.2 本文改進(jìn)方法與傳統(tǒng)MS算法的時(shí)間比較

通過(guò)實(shí)驗(yàn)比較傳統(tǒng)MS方法與本文MS改進(jìn)方法的算法時(shí)間。采用圖2(a)所示Quickbird圖像，計(jì)算時(shí)間見(jiàn)表1，從表1可以看出：本文所采用加速策略是有效的。

2.3 分割效果討論

為了評(píng)價(jià)分割的效果，本文對(duì)圖2(a)所示的Quickbird影像進(jìn)行定量對(duì)比實(shí)驗(yàn)。Polak等[13]提出了一種基于對(duì)象一致性誤差(EOC)的客觀有效的圖像分割評(píng)判標(biāo)準(zhǔn)。EOC(Ig，Is)∈[0，1](其中Ig為地面參考對(duì)象，Is為實(shí)際分割對(duì)象)。當(dāng)EOC=0時(shí)，表示分割完全與地面參考對(duì)象一致；EOC越小，表示分割效果越好。實(shí)驗(yàn)所用圖2(a)的研究范圍較小，圖像空間分辨率較高，且各類地物相對(duì)集中，通過(guò)目視判讀得到實(shí)驗(yàn)參考圖(圖3(b))。將地物分為房屋、林地、水面、農(nóng)田、綠地以及陰影6類。表2所示為本文方法與EDISON軟件分割的EOC，可以看出：采用本文方法所得房屋、林地的分割效果明顯優(yōu)于EDISON的分割效果，水面、農(nóng)田、綠地和陰影4類地物的分割效果相當(dāng)。

圖3 模點(diǎn)融合后分割結(jié)果Fig.3 Segmentation results after mode merging

表1 傳統(tǒng)MS與改進(jìn)MS計(jì)算時(shí)間比較Table 1 Comparative time of traditional MS and improved MS s

表2 改進(jìn)分割方法與EDISON分割方法的Eoc對(duì)比Table 2 Comparison of Eoc between presented improved approach and EDISON

3 結(jié)論

(1) 對(duì)常用的高分辨率遙感影像如 SPOT5和Quickbird的MS模點(diǎn)檢測(cè)的空間帶寬的選擇給出了參考值，其他的遙感影像根據(jù)其空間分辨率參照這2類影像來(lái)確定。

(2) 遙感影像不同波段分別用公式計(jì)算值域帶寬，以實(shí)現(xiàn)自適應(yīng)的MS遙感影像分割。

(3) 模點(diǎn)檢測(cè)時(shí)采用 4個(gè)有效的加速措施，減少了MS迭代時(shí)間；同時(shí)，將遙感影像的多維特征分組，模點(diǎn)檢測(cè)時(shí)采用空間特征和多光譜特征，模點(diǎn)全局融合時(shí)采用紋理特征，這樣加快了MS迭代速度。

(4) 采用基于全局的模點(diǎn)融合提高遙感影像分割精度。用 MATLAB實(shí)現(xiàn)了本文改進(jìn)的遙感圖像分割方法，實(shí)驗(yàn)結(jié)果與由傳統(tǒng)算法所得結(jié)果相比，本文方法具有分割效果好等優(yōu)點(diǎn)。

[1]HONG Yi-ping, YU Jian-qiang, ZHAO Dong-bin. Improved mean shift segmentation approach for natural images[J]. Applied Mathematics and Computation, 2007, 185(2): 940?952.

[2]van der Sande C J, de Jong S M, de Roo A P J. A segmentation and classification approach of IKONOS-2 imagery for land cover mapping to assist flood risk and flood damage assessment[J]. International Journal of Applied Earth Observation and Geoinformation, 2003, 4: 217?229.

[3]楊朝云, 陳光儒, 呂嫦艷, 等. 一種改進(jìn)的高分辨率遙感影像分割方法及應(yīng)用[J]. 電腦信息技術(shù), 2010, 18(4): 45?48.YANG Zhao-yun, CHEN Guang-ru, Lü Chang-yan, et al. An improved high resolution RS images partition way and its applying[J]. Computer and Information Technology, 2010, 18(4):45?48.

[4]Comaniciu D, Meer P. Mean-Shift: A robust approach towards feature space analysis[J]. IEEE Transaction on Pattern Analysis and Machine Intelligence, 2002, 24(5): 603?619.

[5]Comaniciu D, Ramesh V, Meer P. The variable bandwidth mean shift and data-driven scale selection[C]//Proceedings of Eighth IEEE International Conference on Computer Vision. Vancouver,United States, 2001: 438?445.

[6]王兆虎, 劉芳, 焦李成. 一種基于視覺(jué)特性的遙感圖像分割[J]. 計(jì)算機(jī)學(xué)報(bào), 2005, 28(10): 1687?1691.WANG Zhao-hu，LIU Fang，JIAO Li-cheng. A remote sensing image segmentation algorithm based on vision information mean shift algorithm accelerated mean shift algorithm[J]. Chinese Journal of Computers, 2005, 28(10): 1687?1691.

[7]Comaniciu D. An algorithm for data-driven bandwidth selection[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2003, 25(2): 281?288.

[8]Comaniciu D, Ramesh V, Meer P. Real-time tracking of non-rigid objects using mean shift[C]//Proceedings of IEEE Conference on Computer Vision and Pattern Recognition. Hilton Head Island, United States, 2000: 142?149.

[9]Georgescu B, Shimshoni I, Meer P. Mean shift based clustering in high dimensions: a texture classification example[C]//Proceedings of IEEE International Conference on Computer Vision. Nice, France, 2003: 456?463.

[10]Carreira-Perpinan M A. Acceleration strategies for Gaussian mean—shift image segmentation[C]//Proceedings of the 2006 IEEE Computer Society Conference on Computer Vision and Pattern Recognition. New York: IEECS Press, 2006: 543?549.

[11]Kaftan J N. Mean shift segmentation evaluation of optimization techniques[C]//International Conference on Computer Vision Theory and Application. Madeira, Portugal, 2008: 365?374.

[12]Carreira-Perpinan M A. Gaussian mean shift is an EM algorithm[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2007, 29(5): 767?776.

[13]Polak M, ZHANG Hong, PI Ming-hong. An evaluation metric for image segmentation of multiple objects[J]. Image and Vision Computing, 2009, 27: 1223?1227.