時序ndvi數(shù)據(jù)集重建方法綜述

1時序ndvi數(shù)據(jù)的重建算法植被直接反映了自然環(huán)境的條件，是氣候和人文因素影響環(huán)境的敏感指標。因此，它在地球系統(tǒng)中發(fā)揮著重要作用。根據(jù)植被的光譜特性,將衛(wèi)星可見光和近紅外波段進行組合,形成多種植被指數(shù),可以簡單、有效地度量地表植被狀況。其中NDVI(歸一化植被指數(shù))作為植被生長狀態(tài)和植被覆蓋度的最佳指示因子,其時序變化反映植被的生長變化過程,已被廣泛的應用在全球與區(qū)域環(huán)境變化研究、植被動態(tài)變化監(jiān)測、葉面積指數(shù)等植物生物物理量參數(shù)反演以及凈第一性植被生產力研究等許多方面。目前常用的時序NDVI數(shù)據(jù)集是10d最大值合成的產品,例如:PAL(PathfinderAVHRRLand)、GIMMS(GlobalInventoryModelingandMappingStudies)、SPOT_VEGETATION產品等,這些數(shù)據(jù)集雖然經過嚴格的預處理,但由于受到云污染、大氣的不確定性及地表的雙向反射等影響,殘留很多噪聲。因此,需要通過一系列方法降低時序NDVI數(shù)據(jù)中的噪聲水平,為各種研究提供更為可靠的時序數(shù)據(jù)集。迄今為止,國內外學者相繼提出了很多時序NDVI數(shù)據(jù)集的重建算法,例如:最佳指數(shù)斜率提取法、均值迭代濾波法、傅立葉變換法、Savitzky-Golay濾波法、非對稱高斯擬合法、小波分析法以及時間窗內的線性內插法等。這些算法都很好地降低了序列中的噪聲水平,但也存在各自的局限性,因此深入研究時序NDVI數(shù)據(jù)的重建算法,具有十分重要的現(xiàn)實意義。本文介紹了其中6種方法,并以張掖地區(qū)2007年和2008年10d最大值合成的SPOT/VGT的時序NDVI數(shù)據(jù)為例,運行這些方法,對重建結果進行了比較和評價;最后用這些方法對人為的噪聲序列進行重建,進而對重建方法的優(yōu)缺點進行評價。2數(shù)據(jù)的突變點不符合植被生長的漸進過程時序數(shù)據(jù)的重建算法都基于以下兩個假設:①NDVI的時序變化對應著植被的生長和衰老過程;②云和不利的大氣條件往往會降低NDVI的值,時序數(shù)據(jù)中的突變點不符合植被生長的漸進過程,應作為噪聲點修正。最大值合成法(MVC)是目前絕大多數(shù)NDVI數(shù)據(jù)集所采用的合成方法,該方法是在指定的合成周期內選取最大值作為合成數(shù)據(jù)。本文采用的實例數(shù)據(jù)就是周期為10d的最大值合成數(shù)據(jù)。通過MVC方法可以在一定程度上去除云和大氣的干擾,但是,如果一個周期時間內所有圖像都有云,就無法得到該周期的無云圖像。因此,僅僅使用最大值合成法是不夠的,還必須在最大值合成后做進一步的處理。2.1bise算法檢測Viovy提出的最佳指數(shù)斜率提取法(BISE)基于植被變化的可預測性,有別于變化頻繁的視場角和云,植被變化較為穩(wěn)定。該算法是在設定的滑動周期內,向前搜索,如果下一個點的值高于第一個點,則接受該點;如果低于第一個點,接受這個點必須要滿足一個條件:在這個周期內沒有一個點的值大于這個低值點與前一個高值點之差的20%;最后,忽略與前一個點相比增大0.1的點。BISE算法適用于每天的NDVI數(shù)據(jù),因此需要對該算法進行調整才能有效地應用于10d合成的數(shù)據(jù)。對于異常高值的檢測,由于10d合成的數(shù)據(jù),連續(xù)兩個點的時間是從1d跨越到20d,所以一個簡單的閾值檢測并不恰當。另外,當時間序列呈長期遞減趨勢時,BISE算法處理效果不理想。于是Lovell等提出了一種改進的最佳指數(shù)斜率提取算法(ModifiedBISE),提出需要考慮數(shù)據(jù)在局部的變化率,該改進的算法很好地降低了時序數(shù)據(jù)中的噪聲水平。2.2閾值提取和迭代均值迭代濾波法是一種基于濾波的平滑方法,具體處理過程為:在NDVI序列中尋找一個點,其前后相鄰兩個點的平均值與該點的值的差值最大,且該差值超過閾值。如果存在滿足該條件的點,就用前后兩點的平均值替換該點的值。迭代進行,直到NDVI序列中不再存在滿足條件的點。在實際運用中,該閾值可以取NDVI多年平均值的某一個特定的比例。均值迭代濾波法具有算法實現(xiàn)簡單,平滑效果好的優(yōu)點,能夠有效地濾去NDVI序列中的高頻變化。2.3時域信號轉換模塊t傅立葉變換是處理周期性時間序列圖像十分有效的分析方法。傅立葉原理表明:任何連續(xù)測量的時序或信號,都可以表示為不同頻率的正弦波信號的無限疊加。如果輸入的f(x)為離散的序列,其計算公式為:F(u)=1ΝΝ-1∑x=0f(x)e-i2πΝux(1)F(u)=1N∑x=0N?1f(x)e?i2πNux(1)其中:f(x)為時域中的序列,F(u)為傅里葉變換,N為序列長度,x為時域的序列號,u為頻域的序列號(u=0,1,…,N-1);u/N為對應分量正弦曲線的頻率。傅立葉變換是將原來難以處理的時域信號轉換成易于分析的頻域信號,最后再利用傅立葉逆變換將這些頻域信號轉換成時域信號。頻域中的低頻的部分代表背景,而高頻的部分代表隨機成分,即噪聲。去除頻域中的高頻部分就代表去除了原NDVI序列中的噪聲部分,從而得到平滑的NDVI序列。2.4s-g濾波Savizky-Golay濾波是一種基于最小二乘的卷積算法,它是用一定長度的過濾器和待處理數(shù)據(jù)作卷積,對待處理的數(shù)據(jù)作加權多項式擬和,擬合的目標是求得最小均方根誤差,而一些遠離大多數(shù)點的邊沿點不會參與擬合。這樣擬合時過于偏離正常生長趨勢線的噪聲部分會被丟棄。S-G濾波的公式為:Y*j=i=m∑i=-mCi×Yj+1Ν(2)Y?j=∑i=?mi=mCi×Yj+1N(2)其中:Y是原始的NDVI數(shù)值,Y*是平滑過后的NDVI值,Ci是從過濾器首部開始第i個NDVI值的權值,N是過濾器的長度,同時也可以表示為(2m+1)。因此S-G濾波需要人為設定兩個參數(shù):第一個是m,第二個是多項式擬合的階數(shù)。通常m值越大結果越平滑,被平滑的峰谷值也就越多;多項式擬合的階數(shù)一般設在2～4之內,較低的次數(shù)可以得到更平滑的結果,但會引進誤差,較高的次數(shù)可以降低這個誤差,但是可能過于擬合得到一個更多噪聲的結果。2.5a為高斯函數(shù)非對稱高斯函數(shù)擬合法是基于不對稱高斯函數(shù)的非線性最小二乘擬合算法,擬合函數(shù)的基本形式為:f(t)=f(t;c1,c2,a1,…,a5)=c1+c2g(t;a1,…,a5)(3)其中:g(t;a1,…,a5)為高斯函數(shù)。g(t;a1,?,a5)={exp[-(t-a1a2)a3],ift>a1exp[-(a1-ta4)a5],ift<a1(4)a1,?,a5)=?????exp[?(t?a1a2)a3],ift>a1exp[?(a1?ta4)a5],ift<a1(4)其中:c1和c2控制曲線的基準和幅度;a1決定最大值或最小值的位置;a2、a3、a4、a5控制曲線的寬度和陡峭度。當進行完第一次擬合后,如果序列中的NDVI值低于擬合值,可以降低這些點的權重,再進行一次擬合,從而使擬合的曲線更加接近NDVI數(shù)據(jù)的上包絡線。2.6傅立葉序列的計算HANTS是基于傅里葉變換改進的一種有效的算法,它將時間序列分解為平均值和有限個不同頻率正弦函數(shù)(諧波),其表達式為:yi=A0+m∑j=1Ajsin(ωji+θj)yi=A0+∑j=1mAjsin(ωji+θj);i=0,1,…,N(5)其中:A0為諧波的余項,等于序列的平均值;Aj為各諧波的振幅;ωj=2jπ/N為各諧波的頻率,N為序列的長度;θj為各諧波的初相位;m為諧波的個數(shù),等于N-1。這些正弦函數(shù)疊加構成的數(shù)據(jù)序列,即為傅立葉序列。傅立葉序列也可表示為一個正弦函數(shù)和余弦函數(shù)之和,表達式為:yi=a0+m∑j=1(ajcosωji+bjsinωji)(6)其中:a0為序列的平均值,a、b為傅立葉參數(shù)。HANTS算法是根據(jù)幾個最為顯著頻率的傅立葉序列函數(shù)進行最小二乘擬合的過程。3ndvi值的定義本文采用的實例數(shù)據(jù)是SPOT/VEGETATION(VGT)產品。目前,VGT產品主要分為兩類:P產品,提供的輔助數(shù)據(jù)可讓使用者根據(jù)自己的算法來處理原始大氣反射值;S產品,該產品經過幾何校正和大氣校正,提供1天綜合產品(S1)和10d綜合產品(S10)。本文選用的是張掖地區(qū)(97°20′～102°12′E與37°28′～39°57′N)2007年1月到2008年12月的S10產品,其空間分辨率為1km。真實的NDVI值范圍是-1.0～1.0,S10產品是0～255的DN(DigitalNumbers)數(shù)據(jù),從DN值轉化成NDVI值的關系式為:NDVI=0.004×DN-0.1(7)試驗區(qū)范圍內主要的土地類型是:耕地、林地、草地、荒漠和水體。圖1是選取一個草地像元的原始NDVI序列,序列中存在明顯的突降和突增的點。3.1不同方法重建的ndvi時間序列MVI算法唯一的參數(shù)是閾值,閾值的選擇影響重建的效果。作者對比了分別取多年平均值的10%、15%和20%作為閾值的重建結果,發(fā)現(xiàn)取20%無法修正所有的噪聲點;選10%時,曲線被過度平滑,序列中大部分的點都被改變;15%的效果比較好,因此本文中MVI方法的閾值取NDVI多年平均值的15%。改進的BISE算法對滑動窗口的大小和閾值的選取比較敏感,窗口過大會損失一些細節(jié)信息,過小則仍包含很多噪聲。本文分別選用3、4、5和6作為周期,通過對比發(fā)現(xiàn),當周期為5時,修正的噪聲點最多,重建的效果最好。文中傅立葉變換算法使用快速傅立葉變換(FFT),FFT的一個參數(shù)是濾波的截止頻率,選取的截止頻率越小,擬合曲線越光滑,但擬合的諧波數(shù)越少,和原始點的偏離越大。Savitzky-Golay濾波算法有兩個限制因素,可以通過觀察NDVI序列確定。第一個參數(shù)是濾波窗口的一半寬度m,通常m值越大,容易忽略NDVI時序中的一些重要變化信息,而過小的m值,會造成“過于平滑”而帶來更多的噪聲,因此m值范圍一般取4～7;另一個參數(shù)是多項式最小二乘擬合的次數(shù)d,基本上取2～4之內。一系列試驗結果對比發(fā)現(xiàn),m和d取4×3的組合適合本研究區(qū)域。非對稱高斯函數(shù)擬合方法首先需要用滑動窗口均值方法對NDVI曲線進行初步平滑,因為運用傅立葉分析確定生長周期數(shù)時需要基于低噪聲的數(shù)據(jù)。最后用擬合函數(shù)對每一個周期進行一次加權最小二乘擬合。對于NDVI值低于擬合函數(shù)的點,降低其權重,進行第二次擬合。HANTS算法擬合過程需要設定3個參數(shù):頻率個數(shù)n、擬合偏差閾值和額外的數(shù)據(jù)點數(shù)m。n值越大擬合偏差越小,但同時也會使擬合曲線出現(xiàn)不規(guī)則的起伏引入更多的噪聲,本文經過不同n值的對比試驗發(fā)現(xiàn)n取7效果較好。擬合偏差是指原始NDVI值低于擬合值時,兩者差值的絕對值。當擬合偏差超過這個閾值,就把該點剔除并賦零值。做最小二乘擬合時,至少需要2n-1個點才能完成擬合函數(shù)的系數(shù)計算,而用于擬合的點數(shù)越多,效果越好,m用于指定除2n-1個必須點數(shù)以外所使用的數(shù)據(jù)點數(shù)。圖2是選取一些典型例子的原始的和經過各個方法重建之后的NDVI時間序列對比圖(6個方法分別在左右兩列圖中顯示)。分別選取了耕地、林地、草地和沙地4種不同的土地類型。從圖中可以明顯的看出,FFT方法、HANTS算法和非對稱高斯函數(shù)(AG)擬合法重建的曲線比較光滑,但其中FFT重建的結果與原始序列有較大的偏離,這是由于FFT算法是將曲線中高頻部分直接去除,只留下低頻部分,損失了一部分信息。并且原始的時間序列誤差越大、越不規(guī)則,FFT重建的結果偏離原始序列越大。因此,FFT方法給原始數(shù)據(jù)帶來了較大的誤差。HANTS是FFT的改進算法,從左邊圖中可以明顯的看出HANTS算法的重建效果比FFT方法的效果好,它與原始曲線的趨勢大致相同,但改變了序列中大部分點的值。AG擬合法對生長周期明顯的植被類型重建效果較好,而對生長周期不明顯的植被類型處理效果不好,例如圖2右列圖中的沙地圖,AG方法處理的結果較差。MVI、改進的BISE和S-G濾波3種方法重建的結果基本相同,序列中大部分噪聲都被識別和修正,重建效果較好。MVI方法識別了序列中大部分的異常低值點和高低值,但對連續(xù)的噪聲點判斷不夠,例如耕地圖,削弱了正常的峰值點;改進的BISE算法和S-G濾波算法重建的結果都是原始序列的上包絡線,因此無法糾正序列中的異常高值點,例如耕地圖和林地圖中異常的高值點都沒有識別。另外,HANTS算法、改進的BISE算法和FFT算法需要根據(jù)人為經驗不停地嘗試不同的閾值或者周期的選擇,造成一定的困難。在算法運行速度上,MVI算法和改進的BISE算法很快,而AG擬合法和FFT算法需要耗費較長的時間。3.2噪聲序列的重建Jacknife方法可以用來評價NDVI序列的重建效果。其基本原理是:在一個時間序列中,隨機選取一些點,人為加入噪聲,噪聲可以取原值的10%、20%或者更大。經過重建之后,可以直觀的看出各個重建方法的效果。圖3是一組人為噪聲序列經過6種方法重建之后的對比圖,它是在一個假想的光滑曲線上人為放置一些突降點和異常高值點。通過Jacknife方法可以更直觀地看出各個方法的效果,從圖3(a)可以明顯的看出,各個方法都很好的降低了序列中的噪聲水平,其中均值迭代濾波法(MVI)重建后的曲線與真實的NDVI幾乎重疊,能識別序列中大部分的異常點,只有連續(xù)的異常高值點無法識別;改進的BISE算法和S-G濾波法重建的結果相似,噪聲序列中突降點都被識別,但是異常高值點無法判斷;HANTS算法和AG擬合法重建得到的曲線比較光滑,且與假想的NDVI曲線趨勢大致相同,但兩者與假想的曲線都有一定的偏差,改變了大部分的NDVI值。FFT重建得到的曲線也比較光滑,但與假想序列偏離太大,引入了更多的噪聲。4不同算法的對比研究NDVI時間序列數(shù)據(jù)由于云和大氣條件等造成污染,是全球環(huán)境變化等研究面臨的一個重要問題。因此近年來學者設計了個各種算法重建高質量的時序NDVI數(shù)據(jù)。本文介