張 江,袁旻舒,張 婧,李函微,王潔儀,張 賢,鞠佩君,蔣海波,陳 槐,朱求安,3,4,*
1 西北農林科技大學林學院,生態預測與全球變化研究中心,楊凌 712100 2 中國科學院成都生物研究所山地生態恢復與生物資源利用重點實驗室,成都 610041 3 河海大學水文水資源學院,南京 210098 4 國家地球系統科學數據中心,國家科技基礎條件平臺,北京 100101
氣候變化和人類活動對陸地生態系統的影響一直是全球氣候變化研究中的熱點[1-2]。草地生態系統約占陸地總面積的1/5,是陸地生態系統重要組成部分,在調節氣候、水土保持、防風固沙等方面發揮著重要作用[3-4]。有研究表明,高海拔、高緯度地區草地植被動態對氣候變化及人類活動的反應較為敏感[5-6],主要表現為植被生長狀態及分布范圍的變化[7]。植被的分布由于受到氣候模式的制約,在氣候暖化過程中植被可能會向兩極高緯度地區及高海拔地區移動,進而改變物種的分布范圍[8]。與此同時,人口遷移、畜牧業發展及土地利用變化等人類活動的增強也會影響植被的生長狀態及分布范圍[9]。
青藏高原是全球海拔最高的獨特地域單元[10],平均海拔超過4000 m,素有“世界第三極”之稱,屬于氣候變化敏感區和生態脆弱帶[11]。草地是青藏高原主要的植被類型,總面積約為1.59×106km2。因此,青藏高原被認為是研究草地植被動態對氣候變化及人類活動響應機制的理想場所[12]。一些研究表明,在過去的三十多年,隨著氣候變化和人類活動的加劇,青藏高原草地生態系統的結構和功能發生了巨大變化,草地呈現出不同程度的退化且日趨惡化,已嚴重阻礙了當地社會經濟的可持續發展[13-14]。也有一些研究表明,青藏高原草地退化趨勢已有所緩解,呈現出“整體改善、局部退化”的趨勢[7, 15-18]。在青藏高原草地植被變化的研究中,氣溫和降水被認為是導致草地植被變化的主要氣候因子,而放牧及人口數量變化引起的其他人類活動(如土地利用變化、植樹造林等)則被認為是導致草地植被變化的主要人為因子[19-20]。
歸一化植被指數(Normalized difference vegetation index, NDVI)因能有效地反映植被覆蓋程度、生長狀況、生物量及凈初級生產力等植被參數,被廣泛應用于大尺度植被動態研究中[21-22]。例如卓嘎等[23]利用NDVI數據分析了西藏地區植被的分布狀況及變化趨勢;董璐等[24]利用NDVI數據研究了新疆地區植被生長對氣溫的響應;李鵬[25]利用NDVI數據研究了青藏高原植被物候對極端氣候的響應。目前在區域尺度上,關于氣候變化對青藏高原植被動態影響的研究較多,例如楊元合[12]、Gao等[26]、劉軍會等[16]利用NDVI數據研究分析了青藏高原草地植被動態與氣候因子氣溫、降水之間的關系。而關于人類活動對青藏高原草地植被動態影響的研究多為樣點實驗,例如周華坤等[27]、宋磊等[28]在青藏高原高寒草原上進行放牧試驗來研究放牧強度對植被群落的高度、蓋度、生物量及多樣性的影響,在區域尺度上的研究相對較少。因此在區域尺度上綜合分析氣候變化、人類活動對青藏高原高寒草地植被變化影響的研究亟待開展。
本研究基于32年(1982—2013年)青藏高原143個縣區生長季草地的NDVI、平均氣溫、降水、夏季牧場放牧強度及人口數量等數據,采用趨勢分析法和面板數據模型(Panel data model)在區域尺度上研究了青藏高原高寒草地NDVI動態變化及其對氣候變化和人類活動的響應,以期為青藏高原地區生態環境保護和區域經濟發展提供理論參考依據。
青藏高原位于中國西南部邊緣(75—105°E, 27.5—37.5°N),總面積約為2.57×106km2,境內有河流、冰川、湖泊等多種地貌,地形復雜多變,亦有“世界屋脊”之稱,區域氣候易受東亞季風、印度季風、西南季風及西風環流系統的影響,水熱分布差異明顯[10-11]。本研究基于青藏高原境內的143個縣級行政單元,主要包括青海省、西藏自治區及四川省部分地區(圖1)。

圖1 研究區域
1.2.1青藏高原草地空間分布數據
植被空間分布數據來源于歐洲航空局(European Space Agency, ESA(http://www.esa.int/ESA)),時間跨度為1992—2013年,空間分辨率為300 m。將草地植被類型從數據集中分離得到1992—2013年青藏高原草地分布數據。本研究的時間范圍為1982—2013年,由于缺少1982—1991年植被空間分布數據,因此1982—1991年的草地分布數據采用1992年青藏高原的草地分布數據補充完整。將獲得的草地分布數據與各縣區疊加,得到各縣區草地分布,進而統計得到1982—2013年各縣區的草地面積變化序列。
1.2.2NDVI數據
本研究選用美國航空航天局(National Aeronautics and Space Administration, NASA)全球監測與模型研究組(Global Inventory Monitoring and Modelling Studies, GIMMS)提供的1982年1月至2013年12月,每十五天合成的空間分辨率為8 km的NDVI3g(GIMMS NDVI Version3(https://ecocast.arc.nasa.gov/data/pub/gimms/))數據集來描述青藏高原高寒草地的生長狀況及覆蓋程度。該數據集消除了火山爆發、太陽高度角、云覆蓋等因素的影響,已被廣泛應用于植被動態的研究中[29-30]。每月的NDVI采用了最大值合成法(MVC),該方法可以進一步消除云霧、大氣等因素的干擾,提高數據的精度[31]。青藏高原冬季時期草地大部分被積雪覆蓋,因此本研究使用生長季平均NDVI來反映草地生長狀態,以消除積雪覆蓋對NDVI的影響,生長季為每年5月至9月[12]。將獲得的青藏高原生長季平均NDVI數據與青藏高原各縣區草地分布圖進行疊加,得到各縣區草地NDVI的變化序列。
1.2.3氣候數據
為分析氣候變化對青藏高原高寒草地NDVI變化的影響,研究選取了對植被生長影響最為密切的生長季平均溫度和生長季降水作為影響草地NDVI變化主要的氣候因子。基于全國標準氣象站的數據(中國氣象局(http://www.cma.gov.cn/)),利用薄板樣條插值法[32]對氣溫和降水進行插值,生成了全國氣溫和降水數據,分辨率為10 km。該數據集在其他研究中得到了驗證與應用[33-34]。從上述全國氣象數據中提取青藏高原生長季平均氣溫和生長季降水數據,并將其分別與青藏高原各縣區草地分布數據進行疊加,統計得到1982—2013年各縣區草地分布區的生長季平均氣溫與降水量的變化序列。
1.2.4人類活動數據
為了分析人類活動對青藏高原高寒草地NDVI變化的影響,研究選取了影響高寒草地植被變化人類活動因子中最重要的兩個指標,即反映畜牧業發展規模、直接影響草地生長狀態的放牧強度以及代表城鎮規模、直接影響放牧強度變化的人口數量。1982—2013年各縣區畜牧種類及數量和人口數據來源于各省、市、縣(志)統計年鑒(國家數字電子圖書館(http://www.nlc.cn/))。
在計算放牧強度時,為消除不同畜牧種類及數量之間的差異,均以綿羊為基準單位進行羊單位轉換,轉換比值參考《草業科學概論》[35]。為了與生長季草地NDVI、氣溫、降水數據相對應,研究計算了各縣區夏季牧場的放牧強度,其計算公式如下:
(1)
式中,GIi、Ui和Ai分別為第i年該縣區的放牧強度(羊單位/hm2)、羊單位(只/個)和草地面積(hm2);R和S分別為可利用草地面積占比和夏季牧場占比。其中各縣區可利用草地面積、夏季牧場占比數據均來源于各省、市、縣(志)統計年鑒。
1.3.1放牧強度分區
為研究不同放牧強度對草地NDVI的影響,本研究將各縣區的放牧強度擬劃分為輕、中、重度三個放牧等級區。由于青藏高原各縣區草地載畜能力不同,放牧強度的閾值劃分(即單位面積上牲畜的數量)也存在差異。通過文獻收集,對各研究點放牧強度(輕、中、重度)的定義閾值進行了統計(表1),研究點的分布如圖2。將計算的各縣放牧強度的均值與表1中放牧強度的閾值進行對比,對各縣區進行等級劃分。若同一個縣區存在多個放牧實驗,放牧強度的閾值取其均值,區域內缺少放牧實驗時,則參考鄰近縣區。

表1 放牧實驗中放牧強度閾值和數據來源

圖2 青藏高原輕、中、重度放牧分區示意圖
1.3.2趨勢分析
草地NDVI、氣候因子及人為因子在青藏高原及各縣區的變化趨勢采用趨勢分析法。計算公式如下:
(2)
式中,k為斜率;Yi為第i年各變量值;i為年變量;n為年數。
在空間尺度中,根據斜率k及P值大小,將各變量變化趨勢劃分為-3、-2、-1、+1、+2和+3六個等級,分別代表極顯著下降(k<0,P≤0.01)、顯著下降(k<0, 0.01
0.05)、上升不顯著(k>0,P>0.05)、顯著上升(k>0, 0.01
0,P≤0.01)。
1.3.3面板數據模型
面板數據Panel Data即空間和時間同時擴展的混合數據集合,是指由變量Y關于N個不同對象在T個觀測時期所得到二維結構數據。其優勢在于,面板數據能同時反映某一時期各個個體數據的規律和每個個體隨時間變化的規律。與截面數據模型相比,面板數據模型可以控制不可觀測變量所導致的最小二乘法(OLS)估計的偏差,使得模型設定更加合理、模型參數的樣本估計量更加準確;與時間序列模型相比,該模型擴大了樣本信息、降低了變量間的共線性,提高了估計量的有效性[72]。面板分析選取了直接影響草地NDVI變化的氣候因子氣溫、降水及人類活動因子放牧強度。人口數量作為間接影響草地NDVI的因子也被引入模型(主要通過影響放牧數量、退耕還草等方式間接影響草地NDVI)。面板數據模型設定如下:
NDVIit=αi+βi1Tempit+βi2Precit+βi3GIit+βi4Popuit+μit
(3)
式中,NDVIit、Tempit、Precit、GIit和Popuit分別表示第i個縣區在第t年生長季的草地NDVI、生長季平均氣溫(℃)、生長季降水(mm)、夏季牧場的放牧強度(羊單位/hm2)及人口(百萬);αi、βi1、βi2、βi3和βi4分別表示第i個縣區在第t年的截距項及4個自變量的系數;i為縣區個數;t為不同年份;μit為隨機擾動項。
面板數據模型分析中,采用z-score標準化法對各變量進行標準化處理,求得各變量的相對貢獻;采用未標準化的數據分析各變量單位變化量對NDVI的影響。為避免面板模型出現偽回歸現象,采用ADF單位根檢驗法[73]對面板數據進行單位根檢驗。面板數據有混合效應、固定效應和隨機效應3種模型[74],模型的確立采用F檢驗和Hausman檢驗[75-76]。數據分析表明,ADF單位根檢驗結果拒絕其零假設(P<0.01)表明數據在時間序列上整體呈現平穩性不存在單位根,F檢驗和Hausman檢驗結果均拒絕零假設(P<0.01),因此本研究面板分析在時間和截面上設定為固定效應模型。
趨勢分析結果顯示,1982—2013年青藏高原生長季草地NDVI的變化范圍在0.431—0.471之間,草地NDVI整體上呈增長趨勢,但未達顯著水平(P>0.05)。青藏高原大部分縣區生長季草地NDVI呈增長趨勢,少部分縣區呈下降趨勢,各縣區的變化趨勢無明顯空間異質性(圖3)。
32年來,青藏高原生長季平均氣溫呈顯著上升趨勢,增長率約為每10年0.420℃。各縣區生長季平均氣溫增長趨勢顯著,無明顯的空間異質性。青藏高原生長季降水量整體上呈增長趨勢,最大、最小降水量分別為481.179 mm和381.635 mm。各縣區生長季降水量的變化趨勢呈明顯的空間異質性,青藏高原西部和北部降水量顯著增加,而東部和南部降水量增長趨勢不明顯且部分縣區呈下降趨勢(圖3)。
1982—2013年青藏高原放牧數量逐年上升,夏季牧場放牧強度增長趨勢顯著,每年每公頃增長約0.014個羊單位。各縣區夏季牧場放牧強度的變化趨勢存在一定的空間異質性,其中青藏高原中部和北部地區放牧強度呈下降趨勢,而東南部地區放牧強度呈顯著增長趨勢。32年來青藏高原人口增長趨勢顯著,總人口由1982年的817萬增長至2013年的1200萬,每年大約增長13萬人。除部分縣區人口呈現下降趨勢外,其他各縣區人口增長趨勢顯著,且無明顯空間異質性(圖3)。

圖3 1982—2013年青藏高原生長季草地NDVI、氣溫、降水、放牧強度及人口的年際變化與趨勢變化圖
分別對不同的放牧強度分區及青藏高原進行了面板模型分析,4個分析情形下模型結果均具有較高的擬合優度,青藏高原全區以及輕度、中度和重度放牧區的模型決定系數分別為0.56、0.77、0.48和0.41(表2)。在輕度放牧區,生長季草地NDVI與氣溫、降水及放牧強度之間均呈極顯著正相關,與人口之間相關性不顯著。氣溫、降水及放牧強度每增加1單位,草地NDVI相應增加0.0047、0.0007和0.0421,相對貢獻系數分別為0.0576、0.6751和0.2927。在中度放牧區,生長季草地NDVI與氣溫、降水及放牧強度之間呈極顯著正相關,與人口之間相關性仍不顯著。放牧強度對草地NDVI變化的貢獻系數減小至0.0073,而氣溫的貢獻系數上升至0.0234。在重度放牧區,生長季草地NDVI與氣溫、降水呈極顯著正相關,與放牧強度之間則呈極顯著負相關,相對貢獻系數分別為0.2635、0.6889和-0.1525,同時放牧強度每增加1單位,草地NDVI相應減少0.0016。全區分析結果表明,生長季草地NDVI與氣溫、降水之間呈極顯著正相關,與放牧強度呈負相關,但未達到顯著水平(P>0.05)。放牧強度每增加1個單位,草地NDVI減少0.0001。降水的相對貢獻系數最大,溫度次之,分別為0.6963和0.2810。

表2 面板數據模型估計結果
1982—2013年,生長季平均氣溫在不同放牧區中均呈顯著增長趨勢(P<0.01),生長季降水量在輕度放牧區呈顯著增長趨勢(P<0.01),但在中度、重度及全區中,生長季降水增長趨勢未達到顯著水平(P>0.05)(表3)。32年來,輕度放牧區生長季氣溫與降水量分別增加了約1.632℃和43.488 mm,兩者的增加使草地NDVI分別增加了0.0077和0.0304;放牧強度減少了0.096,導致草地NDVI變化了0.0054。在中度放牧區中,生長季平均氣溫、降水及放牧強度的增加分別導致草地NDVI增加了0.0300、0.0271和0.0014。在重度放牧區中,溫度與降水的增加使草地NDVI分別增加了0.0222和0.0238,但放牧強度的增加導致草地NDVI減少了約0.0015。對于整個青藏高原地區,1982—2013年青藏高原草地生長季平均氣溫、降水量及放牧強度的增長率分別為每年0.042℃、1.056 mm和0.014個羊單位/hm2(圖3和表3),溫度和降水的增加使草地NDVI分別增加了約0.0317和0.0270,放牧強度的增加則導致草地NDVI減少了約0.00005。

表3 1982—2013年氣溫、降水、放牧強度及草地NDVI的變化量
青藏高原作為氣候變化的敏感區對氣候變暖已表現出諸多響應。近幾十年來,青藏高原氣溫顯著上升,降水量略微增加,總體呈現出“暖濕化”趨勢[77-79]。本研究結果與其類似,青藏高原生長季平均氣溫、生長季降水量整體上均呈增長趨勢。東部及南部地區氣溫增長趨勢顯著,降水增長趨勢不明顯,區域氣候呈“暖干化”趨勢;青藏高原西部和北部地區氣溫、降水增長趨勢顯著,區域氣候呈“暖濕化”趨勢(圖3)。大尺度環流系統會影響水汽的分布[80],在過去的幾十年青藏高原地區西風、東亞季風、印度季風的增強可能改變了青藏高原地區的水汽分布[81-82],導致了青藏高原東南部和西北部區域氣候差異的出現。
青藏高原是我國重要的草地畜牧業基地[83],放牧是其重要的生產方式之一[7]。近幾十年來,青藏高原人口顯著增加[84]。大規模人口對物質資源(奶、肉、毛皮制品)的需求量不斷加大使得放牧數量也隨之增加,進而導致放牧強度逐年上升[85]。但青藏高原地域廣闊各縣區人口分布、經濟發展水平、區域政策及放牧管理方式的差異可能導致了放牧強度變化趨勢存在一定的空間異質性。青藏高原東部、南部縣區人口、牲畜數量較多且草地分布較少,多屬于重度放牧區,而西北部人口、牲畜數量較少且草地分布較廣,多屬于輕度放牧區。
1982—2013年,青藏高原生長季草地NDVI整體呈增長趨勢,草地NDVI除在部分縣區呈下降趨勢外,均呈增長趨勢(圖3)。劉軍同等[16]利用NDVI研究青藏高原植被覆蓋變化時發現青藏高原植被覆蓋度呈“整體升高、局部退化”的趨勢,且改善區域面積遠大于退化區域面積。本研究結果與前人研究結果基本一致,青藏高原高寒草地生態系統呈現出“整體改善,局部退化”狀態[16, 18, 86-87]。
植被的生長狀態及分布會受到溫度、降水等氣候因子的影響[88]。氣候條件的改變會使植被的種類、分布發生相應的改變以適應新的氣候特點[7]。一些研究表明,氣溫的升高會引起植物生長期提前、生長季的加速以及生長期延長,進而會導致NDVI的增加[89-90]。此外,溫度和降水可以通過生物地球化學反饋影響土壤養分及水分的有效性,進而影響植被的生長狀況[91-92]。在本研究中,面板分析結果均表明生長季平均氣溫、降水與草地NDVI呈極顯著正相關(表2),且趨勢分析顯示1982—2013年青藏高原生長季平均氣溫、降水量均呈增長趨勢(圖3)。生長季溫度和降水量的增加可能提高了土壤養分和水分的有效性,加速了草地的生長進而導致了草地NDVI的增加[12]。魏彥強等[7]利用NDVI數據在研究青藏高原植被帶變化對氣溫、降水變化的響應中也發現“暖溫化”對植被狀況的改善具有明顯的促進作用。
植被覆蓋變化除受氣候變化影響外,還與人類活動息息相關[93]。草地植被變化極易受放牧、土地利用變化等人類活動的影響。一些研究表明,放牧強度不斷增加會導致植被的多樣性和生物量顯著降低[69]。許岳飛等[40]以青藏高原高山嵩草草甸為對象得出隨著放牧牦牛強度的增加,草甸群落的蓋度、群落物種數、物種多樣性均顯著降低;孫大帥[53]對青藏高原東部高寒草甸植被的調查研究表明隨著放牧強度的增加,植被生物量和蓋度顯著降低;周興民等[41]在海北高寒草甸生態系統站長期放牧試驗表明,重度放牧組中牲畜對植被的過度啃食和踐踏作用導致了優良牧草占比大幅度下降,植被地上地下生物量也隨放牧強度的增大而線性下降,與此同時,放牧踐踏的加強導致了土壤板結、土壤有機質和全氮含量的下降。但也有研究表明,草原植物的凈初級生產力、生物多樣性及生物量并不是在不利用的情況下最高,而是在適牧條件下最高,即所謂的“優化放牧理論”和“中度干擾理論”[94]。植物生物量及多樣性隨著放牧強度的增加會呈現先上升后下降的趨勢[95]。例如澤讓東科等[96]關于牦牛放牧對青藏高原草地影響的研究發現,隨著放牧強度的增加,群落生物量及草地植被生物量會呈現先上升后下降的趨勢,適度放牧可以刺激草地植物的生長和提高草地生物多樣性。1982—2013年青藏高原各縣區夏季牧場放牧強度增長趨勢存在一定空間異質性,這可能導致了放牧強度與生長季草地NDVI在不同的區域中呈現出不同的相關關系。在輕、中度放牧區中草地NDVI與放牧強度呈極顯著正相關,放牧強度的適當增加可能刺激了草地植物的生長、增加了植物的生物量,進而促進了草地NDVI的增加。在重度放牧區中草地NDVI與放牧強度呈極顯著負相關,長期高強度的放牧可能使得草地的生物量下降、多樣性降低,加上放牧踩踏作用的增強可能導致了重度放牧區草地NDVI的減少。結果與“中度干擾理論”及“放牧優化理論”相符,即中度或適度放牧能刺激植物產量、改善植被生長狀態、增加植物的地上生物量。在中度或適度放牧干擾下,草地生態系統中植被物種的種間關系也發生改變,競爭優勢種在生長一定程度上因牲畜啃食與踐踏作用而受到抑制,種間競爭作用減弱,競爭優勢種和放牧忍受種共存將可能導致群落內的植物多樣性與生物量的提高。但隨著放牧強度的持續增加,干擾作用的增強使得耐受性較差的物種無法生存,最終導致了物種多樣性的降低、生物量的減少[53]。
青藏高原植被的變化受到水資源、地形地貌、氣候及人類活動等因素的影響[10, 18]。一些研究認為氣候變化在青藏高原植被變化中起主導作用[97-98]。也有研究認為過度放牧、土地利用變化、濫挖濫采等人類活動才是導致青藏高原植被變化的主導因子[27, 84]。在輕度放牧區中,氣溫、降水與NDVI呈顯著正相關,雖然氣溫與降水每增加1個單位對草地NDVI的貢獻小于放牧強度,但32年來氣溫與降雨增長趨勢顯著,其增加量所導致草地NDVI的增加量均大于放牧強度導致草地NDVI的增加量。在長時間序列上,降水對草地NDVI的促進作用更強,主導了青藏高原生長季草地NDVI的變化。在中度放牧區中,氣溫增長趨勢顯著,可能導致了溫度對草地NDVI促進作用的增強,而降雨增長趨勢不顯著,可能導致了降水對草地NDVI促進作用的減弱。在重度放牧區及整個青藏高原地區,盡管放牧強度的增加會導致草地NDVI的減少,但32年來降雨、溫度增加所導致草地NDVI的增長量均大于放牧強度導致NDVI的減少量,因此青藏高原生長季草地NDVI呈現出上升趨勢。氣候因子氣溫與降水的增加均會促進草地植被增長,人類活動放牧強度的大小對植被生長則具有雙重影響,兩者長期的共同作用導致了青藏高原高寒草地植被的變化。
本研究以青藏高原143個縣區32年的植被、氣候、放牧等數據為基礎,利用趨勢分析、面板數據模型對青藏高原生長季草地NDVI、氣候及人類活動因子綜合分析得到以下結論:
(1)1982—2013年,青藏高原生長季平均氣溫與降水量整體增加,區域氣候呈現“暖濕化”趨勢。
(2)青藏高原生長季草地NDVI整體呈緩慢增長趨勢,草地生長狀態呈現出“整體改善、局部退化”趨勢。
(3)氣候因子主導了青藏高原生長季草地NDVI的變化。氣候變暖、降水量的增加改善了青藏高原草地植被的生長狀態,促進了草地NDVI的增加。放牧強度的不斷增加會導致草地NDVI的減少,但其導致草地NDVI減少量小于氣候變化導致草地NDVI的增加量,因此青藏高原草地NDVI整體上呈現出增長趨勢。