基于地形起伏度的冰湖潰決隱患研究
——以希夏邦馬峰東部為例

賀鵬， 童立強， 郭兆成， 涂杰楠， 王根厚

(1.中國地質大學(北京)地球科學與資源學院，北京 100083；2.中國自然資源航空物探遙感中心，北京 100083)

0 引言

在全球氣候變暖背景下，近年來青藏高原暖濕化特征顯著，冰川活動性增強，冰湖潰決事件頻發[1]。西藏地區是我國冰湖潰決災害高發區，以喜馬拉雅山區尤甚[2]。自20世紀30年代以來，該區有文獻記載的潰決冰湖33處，冰湖潰決事件37次，對我國西藏南亞大通道和川藏鐵路等重要交通干線建設構成嚴重威脅[3-6]。冰湖潰決的形成機制復雜，誘發成因多樣，其中冰崩、冰滑坡及相關組合誘因占70%[7]。冰崩、冰滑坡等災害多發于斜坡處，由冰川斷裂崩解形成[8]。目前，冰湖潰決危險性判別及評價流程已形成較為完整的體系[9]，其中斜坡地質災害敏感性是定量評價過程中的重要指標之一[10]，而地形起伏度又是實現準確評價的前提[11]。

作為描述區域地形特征的宏觀指標，地形起伏度是指某一確定面積中最高點與最低點海拔的差值[12]，可以真實反映地表起伏特征變化，與斜坡地質災害，特別是冰崩災害的發生顯著相關[13]。近年來，基于地貌發育的基本理論，前人在不同地區和尺度下展開研究： 涂漢明等[12]利用國家數字地形數據庫，基于最大高差法得出中國地形起伏度最佳統計單元為21 km2； 劉愛利等[14]基于中國1∶100萬數字高程模型(digital elevation model,DEM)，探討DEM定量化、自動化的地貌分類應用，得出在1∶100萬比例尺下最佳統計單元為6.5 km2； 韓海輝等[15]應用均值變點法提取青藏高原地區地形起伏特征并分級分析，認為最佳統計窗口為1.17 km2； 張軍等[16]以新疆1∶25萬DEM數據為基礎，利用鄰域分析及均值變點法得出2.56 km2的窗口尺度為最佳統計單元； 畢曉玲等[17]基于人工判讀法得出四川省地形起伏度的最佳統計單元為2.25 km2，并進一步探討了其在區域地質災害評價中的應用； 千琳勃等[18]利用均值變點法獲取了六盤山區地形起伏度，確認最佳統計單元為1.96 km2，并分析了研究區地形空間變化特征； 張競等[19]針對京津冀地區不同地貌樣本計算了最佳分析窗口，認為地貌起伏在400 m內以4.64 km2為佳，而起伏度大于400 m的以5.35 km2更為適合； 王巖等[20]運用均值變點法，在青海德貴地區開展1∶5萬比例尺下地形起伏度研究，論證了最佳分析窗口為0.16 km2。以往研究中以中低尺度下全國及西部地區的居多，以均值變點法的應用頻率最高、效果最好[21]。然而，在我國藏南冰湖集中區1∶5萬比例尺下的研究較少，將地形起伏度應用于冰湖潰決隱患識別方面的分析案例還較少。在當前冰湖潰決災害頻發的嚴峻形勢下，亟須開展研究，以豐富并完善我國冰湖潰決隱患分析評價方法。本文以希夏邦馬峰東部冰湖集中區為研究區，以國產高分二號(GF-2)衛星影像及ASTER GDEM V3為數據源，利用均值變點法確定地形起伏度最佳提取單元，結合冰湖分布特征開展潰決隱患相關性分析。

1 研究區概況

研究區位于聶拉木縣以北，地處喜馬拉雅山中段希夏邦馬峰東南部區域，在行政區劃上隸屬西藏自治區日喀則市聶拉木縣，范圍涵蓋波絨鄉、門布鄉、亞來鄉及聶拉木鎮等4個鄉鎮，總面積約730 km2(圖1)。地理坐標介于E28.17°～28.5°，N85.75°～86.00°之間。

圖1 研究區遙感影像

受青藏高原新構造運動影響，形成整體以高落差基巖山地夾持縱深河谷的地貌特征，區內冰湖分布廣泛、重力地貌發育、各類地質災害頻發[22]。受近年暖濕化氣候影響，喜馬拉雅山區冰川退縮明顯，冰湖顯著擴張，潰決風險加劇，潰后引發鏈式次生災害事件也隨之增多[23]。據前人研究，希夏邦馬峰東部及周邊地區存在多個高潰決風險冰湖，其中嘉龍錯、次仁瑪錯和帕曲錯等已發生潰決，對下游居民的生命財產安全造成巨大損失[24-25]。

2 研究數據與方法

2.1 冰湖信息提取

采用人機交互解譯方法，基于國產GF-2衛星正射級產品數據，并配合圖像增強提取處理等方法，實現對區內冰湖信息的準確獲取。數據由中國自然資源航空物探遙感中心提供。

2.2 地形起伏度計算

2.2.1 地形數據

本文采用ASTER GDEM V3版地形數據。該數據于2019年8月發布，覆蓋了從N83°～S83°的所有陸地表面，其全球空間分辨率為30 m，垂直分辨率約為10 m，較V2版主要減少了高程值空白區域及水域數值異常[26]。數據來源于美國國家航空航天局(https: //earthdata.nasa.gov/)。

2.2.2 均值變點法

地形起伏度一般基于DEM數據通過窗口分析法求取，其關鍵是確定最佳統計窗口大小。以某一起伏地形上任意點為中心計算地形起伏度，隨著統計窗口尺度的擴大，相對高程差與地形起伏度的變化率呈對數曲線特征，其中必定存在一個變化率由大變小(即曲線由陡轉緩)的“拐點”，該點有且唯一，對應的統計窗口大小即適宜計算尺度。本文通過均值變點法確定研究區最佳統計單元，具有客觀性，相關計算原理參見相關文獻[15-20]，在此不再贅述。

3 適宜計算尺度確定

3.1 主要步驟

以n×n的矩形作為分析窗口，對希夏邦馬峰東部冰湖集中區進行窗口尺度遞增的起伏度計算，規定起止矩形窗口尺度分別為3×3個像元和65×65個像元，單次增幅為2個像元(表1)。

表1 地形起伏度分析窗口尺度設置

基于ArcGIS軟件的Spatial Analyst模塊下neighborhood statistics工具對區內DEM柵格數據進行不同尺度窗口下的遍歷計算，得到各窗口的柵格最大及最小值； 利用Raster Calculator工具計算起伏度值，并將分析窗口大小與該尺度下起伏度值進行關聯(表1)； 運用Band Collection Statistics工具進行柵格數據統計分析，并將結果屬性表輸出為表格； 利用VBA編程采用均值變點法計算統計量，其計算公式為：

(1)

(2)

將表1中各項數值代入式(1)—(2)中，實現對窗口面積與起伏度的曲線擬合與均值變點法中S-Sj的計算。變點的存在會使S和Sj的差距增大，S-Sj的最大值對應的分析窗口大小即為最佳分析窗口。

3.2 地形起伏度擬合曲線

前人研究表明地勢起伏度與分析窗口大小的變化曲線為邏輯斯蒂曲線[12]，因此，對地形起伏度和統計單元面積(分析窗口)進行曲線擬合，如圖2所示。

圖2 窗口大小與地形起伏度擬合曲線

在圖2的擬合曲線中，決定系數R2為0.950 5，擬合度良好，通過統計學檢驗。地形起伏度隨統計單元尺度的擴大而增加，變化曲線表現為對數函數或冪函數特征。統計單元面積小于32.49×104m2(窗口尺度19×19)時，其增長趨勢迅速； 統計單元面積在32.49×104～47.61×104m2區間(窗口尺度在19×19～23×23)時，地形起伏度的增速由急轉緩； 而統計單元面積超過47.61×104m2(窗口尺度23×23)之后，地形起伏度的增速放緩并趨于平穩。

3.3 最佳窗口尺度確定

如圖2所示，擬合曲線內存在一個增速逐漸趨緩的點(非數學拐點)，即所謂“變點”，該點對應的統計面積即為最佳分析窗口尺度。將起伏度計算結果逐一代入統計量公式，得到非線性系統數列樣本j的S與S-Sj值，其中S-Sj最大值所對應的分析窗口大小，即所求的適宜計算尺度，計算結果見表2。

表2 均值變點法分析統計結果

根據表2數據，構建S-Sj與非線性系統數列樣本j的變化擬合曲線(圖3)，可以看到該曲線為近似拋物線型，在j=10時出現S-Sj的最大值，可知此處為所求的變點，即最佳分析窗口的尺度，其窗口尺度為21×21，統計單元面積為39.69×104m2。

圖3 S-Sj值的變化趨勢

4 地形起伏度對冰湖潰決的影響分析

4.1 基于遙感的冰湖提取結果

基于GF-2衛星數據的冰湖信息提取結果顯示(圖4)，區內現存冰湖多達1 020個，總面積17.67 km2。從規模上來看，可分為大型(≥1 km2)、中型(0.1～1 km2)及小型(≤0.1 km2)3類； 就類型而言，包括冰蝕湖、冰磧阻塞湖、冰斗湖及冰川槽谷湖4類。區內冰蝕湖僅3處，分別是位于希夏邦馬峰東部地區的崗西錯、貢錯及嘎龍錯，其總面積達12.21 km2，是全區僅有的3個大型冰湖； 冰磧阻塞湖18處，主要是分布在希夏邦馬峰東南部地區的嘉龍錯和達熱措等冰湖，其中中型冰湖8處、小型冰湖10處，總面積為3.12 km2； 冰斗湖30處，均為小型冰湖，以中南部地區分布居多，總面積0.41 km2； 冰川槽谷湖969處，規模上除2處中型冰湖外其余均為小型，多分布于冰川侵蝕作用形成的谷地內。

圖4 冰湖提取結果及分類

4.2 地形起伏度分級

依據上述統計分析結果，采用21×21的網格大小(39.69×104m2)對研究區內地形起伏度進行計算。為便于后續分析，進一步對區內起伏度進行分級，具體分為6個等級： 起伏度在[0,30) m為平原，[30,70) m為臺地，[70,200) m為丘陵，[200,500) m為小起伏山地,[500,1 000) m為中起伏山地； [1 000,2 500) m為大起伏山地。

結果顯示(圖5)，區內平原分布相對較少，整體以小起伏山地為主，在拉布日峰至俄熱村之間、希夏邦馬峰東南側多有分布； 丘陵及中起伏山地所占比例相近，前者多見于北部納朵熱至大來村、崗西錯東部等地，后者主要分布于中部及南部冰湖集中分布區域； 大起伏山地相對較少，主要位于中西部的希夏邦馬峰及西南部的多爾雷山附近。

圖5 研究區地形起伏度分級

4.3 地形起伏度對冰湖潰決的影響分析

如前文所述，地形起伏度是定量化分析斜坡地質災害敏感性的重要因子之一，更是冰湖潰決危險性評價體系中不可或缺的指標，對研究全區冰湖潰決隱患具有重要的指示意義。本文以冰湖周邊地形起伏組合特征為依據，疊加冰湖空間分布特征，以分析區內地形起伏度對冰湖潰決的影響。

4.3.1 地形起伏度與冰湖潰決的關系

通過地形起伏度與冰湖空間分布特征的疊加分析發現就冰湖后緣及兩側起伏度特征而言，主要可分為單一起伏等級型、多起伏等級組合型及跨起伏等級型3類冰湖(以下簡稱單一型、組合型和跨級型)。

單一型冰湖周緣起伏度類型一致，多發育于冰川槽谷區，湖泊類型上以小型冰磧湖為主。區內地形起伏度較小、地勢變化緩、穩定性高，不具備潰決條件。從影像特征及實地驗證結果看，沖堆普流域(希夏邦馬峰與多爾雷山之間)的典型冰湖(圖6(a))特征與分析結果一致，作為冰川退縮消融過程中的殘留湖泊，多處于同一地形起伏等級(圖6(b))，周邊地勢平緩(圖6(c))，冰湖間多由常年地表徑流相連，由于湖面較淺、流動性好，整體潰決隱患極低。

(a) 單一型冰湖影像 (b) 單一型冰湖起伏特征 (c) 單一型冰湖實地調查照片

(d) 組合型冰湖影像 (e) 組合型冰湖起伏特征 (f) 組合型冰湖實地調查照片

(g) 跨級型冰湖影像 (h) 跨級型冰湖起伏特征 (i) 跨級型冰湖實地調查照片

圖6 典型冰湖與地形起伏度疊加分析圖及實地照片

組合型冰湖周邊起伏度類型具顯著差異，湖泊后緣多與冰川相接，湖泊類型以中-小型冰斗湖和冰磧阻塞湖為主。此類湖泊所處地質環境條件相對復雜，周邊各地勢起伏等級單元均有分布，后緣冰川易崩落入湖，形成涌浪可能誘發潰決，潰決概率較高，以嘎龍錯較為典型。遙感影像及實地驗證結果顯示，嘎龍錯周邊地形特征(圖6(d))與分析結果相符，屬組合型(圖6(e))，冰湖前緣終磧堤兩側及后緣部分區域地形高陡(圖6(f))，局部具備冰滑坡、冰崩等成災條件，可能引發冰湖潰決災害。

跨級型冰湖周邊起伏度的等級變化不連續，多為冰蝕湖、冰磧湖，規模均在中型以下，以多爾雷山東北部的嘉龍錯最為典型。從遙感影像(圖6(g))及實地調查結果看，該湖東側地勢較緩，而西北側起伏度變化明顯(圖6(h))，與分析結果一致。冰湖側磧壟與后緣相接處可見顯著地形起伏度跨級界線(圖6(i)紅線標識位置)，湖面與后緣冰川具有一定高差，上部冰川易失穩發生大規?；?，進而誘發冰湖潰決災害，因此這類冰湖發生潰決的可能性極高。

4.3.2 對冰湖潰決隱患的評價

基于上述分析，可進一步依據地形起伏度對冰湖潰決隱患進行分級評價(表3),其中單一型冰湖規模不一、以小型為主，由于周邊地形平緩，不具備成災條件，故整體潰決隱患低； 組合型冰湖規模以小型居多，具一定成災條件，潰決隱患中等； 跨級型冰湖規模適中，后緣冰川所處地形落差較大，具有顯著的成災條件，潰決隱患高。區內屬低潰決隱患的冰湖993處，占比97.35%，在規模上涵蓋小、中、大3種類型冰湖，較典型的有大型冰湖貢錯及中型冰湖達熱措等； 中潰決隱患的冰湖18處，占比1.77%，規模上包括小型及大型2類，以大型冰湖崗西錯、嘎龍錯等最為典型； 高潰決隱患的冰湖9處，占比0.88%，包括小型及中型2種規模的冰湖，其中屬中型冰湖嘉龍錯潰決隱患最高(圖7)。

表3 冰湖潰決隱患評價結果

圖7 地形起伏度對冰湖潰決隱患的影響分析

從評價結果中不同隱患等級冰湖的數量來看，全區近98%的冰湖潰決隱患極低，僅有不足2%的冰湖存在潰決的可能； 但就湖泊規模而言，中-高隱患的冰湖面積達11.37 km2，占到全區冰湖總面積的64.31%，其對下游河道兩側牧民以及聶拉木縣城構成了極大的威脅。

值得注意的是，冰湖潰決造成的損失與其規模大小并非絕對正相關。例如2014年那隆藏布支溝發生的冰湖潰決事件[6]，該冰湖潰前面積僅0.24 km2，按本文湖泊規模劃分方案屬中型，由于冰湖后緣冰川與常年積雪區邊坡高陡，在連續降雨條件下冰川發生崩塌后沖入湖中，進而誘發潰決； 在2016年7月聶拉木縣樟木鎮上游發生的次仁瑪錯右側支溝冰湖潰決事件中[6]，冰湖潰決前面積僅0.01 km2，屬小型規模，但潰后形成洪水導致樟木口岸被淹，造成較大經濟損失。

因此，本次依據地形起伏度評價出的9處高潰決隱患冰湖雖規模偏小，但在當前全球氣候暖濕化的背景下，由冰川退縮消融、冰崩/冰滑坡、埋藏冰融化、強降水等誘發潰決的可能性極高，潛在威脅較大。特別是位于多爾雷山東北處沖堆普流域內的嘉龍錯曾發生過2次潰決事件[27]，再次發生災害的概率較高，需重點關注。

5 結論

本文基于GIS技術對西藏聶拉木縣境內希夏邦馬峰東南部冰湖集中分布區域開展地形起伏度研究，從冰湖地形起伏度組合特征的角度出發，對區內冰湖潰決隱患進行分析評價，得出如下結論：

1)利用均值變點法實現區內地形起伏度的提取，經有效計算及準確性檢驗，該區在1∶5萬比例尺下最佳采樣單元為21×21(采樣單元面積為0.39 km2)，提取結果可準確反映地形變化特征。

2)研究區整體地勢以小起伏山地為主，中起伏山地及丘陵次之，大起伏山地、臺地及平原分布較少，多見高落差基巖山地夾持縱深河谷的地貌特征。

3)區內冰湖周邊地形起伏度特征可分為單一起伏等級型、多起伏等級組合型、跨起伏等級型3種組合類型，對應潰決隱患等級依次升高，以跨起伏等級型冰湖潰決隱患最大。

4)全區1 020個冰湖中，潰決隱患為低、中、高級別的數量分別占總體的97.35%，1.77%和0.88%，中-高隱患的冰湖面積占到全區總體的64.31%，高潰決隱患冰湖中以嘉龍錯特征最為典型，對下游河道兩側牧民以及聶拉木縣城構成威脅。

5)結果顯示，以地形起伏度作為冰湖潰決隱患評價分析指標，對完善現有冰湖潰決危險性判別方法和體系具有實際意義，可有效提高評價結果精度，為同類地區開展相關研究提供理論參考。