基于DTW距離的時序NDVI數(shù)據(jù)植被信息提取——以秦巴山區(qū)為例

韓曉勇，韓　玲，陳魯皖

(長安大學(xué) 地質(zhì)工程與測繪學(xué)院，陜西 西安 710054)

?

基于DTW距離的時序NDVI數(shù)據(jù)植被信息提取
——以秦巴山區(qū)為例

韓曉勇，韓玲，陳魯皖

(長安大學(xué) 地質(zhì)工程與測繪學(xué)院，陜西 西安 710054)

摘要：秦巴山區(qū)是我國重要的生態(tài)屏障，對該區(qū)的植被信息提取開展研究，可為區(qū)內(nèi)生態(tài)服務(wù)功能及自然資源開發(fā)利用提供基礎(chǔ)數(shù)據(jù)。通過加窗處理改進DTW距離相似性算法，結(jié)合臨近度模糊分類方法對2005—2014年的MODIS NDVI時序數(shù)據(jù)進行植被信息提取。首先利用S-G濾波對MODIS NDVI時序數(shù)據(jù)進行重建；再利用2013年的采樣數(shù)據(jù)構(gòu)建各類植被的標(biāo)準(zhǔn)NDVI時序曲線，逐像元計算與標(biāo)準(zhǔn)NDVI時序曲線的加窗DTW距離，利用臨近度模糊分類實現(xiàn)植被信息提取；最后驗證提取精度。結(jié)果表明，算法具有較高的運行效率，可避免錯誤匹配，以較高的精度(總體精度83.8%，kappa系數(shù)0.77)實現(xiàn)長時間序列的植被信息提取。

關(guān)鍵詞：NDVI；時間序列；加窗DTW；臨近度；模糊分類

植被是生物圈的主要構(gòu)成部分，是陸地生態(tài)系統(tǒng)的主體，其生長和分布受環(huán)境制約，可作為氣候變化的指示因子[1]。衛(wèi)星遙感從不同時間、空間尺度和不同光譜波段進行對地觀測，獲取大量的觀測數(shù)據(jù)，在資源調(diào)查與植被變化監(jiān)測等方面得到廣泛應(yīng)用。植被指數(shù)是植被監(jiān)測的最佳指示因子之一，其中以歸一化植被指數(shù)(NDVI) 的應(yīng)用最為廣泛[2]。MODIS NDVI產(chǎn)品數(shù)據(jù)更以其高時間分辨率，較長時間跨度和容易獲取的優(yōu)勢被廣泛應(yīng)用于植被信息提取及變化監(jiān)測[3-4]。不同植被類型在其生長周期內(nèi)擁有不同NDVI時序特征，例如落葉林表現(xiàn)為振幅較大的單峰曲線，一年兩熟的農(nóng)田則為雙峰曲線，可利用植被自有時序特征進行植被信息的識別[5]。Geerken利用傅里葉濾波進行NDVI時序重建，使用相似性度量方法提取草地植被類型[6]。Wardlow通過MODIS NDVI時序決策樹來提取美國中部農(nóng)作物信息，總體精度達84%[7]。那曉東等將離散的傅里葉變換應(yīng)用于NDVI時序重建，采用傅里葉組分的相似度方法提取三江平原濕地植被信息，總體精度達79.67%，Kappa系數(shù)為0.752 5[8]。Peng利用水稻灌水期NDVI時序特征提取湖南省水稻種植面積[9]。管續(xù)棟對MODIS NDVI時序去噪后，引入動態(tài)時間規(guī)整(DTW)距離相似性方法，提取泰國東北部地區(qū)單季稻、雙季稻面積，總體精度為83.38%[10]。

秦巴山區(qū)是我國重要的生態(tài)屏障，是南水北調(diào)中線工程的主要水源地，對該區(qū)域的植被信息實現(xiàn)高精度提取，可為其生態(tài)服務(wù)功能定位及對自然資源的合理開發(fā)利用提供基礎(chǔ)數(shù)據(jù)。本文以2013年的野外調(diào)查數(shù)據(jù)為基礎(chǔ)，利用時序曲線相似性評價的方法提取研究區(qū)內(nèi)2005—2014年的植被特征，利用2014年的驗證數(shù)據(jù)評價提取精度。

1研究區(qū)概況

選擇陜西省內(nèi)的秦巴山區(qū)為研究區(qū)，主要包括商洛市、漢中市和安康市，由于西安市、寶雞市南部也屬于秦巴山區(qū)，故將兩市也納入研究區(qū)(見圖1)，研究區(qū)總面積98 713 km2。秦嶺是黃河、長江兩大水系的分水嶺和我國南北氣候分界線，北坡為暖溫帶氣候，南坡為北亞熱帶氣候，在陜西境內(nèi)連綿約500 km，海拔多在1 000 m以上，主峰太白山海拔3 767 m。巴山是嘉陵江和漢江的分水嶺，北麓為北亞熱帶氣候，南麓為中亞熱帶，海拔平均1 500~2 500 m。研究區(qū)是陜西省最大最主要的林區(qū)，屬于暖溫帶落葉闊葉林和北亞熱帶常綠落葉闊葉混交林地帶。秦嶺與巴山之間，分布著漢中盆地、西鄉(xiāng)盆地、安康盆地、商丹盆地和洛南盆地，盆地內(nèi)耕地集中，農(nóng)業(yè)生產(chǎn)水平較高，是陜南重要的農(nóng)業(yè)生產(chǎn)區(qū)和人口聚居區(qū)。此外研究區(qū)還包括部分渭河平原和渭北旱塬丘陵溝壑區(qū)。秦巴山區(qū)氣候具有明顯的過渡特征，植被類型豐富，按照國家土地利用現(xiàn)狀分類方法，結(jié)合遙感數(shù)據(jù)的分辨率特征，本文以水田、旱地(含水澆地)、有林地和灌木林4個二級類和1個草地一級類進行特征提取。

2數(shù)據(jù)及預(yù)處理

2.1MODIS NDVI數(shù)據(jù)

從EOSDIS(http://reverb.echo.nasa.gov/)下載研究區(qū)2005—2014年的MODIS植被指數(shù)產(chǎn)品(MOD13A2)，空間分辨率1 km，16 d一期數(shù)據(jù)，10年共230期數(shù)據(jù)。使用MRT工具將每期HDF文件中NDVI和NDVI 質(zhì)量控制文件進行拼接、裁剪，統(tǒng)一轉(zhuǎn)換到WGS84 UTM N49投影帶，并利用質(zhì)量控制文件剔除質(zhì)量差的像元。

2.2野外采樣與驗證數(shù)據(jù)

2013年4月至7月及9月對研究區(qū)內(nèi)植被情況進行野外調(diào)查，共調(diào)查樣點173個，為彌補野外調(diào)查樣點的不足，本文利用2013年、2014年World view的高清影像，結(jié)合野外實測數(shù)據(jù)的判讀特征確定130個精度驗證點，兩類樣點的空間分布見圖1。

圖1　研究區(qū)地理位置及樣點數(shù)據(jù)分布

3研究方法

3.1S-G濾波算法原理

理論上植被的NDVI時序曲線應(yīng)是連續(xù)平滑的，但由于傳感器自身性能、太陽光照角度、觀測視角及云、大氣氣溶膠等觀測條件和隨機干擾因素的影響，導(dǎo)致NDVI時序曲線呈鋸齒狀的不規(guī)則波動，反映季節(jié)變化趨勢不明顯，限制NDVI時序在植被覆蓋變化分析和信息提取等領(lǐng)域的應(yīng)用，為此發(fā)展一系列NDVI時序重建的方法[11]。選用較常用的Savitzky-Golay(S-G)濾波方法對MODIS NDVI產(chǎn)品進行數(shù)據(jù)重建，S-G濾波是一種移動窗口的加權(quán)平均算法，其設(shè)計思想是利用n階多項式來擬合滑動窗口內(nèi)時序數(shù)據(jù)，多項式系數(shù)使用最小二乘法得出，NDVI時序數(shù)據(jù)集的S-G濾波過程可描述為

(1)

式中： Y*為S-G時序重建數(shù)據(jù);Yj+i代表原始時序數(shù)據(jù);Ci為濾波系數(shù);N為歸一化因子;m為滑動窗口寬度。

通過不同參數(shù)設(shè)置對比，選用二階多項式，滑動窗口寬度(m=2)為5，歸一化因子N為35，這既保證擬合效果，又避免過度擬合。由于滑動窗口寬度為5，濾波處理后重建的時序的長度減少為19，即損失最初兩期(1月份)和最后兩期(12月份)的數(shù)據(jù)，這兩個月份對于植被特征標(biāo)志意義不大，故本文只用濾波后每年19期重建數(shù)據(jù)進行相似性計算。

3.2加窗DTW距離

利用2013年調(diào)查數(shù)據(jù)的NDVI均值曲線作為各植被類型標(biāo)準(zhǔn)生長曲線，再逐年計算各像元時序曲線與標(biāo)準(zhǔn)生長曲線的DTW距離。年際間植被生長的水熱條件、光照條件不盡相同，從而導(dǎo)致植被每年的NDVI時序曲線發(fā)生不同程度的扭曲，因此DTW距離更適合該類時序曲線相似性評價。

(2)

(3)

邊界條件：w1=dmn，wK=d11，該路徑的起止元素為距離矩陣的斜對角線的兩端元素；連續(xù)性：若ws=dr,c，ws-1=dr′c′，0≤r-r′≤1 且0≤c-c′≤1，保證路徑中相鄰元素的連續(xù)性；有界性：max(m,n) ≤K≤m+n-1，路徑所經(jīng)過的矩陣元素個數(shù)K存在上限和下限。

3.3臨近度模糊分類

從模糊分類的角度，遙感數(shù)據(jù)的混合像元即不完全屬于或完全不屬于某個類別，使用貼近度來描述模糊集與標(biāo)準(zhǔn)模糊集靠近程度，將混合像元貼近度最高的類別作為該像元的歸屬類別。模糊集定義為待識別的像元與各植被類型參考時序曲線的DTW距離，標(biāo)準(zhǔn)模糊集是各類別參考時序曲線間的DTW距離，式(4)~式(6)為貼近度計算過程。

(4)

(5)

(6)

4結(jié)果及精度評價4.1標(biāo)準(zhǔn)時序曲線

使用迭代處理消除樣點數(shù)據(jù)受混合像元及植被類型變化的影響。首先提取S-G 濾波后的NDVI樣點數(shù)據(jù)，逐期計算各類均值作為初始參考曲線，再計算各樣點與初始參考曲線的DTW距離，將DTW距離大于2倍標(biāo)準(zhǔn)差的樣點進行剔除，將剩余樣本點計算新的參考曲線，重復(fù)迭代過程，迭代計算的控制條件是前后兩次參考曲線的DTW距離變化小于0.000 1，最終確定各類標(biāo)準(zhǔn)時序曲線(見圖2)。圖2中水田和旱地NDVI值的年內(nèi)變化均為雙峰結(jié)構(gòu)，這與一年兩熟的生產(chǎn)方式一致；水田的NDVI上升時間最早，結(jié)束較晚，生長期比旱地區(qū)長。草地、灌木林與有林地均為單峰曲線，65—129日三者的NDVI增長較快，此時為展葉期到快速生長階段；129—273日NDVI處于平穩(wěn)高值區(qū)，植被穩(wěn)定生長；273日以后NDVI快速衰落，植被進入休眠階段。

圖2　各地類標(biāo)準(zhǔn)時序曲線

4.2加窗DTW距離算法效率

加窗DTW距離算法的設(shè)計目的是減小計算量，本文使用IDL語言將傳統(tǒng)DTW算法和加窗

DTW算法進行對比測試，在CPU主頻3.3 G HZ、內(nèi)存4 G的硬件配置下，計算研究區(qū)383行518列每年19期數(shù)據(jù)的平均耗時分別為：傳統(tǒng)DTW算法42.815 73 s，加窗DTW算法25.744 16 s，加窗DTW算法可提高39.87%的計算效率。另外，加窗算法能有效避免錯誤匹配，本文對水田和旱地樣點數(shù)據(jù)進行錯誤匹配測試，加窗算法可避免18.6%的錯誤匹配。

4.3模糊分類

表1為各類參考時序曲線之間的DTW距離，表中的每列(或每行)是對應(yīng)地類的標(biāo)準(zhǔn)模糊集?？梢娝?、旱地與草地、有林地、灌木林的DTW距離較大，說明水田、旱地與另外3類之間較易區(qū)分。灌木林與草地和有林地的DTW距離較小，說明三者之間易發(fā)生錯分，尤其是灌木林與有林地的DTW距離為0.013 8，二者之間錯分的可能性最大。

表1　標(biāo)準(zhǔn)模糊集

按照貼近度最大原則對研究區(qū)2005—2014年的植被類型進行特征提取，結(jié)果見圖3。圖中各地類信息的空間分布合理，水田主要分布在漢中、西鄉(xiāng)和安康盆地，旱地分布在商丹和洛南盆地以及渭河平原區(qū)，有林地主要分布在秦嶺與巴山山區(qū)，灌木林和草地主要分布在耕地與有林地的過渡地帶。表2為各植被類型面積統(tǒng)計，2005—2014年有林地面積緩慢增加了約2 500 km2，旱地面積減少較多，尤其在2011年陜西實施的“避災(zāi)扶貧移民搬遷工程”后旱地面積顯著減少，減少的主要是山區(qū)的坡耕地，其它3類的面積變化較小。

表2　各植被類型面積統(tǒng)計　km2

4.4分類精度

使用130個驗證點對提取結(jié)果進行精度分析，得到分類混淆矩陣見表3。本文算法總體分類精度為83.8%，Kappa系數(shù)0.77，說明利用DTW距離與模糊分類相結(jié)合的方法可以較高精度提取各類植被。水田與旱地的制圖精度較低，二者的漏分誤差較高，圖2中旱地與水田的分布特征也顯示該特點。草地與灌木林的用戶精度較低，說明兩類的錯

表3　分類精度

分誤差較大，圖3中可見灌木林與草地主要分布在有林地與耕地的過渡地帶，這類地物在遙感影像中主要以混合像元形式存在，從而導(dǎo)致它們的用戶精度較低。

利用2006—2014年的《陜西省統(tǒng)計年鑒》中各市(區(qū))農(nóng)作物播種面積對旱地和水田分類結(jié)果進行檢驗，旱地和水田面積的平均相對誤差分別為15.3%和14.5%，最大相對誤差分別為23.8%和 19.7%，分類精度整體較高。

5結(jié)論

本文使用2005—2014年的NDVI數(shù)據(jù)進行植被信息提取，時間跨度較大，基于DTW距離的相似性分析有效解決年度間氣像因子差異引起的曲線偏移，并通過加窗處理有效提高DTW算法效率，減少錯誤匹配，結(jié)合臨近度模糊分類能夠避免直接選取閾值造成的因植被物候期差異造成的植被誤分類，且植被信息提取結(jié)果總體精度較高。

然而在長時序提取過程中仍然存在一些問題有待進一步研究。首先提取結(jié)果受到MODIS NDVI產(chǎn)品數(shù)據(jù)質(zhì)量的制約。例如，研究區(qū)內(nèi)2009年樣點的NDVI數(shù)據(jù)S-G濾波后，DTW距離全部超過2倍標(biāo)準(zhǔn)差，最終2009年的信息提取結(jié)果嚴(yán)重失真，這應(yīng)該是與2009年存在11期大面積異常數(shù)據(jù)有關(guān)，這意味著若年度數(shù)據(jù)存在近一半的噪聲數(shù)據(jù)時要設(shè)計其他方法才能保證特征提取的質(zhì)量；其次是混合像元影響。混合像元的類別只代表像元內(nèi)分布最廣的植被類型，例如在西安東北部土地利用方式多為旱地與果園混雜分布，而果園的NDVI特征與有林地相近，從而產(chǎn)生在某些年度有林地-旱地之間反復(fù)的情況，從而降低信息提取精度，可通過結(jié)合作物物候和中高分辨率遙感影像進行混合像元分解，以提高提取精度實現(xiàn)植被類型進一步細(xì)分。

參考文獻：

[1]任園園,張哲,侯欽磊,等.中國大陸植被覆蓋對降水與溫度變化的響應(yīng)[J].水土保持學(xué)報，2013,33(2)：78-82.

[2]陳燕麗，羅永明，莫偉華，等.MODIS NDVI 與 MODIS EVI 對氣候因子響應(yīng)差異[J].自然資源學(xué)報，2014,29(10):1802-1809.

[3]陳燕麗，龍步菊，潘學(xué)標(biāo),等.基于MODIS NDVI和氣候信息的草原植被變化監(jiān)測[J].應(yīng)用氣象學(xué)報，2010,21(2):229-236.

[4]楊斌，王金生.基于GIS的丘陵區(qū)耕地景觀格局時空演變特征分析[J].測繪工程，2014,23(9):1-4.

[5]李文葉,姜魯光,李鵬.2001—2010年鄱陽湖圩區(qū)水稻多熟種植時空格局變化[J].資源科學(xué)，2014，36(4)：809-816.

[6]GEERKEN R,ZAITCHIK B,EVANS J P.Classifying rangeland vegetation type and coverage from NDVI time series using Fourier Filtered Cycle Similarity[J].International Journal of Remote Sensing,2005,26(24):5535-5554.

[7]WARDLOW B D,EGBERT S L,KASTENS J H.Analysis of time-series MODIS 250m vegetation index data for crop classification in the US Central Great Plains[J].Remote Sensing of Environment,2007,108(3):290-310.

[8]那曉東,張樹清,李曉峰,等.MODIS NDVI時間序列在三江平原濕地植被信息提取中的應(yīng)用[J].濕地科學(xué),2007,5(3):227-236.

[9]PENG D，HUETE A R，HUANG J,et al.Detection and estimation of mixed paddy rice cropping patterns with MODIS data[J].International Journal of Applied Earth Observation and Geoinformation，2011，13(1)：13-23.

[10] 管續(xù)棟,黃翀,劉高煥,等.基于DTW距離的時序相似性方法提取水稻遙感信息：以泰國為例[J].資源科學(xué),2014,36(2):267-272.

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

[責(zé)任編輯：張德福]

Extraction of vegetation information using adding windows DTW distance with NDVI time series data

HAN Xiaoyong,HAN Ling，CHEN Luwan

(School of Geology Engineering and Geomatics，Chang’an University，Xi’an 710054，China)

Abstract：Qinling-Bashan Mountain area is an ecological barrier.The vegetation information extraction research can provide base data for ecological service function and natural resources exploitation.DTW distance similarity measure is improved by adding windows，and vegetation information is extracted by proximity-fuzzy classification from MODIS NDVI in the years of 2005 to 2014.S-G filter method is first applied to weakening the noise and reconstructing the NDVI time series.Each type of vegetation standard NDVI time series curves is generated from the 2013 sample data.For each to-be-classified pixel，a quantitative similarity between its time-curve and standard time series curves is calculated using adding windows DTW，then vegetation information is extracted in research area.Finally，validate data is used to verify the extraction accuracy.The result shows that the algorithm has high efficient and can avoid mismatch，and can realize long time serial vegetation information extraction with high accuracy (overall accuracy:83.8%，Kappa coefficient:0.77).

Key words：NDVI；time series；adding windows DTW；proximity；fuzzy classification

中圖分類號：TP79

文獻標(biāo)識碼：A

文章編號：1006-7949(2016)03-0011-06

作者簡介：韓曉勇(1981-)，男，博士研究生.

基金項目：國家自然科學(xué)基金資助項目(41201099)

收稿日期：2015-06-11