覃巧婷,陳建軍,2*,楊艷萍,趙曉宇,周國清,2,尤號田,2,韓小文,2 (.桂林理工大學測繪地理信息學院,廣西 桂林54004;2.桂林理工大學廣西空間信息與測繪重點實驗室,廣西 桂林 54004)
黃河源區位于青藏高原腹地,是黃河流域的主要產流區、水源涵養區,也是我國重要的生態屏障[1-2].黃河源區植被以高寒草甸、高寒草原和高寒灌叢為主,其草地生態系統是長期自然演化的結果,生態環境十分脆弱,容易因外界的干擾和破壞而發生變化,恢復難度極大且恢復過程緩慢[3].近幾十年來,在氣候變化和人類活動的影響下,黃河源區的生態環境發生了顯著的變化,主要表現為高寒草地嚴重退化[4]、土地沙漠化[5]、生物多樣性和數量銳減[6]、水土流失加劇[7]等,引起了眾多學者和決策部門的廣泛關注[8-10].作為青藏高原陸地生態系統的重要組成部分,黃河源區植被對維護黃河流域生態環境安全及可持續發展具有非常重要的作用[11].因此,研究黃河源區植被的變化及其驅動機制具有重要意義.
利用遙感衛星數據監測植被的生長有助于更好地了解植被信息.當前用于度量植被變化的植被指數有很多,例如歸一化植被指數(NDVI)、增強型植被指數(EVI)、土壤調節植被指數(SAVI)、差值植被指數(DVI)等,但NDVI是應用最廣泛的植被指數[12].植被變化監測的方法主要有兩種,一種是站點長期監測,一種是遙感監測.站點監測精度較高,但空間代表性有限;遙感監測能夠較精確地反映地面上植被的變化情況,是當前研究區域植被變化的主要技術手段[13].近年來,許多學者利用 NDVI時間序列從月尺度[14]、年尺度[15]、像元尺度和區域尺度[16]等不同方面對植被變化進行動態監測,認為全球氣候變化導致植被發生了顯著的變化.影響植被變化的氣候因子主要有氣溫、降水和輻射等,而這些氣候變化主要是通過改變植被生長所需要的能量和水分,影響碳的積累過程、水循環過程和土壤有機碳的分解轉化過程,從而影響植被的生長和分布格局[17].目前關于植被與氣候變化的研究相對較多,大多數學者采用氣溫、降水來研究植被對氣候變化的敏感性[18-20],指出與人類活動、地形因素等相比,氣候變化對植被的生長起決定性作用.以往的研究大多討論了整體植被與氣候變量之間的相關性,但氣溫、降水對植被的影響存在空間異質性,且不同植被類型對氣候變化的響應特征存在差異[21].同時,相關研究表明植被在垂直方向上的變化是不同的,高海拔地區氣候變化對植被生態系統的影響更大[22].黃河源區位于青藏高原東北部,其平均海拔為4000m左右,地形起伏較大,地理條件復雜[23].而目前關于不同地形條件下氣候驅動影響的植被分布及演化的研究較少,因此,綜合考慮地形因子和氣候因子對不同植被類型的影響,對于更好地了解黃河源區植被動態具有重要意義.
本文基于 2000~2019年的時序 MODIS NDVI數據,采用一元線性回歸及趨勢分析等方法,分析黃河源區植被的時空變化特征;同時,結合地形數據和氣候數據,利用相關分析和統計分析等方法,探討植被對地形和氣候因子的響應機制,以期為黃河源區生態環境保護提供理論依據.
黃河源區位于青藏高原東北部,東經 95.5°~103.5°,北緯 31.5°~36.5°,涉及青海、四川、甘肅 3 省的6個州、18個縣,面積約為12.2萬km2[24-25],占黃河流域總面積的 16%[26].黃河源區海拔約在 3000m以上,地勢西高東低,地貌復雜,分布有高山、盆地、峽谷、草原、沙漠和眾多的湖泊、沼澤、冰川及多年凍土等.該區屬于典型的高原大陸性高寒氣候,冷熱兩季交替,干濕季節分明,日照時間長,輻射強烈[27].黃河源區多年平均氣溫為-4.0~5℃,自西北向東南逐漸降低,海拔越高氣溫越低;多年平均降水量為220~780mm,降水年內分配不均,年際變化大,表現為冬干秋旱、夏秋降水集中[28].區內高寒植被分布廣泛,主要有高寒草甸、高寒草原、高寒灌叢和高寒沼澤等,其中高寒草甸和高寒草原約占總面積的 80%.黃河源區主要以畜牧業為主,產業結構單一,生產力不發達[29].
1.2.1 NDVI數據 本研究使用的 MODIS NDVI數據是中分辨率成像光譜儀植被指數產品MOD13Q1,下載自美國國家航空航天局(NASA)(https://lads web.modaps.eosdis.nasa.gov/search/),空間分辨率為 250m,時間分辨率為 16d.采用 MODIS重投影工具 (MRT)對原始數據進行拼接、重采樣、格式轉換等處理,并利用黃河源區矢量邊界對處理后的數據進行裁剪.為了消除云層、大氣和太陽高度角等的影響,將每年 NDVI數據采用最大值合成法進行合成.
1.2.2 氣象數據 采用氣溫和降水作為氣象資料,氣象數據來自歐洲中期天氣預報中心(ECMWF)第三代再分析資料 ERA-Interim (https://apps.ecmwf.int/datasets/data/interim-full-daily/levtype=sfc/),空間分辨率為0.5°.ERA-Interim提供了自1979年以來的再分析資料,并實時更新.相比第二代產品 ERA-40數據同化模型,ERA-Interim采用了四維變分分析(4D-Var),并結合改進的濕度分析,衛星數據誤差校正等技術,實現了再分析資料質量的提升[30].采用克里金插值法對氣象數據進行空間插值,將氣象數據重采樣為空間分辨率為250m的柵格數據,實現氣象數據與NDVI數據空間分辨率的統一.
1.2.3 輔助數據 DEM 數據來源于地理空間數據云(http://www.gscloud.cn/),空間分辨率為 30m,地理坐標為GCS_WGS_1984,對DEM數據進行拼接、裁剪、重投影的處理,獲得黃河源區的高程數據,并提取坡度數據,研究黃河源區植被隨高程、坡度的變化特征.植被類型數據來源于中國科學院資源環境科學數據中心(http://www.resdc.cn/),其植被類型主要包括針葉林、闊葉林、灌叢、荒漠、草原、草甸、沼澤、高山植被、栽培植被等,本文僅選取主要的3種植被類型(高寒草甸、高寒草原、高寒灌叢)進行分析,其空間分布如圖1所示.

圖1 黃河源區主要植被類型空間分布Fig.1 Spatial distribution of main vegetation types in the source region of Yellow River
1.3.1 最大值合成法 采用最大值合成法獲取黃河源區年NDVI的最大值,從而去除殘云、殘霧、大氣等的影響.計算方法如下所示:

式中:i為年序號,取值范圍為2000~2019;NDVIi為第i年的最大 NDVI值;NDVI1、NDVI2、NDVI3、NDVI4分別表示每年7月份上半月和下半月、8月份上半月和下半月的NDVI數據.
1.3.2 NDVI變化趨勢 采用線性回歸趨勢分析方法計算植被 NDVI的變化趨勢.多年回歸方程中的趨勢斜率代表年際變化,通常利用最小二乘法求解.計算公式如下:

式中:θslope為 NDVI的回歸斜率;n為研究時間段累計年數,本文中n=20;i代表年序號;NDVIi為第i年平均 NDVI值.若θslope>0,表明 NDVI在 n 年間呈增加趨勢;θslope<0,表明 NDVI在 n年間為下降趨勢;θslope=0,表明NDVI在n年間無明顯變化[31].
1.3.3 NDVI與氣溫、降水的偏相關分析 偏相關分析是探索植被生長與單一氣候因子之間的聯系,同時剔除其他氣候因子的影響,被廣泛用于檢測影響植被生長的主要氣候驅動因子[32].計算公式如下:

式中:x為 NDVI;y/z為氣溫/降水,單位為℃/mm;rxy,z為將降水/氣溫固定后氣溫/降水與 NDVI的偏相關系數.若 r>0,表明氣溫/降水與 NDVI呈正相關;r<0表明氣溫/降水與 NDVI呈負相關;r=0表明氣溫/降水與NDVI不存在線性相關性.rxy、rxz、ryz分別代表NDVI與氣溫、NDVI與降水、氣溫與降水的相關系數.其中,相關系數的計算公式如下:

式中:x為 NDVI;y為氣溫/降水;rxy為氣溫/降水與NDVI之間的相關系數;xi與yi分別為NDVI、氣溫/降水第 i年的值;ˉx和ˉy分別表示 NDVI、氣溫/降水兩變量n年的平均值;n為樣本數,本文中n=20.
1.3.4 NDVI與氣溫、降水的復相關分析 復相關分析是研究一個變量與多個變量之間的相關程度.植被是受到多種因素影響的,因此運用復相關分析計算出NDVI與氣溫、降水的復相關系數,以此分析NDVI受兩者共同作用的影響有多大.復相關系數計算公式如下所示:

式中:x為NDVI;y為氣溫;z為降水;rx,yz為NDVI與氣溫、降水的復相關系數.
復相關系數的顯著性檢驗采用 F檢驗法,其統計量計算公式如下所示.

式中:n為樣本數,本文中n=20;k為自變量個數.
1.3.5 NDVI與高程、坡度的關系 為了探究NDVI隨海拔梯度的分布和變化特征,將黃河源區分為 37個 100m(2564~6295)的高程帶,統計每一高程帶的平均 NDVI值和所占的像元比例,其中整個研究區共有 156056034個像元.同理,將黃河源區分為30 個 2°(0~60°)和 1 個 30°(60~90°)的坡度帶,統計不同坡度帶所占的像元比例和 NDVI值,分析 NDVI隨坡度的變化特征.
2.1.1 植被NDVI的時間變化特征 如圖 2所示,2000~2019年黃河源區植被NDVI呈波浪式的增加趨勢,線性增長率為 0.0014/a.在 2003年植被 NDVI出現最小值,在2018年植被NDVI出現最大值.黃河源區各植被類型NDVI的年際變化與植被整體基本一致.其中高寒草原的增加速率最大,為 0.0027/a;高寒草甸次之,為 0.0013/a;高寒灌叢的增加速率最小,為0.0008/a.表明近20a來研究區高寒草原的改善趨勢更為明顯.在數值上,3種植被類型NDVI值從大到小依次為高寒灌叢、高寒草甸、高寒草原.其中高寒草甸與植被整體的變化趨勢最為接近.

圖2 黃河源區植被NDVI的年際變化Fig.2 Annual variation of vegetation NDVI in the source region of the Yellow River
2.1.2 植被 NDVI的空間分布格局 將黃河源區植被 NDVI劃分為<0.2、0.2~0.4、0.4~0.6、0.6~0.8、>0.8 5個等級.如圖3所示,2000~2019年黃河源區平均NDVI為0.65,表明多年植被NDVI處于中高水平,但不同地區差異較大,呈現出由東南向西北遞減的變化趨勢.其中中高植被覆蓋區主要分布在東南部的紅原縣、若爾蓋縣、瑪曲縣、久治縣等區域;低植被覆蓋區主要分布在西北部的瑪多縣、曲麻萊縣的部分區域.

圖3 黃河源區植被NDVI的空間分布Fig.3 Spatial distribution of vegetation NDVI in the source region of the Yellow River
2.1.3 植被NDVI的變化趨勢 對黃河源區2000~2019年植被NDVI做一元線性回歸分析,根據回歸顯著性檢驗結果(P<0.01),將黃河源區植被 NDVI變化情況劃分為顯著減少、無顯著變化和顯著增加3類.
如圖4所示,總體上植被NDVI顯著增加的面積大于植被 NDVI顯著減少的面積.分析其原因,一方面是近幾十年來實施的生態保護與建設工程取得的成效,例如圍封、退牧還草、濕地保護等;另一方面是降水增多、氣溫升高等使得植被生長的氣候環境變好.其中植被 NDVI顯著增加的面積占源區總面積的 20.42%,主要呈塊狀分布在黃河源區北部的興海縣、同德縣和澤庫縣的部分區域以及扎陵湖和鄂陵湖周邊地區;植被 NDVI顯著減少的面積占研究區總面積的 3.14%,主要零星分布在東南部曲麻萊縣和南部達日縣和久治縣等部分區域;植被NDVI基本不變的面積占研究區總面積的 76.44%,分布于黃河源區大部分區域.
2.2.1 NDVI隨海拔的變化特征 黃河源區地勢西高東低,最高海拔達到6295m,高海拔地區主要分布在黃河源區西北部,低海拔地區主要分布在東南部.如圖 5所示,黃河源區的高程主要分布在3400~4700m之間,其所占比例達到整個研究區的92.26%.植被 NDVI隨著高程帶的不同而存在一定的差異,在高程小于3500m范圍內,植被NDVI隨著高程的增加而增加;高程在 3500~4100m 之間,植被 NDVI處于一個較平緩的趨勢,之后隨著高程的增加迅速下降;高程在 5300m 以上時,植被 NDVI極小,接近于0.

圖5 黃河源區植被NDVI隨海拔的變化關系Fig.5 The relationship of vegetation NDVI with altitude in the source region of the Yellow River
2.2.2 NDVI隨坡度的變化特征 由圖6可知,研究區內坡度小于20°的區域所占比例較大,占研究區總面積的 72.42%,主要分布在黃河源區的西北部和東南部.植被 NDVI隨著坡度的增加呈現出單峰型曲線,在 24~26°坡度帶達到峰值.當坡度小于 24°時,植被 NDVI隨著坡度的增加而逐漸增大;當坡度大于26°時,植被NDVI隨著坡度的增大而逐漸減小.

圖6 黃河源區植被NDVI隨坡度的變化關系Fig.6 The relationship of vegetation NDVI with slope in the source region of the Yellow River
2.3.1 植被NDVI與氣溫、降水的時間相關性 植被NDVI與氣溫、降水時間序列的相關性采用皮爾遜系數來表示,如表 1所示,黃河源區植被 NDVI與氣溫、降水的相關性均顯著,相關系數分別為0.5109和0.5048,均通過了0.05的顯著性檢驗,表明氣溫、降水對植被的生長起促進作用.各植被類型與氣溫、降水均表現出正相關關系,其中氣溫與高寒草甸的相關性最好,且通過了0.05的顯著性檢驗;降水與高寒灌叢的相關性最好,高寒草甸次之;高寒草原與氣溫、降水的相關性相對較差,其與降水的相關系數通過了 0.05的顯著性檢驗,與氣溫的相關系數不具有統計上的顯著性.

表1 黃河源區植被NDVI與氣溫、降水的相關系數Table 1 The correlation coefficient between vegetation NDVI and temperature and precipitation in the source region of the Yellow River
2.3.2 植被NDVI與氣溫、降水的空間相關性 由于不同地區特殊的自然條件,植被對氣溫、降水的響應程度各不相同.由圖 7可知,黃河源區植被 NDVI與氣溫、降水在空間上呈現出以正相關為主的趨勢.NDVI與氣溫、降水的最大偏相關系數分別為0.922和0.907.植被NDVI與氣溫呈現正相關的區域主要分布在黃河源區西北部,表明在這些區域氣溫是影響植被變化的重要因素;與氣溫呈現負相關的區域主要分布在黃河源區東南部.植被 NDVI與降水呈現正相關的區域主要分布在扎陵湖、鄂陵湖西北部及龍羊峽水庫周邊,與降水呈現負相關的區域主要分布在研究區中部和東南部部分地區.

圖7 黃河源區植被NDVI與氣溫和降水的偏相關系數Fig.7 Partial correlation coefficient between vegetation NDVI and air temperature and precipitation in the source region of the Yellow River
除了單一的氣候因子對植被的影響之外,植被的變化還受到多種因素共同作用的影響.利用復相關公式對氣溫、降水進行計算并用F檢驗法進行顯著性檢驗,分析氣溫、降水共同作用下植被的變化.由圖8可知,黃河源區植被NDVI與氣溫、降水的復相關系數介于0~0.93之間.整體上,植被NDVI與氣候因子復相關性較強的區域主要集中在黃河源區西北部的曲麻萊縣、瑪多縣、達日縣;復相關性較弱的區域分布在龍羊峽水庫周圍及若爾蓋縣等海拔較低的區域.

圖8 黃河源區植被NDVI與氣溫、降水的復相關系數Fig.8 The multiple correlation coefficient between vegetation NDVI and temperature and precipitation in the source region of the Yellow River
本文基于MODIS NDVI數據,采用最大值合成法、趨勢分析法研究黃河源區植被的時空變化特征.研究結果表明,2000~2019年黃河源區植被NDVI呈增長趨勢,與劉啟興[24]、Yang[33]等的研究結論相一致,表明近年來黃河源區生態環境在不斷改善.植被的變化是氣候因子和人類活動共同作用的結果[34],對于生態系統較為脆弱的敏感區域,氣候的突變或人類活動的加劇,都會使植被的變化出現較大的波動[35].黃河源區地處高海拔地區,是我國生態環境脆弱的地區之一.近年來由于人類活動加劇和全球氣候變暖,區域內植被變化在不同年份間出現了較大的波動性,整體上植被呈增加趨勢.但植被的改善主要依賴于局部地區的增長,而不是整個區域.其中植被 NDVI增加的區域主要分布在研究區北部,例如澤庫縣、興海縣和瑪多縣等;植被NDVI減少的區域主要分布在研究區南部和東南部,例如紅原縣、若爾蓋縣、達日縣和甘德縣等.除了氣候因素,植被變化與人類活動也有著不可分割的聯系.黃河源區在過去很長一段時間一直處于嚴重超載放牧的狀態,過度放牧是引起該區草地生態系統退化的主要原因
[36].過度放牧會破壞黃河源區脆弱的草地生態系統,降低其保水性能,從而導致草場退化,地表裸露,出現沙礫地.近年來黃河源區植被的好轉可能與多種生態保護與建設工程有關,例如三江源自然保護區、退耕還草、生態移民等.
氣候變化對草地生態系統的影響主要是通過降水和溫度的變化,改變土壤濕度和溫度,同時影響光合作用和植物呼吸作用,進而調節草地的生長和生態系統生產力[37-38].在全球變暖的大背景下,對植被變化與氣候因子的研究呈現多樣性,尤其是在年尺度上.通過對黃河源區過去20a氣溫和降水的變化曲線進行分析(圖 9),發現近 20a間研究區氣溫和降水整體上均呈增加趨勢,線性增長率分別為 0.0344℃/a和 0.4949mm/a,其中氣溫增長速率與馮曉莉[39]、劉光生[40]等的研究結果一致,降水增長率低于兩位作者的研究結果,可能是研究區域范圍以及研究時間段不一致導致的結果.植被 NDVI隨氣溫和降水上下波動,與降水相比,NDVI與氣溫的變化曲線走勢較為一致,表明植被對氣溫的響應更強烈.

圖9 黃河源區植被NDVI與氣溫和降水的變化曲線Fig.9 Variation curve of NDVI, temperature and precipitation in the source region of the Yellow River
經過計算可知,植被NDVI與氣溫、降水均呈線性正相關關系,與氣溫、降水的相關系數均通過了0.05的顯著性檢驗(表 1).不同植被類型對氣候變化的響應存在差異,其中高寒草甸受到氣溫的影響較大,高寒灌叢對降水的敏感性強于氣溫.高寒草甸主要分布在較高海拔地區,氣溫較低,降水量較大,植被的生長主要受氣溫的影響,所以高寒草甸對氣溫的敏感性更強.高寒灌叢主要分布在較低海拔地區,氣溫高,蒸發量大,可能由于水分條件的限制,植被的生長受到影響,所以其生長更依賴于降水[41].結合相關分析,發現近20a來黃河源區植被NDVI與氣溫、降水均呈現出以正相關為主的趨勢,但與降水相比,氣溫對黃河源區植被變化的影響更為顯著,這與相關研究結果一致[42-43].在一定程度上,溫度的升高可以增強植物的光合強度,延長植物的生長期[44-45],從而促進植被的生長.因此,黃河源區植被覆蓋的提高可能是氣候變暖的結果.
相關研究表明,海拔高度對植被覆蓋有一定的影響,植被覆蓋與溫度的關系隨著海拔的變化具有不同的特征[46-47].在本文的研究結果中,植被 NDVI隨著海拔高度的增加呈現出單峰型曲線,峰值出現在 2700~2800m 和 3500~3600m.植被 NDVI較大的區域主要集中在3100~4800m的中高海拔地區.在高程<3600m 的范圍內,植被 NDVI與高程呈正相關,NDVI隨著海拔高度的增加而增加.海拔在3600~5300m的范圍內,植被NDVI隨著高程的增加呈現出減小的趨勢;當海拔>5300m時,植被NDVI趨于穩定,且數值較低,這可能反映了氣候變化條件下高寒植被對環境的響應.分析其原因,可能是在低海拔地區,地勢相對平緩,水熱條件優越,人類活動較為頻繁,嚴重影響了植被的恢復和生長,使植被覆蓋降低.而高海拔地區人口稀少,人類活動相對較少,從而黃河源區植被呈現出增長的趨勢,并達到一個峰值.而后,隨著海拔的進一步升高,地勢陡峭,降水量增加,氣溫下降,光照強烈,風速大,蒸散迅速,植被的生長受到限制,因此植被呈減少的趨勢.
坡度是影響地表物質流動的重要因子,直接影響物質交流與能量轉換的方式和程度,改變土壤的基本屬性,從而影響植被的分布態勢[48].本文的研究結果表明,隨著坡度的增大,植被 NDVI呈現出先增大后減小的變化趨勢,且在 24~26°坡度帶 NDVI達到最大值(圖6).相對于坡度較大的區域,坡度平緩更有利于植被的生長,但黃河源區大部分坡度相對較緩的地區受到的人為干擾相對較多,所以水熱條件較好的緩坡地區其植被覆蓋反而較低.坡度影響地表徑流和排水狀況,所以坡度越高的區域自然降水就越容易流失,土壤礦質營養等稀薄,植被生長差,植被覆蓋也較低.
植被的動態演變受到氣候因素、地形因素、土壤因素、人類活動等多個方面因素的共同影響,作用機制復雜.本文由于缺少相關資料,只考慮了地形因子和氣候因子對黃河源區植被變化的影響,尚未考慮到其他的因子,若能把影響植被變化的其他因子共同進行研究,則更能深刻地揭示各種因素對植被生長的作用機理.
此外,目前NDVI時間序列數據主要有NOAA/AVHRR NDVI、MODIS NDVI、SPOT/VGT NDVI等,而研究植被變化及其影響機制的方法也有許多種,例如趨勢分析法、相關分析法、殘差分析法、Hust指數法、統計分析法等.本研究使用的數據為MODIS NDVI,采用的方法為趨勢分析法、偏相關分析法、復相關分析法,研究結果因使用的數據、方法不同可能存在一些差異.
4.1 黃河源區植被 NDVI整體呈增加趨勢,線性增長率為0.0014/a,波動性較大.多年平均NDVI為0.65,但空間差異顯著,呈現出由東南向西北遞減的空間分布格局.NDVI顯著增加的區域主要呈塊狀分布在黃河源區北部及扎陵湖、鄂陵湖周邊地區,NDVI顯著減少的區域主要零星分布在研究區東南部及南部的部分區域.
4.2 高海拔地區主要分布在黃河源區西北部,低海拔地區主要分布在東南部.隨著高程的增加,植被NDVI呈現出先增加后減少的趨勢,但在 3500~4100m區間 NDVI變化不顯著.隨著坡度的增大,植被 NDVI呈現出單峰型曲線,并在 24~26°坡度帶達到最大值.
4.3 2000~2019年氣溫和降水整體上均呈上升趨勢.NDVI與氣溫、降水的相關性均通過了0.05的顯著性檢驗,表明氣溫、降水對植被的生長起促進作用.與降水相比,NDVI對氣溫的敏感性更強.