怒江-薩爾溫江流域植被覆蓋時空變化趨勢及驅動力

鐘旭珍,王金亮,鄧云程,李 杰,吳瑞娟,董品亮

1 云南師范大學地理學部,昆明 650500 2 內江師范學院地理與資源科學學院,內江 641100 3 云南師范大學-西南聯合研究生院,昆明 650500 4 云南省高校資源與環境遙感重點實驗室, 昆明 650500 5 云南省地理空間信息工程技術研究中心, 昆明 650500 6 美國北德克薩斯大學地理系, 丹頓 76203

近百年來全球氣溫持續上升,根據IPCC發布的《氣候變化2021:公眾摘要》報告顯示:與19世紀末(工業革命前)相比,2011年至2020年地球表面的平均氣溫升高了1.1℃,并且比過去12.5萬年的任何時候都高[1]。歷年的IPCC發布的報告均認為:近百年來,人類過度的開荒及大規模活動,導致所排放的溫室氣體有增無減,這是造成全球變暖最主要的原因[2]。人類活動不僅導致了全球氣候的變暖,還誘發了如極端天氣、植被退化、水土流失等生態安全問題。生態環境是人民生存的基礎條件,人類可持續發展的保障。而植被作為生態系統主要的組成部分,對區域生態系統環境變化有著重要指示,對全球碳循環和氣候調節具有非常重要的作用。植被覆蓋度(FVC)是指植被在地面的垂直投影面積占總面積的百分比,它量化了植被的茂密程度,反映了植被的生長態勢,是描述生態系統變化的重要基礎數據[3]。怒江-薩爾溫江是東南亞最長的自由流動河流和最重要的跨界河流之一,在中國-中南半島經濟走廊建設和中緬經濟走廊建設和生態保護中發揮著不可替代的作用。在過去的50年里,該流域冬季和春季氣溫呈上升趨勢,氣溫上升的幅度隨著海拔高度的升高而增加[4]。然而,目前對該流域關注不夠,尚沒看到其植被覆蓋研究報道。因此,研究分析其植被覆蓋時空演變,并探討其驅動力,可為其生態環境保護和社會經濟發展提供重要指導作用[5]。

估算植被覆蓋度(FVC)是植被覆蓋相關研究的重要基礎,有多種估算方法,如地表實測法、經驗模型法、植被指數法、決策樹分類法等,部分方法存在對實測數據依賴性強、計算繁瑣、需要大量樣本數據等缺點,特別對于偏遠地區來說很難實現[6],而遙感技術的快速發展以及遙感大數據時代的到來,使得在大時空尺度上精確監測偏遠地區的植被行為成為可能[7],其中MODIS NDVI產品被大量學者用于植被覆蓋度的計算及其動態變化的定量分析,而像元二分模型法是常用的FVC估算方法,其操作較為簡單,可直接利用NDVI數據來進行計算[8-10]。目前針對植被覆蓋演變規律的研究主要是運用一些常規的線性趨勢分析方法,如線性回歸分析方法[11-12]、相關性分析法[13-14]、Theil-Sen斜率與Mann-Kendall檢驗結合的趨勢分析法[15],很少有學者關注其非線性變化特征。然而,植被覆蓋的時空分布與變化受多方面的因素影響,使得其演變規律也具有復雜性,并不僅僅呈現為簡單的線性變化,而BFAST法可用于分析不同類型的時間序列,對長時間序列過程中趨勢的突變檢測具有很好的應用效果,已有學者將其用于植被的非線性特征[15-17],以及氣溫變化[18]、氣溶膠變化[19]等領域研究中,并顯示了其較好的適用性[20-21]。對于南北跨度極大的怒江-薩爾溫江流域,其不同緯度區域的植被變化特征及原因也會存在差異,因此有必要運用BFAST法來更加準確地探究其植被變化的非線性特征及其未來發展趨勢。

對于植被未來可持續性的研究,Hurst指數法是廣大學者常用的方法[22-23]。耦合BFAST法與Hurst指數可以實現植被未來趨勢的更加精確表達。驅動植被變化的因素,有人為因素和自然因素,以往學者大多采用偏相關分析、回歸分析等方法研究植被覆蓋度變化的影響因素[24-25],這些方法假設了植被與各環境因子存在線性關系。然而植被在環境中的生長過程影響因素復雜,可能不存在嚴格的線性關系[26]。而由王勁峰等學者提出的地理探測器模型(Geographical Detector,簡稱Geodetector),不僅可用于探測空間分異性并識別其主要的影響因子,還能實現不同因子間的交互作用探測[27]。此方法無線性假設,可以避免多重共線性問題,在計算驅動因子的貢獻率時結果較為準確[28]。其改進的模型OPGD模型(Optimal Parameters-Based Geographical Detector,OPGD)加入了參數優化模塊,可根據自變量的特征進行多種方式和多種斷點數目的離散化處理,選取出最優的離散化方式和最優的斷點數目,提高各類別間的顯著差異性[29]。因此,研究采用該方法進行植被非線性特征的驅動因子探測。

基于以上分析可知,盡管已經開展了許多關于植被時空演變與驅動因素的研究,但目前還沒有學者采用BFAST法研究其歷史非線性趨勢的同時結合Hurst指數表達其未來更加精確的可持續性,也就是考慮了植被非線性特征的未來趨勢預測。另外,怒江-薩爾溫江是東南亞一條重要的國際河流,對南亞、東南亞的氣候調節具有重要作用,但由于其南北跨度大,涉及中國、緬甸和泰國,對其進行相關研究尤其數據收集較為困難,目前還沒有學者對該流域進行植被覆蓋的相關研究。因此,研究以該流域為研究對象,基于MODIS NDVI數據,在遙感與GIS技術的支持下,采用像元二分模型估算其近22年的植被覆蓋度,并借助 BFAST01模型及 Hurst 指數等,對區域內植被覆蓋的時空分布、非線性變化趨勢和未來可持續性進行分析,并結合坡度、海拔、降雨、氣溫、人口、GDP、土地類型等自然和人類活動因子,引入OPGD對影響研究區植被覆蓋度的主要因素進行定量探究,以期為研究區的生態環境保護和可持續經濟發展提供科學依據。

1 研究區概況及數據來源

1.1 研究區概況

怒江-薩爾溫江發源于青藏高原上的唐古拉山脈,經中國云南(保山、臨滄)流入緬甸,最后注入印度洋的安達曼海,經過中國、緬甸和泰國(圖1)。該流域形狀極為狹長,河流全長約3200 km,流域總面積為32.5萬km2,它在中國的部分被稱為怒江,其余的部分稱為薩爾溫江,是東南亞最長的自由流動河流,也是最重要的跨境河流之一。楊帆等的研究中,將嘉玉橋以上劃分為上游,嘉玉橋到中緬邊界為中游,中國邊界以南為下游[30]。由于該流域南北跨度和海拔高差都很大,橫跨多個氣候帶和生態區,其氣候多樣,地形復雜,不同區域的地形、氣候、土地利用類型都存在很大的空間差異性,本研究根據流域海拔特征并結合行政區劃界限,將云南和西藏交界線以北區域劃分為上游,屬于青藏高原,為典型的高原山地氣候區,氣候變化較為敏感;云南西藏交界線到中緬邊界劃分為中游,主要為高山峽谷區,即狹窄又陡峭,氣候垂直差異明顯,流經三江并流保護區;中國邊界以南為下游,海拔相對北部較低,地勢較為平坦,屬于熱帶季風氣候,以此分區作為后期相關研究和分析的基礎。

1.2 數據源及預處理

研究數據主要包括MODIS NDVI 數據,地形數據,土地利用數據,氣象數據以及人類活動數據,研究區流域邊界和行政區劃數據等。詳述如下:

1.2.1NDVI 數據

研究使用2000-2021年MODIS NDVI數據,通過GEE(Google Earth Engine,GEE)遙感大數據云平臺進行調用、預處理及下載,數據原網站為USGS (https://lpdaac.usgs.gov/products/),其時間分辨率為16d,空間分辨率為250m。在GEE平臺中對其進行鑲嵌、投影、轉換和裁剪等預處理,為了消除云、大氣效應、掃描角度和太陽天頂角等對影像質量的影響,研究采用了最大值合成法MVC(Maximum Value Composite)來進行NDVI數據的合成[15],該方法被廣泛用于植被動態研究的預處理中[31]。借助ENVI IDL平臺,選用像元二分模型(DPM)來計算研究區月度及年度FVC,用于后續植被覆蓋的相關分析。

1.2.2地形數據

研究使用的地形數據為2000年采集的SRTM Digital Elevation Data Version 4,分辨率為90m,數據生產于 NASA (https://srtm.csi.cgiar.org/),通過GEE平臺對其進行拼接、裁剪和下載,并通過ArcGIS 10.2軟件提取如海拔、坡度、坡向等地形因素作為后續植被覆蓋驅動力分析的指標。

1.2.3LULC 數據

研究使用的土地利用數據為ESA WorldCover 10m v100,來自esa (https://esa-worldcover.org/en),也是通過GEE平臺對其進行投影轉換和裁剪等預處理和下載。該產品將土地利用類型分為了11類,比較詳細的表示了不同地表覆蓋類型,我們將其作為植被覆蓋的人類活動驅動因子之一。

1.2.4氣象數據

研究使用的氣象數據主要包括2000-2021年的氣溫、降水、蒸散發以及氣候類型。其中氣溫數據由ERA5 提供的第五代 ECMWF 全球氣候大氣再分析數據,通過GEE平臺進行預處理與下載。我們選用的是mean_2m_air_temperature,單位為K;降水量數據選用的是Climate Hazards Group InfraRed Precipitation With Station Data 2.0 (https://chc.ucsb.edu/data/chirps), 單位為mm,分辨率為0.05°。蒸散發數據為GLEAM ET產品,空間分辨率為0.25°,下載于gleam (https://www.gleam.eu/#downloads)。蒸散發可以反映植被的干旱情況,對植被生長狀況有一定影響,因此研究選取該指標作為驅動因子之一。氣候類型數據為全球柯本氣候類型空間分布數據集,空間分辨率為0.1°,下載于國家地球系統科學數據中心(http://www.geodata.cn/),將其作為植被覆蓋驅動力分析因子之一。

1.2.5人類活動等其他數據

人口格網數據來源于world pop (https://hub.worldpop.org/geodata/listing?id=64)。單位為people/km2,空間分辨率為1km。全球GDP數據來源于GitHub (https://github.com/Nowosad/global_population_and_gdp), 分辨率為0.5° × 0.5°。中國區域人口、GDP數據來源于中國科學院資源環境科學數據中心(http://www.resdc.cn)。其他數據如研究區行政區劃,流域邊界等分別來源于中國科學院資源環境科學數據中心 (https://www.resdc.cn/data.aspx?DATAID=259) 和HydroSHEDS 網站 (https://www.hydrosheds.org/)。數據來源如表1,驅動因子指標如圖2。

表1 研究數據來源與描述Table 1 Source and description of research data

圖2 FVC驅動因子空間分布圖Fig.2 Spatial distribution of FVC driving factorsET:蒸散發 Evapotranspiration;FVC:植被覆蓋度 Fractional vegetation cover

2 研究方法

2.1 像元二分模型

像元二分模型是一種簡單實用的遙感估算模型,對綠度信息具有很好的敏感性,且對不同氣候區、植被類型等的植被信息具有很好的識別效果,被廣泛用于計算植被覆蓋度。其計算植被覆蓋度的公式如下:

(1)

式中,NDVIsoil和NDVIveg表示全部為裸土和全部為植被的像元NDVI值。理論上,NDVIsoil為0,NDVIveg為1。然而,實際中植被覆蓋的影響因素復雜,NDVIsoil和NDVIveg的值并不等于0和1,而是受圖像本身質量的影響。參考相關研究[23,32],根據MODIS NDVI時間序列圖像,對NDVI圖像進行頻率直方圖統計,去除異常值以及歸一化處理,并取累計頻率在2%和98%時的NDVI值分別作為NDVIsoil和NDVIveg參數,然后按照公式(1)計算研究區FVC。

2.2 Hurst指數

Hurst 指數是定量描述時間序列長程依賴性的有效方法[33]。近年來在植被覆蓋變化的未來可持續性研究中得到廣泛應用[34-35]。研究利用matlab計算出H值即Hurst指數,來判斷植被覆蓋時間序列是否存在持續性。H值介于0-1間,參考相關研究成果[23,36],其分類有三種:0

2.3 BFAST模型

BFAST是一種迭代算法,它利用分段線性趨勢和季節性模型將時間序列分解為趨勢、季節性和殘差組分,并能檢測趨勢和季節性成分的突變[7,37]。BFAST假設非線性可以用分段線性模型來近似,因此,當一個序列可以更好地用多個線性段表示,而不是單一的單調趨勢時,就可以使用結構變化測試確定一個斷點并確定日期[26]。相比其他突變檢測方法,在時間序列內,BFAST受季節差異和噪聲的影響較小,可以更快地檢測到變化[19]。Verbesselt等對該方法進行了比較完整的描述和介紹[38-39]。在數學上,BFAST的趨勢分量和季節分量通過以下分解模型得到:

Yt=Tt+St+et,t=1,…,n

(2)

式中,Yt為t時觀測到的值,Tt為趨勢組分,St為季節性組分,et為殘差組分。其中趨勢組分可用下式表示:

Tt=ai+bit(t=1,2,…,n)

(3)

季節性成分可用下式表示:

(4)

式中,ai和bi為趨勢項系數,γi為振幅,δi為分段數量,f為頻率,其中振幅、分段數量為未知項,而頻率為已知項。

季節組分和趨勢組分中突變點的識別需要確定突變點的數量和突變點在時間序列中的位置。此處利用最小二乘移動求和從季節組分和趨勢組分中檢驗判斷是否存在突變點[40]。首先,采用普通最小二乘方法(Ordinary Least Squares, OLS)在一個不斷增加的移動窗口內進行檢驗判斷是否存在斷點,如果根據統計檢驗,窗口內的普通最小二乘殘差的移動和很大,則被描述為一個斷點。接著,通過最小化所有窗口上斷點之間的最小二乘移動總和(OLS-MOSUM),以此來確定斷點的位置[19,38]。

根據BFAST原理,本研究采用了BFAST的改進模型BFAST01。BFAST01模型相比BFAST模型的優勢在于,它同時考慮了季節模型和趨勢模型,只檢測時間序列中最具影響的趨勢變化,而不是多個較小的變化,也就是說BFAST01只檢測0個或者1個斷點,其為評估長期季節性和年際變化的主要轉折點提供了一種有效的、全面的變化檢測方法[26,41]。本研究利用R環境,從GitHub加載最新的bfast包(https://github.com/),完成所有BFAST01的統計分析[42],提取出植被覆蓋的趨勢顯著性、變化強度、趨勢斷點發生的時間和次數、突變類型、突變年份等。參考相關研究[15,19,26],將BFAST01檢測到的植被覆蓋變化趨勢類別分為8類,如表2所示。

表2 BFAST01檢測的FVC變化趨勢類型Table 2 Types of FVC change trends detected by BFAST01

2.4 地理探測器

地理探測器Geodetector,是一種統計工具,可以探測地理現象的空間分異性并揭示其背后驅動力,既適用于點數據,也適用于面數據,對于面數據一般要先進行重分類和離散化處理[27]。其包括4個探測器,我們選用因子探測器和交互探測器來探測植被覆蓋變化的影響因素。

因子探測器,用于探測因變量Y的空間分異性,及某個因素X對因變量Y空間分布的解釋力,其大小可用q值來度量(0≤q≤1),本研究中,若某影響因子的q值越接近于1,說明該因子對FVC的解釋力越強,反之越弱[27,43-44]。

交互探測器可以用于識別不同因子之間交互對植被覆蓋度變化的影響是增強還是減弱,或是獨立[27]。其交互方式見表3。

表3 地理探測器交互作用方式Table 3 Interaction modes of geodetectors

3 結果與分析

3.1 怒江-薩爾溫江流域植被覆蓋時空特征分析

3.1.1空間分布特征

通過像元二分模型計算得到FVC,并利用GIS空間分析可視化工具得到研究區2000-2021年平均FVC空間分布圖(圖3),對其進行分析可知,怒江-薩爾溫江流域低植被覆蓋度(0-0.2)、較低植被覆蓋度(0.2-0.4)、中植被覆蓋度(0.4-0.6)、較高植被覆蓋度(0.6-0.8)、高植被覆蓋度(0.8-1)占比依次為:5.41%,5.66%,10.41%,21.46%,57.06%,FVC分布具有明顯的空間異質性,從圖3可以看出,下游和中游地區植被覆蓋情況明顯優于上游地區,總體表現為南高北低。整個流域平均FVC值為0.73,高植被覆蓋度和較高植被覆蓋度地區面積占比達到78.52%,說明該流域植被覆蓋情況較好。結合研究區土地利用類型分析發現,低植被覆蓋區域主要分布在海拔較高的青藏高原,主要土地利用類型為苔蘚、裸地、灌木等;高植被覆蓋區域主要在中游和下游地區,土地利用類型以林地、草地和耕地為主。

圖3 怒江-薩爾溫江流域2000-2021年平均FVC空間分布圖Fig.3 Spatial distribution of the average FVC in the Nujiang-Salween River Basin from 2000 to 2021

3.1.2時間變化特征

通過統計研究區年、月、季節的FVC均值,并制作FVC變化折線圖(圖4)。從圖4的年際變化圖可以看出,怒江-薩爾溫江流域2000-2021年植被覆蓋整體呈波動上升趨勢,FVC最大值出現在2021年為0.750,最小值出現在2015年為0.722,在2007年和2015年出現了明顯的低值拐點,但2016年以來,植被生長狀況出現明顯上升趨勢。從研究區連續22年植被覆蓋季節變化圖可以看出怒江-薩爾溫江流域植被變化具有明顯的季節規律性,這與研究區氣候有關。而FVC月際變化圖表明,3月份FVC值最低,3-9月份FVC值隨時間逐漸上升,在9月份的時候FVC值達到最高,隨后呈下降趨勢;整體上夏秋季植被覆蓋狀況較春冬季植被覆蓋狀況好,不同季節FVC隨著時間變化呈現較小的波動,春季和夏季波動較秋季和冬季大,但整體較為穩定。

圖4 怒江-薩爾溫江流域FVC時間變化Fig.4 Time variation of FVC in the Nujiang-Salween River BasinFVC:植被覆蓋度Fractional vegetation cover

3.2 怒江-薩爾溫江流域植被覆蓋非線性趨勢特征

運用BFAST01算法對怒江-薩爾溫江流域2000-2021年FVC進行突變檢測,得到FVC的非線性趨勢特征(圖5)。從圖中可知,研究區分為了8種非線性植被覆蓋趨勢類型,其中單調遞增占34.63%,是占比最多的,在流域上、中、下游各區域均有分布;單調遞減類型占比12.22%,主要分布在流域下游區域;“單調遞增+”占比2.3%,主要分布在中緬邊界;“單調遞減-”占比最少為0.63%;“中斷-+”占比22.91%,排名第二,主要分布在流域最北端和云南以及緬甸北部;“中斷+-”占比8.19%,主要分布在流域最南角,云南和北部也有些許分布;“反轉+-”占比7.73%,“反轉-+”占比11.39%;通過分析發現,總體呈增加趨勢類型的占比(71.24%)大于總體呈減少趨勢類型的占比(28.76%),表明研究區植被發展趨勢向好。從顯著性檢測結果來看,變化呈現顯著的占比為44.29%,不顯著的占比為55.71%;突變時間檢測表明,發生突變的時間主要集中在2003-2018年,2003年和2018年發生突變的數量較多,其他年份較為均勻,“中斷-+”發生的突變最多,發生突變占比較少的類型是:“單調遞增+”和“單調遞減-”。需要提醒的是,單調遞增和單調遞減沒有發生斷點,因此圖5(突變時間)中空白區域表示沒有發生突變的區域。

圖5 BFAST01檢測結果Fig.5 BFAST01 detection results

為了更加詳細的描述植被覆蓋突變情況,研究還對FVC變化的最大突變強度、最大突變強度發生的時間、首次季節性斷點發生的時間以及斷點發生次數等信息進行了檢測(圖6)。通過分析發現整個流域突變發生最大強度的地區主要在云南以及青藏高原北部,最大突變強度范圍介于-1.2078-1.0287之間,最大突變強度時間主要在2002-2019年,2006年發生最大強度突變的比例是最大的,2013年、2015年和2018年最大突變強度數量也較多,上游地區最大突變強度斷點發生時間變化趨勢與此相似,說明上游地區在全流域FVC的變化趨勢中具有重要作用;而中游和下游區域最大突變強度斷點發生時間變化趨勢較上游地區波動平穩,中游地區突變數量較多的年份為2006年、2010年和2015年,下游地區在2011年最大強度突變數量最多,占下游地區突變總面積的9.90%。各區域最大強度突變面積從2018年開始大面積減少。此外,根據研究區FVC時間變化可知,其植被變化具有很強的季節性,從圖6可以看出,首次發生季節性斷點的時間主要集中在2010-2014年。

此外,怒江-薩爾溫江流域最大突變強度斷點發生突變的次數分別有1次、2次和3次(圖7),表示某點在研究時間段內植被發生非線性突變次數的空間分布,1次的占比最多為77.02%,3次的最少為4.18%,2次的區域主要分布在中緬邊界和最北部。上游地區發生1次斷點的比重最大,占上游地區發生斷點次數總面積的85.72%,而發生3次斷點的占比僅占本區斷點總面積的0.91%;下游和中游地區各斷點次數占本區域斷點總次數面積比重相差不大,也是以1次為主,占比略小于上游地區,而發生3次斷點的次數占比均在7%左右,略大于上游地區。

圖7 植被覆蓋突變頻率分布Fig.7 Frequency distribution of vegetation cover mutations

3.3 怒江-薩爾溫江流域植被未來可持續性

為了表示植被覆蓋的未來可持續性,我們計算了怒江-薩爾溫江流域FVC 的 Hurst 指數,通過統計分析發現,Hurst 指數小于0.5的區域占總面積的2.76%,這些區域具有反持續性,表示過去增加的趨勢將來總體上減少,反之亦然;Hurst 指數大于0.5的區域占比為 94.89%,表明研究區植被覆蓋的正向持續性較強;Hurst 指數為0.5的區域占2.35%,表明未來植被覆蓋發展趨勢不明確。

為了進一步揭示植被覆蓋的非線性變化趨勢及其持續性,將BFAST變化趨勢結果與 Hurst 指數結果進行疊加,得到非線性變化趨勢與持續性的耦合信息(圖8和表4)。耦合結果總共為17種情形,其中“持續性&單調遞增”耦合類型含義表示未來植被與過去發展趨勢一樣,將持續增加;“反持續性&單調遞增”表示未來植被與過去發展趨勢相反,將從增加變為減少,其它情形解釋相似,比較詳細的展示了植被覆蓋的未來發展趨勢。未來將呈持續改善情景的面積占比達到68.99%,持續改善的趨勢類型分別為:“持續性&反轉-+”、“持續性&中斷-+”、“持續性&單調增加+”、“持續性&單調遞增”、“反持續性&單調遞減”、“反持續性&單調遞減-”、“反持續性&中斷+-”和“反持續性&反轉+-”,這些區域占據了流域大部分地區;未來發展將呈持續退化的情景所占比例為29.09%,這些趨勢類型分別為:“反持續性&反轉-+”、“反持續性&中斷-+”、“反持續性&單調遞增+”、“反持續性&單調遞增”、“持續性&單調遞減”、“持續性&單調遞減-”、“持續性&中斷+-”和“持續性&反轉+-”,這些區域主要分布在流域中游和下游;剩下1.92%的區域未來變化趨勢無法確定,主要分布在流域上游區域。持續退化和未來變化趨勢無法確定的區域的植被變化狀況需要研究人員繼續關注。

表4 BFAST01和Hurst指數耦合類型統計Table 4 Statistics of BFAST01 and Hurst coupling types

圖8 植被覆蓋未來發展趨勢Fig.8 Future development trend of vegetation coverage

3.4 怒江-薩爾溫江流域植被覆蓋驅動因子探測

研究選取了海拔、坡度、坡向、蒸散發、氣溫、降水、氣候區劃等自然因素和GDP、人口、土地利用等人為因素作為植被覆蓋的影響因子,分別為X1-X10,植被覆蓋度FVC為Y,作為地理探測器的數據輸入,空間分辨率統一重采樣為250m。由于這些因子及FVC都是連續變量,而地理探測器的輸入數據要求為類別數據,則根據王勁峰等提出的數據離散化方法[27],利用OPGD在R語言環境中進行參數優化,為了使結果更加符合實際,研究同時對離散化進行了人為調試,最終得到比較理想的結果,具體離散化方法和類別見表5。并利用ArcGIS創建漁網,通過多次調試設置采樣間隔為3km,運用Spatial Analyst-提取分析-采樣工具,輸入柵格為所有重分類的X和Y,采樣點為研究區漁網點,將自變量和因變量的值提取到點,作為地理探測器輸入數據。

表5 地理探測器因子離散化方法與類別Table 5 Factor discretization methods and categories of geographic detectors

利用因子探測和交互探測器,得到全流域以及各子流域植被覆蓋驅動力探測結果(表6和圖9)。除中游地區的坡向因子p值為0.0791外,其余因子p值均為0,通過顯著性檢驗。根據表6可知,全流域中海拔的q值是最大的為0.5761,接著影響較大的是氣溫、降水、蒸散發等氣象因子,其次是土地利用類型、GDP、人口等人為因子,人為因子中以土地利用影響最大。各子流域因子探測結果顯示同一個影響因子在不同區域對植被覆蓋的影響具有差異性,縱向上分析,如海拔和土地利用因子對上游和中游區域的影響大于下游區域;而人口、GDP、蒸散發等因子對下游區域的影響則大于上游和中游區域。橫向上分析,上游區域影響最大的是土地利用因子,其次是海拔,而人口、GDP等人類活動因子對于較為偏遠的上游地區影響較小。對于中游區域,對植被影響最大的是海拔,其次是土地利用,氣溫、降水、蒸散發等也是影響中游區域植被覆蓋變化的重要因素,自然因素對中游區域植被的影響大于人為因素。與中游區域相反,下游區域的人口和GDP是對植被影響最大的因素,其次是蒸散發,土地利用和降水影響也較大,而影響最低的是海拔。

表6 怒江-薩爾溫江流域各區域影響因子q值Table 6 Impact factor q values in the Nujiang Salween River Basin

圖9 怒江-薩爾溫江流域各區域影響因子交互作用探測結果Fig.9 Detection results of the interaction of influencing factors in the Nujiang-Salween River Basin

交互探測結果顯示(圖9),任何兩個因子交互影響都會增強其影響力,對于全流域來說,發現其它每個因子和海拔交互都達到最大影響值,說明海拔對植被覆蓋有著很強的影響力,這與整個研究區海拔差異顯著具有關聯性。上游區域和中游區域與全區相似,各因子與海拔的交互作用對植被覆蓋的影響增強作用是最顯著的。對于下游區域來說,各因子交互作用具有差異性,其中氣溫、降水、蒸散發與坡度的交互作用對植被的影響大于它們與其它因子的交互作用;氣候類型、土地利用、人口與蒸散發的交互作用對植被的影響大于它們與其它因子的交互作用;說明坡度和蒸散發對下游區域的植被覆蓋影響具有很大的可利用空間;而GDP與人口的交互作用也大于GDP與其它任何因子的交互作用,成為所有交互作用中影響最大的,這與因子探測的結果具有一致性。

4 討論

4.1 怒江-薩爾溫江流域植被覆蓋時空總體特征

怒江-薩爾溫江流域區位特殊,上、中、下游地區地形海拔和緯度等差異明顯,植被覆蓋空間分布特征也形成了鮮明的對比。從全區域來看,高植被覆蓋區域占比大于50%,而其中上游地區高植被覆蓋區域占比僅占7%左右,平均FVC值為0.53;中游和下游區域高植被覆蓋區域占比達到85%和90%以上,平均FVC值分別為0.86和0.88,說明中游和下游地區拉高了整個區域的植被覆蓋水平,這與研究區地形地貌有關[45-46],這將在驅動力分析中進行探討。

從FVC時間演變特征來看,怒江-薩爾溫江流域植被覆蓋在2015年以前增長趨勢不明顯,但2016年開始發展趨勢向好。這與李杰、王春雅等的研究結論相似[23,47]。但在2007年、2015年和2019年出現了明顯的低值拐點,結合圖4與圖10,發現上游地區FVC年際變化也在2007年和2015年發生了低值拐點,變化趨勢與全區FVC非常相似;而下游和中游區域FVC年際變化增長趨勢明顯,其中下游地區在2004、2012、2018等年份出現低值拐點,查詢資料發現,2004年中緬邊境發生了一場大森林火災,2006年的颶風“馬拉”和2007年的臺風“利奇馬”,2012年和2013年緬甸發生了嚴重干旱,2014年、2015年的厄爾尼諾現象、極端干旱事件等,此外,據聯合國人類事務協調辦公室稱,緬甸在2015年經歷了幾十年來最嚴重的洪水,同時2015年“4·25”西藏也發生了地震災害,以及2018、2019年的極端氣候自然災害如印度西南部的特大洪災等[48],這些事件均嚴重威脅到當地植被的生長狀態,可能是導致FVC出現低值拐點的原因。對于下游地區來說,有研究表明森林火災、農業轉型、公路建設增加等是導致其植被減少的主要原因[49]。而研究區整體植被變化趨勢向好與生態工程的實施和圍欄禁牧政策以及其他人類活動等有關[47]。

圖10 怒江-薩爾溫江流域各流域FVC年際變化Fig.10 Interannual variation of FVC in each sub-basin of the Nujiang-Salween River Basin

4.2 怒江-薩爾溫江流域植被覆蓋非線性特征

從研究區FVC年、月、季節的變化分析以及相關研究表明,植被變化往往受多種因素的影響使其很難呈現嚴格的線性變化[15]。因此,我們在進行研究時應該考慮其趨勢的非線性特征及其影響因素,特別是大尺度長時間序列的演變分析。研究運用BFAST01方法來分析FVC的演變趨勢以及非線性突變,BFAST方法可以很好地適用于突變檢測[50]。在3.2小節中,我們詳細分析了怒江-薩爾溫江流域FVC變化的突變類型、突變時間、突變顯著性與突變次數等特征,將研究區植被變化詳細分為了8種突變類型,每種突變類型的含義見表2。為了進一步分析各突變類型空間分布特征及驗證分類的準確性,我們結合 Google Earth time-series images來進行突變前后年份影像的比對(圖11)。

圖11 突變年前后FVC及Google Earth imagesFig.11 FVC and Google Earth images before and after the mutation year圖中圓圈4、5、6、8分別代表對應的突變類型;黑色矩形框所指圖像代表突變年份的Google Earth images和FVC圖像;紅色矩形框圖像代表突變年份前后的Google Earth images和FVC圖像

從上、中、下游區域分別選擇了某局部區域進行制圖比對,上游區域選擇的代表像元突變類型為 “單遞減-”,突變時間為2017年,從Google Earth time-series images上我們可以明顯看到該區域范圍突變時間前后土地利用從草地和稀疏植被逐漸轉變為建設用地,FVC圖像也可以看出其顏色逐漸變深,也就是FVC值在減小,這就說明了該區域植被覆蓋在減少,這和突變類型“單調遞減-”的特征是相符的。中游區域選擇的代表像元突變類型為 “中斷+-”,該突變類型有1個明顯突變點,且斷點處值突然增大,即顯著正中斷,然后顯著減小,從Google Earth time-series images和FVC圖像比對來看,也符合其突變特征。下游區域選擇的兩個突變類型代表分別為“中斷-+”和“反轉-+”,“中斷-+”特征與“中斷+- 相反,從圖像特征上可以看出該局部區域植被覆蓋在斷點2012年處突然減小,即負中斷,然后增加;而“反轉-+”表示趨勢從顯著減小轉換為顯著增加,即在突變點處FVC顯著減小,也符合Google Earth time-series images和FVC圖像的比對特征。因此,通過圖像比對驗證,發現用BFAST模型進行長時間序列植被非線性趨勢檢測是適用的,結果也具有一定可靠性,研究可為相關生態環境保護措施的制定提供科學參考依據。

4.3 怒江-薩爾溫江流域植被覆蓋驅動力

為了更有針對性地分析植被覆蓋的影響因素,考慮研究區地形地貌的差異性,研究選擇了適用于具有非線性特征分析的基于最優參數的地理探測器(OPGD),對研究區植被覆蓋影響因素進行分區探討,根據探測結果進行更有針對性的原因分析。探測結果表明,對整個區域來說,海拔的影響是最大的,原因在于研究區海拔南北差異極大,對整個區域的地形地貌起著宏觀控制作用,進而影響局部地區的水熱條件,對植被空間分異產生重要影響。此外,影響較大的是氣溫、降水、蒸散發等氣象因子,這與許多研究以及現實情況也是相符的[51-53],植被受氣候影響較大,特別是怒江-薩爾溫江流域跨越了多個不同氣候區,氣候差異大,因此氣候因子對其影響較顯著;人為因子對全區植被的影響程度比自然因子較低,其中以土地利用影響最大,原因在于研究區地形偏遠,人口經濟發展較緩慢,特別是高海拔地區和高山峽谷區,人類活動很難影響到,而土地利用作為能較好反映人類活動的因素,能很好的表現植被覆蓋的變化,同時也有很大的可利用空間來增強或者減弱植被覆蓋。

對各子流域來說,植被覆蓋的影響因素具有差異性。上游區域土地利用對植被覆蓋的影響大于海拔,原因在于區域內部大部分屬于高原,高原內部地勢起伏不大,因此海拔以及坡度、坡向等對植被的影響小于土地利用,人類活動通過改變土地利用結構,特別是耕地面積、耕作類型等從而影響植被覆蓋;在剩下的因子中降水是影響較大的,原因在于上游地區季風難以到達,降水又是植被生長的重要因素,因此,降水對植被影響較大。由于其較為偏遠,人為因素影響小于自然因素。中游區域,由于其大部分屬于高山峽谷區,海拔差異顯著,海拔成為了植被覆蓋影響最大的因素,而氣溫、降水、蒸散發等受海拔差異的影響,對植被覆蓋影響也較大。下游區域,由于以平原為主,地形較為平坦,人類活動頻繁,因此人口和GDP是對植被影響最大的因素,由于緯度低,氣溫高,蒸散發也成為了植被變化的重要影響因素,與上游和中游地區不一樣的是海拔對植被影響微弱。

綜上可知,由于怒江-薩爾溫江流域地理特征的復雜性,在不同環境條件下,各環境因子對植被覆蓋的影響也存在顯著的差異性。因此,建議綜合考慮不同環境條件下FVC的空間分布與各影響因子的空間關聯性,明確不同區域的環境制約因子,因地制宜制定生態修復措施,規劃流域的經濟建設。對于上游地區來說,植被覆蓋較低,未來植被呈退化趨勢的占比較多為9.03%,應防止植被覆蓋率進一步下降加強生態恢復,如優化土地利用結構、發展現代生態農業、更改耕作類型、進行輪作和休耕管理、退耕還林、建設自然保護區和實施生態移民、以及爭取國家政策的扶持等,注重社會經濟和生態環境的共同發展;對于中游地區,自然環境因子對植被覆蓋影響更大,其經過三江源等重要生態源地區,應該確保生態環境不被破壞;對于下游地區,植被覆蓋度較高,但影響植被覆蓋的主要是社會經濟因素,且未來植被呈退化趨勢的占比是最多的為17.39%,應切實保護好現有資源,維持生態系統的穩定性,并發展高效、多元化、可持續的生態經濟。植被動態是一個復雜的過程,還應加強對極端事件的預防和控制,以防止氣候變化造成植被的進一步退化。總之,應充分發揮人類活動對于生態環境的積極作用,找準生態保護與經濟發展的平衡點,促進生態與經濟的協調可持續發展。

4.4 局限和展望

研究參考楊帆等[30]的成果重新將怒江-薩爾溫江流域劃分為了上游、中游和下游區域進行植被趨勢分析與驅動力探討,我們主要考慮的是海拔和行政區劃,區域劃分較為主觀,這是未來研究需要改善的地方,比如利用具體的氣候區、土壤區、地貌區、植被區等進行分析。此外,研究在進行驅動力指標選取的時候,由于跨境各個國家的數據統計標準不一,精度也有差異,使得部分數據獲取較為困難,特別是人為因素的數量和質量,這可能會影響結果的精度,因此,未來應該考慮在分區探討的同時對不同區域建立不同的指標體系,以提高數據的精度和結果的準確性。

5 結論

研究基于怒江-薩爾溫江流域2000-2021年MODIS NDVI數據,計算出其植被覆蓋度FVC,探討了研究區近22年植被覆蓋的時空趨勢變化、非線性特征、未來可持續性以及驅動力。主要結論為:

(1)近22年,怒江-薩爾溫江流域植被覆蓋總體發展趨勢向好,2007年和2015年出現較大突變,多年平均FVC值為0.73。受水熱條件、地形地貌等影響,植被分布具有明顯的空間異質性,下游和中游區域植被覆蓋明顯優于上游區域,總體表現為南高北低,植被覆蓋狀況較好。(2)BFAST趨勢表明,近22年怒江-薩爾溫江流域植被覆蓋改善和退化的區域面積占比分別為71.24%、28.76%,改善的區域遠大于退化的區域,說明研究區植被得到較好的保護。Hurst指數顯示,未來植被相比以前將得到改善和將有所退化的區域占比分別為94.89%、2.76%。BFAST與Hurst二者的耦合表達了植被覆蓋的未來非線性趨勢,也是改善的區域大于退化的區域,占比分別為68.99%和29.09%。應特別關注未來發展趨勢退化的區域。(3)地理探測器結果表明各區域植被覆蓋影響因素具有差異性,針對全區域海拔的影響是最大的,進一步說明了地形與研究區植被具有密切關系。而針對各子流域,土地利用對上游地區植被覆蓋的影響比中游和下游區域顯著,中游地區海拔的影響力最大,下游地區以人口、GDP等人為因素影響為主。因此,要合理根據研究區不同區域的環境制約因子等制定合理的植被保護措施。