一種多時(shí)相遙感圖像多目標(biāo)檢測算法*

席彥新，計(jì)璐艷，耿修瑞

(1 中國科學(xué)院空天信息創(chuàng)新研究院 中國科學(xué)院空間信息處理與應(yīng)用系統(tǒng)技術(shù)重點(diǎn)實(shí)驗(yàn)室， 北京 100094；2 中國科學(xué)院大學(xué)， 北京 100049) (2019年10月29日收稿； 2020年3月27日收修改稿)

在遙感領(lǐng)域，目標(biāo)檢測一直是重要的研究方向，其在軍用及民用項(xiàng)目上，均有著重要作用，如偽裝目標(biāo)檢測、災(zāi)情控制、土地規(guī)劃、城區(qū)監(jiān)控以及目標(biāo)跟蹤等[1-3]。在高/多光譜圖像目標(biāo)檢測中，目標(biāo)和背景的光譜差異，是目標(biāo)能夠被有效檢測的物理依據(jù)。

理論上來看，遙感圖像目標(biāo)檢測實(shí)際上是一個(gè)二分問題：將像元標(biāo)記為目標(biāo)或者背景。經(jīng)過多年的發(fā)展，該領(lǐng)域提出了很多關(guān)于高/多光譜圖像目標(biāo)檢測算法[4]，其中大部分都是在壓制背景信息的基礎(chǔ)上，突出目標(biāo)，如1)基于概率統(tǒng)計(jì)模型的檢測方法：匹配濾波(matched filter，MF)算法[4]、自適應(yīng)余弦估計(jì)(adaptive cosine estimator，ACE)算法[5]、約束能量最小化(constrained energy minimization, CEM)算法[6]和“慧眼”(clever eye，CE)算法[7]；2)基于光譜間相似性度量：光譜角度填圖(spectral angle mapping，SAM)算法[8]；3)基于子空間模型的檢測算法：正交子空間投影(orthogonal subspace projection，OSP)算法[9-10]。

雖然以上這些方法取得了令人矚目的成就并得到廣泛應(yīng)用，但都是基于單時(shí)相的遙感圖像，沒有有效利用遙感圖像時(shí)間維信息。而目前研究表明，多時(shí)相遙感數(shù)據(jù)對地物目標(biāo)分類和探測能夠起到十分重要的作用，如基于多時(shí)相紅外圖像的淺層地下目標(biāo)探測應(yīng)用[11]；使用多時(shí)相SPOT-4 植被傳感器數(shù)據(jù)的中國東北森林的分類[12]；基于MODIS NDVI數(shù)據(jù)的土地利用/覆蓋分類及變化檢測[13]；基于多時(shí)相MODIS圖像的中國南部的水稻田分布檢測[14]和基于多時(shí)相遙感數(shù)據(jù)及景觀指標(biāo)的地表分類和變化檢測[15]等。在高/多光譜目標(biāo)探測領(lǐng)域，Geng 等[16]提出基于多重線性函數(shù)的多時(shí)相目標(biāo)檢測方法——濾波張量分析(filter tensor analysis，F(xiàn)TA)，其利用不同時(shí)相的遙感圖像作為輸入來檢測目標(biāo)，提高檢測精度。從機(jī)理上講，多時(shí)相的目標(biāo)檢測方法，因?yàn)槭褂枚鄠€(gè)時(shí)相的光譜特征，所以相當(dāng)于對目標(biāo)特征維度的增長。多時(shí)相的數(shù)據(jù)相比于單時(shí)相，因?yàn)槠涮卣骶S度的增加，所以更具有可分性，因此對于一些在單個(gè)時(shí)相上難以區(qū)分的光譜，通過其加長的特征維度，可以實(shí)現(xiàn)更準(zhǔn)確的分類。所以本質(zhì)上多時(shí)相目標(biāo)檢測方法就是通過增加目標(biāo)的特征維度以更好地區(qū)分目標(biāo)和背景。

上述目標(biāo)檢測算法均只能檢測單個(gè)類型的目標(biāo)(簡稱為“單目標(biāo)”)，但如果單個(gè)目標(biāo)的光譜出現(xiàn)變化、同種目標(biāo)的光譜不完全一致或多個(gè)目標(biāo)光譜差異較大時(shí)，單目標(biāo)檢測算法的效果將會(huì)變差，在實(shí)際應(yīng)用中效果不佳。為在一幅遙感圖像中同時(shí)檢測多個(gè)目標(biāo)光譜，Ren等[17]在CEM的基礎(chǔ)上，提出多目標(biāo)約束能量最小化(multiple target CEM，MTCEM)，使得CEM的輸入目標(biāo)增加為多個(gè)；此外，多目標(biāo)探測算法還包括加和約束能量最小化(sum CEM, SCEM)[17]和勝者全贏約束能量最小化(winner-take-all CEM，WTACEM)算法[17]。在此基礎(chǔ)上，尹繼豪等[18]對上述3種多目標(biāo)檢測算法進(jìn)行改進(jìn)，修正原有算法中的自相關(guān)矩陣，減少目標(biāo)向量對自相關(guān)矩陣的影響，如修正的多目標(biāo)約束能量最小化(revised MTCEM，RMTCEM)算法使用修正的自相關(guān)矩陣代替原有的自相關(guān)矩陣。但從時(shí)相角度來看，以上算法都是基于單時(shí)相遙感圖像的算法，無法解決多時(shí)相問題。

因此本文在FTA的基礎(chǔ)上，通過引入多目標(biāo)形式，從FTA推導(dǎo)出多時(shí)相高/多光譜數(shù)據(jù)的多目標(biāo)檢測算法，以解決高/多光譜的多時(shí)相、多目標(biāo)檢測問題。

1 背景介紹

本節(jié)介紹基于多時(shí)相數(shù)據(jù)的目標(biāo)檢測算法FTA[16]。

(1)

(2)

其中vec( )是一個(gè)將L1×L2×…×LM維的張量轉(zhuǎn)化為L1×L2×…×LM維向量的運(yùn)算符，同理可得

W×1d(1)×2d(2)…×Md(M)=

(d(M)?d(M-1)?…?d(1))Tvec(W),

(3)

(4)

可得FTA濾波算子的表達(dá)式

(5)

2 多時(shí)相多目標(biāo)檢測算法

上一節(jié)介紹的FTA在相關(guān)性比較低的多時(shí)相圖像中單目標(biāo)檢測效果良好，但是實(shí)際應(yīng)用中，目標(biāo)光譜通常會(huì)有變化。在大多數(shù)高/多光譜圖像中，存在多個(gè)目標(biāo)情況或者目標(biāo)本身的光譜曲線會(huì)因?yàn)閭鞲衅骰蚬庹盏劝l(fā)生變化，不同區(qū)域的同一個(gè)目標(biāo)會(huì)有不同的光譜。但是FTA無法適用于上述情況，因此需要多目標(biāo)的檢測算法來解決這個(gè)問題。本節(jié)以多目標(biāo)替換單目標(biāo)的方式，提出多時(shí)相多目標(biāo)檢測算法(multiple target FTA，MTFTA)。

(6)

由于

(7)

(8)

(9)

可以推出

(10)

從而得到MTFTA的輸出為

(11)

3 實(shí)驗(yàn)與分析

本實(shí)驗(yàn)共分為2個(gè)部分，包括一份模擬數(shù)據(jù)和一份真實(shí)遙感數(shù)據(jù)。真實(shí)遙感數(shù)據(jù)使用多時(shí)相的Landsat8 OLI (operational land imager) 數(shù)據(jù)。本實(shí)驗(yàn)中Landsat8數(shù)據(jù)有7個(gè)波段，包括可見光部分(深藍(lán)、藍(lán)、綠、紅)、1個(gè)近紅外波段和2個(gè)短波紅外波段。本節(jié)將本文提出的MTFTA算法檢測結(jié)果與FTA、MTCEM、RMTCEM和SCEM等算法相比較，其中在使用單時(shí)相算法MTCEM、RMTCEM及SCEM時(shí)，將多個(gè)時(shí)相的模擬數(shù)據(jù)組合在一起，對整體進(jìn)行目標(biāo)探測。為進(jìn)行定量化比較，繪出算法的受試者工作特性曲線(receiver operating characteristic curve, ROC curve)，并且計(jì)算ROC曲線下的面積(area under curve, AUC)；此外還比較總體精度(overall accuracy, OA)和F-score。其中，ROC曲線和AUC值可直接從算法的輸出圖像獲得，而OA和F-score則需要從算法檢測結(jié)果的二分類圖像獲得。在進(jìn)行二分類時(shí)，需使用閾值，本文使用Youden 法則[20]確定閾值。OA和F-score都可以從表1混淆矩陣計(jì)算出。

表1 二分類的混淆矩陣Table 1 Confusion matrix of binary classification

(12)

(13)

3.1 模擬數(shù)據(jù)實(shí)驗(yàn)

在模擬實(shí)驗(yàn)中，模擬3個(gè)時(shí)相的多光譜圖像，每個(gè)圖像有7個(gè)波段，每個(gè)時(shí)相的圖像大小為200像元×200像元，其中設(shè)置4個(gè)目標(biāo)，分別為40像元×40像元，添加的目標(biāo)光譜曲線分別為圖1(a)～1(c)所示。對圖像的上半部分(背景1，100像元×200像元)和下半部分(背景2，100像元×200像元)添加不同的背景光譜(背景光譜數(shù)據(jù)見圖1(d))。為方便表述，將左上角目標(biāo)稱作目標(biāo)1，順時(shí)針方向標(biāo)記目標(biāo)2、3、4，如圖1(e)所示。為了模擬自然界可能出現(xiàn)的異物同譜現(xiàn)象，在3個(gè)時(shí)相中，分別人為設(shè)置目標(biāo)與背景光譜相似的情況：在第1時(shí)相模擬數(shù)據(jù)(圖1(f))中，目標(biāo)1、3分別與背景1、2光譜相似；在第2時(shí)相數(shù)據(jù)中(圖1(g))，令目標(biāo)1和背景2光譜相似；在第3時(shí)相中(圖1(h))，目標(biāo)3的光譜與背景1光譜完全一樣。之后對目標(biāo)和背景光譜添加信噪比為13 dB，均值為0的高斯白噪聲。模擬數(shù)據(jù)的目標(biāo)和背景光譜(光譜數(shù)據(jù)值是從2014年的Landsat8 (path=029, row=031)圖像隨機(jī)挑選的植被或裸土的光譜)、真值圖及假彩色圖像(已加高斯白噪聲)如圖1所示。

圖1 模擬數(shù)據(jù)光譜、真值圖及假彩色圖像(紅:波段5, 綠:波段4, 藍(lán):波段3)Fig.1 Spectral signatures, ground truth images, and false-color images (R:band5, G:band4, B:band3) of the simulation data

首先將目標(biāo)1、2、3和4作為目標(biāo)光譜分別使用FTA進(jìn)行檢測，分別稱作FTA1、FTA2、FTA3和FTA4。之后將4個(gè)目標(biāo)的光譜向量作為多目標(biāo)輸入，使用 MTCEM、RMTCEM、SCEM與MTFTA進(jìn)行檢測(在使用MTCEM、RMTCEM與SCEM時(shí)，將3個(gè)時(shí)相的模擬數(shù)據(jù)組合在一起，形成200像元×200像元×21波段的圖像進(jìn)行目標(biāo)探測)。實(shí)驗(yàn)結(jié)果顯示如圖2。

從檢測結(jié)果來看：首先，單目標(biāo)的檢測算法FTA(圖2(a)～2(d))只能檢測對應(yīng)的單個(gè)目標(biāo)，因此FTA只具備單目標(biāo)檢測能力；其次，MTCEM、RMTCEM和SCEM在多目標(biāo)檢測方面均要好于FTA(圖2(e)～2(g))，因?yàn)樗麄兙昧怂心繕?biāo)的光譜；但是在背景壓制方面效果較差，而且對噪聲敏感。從視覺角度來看：MTFTA(圖2 (h))對背景的壓制要好于前3種單時(shí)相算法，且對噪聲敏感程度下降，檢測結(jié)果更接近于真值圖。

為了定量化比較各算法性能，圖3給出所有算法的ROC曲線，可以看出，MTFTA的ROC曲線要明顯高于其余算法，表明檢測性能最佳。此外，還計(jì)算了相應(yīng)的AUC、OA和F-score(表2)，同樣表明MTFTA在模擬數(shù)據(jù)上的檢測效果要明顯好于其他算法。

圖3 模擬數(shù)據(jù)檢測ROC曲線Fig.3 ROC curves of detection results on the simulation data

3.2 真實(shí)數(shù)據(jù)實(shí)驗(yàn)

此部分實(shí)際數(shù)據(jù)包含3景Landsat8圖像，圖像大小均為137像元×300像元×7波段，3景圖像分別為2014年的第148天(春天)、第180天(夏天)和第308天(冬天)(假彩色圖像見圖4(a)～4(c))。研究區(qū)域位于美國內(nèi)布拉斯卡州，主要包含6種地物，分別是水體(池塘)、大豆、玉米、草地、冬小麥和苜蓿。真值圖來源于2014年美國農(nóng)業(yè)部網(wǎng)站上的農(nóng)田數(shù)據(jù)層(cropland data layer, CDL)(https:∥nassgeodata.gmu.edu/CropScape/)，如圖4(d)所示(藍(lán)色為水體，棕色為冬小麥田，粉色為苜蓿田，綠色為大豆田，淺黃色為草地，黃色為玉米田)。

本實(shí)驗(yàn)選取圖4(d)左上角藍(lán)色的水體(目標(biāo)1)、右下角棕色的冬小麥田(目標(biāo)2)以及粉色的苜蓿田(目標(biāo)3)作為檢測目標(biāo)。根據(jù)圖4(e)所示的光譜曲線來看，水體的光譜曲線在3個(gè)時(shí)相中變化不大；冬小麥田也一直都是植被；而苜蓿田在第1時(shí)相中是裸土，在后2個(gè)時(shí)相中展現(xiàn)了植被的光譜特征。與模擬實(shí)驗(yàn)類似，對目標(biāo)1、2和3分別使用FTA算法，結(jié)果稱作FTA1、FTA2 和FTA3；對3個(gè)目標(biāo)同時(shí)應(yīng)用MTCEM、RMTCEM、SCEM和MTFTA，檢測結(jié)果如圖5所示。

圖4 Landsat8農(nóng)田數(shù)據(jù)假彩色圖像(紅:波段5, 綠:波段4, 藍(lán):波段3)，真值圖及目標(biāo)光譜曲線Fig.4 False color images (R:band5, G:band4, B:band3), ground truth image, and target signatures of the farmland data of Landsat8

從檢測結(jié)果(圖5)可以看出，F(xiàn)TA可以很好地將某一個(gè)目標(biāo)單獨(dú)檢測出來，但是僅對單個(gè)目標(biāo)有效(圖5 (a)～5(c))。而將3個(gè)時(shí)相的影像組合在一起，使用MTCEM、RMTCEM和SCEM檢測的效果不佳(圖5 (d)～5(f))，原因在于：有些背景像元與目標(biāo)在某個(gè)時(shí)相內(nèi)具有相似的光譜特征。比如，圖5 (g)中紅框內(nèi)部的很多背景像元，在圖5 (d)～5(f)的檢測結(jié)果中被錯(cuò)檢測為目標(biāo)。而使用基于濾波張量分析的MTFTA算法時(shí)，可以較好地壓制圖5 (g)中紅框內(nèi)部的大部分背景像元，原因是MTFTA能夠更加有效地使用多時(shí)相信息,從而其檢測結(jié)果(圖5 (g))效果最好。

圖5 Landsat8農(nóng)田數(shù)據(jù)檢測結(jié)果圖Fig.5 Detection results on farmland data of Landsat8

以上檢測方法對應(yīng)的ROC曲線如圖6所示，可以看出MTFTA算法的ROC曲線表現(xiàn)也是最好的，相對應(yīng)的AUC、OA和F-score(表3)也最高。綜上，MTFTA在農(nóng)田的檢測效果是最佳的。

表3 Landsat8農(nóng)田數(shù)據(jù)檢測結(jié)果AUC, OA和F-scoreTable 3 The AUC, OA, and F-score of detection results on the farmland data of Landsat8

圖6 Landsat8農(nóng)田數(shù)據(jù)檢測ROC曲線(橫軸使用對數(shù)坐標(biāo))Fig.6 ROC curves of detection results on the farmland data of Landsat8 (The x-axis is plotted in logarithmic coordinates)

4 分析與討論

本節(jié)主要討論多時(shí)相多目標(biāo)檢測算法MTFTA的時(shí)間復(fù)雜度。

MTFTA的計(jì)算復(fù)雜度包含3部分：

其中：L=L1×L2×…×LM，因此總共的時(shí)間復(fù)雜度是O(L2N+L3+(L2q+Lq2))。在L>>q時(shí)，時(shí)間復(fù)雜度將簡化為O(L2(N+q)+L3)，也就是((L1×L2×…×LM)2(N+q)+(L1×L2×…×LM)3)。因此可以知道，MTFTA的時(shí)間復(fù)雜度隨著N和L1×L2×…×LM增加而增加。尤其是當(dāng)時(shí)相數(shù)M或者相應(yīng)的波段數(shù)L1,L2,…,LM很大的時(shí)候，計(jì)算復(fù)雜度將會(huì)變得很大。這里還給出FTA、MTCEM、RMTCEM和SCEM的時(shí)間復(fù)雜度，如表4所示。

表4 FTA, MTCEM, RMTCEM, SCEM及MTFTA的計(jì)算復(fù)雜度Table 4 Computational complexity of FTA, MTCEM, RMTCEM, SCEM, and MTFTA

根據(jù)表4，可以得到各個(gè)算法的計(jì)算復(fù)雜度隨著單個(gè)時(shí)相的波段數(shù)La(此處令L1=L2=…=LM=La)、總像元數(shù)N、時(shí)相數(shù)M變化的曲線如圖7所示。從圖7可以看出，MTFTA和FTA的計(jì)算復(fù)雜度基本相當(dāng)，兩者的計(jì)算復(fù)雜度在大部分情況下都大于MTCEM、SCEM和RMTCEM，并且隨著波段數(shù)La、像元數(shù)N和時(shí)相數(shù)M的增長，增長幅度遠(yuǎn)大于MTCEM、RMTCEM和SCEM。

圖7 算法的浮點(diǎn)運(yùn)算次數(shù)與波段數(shù)La、像元數(shù)N和時(shí)相數(shù)M的關(guān)系(縱軸使用對數(shù)坐標(biāo)系)Fig.7 Number of flops varying with number of bands La, number of pixels N, and number of time phases M (The y-axis is plotted in logarithmic coordinates)

5 結(jié) 論

隨著遙感技術(shù)的快速發(fā)展，遙感圖像的獲取變得極為方便，圖像的質(zhì)量也愈來愈好，空間分辨率、時(shí)間分辨率及光譜分辨率得到很大的提升。但是作為遙感圖像處理領(lǐng)域的一個(gè)很重要的分支，遙感圖像目標(biāo)檢測算法要么集中于使用單時(shí)相的圖像，對于多時(shí)相的信息利用不夠充分；要么僅限于單目標(biāo)的多時(shí)相探測，如FTA，對目標(biāo)光譜的一致性要求較高。本文通過對限制條件進(jìn)行改進(jìn)，令算法對每個(gè)目標(biāo)的輸出能量都限制為1的前提下，得到基于FTA的多時(shí)相遙感圖像的多目標(biāo)檢測算法，在模擬多光譜圖像及真實(shí)多光譜圖像上檢測效果良好。但是與CEM和FTA一樣，本文的MTFTA算法目前只適用于目標(biāo)在圖像中的比例比較小或者目標(biāo)光譜變化比較小的情況。另外，由于MTFTA的計(jì)算復(fù)雜度較高，所以當(dāng)遙感圖像的波段數(shù)或者時(shí)相數(shù)很大時(shí)，會(huì)出現(xiàn)耗時(shí)較多以及占用內(nèi)存增大的情況，因此未來的工作將致力于提高M(jìn)TFTA的計(jì)算效率和存儲效率。