吳勇鋒,江 洪
(福州大學數(shù)字中國研究院,空間數(shù)據(jù)挖掘與信息共享教育部重點實驗室,衛(wèi)星空間信息技術(shù)綜合應用國家地方聯(lián)合工程研究中心,福建 福州 350108)
地形陰影是地形影響的一種,它的存在嚴重干擾山區(qū)植被信息提取精度,給山區(qū)土地覆被解譯以及生態(tài)參量的遙感反演帶來巨大挑戰(zhàn)[1].目前,地形陰影校正方法主要可以概括為3類:1)基于山地輻射傳輸模型方法同時考慮太陽直射輻射、大氣散射輻射以及周圍地形反射輻射,效果較好,如Proy校正模型[2]、Sandmeier模型[3]和光學遙感反射率計算模型[4],但模型參數(shù)較多且計算過程復雜,推廣難度大;2)基于DEM數(shù)據(jù)的經(jīng)驗校正方法以反射率與cosi的統(tǒng)計關(guān)系為基礎(chǔ)進行經(jīng)驗校正,包括統(tǒng)計-經(jīng)驗模型(如Teil?let-回歸校正[5]、b校正[6]、VECA模型[7])、歸一化模型(如二階校正模型[8]、坡度匹配模型[9])、朗伯體反射率模型(如余弦校正模型、C校正模型[3]、SCS模型[10]、SCS+C模型[11]等)以及非朗伯體反射率模型(如Minnaert模型[12]),但該類方法只對本影有效,對落影校正效果不佳;3)基于波段組合優(yōu)化計算模型的方法通過構(gòu)建特殊的植被指數(shù)來獲取消除地形影響后的植被信息,如波段比模型[13]、地形調(diào)節(jié)植被指數(shù)(Topography Adjusted Vegetation Index,TAVI)[14-15]、歸一化差值山地植被指數(shù)(Normalized Difference Mountain Vegeta?tion Index,NDMVI)[16-17]、陰影消除植被指數(shù)(Shadow Eliminated Vegetation Index,SEVI)[18-21]、植被區(qū)分陰影消除植被指數(shù)(Vegetation Distinguished and Shadow Eliminated Vegetation Index,VDSEVI)[22]等,但該類模型只能獲取落影的指數(shù)信息,缺乏多光譜信息.多光譜信息是遙感影像分類的基礎(chǔ),對山區(qū)土地利用分類至關(guān)重要.
為了獲取消除落影影響的多波段反射率信息,提升山區(qū)土地利用分類精度,筆者以福建省連江縣內(nèi)某山地丘陵區(qū)域的Landsat8 OLI影像為例,提出一種地形落影校正方法.首先采用地表反射率計算的SE?VI作為陰影校正信息源,通過構(gòu)建光照區(qū)SCS+C校正后紅綠藍波段反射率與SEVI間的隨機森林回歸模型,進而利用落影區(qū)SEVI校正地形落影,從而獲取消除地形落影影響的紅綠藍波段地表反射率的光譜信息.
1.1研究區(qū)概況研究區(qū)位于福建省福州市東部(圖1),地理坐標范圍介于119°25′46″E~119°30′55″E,26°7′21″N~26°12′10″N之間,總面積為72 km2,最小坡度為0°,最大坡度大于60°,平均坡度21.7°,標準方差10.37°.該地區(qū)以山地為主,地形崎嶇復雜,區(qū)內(nèi)植被以常綠闊葉林和竹林為主.
圖1 研究區(qū)地理位置
1.2數(shù)據(jù)源及預處理本研究采用Landsat8 OLI衛(wèi)星影像、ASTERGDEM V2高程數(shù)據(jù)2種數(shù)據(jù)源,如表1所示.其中Landsat8 OLI衛(wèi)星影像過境日期為2019年12月11日,行列號為119/042,太陽高度角36.98°,太陽方位角155.56°.
表1 數(shù)據(jù)源
為消除大氣和光照等因素對地物反射的影響,利用Landsat8 OLI衛(wèi)星影像數(shù)據(jù)元文件提供的輻射定標參數(shù)將遙感影像的DN值轉(zhuǎn)換為輻射亮度值.采用ENVI5.3中的FLAASH大氣校正模塊對影像的輻射亮度進行大氣校正,從而獲得更接近真實地物的地表反射率數(shù)據(jù),并對DEM數(shù)據(jù)進行投影變換,使其與Landsat8 OLI衛(wèi)星影像具有相同的投影坐標系統(tǒng).
為了實現(xiàn)山區(qū)遙感影像地形落影校正,首先利用隨機森林分類方法將研究區(qū)分為陰影區(qū)、光照區(qū)和水體,利用DEM數(shù)據(jù)和影像角度信息相結(jié)合提取本影,落影則由陰影區(qū)和本影的差運算得到;其次,利用研究區(qū)地表反射率數(shù)據(jù)進行SCS+C校正與SEVI計算,以光照區(qū)為掩膜對SCS+C校正結(jié)果和SEVI進行隨機采樣,生成樣本集,進而利用訓練樣本集生成隨機森林回歸模型,并結(jié)合測試樣本集完成隨機森林回歸模型精度評價.最后,結(jié)合落影區(qū)域的SEVI和隨機森林回歸模型,實現(xiàn)地形落影校正,從而獲取消除地形落影影響的地表反射率信息,之后對地形落影校正效果進行評價.研究的技術(shù)流程圖如圖2所示.
圖2 地形落影校正技術(shù)流程圖
2.1落影提取山區(qū)陰影依據(jù)其形成原因可分為本影和落影,如圖3所示.
圖3 地形陰影示意圖
本影可根據(jù)公式計算,落影由總體陰影和本影的差值求得.
其中,Sself表示陰影檢測結(jié)果,若為1則為本影,若為0則為非本影;σ表示地形坡度角,π表示180°,ω表示太陽方位角,β表示地形坡向角,γ表示太陽高度角.
2.2落影校正
2.2.1SCS+C校正SCS+C校正[11]計算公式
其中,LSCS+C表示SCS+C校正之后的像元值,LT表示SCS+C校正之前的像元值,θ表示太陽天頂角,σ表示地形坡度角,c表示地形校正參數(shù),i表示太陽入射角.
cosi的計算公式
其中,i表示太陽入射角,σ表示地形坡度角,θ表示太陽天頂角,β表示地形坡向角,ω表示太陽方位角.
2.2.2SEVI計算SEVI[18]計算公式
其中,RVI表示比值植被指數(shù),f(Δ)表示地形調(diào)節(jié)因子,SVI表示陰影植被指數(shù),Bnir為近紅外波段反射率,Br為紅光波段反射率.
f(Δ)的計算步驟:
步驟1選擇具有明顯地形陰影效應的樣區(qū),確保陰陽坡陰影對等;
步驟2使f(Δ)從0開始,以0.001為間隔進行循環(huán)疊加,同時計算SEVI與RVI的相關(guān)系數(shù)r1以及SEVI與SVI的相關(guān)系數(shù)r2;
步驟3當r1和r2的差值無限接近時,得到最優(yōu)的f(Δ).
其中,r1為SEVI和RVI的相關(guān)系數(shù),r2為SEVI和SVI的相關(guān)系數(shù),n為參與f(Δ)計算的影像像元數(shù),x為SE?VI,y1為RVI,y2為SVI.
2.2.3隨機森林回歸校正
1)在光照區(qū)隨機選取樣本點提取SEVI以及SCS+C校正后紅綠藍波段地表反射率數(shù)據(jù),將SEVI數(shù)據(jù)作為自變量,將與其對應的SCS+C校正后紅綠藍波段地表反射率數(shù)據(jù)作為因變量,制作原始樣本集.
2)使用Python3.7的Sklearn模塊將原始樣本集隨機分為訓練樣本集和測試樣本集,其中訓練樣本和測試樣本的比例為7∶3.在訓練樣本集中通過Bootstrap重抽樣得到K個與訓練樣本集相等的訓練樣本并生成K棵回歸樹,組成隨機森林回歸模型[23](見式9).在Bootstrap重抽樣時,未被抽取的概率為當N→∞時,表明每次有36.8%的數(shù)據(jù)未被抽取,這些數(shù)據(jù)被稱為袋外數(shù)據(jù)(Out of Bag,OOB)[24],可用于估計預測誤差.
其中,x為輸入向量,θk為獨立同分布隨機變量,K為回歸樹的個數(shù),N為訓練樣本集樣本數(shù).
3)回歸樹生長過程中,每個分裂節(jié)點隨機抽取所有變量中的M個變量作為當前節(jié)點的特征子集(一般選取總變量的1/3,本文中M為1)進行分裂,并使每棵回歸樹得到最大限度生長,分裂過程中不需要剪枝.
4)將所有回歸樹的預測平均值作為最終預測結(jié)果[25].其最終分類決策如下
其中,H(x)表示K棵回歸樹得到的預測值的平均值,hk(x)表示在訓練集上訓練出的第k棵回歸樹得到的估計值.
5)使用Python3.7的bayes_opt模塊的貝葉斯優(yōu)化(Bayesian Optimization)調(diào)參方法尋找使隨機森林回歸模型表現(xiàn)最優(yōu)的超參數(shù),并用最優(yōu)超參數(shù)進行模型訓練,分別得到紅綠藍波段SCS+C校正后地表反射率與SEVI的隨機森林回歸模型,結(jié)合落影區(qū)域的SEVI數(shù)據(jù)即可預測出落影區(qū)紅綠藍波段的地表反射率數(shù)據(jù).
2.3落影校正效果評價模型評價是建模過程中的核心工作之一,通過計算隨機森林回歸模型的決定系數(shù)(公式(11))來評價模型的預測能力.
其中,yi為實際觀測值為模型預估值為樣本平均數(shù),n為樣本數(shù).
紅綠藍波段隨機森林回歸模型的預測精度(決定系數(shù))如表2所示.
表2 Landsat8 OLI影像各波段隨機森林回歸模型的預測精度
采用目視分析、統(tǒng)計特征分析、光譜特征分析、相對誤差分析和反射率與cosi相關(guān)性分析等方法進行地形校正效果評價.相對誤差絕對值(Absolute Relative Error,ARE)計算公式
其中,ARE表示相對誤差絕對值,Bshadow表示陰影區(qū)域內(nèi)光譜波段的平均值,Bsunlit表示相鄰非陰影陽坡光譜波段的平均值.
3.1目視分析原始影像、SCS+C校正后影像地以及本文方法校正后真彩色影像如圖4.對比影像校正效果可知,原始影像中地形陰影效應十分明顯,陰影區(qū)地表反射率偏低,光譜信息缺失,如圖4a所示.經(jīng)SCS+C校正后,地形陰影效應削弱,本影區(qū)域光譜信息得以恢復,但仍有少部分落影區(qū)域光譜信息難以恢復,如圖4b所示.經(jīng)本文方法校正后,地形陰影效應基本消除,落影區(qū)光譜信息進一步恢復,光照區(qū)和陰影區(qū)的光譜差異減小,如圖4c所示.
圖4 地形校正前后效果對比
3.2統(tǒng)計特征分析為了定量分析地形校正方法的校正效果,分別統(tǒng)計地形校正前后影像均值、標準差以及變異系數(shù)(Coefficient of Variation,CV)[26],并進行對比分析.影像標準差反映的是影像各波段的離散程度,地形校正后,各波段影像像元間的差異縮小.若校正后影像標準差小于原始影像,表明該校正方法具有一定的校正效果.影像變異系數(shù)是指其標準差和平均值之比,地形校正后變異系數(shù)越小,表明地形校正效果越好.由地形校正前后研究區(qū)影像各波段均值、標準差及變異系數(shù)可知,SCS+C校正后,影像各波段標準差及變異系數(shù)相對于原始影像均有所下降,說明SCS+C校正都有一定的校正效果.本文方法校正后,影像各波段標準差及變異系數(shù)均小于SCS+C校正,說明本文方法彌補了SCS+C校正的不足,使影像像元間的差異進一步縮小,如表3所示.
表3 地形校正前后影像統(tǒng)計參數(shù)對比
圖5 局部地形校正前后效果對比
3.3光譜特征分析依據(jù)陰影檢測結(jié)果,選取30組本影、落影及相鄰非陰影陽坡樣本,并分別統(tǒng)計30組樣本在原始影像、SCS+C校正以及本文方法校正后紅綠藍三波段的光譜均值.統(tǒng)計結(jié)果顯示,原始影像的紅綠藍波段光譜值在本影和落影處均低于相鄰非陰影陽坡.SCS+C校正后本影處的紅綠藍波段光譜值恢復到陽坡水平,而落影依舊低于陽坡,證明SCS+C校正對本影的校正效果良好,但對落影的校正效果欠佳.經(jīng)本文方法校正后,落影處的紅綠藍波段光譜值校正到陽坡水平,證明本文方法能夠?qū)β溆斑M行校正.
3.4本影、落影相對誤差分析地形校正前后,紅綠藍波段對應的本影、落影相對誤差絕對值統(tǒng)計結(jié)果表明,如表4所示:
表4 紅綠藍波段本影落影區(qū)域相對誤差絕對值
1)未經(jīng)地形校正的各光譜波段本影、落影與其相鄰陽坡的相對誤差較大,相對誤差絕對值在30%~50%左右.
2)SCS+C校正后,紅綠藍波段本影與其相鄰陽坡的相對誤差大幅降低,分別從40.40%、43.43%、29.28%降至8.75%、13.35%、6.28%,而落影與其相鄰陽坡的相對誤差降幅相對較小,分別從48.76%、51.30%、38.50%降至38.25%、41.23%、29.84%,表明SCS+C校正對本影校正效果較好,而對落影校正效果欠佳.
3)本文方法校正后,紅綠藍波段落影與其相鄰陽坡的相對誤差大幅降低,分別從48.76%、51.30%、38.50%降至0.43%、1.82%、1.91%,表明本文方法在落影區(qū)域能夠取得較好的校正效果.
圖6 30組樣本紅綠藍波段折線圖
3.5反射率與cosi相關(guān)性分析由于地形效應影響,原始遙感影像地表反射率與cosi之間存在一定相關(guān)性.地形校正后,此相關(guān)性將被削弱甚至消除,故地形校正前后兩者決定系數(shù)的變化被廣泛用于檢驗地形校正的效果[21].地形校正前后各波段反射率與cosi的斜率及決定系數(shù)如表5所示.影像校正前,各波段像元值與其cosi都存在一定相關(guān)性.SCS+C校正后,研究區(qū)影像各波段斜率和決定系數(shù)均顯著降低,表明SCS+C校正方法能夠削弱地形陰影效應.本文方法校正后,斜率和決定系數(shù)的下降幅度都高于SCS+C校正,表明本文方法能夠彌補SCS+C校正的不足,使整體校正效果更優(yōu).
表5 地表反射率與光照系數(shù)線性相關(guān)性分析
結(jié)合SCS+C模型和SEVI指數(shù)完成山區(qū)植被地形落影可見光波段校正,并對本文方法及傳統(tǒng)地形校正方法進行比較.未經(jīng)地形校正的原始影像地形陰影效應明顯,本影和落影與相鄰陽坡的相對誤差在30%~50%左右.因此,利用遙感影像研究崎嶇山區(qū)植被時,要考慮地形陰影效應.SCS+C等傳統(tǒng)地形校正對本影的校正效果較好,但對落影校正效果欠佳.本文方法對山區(qū)落影區(qū)域植被可見光多光譜信息校正效果明顯,彌補了傳統(tǒng)地形校正方法對落影校正失效的缺陷和SEVI只能獲取消除地形陰影干擾的單一指數(shù)信息的不足.本文地形落影校正方法目前對可見光波段效果較好,對其他波段的適用性還有待進一步研究.