基于秦嶺樣區(qū)的四種時(shí)序EVI函數(shù)擬合方法對(duì)比研究

劉亞南，肖　飛，杜　耘

1 中國(guó)科學(xué)院測(cè)量與地球物理研究所，武漢　430077 2 湖北省環(huán)境與災(zāi)害監(jiān)測(cè)評(píng)估重點(diǎn)實(shí)驗(yàn)室，武漢　430077 3 中國(guó)科學(xué)院大學(xué)，北京　100049

?

基于秦嶺樣區(qū)的四種時(shí)序EVI函數(shù)擬合方法對(duì)比研究

劉亞南1,2,3，肖飛1,2,*，杜耘1,2

1 中國(guó)科學(xué)院測(cè)量與地球物理研究所，武漢430077 2 湖北省環(huán)境與災(zāi)害監(jiān)測(cè)評(píng)估重點(diǎn)實(shí)驗(yàn)室，武漢430077 3 中國(guó)科學(xué)院大學(xué)，北京100049

函數(shù)曲線擬合方法是植被指數(shù)時(shí)間序列重建的一個(gè)重要方法，已經(jīng)廣泛應(yīng)用于森林面積動(dòng)態(tài)變化監(jiān)測(cè)、農(nóng)作物估產(chǎn)、遙感物候信息提取、生態(tài)系統(tǒng)碳循環(huán)研究等領(lǐng)域?；谇貛X樣區(qū)多年MODIS EVI遙感數(shù)據(jù)及其質(zhì)量控制數(shù)據(jù)，探討并改進(jìn)了時(shí)序EVI重建過(guò)程中噪聲點(diǎn)優(yōu)化和對(duì)原始高質(zhì)量數(shù)據(jù)保真能力的評(píng)價(jià)方法；在此基礎(chǔ)上，比較了常用的非對(duì)稱性高斯函數(shù)擬合法(AG)、雙Logistic函數(shù)擬合法(DL)和單Logistic函數(shù)擬合法(SL)?；赟L方法，調(diào)整了模型形式并重新定義d的參數(shù)意義，提出了最值優(yōu)化單Logistic函數(shù)擬合法(MSL)，并與其他3種方法進(jìn)行對(duì)比。結(jié)果表明；在噪聲點(diǎn)優(yōu)化及保留原始高質(zhì)量數(shù)據(jù)方面，AG方法和DL方法二者整體差別不大，而在部分像元的處理上AG方法表現(xiàn)出更好的擬合效果；MSL方法和SL方法相比于AG方法和DL方法其效果更為突出；在地形氣候復(fù)雜，植被指數(shù)噪聲較多的山區(qū)，MSL方法表現(xiàn)出更好的適用性。

MODIS EVI；曲線擬合；時(shí)間序列重建；Logistic；秦嶺

植被指數(shù)數(shù)據(jù)已經(jīng)廣泛應(yīng)用于農(nóng)情監(jiān)測(cè)、物候信息提取、碳通量估算、陸地生態(tài)系統(tǒng)對(duì)氣候變化的響應(yīng)等研究領(lǐng)域[1- 3]。然而，遙感獲取的植被指數(shù)時(shí)間序列曲線包含很多噪聲點(diǎn)，往往對(duì)研究精度造成很大影響[4- 5]。為有效濾除噪聲，學(xué)者發(fā)展了一系列植被指數(shù)時(shí)間序列重建方法。主要包括：NDVI閾值法[6- 7]、后向移動(dòng)平均法[8]、最大上升速率判斷方法[9]、曲線擬合法[10]和經(jīng)驗(yàn)回歸方程法[11]，動(dòng)態(tài)權(quán)重濾波法[12]。其中，函數(shù)曲線擬合法針對(duì)每一像元進(jìn)行單獨(dú)處理，利用函數(shù)方程對(duì)植被指數(shù)時(shí)間序列曲線進(jìn)行重建，不需設(shè)定閾值或經(jīng)驗(yàn)系數(shù)，故在不同環(huán)境條件中得到廣泛應(yīng)用[13- 16]。有學(xué)者研究發(fā)現(xiàn)函數(shù)曲線擬合方法總體上優(yōu)于濾波方法[17- 18]。

函數(shù)曲線擬合最常用的方法包括：非對(duì)稱性高斯函數(shù)擬合法(AG)、雙Logistic函數(shù)擬合法(DL)和單Logistic函數(shù)擬合法(SL)。Jonsson等利用非對(duì)稱高斯函數(shù)擬合法對(duì)西非地區(qū)時(shí)序NDVI數(shù)據(jù)進(jìn)行重建，發(fā)現(xiàn)其效果要優(yōu)于濾波方法[17]。雙 Logistic 函數(shù)擬合法(DL)對(duì)于生長(zhǎng)季較短的高緯度地區(qū)具有很好的適用性，但是對(duì)中緯度地區(qū)的植被指數(shù)的重建效果需要進(jìn)一步驗(yàn)證[15]。單Logistic函數(shù)擬合方法針對(duì)各個(gè)生長(zhǎng)季的不同特點(diǎn)分段擬合，對(duì)于生長(zhǎng)季較長(zhǎng)或年內(nèi)具有多個(gè)生長(zhǎng)季的地區(qū)也能適用[10]。

不同函數(shù)曲線擬合方法對(duì)不同地理環(huán)境有不同的適應(yīng)性。尤其山地區(qū)域，因其地形復(fù)雜、數(shù)據(jù)噪聲大，不同函數(shù)曲線擬合方法往往結(jié)果差異較大。因而，在區(qū)域研究中，通常需要對(duì)比選擇一種適用于區(qū)域植被指數(shù)波動(dòng)特征的曲線擬合方法。

目前，擬合方法的對(duì)比研究大多從直觀視覺(jué)角度出發(fā)來(lái)判別噪聲點(diǎn)及其擬合重建效果。傳統(tǒng)的求算擬合相關(guān)系數(shù)和回歸估計(jì)標(biāo)準(zhǔn)差(RMSE)的方法，只能表現(xiàn)對(duì)原始數(shù)據(jù)的擬合程度，無(wú)法定量表征噪聲點(diǎn)的影響程度，因而也無(wú)法充分說(shuō)明重建后的曲線保留原始高質(zhì)量數(shù)據(jù)的能力[18- 21]。

本文基于秦嶺樣區(qū)MODIS EVI多年遙感數(shù)據(jù)及其質(zhì)量評(píng)價(jià)數(shù)據(jù)，探討并發(fā)展了時(shí)序NDVI重建過(guò)程中噪聲點(diǎn)優(yōu)化和原始高質(zhì)量數(shù)據(jù)保真能力的評(píng)價(jià)方法；對(duì)比了常用的 AG方法、DL方法和SL方法的擬合效果，并將上述三種方法和我們基于山區(qū)植被指數(shù)波動(dòng)特征改進(jìn)的MSL方法做了比較研究。

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

本文選取了均勻分布于秦嶺的10個(gè)樣區(qū)進(jìn)行分析。秦嶺為濕潤(rùn)地區(qū)和半濕潤(rùn)地區(qū)的分界線、地形起伏較大、植被時(shí)空分異明顯，以其作為樣區(qū)具有較好典型性和代表性。樣區(qū)位置如圖1所示。

研究使用的實(shí)例數(shù)據(jù)選擇EOS/Terra衛(wèi)星的MODIS MOD13Q1產(chǎn)品，其中包括基于MVC方法16d合成的250m分辨率EVI及其質(zhì)量控制數(shù)據(jù)。與EVI相比，NDVI 在植被生長(zhǎng)旺盛期容易達(dá)到飽和[20]，EVI比NDVI的動(dòng)態(tài)變幅更大[21]，能更真實(shí)地反映植被的生長(zhǎng)變化過(guò)程。數(shù)據(jù)來(lái)源于中國(guó)科學(xué)院計(jì)算機(jī)網(wǎng)絡(luò)信息中心地理空間數(shù)據(jù)云(http://www.gscloud.cn)。該產(chǎn)品經(jīng)過(guò)幾何校正和大氣校正。EVI值的標(biāo)準(zhǔn)范圍為-1.0—1.0。MOD13Q1 EVI產(chǎn)品是-3000—10000的DN值，-3000為填充值，從DN值轉(zhuǎn)化成 EVI 值的關(guān)系式為:

NDVI=0.0001·DN

(1)

秦嶺樣區(qū)是利用MOD13Q1的兩景圖像通過(guò)MRT拼接裁剪和投影變換得到。此外，還需從MOD13Q1的兩景EVI產(chǎn)品中提取 data pixel reliability數(shù)據(jù)集，其空間、時(shí)間分辨率與EVI數(shù)據(jù)集匹配。該數(shù)據(jù)對(duì)區(qū)域每個(gè)像元EVI數(shù)據(jù)質(zhì)量進(jìn)行說(shuō)明，分為5個(gè)不同等級(jí)，其中DN值等級(jí)為0的數(shù)據(jù)質(zhì)量最高。本研究中將DN值等級(jí)為0和1的符合質(zhì)量要求的像元作為擬合的原始數(shù)據(jù)，其他則全部設(shè)為NAN不參與計(jì)算。所以可能出現(xiàn)某些期EVI數(shù)據(jù)是空值的情況。樣區(qū)2011—2013的原始EVI數(shù)據(jù)，共8328個(gè)像元(2776×3年)。

圖1　秦嶺樣區(qū)位置分布圖Fig.1　Distribution of the sample sites in the Qinling Mountains

2　4種函數(shù)曲線擬合方法

2.1常用函數(shù)曲線擬合方法介紹

非對(duì)稱性高斯函數(shù)擬合法(AG)、雙Logistic函數(shù)擬合法(DL)和單Logistic函數(shù)擬合法(SL)是目前最常用的函數(shù)曲線擬合方法。

非對(duì)稱高斯函數(shù)擬合方法(AG)是一個(gè)從局部擬合到整體擬合的方法過(guò)程, 使用分段高斯函數(shù)來(lái)模擬植被生長(zhǎng)過(guò)程, 最后通過(guò)平滑連接各高斯擬合曲線實(shí)現(xiàn)時(shí)間序列重建[17]。其過(guò)程大致可分為3 步驟：區(qū)間提取、局部擬合和整體連接[20]。

Pieter S.A.等，于2006年研究高緯度地區(qū)植被物候提取問(wèn)題提出雙Logistic曲線擬合法[13]。首先，將整個(gè)時(shí)間序列中時(shí)間點(diǎn)對(duì)應(yīng)的值按極大或極小值分成多個(gè)區(qū)間, 分別對(duì)該區(qū)間進(jìn)行雙 Logistic 函數(shù)局部擬合，最后再和AG一樣，進(jìn)行整體連接。

張曉陽(yáng)[10]等提出了利用單Logistic函數(shù)擬合方法。單Logistic函數(shù)是一種分段式Logistic函數(shù)擬合的方法，利用擬合曲線曲率變化的特點(diǎn)，確定EVI時(shí)序曲線上植物各物候轉(zhuǎn)換期。具體形式為：

(2)

式中，y(t)是t時(shí)刻的EVI值，a、b為擬合參數(shù)，d為植被指數(shù)初始背景值，c+d為最大值。此方法針對(duì)各個(gè)生長(zhǎng)季的不同特點(diǎn)分段擬合，對(duì)于生長(zhǎng)季較長(zhǎng)或年內(nèi)具有多個(gè)生長(zhǎng)季的地區(qū)也能適用，而且有十分廣泛的應(yīng)用[13-14,22-23]。Logistic模型方法同時(shí)被NASA(National Aeronautics and Space Administration)作為計(jì)算其MODIS全球物候數(shù)據(jù)產(chǎn)品MGLCD(MODIS Global Land Cover Dynamics Product)的方法之一。

2.2最值優(yōu)化單Logistic函數(shù)擬合法(MSL)

原單Logistic模型中，參數(shù)d作為植被指數(shù)初始背景值，其通常通過(guò)植被指數(shù)時(shí)間數(shù)據(jù)序列直接計(jì)算得出，一般為植被指數(shù)時(shí)間序列中的最小值[21]。然而，由于植被指數(shù)最小值處于植被休眠期內(nèi)(生長(zhǎng)季初始或生長(zhǎng)季末期)，此時(shí)地表植被覆蓋度很低，地表下墊面背景對(duì)數(shù)據(jù)干擾較大，其時(shí)空變化往往導(dǎo)致植被指數(shù)最小值產(chǎn)生較大波動(dòng)，這是導(dǎo)致模型擬合精度降低的原因之一。在分析了秦嶺地區(qū)多年EVI最值波動(dòng)特征的基礎(chǔ)上，本文發(fā)現(xiàn)其最大值相對(duì)其最小值更為穩(wěn)定。

基于上面的原因，本文對(duì)模型形式進(jìn)行改動(dòng)，并重新定義d的參數(shù)意義。

(3)

式中,d為最大EVI值；a、b、c是待擬合參數(shù)，其意義同前。由于EVI數(shù)據(jù)序列最大值相對(duì)較為穩(wěn)定，使得上述改進(jìn)的模型形式希望能達(dá)到較好的擬合效果。

3　基于秦嶺樣區(qū)4種方法的擬合效果對(duì)比分析

3.1時(shí)序植被指數(shù)擬合效果評(píng)價(jià)方法探討

對(duì)呈鋸齒狀波動(dòng)的植被指數(shù)平滑重建是為了最大程度消除噪聲點(diǎn)的影響，同時(shí)保留高質(zhì)量原始數(shù)據(jù)，還原植被生長(zhǎng)的季節(jié)性變化特征。然而不同的重建方法差別很大，如何客觀定量的比較不同方法的優(yōu)劣是十分重要的問(wèn)題。上包絡(luò)線分析法、擬合曲線特征直觀比較法、擬合指標(biāo)評(píng)價(jià)法和相關(guān)統(tǒng)計(jì)量分析法是目前常用的幾種評(píng)價(jià)方法。在實(shí)際研究中往往需要幾種方法綜合使用，如Jonsson等參考原始時(shí)序NDVI數(shù)據(jù)上包絡(luò)線，同時(shí)比較不同擬合曲線形態(tài)特征得出AG方法總體上要優(yōu)于兩種濾波方法[17]。上包絡(luò)線分析法及擬合曲線特征直觀比較法往往從直觀的視覺(jué)角度出發(fā)來(lái)衡量不同擬合方法對(duì)噪聲點(diǎn)的擬合重建效果，主觀干擾較強(qiáng)，因而有些學(xué)者輔助其他方法來(lái)進(jìn)行定量分析。如宋春橋等除了比較擬合曲線和原始數(shù)據(jù)的上包絡(luò)線之外，還利用擬合的相關(guān)系數(shù)和回歸估計(jì)標(biāo)準(zhǔn)差對(duì)AG方法、DL方法和Savitzky-Golay濾波方法進(jìn)行了比較研究[18]。曹云鋒等引入NDVI均值、NDVI平均絕對(duì)誤差和相對(duì)誤差等統(tǒng)計(jì)量來(lái)比較分析3種濾波算法[21]。朱文泉等，基于模擬的真實(shí)的EVI時(shí)序曲線再混進(jìn)不同程度的噪聲來(lái)評(píng)價(jià)不同擬合方法抗噪聲干擾的能力[12]。

通過(guò)比較原始數(shù)據(jù)上包絡(luò)線和不同擬合曲線的特征，能夠準(zhǔn)確判斷不同擬合曲線對(duì)原始數(shù)據(jù)的擬合重建效果，對(duì)噪聲點(diǎn)的敏感性及其優(yōu)化噪聲點(diǎn)的能力。然而，這種比較只局限于個(gè)別像元的對(duì)比，缺乏客觀的定量指標(biāo)來(lái)衡量不同擬合方法對(duì)整個(gè)研究區(qū)域所有像元的擬合能力；同時(shí)此方法對(duì)噪聲點(diǎn)的正確判讀依賴很大，容易受到研究者主觀因素的影響，將個(gè)別突降的低值判定為噪聲點(diǎn)，卻忽略了連續(xù)噪聲點(diǎn)的情況。有些學(xué)者選用了相關(guān)系數(shù)和回歸估計(jì)標(biāo)準(zhǔn)差(RMSE)等指標(biāo)來(lái)對(duì)不同曲線擬合效果進(jìn)行定量分析[14],相關(guān)系數(shù)表明兩組數(shù)據(jù)之間的相關(guān)強(qiáng)度；而RMSE可以表明兩組數(shù)據(jù)之間平均差異程度。相對(duì)于NDVI均值、NDVI平均絕對(duì)誤差和相對(duì)誤差等統(tǒng)計(jì)量，相關(guān)系數(shù)和RMSE能更好的衡量不同擬合曲線的擬合效果。通過(guò)對(duì)不同方法擬合區(qū)域每個(gè)像元得到的相關(guān)系數(shù)和RMSE再分別求取均值，或者進(jìn)行差值比較，可以從整體上評(píng)價(jià)不同方法對(duì)研究區(qū)的適用性。然而，如果只簡(jiǎn)單將擬合前后的數(shù)據(jù)進(jìn)行組合，求算得到的相關(guān)系數(shù)和RMSE則不能突出擬合后的數(shù)據(jù)與高質(zhì)量原始數(shù)據(jù)之間的相關(guān)程度，無(wú)法充分說(shuō)明重建后的曲線保留原始高質(zhì)量數(shù)據(jù)的能力。

本文利用質(zhì)量控制數(shù)據(jù)定位噪聲點(diǎn)(質(zhì)量級(jí)別為1)，而不是主觀判定突降的低值為噪聲點(diǎn)，能相對(duì)客觀的比較不同擬合曲線對(duì)噪聲的重建效果。在此基礎(chǔ)上，用去除噪聲點(diǎn)的高質(zhì)量原始數(shù)據(jù)(質(zhì)量級(jí)別為最高0)和與其相對(duì)應(yīng)的擬合重建后的EVI數(shù)據(jù)進(jìn)行組合求算相關(guān)系數(shù)和回歸估計(jì)標(biāo)準(zhǔn)差(RMSE)，因而求算出的相關(guān)系數(shù)和RMSE能夠充分說(shuō)明重建后的曲線與原始高質(zhì)量數(shù)據(jù)的相關(guān)性及差異性。

3.24種方法擬合效果的直觀比較

基于上述方法，本文從秦嶺10個(gè)樣區(qū)中抽取4種方法擬合效果差別較大的6個(gè)典型像元對(duì)EVI時(shí)間序列重建前后的曲線進(jìn)行對(duì)比分析(圖2)。

圖2　2011—2013年6個(gè)生長(zhǎng)季的典型像元EVI(Enhanced Vegetation Index) 時(shí)間序列重建前后曲線對(duì)比Fig.2　Comparison of the original and reconstructed MODIS EVI time-series during 2011—2013

TIME表示影像期號(hào)，MOD13Q1數(shù)據(jù)是16d合成，一年23期。圖2中可以看出，4種法都不同程度的去除了鋸齒波動(dòng)，重建后的EVI曲線更加平滑。4種方法都會(huì)受到連續(xù)噪聲點(diǎn)(質(zhì)量級(jí)別為1的BAD EVI)的影響，如2012上半年的1、2、3、5期和2013上半年年的3、4、5期。

總體上，AG方法和DL方法的擬合曲線比較接近，但都對(duì)噪聲點(diǎn)比較敏感，曲線形態(tài)容易受到噪聲點(diǎn)的影響發(fā)生變化。這種影響大多發(fā)生在生長(zhǎng)季的初期和末期。但是在生長(zhǎng)季的中期也會(huì)受到影響。如2012下半年的第18(9月份)期和2013上半年的第14期(7月份)。MSL方法和SL方法的擬合曲線有較好的一致性；相比于AG方法和DL方法，MSL方法和SL方法對(duì)高質(zhì)量EVI原始數(shù)據(jù)(質(zhì)量級(jí)別為0的 EVI數(shù)據(jù))有更好的重建效果。

3.34種方法擬合效果的定量分析

本文利用4種方法擬合3年6個(gè)生長(zhǎng)季EVI數(shù)據(jù)(6×2776個(gè)像元)，并用去除噪聲點(diǎn)的高質(zhì)量原始數(shù)據(jù)(質(zhì)量級(jí)別為最高0)和與其相對(duì)應(yīng)的擬合重建后的數(shù)據(jù)進(jìn)行組合求算相關(guān)系數(shù)和RMSE，對(duì)比其各生長(zhǎng)季的整體平均值，以此反映4種方法對(duì)不同生長(zhǎng)季EVI 重建保真性能力(表1,表2)。

表1和表2中可以看出，MSL方法擬合的相關(guān)系數(shù)各期都高于其他3種方法，MSL方法擬合的RMSE都各期要小于其他3種方法，其中2011年上半年和2012年，這種效果比較明顯，反映出MSL方法保留高質(zhì)量原始數(shù)據(jù)的能力更強(qiáng)。同時(shí)說(shuō)明MSL方法擬合值與原始值之間的平均差異程度更小，代表性更好。SL方法的表現(xiàn)比MSL方法稍差，但是優(yōu)于AG方法和DL方法。除了在2011和2012上半年，AG方法與DL方法差別不是很大。

表1　4種方法擬合2011—2013年EVI數(shù)據(jù)的相關(guān)系數(shù)的總體均值

非對(duì)稱性高斯函數(shù)擬合法(AG)、雙Logistic函數(shù)擬合法(DL)和單Logistic函數(shù)擬合法(SL)

表2　4種方法擬合2011—2013年EVI數(shù)據(jù)的RMSE的總體均值

本文統(tǒng)計(jì)了4種方法擬合相關(guān)系數(shù)差值大于0的像元點(diǎn)個(gè)數(shù)占該生長(zhǎng)季總像元數(shù)(2776)的比例(圖3)。

圖3　4方法擬合的相關(guān)系數(shù)差值大于0的像元的百分比　Fig.3　The percentage of pixels which correlation coefficient difference greater than zero among the four methods

從圖3可以看出，MSL與SL相比，在2011年上半年和2013上半年，相關(guān)系數(shù)差值大于0的像元占該期總像元數(shù)的比例都大于50%；MSL方法與AG方法和DL方法相比，相關(guān)系數(shù)差值大于0的像元占該期總像元數(shù)的比例都不小于50%，其中2011年上半年和2013年達(dá)到70%以上。SL方法與DL方法的相關(guān)系數(shù)差值大于0的像元占該期總像元數(shù)的比例都大于50%；SL與AG方法的比較中，2011上半年，2012下半年和2013上半年，這種比例都超過(guò)60%。綜上所述，相對(duì)于其他3種方法，MSL方法在對(duì)絕大多數(shù)像元的擬合上，都表現(xiàn)出更好的效果。MSL方法和SL方法比另外兩種方法保留高質(zhì)量原始數(shù)據(jù)的能力更強(qiáng)。

經(jīng)分析3年6個(gè)生長(zhǎng)季的EVI實(shí)驗(yàn)數(shù)據(jù)(6×2776個(gè)像元)發(fā)現(xiàn)，噪聲點(diǎn)數(shù)目小于3個(gè)的象元比例為17%，小于等于4個(gè)的為82%，為了進(jìn)一步分析4種方法對(duì)不同程度噪聲的抗干擾能力，本文分別隨機(jī)從實(shí)驗(yàn)數(shù)據(jù)中提取噪聲點(diǎn)數(shù)目小于3個(gè)，3—4個(gè)和大于4個(gè)的像元各500個(gè)。利用4種方法分別擬合這3種EVI數(shù)據(jù)得到RMSE，并對(duì)比分析其均值(表3)。

從表3可以看出，MSL和SL比DL和AG抗噪聲干擾的能力更強(qiáng)。對(duì)包含較多噪聲的像元的擬合，MSL比SL、DL和AG表現(xiàn)出更好的擬合效果，RMSE總體均值明顯低于其它3種方法。MSL擬合方法抗噪聲干擾的能力最強(qiáng)。

表3　4種方法擬合EVI數(shù)據(jù)的RMSE均值

3.4分析與討論

函數(shù)曲線擬合法針對(duì)每一像元進(jìn)行單獨(dú)處理，利用函數(shù)方程對(duì)植被指數(shù)時(shí)間序列曲線進(jìn)行重建，不需設(shè)定閾值或經(jīng)驗(yàn)系數(shù)。而其他方法，如NDVI閾值法、后向移動(dòng)平均法、最大上升速率判斷方法都需要先設(shè)定一個(gè)合理的閾值或經(jīng)驗(yàn)系數(shù)。地形復(fù)雜的山區(qū)往往缺乏足夠的物候觀測(cè)數(shù)據(jù)來(lái)確定閾值或經(jīng)驗(yàn)系數(shù)，故本文重點(diǎn)對(duì)比分析了常用的函數(shù)曲線擬合方法。然而，函數(shù)擬合方法會(huì)受其限定的函數(shù)曲線形態(tài)的影響，與其他方法相比，也表現(xiàn)出一定的局限性。如在人為干預(yù)較大的農(nóng)耕地區(qū)，植被指數(shù)曲線往往呈現(xiàn)多鋒，動(dòng)態(tài)濾波方法則能達(dá)到更好的效果。Logistic函數(shù)曲線的特點(diǎn)是開(kāi)始增長(zhǎng)緩慢，而在以后的某一范圍內(nèi)迅速增長(zhǎng)，達(dá)到某限度后，增長(zhǎng)又緩慢下來(lái)。曲線略呈拉長(zhǎng)的“S”型。這種特點(diǎn)符合絕大多數(shù)自然植被的生長(zhǎng)規(guī)律，在函數(shù)曲線擬合方法中最為常用。AG和DL模型都是以整個(gè)生長(zhǎng)季為擬合區(qū)間的(往往為1a)，其曲線形態(tài)近似為“幾“字形。其優(yōu)點(diǎn)在于，在缺乏有效觀測(cè)數(shù)據(jù)的情況下，也能大致模擬出植被整個(gè)生長(zhǎng)季的生長(zhǎng)特征。在生長(zhǎng)季較短的高緯度地區(qū)，一年內(nèi)植被長(zhǎng)時(shí)間處于休眠期，植被生長(zhǎng)期非常短，植被指數(shù)的振幅很小，加上雪蓋云掩等影響，往往缺乏有效觀測(cè)數(shù)據(jù)來(lái)分段模擬植被生長(zhǎng)期或者衰老期的變化特征。在這種情況下，與MSL和SL相比，AG和DL就有一定優(yōu)勢(shì)。一年中，各個(gè)生長(zhǎng)季EVI波動(dòng)特征有不同特點(diǎn)，很多象元上下半年的EVI曲線呈現(xiàn)非對(duì)稱性。分段擬合能較合理的還原植被生長(zhǎng)的季節(jié)性變化特征，而且有更大的靈活性。在地形氣候復(fù)雜的中緯度山區(qū)，植被上下半年的自然生長(zhǎng)規(guī)律都呈現(xiàn)出較好的“S”型曲線形態(tài)，與AG和DL相比，SL和MSL的擬合效果更佳。與其他3種函數(shù)曲線擬合方法相比，MSL擬合方法有更強(qiáng)的抗噪聲干擾的能力。

本文所采用的幾種函數(shù)曲線擬合方法都屬于最小二乘非線性擬合方法，非線性回歸過(guò)程可以獲得模型參數(shù)的最小二乘無(wú)偏性估計(jì)。所采用的算法分為3種：Trust-Region算法、Gauss-Newton算法、和Levenberg-Marquardt算法。利用這此算法進(jìn)行迭代時(shí),都必須先給出模型中各參數(shù)的初始值。初始值選擇不當(dāng)會(huì)對(duì)擬合結(jié)果造成很大影響。本文選用的是最為常用的Levenberg-Marquardt算法，它是利用梯度求最值的算法，屬于“爬山”法的一種，同時(shí)具有Trust-Region算法和Gauss-Newton算法的優(yōu)點(diǎn)。在參數(shù)初始值的確定上，本文先基于模型參數(shù)的生態(tài)學(xué)意義限定某參數(shù)的取值范圍并利用三點(diǎn)法進(jìn)行估算(如logistic函數(shù)中，c為環(huán)境負(fù)載力或容納量，其值大于0)，然后利用線性化回歸求出模型其他參數(shù)的估計(jì)值并且進(jìn)行顯著性檢驗(yàn)，再以此估計(jì)值為初始值進(jìn)行非線性回歸。這樣可以最少的迭代次數(shù)獲得參數(shù)滿意的擬合精度。還有很多可以對(duì)參數(shù)初始值進(jìn)行估算的方法，如四點(diǎn)法、拐點(diǎn)法、遺傳算法等。有關(guān)非線性函數(shù)曲線擬合參數(shù)初始值的估算目前還沒(méi)有一種公認(rèn)的最優(yōu)算法。

4　結(jié)論

本文基于質(zhì)量控制數(shù)據(jù)，改進(jìn)了時(shí)序EVI數(shù)據(jù)重建效果的定量評(píng)價(jià)方法，能相對(duì)客觀的比較不同擬合曲線對(duì)噪聲點(diǎn)的重建效果及其保持原始高質(zhì)量數(shù)據(jù)的能力?；诟倪M(jìn)的評(píng)價(jià)方法，本文以秦嶺均勻分布的10個(gè)樣區(qū)為研究區(qū)，對(duì)比了非對(duì)稱性高斯函數(shù)擬合法(AG)、雙Logistic函數(shù)擬合法(DL)、單Logistic函數(shù)擬合法(SL)和最值優(yōu)化單Logistic函數(shù)擬合法(MSL)，結(jié)果表明：與AG和DL相比，MSL和SL在噪聲點(diǎn)優(yōu)化和保留原始高質(zhì)量數(shù)據(jù)的能力方面更為突出，同時(shí)擬合值與原始值之間的平均差異程度更小，代表性更好。AG和DL的擬合曲線比較接近，但都對(duì)噪聲數(shù)據(jù)比較敏感，曲線形態(tài)容易受到噪聲點(diǎn)的影響發(fā)生變化。AG相對(duì)于DL，在部分像元的處理上擬合效果更好。MSL比SL、DL和AG表現(xiàn)出更好的擬合效果，在噪聲點(diǎn)較多的情況下這種效果更為突出。4種方法擬合效果的優(yōu)劣依次是MSL>SL>AG>DL。

[1]李正國(guó), 楊鵬, 周清波, 王仰麟, 吳文斌, 張莉, 張小飛. 基于時(shí)序植被指數(shù)的華北地區(qū)作物物候期/種植制度的時(shí)空格局特征. 生態(tài)學(xué)報(bào), 2009, 29(11): 6216- 6226.

[2]何月, 樊高峰, 張小偉, 柳苗, 高大偉. 浙江省植被NDVI動(dòng)態(tài)及其對(duì)氣候的響應(yīng). 生態(tài)學(xué)報(bào), 2012, 32(14): 4352- 4362.

[3]Zhu W Q, Chen G S, Jiang N, Liu J H, Mou M J. Estimating carbon flux phenology with satellite-derived land surface phenology and climate drivers for different biomes: A synthesis of AmeriFlux observations. Plos One, 2013, 8(12): e84990.

[4]Josef Cihlar, Hung Ly, Li Z Q, Chen J, Hartley Pokrant, Huang F T. Multitemporal, multichannel AVHRR data sets for land biosphere studies——artifacts and corrections. Remote Sensing of Environment, 1997, 60(1): 35- 57.

[5]李杭燕, 頡耀文, 馬明國(guó). 時(shí)序NDVI數(shù)據(jù)集重建方法評(píng)價(jià)與實(shí)例研究. 遙感技術(shù)與應(yīng)用, 2009, 24(5): 596- 602.

[6]Daniel Lloyd. A phenological classification of terrestrial vegetation cover using shortwave vegetation index imagery. International Journal of Remote Sensing, 1990, 11(12): 2269- 2279.

[7]Michael A White, Peter E Thornton, Steven W Running. A continental phenology model for monitoring vegetation responses to interannual climatic variability. Global Biogeochemical Cycles, 1997, 11(2): 217- 234.

[8]Bradley C Reed, Jesslyn F Brown, Darrel VanderZee, Thomas R Loveland, James W Merchant, Donald O Ohlen. Measuring Phenological Variability from Satellite Imagery. Journal of Vegetation Science, 1994, 5(5): 703- 714.

[9]J?rg D Kaduk, Martin Heimann. A prognostic phenology scheme for global terrestrial carbon cycle models. Climate Research, 1996, 6(1): 1- 19.

[10]Zhang X Y, Mark A Friedl, Crystal B Schaaf, Alan H Strahler, John C F Hodges, Gao F, Bradley C Reed, Alfredo Huete. Monitoring vegetation phenology using MODIS. Remote Sensing of Environment, 2003, 84(3): 471- 475.

[11]Zhu W Q, Pan Y Z, He H, Wang L L, Mou M J, Liu J H. A changing-weight filter method for reconstructing a high-quality NDVI time series to preserve the integrity of vegetation phenology. IEEE Transactions on Geoscience and Remote Sensing, 2012, 50(4): 1085- 1094.

[12]Moulin S, Kergoat L, Viovy N, Dedieu G. Global-scale assessment of vegetation phenology using NOAA/AVHRR satellite measurements. Journal of Climate, 1997, 10(6): 1154- 1170.

[13]Pieter S A Beck, Clement Atzberger, Kjell Arild H?gda, Bernt Johansen, Andrew K Skidmore. Improved monitoring of vegetation dynamics at very high latitudes: a new method using MODIS NDVI. Remote Sensing of Environment, 2005, 100(3): 321- 334.

[14]吳文斌, 楊鵬, 唐華俊, 周清波, Shibasaki Ryosuke, 張莉, 唐鵬欽. 過(guò)去20年中國(guó)耕地生長(zhǎng)季起始期的時(shí)空變化. 生態(tài)學(xué)報(bào), 2009, 29(4): 1777- 1786.

[15]苗翠翠, 江南, 彭世揆, 呂恒, 李揚(yáng), 張瑜, 王妮, 李軍. 基于NDVI時(shí)序數(shù)據(jù)的水稻種植面積遙感監(jiān)測(cè)分析——以江蘇省為例. 地球信息科學(xué)學(xué)報(bào), 2011, 13(2): 273- 280.

[16]Aziz Habib, Chen X L, Gong J Y, Wang H Y, Zhang L. Analysis of China vegetation dynamics using NOAA-AVHRR data from 1982 to 2001. Geo-spatial Information Science, 2009, 12(2): 146- 153.

[17]Jonsson P, Eklundh L. Seasonality extraction by function fitting to time-series of satellite sensor data. IEEE Transactions on Geoscience and Remote Sensing, 2002, 40: 1824- 1832.

[18]宋春橋, 游松財(cái), 柯靈紅, 劉高煥. 藏北地區(qū)三種時(shí)序NDVI重建方法與應(yīng)用分析. 地球信息科學(xué)學(xué)報(bào), 2011, 13(1): 133- 143.

[19]Jennifer N Hird, Gregory J McDermid. Noise reduction of NDVI time series: an empirical comparison of selected techniques. Remote Sensing of Environment, 2009, 113(1): 248- 258.

[20]李儒, 張霞, 劉波, 張兵. 遙感時(shí)間序列數(shù)據(jù)濾波重建算法發(fā)展綜述. 遙感學(xué)報(bào), 2009, 13(2): 335- 341.

[21]曹云鋒, 王正興, 鄧芳萍. 3種濾波算法對(duì)NDVI高質(zhì)量數(shù)據(jù)保真性研究. 遙感技術(shù)與應(yīng)用, 2010, 25(1): 118- 125.

[22]李紅軍, 鄭力, 雷玉平, 李春強(qiáng), 周戡. 基于EOS/MODIS數(shù)據(jù)的NDVI與EVI比較研究. 地理科學(xué)進(jìn)展, 2007, 26(1): 26- 32.

[23]劉珺, 田慶久, 黃彥, 杜靈通. 黃淮海夏玉米物候期遙感監(jiān)測(cè)研究. 遙感信息, 2013, 28(3): 85- 90.

[24]張明偉. 基于MODIS數(shù)據(jù)的作物物候期監(jiān)測(cè)及作物類型識(shí)別模式研究[D]. 武漢: 華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院, 200

Analysis of four time series EVI data reconstruction methods

LIU Yanan1,2,3, XIAO Fei1,2,*, DU Yun1,2

1InstituteofGeodesyandGeophysics,ChineseAcademyofSciences,Wuhan430077,China2KeyLaboratoryforEnvironmentandDisasterMonitoringandEvaluation,Wuhan430077,China3UniversityofChineseAcademyofSciences,Beijing100049,China

Time-series vegetation index data are contaminated with residual noise and cannot be used in land-cover change detection and crop yield estimation directly. To remove noise effectively, researchers have developed a series of methods for vegetation-index data reconstruction. Function curve-fitting methods are popular in the reconstruction of time-series vegetation index data and have been widely applied in many fields. Different function curve-fitting methods have specific adaptabilities to different geographical environments. In practice, researchers usually have to compare many function curve-fitting methods and select the most suitable one according to the characteristics of regional time-series vegetation index curve fluctuation. Therefore, the means of comparing different function curve fitting methods objectively and quantitatively is very important. Based on ten-year MODIS EVI data of evenly distributed sample areas and its quality control data from the Qinling Mountains, the evaluation method for EVI time-series data reconstruction was discussed and developed in this study. The new evaluation method can compare different function curve fitting methods objectively and quantitatively on two important aspects. One is the function curve-fitting effect under the disturbance of noise points, and another is the ability of retaining original high-quality data. In this study, we used EVI time series data of the sampling area in the Qinling Mountains to analyze the stability of the maximums and minimums of the EVI curves and found that the maximums are more stable in the EVI time series data than the minimums. Then we modified the form of the single logistic model on the basis of the above analysis. Finally, the Maximum optimization Logistic function fitting method (MSL) was proposed to improve the accuracy of the EVI time series reconstruction with large noise in complex mountains. In this study, a new evaluation method was used to compare the Asymmetry Gauss function fitting method (AG), Double Logistic function fitting method (DL), and Single Logistic function fitting method (SL) with the Maximum optimization Logistic function fitting method (MSL). The results show the following: (1) For the function curve fitting effect under the disturbance of the noise points while maintaining original high-quality data, AG showed better results in the treatment of several pixels. (2) Compared to the AG and DL, the fitting effect of SL and MSL is more significant. They not only undisturbed by the noise points but also have a stronger ability to maintain the original high-quality data than AG and DL. (3) Compared to the other methods, MSL was found to be more applicable for EVI time-series data reconstruction in complex mountains with large noise.

MODIS EVI; logistic curve fitting; time-series reconstruction; phenology; Qinling

國(guó)家自然科學(xué)基金資助項(xiàng)目(41271125)；國(guó)家重點(diǎn)基礎(chǔ)研究發(fā)展計(jì)劃項(xiàng)目(2012CB417001)

2015- 01- 07; 網(wǎng)絡(luò)出版日期:2015- 11- 16

Corresponding author.E-mail: xiaof@whigg.ac.cn

10.5846/stxb201501070054

劉亞南，肖飛，杜耘.基于秦嶺樣區(qū)的四種時(shí)序EVI函數(shù)擬合方法對(duì)比研究.生態(tài)學(xué)報(bào),2016,36(15):4672- 4679.

Liu Y N, Xiao F, Du Y.Analysis of four time series EVI data reconstruction methods.Acta Ecologica Sinica,2016,36(15):4672- 4679.