伍良旭, 王 晗, 邵懷勇, 周思彤

(成都理工大學 地球科學學院, 成都 610059)

植被是地球生態(tài)系統(tǒng)主要的組成部分，是連接氣候與土壤的紐帶，在全球物質(zhì)能量交換和氣候變化中具有不可替代的作用[1]。植被動態(tài)變化在高原生態(tài)環(huán)境中具有敏感的指示作用，合理利用區(qū)域資源和有效治理生態(tài)環(huán)境提供依據(jù)，已成為生態(tài)環(huán)境研究的熱點[2]。NDVI是表征植被動態(tài)變化的一個重要指標，能較精確地反映地表植被的生長狀態(tài)及演變特征，在區(qū)域植被退化、濕地景觀格局等領(lǐng)域發(fā)揮著重要的作用[3]。川西高原位于我國青藏高原向川西南山地和四川盆地的過渡地區(qū)，是我國生態(tài)安全戰(zhàn)略格局中青藏高原生態(tài)屏障和黃土高原—川滇生態(tài)屏障的重要組成部分，地形復雜，且受高原背風坡和西南季風等系統(tǒng)的共同影響，是典型的生態(tài)環(huán)境脆弱區(qū)和生態(tài)氣候敏感區(qū)[2]，但目前針對該區(qū)域植被時空格局與氣候因子的相關(guān)研究較少，因此對該地區(qū)的植被時空格局及其對氣候變化響應的研究具有重要意義。

此外，現(xiàn)在的遙感數(shù)據(jù)通常不能同時具備較長的時間序列和較高的空間分辨率，如GIMMS NDVI3g是目前最長時間序列的NDVI數(shù)據(jù)集之一[4]，前人利用其進行了1982—2009年東北地區(qū)[5]、1982—2006年蒙古高原[6]和1982—2013年新疆地區(qū)[7]的植被變化研究，但8 km的空間分辨率應用在區(qū)域植被變化研究中過于粗糙。雖然如今高分辨率輻射計(AVHRR)、地球觀測系統(tǒng)(SPOT)-VEGETATION，MODIS和ENVISAT中分辨率成像光譜儀(MERIS)等衛(wèi)星傳感器所采集的數(shù)據(jù)有較高的時間分辨率，通過合成一定時間內(nèi)的圖像，能夠減少甚至消除云污染[4]，遙感數(shù)據(jù)質(zhì)量不斷提高，使得區(qū)域植被趨勢分析在全球范圍內(nèi)得到廣泛應用[8-10]，但如何提高數(shù)據(jù)的空間分辨率仍是亟待解決的問題。相對于GIMMS NDVI3g數(shù)據(jù)精度更高、更能表達植被時序變化的MODIS NDVI產(chǎn)品，卻因其時間覆蓋范圍(2001年7月至今)無法針對長時間序列進行研究。因此，如何提高早期遙感數(shù)據(jù)的分辨率對于長時間序列的植被變化監(jiān)測研究顯得尤為重要。

針對遙感數(shù)據(jù)時空尺度的融合問題，前人做了大量的研究。如Gao等[11]提出的STARFM模型(spatial and temporal reflection fusion Model)可融合Landsat和MODIS地表反射率，生成具有高時空分辨率的合成反射率數(shù)據(jù)。Singh等[12]利用STARFM融合模型充分提取了農(nóng)作物8 a生長季節(jié)的物候。Bhandari等[13]和Schmidt等[14]評估了STARFM合成長期時間序列數(shù)據(jù)的能力，并獲取不同植被的物候特征進行研究。Meng等[15]和Zhu等[16]對STARFM算法改進，分別提出了時空融合自適應植被指數(shù)模型(STAVFM)和增強型自適應時空融合模型ESTARFM，以更準確地表征地表的時空變化。由于在地形復雜地區(qū)，ESTARFM的融合效果比STARFM較好，使得ESTARFM算法被頻繁用于生成高時空分辨率的合成遙感圖像[17-18]。然而，目前較少有學者針對川西高原使用ESTARFM模型將GIMMS NDVI3g與MODIS NDVI進行時空融合并開展1995年至今植被時空格局的相關(guān)研究，為此本文基于ESTARFM模型開展川西高原GIMMS NDVI3g和MODIS NDVI數(shù)據(jù)的時空融合，這對開展其長期的植被變化監(jiān)測具有重要意義。

本研究嘗試基于ESTARFM模型建立川西高原GIMMS NDVI3g與MODIS NDVI數(shù)據(jù)的函數(shù)關(guān)系，以提高GIMMS NDVI3g數(shù)據(jù)的空間分辨率。在此基礎(chǔ)上分析川西高原長時間的植被變化特征及其對氣候變化的響應，為當?shù)氐沫h(huán)境保護決策提供基礎(chǔ)數(shù)據(jù)和科學參考。

1 研究區(qū)概況

川西高原(27.59°—34.31°N，97.36°—104.62°E)地屬青藏高原與四川過渡地帶，覆蓋阿壩州和甘孜州等高原藏區(qū)，區(qū)域面積約為25.38萬km2。研究區(qū)內(nèi)地勢高，地形起伏大，平均海拔約4 000 m。境內(nèi)屬亞熱帶高原季風氣候，地域差異明顯，多年降水量約435.6～1 584.1 mm。近年來，由于全球氣候變化以及人類活動的增加，當?shù)刂脖桓采w逐漸變化，植被的退化對境內(nèi)生態(tài)系統(tǒng)安全會造成極大威脅[19]，因此，研究分析川西高原植被覆蓋時空分布特征，并掌握主要氣候因素對植被變化的影響非常必要[20]。

2 數(shù)據(jù)與方法

2.1 數(shù)據(jù)來源與處理

(1) 遙感數(shù)據(jù)：MODIS 13Q1植被指數(shù)數(shù)據(jù)集產(chǎn)品，空間分辨率500 m×500 m，時間分辨率30 d。首先進行Savitzky-Golay濾波處理降低異常值等噪聲影響，隨后利用最大值合成法合成2005年和2015年MODIS NDVI數(shù)據(jù)。該方法能減少云霧等大氣效應對數(shù)據(jù)的影響，提高數(shù)據(jù)精度[21]。

1995年、2005年和2015年的GIMMS NDVI3g數(shù)據(jù)空間分辨率為8 km×8 km，時間分辨率為15 d。采用最大值合成法合成3個時間段的年GIMMS NDVI3g數(shù)據(jù)，并將其空間分辨率重采樣為500 m×500 m，為后續(xù)使用ESTARFM模型進行數(shù)據(jù)融合做準備。

(2) 氣象數(shù)據(jù)：年均氣溫、年降水量空間插值數(shù)據(jù)集產(chǎn)品。薄板光順樣條法采用平滑樣條函數(shù)對多變量數(shù)據(jù)進行分析和插值，即使用函數(shù)逼近曲面，是一種較好的氣候要素插值的方法[22]。因此，該數(shù)據(jù)集采用基于薄板光順樣條法開發(fā)的ANUSPLIN插值軟件(Hutchinson，2001)對氣候要素進行插值，其中將氣象站點的經(jīng)度、緯度和海拔高程作為年降水量空間插值的自變量，而年均溫度空間插值則考慮了溫度隨垂直高度的變化[22]。

MODIS NDVI數(shù)據(jù)來源于中國科學院計算機網(wǎng)絡(luò)信息中心地理空間數(shù)據(jù)云平臺(http:∥www.gscloud.cn/)；GIMMS NDVI3g數(shù)據(jù)來源于ECOCAST網(wǎng)站(https:∥ecocast.arc.nasa.gov/data/pub/gimms/3 g.v1/)；氣象數(shù)據(jù)和高程數(shù)據(jù)來源于中國科學院資源環(huán)境科學數(shù)據(jù)中心(http:∥www.resdc.cn)；為提高試驗結(jié)果精度，轉(zhuǎn)換所有數(shù)據(jù)到統(tǒng)一的投影坐標系WGS_1984_Zone_48N。

2.2 研究方法

2.2.1 ESTARFM算法 為生成高時空分辨率的遙感圖像，使用ESTARFM模型融合MODIS NDVI和GIMMS NDVI3g數(shù)據(jù)。ESTARFM模型是在STARFM模型[11]基礎(chǔ)上根據(jù)像元的異質(zhì)性調(diào)整賦權(quán)方法，通過設(shè)置轉(zhuǎn)換系數(shù)來改進預測結(jié)果，異質(zhì)性越大，預測精度越高，可以保存更多的空間特征細節(jié)[23]。

本文根據(jù)第1時期m時刻MODIS NDVI與tp時刻的GIMMS NDVI3g得到融合后tp時刻的GIMMS NDVI3g數(shù)據(jù)，記為mm=(Xw/2,yw/2,tp,B)，或利用第2時期n時刻的觀測數(shù)據(jù)得到融合后的tp時刻GIMMS NDVI3g數(shù)據(jù)，記為mn=(Xw/2,yw/2,tp,B)。通過tm和tn與時刻tp間檢測到的GIMMS NDVI3g變化幅度來計算權(quán)重，如公式(1)所示[16]：

(k=m,n)

由公式(2)計算中心像元的值：

M(xw/2,yw/2,tp,B)=Tm·Mm(xw/2,yw/2,tp,B)+

Tn·Mn(xw/2,yw/2,tp,B)

(2)

由于篇幅限制，具體融合原理請見參考文獻[16]。統(tǒng)計2005年和2015年融合后的GIMMS NDVI3g和MODIS NDVI對應的像元值繪制散點圖(圖1)。融合的GIMMS NDVI3g數(shù)據(jù)與MODIS NDVI數(shù)據(jù)具有高度相似性，各像元值大致相同，兩者基本分布在y=x的趨勢線上(p<0.01)，預測效果較好。

圖1 真實MODIS NDVI與融合后的GIMMS NDVI3g數(shù)據(jù)散點圖

對比2005年和2015年GIMMS NDVI3g、融合后的GIMMS NDVI3g和真實MODIS NDVI數(shù)據(jù)細節(jié)(圖2)可以看出，融合后的GIMMS NDVI3g與真實MODIS NDVI數(shù)據(jù)像元分辨率較好，地物紋理清晰，部分特征地區(qū)的細節(jié)能夠凸現(xiàn)出來，能較好地再現(xiàn)區(qū)域植被覆蓋細節(jié)特征。

注：A為2005年GIMMS NDVI3g；B為2005年融合后的GIMMS NDVI3g；C為2005年MODIS NDVI；D為2015年GIMMS NDVI3g；E為2015年融合后的GIMMS NDVI3g；F為2015年MODIS NDVI。圖2 融合數(shù)據(jù)對比圖

2.2.2 一元線性回歸模型 基于融合處理后的GIMMS NDVI3g數(shù)據(jù)(值域為[-1,1])，使用像元二分模型計算植被覆蓋度，由于篇幅限制，具體的植被覆蓋度計算公式和原理請見參考文獻[24]。

植被覆蓋度變化主要用于分析植被覆蓋度隨時間變化的關(guān)系，他們的關(guān)系滿足一元線性回歸模型的條件，以局部像元時間變化特征表示整個空間的變化規(guī)律[24]，其計算方式見公式(3)：

(3)

式中:n為總監(jiān)測年數(shù);fi為第i年植被覆蓋度;Slope為植被覆蓋度線性擬合斜率。利用兩者的相關(guān)關(guān)系來判斷植被覆蓋度時間變化的顯著性,斜率為負表示植被覆蓋度下降,反之則表示植被覆蓋度上升,擬合斜率絕對值越大,變化越明顯[24-25]。

2.2.3 偏相關(guān)性分析法 地理系統(tǒng)是一個由多要素組成的復雜系統(tǒng)，偏相關(guān)分析可剔除其他變量的影響，只研究一個因素與另一個因素的相關(guān)程度[25]。本文基于像元對川西高原1995年、2000年、2015年植被覆蓋度與年降水量、年均氣溫兩個氣候因素進行偏相關(guān)性研究，見公式(4)[25]：

(4)

式中:rab,c表示控制變量c的值不變,變量a，b間的偏相關(guān)系數(shù)。rab，rbc，rbc分別表示變量a和b，b和c，a和c的簡單相關(guān)系數(shù)。其中a，b和c分別為植被覆蓋度、年均溫度和年降水量。rab,c>0時表示兩變量呈正相關(guān),rab,c<0時表示兩變量呈負相關(guān),rab,c的絕對值越大,兩變量的相關(guān)性越大。

簡單相關(guān)系數(shù)計算見公式(5)[25]：

(5)

3 結(jié)果與分析

3.1 植被覆蓋特征分析

3.1.1 植被覆蓋空間分布特征 參照前人研究[26]將川西高原植被覆蓋度分為極低覆蓋度(0～0.2)、低覆蓋度(0.2～0.4)、中覆蓋度(0.4～0.6)、中高覆蓋度(0.6～0.8)和高覆蓋度(0.8～1.0)5個等級，得到川西高原1995年、2005年和2015年植被覆蓋度空間分布圖(圖3)。

圖3 不同等級植被覆蓋度空間分布

分析可得，1995年、2005年和2015年研究區(qū)極低、低和中覆蓋度區(qū)域主要分布在西部高海拔的白玉、巴塘、理塘和東部低海拔的瀘定、康定、九龍、理縣等縣區(qū)，南部的稻城、和西北部的瑪多縣主要為中高覆蓋度，其余地區(qū)多屬于高覆蓋度區(qū)域。西部地區(qū)海拔較高，多為高山植被，西部德格、甘孜一帶和東部瀘定縣區(qū)有荒漠和草甸區(qū)，其余地區(qū)水熱條件良好，主要為針葉林和闊葉林。2005年川西高原高覆蓋度區(qū)域面積相對1995年有所減少，中植被覆蓋和中高植被覆蓋區(qū)域面積增加。2015年高植被覆蓋區(qū)域覆蓋面積繼續(xù)降低，但趨勢減緩，中、中高植被覆蓋區(qū)域面積增加，且后者相對前者增加明顯。

川西高原整體的植被覆蓋較高，呈現(xiàn)較為明顯的區(qū)域差異。高植被覆蓋區(qū)主要分布在北部和南部地區(qū)，低值區(qū)分布在中部地區(qū)，此外沿海拔梯度，植被覆蓋呈先升高后降低的垂直分布特點。

3.1.2 植被覆蓋變化趨勢特征 為研究1995—2015年川西高原植被覆蓋的時間變化規(guī)律，對植被覆蓋與時間進行一元線性回歸分析(圖4)，根據(jù)趨勢變化的方向和程度，劃分為5個等級(表1)[24]。植被覆蓋呈增長趨勢的面積占17.63%，其中明顯增加面積占15.53%，輕度增加面積占2.10%，主要分布在高海拔地區(qū)的瑪多、石渠、白玉等縣區(qū)，以瑪多縣分布尤為密集。植被覆蓋呈減少趨勢的面積占27.36%，其中明顯減少面積占18.29%，輕度減少面積占9.07%，減少區(qū)域廣泛分布在各個區(qū)域，以理塘、康定縣區(qū)分布較廣，結(jié)合前人研究分析認為[3]，理塘縣自20世紀80年代開始，因過度挖采、放牧等人為因素、鼠蟲害、氣候變化、暴雨沖刷等多種因素綜合并長期作用使理塘縣草原退化、沙源裸露導致植被覆蓋減少；康定縣由于亂砍濫伐，開荒種地，旅游開發(fā)等人為干擾因素和自然災害等多因素引起大面積土地退化，山坡地水土流失十分嚴重，植被減少趨勢明顯。植被覆蓋基本穩(wěn)定的區(qū)域面積占55.01%，分布廣泛，主要在研究區(qū)的東北部的阿壩、若爾蓋縣區(qū)，該地區(qū)主要為沼澤、草甸等，植被生長狀況穩(wěn)定，植被變化不明顯?？傮w而言，1995—2015年川西高原植被覆蓋呈增長趨勢面積占17.63%，呈增長趨勢面積占27.36%，基本穩(wěn)定面積占55.01%，植被整體變化趨勢穩(wěn)定。

圖4 川西高原1995-2015年植被覆蓋變化趨勢

表1 1995-2015年川西高原植被覆蓋變化面積及比例

3.2 植被對氣候變化的響應

3.2.1 年均氣溫與年降水量空間分布特征 植被的變化受到多方面因素的影響，其中年降水量和年均氣溫是對植被變化影響較大的兩個氣候因素[25]。計算1995年、2005年和2015年的年均氣溫和年降水量(圖5)。

圖5 1995-2015年年均氣溫和年降水量分布

分析可得，1995年、2005年和2015年川西高原年均氣溫為3.3℃，最高溫度為17.9℃，位于最南部的鹽源縣，最低溫度為-13.1℃，位于西北部的瑪多縣；年均降水量為848.6 mm，最大降水量為1 584.1 mm，位于西南部的瀘定縣；最小降水量為435.7 mm，位于西北部的瑪多縣。研究區(qū)內(nèi)年均氣溫和年降水量的分布具有空間異質(zhì)性，瑪多縣所在的西北部和瀘定縣所在的東南部區(qū)域氣溫整體低于中部地區(qū)，瀘定縣所在的東南地區(qū)，位于青藏高原和四川盆地過渡區(qū)，受東南季風和西南季風影響，降水充沛，氣溫不高；瑪多縣地區(qū)靠近青藏高原，從研究區(qū)中部到西北部地勢逐漸增高，氣溫逐漸下降。

3.2.2 植被覆蓋對氣候變化的響應 計算1995—2015年年均氣溫和年降水量與植被覆蓋的偏相關(guān)系數(shù)(圖6A—B)。年植被覆蓋與年均氣溫的偏相關(guān)系數(shù)變化范圍為[-0.71，0.71]，年植被覆蓋與年降水量的偏相關(guān)系數(shù)變化范圍為[-0.71，0.71]。統(tǒng)計分析可得年植被覆蓋與年均氣溫呈正相關(guān)的區(qū)域面積占14.9%，年植被覆蓋與年降水量呈正相關(guān)的區(qū)域面積占35.4%。根據(jù)年均氣溫和年均降水量與年植被覆蓋的偏相關(guān)系數(shù)方向，對研究區(qū)域進行分區(qū)討論[27]，年植被覆蓋與年均氣溫負相關(guān)，與年降水量負相關(guān)，歸于分區(qū)一；年植被覆蓋與年均氣溫負相關(guān)，與年降水量正相關(guān)，歸于分區(qū)二；年植被覆蓋與年均氣溫正相關(guān)，與年降水量負相關(guān)，歸于分區(qū)三；年植被覆蓋與年均氣溫正相關(guān)，與年降水量正相關(guān)，歸于分區(qū)四；年植被覆蓋與年均氣溫和年降水量的偏相關(guān)系數(shù)較低(都為0)，歸于分區(qū)五(圖6C)。

圖6 1995-2015年年均氣溫和年降水量與植被覆蓋度偏相關(guān)分析及偏相關(guān)性分區(qū)

分區(qū)一占總面積的6.4%，主要分布在石渠縣、甘孜縣和色達縣的北部，植被類型主要為高寒草甸，分析認為該區(qū)域鼠災泛濫和人為過度放牧導致植被覆蓋稀疏，氣溫和降水反而加速了水土流失，抑制植被生長。分區(qū)二占總面積的30.8%，主要分布在西部的白玉縣、理塘縣和東北部的黑水縣等地，植被類型主要為草甸、常綠闊葉林、針葉林和灌叢，地形復雜，氣溫較高而降水量低，水土蒸發(fā)量大，水土流失較嚴重[28]，氣溫升高減弱其光合作用，抑制植被的生長，適當?shù)慕邓稍黾油寥篮浚龠M植被生長。分區(qū)三占總面積的10.3%，主要分布在馬爾康縣、道孚縣、九龍縣一帶，該區(qū)域氣溫相對較低而降水量高，降水持續(xù)增加引發(fā)泥石流、滑坡等地質(zhì)災害，并使植被光照條件下降，減弱其光合作用，抑制植被生長，從而使植被覆蓋降低[29-30]。分區(qū)四占總面積的4.6%，主要分布在西北部的瑪多縣及東南部的汶川縣和康定縣地區(qū)。汶川地區(qū)自2008年汶川地震后，該區(qū)域的地質(zhì)活動較為頻繁，生態(tài)系統(tǒng)也在逐漸恢復，氣溫和降水作為植被生長的兩個重要因素，促進該區(qū)域生態(tài)重建[31]?，敹嗫h位于高海拔地區(qū)，氣溫低、降水量少、人類活動影響程度小，氣溫和降水成為植被覆蓋變化的主要誘因[32]。分區(qū)五占總面積的43.8%，主要分布在若爾蓋縣和紅原縣等東北部高寒高緯度地區(qū)，土地沙化嚴重，植被變化的主導驅(qū)動因素是人類活動，因此，該區(qū)域的年植被覆蓋變化與年均氣溫和年降水量的相關(guān)性不大[33]。

在研究區(qū)不同的區(qū)域，年植被覆蓋與年均氣溫的偏相關(guān)性和年植被覆蓋與年均降水量的偏相關(guān)性的方向與程度不同，整體而言，植被覆蓋對降水的變化相較于氣溫更加敏感，此外研究區(qū)內(nèi)年植被覆蓋與年降水量呈負相關(guān)的區(qū)域，往往與年均氣溫呈正相關(guān)，這與仙巍等[20]研究川西高山高原過渡帶植被對氣候變化響應的結(jié)果相符。

4 結(jié) 論

(1) 運用ESTARFM模型融合的2005年和2015年的GIMMS NDVI3g數(shù)據(jù)與對應的真實MODIS NDVI數(shù)據(jù)吻合度較高，能凸顯地區(qū)特征，融合結(jié)果可用于地形復雜的川西高原植被變化的研究，為早期的大尺度植被變化研究提供了一種新的思路。

(2) 1995年、2005年和2015年川西高原植被覆蓋結(jié)構(gòu)整體較好，空間上呈現(xiàn)沿海拔先升高后降低的垂直分布特點，時間上植被變化趨勢穩(wěn)定。植被覆蓋度較低的區(qū)域主要分布在低海拔區(qū)，近20 a來研究區(qū)植被覆蓋度整體變化趨勢比較穩(wěn)定，其中以高覆蓋度區(qū)域面積的逐年減少和較低覆蓋度區(qū)域面積的逐年增加趨勢較為明顯。

(3) 1995年、2005年和2015年川西高原年均氣溫和年降水量存在空間異質(zhì)性，氣溫空間差異表現(xiàn)為研究區(qū)北部氣溫低于南部、西部氣溫低于東部，降水空間差異表現(xiàn)為研究區(qū)西部降水量低于東部。

(4) 年植被覆蓋與年均氣溫呈正相關(guān)的區(qū)域面積占14.9%，年植被覆蓋與年均降水量呈正相關(guān)的區(qū)域面積占35.4%，此外研究區(qū)內(nèi)年植被覆蓋與年降水量呈負相關(guān)的區(qū)域，往往年植被覆蓋與年均氣溫呈正相關(guān)，植被覆蓋度對降水的變化相較于氣溫的變化更加敏感。

本研究中使用的數(shù)據(jù)時間間隔較長，可能會忽略一些細節(jié)。此外，由于數(shù)據(jù)獲取等原因，僅選取了年均氣溫和年降水量兩個主要的氣候要素與植被進行響應分析，未考慮地表因素和人為因素等其他因素對植被的響應。未來的研究將會考慮更短的時間間隔和更多因素對植被變化的影響研究。