STARFM算法生成濕地類型TM反射率數據的應用評價

趙艷麗　李大成　賈琇明　崔鵬燕

(太原理工大學測繪科學與技術系　山西 太原 030024)

?

STARFM算法生成濕地類型TM反射率數據的應用評價

趙艷麗李大成賈琇明*崔鵬燕

(太原理工大學測繪科學與技術系山西 太原 030024)

當前數據獲取條件下，很難直接獲得兼具高時間與高空間分辨率的多光譜遙感數據,提出利用STARFM(SpatialandTemporalAdaptiveReflectionFusionModel)算法來合成高時間序列的高空間分辨率數據。該算法在我國地理區域的適用性與預測精度驗證等工作尚未充分展開。為此，以內蒙古呼倫湖濕地自然保護區為研究樣區，并借助于Landsat-5TM(ThematicMapper)與高時序MODIS反射率產品，利用STARFM算法生成具有高時序特征的TM數據，進而將其與真實TM數據進行對比驗證分析。結果表明：STARFM算法能夠在空間上保持一定預測精度的條件下，對濕地區域內不同地物類別隨時相的變化特征具有較好的預測能力，尤其適用于對反射特征隨時相變化較小的濕地區域進行時空擬合或數據預測研究。

呼倫湖濕地STARFM高時空分辨率Landsat-5TM

0　引　言

在當前遙感應用條件下，獲取和利用高分辨率遙感數據具有重要的應用和發展前景，例如碳循環研究、土地覆蓋和土地利用研究[1]，以及大區域制圖[2]等均需要幾十米級的地表覆蓋數據。美國陸地資源衛星Landsat的TM/ETM+(EnhancedThematicMapperPlus)數據(30米)已廣泛應用于土地覆蓋、土地覆蓋變化監測[3]以及生態系統動態變化監測中[4]。但較長的重訪周期(16天)以及大氣條件和云污染等因素的影響，嚴重制約了Landsat數據在快速地表變化監測方面的應用[5]。迄今，雖然很多具有高時間分辨率的遙感平臺也已投入運行，但其在應用上的不足之處是空間分辨率很低，如MODIS(Moderate-resolutionImagingSpectroradiometer)數據的空間分辨率最高為250米。所以，如何借助這類數據在高時序特征上的優勢來彌補LandsatTM/ETM+等中-高空間分辨率數據在捕捉時相特征能力上的不足就成為解決上述問題的有力途徑。

時空融合的思想最早由Gao等在2006年提出，他們利用所發展的STARFM融合算法來生成高時間序列的地表反射率數據，并利用模擬數據與真實數據進行了初步驗證[5]。該算法在局部尺度下的預測精度高，并得到了廣泛的應用推廣及驗證[6-9]。例如，Schmidt等對利用STARFM算法生成長期密集時間序列數據集進行了分析，取得了很好的預測結果[10]；黃登成等則驗證了STARFM算法在提高NPP估算精度中的應用潛力[11]；Hilker等[12]在進行深林干擾的制圖研究時，基于STARFM算法提出一種STAARCH(SpatialTemporalAdaptiveAlgorithmformappingReflectanceChange)時空融合算法，使用該算法可以獲取地表變化區域和發生變化的具體時相，通過選擇最佳時相的數據進行融合來提高融合精度。但該方法在確定變化區域及其對應時相時需要借助多時序的高空間分辨率數據以及相應的低空間分辨率輔助數據，這極大限制了其應用范圍。為提高STARFM算法對于地表反射率隨時相變化較大的地物(如植被等)的預測精度，ZHU等[13]在STARFM算法的基礎上提出了ESTARFM算法，當研究區有兩幅以上清晰無云的高空間分辨率數據覆蓋時，它具有比STARFM算法更為精確的融合能力[14]。然而，因該算法在推導過程中假設高空間分辨率圖像內某一地物類別在兩個特定時刻的反射率變化量是一個只與時相差異有關的增量，所以該算法更適用于那些起始融合時相與預測時相均分布在一個季相特征呈單向、線性變化的時段內的情形。

為明確處理二向性反射BRDF(Bidirectionalreflectancedistributionfunction)問題，Roy等[15]在2008年提出了一種基于半物理模型的時空融合算法，該算法在區域尺度具備一定的應用潛力；為解決該算法在計算復雜度、尺度差異等方面的限制，Li等[16]進一步提出了一種擴展的融合策略，該策略對處理大區域數據的時空融合問題具有效率高且光譜連續等優勢，并通過引入像元的鄰域光譜信息消除了斑塊效應的影響[17]。但相比于STARFM算法，上述基于半物理模型的時空融合算法在預測精度以及模型的推廣上都有待進一步的應用驗證。

本文即以內蒙古呼倫湖濕地自然保護區為例，以Landsat-5TM與MODIS反射率數據為基礎研究數據，通過利用STARFM算法生成研究區高時序的TM反射率數據，并分別從空間與時間維度上分析STARFM算法的預測精度及其對于時相變化特征的捕捉能力。

1　研究區概況及數據預處理

1.1研究區概況

呼倫湖是我國第四大淡水湖，以其為核心的濕地自然保護區位于內蒙古呼倫貝爾市轄區內(圖1中虛線框)，地理位置為東經116°50′10″～118°10′10″，北緯47°45′50″～49°20′20″，總面積達400 000公頃，平均海拔約545米。保護區內以沙地、鹽堿地、湖泊、草地、濕地等分布為主，年平均降雨量約239毫米，主要集中于6至8月份。本文所選研究樣區處于呼倫湖東側的烏爾遜河上，其位置如圖1中黑色方框所示。

圖1　研究區地理位置

1.2數據預處理

本文所用實驗數據為Landsat-5TM影像(空間分辨率30米，回訪周期為16天)和MODIS每8天的500米地表反射率產品(MOD09A1)，并搜集了研究區夏季(5月至9月)的TM數據4幅和MODIS數據20幅 (如圖2所示，粗線表示起始融合時相下的Landsat-5TM)。其中，TM數據的軌道號為124/26(列/行)，過境日期分別為2009年的6月25日、7月27日、8月12日、9月29日；MODIS數據的軌道號為25/4(列/行)，起終日期為2009年的5月1日和9月30日。在進行融合處理前，需對Landsat-5TM數據和MODIS反射率數據進行相應的數據預處理。

Landsat-5TM數據首先需進行輻射校正，從而將DN值轉化為可用于大氣校正處理的輻射亮度值。本文利用ENVI的FLAASH大氣校正模塊對經輻射校正處理后的TM數據進行大氣校正，進而將輻射亮度值轉換為地表反射率；為與TM數據保持一致，MODIS數據需利用相應的重投影工具(MODISReprojectionTool，MRT)將反射率數據(MOD09A1)的投影類型轉換為UTM-WGS84投影坐標系，并采用最鄰近法重采樣到30米分辨率。經上述處理后的TM數據和MODIS數據還需進行精確的幾何配準，從而使它們能夠在空間上精確匹配。

圖2　Landsat-5 TM和MODIS數據的獲取日期

2　方　法

2.1STARFM算法理論基礎

當影像間配準誤差和大氣校正誤差影響較小時，某同質像元(xi,yj)在TM影像上t0時刻下的反射率L(xi,yj,t0)可以由其在MODIS影像上同時刻下的反射率M(xi,yj,t0)表示為[5]：

L(xi,yj,t0)=M(xi,yj,t0)+ε0

(1)

其中，εk為由傳感器系統誤差造成Landsat影像與MODIS影像的光譜差異值。

同理，該像元在時刻tk下的TM反射率數據與相應時刻的MODIS反射率數據間存在以下關系：

L(xi,yj,tk)=M(xi,yj,tk)+εk

(2)

其中，L(xi,yj,tk)和M(xi,yj,tk)分別為(xi,yj)位置上Landsat和MODIS的反射率值。

若Landsat影像與MODIS影像反射率差值不隨時間發生變化，則有:

εk=ε0

(3)

因此，該像元在tk時刻下TM影像上的反射率值可表達為：

L(xi,yj,t0)=M(xi,yj,t0)+L(xi,yj,tk)-M(xi,yj,tk)

(4)

然而，由于實際中地表類型的復雜性，它們通常不是同質像元，因而無法滿足上述條件下的推導過程。為了解決這一問題，STARFM算法通過借助目標像元周圍的鄰近像元的光譜信息來構建一個充分考慮了光譜、時間和空間等因素的權重函數Wijk來計算目標像元的反射率值。此權重函數決定了周圍鄰近像元對目標(中心)像元光譜信息的貢獻大小，由光譜差異Sijk、時間差異Tijk以及空間差異Dijk三個因子組成：

(5)

(6)

(7)

其中：A是特定常數。上述三個參量的值越小，則該鄰近像元對中心像元的反射率貢獻越大，所被賦予的權重系數就越高。權重函數Wijk的計算公式為：

(8)

目標像元的反射率值的最終計算公式如下：

(M(xi,yj,t0)+L(xi,yj,tk)-M(xi,yj,tk))

(9)

其中：W為局部移動窗口的大小，(xw/2,yw/2)為移動窗的中心像元位置，t0為預測時相，L(xw/2,yw/2,t0)為t0時刻中心像元的反射率值，為保證鄰近像元光譜信息的有效性，只有在移動窗內且與中心像元光譜相似的清晰像元才參與計算過程。

該融合模型在局部尺度下的預測精度高，本文主要分析和探討該模型在本研究區對不同隨著時相變化的地物特征的預測能力。

2.2時空融合策略與評價方法

本文選取2009年8月12日的Landsat-5TM數據作為STARFM算法在起始融合時相下的數據，并選取與該時相最為接近的2009年8月13日的MODIS數據為起始融合時相下的MODIS數據，然后將2009年5月1日至9月30日間其他時刻的MODIS數據依次作為預測時相下的MODIS融合數據來生成相應時相下的TM反射率數據(共19幅)。

對于本實驗所需生成的19個時相的預測結果中，6月26日、7月28日以及9月30日的TM反射率數據預測結果分別與6月25日、7月27日以及9月29日下真實的TM反射率數據進行融合精度的對比分析，由于這兩組數據在時相上僅相差一天，因而可以認為兩組數據間的反射率差異能夠反映出預測數據與真實數據的反射率差異。本文首先利用基于反射率的散點分布圖來分析上述兩組數據分別在三個時相下反射率數據的相關性，然后進一步對兩組數據中不同地物類別的反射率數據的相關性以及統計特征進行對比分析；最后分析短時期內(一個季相)不同地物類別在所預測的多時序TM數據中的NDVI變化特征與MODIS中NDVI數據的一致程度，從而驗證不同地物隨時相的變化特征。

考慮到研究區的分類精度會直接影響到本文后續的實驗分析及評價過程，因而需選取分類精度高，并且處理性能穩定的監督分類方法來進行分類處理。在充分了解研究區地物分布特征的基礎上，本文采用最大似然分類方法來對研究區在起始融合時相(2009年8月12日)下的Landsat-5TM數據進行分類處理。這里將研究區域按照草地、濕地、鹽堿地、水體、沙地等共5種地物類別進行劃分，其中每種類別至少選取10個以上的訓練樣本作為分類依據，樣本可分離度均大于1.9。將分類結果中的小圖斑剔除后，用分類結果疊加方法目視判斷分類結果的準確性，將部分不準確的分類結果合并到相應類中。

對于本文所搜集的4個時相的TM數據， 2009/07/27、2009/08/12以及2009/09/29三個時相下的影像上有水體的存在，而2009年6月25日的TM影像上卻不存在水體類別。因此，雖然上述4幅TM數據基本分布于一個季相范圍內，但相對于其他三個時相，研究區在6月25日下的地物反射特征存在顯著變化。

3　實驗分析與討論

本文擬分別從空間預測精度和時序數據擬合程度上分析STARFM算法的預測精度及其對不同地物類別隨時相變化特征的捕捉能力。評價參數采用相關系數r2來分析實驗結果的預測一致性。相關系數([0,1])表征著預測結果與真實值的接近程度，其值越大，相關程度越高，數據間的一致性就越好。

圖3展示了三個預測時相下的融合結果與真實TM反射率數據在近紅、紅和綠三個波段上的散點分布圖。其中，橫坐標為真實的TM反射率數據，縱坐標為STARFM算法的預測結果；虛線為理想擬合線，實線為預測結果與真實數據的擬合線(其直線表達式已在圖中標出)，預測結果與真實反射率數據間的相關系數r2標于圖中左上方。由圖可知，三個時相在近紅外、紅、綠這三個波段的相關系數均在0.7以上(r2>0.73)，尤其是在紅波段和綠波段上的相關系數均在0.8以上，相關程度很高。若對比每個波段在不同時相下的相關性，可發現：9月30日、6月26日以及7月28日下的預測結果與相應真實數據間的相關程度依次增高。同時，對于上述任一時相而言，紅波段和綠波段上的相關程度均比近紅外波段上的相關程度高。由于7月28日的預測結果在時相上與起始融合時相(8月12日)相差較小，地表覆蓋變化很小，因而其相關程度也較其他兩個時相高。此外，研究區域在6月和7月間正處于植被生長旺盛期，而9月30日下的植被生長則趨于衰落階段，因而6月26日的相關性較9月30日的相關性要好。

圖3　真實Landsat-5 TM影像與其預測結果的相關性

表1展示了不同地物類別的真實反射率與其預測結果間的相關性統計信息。可以發現，在研究區的不同地物類別中，鹽堿地、濕地以及沙地的真實反射率與其預測結果的相關程度較高，而每一種地類在三個目標時相(9月30日、6月26日以及7月28日)下的相關程度依次增高，且在紅波段和綠波段上的相關系數整體上要高于近紅外波段。

表1　濕地區域不同地物類別預測結果的相關性統計

此外，在2009年6月26日的預測結果中，水體的預測精度最差(相關性最低)。參照6月25日的TM真實數據可知，該時相下研究區內的水體已接近干涸，地物類別發生了明顯轉變，而其它兩個時相(9月30日和7月28日)下的水體分布特征則與起始融合時相一致。不同地物類別在不同預測時相下的融合結果與相應的TM真實圖像的差異數據的統計信息(均值與標準差)如圖4所示。橫坐標表示不同地物類別，縱坐標為預測結果影像與真實反射率數據的差異大小。十字符、米字符、三角符分別代表近紅外、紅、綠三個波段，實線和虛線分別表示均值和標準差的大小。可見，每種地物類別在所有時相和波段上的差異大小均在零值附近波動，其中沙地的統計特征最接近零值，而水體的擬合結果最差。總體上看，所有類別的差異程度從9月30日、6月26日至7月28日基本上依次降低，即融合精度依次增大，即與圖5的結論一致。

圖5為不同地物類別(鹽堿地-saline，濕地-wetland，沙地-sandyland，草地-grass，水體-water)在預測的時序數據上及其在真實MODIS數據上的NDVI曲線。易發現，所預測的研究區各地類的NDVI變化曲線與MODIS數據的NDVI參考曲線基本一致。具體而言，5月到6月間，植被生長相對緩慢，兩者的NDVI曲線較為接近；7月至8月間，植被處于快速生長時期，因而兩者的NDVI數據的差異逐漸增大，且對于草地的而言，這種差異最為突出；而在九月底時，植被生長季結束，兩者在NDVI數據上的差異程度也相應降低。

總體而言，STARFM算法能夠在保持空間預測精度的基礎上，對不同地物類別隨時相的變化特征有較好的擬合能力；但對于那些反射特征隨時相變化較大的地物，其預測精度稍差。

圖4　STARFM對濕地區域不同地物在不同時相下的反射率預測精度

圖5　STARFM對濕地區域不同地物在一個季相下的NDVI預測曲線

4　結　語

本文采用STARFM算法來預測濕地研究區下多時序的TM反射率數據，并從單時相數據的空間預測精度與多時相數據的時序變化特征上對預測結果進行充分的應用分析與評價。依據本文的實驗及分析結果，可初步得到以下結論：

(1) 在一個季相分布內，起始融合時相與預測時相的差異越小，STARFM算法對于TM反射率數據的預測精度就越高(如實驗中2009年7月28日的預測結果)。而研究區在兩個時相間的地物類別發生改變會顯著影響該算法的融合效果，如本實驗中因季相特征變化引起的部分湖泊干涸等。

(2)STARFM算法對于濕地區域不同地物類別的預測精度也不相同，對于反射特征隨時相變化較小的地物，其預測精度較高；反之，其預測效果稍差。因此，該算法適用于地物反射特征隨時相變化較小的地理區域進行時空融合研究。

(3)STARFM算法對一個季相條件下濕地區域不同地物類別隨時相的變化或分布特征具有很好的擬合能力，其本質上是能夠很好地保持融合前后各波段反射率數據間的相對關系，因而該算法在模擬或預測地物反射的變化特征或規律上具有一定的應用潛力。

[1]IGBP,HDP.Land-useandland-coverchange：science/researchplan[R].IGBPReportNo.35andHDPReportNo.7.Stockholm, 1995.

[2]CihlarJ.Landcovermappingoflargeareasfromsatellites:statusandresearchpriorities[J].InternationalJournalofRemoteSensing, 2000,21(6-7):1093-1114.

[3]WulderMA,WhiteJC,GowardSN,etal.Landsatcontinuity:Issuesandopportunitiesforlandcovermonitoring[J].RemoteSensingofEnvironment, 2008,112(3):955-969.

[4]HealeySP,CohenWB,YangZQ,etal.ComparisonofTasseledCap-basedLandsatdatastructuresforuseinforestdisturbancedetection[J].RemoteSensingofEnvironment, 2005,97(3):301-310.

[5]GaoF,MasekJ,SchwallerM,etal.OntheblendingoftheLandsatandMODISsurfacereflectance:PredictingdailyLandsatsurfacereflectance[J].GeoscienceandRemoteSensing,IEEETransactionson, 2006,44(8):2207-2218.

[6]HilkerT,WulderMA,CoopsNC,etal.GenerationofdensetimeseriessyntheticLandsatdatathroughdatablendingwithMODISusingaspatialandtemporaladaptivereflectancefusionmodel[J].RemoteSensingofEnvironment, 2009,113(9):1988-1999.

[7]LiuF,WangZ.SyntheticLandsatdatathroughdataassimilationforwinterwheatyieldestimation[C]//Geoinformatics, 2010 18thInternationalConferenceon.IEEE, 2010:1-6.

[8]WattsJD,PowellSL,LawrenceRL,etal.Improvedclassificationofconservationtillageadoptionusinghightemporalandsyntheticsatelliteimagery[J].RemoteSensingofEnvironment,2011,115(1):66-75.

[9]SinghD.Evaluationoflong-termNDVItimeseriesderivedfromLandsatdatathroughblendingwithMODISdata[J].Atmósfera,2012,25(1):43-63.

[10]SchmidtM,UdelhovenT,GillT,etal.LongtermdatafusionforadensetimeseriesanalysiswithMODISandLandsatimageryinanAustralianSavanna[J].JournalofAppliedRemoteSensing,2012,6(1):063512-1-063512-18.

[11] 黃登成,張麗,尹曉利,等.數據融合技術在提高NPP估算精度中的應用[J/OL].計算機工程與應用,2013.http://www.cnki.net/kcms/detail/11.2127.TP.20130418.1618.016.html.

[12]HilkerT,WulderMA,CoopsNC,etal.Anewdatafusionmodelforhighspatial-andtemporal-resolutionmappingofforestdisturbancebasedonLandsatandMODIS[J].RemoteSensingofEnvironment, 2009,113(8):1613-1627.

[13]ZhuX,ChenJ,GaoF,etal.Anenhancedspatialandtemporaladaptivereflectancefusionmodelforcomplexheterogeneousregions[J].RemoteSensingofEnvironment, 2010,114(11):2610-2623.

[14] 朱孝林.高時空分辨率遙感數據生成技術研究[D].北京師范大學，2010.

[15]RoyDP,JuJ,LewisP,etal.Multi-temporalMODIS-Landsatdatafusionforrelativeradiometricnormalization,gapfilling,andpredictionofLandsatdata[J].RemoteSensingofEnvironment, 2008,112(6):3112-3130.

[16]LiDC,TangP,HuCM,etal.Spatialtemporalfusionalgorithmbasedonanextendedsemi-physicalmodelanditspreliminaryapplication[J].JournalofRemoteSensing, 2014,18(2):307-319.

[17] 李大成.時空融合方法生成高時空分辨率區域地表反射率數據集的技術研究[D].中國科學院大學，2013.

APPLICATIONEVALUATIONOFSTARFMALGORITHMINGENERATINGWETLAND-TYPETMREFLECTANCEDATA

ZhaoYanliLiDachengJiaXiuming*CuiPengyan

(Department of Mapping Science and Technology, Taiyuan University of Technology, Taiyuan 030024, Shanxi, China)

Multispectralremotesensingdatawithbothhightemporalandhighspatialresolutionscanhardlybeobtaineddirectlyundercurrentdataacquisitionconditions.WeproposedtousSTARFM(spatialandtemporaladaptivereflectancefusionmodel)tosynthesisethehightimeseriesdatawithhighspatialresolution.Nevertheless,theapplicabilityofthismodelforthegeographicareasinChinaanditspredictionprecisionvalidationhavenotyetbeensubstantiallydone.Accordingly,theHulunWetlandNatureReserveinInnerMongoliaisselectedasthestudyarea,andwiththehelpofLandsat-5TMdataandhightimeseriesMODISreflectanceproducts,weusedSTARFMalgorithmtogeneratetheTMdatawithhightimeseriesfeature,andthenmadethecomparativevalidationanalysisonitwiththeactualTMreflectancedata.ResultshowedthattheSTARFMalgorithmhasahigherpredictioncapabilityonthephasefeaturesofdifferentterraincategoriesinwetlandareavaryingalongwithtimephasesundertheconditionofpreservingcertainpredictionaccuracyinspace,itisparticularlysuitableforthespatiotemporalfittingordatapredictionresearchinthosewetlandareaswherethereflectioncharacteristicschangelittlealongwiththetimephases.

HulunwetlandSTARFMHighspatiotemporalresolutionLandsat-5TM

2014-07-24。國家高技術研究發展計劃項目(2009AA122002)。趙艷麗，碩士，主研領域：遙感信息處理以及遙感與地理信息系統的一體化集成。李大成，講師。賈琇明，教授。崔鵬燕，碩士。

TP751

ADOI:10.3969/j.issn.1000-386x.2016.03.063