基于Sentinel-2和Landsat8 OLI數(shù)據(jù)融合的土地利用分類研究

趙亞杰,王立輝,孔祥兵,陰海明,池 泓,黃進(jìn)良

(1.水利部黃土高原水土流失過程與控制重點實驗室,河南 鄭州 450003；2.中國科學(xué)院測量與地球物理研究所環(huán)境與災(zāi)害監(jiān)測評估湖北省重點實驗室，湖北 武漢 430077；3.中國科學(xué)院大學(xué)，北京 100049;4.黃河水利科學(xué)研究院，河南 鄭州 450003)

土地利用變化是全球變化的主要驅(qū)動力之一，是全球變化與可持續(xù)發(fā)展研究中的熱點[1]，通過土地利用分類，既可以了解各種土地覆被類型的基本屬性，進(jìn)行農(nóng)作物種植面積提取、產(chǎn)量預(yù)估[2]，還可以認(rèn)識土地利用的區(qū)域結(jié)構(gòu)和分布特點，為進(jìn)一步研究地域差異奠定基礎(chǔ)[3].遙感是土地利用與土地覆被變化研究技術(shù)體系中的主要組成部分，目前，具有多光譜與多時相特征的遙感技術(shù)已成為獲取土地覆被動態(tài)監(jiān)測數(shù)據(jù)和土地覆被分類的最有效工具[4].常見的遙感影像有空間分辨率較低的MODIS數(shù)據(jù)，也有中空間分辨率的Landsat系列，以及高空間分辨率的高分一號、QuickBird等[5].其中，Landsat8 OLI影像獲取的多光譜影像包含15 m分辨率的全色波段和30 m分辨率的8波段影像，Sentinel-2影像則包括10 m空間分辨率的4波段影像、20 m分辨率的6波段影像和60 m分辨率的3波段影像.不同的影像在土地利用研究中已被廣泛應(yīng)用，李靜等[6]利用Landsat影像對大棚菜地信息進(jìn)行提取并取得較好的地物識別效果，李曉慧等[7]利用Landsat8 OLI時序遙感數(shù)據(jù)結(jié)合地物細(xì)微光譜特征對農(nóng)作物進(jìn)行分類，并獲得較高的總體分類精度，張衛(wèi)春等[8]基于Sentinel-2影像和地形數(shù)據(jù)在低山丘陵區(qū)采用隨機(jī)森林方法獲得的土地利用分類結(jié)果優(yōu)于傳統(tǒng)分類方式所獲得結(jié)果，為地形復(fù)雜區(qū)域的土地利用分類提供了參考，楊斌等[9]同樣使用Sentinel-2A衛(wèi)星數(shù)據(jù)提取了四川黑水河流域的干旱河谷，結(jié)果表明Sentinel-2A數(shù)據(jù)豐富的光譜特征能夠較好地滿足國土資源普查和生態(tài)環(huán)境變遷分析.但對于Sentinel-2和Landsat數(shù)據(jù)融合及其在土地利用中的應(yīng)用等方面，國內(nèi)的研究很少，國外學(xué)者如Gemine et al[10]也只是針對多種影像融合模型進(jìn)行了質(zhì)量評價.

鑒于此，本研究以山東省平邑縣為研究區(qū)，采用Sentinel-2和Landsat8 OLI遙感數(shù)據(jù)作為研究區(qū)數(shù)據(jù)源，探討在基于亮度的平滑濾波(smoothing filter-based intensity modulation,SFIM)、高通濾波(high-pass filtering,HPF)、小波變換(a trous wavelet transform,ATWT) 等3種不同融合模型下，Sentinel-2的10 m高分辨率波段與Landsat8 OLI 30 m波段數(shù)據(jù)的融合效果，并以最佳的融合波段與融合算法，獲取整個研究區(qū)的融合影像，采用隨機(jī)森林機(jī)器學(xué)習(xí)算法分別對Sentinel-2原始影像、融合獲得影像和Landsat8 OLI影像進(jìn)行分類，從分類精度上進(jìn)一步探討該融合方法在土地利用分類研究中的適用性.

1 研究區(qū)概況及數(shù)據(jù)源

1.1 研究區(qū)概況

圖1 研究區(qū)2017年Landsat8真彩色合成影像Fig.1 Landsat 8 true color composite diagram of the study area in 2017

研究區(qū)平邑縣位于山東省中南部、臨沂市西陲(117°25′—117°56′E,35°07′—35°43′N)，東西最大橫距47.35 km，南北最大縱距66.75 km，總面積約1 824.97 km2.研究區(qū)屬季風(fēng)區(qū)域大陸性氣候，具有冬季寒冷、夏季炎熱、光照充足、無霜期長的特點，同時具有明顯的山區(qū)特征，地勢南北高，中間低，略向東南傾斜，中部谷地、陸地土質(zhì)肥沃，田地分布較多，境內(nèi)土地覆被類型多樣.研究區(qū)影像如圖1所示.

1.2 數(shù)據(jù)來源及預(yù)處理

遙感數(shù)據(jù)包括Sentinel-2 L1C級別數(shù)據(jù)的藍(lán)、綠、紅、近紅外4個10 m空間分辨率波段，Landsat8 OLI數(shù)據(jù)的對應(yīng)4個波段，成像時間均為2017年9月14日.針對Landsat8數(shù)據(jù)，利用ENVI進(jìn)行輻射定標(biāo)和大氣校正等預(yù)處理.Sentinel-2數(shù)據(jù)則利用哨兵數(shù)據(jù)應(yīng)用平臺(sentinel application platform,SNAP)中的sen2cor模塊進(jìn)行大氣校正、產(chǎn)品格式轉(zhuǎn)換等預(yù)處理.對兩種預(yù)處理后的影像，再利用ENVI軟件裁剪出圖1紅色區(qū)域的子研究區(qū)，以便后續(xù)研究的應(yīng)用.

2 研究方法

首先在研究區(qū)內(nèi)隨機(jī)選取面積約為5.76 km2的子研究區(qū)，以SFIM、HPF、ATWT 3種不同融合模型為基礎(chǔ)方法，利用Sentinel-2的4個10 m高分辨率波段分別與Landsat8 OLI 30 m波段數(shù)據(jù)進(jìn)行融合，從而得到最佳的融合波段與融合算法，在此基礎(chǔ)上獲取整個平邑縣研究區(qū)的融合影像，根據(jù)2017年8月份平邑縣土地調(diào)查數(shù)據(jù)和谷歌地球高分辨率影像獲取的樣本，參考基于遙感技術(shù)的全國生態(tài)系統(tǒng)分類體系[11]，并采用隨機(jī)森林機(jī)器學(xué)習(xí)算法分別對Sentinel-2原始影像、融合獲得影像和Landsat8 OLI影像進(jìn)行分類，最后進(jìn)行分類精度評定.

2.1 數(shù)據(jù)融合算法

目前，已經(jīng)有多種圖像融合算法在現(xiàn)實中廣泛應(yīng)用，本研究選取SFIM，HPF和ATWT 3種算法作為融合基礎(chǔ).SFIM算法是一種基于空間域濾波變換的融合算法[12]，該算法簡單易用,融合影像的光譜幾乎不失真,具有較高的光譜保真度,但在高頻空間信息融合上略遜色于Brovey變換方法[13].SFIM變換的融合算法公式：

式中,BS為用該算法生成的融合影像；i為波段值，j和k分別為影像的行、列值；Bl為低分辨率影像，在此代表Landsat的4個波段(30 m)；Bh為高分辨率影像，即為Sentinel-2波段(10 m)；Bm為模擬低分辨率影像，它可以通過對10 m波段進(jìn)行低通均值濾波處理來獲得.

HPF是一種基于高通濾波的影像融合算法[14]，HPF融合算法能夠較好地提取全色影像細(xì)節(jié)信息，并把全色影像的細(xì)節(jié)信息直接疊加到多光譜影像上[15].融合后影像的邊界清晰度增強(qiáng)，可更好地保留光譜特征[16].HPF算法的定義公式：

HPK=Wa×PANKH+Wb×MULKL

式中，Wa、Wb為權(quán)，且滿足Wa+Wb=1，PANKH為Sentinel-2高空間分辨率影像波段進(jìn)行高通濾波的結(jié)果，MULKL為Landsat8影像進(jìn)行低通濾波的結(jié)果，HPK為銳化后的輸出影像.

ATWT是基于小波變化的圖像融合算法[17]，該算法利用小波變換將源圖像分解到不同頻段的特征域中，然后在不同的特征域中采用不同的融合規(guī)則和融合算子進(jìn)行圖像融合處理，從而形成新的小波金字塔結(jié)構(gòu)，最后利用小波反變換技術(shù)獲得融合圖像，所以小波變化算法能夠?qū)D像細(xì)節(jié)及邊緣信息進(jìn)行較好的體現(xiàn)[18].

2.2 隨機(jī)森林分類

隨機(jī)森林是屬于集成學(xué)習(xí)的一種組合分類算法，集成學(xué)習(xí)的核心思想就是將若干個弱(基)分類器組合起來，得到一個分類性能顯著優(yōu)越的強(qiáng)分類器.如果各弱分類器之前沒有強(qiáng)依賴關(guān)系、可并行生成，就可以使用隨機(jī)森林算法.隨機(jī)森林利用自主抽樣法從原數(shù)據(jù)集中有放回地抽取多個樣本，對抽取的樣本先用弱分類器—決策樹進(jìn)行訓(xùn)練，然后把這些決策樹組合在一起，最終的分類或預(yù)測結(jié)果由多顆樹分類器通過投票得出[19].隨機(jī)森林算法較強(qiáng)的抗干擾能力，表現(xiàn)性能的優(yōu)異，相較于人工神經(jīng)網(wǎng)絡(luò)、決策樹、支持向量機(jī)等機(jī)器學(xué)習(xí)方法，在分類精度和分類效率方面更具優(yōu)勢[20],因此，越來越多的學(xué)者在進(jìn)行遙感影像分類研究時采用隨機(jī)森林機(jī)器學(xué)習(xí)算法[21-22].

3 結(jié)果分析

3.1 融合結(jié)果分析

首先在子研究區(qū)中，利用SFIM、HPF和ATWT 3種算法，分別以Sentinel-2影像的藍(lán)、綠、紅、近紅外波段作為高分辨率影像參與到融合中，采用不同算法時各波段的融合效果如圖2、圖3、圖4所示.從圖中的視覺效果來看，融合影像能夠?qū)entinel-2影像波段的空間信息和Landsat8多光譜影像的色彩信息集成在一起，相比于原始的多光譜影像，融合影像在保留其波段信息前提下?lián)碛辛薙entinel-2影像的細(xì)節(jié)信息，同時在色彩方面，融合結(jié)果與Landsat8 OLI影像的色彩差距并不明顯.但在空間細(xì)節(jié)方面，利用SFIM和HPF算法獲得的融合影像中，道路及人工表面等地類的紋理模糊感較強(qiáng)，邊界信息模糊，而圖4的融合影像表明ATWT算法對細(xì)節(jié)的處理明顯強(qiáng)于前兩種算法，融合的結(jié)果也更加精確.以Sentinel-2不同波段作為高分辨率融合影像時，融合效果也不同，可以明顯注意到，不論采用哪種算法，近紅外波段的融合影像都丟失了很多細(xì)節(jié)信息，目視效果較差.藍(lán)、綠、紅、近紅外4個波段的融合效果采用相關(guān)系數(shù)(correlation coefficient,CC)、光譜角制圖(spectral angle mapper,SAM)和相對整體維數(shù)綜合誤差(relative global-dimensional synthesis error,ERGAS)3個指標(biāo)作定量比較分析.相關(guān)系數(shù)能夠反映兩幅影像間灰度信息的相關(guān)程度，通過計算融合前后影像間的灰度相關(guān)系數(shù)，可以判斷影像融合后的光譜保持程度，若相關(guān)系數(shù)趨向于1，則說明光譜保持性好，反之則差；SAM是把圖像中的每個像元的光譜視為一個高維向量，通過計算兩向量間的夾角來度量光譜間的相似性，夾角越小，兩光譜越相似，SAM理想值為0[23]；ERGAS主要評價在光譜范圍內(nèi)所有融合波段的光譜質(zhì)量，考慮光譜變化的整體情況.它的值越接近0，表明在光譜范圍內(nèi)，融合圖像的光譜質(zhì)量越好[24].

a:Landsat8 OLI影像(30 m);b:Sentinel-2影像(10 m);c:藍(lán)波段融合結(jié)果(10 m);d:綠波段融合結(jié)果(10 m);e:紅波段融合結(jié)果(10 m);f:近紅外波段融合結(jié)果(10 m).圖2 SFIM融合結(jié)果Fig.2 Results of SFIM-based fusion

表1列出了在不同算法下各波段的融合影像質(zhì)量評價結(jié)果，從表中可以看出，ATWT算法的融合效果總體優(yōu)于SFIM和HPF算法，特別是以Sentinel-2藍(lán)、綠、紅3個高分辨率波段進(jìn)行融合所獲得影像的相關(guān)系數(shù)最高，且SAM和ERGAS值最小，但以近紅外波段進(jìn)行數(shù)據(jù)融合后，3種算法獲得融合影像各波段的相關(guān)系數(shù)較低，SAM值和ERGAS值高于其他波段融合結(jié)果，總體符合目視效果.同時，相比于藍(lán)、綠、近紅外波段，在利用Sentinel-2的紅波段作為融合時高分辨率參考影像時，融合影像的藍(lán)、綠、紅3個波段的相關(guān)系數(shù)最高，分別達(dá)到0.850 405、0.851 331、0.838 229，SAM最低，為5.644 6.

a:Landsat8 OLI影像(30 m);b:Sentinel-2影像(10 m);c:藍(lán)波段融合結(jié)果(10 m);d:綠波段融合結(jié)果(10 m);e:紅波段融合結(jié)果(10 m);f:近紅外波段融合結(jié)果(10 m).圖3 HPF融合結(jié)果Fig.3 Results of HPF-based fusion

3.2 土地利用分類結(jié)果對比

圖5 Landsat8影像分類結(jié)果Fig.5 Classification result based on Landsat 8 image

結(jié)合研究區(qū)主要地類分布，分為人工表面、水體、耕地、林地和裸地等5個類別，影像分類結(jié)果如圖5、圖6和圖7所示.通過總體分類精度、Kappa系數(shù)以及各地類的制圖精度、用戶精度，對分類結(jié)果進(jìn)行定量分析，結(jié)果如表2所示.考慮到分類所用影像空間分辨率的不同，總體而言，3種分類結(jié)果的總體精度和Kappa系數(shù)之間的關(guān)系為：Sentinel-2影像>ATWT融合影像>Landsat影像.ATWT融合影像的總體分類精度及Kappa系數(shù)為88.765 7%和0.869 0，低于Sentinel-2的89.693 9%和0.886 0，差距很小，但比Landsat8影像的總體分類精度和Kappa系數(shù)提高了6.101 4%和0.056 7.結(jié)果表明“紅波段+ATWT”融合方法能較好的保留原始影像的光譜信息，雖然在土地利用分類中只能獲得接近原始影像分類精度的結(jié)果，但該融合方法能為解決Sentinel-2影像中云污染所導(dǎo)致的影像缺失問題提供理論依據(jù)，從而為后續(xù)融合Sentine-2和Landsat8 OLI數(shù)據(jù)以獲得完整的10 m高分辨率時間序列影像提供了一種有效的方法.

4 結(jié)論

本研究結(jié)果表明：(1)利用ATWT算法并以Sentinel-2影像的紅波段為高分辨率融合參考影像時，融合效果最好.所獲得融合影像的B2、B3、B4波段的相關(guān)系數(shù)最高，分別為0.850 405、0.851 331、0.838 229，同時SAM值最低，為5.644 6.

(2)采用隨機(jī)森林對整個研究區(qū)進(jìn)行土地利用分類后，融合影像的總體分類精度為88.765 7%，Kappa系數(shù)為0.869 0，略低于Sentinel-2影像分類精度，其總體分類精度為89.693 9%，Kappa系數(shù)為0.886 0，考慮到空間分辨率的差異，Landsat8影像的總體分類精度為82.664 3%，Kappa系數(shù)為0.812 3，遠(yuǎn)低于之前兩種數(shù)據(jù)的精度.結(jié)果表明“紅波段+ATWT”融合影像在土地利用分類中可以獲得接近原始影像分類精度的結(jié)果，能較好的保留原始影像的光譜信息并應(yīng)用于較大區(qū)域的土地利用精細(xì)提取.

圖6 Sentinel-2原始影像分類結(jié)果
Fig.6 Classification result based on Sentinel-2 image

圖7 ATWT融合影像分類結(jié)果
Fig.7 Classification result based on ATWT fusion image

表2 不同影像分類精度評價表Table 2 Accuracy evaluation of images from different sources

(3)圖像融合時選取的Sentinel-2數(shù)據(jù)與Landsat8數(shù)據(jù)日期一致，而實際應(yīng)用中，往往難以獲取同一時間的影像，故不同時間Sentinel-2影像與Landsat影像融合在土地利用研究中還有待進(jìn)一步探討，以期能得到連續(xù)高分辨率時序影像，更好應(yīng)用于地面種植作物復(fù)雜區(qū)域的遙感制圖.