青藏高原NDVI變化趨勢及其對氣候的響應

孟 夢, 牛 錚, 馬 超, 田海峰, 裴 杰

(1.中國科學院 遙感與數字地球研究所, 遙感科學國家重點實驗室, 北京 100101; 2.中國科學院大學, 北京 100049; 3.河南理工大學 測繪與國土信息工程學院, 河南 焦作 454000)

青藏高原總面積約250萬km2，約占中國陸地總面積的1/4，平均海拔在4 000 m以上，是中國面積最大、世界上海拔最高的高原，被稱為世界的“第三極”[1]。其特殊的中低緯度地理位置、地貌條件以及巨大的海拔高度，形成了濕潤、干旱、熱帶和寒帶等多種氣候類型，幾乎囊括了東南亞所有類型的生態系統[2]。青藏高原具有遼闊的生態脆弱區，生態系統敏感而脆弱，且受人類影響較小，這也成為其受廣大生態學研究者鐘愛的一個重要原因[3]。

近年來，在全球變化和人為因素的共同作用下，青藏高原生態環境產生了一系列負面變化，具體表現為濕地草地面積減小、凍土退化、土地沙漠化嚴重、水域減少、多種動植物瀕臨滅絕以及自然災害增多等現象[4]。青藏高原作為中國乃至亞洲地區的生態安全屏障，其“一舉一動”都備受社會和學術界關注，因此科學有效的評判青藏高原生態環境變化以及探討其變化規律成為許多學者研究的重要課題[5]。植被作為高原生態系統中重要的組成部分，植被覆蓋變化在一定程度上能夠反映生態環境的變化，因此植被覆蓋變化研究對青藏高原生態環境變化規律的探索具有深刻意義。

歸一化植被指數(Normalized Difference Vegetation Index，NDVI)對植被的生長量和長勢十分敏感，是反映地表植被生長狀況的重要指標，因而在定性和定量評價植被覆蓋及其生長活力時被作為重要的信息源[6-7]。其定義為近紅外波段與可見光紅外波段波譜數值之差與兩波段波譜數值之和的比值，計算公式為NDVI=(NIR-IR)/(NIR+IR)，NDVI能較好的反映植被生長狀況，且可以消除地形、部分太陽角以及云等因素的影響[8-9]，其次，NDVI對干旱和半干旱區地區的低值被覆蓋信息十分敏感[10]，因此在低值被覆蓋區用NDVI反映植被信息得到了廣泛應用。

近年來，很多學者利用NDVI研究青藏高原植被的動態變化及植被與溫度、降水等氣候因子的關系。李艷芳等[11]基于SPOT-Vegetation植被指數數據，得出了青藏高原NDVI總體呈上升趨勢，季節變化趨勢不明顯的結論。向波等[12]基于全球歸一化植被指數(Glob-NDVI)發現近年來青藏高原植被向好的方面發展的現象。楊元合等[13]利用NOAA-AVHRR NDVI數據，發現青藏高原植被指數在生長季階段顯著增加，春季對NDVI增長率貢獻最大，且不同植被類型在不同的生長季NDVI增長率也有所區別。王青霞等[14]基于1982—2006年GIMMS NDVI數據，研究表明青藏高原植被與溫度降水等氣候因子具有相關性，不同等級的植被覆蓋度對溫度降水的響應程度是不同的，植被生長較好的區域對氣候響應更為明顯。楊保等[15]利用冰芯中的花粉以及氧同位素記錄，恢復了青藏高原過去500多年的植被演變與氣候變化趨勢，研究表明當高原為暖干氣候時，植被生長范圍擴大，高原氣候變冷變濕則植被生長范圍縮小。陸晴等[16]基于GIMMS NDVI和地面氣象站臺觀測數據，發現青藏高原植被與溫度降水的相關關系具有空間差異性，在溫度較高、熱量條件較好的高原東北部地區，降水為植被NDVI變化的主導因素；在降水充足的高原東中部地區，溫度為植被NDVI變化的主導因子。

盡管前人做了大量的研究，但多以行政邊界劃分研究區域，并且研究區中包括了西藏自治區東南部、四川省西南部、云南省西北部和青海省中部植被十分茂盛的地區，且多數研究所用遙感數據和氣象數據時序短，因此不能突出地反映青藏高原高原植被的時空變化。鑒于此，本文以地理界線劃分研究區，突破人為劃定的行政邊界，選取青藏高原海拔4 000 m以上的生態極端脆弱敏感區作為本次的研究區域，此外，選用長時間序列1982—2012年GIMMS NDVI數據、1982—2012年溫度降水數據為數據源，結合多種統計方法，分析青藏高原植被NDVI的變化趨勢及NDVI對氣候的響應，旨在揭示全球變化背景下青藏高原植被變化趨勢以及氣候對高原植被生長的影響機制。

1　研究區與數據

1.1　研究區概況

青藏高原雄踞亞洲東南，西起帕米爾高原，東到橫斷山脈，北以昆侖山、阿爾金山和祁連山為鄰，南抵喜馬拉雅山脈，西南邊緣部分分屬印度、巴基斯坦、尼泊爾、錫金、不丹及緬甸等國，其地理位置大致為27°17′—40°19′N，68°52′—102°43′E(圖1)。青藏高原氣溫隨緯度和海拔的升高而降低，氣溫日差大，年平均氣溫和降水量由東南向西北遞減。該區域的氣候變化不僅對北半球有巨大的影響，甚至對全球氣候變化也具有顯著的超前性和敏感性。

1.2　數據及數據預處理

1.2.1數 據文中所用主要數據：GIMMS NDVI數據，來源于美國國家航空航天局(https:∥ecocast.arc.nasa.gov/data/pub/)，時間分辨率為15 d，空間分辨率為8 km×8 km，時間跨度為1982—2012年，31 a的數據共包括了744期半月合成圖像。GIMMS NDVI數據經過輻射定標、大氣校正、除云、除壞線等處理，是目前研究干旱、半干旱區植被變化質量比較可靠的NDVI數據集[17]。氣象數據來源于中國氣象科學數據共享服務網(http:∥cdc.nmic.cn)提供的氣候資料年值數據集。高程數據來源于GTOP30項目所提供的DEM數據，分辨率為90 m。

1.2.2數據預處理首先利用青藏高原矢量數據對GIMMS NDVI數據批量裁剪，其次采用最大值合成法(Maximum Value Composite，MVC)獲取青藏高原月度、年度NDVI數據，最后對得到的月度、年度NDVI值進行計算，獲取研究區月度、年度NDVI均值并統計入Excel表格。

位于青藏高原地區的氣象站點數目極少，本文選取該數據集青藏高原地區102個地面氣象站臺的氣象數據(1982—2012年)進行克里金插值，獲取分辨率為8 km×8 km的年平均氣溫和年平均降水量圖像，從而獲取青藏高原地區年氣溫、年降水量信息。

使用ArcGIS軟件中的雙線性插值法將DEM數據重采樣為和NDVI相同的分辨率，即8 km×8 km，其次在Global Mapper中生成海拔為4 000 m的等高線，確定本次研究區范圍，即青藏高原海拔高于4 000 m的區域。

圖1研究區域及氣象站點分布

2　研究方法

(1) 最大值合成法。為了研究青藏高原地區植被在空間上的整體變化特征，采用最大值合成法獲取研究區月NDVI值和年NDVI值，最大值合成法可以有效減少大氣中氣溶膠、云陰影、太陽高度角等因素的影響[18]。計算公式如下：

(1)

式中：MNDVIi為第i月NDVI最大值；i=1,2，…，12為月序號；NDVI1為第i月上半月NDVI值；NDVI2為第i月下半月NDVI值。以求取的月NDVI值為基礎求取年NDVI。

(2) 植被變化分析。結合植被綠度變化率(Greenness Rate of Change,GRC)[19]分析青藏高原植被變化趨勢。GRC定義為在一定時間范圍內采用最小二乘法擬合年度NDVI均值的斜率，它能夠反映出每個柵格的變化趨勢[20]。計算公式如下：

(2)

式中：變量i為年份序號，取值范圍為1～31；n為研究的時間序列長度，n=31；YMNDVIi為第i年的年度NDVI平均值。GRC>0說明NDVI在31 a間的變化趨勢是增加的，反之則是減少的[21]。

(3) 植被覆蓋度分析。計算公式如下[22]：

(3)

(4) NDVI變化率分析。由下式計算：

(4)

式中：直線斜率表示31 a(1982—2012)間NDVI年平均值與年份之間求直線回歸所得直線的斜率，均值為31 a的年度NDVI平均值，該變化率在數值上相當于初期值除以研究期間年平均NDVI的末期值與初期值之差[22]。

(5) 相關性分析。此外，考慮到植被生長與地形有密切的關系，借助GTOP 30DEM數據，分析了青藏高原歸一化植被指數隨地形的變化情況。

3　結果與分析

3.1　青藏高原NDVI時間序列變化特征

圖2　青藏高原月度NDVI變化趨勢

3.1.2年度NDVI變化將利用最大值合成法獲取的的年NDVI進行計算得到年NDVI均值，將其與年份線性回歸，得到青藏高原31 a NDVI均值變化趨勢(圖3)。可知31 a來青藏高原總體年NDVI均值呈緩慢增加趨勢(slope=0.000 2)，利用公式(4)可計算得到青藏高原31 a年NDVI變化率為0.029，表明該研究區植被向好的方向發展，與王敏等[25]得出的結論一致。

圖3　青藏高原年度NDVI變化趨勢

植被覆蓋度是單位面積內植被(包括葉、莖、枝)垂直投影面積所占的百分比[26]，它是陸面過程模式中最重要的參數[27]。植被覆蓋度不僅可以表示植被進行光合作用面積的大小，它在一定程度上對植被的茂密程度、生長態勢以及區域生態系統環境變化也有所指示[28]。由公式(3)計算得到青藏高原31 a植被覆蓋度(圖4)，將其與年份線性回歸，可知青藏高原31 a間植被覆蓋度呈0.000 3/a 的速率緩慢增長，其中在1988年、1994年、2000年、2010年有大幅度的躍升。認為全球變化對青藏高原起到了正面影響，植被長勢逐漸向好的方向發展，植被覆蓋度相應增加[29]。

3.2　青藏高原NDVI對氣候因子的響應

陸氣與氣候間相互作用的變化，通過降水機制改變和氣溫變化等對植被的生長和分布產生重要影響。圖5A、圖5B分別展示了青藏高原1982—2012年31 a年平均氣溫以及年平均降水量的變化趨勢。

圖4　青藏高原31 a植被覆蓋度變化趨勢

圖5A展示了青藏高原1982—2012年共31 a的年平均降水量變化趨勢，年最低降水量為255.2 mm(1994年)，最高降水量為343.7 mm(2010年)，31 a平均降水量為309.9 mm。1982—1995年，年降水量呈下降趨勢，1995—2012年，年降水量呈上升趨勢，31 a青藏高原地區總體降水量呈逐漸增加趨勢，變化率為0.223 mm/10 a。

由圖5B可以看出，近31 a青藏高原海拔高于4 000 m的區域平均氣溫呈明顯增加趨勢，年最低氣溫為4.82℃(1983年)，最高氣溫為6.99℃(2009年)，增長速率為0.018℃/10 a。氣溫增長速率在1995年前后存在轉折點，1982—1995年年平均氣溫增長速率為0.037℃/10 a，1995—2012年年平均氣溫增長速率為0.029℃/10 a，宋辭等[30]認為青藏高原氣溫發生突變，整體呈波動上升的趨勢。

圖5　青藏高原31 a降水、溫度變化趨勢

對青藏高原地區NDVI與降水、溫度做相關性分析，見圖6A,6B所示。研究區年均NDVI與降水、溫度的皮爾遜相關系數分別為0.183，0.281，均呈現不顯著的弱相關。相關系數表明該地區溫度對植被的影響大于降水對植被的影響。此外，對降水、溫度與NDVI做超前滯后分析，發現降水對研究區NDVI有三年的滯后效應，與劉雅勤等[31]的研究結果保持一致，當降水向后移動三年時，降水與NDVI的相關系數達到最大值0.387。溫度對植被指數滯后性影響不明顯。

由圖6A，6B可以看出，青藏高原地區年NDVI均值在2000—2008年呈現明顯下降趨勢。對此，分別對2000—2008年的年NDVI均值、年降水量、年氣溫與年份做一元線性回歸分析(圖7)。由圖7可以看出，2000—2008年年降水量呈下降趨勢，年氣溫呈升高趨勢，溫度的升高及降水量的減少有可能打破了青藏高原地區降水量與蒸發量之間的平衡狀態，使青藏高原干旱地區更加干旱。已有研究表明[32]，干旱對青藏高原高山草甸植被的生長具有抑制作用，加劇了植被的退化，因此對總體區域植被NDVI產生了重要影響。

此外，在2008—2012年年NDVI均值呈增長趨勢，對于這一變化趨勢，考慮到時間序列短、所用氣象數據和NDVI數據集分辨率低、以及青藏高原自然環境復雜等因素，因此尚不能對此情況做出合理的解釋。

圖6　青藏高原地區NDVI、降水相關曲線NDVI(A)、溫度相關曲線(B)

圖7　青藏高原地區(2000-2008)年NDVI、降水量、氣溫變化趨勢

3.3　青藏高原NDVI與海拔的相關性分析

高原植被在全球變化中發揮著重要的指示作用，因此備受關注。高原植被不僅對溫度、降水等氣候因子響應敏感，而且受海拔高度的影響也很大[33]。在研究區2000年的NDVI影像中均勻選取100個樣本點，分別統計它們所對應的NDVI值和高程值，將NDVI值和高程值映射到以NDVI值為縱坐標，海拔高度為橫坐標的的坐標系中(圖8)。結果表明整體NDVI與海拔高度呈現相關程度不是很強的負相關，二者皮爾遜相關系數為-0.309，究其原因可能是植被類型過多，且植被指數受植被覆蓋度、海拔高度、氣候因子等綜合影響，因此使海拔高度的影響變的不是太明顯[34]。

在ArcGIS軟件中用2000年的NDVI影像生成間距為0.1的等值線(圖9)，并將等值線與青藏高原DEM影響疊加，可知在青藏高原中部、東南部、以及西北部海拔較低的區域NDVI等值線分布比較密集，印證了NDVI與海拔相關的事實。

圖8　青藏高原地區NDVI與DEM相關曲線

圖9　青藏高原地區NDVI等值線與DEM疊加

4　結 論

(1) 在全球變化影響下，1982—2012年各研究區NDVI值呈現以年為周期的變化特點，年NDVI值表明，植被的生長旺季處于每年的7月、8月、9月份，且年NDVI值總體呈增長趨勢，年NDVI變化率為0.029，表明青藏高原地區植被生長狀況發展良好。

(2) 青藏高原地區1982—2012年植被覆蓋度呈0.000 3/10 a的速率增長，認為全球變化對青藏高原起到了正面影響，植被長勢逐漸向好的方向發展，植被覆蓋度相應增加。

(3) 青藏高原地區溫度、降水都呈增長趨勢，在1995年之前，年降水量呈下降趨勢，1995年之后，年降水量呈增長趨勢。溫度在1995年之后增長趨勢比1995年之前慢，總體增長趨勢為前高后低。

(4) 研究區NDVI與溫度降水的皮爾遜相關系數分別為0.281，0.183，均呈現不顯著的弱相關，降水對歸一化植被指數的影響有三年的滯后效應，溫度對植被指數滯后性影響不明顯。

(5) 植被生長受海拔高度影響明顯，NDVI與海拔高度皮爾遜相關系數為-0.309，青藏高原中部、東南部、以及西北部海拔較低的區域植被長勢相對較好，且NDVI值相對較高。

此外，青藏高原地區地面氣象站臺分布稀疏，即使我們采用了國際相對認可的克里金插值方法，但氣象站臺稀疏地區插值數據的不穩定性仍然是存在的，這種不穩定性對于后期的趨勢分析會產生影響，因此，采用更多的氣象觀測數據來提高青藏高原地區氣候插值的結果精度是十分必要的。

參考文獻：

[1]鄭度,趙東升.青藏高原的自然環境特征[J].科技導報,2017,35(6):13-22.

[2]莫申國,張百平,程維明,等.青藏高原的主要環境效應[J].地理科學進展,2004,23(2):88-96.

[3]梁四海,陳江,金曉媚,等.近21年青藏高原植被覆蓋變化規律[J].地球科學進展,2007,22(1):33-40.

[4]姚檀棟,朱立平.青藏高原環境變化對全球變化的響應及其適應對策[J].地球科學進展,2006,21(5):459-464.

[5]孫鴻烈,鄭度,姚檀棟,等.青藏高原國家生態安全屏障保護與建設[J].地理學報,2012,67(1):3-12.

[6]田慶久,閔祥軍.植被指數研究進展[J].地球科學進展,1998,13(4):327-333.

[7]孫紅雨,王長耀,牛錚,等.中國地表植被覆蓋變化及其與氣候因子關系：基于NOAA時間序列數據分析[J].遙感學報,1998,2(3):204-210.

[8]韓瑞梅,劉培,馬超,等.鄂爾多斯植被的NDVI3g動態及氣候響應[J].水土保持通報,2016,36(5):28-33.

[9]李春暉,楊志峰.黃河流域NDVI時空變化及其與降水/徑流關系[J].地理研究,2004,23(6):753-759.

[10]刁鳴軍,夏朝宗.1982—2013年準噶爾盆地植被長勢變化分析[J].林業資源管理,2016(5):39-46.

[11]李艷芳,孫建.青藏高原NDVI時空變化特征研究(1982—2008)[J].云南農業大學學報,2015,30(5):790-798.

[12]向波,繆啟龍,高慶先.青藏高原氣候變化與植被指數的關系研究[J].高原山地氣象研究,2001,21(1):29-36.

[13]楊元合,樸世龍.青藏高原草地植被覆蓋變化及其與氣候因子的關系[J].植物生態學報,2006,30(1):1-8.

[14]王青霞,呂世華,鮑艷,等.青藏高原不同時間尺度植被變化特征及其與氣候因子的關系分析[J].高原氣象,2014,33(2):301-312.

[15]楊保,唐領余,李春海,等.近550年來青藏高原中部植被演化與氣候變化研究[J].科學通報,2010,55(21):2139-2147.

[16]陸晴,吳紹洪,趙東升.1982—2013年青藏高原高寒草地覆蓋變化及與氣候之間的關系[J].地理科學,2017,37(2):292-300.

[17]張戈麗,徐興良,周才平,等.近30年來呼倫貝爾地區草地植被變化對氣候變化的響應[J].地理學報,2011,66(1):47-58.

[18]郭鵬,徐麗萍.基于GIMMS—NDVI的新疆植被覆蓋時空變化[J].水土保持研究,2014,21(4):97-100.

[19]陳志彪,朱鶴健.不同水土流失治理模式下的土壤理化特征[J].福建師大學報:自然科學版,2006,22(4):5-9.

[20]仝莉棉,曾彪,王鑫.山西省不同生態區NDVI時空變化及其影響因素[J].水土保持研究,2016,23(3):71-76.

[21]張戈麗,歐陽華,張憲洲,等.基于生態地理分區的青藏高原植被覆被變化及其對氣候變化的響應[J].地理研究,2010,29(11):2004-2016.

[22]方精云,樸世龍,賀金生,等.近20年來中國植被活動在增強[J].中國科學,2003,33(6):554-565.

[23]劉雙俞,張麗,王翠珍,等.基于MODIS數據的青藏高原植被物候變化趨勢研究(2000年—2010年)[J].遙感信息,2014，29(6):25-30.

[24]李娜,謝虹.青藏高原NDVI對ENSO事件的響應[J].遙感技術與應用,2014,29(5):818-822.

[25]王敏,周才平,吳良,等.近10a青藏高原干濕狀況及其與植被變化的關系研究[J].干旱區地理,2013,36(1):49-56.

[26]郭敏杰,張亭亭,張建軍,等.1982—2006年黃土高原地區植被覆蓋度對氣候變化的響應[J].水土保持研究,2014,21(5):35-40.

[27]何建軍,余曄,陳晉北,等.植被覆蓋度對蘭州地區氣象場影響的模擬研究[J].高原氣象,2012,31(6):1611-1621.

[28]穆少杰,李建龍,陳奕兆,等.2001—2010年內蒙古植被覆蓋度時空變化特征[J].地理學報,2012,67(9):1255-1268.

[29]張憲洲,楊永平,樸世龍,等.青藏高原生態變化[J].科學通報,2015,60(32):3048-3056.

[30]宋辭,裴韜,周成虎.1960年以來青藏高原氣溫變化研究進展[J].地理科學進展,2012,31(11):1503-1509.

[31]劉雅勤,范廣洲,周定文,等.青藏高原冬、春植被歸一化指數變化特征及其與高原夏季降水的聯系[J].氣象學報,2007,65(6):959-967.

[32]王敏,周才平,吳良,等.近10a青藏高原干濕狀況及其與植被變化的關系研究[J].干旱區地理,2013,36(1):49-56.

[33]魏偉,趙軍,王旭峰.天祝高寒草原區NDVI, DEM與地表覆蓋的空間關系[J].干旱區研究,2008,25(3):394-401.

[34]萬紅,李希燦,萬劍華,等.青藏高原地區不同植被類型MODIS植被指數與海拔高度的關系研究[J].干旱區資源與環境,2013,27(11):39-44.