宗加權, 白淑英, 馮朝陽, 蒲 陽

(1.南京信息工程大學, 南京 210044 ; 2.中國環(huán)境科學研究院,國家環(huán)境保護區(qū)域生態(tài)過程與功能評估重點實驗室, 北京 100012)

植被是陸地生態(tài)系統(tǒng)中的重要組成部分，在陸地表層能量交換、水分循環(huán)和生物地球化學循環(huán)過程中扮演著重要角色[1]，植被覆蓋度變化是生態(tài)環(huán)境變化的直接結果，很大程度上表示了生態(tài)環(huán)境總體狀況[2]。相對于傳統(tǒng)實地調查植被覆蓋特征費時費力，通過遙感手段獲取數據，優(yōu)勢有空間范圍覆蓋面積大、時間序列長、方便獲取等。歸一化差值植被指數(NDVI)是目前常用的研究地表植被覆蓋特征的重要指標之一[3]，長時間序列監(jiān)測植被覆蓋特征變化對更好地模擬陸地生態(tài)系統(tǒng)的動態(tài)變化特征具有重要意義[4]。AVHRR NDVI是目前時間序列最長的全球連續(xù)數據集，具有很高的科研價值和實際意義[5-6];MODIS NDVI是從2000年開始的數據集，提高了葉綠素敏感度和空間分辨率、排除了大氣水汽的干擾、改善了合成方法，但時間序列長度有限。探尋兩者間的關系，通過數據擬合構建長時間序列的NDVI數據集，是實現(xiàn)對植被長期觀測的關鍵[7-8]。一些學者基于不同NDVI數據集重構NDVI數據，并結合氣象數據等分析了植被狀況的時空變化特征，如:李凈等[9]基于重構后的數據產品和植被類型數據集，分析了不同時期植被NDVI變化特征，發(fā)現(xiàn)西北地區(qū)NDVI整體呈上升趨勢;李飛等[10]根據MODIS NDVI與AVHRR數據構建了較長時間序列NDVI數據集，采用相關性分析法分析了全國農業(yè)、森林、草地與稀疏植被區(qū)植被活動的變化特征，結果指出近30 a來中國陸地植被活動整體增強，這與方精云等[11]利用AVHRR NDVI研究中國地區(qū)植被活動的結果相同。趙玉萍等[12]利用1982—2003年GIMMS NDVI數據集和逐月氣象資料，用相關分析法研究了藏北高原草地生態(tài)系統(tǒng)NDVI與氣候因子的相關性及滯后性，得到NDVI與月平均氣溫等呈高度正相關(p<0.001)。劉可等[13]利用1982—2012年GIMMS NDVI3g和中國陸地生態(tài)系統(tǒng)類型數據，采用一元線性經驗模態(tài)分解和相關分析方法，研究了近30 a中國各生態(tài)系統(tǒng)NDVI的時空變化特征，結果表明近30 a中國植被活動顯著上升，農田、森林、草地和水體與濕地生態(tài)系統(tǒng)的NDVI總體為非穩(wěn)定上升趨勢。中國陸地生態(tài)系統(tǒng)中的森林、草原、荒漠、水體與濕地和其他生態(tài)系統(tǒng)類型是對土地利用/覆蓋的分類判別匯總，對生態(tài)系統(tǒng)狀況有指示意義。水體與濕地在生態(tài)環(huán)境評價中起到重要作用，其水體與濕地狀況好壞和面積變化影響著生態(tài)環(huán)境狀況，水域面積的減少會導致生態(tài)環(huán)境惡化[14]。水體和濕地因為水源豐富、土壤肥沃，很容易被開墾為良田和進行水產的養(yǎng)殖，水體與濕地生態(tài)系統(tǒng)因此被嚴重破壞[15]，牛振國等[16]研究表明，近30 a，中國濕地類型由濕地轉化為非濕地主要是氣候變化和農業(yè)活動引起的。

以上研究存在一些問題，首先采用的AVHRR數據為一代數據產品，數據時間跨度不長，數據的空間分辨率不高(8 km×8 km);其次大多數研究多從行政區(qū)劃和地理單元探究植被活動特征，未考慮到不同生態(tài)系統(tǒng)間的鑲嵌性，因此難以詳細地探究不同生態(tài)系統(tǒng)間NDVI的差異。對此，本研究基于重構的1981—2018年長時間序列植被指數數據集(GIMMS NDVI3g和MODSI NDVI)，采用分區(qū)統(tǒng)計、NDVI趨勢分析法，分析中國大陸范圍(包含香港、澳門，未包含中國臺灣省)非耕地、非聚落和濕地生態(tài)系統(tǒng)NDVI的空間分布和時間動態(tài)變化特征，試圖為中國生態(tài)環(huán)境保護和可持續(xù)發(fā)展提供決策依據。

1 數據來源及預處理

1.1 數據來源

本研究利用2000年2月—2018年12月MODIS13 A3數據集植被指數產品的月合成數據(https:∥ladsweb.modaps.eosdis.nasa.gov/search/)，空間分辨率為1 km。還有美國國家航天航空局NASA Goddard Space Flight Center的全球監(jiān)測與模型研究組(GIMMS)采用最大值合成法(MVC)合成的15 d、空間分辨率為8 km的最新GIMMS NDVI 3g.v1數據，時間范圍為1981年7月—2015年12月(https:∥ecocast.arc.nasa.gov/data/pub/gimms/3g.v1)。

中國陸地生態(tài)系統(tǒng)數據來源于中國科學院地理科學與資源研究所、中國科學院資源環(huán)境科學數據中心的“資源環(huán)境數據云平臺”。包括1980年、1990年、1995年、2000年、2005年、2010年和2015年，共7時期數據。該數據是在遙感解譯獲取的1∶10萬比例尺土地利用/土地覆蓋數據的基礎上，通過對各生態(tài)系統(tǒng)類型進行辨識和研究，經過分類處理形成多期中國陸地生態(tài)系統(tǒng)類型空間分布數據集，矢量數據柵格化生成的1 km柵格數據。生態(tài)系統(tǒng)類型依據土地利用/土地覆蓋遙感分類系統(tǒng)，歸分為7大生態(tài)系統(tǒng)類型:(1) 農田生態(tài)系統(tǒng)，主要包括水田和旱地;(2) 森林生態(tài)系統(tǒng)，主要包括密林地(有林地)、灌叢、疏林地、其他林地;(3) 草地生態(tài)系統(tǒng)，主要包括高覆蓋度草地、中覆蓋度草地、低覆蓋度草地;(4) 水體與濕地生態(tài)系統(tǒng)，主要包括沼澤地、河渠、湖泊、水庫、冰川與永久積雪、灘地;(5) 聚落生態(tài)系統(tǒng)，主要包括城鎮(zhèn)、農村居民地、工礦;(6) 荒漠生態(tài)系統(tǒng)，主要包括沙地、戈壁、鹽堿地;高寒荒漠;(7) 其他生態(tài)系統(tǒng)，主要包括裸土地和裸巖礫石地。

1.2 數據預處理

采用MRT軟件對MODIS NDVI數據進行拼接，投影轉換，數據格式轉換。采用Matlab軟件對GIMMS NDVI數據進行逐月最大值合成(MVC)，進一步減少云的影響，并降低月內物候循環(huán)的影響[17]，還原真值，格式轉換。由于GIMMS NDVI數據的空間分辨率約為8 km，將GIMMS NDVI由8 km重采樣為1 km的數據，不損失影像信息，以匹配MODIS NDVI。本文采用的地理坐標系是WGS-84,Albers投影坐標系，隨后統(tǒng)一采用IDL程序對兩套數據進行逐月裁切、掩膜、還原真值等操作。

2 研究方法

MODIS NDVI數據被稱為GIMMS NDVI數據的延續(xù)[18]，兩套數據集有一段時間的重疊部分，即2000年2月份-2015年12月份。已有學者[17,19]的研究表明可以利用兩套數據的重疊時間部分，以此建立逐月線性轉換方程，但未能從文章中找出適合本研究區(qū)域的轉換方程，本研究中建立的最佳轉換方程為y=0.970034x+0.016953，擬合方程的判定系數R2為0.88,p<0.01，擬合效果較好。利用轉換方程將1981年7月—2001年12月的GIMMS NDVI進行修正擬合，為了保留MODIS數據1 km的空間分辨率，將GIMMS NDVI數據也重采樣為1 km，并連接MODIS NDVI 2002年1月至2018年12月數據，最終獲得中國大陸地區(qū)的1981年7月至2018年12月長時間序列空間分辨率為1 km的NDVI數據集。

2.1 NDVI時空變化趨勢分析法

2.1.1 NDVI時間變化趨勢分析法 趨勢分析是對一組隨時間變化的變量進行回歸分析，得到其變化趨勢的方法。對多年植被指數數據逐像元NDVI值進行一元線性回歸，得到研究區(qū)內植被隨時間的變化趨勢，即NDVI的年際變化。計算公式為[20-21]:

(1)

式中:k為研究時間段內的年數,取1-n；YNDVIk是第k年某季相的NDVI平均值,如夏季NDVI平均值采用當年6月、7月、8月份NDVI值求得平均；slope>0表示植被覆蓋呈增加趨勢;slope小于0則相反;slope=0表示植被覆蓋沒有明顯的變化。

2.1.2 NDVI空間變化分析法 將大陸地區(qū)分為六大區(qū)域(東北、華北、華東、中南、西北、西南)，分區(qū)統(tǒng)計分析NDVI的總體變化特征，能整體上體現(xiàn)NDVI均值和趨勢變化情況。

2.2 主要生態(tài)系統(tǒng)類型的NDVI時空變化分析法

采用分區(qū)統(tǒng)計的方法來統(tǒng)計不同生態(tài)系統(tǒng)NDVI的變化情況，及其空間分異特征。由于生態(tài)系統(tǒng)數據為多年一幅的柵格數據，為充分利用生態(tài)系統(tǒng)數據，將連續(xù)時間序列的NDVI與最鄰近年的生態(tài)系統(tǒng)數據進行對應統(tǒng)計，表1是分區(qū)統(tǒng)計對應表。

表1 分區(qū)統(tǒng)計年份對應

3 結果與分析

3.1 中國植被NDVI時空變化分析

3.1.1 中國植被NDVI時間變化分析 圖1是整個中國大陸范圍內NDVI大于0值，2000年2月至2015年12月GIMMS NDVI平均值修正前后與MODIS NDVI平均值變化曲線，由圖1可知，原始GIMMS NDVI的周期性變化趨勢與MODIS NDVI相似，每年冬季原始GIMMS NDVI都小于MODIS NDVI，在植被生長季，兩者相差不大，修正后GIMMS NDVI與MODIS NDVI值相接近，具有一致性。值得注意的是MODIS NDVI數據在2000年初和2001年冬季，其值明顯低于修正前后GIMMS NDVI數據，可見MODIS數據在剛獲取的階段，數據質量可能不穩(wěn)定，所以對于中國大陸地區(qū)，在使用MODIS NDVI數據時，可以從2001年7月份開始選用。

圖1 2000年2月至2015年12月MODIS NDVI、原始與修正后GIMMS NDVI時間序列曲線對比

由表2看出，修正后GIMMS NDVI在3種數據中的平均值和最小值都為最大，平均值最大為0.373，最小值最大為0.273;除了2000年和2001年，修正后GIMMS NDVI與MODIS NDVI的年平均值環(huán)比增長率和最大值有相似的變化規(guī)律，且有相同的NDVI線性趨勢(0.000 8/a)，最大值最小僅相差0.002。這說明，除了2000年、2001年，修正后GIMMS NDVI在年尺度上與MODIS NDVI相似的變化規(guī)律。

表2 不同數據年平均值、最大值、最小值和年平均值環(huán)比增長率

從表3中看出，MODIS NDVI和修正后GIMMS NDVI在2月至5月NDVI平均值環(huán)比增長率逐漸增加，而從5月至6月NDVI值月平均環(huán)比增長率有所下降，這是因為冬小麥成熟收割，新播種的作物未生長成熟造成的NDVI值下降現(xiàn)象，原始GIMMS NDVI未能捕捉到這個現(xiàn)象;修正后GIMMS NDVI平均值和最大值與MODIS NDVI有相似的變化規(guī)律，GIMMS NDVI月均值環(huán)比增長率最大為0.266%,MODIS NDVI為0.273%;修正后GIMMS NDVI最小值大于MODIS NDVI,2月份最大差值為0.096，這可能是由于2000年2月份MODIS NDVI數據剛剛獲取，數據質量不穩(wěn)定導致NDVI值偏低，10月份相差值最小為0.008。

表3 不同數據月平均值、最大值、最小值和月平均值環(huán)比增長率

3.1.2 中國植被NDVI空間變化分析 從圖2中可以看出中國大陸夏季的NDVI的空間分布整體呈現(xiàn)東南和東北高、西北低的格局。在沿海的華北平原、長江三角洲和珠江三角洲等經濟發(fā)達地區(qū)呈現(xiàn)低NDVI值的現(xiàn)象，這是由于人類自身發(fā)展的需要，NDVI值隨著城鎮(zhèn)化的發(fā)展而降低;在新疆的阿爾泰山南坡、天山山脈北坡的NDVI值卻反常高，最高可達到0.896，是因為在準格爾盆地西部有缺口，大西洋和北冰洋的水汽能夠影響到兩側的山脈，形成迎風坡，對于植被生長有利;另外，在西北新疆塔里木盆地、準格爾盆地、甘肅西北部、青海西北部、內蒙古西部有沙漠、荒漠戈壁等，西藏自治區(qū)內的高原有高原荒漠和裸巖，氣候惡劣，不利于植物生長的環(huán)境，所以NDVI值普遍偏低，平均值為0.275。華北平原主要有農田，冬小麥收割之后，夏季播種農作物還未完全長成，摻雜土壤背景信息，因此NDVI值較低，平均值為0.471。位于華東地區(qū)的浙江、福建有丘陵，貴州、重慶、湖南和湖北西部有山地適宜森林生長，NDVI值較高，平均值為0.668。內蒙古東北部和東北三省內有大興安嶺、小興安嶺和長白山山脈，夏季森林和農作物生長旺盛NDVI值較高，平均值為0.775。

圖2 2018年夏季中國大陸NDVI平均值

表4 2018年夏季中國大陸六大分區(qū)NDVI均值

由圖3和表5，6可以看出過去38 a中，有49.923%的NDVI穩(wěn)定不變，49.846%的NDVI輕微改善，中國大陸地區(qū)整體植被沒有發(fā)生劇烈變化。西北地區(qū)NDVI趨勢整體無明顯變化(0.000 098/a)，但新疆天山山脈北坡和塔里木河流域處發(fā)現(xiàn)NDVI的變化較為劇烈，最大為0.030 272/a，呈上升趨勢;內蒙古陰山南、甘肅省內祁連山、寧夏賀蘭山等NDVI變化為上升趨勢;在經濟發(fā)達的長江三角洲、珠江三角洲、京津翼地區(qū)NDVI趨勢呈現(xiàn)下降現(xiàn)象。在西藏藏南地區(qū)、內蒙古大興安嶺、東北三省內的小興安嶺和長白山，夏季NDVI值較高，而年際變化趨勢中呈現(xiàn)下降趨勢，說明植被覆蓋濃郁，多年以來的植被覆蓋有所下降。

表5 1981-2018年六大分區(qū)NDVI趨勢均值

表6 1981-2018年NDVI變化趨勢分類級別

圖3 1981-2018年中國大陸夏季NDVI變化趨勢

3.2 主要生態(tài)系統(tǒng)類型的NDVI變化分析

中國陸地生態(tài)系統(tǒng)具體分為農業(yè)、森林、草地、水體與濕地、聚落、荒漠和其他共7種生態(tài)系統(tǒng)類型。農田和聚落受人類活動影響較大，近30 a來農作物總產和單產都有所提高[22]，對于全國范圍內植被NDVI增加是正向作用，但想要更詳細地分析中國大陸的生態(tài)狀況則需要中國陸地生態(tài)系統(tǒng)分類數據的支持。自然生態(tài)系統(tǒng)有森林、草地、水體與濕地、荒漠或者其他生態(tài)系統(tǒng)，這些生態(tài)系統(tǒng)類型受人類影響較小，對生態(tài)系統(tǒng)狀況有指示意義，按表1對上述生態(tài)系統(tǒng)分區(qū)統(tǒng)計其NDVI長時間序列變化，有助于中國生物多樣性保護和開展生態(tài)環(huán)境狀況評價。

森林和農田生態(tài)系統(tǒng)在過去38 a中NDVI呈現(xiàn)波動式上升，上升趨勢分別為0.001 4/a，0.001 1/a，森林生態(tài)系統(tǒng)NDVI最低值為1982年0.692，最高值為2017年0.766。1987年5月6日至6月2日大興安嶺發(fā)生特大森林火災;2000年6月17，18日黑龍江大興安嶺呼中區(qū)呼源、松嶺區(qū)發(fā)生特大森林火災，同年8月6日新疆裕民縣巴旦杏林保護區(qū)發(fā)生重大森林火災;2002年7月28日和2004年6月22日內蒙古大興安嶺北部原始林區(qū)分別發(fā)生特大、重大森林大火;2018年6月1日和2日內蒙古大興安嶺汗馬國家級自然保護區(qū)和北部原始林區(qū)發(fā)生森林火災等，這些森林大火使NDVI值在時間序列上呈現(xiàn)出偏低或者谷。在過去38 a中草地生態(tài)系統(tǒng)的NDVI在0.358至0.402間變化，草地生態(tài)系統(tǒng)NDVI變化較為劇烈，說明NDVI對于不同覆蓋度草地的變化敏感，1995年、1997年錫林浩特市草原分別發(fā)生3次和1次特大火災，該年份NDVI值有所下降，草地生態(tài)環(huán)境狀況變差?；哪鷳B(tài)系統(tǒng)的NDVI在過去38 a中呈下降趨勢(-0.000 6/a)，最低為0.094，最高為0.126，荒漠生態(tài)系統(tǒng)主要包含沙地、戈壁、鹽堿地、高寒荒漠，這些土地類型環(huán)境惡劣，十分脆弱，NDVI能夠反映出其變化狀況，在1981—2001年NDVI幾乎沒有大的變化，處于相對穩(wěn)定的狀態(tài)，而在2002—2018年NDVI值的上升趨勢為0.000 9/a，說明近17 a的生態(tài)環(huán)境呈現(xiàn)向好的態(tài)勢發(fā)展。其他生態(tài)系統(tǒng)的NDVI在過去38 a中呈現(xiàn)下降趨勢(-0.000 7/a)，最低為0.138，最高為0.184，與荒漠生態(tài)系統(tǒng)類似，在1981—2001年，NDVI也變化不大，2002—2018年，其NDVI的上升趨勢與荒漠生態(tài)系統(tǒng)一樣，生態(tài)系統(tǒng)呈現(xiàn)向好態(tài)勢。水體與濕地生態(tài)系統(tǒng)NDVI在過去38 a中呈現(xiàn)下降趨勢(-0.001 7/a)，最低為0.311，最高為0.401，在1981—2001年NDVI呈現(xiàn)較穩(wěn)定狀態(tài)，在2002—2018年中呈現(xiàn)上升趨勢(0.001 6/a)。

圖4 1981-2018年森林和農田生態(tài)系統(tǒng)NDVI變化

4 討 論

本研究利用NDVI探究生態(tài)系統(tǒng)狀況，但可能有所欠缺，另如氣候、地形等多種因素可以加以綜合考慮;趨勢線分析方法能很好地體現(xiàn)NDVI的變化趨勢，但是如何尋找合適的NDVI變化趨勢分級區(qū)間，需要進一步研究。劉可等[13]僅使用2010年一時期的中國陸地生態(tài)系統(tǒng)類型數據統(tǒng)計近30 a GIMMS NDVI3g數據，指出森林、農田、草地和水體與濕地總體NDVI非穩(wěn)定上升，本文利用生態(tài)系統(tǒng)類型數據統(tǒng)計重構后的NDVI數據，指出森林和農田生態(tài)系統(tǒng)NDVI整體呈增強趨勢，草地生態(tài)系統(tǒng)NDVI呈非穩(wěn)定上升趨勢，這可能是因為本文利用的是多期生態(tài)系統(tǒng)類型數據。利用生態(tài)系統(tǒng)分區(qū)統(tǒng)計不同年份的NDVI數據時，雖然盡可能利用生態(tài)系統(tǒng)數據，但最終統(tǒng)計的結果可能與實際的有差別，荒漠、水體與濕地和其他生態(tài)系統(tǒng)的NDVI曲線都在2000—2001年處有下降的陡坡，這可能是因為對GIMMS NDVI數據進行重采樣，仍無法解決其分辨率粗的問題，但若將MODIS NDVI數據進行降分辨率，則會損失大量信息。

5 結 論

(1) 構建數據集時，修正后GIMMS NDVI和MODIS NDVI相似的年份和月份變化規(guī)律;中國大陸地區(qū)MODIS NDVI在2000年、2001年的冬季數據質量不佳，可以從2001年7月開始選用數據。中國大陸地區(qū)夏季NDVI整體分布特征呈現(xiàn)西北低，東南和東北高，在經濟發(fā)達區(qū)域，由于人類經濟社會活動因素出現(xiàn)低NDVI值現(xiàn)象。

(2) 在過去38 a中，對受人類活動影響較小的生態(tài)系統(tǒng)，森林生態(tài)系統(tǒng)植被呈上升趨勢，森林火災對森林的影響在NDVI序列變化能夠反映出來;草地生態(tài)系統(tǒng)呈現(xiàn)穩(wěn)定趨勢，NDVI變化較為劇烈，NDVI序列可以作為監(jiān)測草地生態(tài)狀況變化的指示;荒漠、水體與濕地、其他生態(tài)系統(tǒng)的都呈下降趨勢，三者在1981—2001年NDVI變化都較為穩(wěn)定，且2002年以后NDVI呈現(xiàn)向好態(tài)勢發(fā)展，這與MODIS傳感器分辨率和波段寬度設置優(yōu)于AVHRR傳感器有關。