姚 昆, 何 磊, 白 琳, 羅 涵

(1.西昌學院 資源與環(huán)境學院, 四川 西昌 615000; 2.成都信息工程大學 軟件工程學院,成都 610225; 3.四川省信息化應用支撐軟件工程技術研究中心, 成都 610225)

植被作為陸地生態(tài)系統(tǒng)的重要構成要素之一,是聯(lián)結大氣、土壤和水體等的關鍵紐帶,也是物質能量交換的重要場所[1-2]。近年來,隨著全球變化研究的不斷深入,植被覆蓋變化作為其中一項重要的組成部分,已成為自然資源學、生態(tài)學和地理學等領域關注的一大熱點[3]。因此,掌握區(qū)域長時間序列植被覆蓋變化狀況,并對其驅動因素進行揭示十分必要[4]。

植被NDVI作為定量表征植被生長狀態(tài)的一項重要參數(shù),其常被用于分析區(qū)域植被覆蓋變化,它能較全面描述區(qū)域植被的生長狀態(tài)[5-6]。肖建勇等[7]以MODIS NDVI為基礎,對喀斯特關鍵地帶植被覆蓋時空變化規(guī)律展開分析,并對其驅動因素進行探索;陳超男等[8]基于多源NDVI實現(xiàn)了秦巴山區(qū)1982—2017年植被覆蓋時空變化及與氣候響應規(guī)律的分析,并針對重點地區(qū)植被覆蓋與地形、氣候因子的響應關系展開探討;何亮等[9]基于NDVI結合斜率、顯著性檢驗和Hurst指數(shù)等模型,發(fā)現(xiàn)近20年杭錦旗植被NDVI改善區(qū)面積呈先增加后減少的變化狀態(tài),降水對其變化起正相關驅動作用而氣溫則相反;鄭勇等[10]利用NDVI對川西北高原近20年植被覆蓋變化及與氣候因子的相關性進行探討,結果表明降水對其變化的驅動作用比氣溫更顯著;伍良旭等[11]認為川西北高原1995—2015年植被覆蓋在海拔梯度上呈現(xiàn)顯著的分布差異,相比年均氣溫對其變化的影響,年降雨量的驅動作用更明顯;韓繼沖等[12]對長江上游植被覆蓋變化及與氣候因子的響應關系展開探討,發(fā)現(xiàn)植被NDVI對氣候的響應模式在動態(tài)和靜態(tài)視角下存在顯著差異。以上成果均對植被覆蓋動態(tài)變化規(guī)律及驅動力展開探討,在區(qū)域植被恢復與保護中起到了積極促進效果;然而,它們在進行驅動力分析時,均以相關性為主要研究方法,又側重僅從氣候單系統(tǒng)展開探討;然而,驅動植被覆蓋變化的是因素相對復雜的,各因素的影響作用也存在較顯著差異,若僅采用相關性分析僅從單系統(tǒng)對此展開探討,那么研究的準確性也會存在一定片面性。

地理探測器不僅能實現(xiàn)植被覆蓋空間分異性探索,還能對其背后復雜的驅動力展開系統(tǒng)分析。Peng等[13]以地理探測為工具完成自然因子對四川省植被覆蓋變化的驅動力解析;馬曉妮等[14]基于地理探測器對影響砒砂巖區(qū)植被覆蓋的因子展開探討。

川西北高原區(qū)生態(tài)環(huán)境先天脆弱,加之受過度放牧、草原鼠害和坡耕地開墾等影響,區(qū)內植被覆蓋受損嚴重;爾后,隨著“天然林保護”“退耕還林”和“退牧還草”等工程的實施,地區(qū)植被生態(tài)環(huán)境逐漸有所改善。因此,本文以MODIS NDVI數(shù)據(jù)為基礎,結合Sen氏趨勢、Hurst指數(shù)和地理探測器,完成其近20年植被覆蓋動態(tài)變化及驅動力分析。研究期望該成果能為區(qū)域植被保護與恢復建設提供理論參考。

1 研究區(qū)概況

川西北高原區(qū)位于青藏高原東南緣,包括石渠、德格和甘孜等8個縣,面積約8.2×104km2。氣候以高原山地氣候為主要類型,由于海拔落差大,導致氣候立體變化明顯,其整體特征表現(xiàn)為山地濕冷,光照豐富,降雨量少。該地區(qū)年均氣溫-7.36～9.07℃,年降雨量約447.54～859.93 mm,植被類型主要為山地草甸。在地殼活動頻繁、地形陡峭和氣候條件多樣等共同作用下,區(qū)內自然災害頻發(fā),加之區(qū)域土地生產力低下,植被生長季短,植被遭破壞后不易恢復,這對自然環(huán)境變化也產生了重大影響。

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

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

基本數(shù)據(jù)包括:2000—2020年植被NDVI為MODIS13Q1產品,源于NASA;90 m的DEM源于地理空間數(shù)據(jù)云;2020年土地利用類型(矢量)、植被和土壤類型(矢量)、人口密度和GDP km網格分布數(shù)據(jù),均源于中國科學院地理資源環(huán)境據(jù)中心(https:∥www.resdc.cn/);各氣象站觀測數(shù)據(jù)由國家氣候中心提供。以最大化合成(MVC)法[15]獲取各年份NDVI;以Anusplin4.37為數(shù)據(jù)處理平臺,采用薄盤樣條函數(shù)插值實現(xiàn)氣象數(shù)據(jù)柵格化[16]。

2.2 研究方法

2.2.1 斜率 Sen氏趨勢主要用于分析時間序列數(shù)據(jù)未來發(fā)展傾向,作為一種數(shù)據(jù)統(tǒng)計分析方法,它能較有效降低個別數(shù)據(jù)缺失或形態(tài)異常對計算結果準確性的干擾[17-18]。公式如下:

式中:xi和xj分別代表第i和j年的NDVI;n為時間長度。若Slope>0,表明隨著時間年限增長植被NDVI呈現(xiàn)數(shù)值上升且植被覆蓋也表現(xiàn)為改善的發(fā)展趨勢;相反,若Slope<0,植被NDVI表現(xiàn)為數(shù)值下降且植被覆蓋呈現(xiàn)退化的發(fā)展傾向。

作為一種非線性的數(shù)據(jù)顯著性檢驗方法,Manna-Kendall檢驗不需要樣本數(shù)據(jù)滿足特定分布,還能較有效規(guī)避異常像元的干擾[11,19]。研究將Sen氏趨勢和顯著性檢驗相結合,完成該地區(qū)植被NDVI發(fā)展趨勢顯著性檢驗(表1)。

表1 植被NDVI斜率變化與顯著性檢驗Table 1 Vegetation NDVI slope change and significance test

2.2.2 Hurst指數(shù) 為了解植被覆蓋未來發(fā)展趨勢的持續(xù)狀態(tài),研究引入Hurst指數(shù)模型。它在植被覆蓋變化研究中已被廣泛應用[20-22]。參考已有成果[19]將Hurst數(shù)值分為:0

2.2.3 地理探測器 地理探測器是用于分析地理現(xiàn)象空間分異性并對其背后驅動力進行探索的工具[23],核心思想為:如果某因子對地理現(xiàn)象具有重要影響,那么它們在空間分布上應具有一定吻合性[24-25]。

因子探測器用于分析地理因變量Y的分異性及各驅動因素X對Y的解釋力。

(2)

(3)

交互探測器用于判定,若任意兩個驅動因素X1和X2發(fā)生交互時,該協(xié)同作用對Y分異性的解釋力為增強、減弱或未產生影響[26-27]。核心思想為:先分別計算因素X1和X2對Y的解釋力q(X1)和q(X2);其次,計算X1∩X2對Y的解釋力q(X1∩X2);最后,對q(X1)、q(X2)與q(X1∩X2)進行對比,完成解釋力類別判定(表2)。

表2 交互作用判定依據(jù)Table 2 Basis for determining interaction

風險探測器用于識別兩個影響因素子區(qū)間屬性是否具有顯著差異,并對最適宜植被生長的區(qū)間進行探測[28]。

在參考相關成果[2,13,29-30]基礎上,結合數(shù)據(jù)的可獲得性、代表性和整體性,最終從地形、植被、土壤、氣候和社會經濟等方面共選擇11個指標。其中,氣候和社會經濟要素(年均氣溫、年降雨量、年均日照時長、人口密度和GDP)均取近20年平均值運算。

由于地理探測器在進行驅動力探測時要求各變量均為類型量,土地利用、土壤和植被類型分別參照《中科院土地利用覆蓋分類體系》、《1∶100萬中華人民共和國土壤圖》和《1∶1 000 000中國植被圖集》將數(shù)據(jù)分為6,22,7個一級類別。除去以上3個指標外的其他7個因素,根據(jù)王勁峰等[23]提出的數(shù)據(jù)離散化原則及參考相關成果[2,13,25,29],以“自然間斷點”將各指標分為9類(高程、年均氣溫、年降雨量和年均日照時長)和5類(人口密度和GDP),以“相等間隔”將坡向和坡度分別劃分為10,9類。以3 km×3 km格網為采樣點完成各變量信息采集。

3 結果與分析

3.1 2000-2020年植被NDVI空間分布格局

為較全面掌握該區(qū)植被覆蓋空間分布特征,在參考前人成果[2,10,11]基礎上,將該地區(qū)近20年植被NDVI平均值分為:低植被覆蓋區(qū)(0～0.2)、中低植被覆蓋區(qū)(0.2～0.4)、中植被覆蓋區(qū)(04～0.6)、中高植被覆蓋區(qū)(0.6～0.8)和高植被覆蓋區(qū)(0.8～1)共計5個區(qū)間(圖1)。

圖1 研究區(qū)植被覆蓋

由圖1可知,川西北高原區(qū)植被覆蓋整體表現(xiàn)出較明顯分布差異,具體表現(xiàn)為:東部地區(qū)植被NDVI數(shù)值相對較高,越向西部靠近則NDVI值相對越低。

低植被覆蓋和中低植被覆蓋區(qū)全域面積占比僅約3.22%,結合地形、土地利用類型和植被類型等數(shù)據(jù)分析可知,它們主要分布在海拔約5 000 m以上的地區(qū),這些區(qū)域地勢險峻,主要以裸巖石質地、沙地或裸土地為主要土地生態(tài)景觀,除部分區(qū)域分布有極少量高山稀疏植被外,其余大多地區(qū)無植被生長。中植被覆蓋區(qū)全域面積占比約7.70%,主要集中分布于石渠縣北和西南部地區(qū),此外在德格縣西南部以及色達和壤塘縣交界處等部分地區(qū)也有零散分布;結合相關資料分析可以發(fā)現(xiàn),這些區(qū)域海拔高度介于4 200～5 000 m,主要分布著高寒草甸、灌叢和高山稀疏植被等,土地利用類型則以低覆蓋度草地和灌木林地為主,少部分區(qū)域夾雜著中覆蓋度草地。相比下,中高和高植被覆蓋區(qū)全域分布范圍最廣,面積占比達89.08%,在海拔小于4 200 m的大部分地區(qū)均有分布。

3.2 2000-2020年植被NDVI趨勢變化分析

結合川西北高原區(qū)植被NDVI斜率計算結果(圖2A)分析可知:近20年內,整個高原區(qū)植被NDVI的Slope<0的地區(qū)全域面積占比僅約28.62%,而Slope>0的地區(qū)面積占比卻高達71.38%;其暗示近20年內整個高原區(qū)植被NDVI數(shù)值整體呈現(xiàn)升高的發(fā)展趨勢。

圖2 植被NDVI斜率變化(A)及顯著性檢驗(B)

結合顯著性檢驗(圖2B)分析可知:21年內,植被覆蓋呈極顯著和顯著退化的區(qū)域在整個高原區(qū)分布范圍最窄,全域面積占比僅1.64%,主要集中分布于石渠縣和色達縣中部,此外在阿壩縣、紅原縣和若爾蓋縣交界處也有極小范圍零散分布;顯著和極顯著改善地區(qū)全域面積占比約17.80%,主要集中分布于若爾蓋縣、紅原縣和阿壩縣全境的大部分地區(qū)、石渠縣北和西南部部分區(qū)域,剩余小部分零散分布于其他區(qū)縣;和以上兩個類別相比,植被覆蓋無顯著(改善或退化)變化的區(qū)域全域分布范圍最廣,占比約80.56%,整個高原區(qū)均有分布;以上成果分析表明,在2000—2020年內,整個高原區(qū)植被覆蓋整體發(fā)展趨勢較好。

3.3 植被NDVI的可持續(xù)發(fā)展穩(wěn)定性分析

為了解整個高原區(qū)未來植被覆蓋發(fā)展趨勢的可持續(xù)性,完成該地區(qū)近20年植被NDVI的Hurst指數(shù)計算(圖3)。

圖3 植被NDVI的Hurst指數(shù)(A)及未來變化趨勢預測(B)

結合圖3分析可知:全時段內,高原區(qū)植被NDVI的Hurst指數(shù)介于0.088 7～0.890 7間,Hurst<0.5的區(qū)域約占整個高原區(qū)總面積的74.91%,這暗示著這些地區(qū)植被NDVI在未來的發(fā)展趨勢不具有可持續(xù)性,即當前變化趨勢為改善或退化而將來則轉變?yōu)橥嘶蚋纳啤?/p>

為更詳細了解高原區(qū)植被NDVI發(fā)展趨勢的穩(wěn)定性和可持續(xù)發(fā)展狀態(tài),將Hurst指數(shù)與植被NDVI顯著性檢驗結果進行柵格疊加。植被NDVI未來呈現(xiàn)持續(xù)極顯著和顯著退化的地區(qū)全域面積占比約14.12%,它們主要集中分布于石渠縣(中部、北部和西南部)小部分地區(qū),紅原與阿壩縣交界的大部分范圍以及若爾蓋縣北和東部大部分。結合土地利用類型、地形和人類活動等實際分析可發(fā)現(xiàn):驅動以上地區(qū)植被覆蓋未來發(fā)展趨勢處于持續(xù)退化的原因可能有所差異;石渠縣在全域海拔相對較高、氣溫相對較低且降雨量也相對較少,這些自然條件均不利于植被生長,地類以低覆蓋度草地為主,先天自然環(huán)境的限制加之不合理人類活動(坡耕地過度開墾和高山草甸區(qū)過度放牧等)的干擾都可能驅使域植被覆蓋持續(xù)退化;紅原、阿壩和若爾蓋縣在整個高原區(qū)內地形和氣候條件相比石渠縣雖然較適宜植被生長,但植被類型主要為高山草甸,人類活動和社會經濟發(fā)展方式卻比石渠縣復雜多樣,不合理人類活動對植被生長環(huán)境干擾程度也更相對顯著。極顯著和顯著持續(xù)改善區(qū)全域面積占比約5.31%;相比下,未來發(fā)展狀態(tài)持續(xù)處于無顯著(改善或退化)變化地區(qū)在全域內分布范圍最廣,約占整個地區(qū)總面積的80.58%。

綜上,未來川西北高原區(qū)植被覆蓋變化趨勢處于持續(xù)無顯著變化的地區(qū)仍占據(jù)全域的絕大部分,而明顯退化卻比明顯改善地區(qū)全域總面積占比高出了8.81%。

3.4 植被NDVI驅動因素分析

3.4.1 因子探測分析 研究利用地理探測器完成川西北高原區(qū)植被NDVI空間分布格局驅動因子探測分析(表3)。

表3 植被NDVI因子探測結果Table 3 Vegetation NDVI factor detection results

結合表3分析可知:各驅動因子對植被NDVI空間分布總體格局差異形成的影響力存在較明顯差異;各因子的q值呈現(xiàn)高程>年均氣溫>土壤類型>植被類型>年均日照時長>土地利用類型>年降雨量>坡向>人口密度>坡度>GDP;高程、年均氣溫和土壤類型的q值分別為0.453 6,0.353 3,0.300 7,解釋力均在30%以上,它們3個基本占據(jù)了絕對優(yōu)勢,客觀程度可判定這3個因素是驅動高原區(qū)植被覆蓋空間分布總體格局差異形成的最主要驅動因素;植被類型的q值為0.173 8,解釋力在15%以上,為主要因素;解釋力小于15%,但大于10%的因子中年均日照時長、土地利用類型和年降雨量的q值分別為0.119 6,0.112 9,0.108 2,它們的解釋力基本相同,均為次要因素;坡向的解釋力小于1%,驅動作用非常有限;坡度、人口密度和GDP解釋力均在1%以下,且均未通過顯著性檢驗,說明若僅從靜態(tài)角度分析,它們幾乎未對高原區(qū)植被覆蓋空間分布總體格局差異的形成產生影響作用。

3.4.2 交互探測分析

交互探測結果表明,當各驅動因素產生交互時,其協(xié)同作用為線性或非線性增強兩種類型,不存在減弱或者獨立(圖4)。

圖4 因子交互探測結果

結合圖4分析可知:當高程和其他因素交互時,其協(xié)同作用對植被NDVI總體空間分布格局形成的影響作用均大于單因素的驅動作用,這也間接證明了若僅從單因素驅動作用分析,高程對高原區(qū)植被NDVI空間分布總體格局變化的影響作用最顯著,該結果與因子探測基本吻合。高程∩年降雨量(q=0.5013),是高程與其他因子交互類型中最顯著的類別,這暗示著在特定高程區(qū)間范圍內,若年降雨量也滿足一定要求,那么該地區(qū)植被覆蓋空間布局變化受高程影響更明顯。相似的,年均氣溫∩年降雨量(q=0.5030),解釋力已經達到50.30%,在各因子交互類型中也表現(xiàn)地非常顯著,這表征出特定水熱條件對植被生長起著重要影響作用。此外,人口密度和GDP僅分別與高程、年均氣溫和土壤類型3個因素產生交互時,它們協(xié)同作用的解釋力才會達到30%以上。當人口密度和GDP在與除高程、年均氣溫和土壤類型3個因素外的其他因子交互時,它們的解釋力均在20%以下,這也表明僅從靜態(tài)角度分析,人口密度和GDP對高原區(qū)植被覆蓋空間總體布局的影響作用非常有限。

3.4.3 適宜范圍或類型分析 結合圖5分析可知:各因素適合植被生長的區(qū)間或者類別也存在顯著差異;地形(高程、坡度和坡向)系統(tǒng)最適合植被生長的區(qū)間分別為3 619.00～3 859.00 m,13.67°～20.51°和0°～35.10°,它們各自對應的植被NDVI均值分別為0.819 6,0.759 2,0.770 5;氣候系統(tǒng)(年均氣溫、年降雨量和年均日照時長)在2.55～3.59℃,744.15～794.58 mm和2 255.67～2 339.35 h最適宜植被生長,對應區(qū)間植被NDVI均值分別為0.823 4,0.791 5,0.781 7;林地植被生長狀態(tài)最好,NDVI均值為0.812 9;相似,闊葉林和棕壤也為植被生長的最理想類型,NDVI均值分別為0.838 9,0.846 8;人口密度和GDP數(shù)值越小的地區(qū)越適合植被生長。

圖5 因子風險探測結果

研究還發(fā)現(xiàn),當高程<3 619.00 m和坡度<13.67°時,植被NDVI均值隨著它們數(shù)值增加而升高;相反,若高程>3 859.00 m、坡度>20.51°時,植被NDVI的均值則呈現(xiàn)相反的變化形態(tài),該成果與韓繼沖等[12]采用格網模型完成長江上游植被NDVI與地形因子響應關系分析的結果基本吻合。相似,靜態(tài)視角下植被NDVI與年均氣溫和年降雨量的響應模式也類似拋物線,這也與前人成果基本吻合[12];然而,植被NDVI與年均日照時長的響應關系卻相對復雜,類似為分段函數(shù)。川西北高原區(qū)為四川省重要的水源涵養(yǎng)和林草資源區(qū),隨著系列生態(tài)工程的實施,區(qū)內森林資源得到較有效恢復且草地退化也被較有效遏制。因此,在林地和草地類型上植被NDVI均值相對較高。然而,隨著人口密度和GDP數(shù)值升高,植被NDVI均值卻呈下降變化趨勢,這主要是高原區(qū)社會經濟生產方式相對傳統(tǒng)、簡單,人口密度和GDP數(shù)值越高的地區(qū)人類活動明顯,對植被的擾動干擾也相對越突出。不同土壤和植被類型區(qū)植被生長狀態(tài)也存在較明顯差異,不同類型的植被對光的敏感性、抗旱(寒)能力等都有較明顯不同;各類型的土壤其自身的保水性、可侵蝕性和養(yǎng)分均存在差異,對植被生長也起著促進或抑制的不同作用。

4 討論與結論

4.1 討 論

近20年內,高原區(qū)植被覆蓋整體呈現(xiàn)良好的發(fā)展趨勢,該結果與前人[10,31]基本一致;高程、年均氣溫、土壤和植被類型對植被NDVI分布變化的影響作用非常顯著,這與彭文甫等[2]成果基本一致。

研究區(qū)地處高原,貧瘠的自然條件對區(qū)域植被生長造成了一定限制。分析各因子對植被覆蓋空間分布格局形成與變化的影響作用可發(fā)現(xiàn):雖然地形、植被和土壤類型主導了區(qū)域植被覆蓋空間布局的形成,但這些因素自身穩(wěn)定性也相對較強,20年內幾乎未發(fā)生較明顯改變,對植被覆蓋布局變化產生明顯驅動作用較弱可忽略;同理,結合各氣象站點觀測的年均氣溫和年均日照時長資料分析,這兩個因素的穩(wěn)定性也較強,也未對植被覆蓋布局的動態(tài)變化產生較明顯驅動作用。適宜植被生長的降雨條件也相對苛刻,影響作用也相對有限。

綜上且結合區(qū)域實際分析可知,人類活動對該地區(qū)植被覆蓋空間布局改變的影響作用相對穩(wěn)定且突出,這與鄭勇等[32]研究成果基本一致。一方面,約2000年起,國家對川西北高原區(qū)相繼實施了系列生態(tài)保護工程。結合土地利用類型遙感解譯成果分析可知:約2000年起,國家川西北高原區(qū)相繼實施了系列生態(tài)保護工程,21年內該地區(qū)林地面積共增加了849.07 km2;它也是四川省重要的生態(tài)牧區(qū),草地類型以高寒草甸為主,若它被破壞想要再次進行修復難度極高,隨著眾多環(huán)保工程的實施,區(qū)域草地退化現(xiàn)象被有效遏制;林草地資源的有效保護與恢復較大程度促進了區(qū)域植被覆蓋狀況的好轉。另外一方面,受如過度放牧和超承載力耕作等不合理人類活動的影響,極小部分地區(qū)植被仍然處于被破壞的狀態(tài)。

相比以往成果,本研究采用最大值合成完成各年份植被NDVI數(shù)據(jù)獲取,較有效降低了噪聲對數(shù)據(jù)準確性的干擾;將地理探測器用于驅動力分析,較有效彌補了以往成果僅從氣候或地形等單系統(tǒng)對驅動力進行討論,研究結果存在不詳實的局限;論文不僅較詳細分析了各因素對植被覆蓋空間分布總體格局差異形成的影響作用,還結合區(qū)域實際對驅動其布局變化的因素進行了進一步探討,研究的價值也得到了進一步升華;然而,在進行數(shù)據(jù)進行離散化時,不同的數(shù)據(jù)分類方法也會導致計算結果存在一定差異,因此結合區(qū)域實際嘗試多種不同數(shù)據(jù)分類方式,綜合尋找更恰當?shù)臄?shù)據(jù)離散方法也將有助于進一步提高研究的準確性。

4.2 結 論

(1) 從空間分布變化看,川西北高原區(qū)植被覆蓋整體呈現(xiàn)自西向東逐漸升高的變化,越向東邊發(fā)展其整體的植被覆蓋度越高。

(2) 從時間角度分析:近20年內,僅有19.44%的地區(qū)植被覆蓋呈現(xiàn)出顯著及極顯著退化的發(fā)展趨勢;未來整個高原區(qū)仍有80.56%的地區(qū)植被覆蓋將持續(xù)處于無顯著改變的變化趨勢。

(3) 地理探測結果顯示:高程、年平均氣溫和土壤類型對川西北高原區(qū)植被覆蓋空間分布總體格局差異形成的驅動解釋力均在30%以上;各因素交互產生的協(xié)同作用均高于單因素的驅動效果;各因素適宜植被生長的區(qū)間也各不相同。