基于GEE的祁連山自然保護區景觀格局變化分析

劉寬梅

（蘭州大學 資源環境學院，甘肅 蘭州730000）

景觀動態變化研究具體表現為景觀結構、景觀功能、景觀空間格局隨時間的變化，是景觀生態學的核心問題[1-2]。景觀格局是指在空間上大小、形狀不同、屬性有別的景觀斑塊在景觀空間上的排列[3]。研究景觀格局的動態過程，能掌握其動態變化規律，進一步準確認識研究區面臨的各種生態問題。隨著國內外景觀格局研究成果不斷豐富，大致可以將景觀格局分析方法分為三大類：（1）景觀格局指數法，是由文獻[4]中Turner和Gardner、Li和Reynolds等研究發展的反映景觀結構與空間配置特征的量化指標，能高度濃縮景觀格局信息。Fragstats等軟件的出現讓指數計算更便利，如Yang C Y、Wang J J及Zhuang G S等均通過Fragstats計算相應的景觀格局指數，對不同地區的景觀格局進行分析[5-7]；Song H H等通過景觀指數分析發現景觀格局的不足從而優化景觀格局[8]。（2）空間統計學方法，包含了地統計學方法、空間自相關及空間插值法等[9]。如黃玲[10]、肖飛艷[11]、潘竟虎[12]、Pan D Y[9]等利用轉移矩陣法分別對不同地區景觀格局變化進行了探討。（3）景觀格局模型與模擬法，便于研究者對數據隱含的內在規律進行探索[13-14]。張小松[15]、廖芳均[16]、羅攀[17]、陳頤[18]、羅紫薇[19]等基于CA-Markov模型分別對城市、自然保護區等景觀格局進行預測。從研究內容來看，主要有景觀格局的時空變化及其驅動力因素[20]、景觀空間格局的粒度效應[21-22]、景觀格局的演化過程與特征及景觀格局變化的尺度效應等方面的內容[23-24]。

祁連山作為水源涵養林和生物多樣性保護熱點區，其生態環境問題受到廣泛關注。目前，關于景觀格局的研究多在祁連山東段開展。草地為該區的景觀基質，其優勢度高，連接性好[25-28]。斑塊水平格局特征為：草地形狀結構趨于松散且簡單化，森林和冰雪形狀結構趨于緊密但簡單化，灌木形狀結構趨于松散但復雜化，裸地和水域形狀結構趨于緊密但簡單化[27，29]。景觀水平格局特征為：景觀整體多樣性水平增加，破碎化程度增大[28-29]，景觀均勻度增大，景觀優勢度、分離度減小[27]。對于景觀整體異質性的研究，不同學者研究時段及選取指標不同，會導致結果有所差異，劉晶[27]和劉宇[30]等對祁連山東段1994-2008年的景觀格局進行探討，劉晶等認為該區景觀異質性增大，景觀向著均勻化方向發展；而劉宇等認為該區景觀異質性呈下降趨勢；徐亞男等[29]發現2000-2016年間景觀異質性是增大的。人類活動是影響該區景觀格局變化、生態敏感性變化和景觀脆弱性變化的決定性因素[26-29]。而對于祁連山其他區域，不同地區景觀格局變化不同。祁連山地區最主要的土地利用類型是草地和未利用地[31]；祁連山南坡基質由裸地向草地轉變，景觀破碎化水平較高，放牧活動導致草地和灌叢破碎化，水域形狀差異大，沼澤分布集中成片，居民區和工礦用地分布破碎[32-33]；中祁連山淺山區草地面積最大，受海拔和坡向的影響，各景觀組分空間分布不均衡[34]，喬木林、灌木林和農田面積在增加，疏林、草地和裸地面積減少[35]；祁連山三河源地區景觀格局向多樣化、破碎化趨勢演進，人類活動對景觀格局的干擾在增強[36]；石羊河上游祁連山區草地、工礦居民用地、裸土和戈壁的破碎度增大，耕地和林地的破碎度減小，景觀多樣性略有上升，單一性略微下降[37]；祁連山甘肅區域景觀破碎化程度較高，Ⅰ級生態景觀格局比較穩定，Ⅱ級生態景觀格局受自然和社會環境等因素影響變化較大[38]。

祁連山地區云量分布有明顯的西部向東部逐漸減少的趨勢，各季節總云量空間分布也有明顯的變化，遙感影像數據質量會影響分類結果的好壞，已有研究多基于RS和GIS平臺獲取景觀分類結果，其會受計算機空間存儲和處理環境配置限制，還會受單一數據質量限制，而谷歌地球引擎（Google Earth Engine，簡稱GEE）云平臺利用多個CPU同時對遙感數據進行運算，可實現數據快速處理，而且可以提取研究時段內多景遙感影像合成云量最小的影像，克服受云量影響難以獲得高質量數據的缺陷，提高數據可用性。因此，基于GEE平臺解譯保護區近30年的景觀類型圖，并根據斑塊水平和景觀水平選取景觀格局指數對該區域內景觀格局變化進行分析，以期為保護區管護、生態規劃提供參考。

1 研究區概況

祁連山自然保護區（簡稱保護區）位于北緯36°29′-39°43′，東經96°46′-103°45′之間，位于祁連山中、東段區域范圍內，面積約282.19萬公頃（包括甘肅祁連山國家級自然保護區198.72萬公頃和青海祁連山省級自然保護區83.47萬公頃）。祁連山東西長達1200km，南北寬50-120km，海拔3000-5000m。祁連山深居內陸，具有典型的高原大陸性氣候特征，全年降水量主要集中在5-9月。祁連山氣候要素的垂直地域分異明顯，一般山前低山是荒漠氣候，年降水量約150mm，年均溫度6℃左右；山體中下部屬于半干旱草原氣候，年降水量250-300mm，年均溫2-5℃；山體中上部為半濕潤森林草原氣候，年降水量400-500mm，年均溫0-1℃；山體上部為亞高山和高山冷濕氣候，年降水量約800mm，年均溫-5℃左右。研究區位置示意圖如圖1所示。

圖1 研究區位置示意圖

2 數據來源與方法

2.1 數據來源與處理

（1）研究區邊界：青海祁連山省級自然保護區邊界是基于《青海祁連山自然保護區科學考察集》[39]中的邊界圖，對其進行地理配準矢量化得到，甘肅祁連山國家級自然保護區邊界來源于國家地球系統科學數據中心（http：//www.geodata.cn）。

（2）景觀分類標準體系：由于景觀是具有特定結構、功能和動態特征的不同生態系統鑲嵌組成的異質性區域[40]，因此分類體系以《全國生態系統分類體系》為基礎，參考歐陽志云等人[41]的《基于遙感技術的全國生態系統分類體系》，結合保護區景觀特征，將研究區的景觀類型劃分為八類，分別為濕地、高寒荒漠、冰川/永久積雪、城鎮、農田、森林、草地、灌叢。

（3）遙感數據及預處理：基于GEE平臺調用研究區1989、1995、2000、2005、2010、2015、2019年每年4-10月已經過預處理的Landsat地表反射率產品（Surface Reflectance，SR）（https：//code.earthengine.google.com）；利用多景影像重構及合成目標年份最小云量影像集；計算歸一化建筑指數、歸一化水體指數及歸一化植被指數作為進行隨機森林分類時的特征變量以提高分類精度；結合GEE與Google Earth平臺，采用人工目視判讀得到8類景觀類型各年的總樣本數分別為1193、1092、1100、1023、1101、1110、1161個，其中70%用于分類器訓練，30%用于結果驗證。

2.2 研究方法

2.2.1 基于GEE的隨機森林景觀類型分類

隨機森林分類方法是由多棵決策樹構成的機器學習算法[42]，該方法已被證明比其他方法具有更高的分類精度[43]。設置決策樹數量為100對訓練樣本集進行訓練，并提取每一目標年份的景觀類型信息，以總體精度、Kappa系數作為衡量景觀分類精度的指標，具體各指標計算公式參見相關文獻[44-45]，七個時期的總體精度分別為：90.8%、91.21%、91.78、85.81%、88.90%、83.71%及89.66%，其Kappa系數分別為0.88、0.88、0.87、0.82、0.85、0.79及0.87。

2.2.2 景觀格局指數

景觀格局指數能夠高度濃縮景觀格局信息[4]，而且其計算簡單，現被廣泛使用。根據以往的研究來看，指數間有很高的相關性，若是引入同一類型指數太多，容易造成信息冗余[46]。因此，采用基于統計的景觀指數選擇方法，該方法能解決個人選擇指數的主觀性，降低選取大量景觀指數進行計算的工作量，具有一定的合理性[17]。以“自然保護區景觀格局”為主題檢索詞，在知網上檢索2015－2020年的論文數據（共計76篇）；然后統計論文中出現的景觀指數，最后提取到使頻次為10次及以上的景觀指數17個，對入選指數進行皮爾遜相關分析，最后綜合考慮指數生態學意義、不重復性、常用性及研究需求，斑塊類型水平選擇4類特征指標共8個：（1）優勢度指標：最大斑塊面積指數（LPI）、平均斑塊面積（MPS）；（2）破碎化指標：斑塊密度（PD）、邊緣密度（ED）；（3）形狀特征指標：景觀形狀指數（LSI）、分維數（PAFRAC）；（4）聚集度指標：散布與并列指數（IJI）、聚合度指數（AI）。景觀水平選擇4類特征指標共6個：（1）聚集度指標：聚合度指數（AI）、散布與并列指數（IJI）；（2）優勢度指標：最大斑塊面積指數（LPI）；（3）形狀特征指標：景觀形狀指數（LSI）；（4）景觀多樣性指標：香農多樣性指數（SHDI）、香農均勻度指數（SHEI）。指標含義見表1，指數計算采用Fragstats4.2軟件完成，計算方法參考Fragstats4.2軟件幫助文件（http：//www.umass.edu/landeco/research/fragstats/fragstats.html）。

表1 景觀格局指數意義

3 結果與分析

3.1 景觀類型結構變化

由圖2可知，保護區景觀類型以高寒荒漠和草地為主。1989-2019年間，高寒荒漠面積總體上呈減少趨勢，凈減少面積為3130.45km2，其優勢地位慢慢被草地代替，說明荒漠化環境有所改善；草地面積在各研究時段內波動變化，但總體是上升的，凈增加面積為4223.62km2，說明草地景觀比較脆弱，極易受放牧活動等人類活動的影響；森林面積經歷先減少后增加再減少的趨勢，總體上面積是減少的，凈減少面積為261.09km2，灌叢面積變化呈現先增加后減少的趨勢，總體上面積是減少的，凈減少面積為242.47km2，可以看出植被結構遭到破壞，出現退化現象；冰川/永久積雪面積呈現緩慢減少趨勢，凈減少面積為662.18km2；濕地、城鎮與農田面積占比均較小，面積變化不明顯。

圖2 1989-2019年祁連山自然保護區景觀類型面積

3.2 基于景觀指數的景觀格局分析

3.2.1 類型水平的景觀指數分析

（1）優勢度分析

由表2可知，最大斑塊面積和平均斑塊面積最大值由高寒荒漠和草地主導，2019年之前，高寒荒漠是保護區的基質景觀，但在1989-2019年間，由于封沙育林育草、黑土灘治理、荒漠化治理等環境保護政策和手段的實施，使得高寒荒漠被草地等景觀分割的程度逐漸增加，其最大斑塊面積指數和平均斑塊面積值不斷減小，破碎化程度增大，優勢度降低，對景觀的控制作用正在削弱，而草地的最大斑塊面積指數和平均斑塊面積不斷增加，優勢度則不斷增強，逐漸代替高寒荒漠成為新的優勢景觀，說明荒漠化環境得到改善。其他景觀類型的最大斑塊面積指數值變化相對較小，結合其平均斑塊面積變化，可以發現冰川/永久積雪和森林的優勢度次于高寒荒漠和草地，灌叢和濕地的優勢度較低，城鎮和農田的優勢度最低。

表2 最大斑塊面積指數和平均斑塊面積

（2）破碎度分析

由表3可知，高寒荒漠和草地的斑塊密度值和邊緣密度值比其他景觀類型的值高，說明要素斑塊與其相鄰異質斑塊間的接觸多，邊緣復雜，破碎化程度較其他景觀類型高。1989-2019年間，高寒荒漠和草地的斑塊密度值波動幅度小，但高寒荒漠的邊緣密度值呈下降趨勢，說明其邊緣由復雜向簡單化演變，景觀的破碎化一定程度上在減輕，草地的邊緣密度值呈上升趨勢，其邊緣復雜化，破碎度增加，道路交通及放牧等人類活動導致其破碎度增加；灌叢的邊緣密度值在不斷下降，斑塊密度值先上升后下降，說明該景觀的破碎化程度在由高向低轉變；冰川/永久積雪、農田的邊緣密度值在不斷下降，但冰川/永久積雪的斑塊密度值變化不大，農田的斑塊密度呈下降趨勢，這說明這兩種景觀要素的破碎化程度好轉；森林的斑塊密度和邊緣密度值的波動趨勢一致，均呈“倒N型”分布，破碎化程度好轉；濕地和城鎮的斑塊密度和邊緣密度值均較小且波動程度也小，其破碎化程度較低，這與城鎮呈團聚式發展有關。

表3 斑塊密度和邊緣密度指數

（3）形狀特征分析

由表4可知，1989-2019年間所有景觀類型的形狀指數都較大，說明研究區內各景觀斑塊形狀相對復雜。結合表5可看出，各類景觀類型的形狀指數和周長面積分維數變化規律基本一致。具體來看：草地的形狀指數值和周長面積分維數指數值均較大且變化平穩，說明草地斑塊形狀復雜度較大且較穩定；灌叢、農田、冰川/永久積雪的形狀指數和周長面積分維數值均呈下降趨勢，說明這三類景觀類型的形狀在趨于簡單化；森林和高寒荒漠的形狀指數在增加，周長面積分維數基本不變，說明其形狀在趨于復雜；濕地和城鎮的斑塊形狀指數和周長面積分維數均呈波動變化狀態，但形狀指數是除冰川/永久積雪外最小的，說明這兩類景觀在所有景觀中的形狀最簡單。

表4 景觀形狀指數與周長面積分維數

（4）聚集度分析

從表5可以看出，聚合度指數最大的是高寒荒漠，其次為冰川/永久積雪和草地，城鎮和農田最低。說明高寒荒漠、草地和冰川/永久積雪聚集程度高，趨向連片分布，前兩者的聚集度高與其面積占優勢有關，而冰川/永久積雪則往往在高海拔地區成片分布，而城鎮和農田趨于分散分布，這與保護區內城鎮分布實際情況相符。在1989-2019年間，高寒荒漠、冰川/永久積雪、草地和森林的聚合度指數基本沒發生變化，說明這四類景觀各自聚集度高且穩定；濕地和城鎮的斑塊聚集度呈增加趨勢；而灌叢和農田的聚合度指數呈波動變化，灌叢的斑塊分布由聚集走向分散再走向聚集，2010年以后趨于平緩，而農田的聚集度則是波動下降的。高寒荒漠和草地的散布與并列指數值比較大，說明這兩者與其他景觀斑塊間的連通度較大，這與高寒荒漠和草地在景觀中的優勢地位有關，但在1989-2019年間，它們與其他斑塊的連通度都呈先上升后下降的趨勢；灌叢和農田與其他的景觀斑塊類型的連通度也是下降的；濕地和城鎮與其他斑塊的鑲嵌程度不斷增加；而冰川/永久積雪的散布與并列指數盡管有上升趨勢，但仍然是最低的，說明它與其他斑塊的連通度最小，這跟該景觀類型分布在高海拔地區有關，除高寒荒漠外，其余景觀類型極少分布在冰川/永久積雪所在的高海拔位置。

表5 聚合度指數和散布與并列指數

3.2.2 景觀水平的景觀指數分析

（1）聚集度分析

由圖3可知，散布與并列指數表現出逐漸下降的趨勢，說明各種景觀類型之間空間鄰接分布的復雜程度有所減小，彼此之間相鄰的景觀要素類型由多變少，空間關系趨于簡單化，景觀的連通性變弱。通過聚合度指數變化趨勢可以看出，該值在研究時段內聚合度指數值均較大（均大于88），且有緩慢增長趨勢，說明同類斑塊趨于集中分布，景觀破碎化減小，景觀間的空間異質性逐漸縮小。

圖3 景觀水平上的散布與并列指數和聚合度指數

（2）景觀優勢度及形狀特征分析

由圖4可知，最大斑塊面積指數呈現出先增加再減小再增加的變化趨勢，進一步驗證高寒荒漠景觀優勢度由強變弱，草地景觀優勢度由弱變強的發展過程。1989-1995年的景觀形狀復雜度減小，1995-2005年景觀復雜度呈增加趨勢，2005年以后景觀斑塊形狀在趨于簡單化發展，說明景觀形狀變化較頻繁，人類活動對景觀格局的影響呈現出有目的、有規劃改造的趨勢。

圖4 景觀水平上的最大斑塊面積指數和景觀形狀指數

（3）景觀多樣性分析

由圖5可知，香農多樣性指數和香農均勻度指數值變化趨勢一致，1995-2010年間香農多樣性指數和香農均勻度指數值呈緩慢增長趨勢，表明該時段內景觀多樣性呈增長趨勢，且各景觀類型分布趨向均勻化發展，1989-1995年和2010-2019年間兩個指數為下降趨勢。但從總體來看，香農多樣性指數值在1.14-1.24之間，反映出研究區景觀類型多樣性低，所有景觀類型在研究區內呈現不均衡化分布趨勢，香農均勻度指數值為0.55-0.60之間，說明各景觀類型所占比例差異較大，景觀受優勢景觀支配，進一步驗證了本保護區內高寒荒漠、草地景觀占主導地位的現狀。

圖5 景觀水平上的香農多樣性指數和香農均勻度指數

4 結論

（1）保護區景觀類型以高寒荒漠和草地為主。1989-2019年間，高寒荒漠面積總體上呈減少趨勢，凈減少面積為3130.45km2，其優勢地位慢慢被草地代替，荒漠化環境有所改善；草地面積總體上是波動上升的，凈增加面積為4223.62km2，說明草地景觀比較脆弱，極易受人類活動尤其是放牧活動的影響；森林面積總體上是減少的，凈減少面積為261.09km2，灌叢面積總體上是減少的，凈減少面積為242.47km2，可以看出植被結構遭到破壞，出現退化現象；冰川/永久積雪面積呈現緩慢減少趨勢，凈減少面積為662.18km2；濕地、城鎮與農田面積占比均較小，面積變化不明顯。

（2）斑塊類型水平景觀格局變化：30年間，優勢景觀由高寒荒漠向草地過渡；草地的斑塊密度值和邊緣密度值在各年都是最大的，其破碎度最大且呈增加趨勢；形狀指數值和周長面積分維數值最大值由草地和灌叢主導，兩者的斑塊形狀復雜度較大，草地的斑塊形狀復雜度穩定，而灌叢的斑塊形狀由復雜趨向簡單化發展；高寒荒漠和草地在景觀中的優勢地位使得高寒荒漠和草地的聚合度和連通度均比其他景觀高，兩者的聚合度變化不大，而連通度都呈先上升后下降的趨勢，冰川/永久積雪因為其分布海拔高且成片分布，致使其聚合度跟草地相當但連通度低。

（3）景觀水平格局變化：景觀整體間的連通度在減小，同類景觀趨于集中分布，景觀破碎度在減小；景觀斑塊形狀總體在趨向簡單化發展；景觀類型多樣性低，景觀類型呈現不均衡化分布趨勢，景觀受優勢景觀所支配。

5 討論

有研究表明，草地是祁連山大部分區域內的基質[25-27]，本文中高寒荒漠占比最大的原因是由于本研究區域青海境內大片高寒荒漠分布，這在湯萃文[32]等對祁連山南坡植被景觀格局的研究中得到印證，而農田、城鎮及水體小面積零散分布與的王清忠[35]的研究結果相同。劉晶[27]，蔣剛[33]等學者對祁連山不同地區的研究都表明，人類活動是景觀頻繁變化的動因，甚至使某些區域的植被類型出現嚴重退化現象，不同地區受到的人類活動影響不同，各研究學者選取的景觀格局指數及研究時段也不同，所以對景觀格局變化結果也會有差異。如辛雨瓊[26]等發現東祁連山農牧交錯區1988-2018年景觀呈復雜化和破碎化趨勢，劉晶[28]等發現祁連山東段1994-2008年間景觀多樣性、均勻度、破碎度指數增大，景觀優勢度、分離度減小，但是徐亞男[29]等發現祁連山東段2000-2008年景觀整體多樣性水平降低，異質性和破碎化程度減小，各景觀組分分配不均勻且朝著聚集的趨勢發展，2008-2016年，景觀整體多樣性水平增加，異質性和破碎化程度增大，各景觀組分向著均勻化且分散的趨勢發展。綜上可以看出，即使是同一個大區域，也可能會因為指標選取及小區域內人類活動影響程度不同等因素導致景觀格局結果有差異，因此祁連山這樣的山地景觀格局指標選取還有待進一步探討。此外，研究區內植被存在抵御能力弱及退化現象，仍需對其進行加強管護。