陳啟英，安裕倫*，奚世軍

(1.貴州師范大學地理與環(huán)境科學學院，貴陽 550000；2.貴州省山地資源與環(huán)境遙感重點實驗室，貴陽 550000)

隨著遙感技術的應用發(fā)展，遙感領域進入了更大范圍、更高精度、更深層次的應用階段[1-7]。在喀斯特地形破碎與生態(tài)脆弱等復雜條件下，迫切需要高時空分辨率的植被覆蓋度等數(shù)據(jù)來支撐植被動態(tài)、石漠化綜合防治等相關研究[8-9]。由于技術和預算的限制，單一傳感器難以獲取同時滿足高空間分辨率和高時間分辨率的數(shù)據(jù)[10-11]。為滿足目前研究對高時空分辨率數(shù)據(jù)的需求，進行大范圍、高精度、快速變化的地面信息遙感監(jiān)測研究，并解決衛(wèi)星遙感數(shù)據(jù)獲取能力不足的問題，一些學者提出了一種利用高空間低時間分辨率數(shù)據(jù)和低空間高時間分辨率數(shù)據(jù)融合生成高時空分辨率數(shù)據(jù)的技術，即遙感數(shù)據(jù)時空融合技術[7,12-13]。Gao等[14]綜合像元間的空間距離、光譜差異，以及像元時間上的差異等因子，最早提出了一種時空自適應反射融合模型(STARFM),用于預測Landsat空間分辨率和MODIS時間頻率下的日地表反射率；Hilker等[15]提出了一種映射地表反射率變化的時空自適應融合算法(STAARCH),用于生成合成圖像和檢測變化。Liu等[16]在假設多源遙感數(shù)據(jù)交叉使用可以通過有效的混合技術最小化觀察誤差基礎上，基于數(shù)據(jù)同化理論提出了一種實用的動態(tài)同化STARFM算法(DASTARFM)，用于農(nóng)作物產(chǎn)量估算。Zhu等[10]基于已有的STARFM模型提出了增強型時空自適應反射融合模型(ESTARFM),該算法考慮像元反射率的時間變化特征，結(jié)合光譜解混理論，更好地預測異質(zhì)景觀下的反射率。Wu等[17]綜合考慮空間變化和非線性時間變化信息，提出了一種生成高時空分辨率影像的多源遙感數(shù)據(jù)時空融合方法(STDFM)。Zhu等[18]基于混合的方法，將空間插值和STARFM模型的思想集成，提出了一種靈活的時空數(shù)據(jù)融合模型(FSDAF)，該模型能夠有效地預測高空間分辨率圖像的土地覆蓋類型的變化，并且保留更多的空間細節(jié)。上述時空融合模型研究主要集中于平原地區(qū)，其中STARFM、ESTARFM、FSDAF三種時空融合模型較為常用且較適用于異質(zhì)景觀，但在喀斯特高原復雜條件下的表現(xiàn)有待深入研究。

圖1 研究區(qū)地理空間特征及分類Fig.1 Geographical spatial characteristics and land cover of the study area

針對中國西南喀斯特地區(qū)持久性云雨天氣等原因而導致Landsat系列衛(wèi)星獲取的遙感數(shù)據(jù)缺失的問題，遙感數(shù)據(jù)“時空矛盾”問題，MODIS數(shù)據(jù)空間分辨率較低不適合小尺度應用，喀斯特山區(qū)迫切需要高時空分辨率數(shù)據(jù)來支撐植被動態(tài)、石漠化綜合防治相關研究等問題，采用Landsat8 OLI數(shù)據(jù)和MODIS數(shù)據(jù)，以喀斯特高原的20 km×20 km區(qū)域為實驗區(qū)，生成30 m具有MODIS時間分辨率數(shù)據(jù)，對比研究STARFM、ESTARFM、FSDAF三種時空融合模型在異質(zhì)性較高的喀斯特高原區(qū)的適用性，以期為后續(xù)喀斯特山區(qū)的相關研究提供科學依據(jù)。

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

1.1 研究區(qū)概況

為了對比研究STARFM、ESTARFM、FSDAF三種時空融合模型在喀斯特復雜破碎地形的應用能力，選取貴州省遵義市境內(nèi)20 km×20 km范圍作為研究區(qū)。該區(qū)域跨東經(jīng)107°1′30″～107°12′30″，北緯27°23′00″～27°34′00″，海拔最低為572 m，最高海拔為1 342 m，屬于典型的喀斯特高原地貌區(qū)。該區(qū)域內(nèi)土地覆蓋類型復雜多樣，其主要包括有林地、灌木林地、疏林地、草地、水域、建設用地、水田、旱地等，將其作為實驗區(qū)域可以驗證STARFM、ESTARFM、FSDAF三種模型在喀斯特高原區(qū)的時空融合精度。研究區(qū)地理空間特征和土地覆蓋類型如圖1所示。

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

圖2 研究區(qū)Landsat8和MODIS原始數(shù)據(jù)Fig.2 Original data of Landsat8 and MODIS images of study eare

研究選取Landsat8搭載的OLI傳感器所獲取的多光譜影像作為高空間低時間分辨率數(shù)據(jù)，以Terra衛(wèi)星上搭載的MODIS傳感器獲取的地表反射率產(chǎn)品(MOD09GA)作為低空間高時間分辨率數(shù)據(jù)。利用STARFM、ESTARFM、FSDAF三種時空融合模型，將OLI傳感器和MODIS傳感器獲取的數(shù)據(jù)之間的對應波段地表反射率進行時空融合，生成同時具備高時間、高空間分辨率的地表反射率數(shù)據(jù)。研究選取2015年4月12日和2015年10月21日的OLI數(shù)據(jù)和MOD09GA 數(shù)據(jù)作為基礎數(shù)據(jù)，以2016年2月10日的MOD09GA數(shù)據(jù)為預測時期的高時間低空間分辨率數(shù)據(jù)，預測2016年2月10日的30 m分辨率數(shù)據(jù)，如圖2所示。為了對結(jié)果進行驗證，以2016年2月10日的Landsat8 OLI數(shù)據(jù)與時空融合結(jié)果影像進行對比評價。OLI傳感器獲取的數(shù)據(jù)空間分辨率為30 m，MODIS傳感器獲取的數(shù)據(jù)空間分辨率為250 m和500 m。通過ENVI對OLI數(shù)據(jù)進行輻射定標、大氣校正、裁剪。MOD09GA地表反射率原始數(shù)據(jù)利用MODIS Reprojection Tool (MRT)進行重投影，與OLI數(shù)據(jù)的UTM-WGS84坐標系一致，且采用雙線性插值法將250 m和500 m空間分辨率數(shù)據(jù)重采樣為30 m。

2 時空融合方法及評價指標

2.1 STARFM模型

STARFM模型[14,19]通過一對或兩對ta時刻的中高分辨率數(shù)據(jù)、低分辨率數(shù)據(jù)以及tb時刻的低分辨率數(shù)據(jù)，結(jié)合不同的空間權重融合計算出tb時刻的高分辨率數(shù)據(jù)，其計算公式如下：

Rhigh(xi,yi,ta)-Rlow(xi,yj,ta)]

(1)

(2)

式中：Rhigh(xi,yj,ta)和Rlow(xi,yj,ta)分別為ta時刻時在(xi,yj)處的高空間分辨率和低空間分辨率影像對應的地表反射率值；(xw/2,yw/2)為移動窗口的中心像元；w為搜索窗口的大??；n為不同時刻；W為權重矩陣；Sijk為給定位置的高、低空間分辨率數(shù)據(jù)間的光譜差異；Tijk為給定位置的輸入的低空間分辨率影像和預測時期的低空間分辨率影像的時間距離；Dijk為目標像元和中心像元的空間距離。

2.2 ESTARFM模型

ESTATFM模型[11,20-23]考慮了像元之間的空間和光譜相似性，利用兩對初期數(shù)據(jù)生成預測期的高時空分辨率數(shù)據(jù)，在初期高空間分辨率數(shù)據(jù)中搜索與中心像元光譜相似的像元，根據(jù)其空間和光譜相似性計算每個相似像元得權重和轉(zhuǎn)換系數(shù)，用相似像元計算中心像元從而計算出預測期的高空間分辨率數(shù)據(jù)，計算公式如下：

(3)

t=a,c

(4)

T′t=

(5)

2.3 FSDAF模型

FSDAF模型[18,24-25]，需輸入一景ta時刻高空間分辨率影像Fa與兩景ta、tb時刻低空間分辨率影像Ca與Cb，對ta時刻的高空間分辨率影像進行非監(jiān)督分類，獲取每一個粗像元內(nèi)的每個地物類別的比重fc，根據(jù)光譜線性分解理論，選擇從ta～tb時段地表覆蓋類型未發(fā)生變化的n個像元利用式(6)進行線性混合方程求解，計算得到每個波段的時間變化ΔF(c)。

(6)

式(6)中：l為非監(jiān)督分類得到的分類總數(shù)量，ΔC(xi,yi)=Cb(xi,yi)-Ca(xi,yi)。

(7)

(8)

(9)

E(xi,yi)[1-HI(xij,yij)]

(10)

將CW(xij,yij)歸一化得到W(xij,yij)，最后通過計算得到殘差值r(xij,yij)和ΔF(c)，解算出類別為c時tb和ta時刻高分辨率影像的變化值ΔF(xij,yij)：

ΔF(xij,yij)=mE(xi,yi)W(xij,yij)+ΔF(c)

(11)

通過空間距離權重對所有相似像素的變化信息求和，獲得目標像素的總變化值，并將其與ta時刻的高空間分辨率影像的地表反射率值相加得到tb時刻的最終預測值Fb(xij,yij)，其中c為移動窗口中心像元位置，w為窗口大小，具體計算公式為

(12)

2.4 時空融合定量評價指標

研究采用相關系數(shù)(R)和均方根誤差(RMSE)，以及結(jié)構相似性指數(shù)(SSIM)3個評價指標分析3種時空融合模型的精度。

相關系數(shù)(R)[25-26]，用于反映時空融合的結(jié)果影像與參考影像之間的光譜相似性，R越接近1，光譜相似性越高，其計算公式為

(13)

均方根誤差(RMSE)[25-26]，是用來評價融合影像與真實影像之間的差異，均方根誤差越小，影像融合質(zhì)量越高，其計算公式為

(14)

結(jié)構相似性指數(shù)(SSIM)[25-26]，用來評價融合影像與真實影像之間結(jié)構的相似程度，SSIM越接近于1，則兩景影像之間在結(jié)構上相似度越大，通過式(15)計算而得。

(15)

3 結(jié)果與分析

3.1 STARFM結(jié)果分析

圖3 實驗區(qū)真實影像與預測影像局部放大圖Fig.3 Partial enlargement of real image and predicted image in study area

圖4 STARFM預測反射率與真實反射率散點圖Fig.4 Scatter plots of the actual reflectance and the predicted ones by the STARFM

由實驗區(qū)真實影像與預測影像局部放大圖(圖3)可得，STARFM預測影像較為模糊，層次性較差，失去部分空間細節(jié)，且存在一定的斑塊化現(xiàn)象。圖4分波段繪制了利用STARFM模型融合得到的2016年2月10日的預測反射率與實際反射率之間的關系，分別為Blue、Green、Red、NIR、SWIR1、SWIR2的預測反射率與實際反射率散點分布圖。散點圖中每個波段的預測反射率與實際反射率模型的值均大致分布在等值線1∶1線兩邊，且每個波段對應的R均高于0.6，分別為0.639、0.748、0.713、0.699、0.805、0.754，在SWIR1、SWIR2波段的相關性較高。六個波段的融合影像與實際影像的均方根誤差(RMSE)分別為0.018、0.018、0.021、0.041、0.026、0.022。結(jié)構相似性指數(shù)(SSIM)分別為0.986、0.986、0.980、0.930、0.957、0.968，從整體上看，STARFM模型融合結(jié)果的地表反射率與真實的地表反射率值具有較高相關性。

3.2 ESTARFM結(jié)果分析

從圖3的ESTARFM預測影像可以得出，ESTARFM模型預測的影像空間細節(jié)清晰，層次性明顯，與真實影像差異較小。Blue、Green、Red、NIR、SWIR1、SWIR2等六個波段的預測反射率與實際反射率之間的關系如圖5中的散點圖所示。散點圖中每個波段的預測反射率與實際反射率模型的值均較好的分布在等值線1∶1線兩邊，且六個波段對應的R分別為0.779、0.864、0.843、0.672、0.930、0.883，除Blue、NIR波段外，其余波段的R均高于0.8；RMSE分別為0.017、0.017、0.020、0.043、0.021、0.020；SSIM分別為0.991、0.990、0.986、0.918、0.984、0.984。從整體上看，ESTARFM模型預測影像與真實影像之間具有較高的相關性，較小的誤差和較高的的結(jié)構相似性。

圖5 ESTARFM預測反射率與真實反射率散點圖Fig.5 Scatter plots of the actual reflectance and the predicted ones by the ESTARFM

3.3 FSDAF結(jié)果分析

由圖3的FSDAF預測影像的局部放大圖可得，F(xiàn)SDAF模型預測的影像較好保留了空間細節(jié)，具有較好的層次性，與真實影像具有較高的相似性。圖6中的散點圖顯示了Blue、Green、Red、NIR、SWIR1、SWIR2等六個波段與實際影像對應波段之間的反射率關系。散點圖中每個波段的預測反射率與實際反射率模型的值大致分布在等值線1∶1線兩邊，且每個波段對應的R均高于0.6，分別為0.609、0.704、0.653、0.645、0.742、0604；RMSE分別為0.017、0.018、0.025、0.046、0.030、0.024；SSIM分別為0.983、0.982、0.978、0.912、0.945、0.955。FSDAF預測影像在Green、SWIR1波段獲得與真實影像較大的相關性，分別為0.704、0.742。六個波段的R介于0.6～0.8，在喀斯特高原區(qū)地形復雜條件下，采用非監(jiān)督分類對影像進行分類，獲取地物類別比重對融合精度產(chǎn)生一定的影響，導致FSDAF模型預測影像與真實影像之間未能取得更大的相關性。

3.4 三種方法對比分析

圖3顯示了三種時空融合方法預測的2016年2月10日的Landsat 8影像，結(jié)果顯示ESTARFM模型預測的Landsat 8影像比STARFM和FSDAF預測的更接近真實影像。局部放大影像顯示ESTARFM模型的融合影像與STARFM和FSDAF模型融合影像相比，空間細節(jié)更加清晰，層次性明顯，更好地保留了小地物的形狀。相比之下，STARFM預測的影像較為模糊，層次性較差，失去部分空間細節(jié)，存在斑塊化現(xiàn)象。

圖6 FSDAF預測反射率與真實反射率散點圖Fig.6 Scatter plots diagram of the actual reflectance and the predicted ones by the FSDAF

從STARFM、ESTARFM、FSDAF三種時空數(shù)據(jù)融合模型分別生成的2016年2月10日的30 m影像的Blue、Green、Red、NIR、SWIR1、SWIR2等六個波段的預測反射率與實際反射率之間關系的散點圖來看(圖4～圖6)，ESTARFM模型得到的各波段預測反射率與實際反射率的散點圖更接近1∶1線，且從整體上看，R均大于STARFM、FSDAF模型預測結(jié)果。由此可以說明ESTARFM模型融合得到的影像與實際觀測的影像更相似，精度更高。

三種時空融合模型生成的Blue、Green、Red、SWIR1、SWIR2等六個波段整體上為ESTARFM 模型得到更高的R和SSIM，更小的RMSE(表1)，相較于STARFM和FSDAF模型，ESTARFM模型得到的融合影像R平均提高0.102、0.169；RMSE平均降低0.001、0.003、SSIM平均提高0.008、0.016。STARFM、ESTARFM、FSDAF模型得到的NIR波段的R均偏低，分別為0.699、0.672、0.645，表明預測期NIR波段影像受初期高空間分辨率NIR波段影響較大。由于NIR波段主要體現(xiàn)植被的生長狀況，物候作用會引起地表反射率發(fā)生一定的變化。三種模型得到的預測反射率與真實反射率的相關性在波長較長的SWIR1、SWIR2波段最高。從整體上看，ESTARFM模型預測的地表反射率比STARFM和FSDAF模型精度更高。

表1 三種方法整體質(zhì)量評價結(jié)果 Table 1 Accuracy assessment of three data fusion methods

4 結(jié)論

為了解中國西南喀斯特地區(qū)持久性云雨天氣等原因而導致Landsat系列衛(wèi)星獲取的遙感數(shù)據(jù)缺失的問題，遙感數(shù)據(jù)“時空矛盾”問題，MODIS數(shù)據(jù)空間分辨率較低不適合小尺度應用，喀斯特山區(qū)迫切需要高時空分辨率數(shù)據(jù)來支撐植被動態(tài)、石漠化綜合防治相關研究等問題，采用Landsat 8 OLI數(shù)據(jù)和MODIS數(shù)據(jù)，以喀斯特高原的20 km×20 km區(qū)域為實驗區(qū)，生成30 m具有MODIS時間分辨率數(shù)據(jù)，對比研究STARFM、ESTARFM、FSDAF三種時空融合模型在異質(zhì)性較高的喀斯特高原區(qū)的適用性。主要得出以下結(jié)論。

(1)ESTARFM模型得到的融合影像空間細節(jié)最為清晰，層次性更明顯，F(xiàn)SDAF模型次之，STARFM模型的結(jié)果影像較為模糊，層次性較差，失去部分空間細節(jié)，且存在斑塊化現(xiàn)象。

(2)ESTARFM模型生成的預測期影像與真實影像之間的相關性最高，結(jié)構最為相似，融合效果

最佳。三種時空融合模型中，ESTARFM模型生成的融合影像與真實影像之間相關系數(shù)最高，R整體上高于0.75，且RMSE最低，SSIM最接近于1。

(3)ESTARFM模型利用像元間反射率隨時間的變化趨勢，結(jié)合光譜解混理論，能更好地預測異質(zhì)景觀下的地表反射率，在地形復雜、異質(zhì)性較大的喀斯特高原區(qū)具有較好的適用性。

通過STARFM、ESTARFM、FSDAF三種時空融合模型在喀斯特高原區(qū)的適用性分析，發(fā)現(xiàn)三種模型在NIR波段的R較其余波段低，由于初期高分辨率影像對預測結(jié)結(jié)果影響較大，而NIR波段主要體現(xiàn)植被的生長狀況，物候作用會引起地表反射率發(fā)生一定的變化。STARFM 模型受景觀的斑塊大小影響較大，在地表破碎喀斯特高原區(qū)適應性較差。FSDAF模型得到的預測影像相關性較低，主要是在估算每一類地物時間變化是采用的非監(jiān)督分類會對精度造成一定影響，從而增大了FSDAF模型的誤差。ESTARFM模型利用兩對初期影像計算相似像元權重和轉(zhuǎn)換系數(shù)，對于異質(zhì)性強的喀斯特高原區(qū)具有更好的適應性。下一步將研究采用時空融合方法來生成長時間序列NDVI數(shù)據(jù)集、植被覆蓋度數(shù)據(jù)集，對典型喀斯特小流域植被動態(tài)變化進行研究。