ndvi時(shí)間序列曲線重建方法綜述

ndvi時(shí)間序列曲線重建方法綜述

1ndvi時(shí)間序列重建ndvi（歸一化植被指數(shù)）的時(shí)間序列數(shù)據(jù)能夠準(zhǔn)確反映陸地生態(tài)系統(tǒng)植被的生長狀態(tài)、季節(jié)特征和年際變化特征。該方法在監(jiān)測和模擬世界和區(qū)域生態(tài)環(huán)境變化、植被覆蓋動(dòng)態(tài)變化、植被物候特征的識別和信息提取方面得到廣泛應(yīng)用。目前常用的NDVI時(shí)間序列主要來自NOAAAVHRR、SPOTVEGETATION和TERRA/AQUAMODIS等傳感器。NDVI反映植被生長過程的時(shí)間序列曲線通常是連續(xù)平滑的,但由于傳感器本身性能、數(shù)據(jù)傳輸過程失誤、太陽光照角度、觀測視角、地物雙向性反射以及云、大氣氣溶膠等觀測條件因時(shí)間而異,此外受地表水、冰雪等隨機(jī)干擾因素,得到的觀測值包含很多不可預(yù)測的噪音,得到的NDVI時(shí)間序列呈鋸齒狀的不規(guī)則波動(dòng),反映季節(jié)變化趨勢不明顯,使陸地生態(tài)系統(tǒng)監(jiān)測、植被覆蓋變化分析和信息提取等研究結(jié)論引進(jìn)了諸多干擾。盡管常用的NDVI數(shù)據(jù)集廣泛基于最大值合成法(MVC)來降低這些噪聲,但是一些假高值、合成期內(nèi)連續(xù)的NDVI低值等噪聲殘差還是存在,為此研究發(fā)展了一系列NDVI時(shí)間序列重建的方法。目前NDVI時(shí)間序列曲線重建的主要方法大致歸納為3類:閾值法、基于濾波函數(shù)的平滑法和曲線擬合法。3種方法既可以基于原始的NDVI時(shí)間序列曲線重建,也可以對MVC處理后的數(shù)據(jù)集進(jìn)行再處理。閾值法包括Viovy等提出的最佳指數(shù)斜率提取法(BISE)和Lovell等在此基礎(chǔ)上提出的改進(jìn)算法;基于濾波的平滑方法主要包括Savitzky-Golay(S-G)濾波法、均值迭代濾波法(MVI)、傅立葉變換和時(shí)間序列諧波分析法(HANTS);曲線擬合法以雙Logistic函數(shù)擬合法(D-LCF)和非對稱性高斯函數(shù)擬合法(AGFF)為代表。以上基本涵蓋了當(dāng)前遙感NDVI時(shí)間序列曲線重建的主要方法。目前,對于NDVI數(shù)據(jù)集預(yù)處理的方法除MVC得到大家基本認(rèn)可,廣泛應(yīng)用于NOAA、MODIS、SPOT等遙感數(shù)據(jù)生產(chǎn)合成時(shí)間序列數(shù)據(jù)初級產(chǎn)品,其他方法分別用于國內(nèi)外學(xué)者在全球或不同區(qū)域的研究,尚未達(dá)成共識。研究者在不斷改進(jìn)或者提出新的算法,希望能夠克服現(xiàn)有方法本身的適應(yīng)性或缺陷,有效濾除噪聲、反映植被真實(shí)的生長狀況和物候規(guī)律等。同時(shí),部分國內(nèi)外學(xué)者選擇多種或兩種不同的重建算法進(jìn)行了比較研究,分析各自的優(yōu)劣。比如:Julien等對其提出的IDR(IterativeInterpolationforDataReconstruction)方法與HANTS、D-LCF進(jìn)行實(shí)例樣點(diǎn)對比研究,IDR在多數(shù)樣例中具有最好的效果,但是在強(qiáng)噪聲中對NDVI低值過度擬合;Sakamoto等以日本為研究區(qū),對小波變換和傅立葉變換方法進(jìn)行了比較分析,結(jié)果表明,小波變換總體上優(yōu)于傅立葉變換;Jennifer選擇加拿大阿爾伯塔中西部的洛磯山脈中一個(gè)典型區(qū)域?yàn)槭痉秴^(qū),比較兩種曲線擬合方法(D-L擬合法和AG擬合法)與4種濾波方法(S-G濾波法、MVI、ARMD3-ARMA5濾波及4253H兩次濾波法)對NDVI時(shí)間序列重建效果,認(rèn)為兩種曲線擬合方法總體上優(yōu)于4種濾波方法。在國內(nèi),僅有極少數(shù)相關(guān)的研究,且多針對上述各種方法的原理、優(yōu)缺點(diǎn)進(jìn)行定性分析。如顧娟、李儒對目前主要的幾種NDVI時(shí)間序列方法進(jìn)行概述;吳文斌等研究發(fā)現(xiàn),S-G濾波法與非對稱高斯擬合法對NDVI時(shí)間序列數(shù)據(jù)集擬合重建效果具有區(qū)域差異性。由于不同區(qū)域大氣及地面條件、植被覆蓋類型及作物制度的差異,以及缺乏定量比較手段,本研究以中國西部地區(qū)高寒和干旱環(huán)境具有代表性的藏北草地為例,利用非對稱性高斯函數(shù)擬合、雙Logistic曲線擬合和Savitzky-Golay濾波3種最常用的擬合重建算法對2007～2009年MODIS16d合成的NDVI時(shí)間序列進(jìn)行處理,對比分析擬合結(jié)果曲線與NDVI序列上包絡(luò)曲線擬合效果以及保留高質(zhì)量原始真值程度兩方面的效果,為今后NDVI應(yīng)用研究中選擇時(shí)間序列擬合的方法提供參考。2藏北地理位置及草地分布藏北地區(qū)位于青藏高原腹地,素有“江河源”、“生態(tài)源”、“中華水塔”等美譽(yù),是我國長江、怒江、瀾滄江等主要江河的發(fā)源地,圖1為藏北地理位置及草地類型空間分布。區(qū)域幅員遼闊,緯度位置為30°27′25″～35°39′13″N,83°41′14″～95°10′46″E。由于疆域跨度大,加之多變的垂直地帶性影響,氣候及水熱條件差異較大,所以從東南向西北,植被類型分布與生物氣候帶相適應(yīng),依次出現(xiàn)亞高山疏林—灌叢草甸、高寒草甸、高寒草原、高寒荒漠草原四大基帶及相鄰的過渡亞帶。3基于平臺的ndvi研究使用的實(shí)例數(shù)據(jù)選擇EOS/Terra衛(wèi)星的MODIS產(chǎn)品之一MOD13Q1v005,包括基于MVC方法16d合成的250m分辨率NDVI及其質(zhì)量控制數(shù)據(jù),該產(chǎn)品經(jīng)過幾何校正和大氣校正。該時(shí)間序列數(shù)據(jù)集從美國NASALPDAAC工作組(LANDPROCESSESDISTRIBUTEDACTIVEARCHIVECENTER)網(wǎng)站下載(/lpdaac/products/modis_products_table),時(shí)間覆蓋2007年1月～2009年12月,3a時(shí)間共69期影像。NDVI值標(biāo)準(zhǔn)范圍為-1.0～1.0,MOD13Q1NDVI產(chǎn)品是-3000～10000的DN值,-3000為填充值,從DN值轉(zhuǎn)化成NDVI值的關(guān)系式為:NDVI=0.0001×DN(1)此外,還需從MOD13Q1v005產(chǎn)品中提取datapixelreliability數(shù)據(jù)集,其空間、時(shí)間分辨率與NDVI數(shù)據(jù)集匹配。該數(shù)據(jù)是對區(qū)域每個(gè)像元NDVI數(shù)據(jù)質(zhì)量進(jìn)行說明(如表1),分為5個(gè)不同等級,為該NDVI產(chǎn)品不同應(yīng)用級別提供參考。本研究中,NDVI時(shí)間序列擬合算法中填充NDVI點(diǎn)不參與計(jì)算,DN為0的NDVI點(diǎn)權(quán)值為1,DN為1的點(diǎn)權(quán)值設(shè)為0.8,而DN為2、3的點(diǎn)值權(quán)重為0.2。4學(xué)習(xí)方法4.1ndvi時(shí)間序列提取TIMESAT軟件包是J?nsson與Eklundh共同開發(fā)的用于植被指數(shù)時(shí)間序列數(shù)據(jù)集重建及提取植被生長物候信息的程序包,以便于獲取季節(jié)性植被生長規(guī)律信息等。該軟件平臺以S-G濾波法、AG擬合法及D-L擬合法為核心算法,能對NOAA/NOAAPAL、TERRA或AQUA/MODIS等NDVI時(shí)間序列數(shù)據(jù)進(jìn)行有效地處理提取信息。基于上述3種方法對2007～2009年共23期NDVI時(shí)間序列重建效果對比分析處理主要包括以下4個(gè)步驟:①通過MRT軟件對從NASA網(wǎng)站上獲取的HDF格式MOD13Q1產(chǎn)品進(jìn)行轉(zhuǎn)換,并進(jìn)行區(qū)域拼接、投影轉(zhuǎn)換等預(yù)處理;②基于TIMESAT2.3平臺對原始NDVI時(shí)間序列進(jìn)行平滑處理,其關(guān)鍵步驟為設(shè)置參數(shù),包括NDVI有效值域(Range,-0.05～1.00)、原始NDVI值權(quán)重(Weightvalue,詳見表1)、噪聲去除閾值(Spike,2)、滑動(dòng)窗口大小(W-size,5-6-7)、擬合峰值參數(shù)(Alptitude,5)、迭代次數(shù)(Fittingsteps,3)等,基于S-G濾波法、非對稱高斯函數(shù)擬合法及雙Logistic曲線擬合法產(chǎn)生擬合后的NDVI時(shí)間序列數(shù)據(jù)集;③基于圖1中的草地類型斑塊,各種類型草地隨機(jī)選取60個(gè)左右(共計(jì)300個(gè))的年時(shí)間序列樣點(diǎn),得到各樣點(diǎn)2007～2009年原始及3種算法擬合后的NDVI時(shí)間序列數(shù)據(jù)集;④利用Matlab程序獲取原始NDVI時(shí)間序列曲線的上包絡(luò)線,并計(jì)算3種方法擬合后的NDVI序列對原始高質(zhì)量NDVI值和上包絡(luò)線的擬合效果,評價(jià)3種方法的優(yōu)劣。擬合重建處理中的任何時(shí)刻像元NDVI值依據(jù)datapixelreliability數(shù)據(jù)設(shè)置相應(yīng)的權(quán)重,圖2中像元在3、4、6、10、16、19、20、21期為噪聲點(diǎn)(圓圈標(biāo)識)。圖中3種擬合算法的加權(quán)與未加權(quán)處理的NDVI時(shí)間序列曲線顯示,前者擬合后曲線均有不同程度的提升,因?yàn)槭茉蒲?點(diǎn)10、16)及地面冰雪覆蓋(其他圈出的點(diǎn))的干擾造成NDVI值降低,使重建曲線對上包絡(luò)線擬合不夠,擬合后的NDVI最大值偏低。另外,由于春季返青之前和進(jìn)入枯黃期階段的連續(xù)噪聲,生長季開始階段和結(jié)束階段的擬合曲線都發(fā)生右移,如果未進(jìn)行不同權(quán)重處理,對于植被物候信息(如返青期、生長峰值、枯黃期等)的提取都將引入新的噪聲。本文對3種擬合算法的比較研究都將引入數(shù)據(jù)質(zhì)量權(quán)重系數(shù)進(jìn)行處理。4.23基于自治重建算法的基本思想(1)基于分段高斯函數(shù)的植被生長模擬基于非對稱高斯函數(shù)擬合方法是一個(gè)從局部擬合到整體擬合的過程,使用分段高斯函數(shù)來模擬植被生長過程,最后通過平滑連接各高斯擬合曲線實(shí)現(xiàn)時(shí)間序列重建。其過程大致分為區(qū)間提取、局部擬合和整體連接三步驟。(2)雙邏輯函數(shù)濾波雙Logistic曲線擬合法是由Beck等于2006年提出的一種新算法,與非對稱高斯函數(shù)擬合法類似,雙邏輯函數(shù)濾波也是一種半局部擬合方法。首先將整個(gè)時(shí)間序列中時(shí)間點(diǎn)對應(yīng)的值按極大或極小值分成多個(gè)區(qū)間,分別對該區(qū)間進(jìn)行DoubleLogistic函數(shù)局部擬合,其處理過程與AG方法類似。(3)zki和go爾雅Savitzky-Golay濾波法最初由Savitzky和Golay于1960年提出,是一種移動(dòng)窗口的加權(quán)平均算法,但其加權(quán)系數(shù)不是簡單的常數(shù)窗口,而是通過在滑動(dòng)窗口內(nèi)對給定高階多項(xiàng)式的最小二乘擬合得出。4.3重建前后ndvi時(shí)間序列的特征在基于植被NDVI的各種陸地生態(tài)系統(tǒng)應(yīng)用研究中,通常選取生長季峰期的NDVI值反映全年植被生長狀況。因此有必要比較分析各種算法擬合重建后的NDVI值在7、8月份相對于初始NDVI時(shí)間曲線的上包絡(luò)線的擬合效果。此外,利用NDVI時(shí)間序列提取物候信息及應(yīng)用全年或生長季期的NDVI平均狀況的研究中,選擇一種在去除原始NDVI數(shù)據(jù)集中噪聲的同時(shí),能夠有效地保持質(zhì)量好的NDVI點(diǎn)值的真實(shí)性的算法也十分關(guān)鍵。為了比較AG擬合、D-L擬合與S-G濾波法3種方法對高寒干旱環(huán)境的藏北草地MODIS16d合成的NDVI時(shí)間序列曲線擬合重建的效果,將著重從3個(gè)方面分析:①分析3種算法對整個(gè)區(qū)域NDVI均值在時(shí)間序列上改變的程度,以及擬合的時(shí)間曲線與原始NDVI曲線的相關(guān)系數(shù);②比較各像元擬合后的NDVI時(shí)間序列曲線在成熟期對上包絡(luò)線的擬合效果;③分析全年、生長季內(nèi)、生長季開始和結(jié)束時(shí)段以及生長季峰期重建后的NDVI與高質(zhì)量NDVI原始值的偏離程度。其中,相關(guān)系數(shù)用于反映擬合后的NDVI時(shí)間序列曲線保持原始植被生長季特征的能力,如式(2)所示:r=∑i=123(NDVIpi?NDVIpˉˉˉˉˉˉˉˉˉˉˉ)(NDVIoi?NDVIo)∑i=123(NDVIpi?NDVIpˉˉˉˉˉˉˉˉˉˉˉ)2×∑i=123(NDVIoi?NDVIoˉˉˉˉˉˉˉˉˉˉˉ)2ue001?ue000ue000(2)r=∑i=123(ΝDVΙpi-ΝDVΙpˉ)(ΝDVΙoi-ΝDVΙo)∑i=123(ΝDVΙpi-ΝDVΙpˉ)2×∑i=123(ΝDVΙoi-ΝDVΙoˉ)2(2)其中:NDVIoi、NDVIpi分別是年內(nèi)時(shí)間序列中第i期擬合處理前后的NDVI值,NDVIoˉˉˉˉˉˉˉˉˉˉΝDVΙoˉ、NDVIpˉˉˉˉˉˉˉˉˉˉΝDVΙpˉ分別為處理前后的年內(nèi)NDVI時(shí)間序列的均值。通過計(jì)算每個(gè)樣本像元的擬合重建前后的整年NDVI時(shí)間序列相關(guān)系數(shù),統(tǒng)計(jì)所有像元相關(guān)系數(shù)的平均值,反映3種擬合算法保持原始NDVI時(shí)間曲線整體特征的能力?；谶@3種主流算法的重建前后NDVI時(shí)間序列的擬合效果采用兩種方式表示:一是擬合前后NDVI偏差分布的箱須圖(Box-whiskerPlot)表示法,利用極小值、第一四分位數(shù)、中位數(shù)、第四分位數(shù)和極大值5個(gè)表示兩組NDVI序列殘差分布的離散程度及對稱性等;二是回歸估計(jì)標(biāo)準(zhǔn)差(如式3),描述重建后的NDVI年內(nèi)時(shí)間序列與原始值之間的平均差異程度,表明重建前后NDVI值的代表性強(qiáng)弱,其值越小,擬合值的代表性越強(qiáng)。RMSE=∑i=1N(NDVIpi?NDVIoi)2N??????????????ue001?ue000(3)RΜSE=∑i=1Ν(ΝDVΙpi-ΝDVΙoi)2Ν(3)其中:NDVIoi、NDVIpi分別是年內(nèi)時(shí)間序列中第i期擬合處理前后的NDVI值,N表示質(zhì)量控制文件中標(biāo)識高質(zhì)量(DN=0,1)的NDVI像元樣本數(shù)量。5結(jié)果與分析5.1基本模型的擬合選擇3種算法對2007～2009年全區(qū)300個(gè)樣點(diǎn)像元的NDVI時(shí)間序列曲線擬合重建,均將不同程度地去除各類噪聲,通常情況一定區(qū)域內(nèi)的像元受地面水分或冰雪覆蓋以及大氣云層、氣溶膠影響造成的NDVI值降低的情況遠(yuǎn)多于“假”高值的情形。因此,各種重建方法將使區(qū)域內(nèi)的像元NDVI整體提高,如圖3所示,5條曲線分別為900個(gè)(300×3年)像元NDVI重建前、3種算法擬合后以及像元點(diǎn)原始NDVI序列的上包絡(luò)線取均值的整年時(shí)間序列。其中,3種算法分別擬合重建后的NDVI值在全年23期較原始NDVI時(shí)間序列曲線均不同幅度提升。特別在植被返青(Number8～10)及進(jìn)入枯黃期(Number18～21)時(shí)段,NDVI值重建后的整體水平相對提高11.9%～30.2%。這與實(shí)例區(qū)域的本身特征有關(guān),這兩個(gè)時(shí)段分別是草地覆被NDVI序列發(fā)生突變的“拐點(diǎn)”,同時(shí)也是藏北地區(qū)自然條件發(fā)生轉(zhuǎn)變的時(shí)刻,土壤分別處于開始解凍和發(fā)生凍融作用的階段,受冰雪覆蓋變化影響較大;藏北高原的大部分地區(qū)氣溫、降水以及大氣狀況也是在一年中的5月份和10月份發(fā)生轉(zhuǎn)折。然而,3種擬合重建方法對原始NDVI值的改變程度不一致。D-L曲線擬合與非對稱高斯函數(shù)擬合方法擬合后的結(jié)果比較吻合,僅在生長季峰期D-L擬合法較AG擬合法的NDVI均值略高0.002左右。S-G濾波法的“過度”擬合特征對噪聲比較敏感,其處理的NDVI均值與前兩種方法有較大的差異,在生長季開始和結(jié)束時(shí)間,擬合的曲線相對于原始及前兩種算法重建的NDVI時(shí)間序列分別發(fā)生較大程度的“左偏”和“右移”,這對基于NDVI時(shí)間序列提取植被物候信息(如返青期或枯黃期)將引入新的噪聲。在生長季峰期(Number13～15),3種算法擬合后的NDVI平均值較初始值均有提高,AG函數(shù)與D-L曲線擬合處理結(jié)果甚至超出原始NDVI曲線的上包絡(luò)線,而S-G濾波法重建結(jié)果低于上包絡(luò)線值。圖4所示為900個(gè)全年及生長季內(nèi)(5月9日～10月15日)NDVI時(shí)間序列基于上述3種算法擬合重建前后的相關(guān)系數(shù)。其中非對稱高斯函數(shù)擬合方法重建的時(shí)間序列與植被原始NDVI時(shí)間曲線相關(guān)性最好,全年23期的NDVI點(diǎn)值及生長季內(nèi)NDVI序列的擬合相關(guān)系數(shù)分別為0.791、0.803。D-L擬合法比AG擬合略低,相關(guān)程度分別為0.788、0.798;而S-G濾波法受噪聲影響對整個(gè)時(shí)間序列曲線整體特征的保持性不及前兩種方法。相對于全年的擬合效果,生長季內(nèi)的時(shí)間序列采用任何一種算法重建前后相關(guān)程度都更好,因?yàn)閷?shí)例區(qū)域在生長季以外的時(shí)間內(nèi)地表覆被NDVI受噪聲影響更嚴(yán)重,重建結(jié)果可靠性相對稍低。5.2模型擬合效果比較生長季NDVI最大值(藏北草地時(shí)間為7～8月份)通常反映出一個(gè)地區(qū)植被生長最好的狀況,而NDVI受噪聲污染引起的NDVI值嚴(yán)重偏低將給研究結(jié)論帶來偏差。圖5反映了基于AG擬合、D-L擬合及S-G濾波3種算法重建的NDVI時(shí)間序列在13～15期(7月12日～8月28日)相對于原始NDVI序列上包絡(luò)線的擬合效果比較。圖3中顯示,AG擬合與D-L擬合的結(jié)果超出上包絡(luò)線曲線,S-G濾波法擬合重建結(jié)果低于上包絡(luò)曲線。圖5中不難發(fā)現(xiàn),AG函數(shù)擬合與D-L曲線擬合結(jié)果相近,效果均優(yōu)于S-G濾波法。5.3ag擬合與d-l擬合法處理效果對比圖6所示為900個(gè)(2007～2009年300個(gè)均勻分布于全區(qū)域的樣點(diǎn)像元)全年NDVI時(shí)間序列樣本點(diǎn),基于上述3種算法擬合后的NDVI時(shí)間序列與原始高質(zhì)量點(diǎn)值的偏差的箱須圖,縱坐標(biāo)代表樣本像元NDVI在一年中23期數(shù)據(jù)的重建前后偏差(Residualvalue)分布,橫坐標(biāo)為年內(nèi)NDVI序列的期號(1～23)。正值代表擬合后的值偏高,負(fù)值表示低估了原始NDVI值。由該圖可知,AG擬合與D-L擬合結(jié)果在各時(shí)間點(diǎn)的NDVI偏差程度以及整個(gè)時(shí)間序列的分布模式基本一致,NDVI偏差時(shí)間序列上均為“單峰”,而S-G濾波法擬合的時(shí)間序列相對原始高質(zhì)量NDVI點(diǎn)值的偏差在生長季開始和結(jié)束兩個(gè)時(shí)間段達(dá)到最大,呈現(xiàn)“雙峰”現(xiàn)象。由于生長季NDVI的絕對值大于年內(nèi)開始與結(jié)束時(shí)段,因此AG擬合和D-L擬合法處理后的NDVI偏差在13～15期的合成期內(nèi)取值范圍及離散程度達(dá)到最大,并且D-L擬合法處理結(jié)果殘差取值范圍和離散程度均稍高于AG函數(shù)擬合。S-G濾波法擬合結(jié)果偏差在13～15期范圍低于前兩種方法,但是在植被返青和進(jìn)入休眠時(shí)間段NDVI殘差對水平軸偏離程度更大,即高估擬合了高質(zhì)量的NDVI值,并且分布較離散。圖7表示的是全年、生長季、返青期、成熟期、枯黃期5個(gè)統(tǒng)計(jì)時(shí)段內(nèi)基于3種算法擬合重建前后的NDVI序列回歸標(biāo)準(zhǔn)差。圖中比較分析得出,除了成熟期外其他4個(gè)統(tǒng)計(jì)項(xiàng),AG擬合法與D-L擬合法在保持高質(zhì)量NDVI點(diǎn)值真實(shí)性效果均有優(yōu)勢;而S-G濾波法擬合生長季旺盛期時(shí)段,對高質(zhì)量NDVI樣本點(diǎn)其擬合結(jié)果較前兩種方法更有代表性。5.4基于合成期的ndvi重建效果圖8所示為DOY(從一年中1月1日開始計(jì)算的年序日,Dateofyear)209～224合成期的MODISNDVI重建前后對比,總體上,NDVI時(shí)間序列重建最大程度地去除地面狀況、大氣條件造成NDVI突降和假高情況,從而使區(qū)域內(nèi)NDVI的異質(zhì)性減弱,植被覆蓋的空間連續(xù)性變好,從圖8(a)、(b)中的1、2、3矩形框標(biāo)識均能得到明顯反映,不同區(qū)域的像元值均有提升。6ag與d-l擬合法重建的ndvi時(shí)間序列擬合NDVI時(shí)間序列擬合重建在消除各種噪聲值的同時(shí),也可能使原始NDVI真實(shí)值發(fā)生偏離。由于通常噪聲影響的NDVI值使其低