基于MODIS時序的陜西省植被物候時空變化特征分析

韓紅珠， 白建軍， 張波， 馬高

(1.陜西師范大學地理科學與旅游學院，西安 710062；2.陜西師范大學地理學國家級實驗教學示范中心，西安 710062；3.中國空間技術研究院西安分院，西安 710000)

0 引言

植被物候是指植物受氣候和其他環境因子的影響而出現的以年為周期的自然現象，包括植物的發芽、展葉、開花、葉變色和落葉等，是植物長期適應季節性變化的環境而形成的生長發育節律[1]。因此，植被物候具有反映氣候變化和環境變化的作用，是指示氣候與自然環境變化的重要指標[2-3]。同時，植被物候還是陸面過程模型、動態全球植被模型的重要參數[1,4]，在全球變化研究中具有重要作用。

傳統植被物候觀測主要通過氣象站點的觀測和記錄實現，這種觀測方式具有覆蓋面積小、成本耗費高等缺點，不足以應對大面積的植被物候觀測要求。而衛星遙感數據具有觀測范圍廣、空間連續性好、時間序列長、信息量豐富、數據種類多、經濟效益高等特點，在大尺度、長時間范圍內的土地利用、植被覆蓋和農業資源調查等工作中被廣泛使用，應用衛星數據進行遙感監測也成為監測區域或全球尺度植被動態變化的強大手段[5]。

衛星遙感數據植被指數產品在大尺度上監測和描述植被物候變化趨勢的能力，能夠很好地描述整個生態系統的物候變化[6]。目前，不少專家和學者已經注意到衛星遙感數據植被指數產品應用到植被物候方面研究的潛力，且已有部分研究工作者開展了這方面的研究，在中國溫帶地區[6]、中國北方溫帶地區[7]、青藏高原[8]、中國農牧交錯地帶[9]、東北[10]、華北[11]和浙江[12]等不同尺度的研究區域取得了許多研究成果。但這些成果因氣候類型不同、地形差異大，具有區域局限性而不適于描述其研究區以外的植被物候，尤其對于氣候類型差異顯著的區域有較大的局限性。

陜西省地處中國大陸腹地，緯度跨度大，地形和氣候類型豐富，開展陜西省的植被物候研究、分析其時空變化特征，有助于理解在多種地形和氣候類型影響下的植被物候響應機制，豐富我國植被物候研究的區域范圍，同時為其他區域的植被物候研究提供對比資料。

1 研究區概況與數據源

1.1 研究區概況

陜西省位于中國內陸，地理范圍在E105°29′～111°15′，N31°42′～39°35′之間，地處東部濕潤地區和西部干旱區的交界地帶，因此省內氣候類型多樣，可分為暖/中溫帶半干旱地區(約N37°以北)、暖溫帶半濕潤地區(約N33°～37°)、北亞熱帶濕潤地區(約N33°以南)3大氣候區[13-14]；省內地形地貌類型豐富，自北向南依次劃分為毛烏素沙地陜北部分、黃土高原、關中平原、秦嶺山地、漢江盆地和大巴山山地。

1.2 數據來源和預處理

研究選用MODIS產品之一的MOD13Q1 v006植被指數產品數據集。該數據集包括基于限制視角的最大值合成方法(constrained view angle- maximum value composite，CV-MVC)合成的16 d時間分辨率、250 m空間分辨率的歸一化植被指數(normalized difference vegetation index，NDVI)和增強型植被指數等數據。與MOD13Q1 v005數據相比，v006數據在植被指數的運算中使用8 d預合成地表反射率數據而非天數據，改進了NDVI產品在生產過程中的誤差，提高了數據精度。產品已經過幾何糾正和大氣校正。該時間序列數據集從美國國家航空航天局(National Aeronautics and Space Administration，NASA)陸地過程分布式主動存檔中心(Land Processes Distributed Active Archive Center, LP DAAC)網站下載(https: //lpdaac.usgs.gov/lpdaac/products/modis products table)，時間覆蓋2001年1月—2016年12月，16 a間共368期影像，用于生產陜西省遙感植被物候數據。另外，利用來源于清華大學基于30 m空間分辨率TM影像的土地覆蓋分類圖(http: //data.ess.tsinghua.edu.cn/)剔除陜西省的非植被覆蓋區。

利用MRT工具從MOD13Q1 v006的368期影像中提取NDVI和數據像元可信度這2個數據層，并進行拼接、格式轉換和重投影。數據像元可信度是對NDVI數據每個像元數值可信度的說明，與NDVI數據的時間和空間分辨率一一對應，取值范圍為-1～3。根據參考文獻[15]，利用數據像元可信度對NDVI數據進行預處理，進一步提高NDVI像元的可信度；然后結合30 m空間分辨率的土地覆蓋分類圖在ArcMap中利用重采樣和掩模提取陜西省的植被覆蓋區，同時剔除NDVI<0.05的像元，盡量保證每個NDVI像元都是地表植被覆蓋區。

2 物候提取

2.1 時間序列重構

目前，用于植被指數數據時間序列重構的方法很多，包括最佳指數斜率提取法、均值迭代濾波法、Hants算法、Savitzky-Golay(S-G)濾波法、非對稱高斯函數擬合法、時間窗的線性內插法、雙函數擬合方法和小波變換等[16-18]。這些方法中既有基于設置閾值的噪聲去除法，也有通過濾波達到平滑效果和以非線性擬合為基礎的算法。每種方法均有其優缺點及適用范圍，例如，最佳指數斜率提取法能有效處理NDVI曲線中突降的點值，但對異常高值缺乏相應對策，且在處理過程中，經驗性要求高[19]；基于傅里葉變換的濾波方法對較為對稱的NDVI時序曲線處理效果較好，但對于存在波動的NDVI時序曲線進行平滑容易引入新的人為噪聲；非對稱高斯函數擬合法對時序曲線進行分段擬合，具有較高的擬合精度，但在實際運算過程中，需要人為確定分段擬合的起始值，明顯受制于工作經驗，且運算復雜耗時，效率較低；而S-G濾波法采用NDVI數據的上包絡線來擬合NDVI時序數列的變化趨勢，保留了NDVI曲線的特征，通過一個迭代的過程可消除偏離正常生長趨勢線的噪聲，而使平滑達到最好的效果[19]。因此，本文選擇S-G濾波法進行NDVI時序重構。

S-G濾波法由Savitzky和Golay于1960年提出，是一種滑動窗口的加權擬合算法。其加權系數通過在滑動窗口內對給定高階多項式進行最小二乘擬合得出。因此，在利用S-G濾波對NDVI時序曲線平滑去噪的過程中，2個參數的設置至關重要，一個是滑動窗口數值，另一個是多項式階數。一般來說，滑動窗口數值越大，參與擬合的值越多，效果越好，但細節容易被過濾掉，反之則平滑效果不明顯；多項式的階數越高，越易出現過擬合，階數越低，平滑效果越好，但會引進誤差[13,19]。已有研究結果表明，滑動窗口取值一般在3～7之間，多項式階數一般在2～4之間比較合適[13,19]。因此，本文經過運算結果比較研究后，設置滑動窗口數值為5，多項式階數為2，進行NDVI時序曲線的S-G濾波，實現NDVI時序曲線重構。

2.2 物候期提取

本研究中提取的植被物候信息主要針對遙感監測，且只考慮一年中只有一個生長期的情況(對于一年中存在2個生長期的地區，提取該地區植被生長期中NDVI值最大的一個生長期)，提取植被生長過程中共同的3個物候期特征指數，分別是植被生長季開始期(the start of season，SOS)、生長季結束期(the end of season，EOS)和生長季長度(the length of season，LOS)。其中，SOS指植被開始增長的日期或者光合作用開始加強的日期；EOS指植被光合作用和綠色葉面積開始迅速下降的日期；LOS指SOS到EOS以天數計的時間跨度[5]。

目前，基于植被指數時序曲線提取遙感物候的方法不少，較常用的有閾值法、最大斜率法、最大比率法和累積頻率法等[6,20-22]。考慮到計算效率以及可操作性，采用最大比率法提取陜西省植被物候期。文中日期均采用儒略日計算法，即將每年的1月1日記為第1天。

3 結果與分析

3.1 植被物候多年均值的空間分布特征

圖1展示了陜西省遙感植被物候2001—2016年間平均值的空間分布。從圖中可以看出，各物候期隨緯度變化形成南北差異，隨地形起伏變化形成區域差異，總體植被物候的空間分布特征與陜西省不同的地形地貌的空間分布具有較好的一致性。因此在本研究分析中，按照陜西省地形地貌的空間分布特征[23-24]將圖1自北向南依次劃分為毛烏素沙地陜北部分、黃土高原、關中平原、秦嶺山地、漢江盆地和大巴山山地進行討論。

(a) SOS (b) EOS (c) LOS

陜西省遙感植被物候空間分布特征顯示出地域分異規律，展現植被生長受氣候和地形影響的一般規律，表現出在不同氣候和地形條件下，植被物候時間的差異。例如，暖/中溫帶半干旱地區(約N37°以北)溫度和降水條件相對較差，SOS和EOS均較晚，LOS較短；暖溫帶半濕潤地區(約N33°～37°)、北亞熱帶濕潤地區(約N33°以南)溫度和降水條件較好，SOS較早，EOS較晚，LOS較長。不同氣候區物候期的差異表明氣候對物候產生影響。這與前人的研究成果相似[25]。而在同一氣候區內，物候期仍表現出明顯的差異，主要表現在地形起伏較大的地區。如暖溫帶半濕潤地區(約N33°～37°)中黃土高原南部、關中平原和秦嶺山地3個地區物候期明顯不同。

陜西省遙感植被物候的空間分布總體上隨緯度和地形發生變化。相對于緯度較低的地區，緯度較高的地區SOS和EOS均較晚，LOS較短。在同一緯度，海拔高的植被SOS晚，EOS早，LOS較短；海拔低的植被SOS早，EOS晚，LOS較長。陜西省遙感植被物候的空間分布特征表現出與其不同地形地貌空間分布較好的一致性。

對于不同的物候期，陜西省遙感植被物候的空間分布差異不同。SOS最明顯，與陜西省不同地形地貌的空間分布一致性最佳。主要原因可能是在植被開始生長時，受到環境因素(溫度、降水等)的作用，而環境因素又與緯度、地形等條件密切相關，因此SOS表現的一致性最佳。具體來說，從高緯度至低緯度SOS日期在逐漸提前，在陜西黃土高原內部SOS時間呈現出明顯的南北差異，黃土高原北部的SOS明顯遲于南部，再往南到達關中平原地區，SOS更加提前；在秦嶺山地，由于海拔上升，SOS總體較遲且山頂遲于山腳；在緯度進一步降低和海拔有所下降的漢江盆地，SOS又較早；在大巴山山地因海拔升高，SOS又較遲。

EOS的空間分布差異不如SOS明顯，與陜西省不同地形地貌的空間分布一致性相對較差。主要原因可能是植被在后期發育以及成熟階段，植被生長不僅受到環境因素影響，而且植被自身生理條件、土壤狀況和病蟲害等條件同樣對植被生長起作用。盡管如此，EOS仍表現出一定的規律性。從高緯度到低緯度，總體上EOS日期在提前，因為相對于高緯度地區，低緯度地區的熱量更充足，植被完成生長所需的積溫時間較短；EOS受地形的影響明顯，在秦嶺山地和大巴山山地可以看出，EOS總體保持著隨海拔升高而時間提前的趨勢。

LOS的空間分布是SOS和EOS兩者的綜合，表現出一定的區域差異特征。且LOS從高緯度到低緯度、從高海拔到低海拔，總體上保持著逐漸延長的趨勢。從陜西省植被物候2001—2016年各年的時間看，SOS最早出現在第80天之前，最晚出現在一年中第162天之后，平均時間在第120天；EOS最早出現在第220天之前，最晚出現在第300天之后，平均時間在第280天；LOS最短為112 d以內，最長在192 d以上，平均時間為160 d。區域差異明顯的物候期平均時間見表1。在氣候和地形的共同影響下，漢江盆地的SOS是陜西省幾種地形地貌類型中最早的，比最晚的黃土高原早45.9 d。以約N35°為界，EOS北部的地貌類型與南部的地貌類型差異較大，在界限內的地貌類型EOS差異均在7 d以內。LOS的區域差異情況與SOS相似，漢江盆地LOS最長，其次是關中平原，最后是黃土高原。

表1 2001—2016年間各區域物候期平均時間

3.2 植被物候年際變化特征分析

圖2顯示了2001—2016年間陜西省遙感物候年際變化特征。

(a) SOS (b) EOS (c) LOS

由圖2可知，2001—2016年間各物候期呈現提前或推遲的波動變化。總體上，SOS呈波動提前的趨勢，變化率約為-0.79 d/a(R2=0.40，P<0.01)，處在1981—1991年期間全球N45°～70°地區的-1～-0.7 d/a[26]的范圍內，與中國溫帶地區1982—1999年間的-0.79 d/a[6]持平，高于全球1981—2003年間的-0.38 d/a[27]，與北半球2000—2008年間的-0.2 d/a[3]的幅度差距較大；EOS呈波動推遲，變化率約為0.50 d/a(R2=0.25，P<0.05)，高于1981—1991年間全球N45°～70°的0.4 d/a[26]、中國溫帶地區1982—1999年間的0.37 d/a[6]及全球1981—2003年間的0.45 d/a[27]，低于北半球2000—2008年間的0.66 d/a[3]；LOS表現為波動延長，變化率約為1.29 d/a(R2=0.37，P<0.05)，高于中國溫帶地區1982—1999年間的1.16 d/a[6]。

以上分析可以看出本文提取的陜西省遙感植被物候時間變化率與前人的研究結果存在差異，主要原因可能是時間尺度和區域范圍不同導致。在3個物候期的年際變化中，SOS的變化趨勢最顯著，其次為LOS，EOS的變化趨勢顯著性最低。對SOS，EOS和LOS分別做相關分析，發現SOS與LOS的相關性比EOS對LOS的相關性更顯著。SOS與LOS相關分析結果為N=16，R2=0.90，P=0; EOS與LOS相關分析結果為N=16，R2=0.85，P=0。說明LOS主要與SOS相關，SOS越早，LOS越長。

3.3 植被物候年際變化的空間分布特征

圖3展示了2001—2016年間陜西省遙感物候變化趨勢的空間分布特征。從圖3可以看出，陜西省植被物候在各物候期均有提前和推遲的趨勢變化，且SOS和LOS的變化幅度較EOS強；在不同的物候期，不同的趨勢變化和變化幅度存在明顯的空間分布差異。根據年變化率k和P值的大小，將陜西省2001—2016年間遙感物候年際變化趨勢分為3種情況：①當k<0且P<0.2時，顯著提前；②當k>0且P<0.2時，顯著推遲；③當P>0.2時,無顯著變化。對于SOS，陜西省大多數地區表現為提前趨勢，且從北到南均有分布，提前幅度在一定區域范圍內表現為從北向南逐漸減弱，如黃土高原北部與南部；具有推遲趨勢的區域主要分布在毛烏素沙地陜北部分以及約N35.5°以南地區。但在SOS的趨勢變化中只有約26.55%(P<0.2)的像元表現出顯著性，其中顯著提前的區域占總面積的20.96%，提前幅度范圍約為1～6 d，且從北向南幅度逐漸減小；顯著推遲的區域占總面積的5.58%，推遲幅度約為2～9 d且從北向南幅度逐漸增加。對于EOS，總共約16.46%(P<0.2)的像元表現出顯著性。大部分區域表現為推遲趨勢，占陜西省總面積的10.1%，推遲幅度約為1～7 d且從北向南幅度逐漸增加。提前趨勢主要出現在陜西省中部和東南部，約占總面積的6.39%，提前幅度約為1～7 d，且從北向南幅度逐漸增加。從LOS的變化來看，總共約23.89%(P<0.2)的像元表現出顯著性。絕大部分地區表現為延長趨勢，約占總面積的18.07%，延長幅度約為2～9 d。縮短趨勢主要集中在毛烏素沙地陜北部分和陜西東南角，約占5.82%，縮短幅度約為3～12 d。

(a) SOS年變化率 (b) EOS年變化率 (c) LOS年變化率

(d) SOS變化顯著水平 (e) EOS變化顯著水平 (f) LOS變化顯著水平

3.4 不確定性分析

遙感植被物候的監測基礎是植被指數時序變化，因此繼承了來源于植被指數本身由于“同物異譜”、“同譜異物”的光譜差異造成的反映地物真實情況的誤差。植被信號的強弱也會影響物候的遙感識別，植被信號越強越有利于物候的遙感識別[28]。因此在植被信號相對較弱的陜西省最北部的毛烏素沙地，植被物候期識別可能存在較大誤差。在地形相對復雜的秦嶺山地和大巴山山地，也因地形(坡度、坡向)、山地氣候(云、雨、霧)對植被反射率的影響而影響遙感植被物候的探測精度。此外，本文選取的遙感數據MODIS NDVI的時間分辨率為16 d，在監測植被物候的實時發生日期上可能存在偏差。

4 結論

以具備多種地形地貌和氣候類型的陜西省為研究區，在采用S-G濾波法對MODIS NDVI數據進行時間序列重構的基礎上，提取了陜西省2001—2016年間的植被物候參數并對結果進行時空特征分析，得出以下結論：

1)陜西省的遙感植被物候空間分布特征與其不同地形地貌的空間分布具有較好的一致性，顯示出地域分異規律，其中，SOS的區域差異分布與地形地貌的空間分布一致性最佳，區域差異最明顯。

2)陜西省2001—2016年間植被物候的SOS平均時間在每年的第120天，EOS平均時間在第280天，LOS平均時間為160 d。

3)2001—2016年間陜西省植被SOS呈提前趨勢，變化率約為-0.79 d/a(R2=0.40，P<0.01)；EOS呈推遲趨勢，變化率約為0.50 d/a(R2=0.25，P<0.05)；LOS呈延長趨勢，變化率約為1.29 d/a(R2=0.37，P<0.05)。

4)陜西省植被LOS主要與SOS相關，SOS越早，LOS越長。

5)在不同的物候時期，陜西省植被的物候變化趨勢空間分布差異較大。SOS顯著提前的區域主要分布在陜西黃土高原地區，顯著推遲的區域主要分布在陜南和毛烏素沙地陜北部分；EOS顯著提前的區域主要分布在關中平原地區，顯著推遲的區域主要分布在黃土高原中部和陜西西部；LOS顯著縮短的區域主要集中在毛烏素沙地陜北部分，顯著延長的區域集中在陜西黃土高原和陜西西部。