不同精度下地表穩定性模型的評價

趙曉東， 王順東， 張泰麗， 吳鑫俊， 丁 茜

(1.大連大學建筑工程學院, 大連 116622； 2.中國地質調查局南京地質調查中心, 南京 210016)

地質災害評價，國外亦稱為滑坡敏感性評價，中國對地質災害評價的劃分包括易發性、危險性、易損性及風險評價等，其中易發性評價是地質災害條件和地質災害發生的空間概率統計分析，注重靜態地質災害敏感條件分析和發生的空間概率預測，是地質災害危險性和風險評估的基礎。易發性評價的核心內容包括地質災害特征、空間密度、易發條件和潛在易發區域預測評價，最終形成研究區地質災害易發性分布圖，其主要目的和用途就是為政府制定土地利用規劃及災害預警提供基礎資料和依據[1]。

目前滑坡災害評價的模型眾多，常用的包括SINMAP、SHALSTAB、CAHSM、TRIGRS、SHETRAN、GEOTOP-FS和SUSHI等[2]。這些評價模型大致可分為確定性模型、統計模型、和非參數模型[3]。確定性模型是以滑坡發生的物理力學為基礎，需要邊坡較為詳細的巖土力學參數，結果會受到數據量和數據精度的影響。相對于其他模型，確定性模型可以較為客觀地反映邊坡的物理變化，其缺點是受到人力、物力等方面的考慮，缺乏較為詳細的巖土力學參數，該方法通常只能適用于單個滑坡或者小范圍區域[4]。在常用模型中，侵蝕地表穩定性模型(stability index mapping，SINMAP)由于耦合了無限斜坡模型和穩態水文模型，利用由穩定狀態水文模型獲取的地形濕度指數、由數字高程模型(digital elevation model，DEM)獲取的坡度、單位寬度匯水面積等數據，結合地質考察資料，可在地理信息系統(GIS)空間分析平臺上建立定量分析模型，獲得地表穩定性分級，實現對研究區域的地表穩定性評價[5]。

SINMAP模型是Pack等[6]在Montgomery和Dietrich提出的淺層滑坡穩定性模型(shallow landsliding stability，SHALSTAB)基礎上的改進模型，考慮了土壤物理參數的不確定性對輸出結果的影響，保留了土壤黏聚力及植物根系對斜坡穩定性的加強作用，可模擬不同降雨條件下斜坡穩定性評價，非常適用于東南沿海地區臺風暴雨型地質災害易發性的評價。

SINMAP模型的輸入參數分為兩類，一類為由DEM轉化而來的地形數據，另一類是巖土物理力學參數。DEM精度會對第一類參數產生直接影響，隨著精度的降低坡度會更加平緩；DEM分辨率越高，地圖越能清晰呈現真實的流域面積和坡度。但DEM精度越高，其獲取的成本就越高。要做到技術參數可行，經濟成本合理就必須考察不同精度DEM對易發性評價模型的影響程度。為了探明DEM精度對應用SINMAP模型進行地質災害易發性評價的影響，本文選用免費全球30 m DEM、日本ALOS衛星的12.5 m和課題組購買的5 m三種DEM精度，兩種日最大降雨量工況，對SINMAP模型的地質災害易發性評價結果進行了對比分析研究。

1 SINMAP模型及應用

SINMAP是由Pack和Tarboton等開發出的基于無限斜坡模型和穩態水文模型的侵蝕地表穩定性模型。在該模型中，假設了滑動面平行于地表且忽略其邊緣作用，采用地表土層穩定的抗滑力與滑動力之比為安全系數，其表達式為

(1)

(2)

(3)

(4)

(5)

式中：Fs為安全系數；C為黏聚力因子；w為地形濕度指數；z為土壤厚度，m；Cr為植物根系產生的黏聚力，N/m2；Cs為土體自身黏聚力，N/m2；θ為地形坡度，(°)；ρs為土壤密度,kg/m3；ρw為水的密度,kg/m3；r為水與土壤密度之比；g為重力加速度,m/s2；φ為內摩擦角,(°)；T為土壤導水系數,m2/d；R為地下水側向穩態補給,mm/d；a為單位匯水面積,m2；h為地下水高度,m；q為流量,m3/d。

穩定性指數(SI)的定義來自安全系數，模型假設研究區域內巖土物理參數呈均勻分布，對于某一地面點來說，當C、tanφ取最小值，R/T取最大值時，此時安全系數Fs最小，該點穩定性處于最差狀態。當C、tanφ取最大值，R/T取最小值時，此時安全系數Fs為最大，該點穩定性處于最佳狀態。當Fs，min>1時，該點無條件穩定，此時穩定性指數SI=FSmin；當Fs,min<1，Fs,max>1時，該點存在發生不穩定的可能，此時穩定性指數SI定義為該地面點處于穩定狀態的概率，即SI=Prob(FS>1)；當Fs,max<1時，該點處于不穩定狀態，此時穩定性指數SI=Prob(Fs>1)=0，根據SI的大小將邊坡的穩定性劃分為5個等級，SINMAP輸出結果分類見表1，關于SINMAP的詳細理論參見文獻[6]。

表1 SINMAP模型輸出值分類Table 1 Output value classification of SINMAP model

中外學者對滑坡穩定性的預測大多應用基于ArcGIS軟件的統計模型，復合了確定性模型、信息量、邏輯回歸、潛勢度等多種模型及方法[7-10]。將SINMAP模型應用于地質災害易發性評價中國也早有先例，朱遠樂等[11]將依據層次分析法得到的滑坡預測易發性評價柵格圖與應用SINMAP模型集成分析得到尾礦庫地表穩定性指數分布圖。姚志永[12]對比了SINMAP和TRIGRS(transient rainfall infiltration and gridbased regional slope-stability model)兩種模型在湖南省辰溪縣暴雨型淺層滑坡易發性評價中的差異，得出在該區域SINMAP模型更為準確。周偉[13]研究了在白龍江流域，邏輯回歸模型在整體上優于SINMAP模型但針對單個滑坡而言后者能更好反映坡體穩定性。莊建琦等[14]研究發現在降雨強度為30 mm的條件下SINMAP模型與其建立的黃土地區淺層滑坡發育危險性評價模型結果最吻合。

中外學者將穩定性模型引入地質災害評價非常廣泛，Rabonza等[15]在菲律賓中部使用SINMAP軟件，將地形、土壤強度、水文參數和DTM數據輸入模型，計算相應的安全系數，利用SINMAP生成的滑坡圖與2002—2014年高分辨率衛星圖像得到的滑坡清查結果高度一致。Michel等[16]對巴西南部山區的達河流域應用SHALSTAB和SINMAP分別進行滑坡易發性制圖，得出該流域SINMAP的模擬結果穩定區域預測面積過大，而SHALSTAB模型模擬結果良好。

前文敘述過SINMAP模型是在SHALSTAB模型基礎上改進而成，SHALSTAB模型的提出者Dietrich等[17]及研究者Ramos等曾指出：模型的模擬精度會隨著輸入參數精度的提高而提高，地形數據分辨率以及土壤物理參數對模型模擬效果的影響尤為顯著，同樣推測DEM精度也會對SINMAP模型的輸出結果產生很大影響。

如何選取量化指標進行模擬結果的對比，不同學者使用了不同的方法，大多數研究人員應用野外獲得的土壤物理參數如內摩擦角、抗剪強度、飽和密度，但這些參數在一定區域內變化較大，而且帶有誤差，Seefelder等[18]提出了一種用于山區滑坡問題分區的方法，該方法有效且最適合于水文盆地，他認為可以用流域內工程地質巖組對區域進行劃分工況。Regiane Mara Sbroglia應用這種劃分方法對巴西某山區地帶進行了研究，模擬結果良好[19]。

2 研究區域數據整理

2.1 研究區域概況

研究數據范圍處于溫州市飛云江流域，以玉壺流域重點研究，研究區域位置見圖1。溫州陸域面積11 784 km2，其中山地面積占78.2%，地勢從西南向東北呈現梯形傾斜，有洞宮、括蒼、雁蕩山脈，泰順的白云尖，海拔1 611 m，為全市最高峰。屬亞熱帶季風氣候區，年降水量1 500～1 900 mm，8—10月是臺汛期，多年平均臺風登陸和受臺風影響3.5次/年。玉壺流域位于溫州市東南地區，隸屬于文成縣，面積103.9 km2。

圖1 研究區域位置Fig.1 Location of study area

受地質構造和第四紀晚期快速海退的控制，山區河流深切、相對高差大、坡面陡峭，坡度多大于30°，坡面穩定性差；組成山體的巖石以侏羅系白堊系火山巖為主，山坡表面和坡腳殘積層一般1～5 m；山區植被覆蓋率高，河谷和緩坡多被開墾成梯田或為居民建設用地，山區人口分散。域內山地多，山體穩定性差，人為工程活動強，臺風暴雨頻繁，是地質災害高易發區和多發區，是浙江省突發性地質災害最嚴重的市。據不完全統計，已查明的地質災害隱患點近2 000處，約占全省的20%以上；受威脅人口30 007人，占全省的21.5%。

2.2 地形數據

5 m數字高程數據為課題組購入數據，使用ArcGIS10.2軟件從DEM當中提取山脊線，繪制流域分布圖，依據中華人民共和國水利行業標準控制流域面積并命名[20]，最終選取飛云江上游玉壺流域作為研究區域。12.5 m DEM數據來源于日本ALOS衛星，30 m DEM取自NASA全球30 m高程數據，因其相鄰兩種精度比值大小相近，便于比較，故均取其原始大小計算。利用ArcGIS軟件Hillshade工具從DEM數據當中提取出玉壺流域山體陰影(圖2)。

2.3 巖土物理力學參數

巖土物理參數來自土體采樣點取樣的室內土工試驗，土體抗剪強度指標黏聚力和內摩擦角采用快剪實驗得到。風化層厚度z來自于對鉆孔土樣的分析。

工程地質巖組由1∶5 萬區域地質圖地層巖性歸類獲得，如圖3所示玉壺流域內主要分布有六種巖組，Hif(軟硬相見塊狀晶、玻屑熔結凝灰巖夾凝灰質粉砂巖)、Hs(軟硬相間的層狀沉凝灰巖、凝灰質粉砂巖夾流紋晶屑玻屑凝灰巖)、Pb(堅硬塊狀玄武巖等基性噴出巖不易崩滑巖組)、Qg[堅硬塊狀花崗(斑)巖等酸性侵入巖、潛火山巖]、Rr[堅硬塊狀流紋(斑)巖等酸性巖]、Scf(軟硬相間凝灰質粉砂巖、沉凝灰巖、粉砂巖、砂礫巖、砂巖)，研究區域內巖組分布見圖3，具體占比見表2。

圖2 玉壺流域山體陰影圖Fig.2 Hillshade map of yuhu basin

圖3 巖組分布Fig.3 Distribution of rock group

表2 巖土參數分布

2.4 SINMAP模型實現

SINMAP模型的實現環境為ArcGIS10.2版本，軟件采用柵格單元計算，將采樣點巖土參數數據整理成覆蓋整個玉壺流域的柵格數據(黏聚力c、內摩擦角φ、風化層厚度z、密度ρ)。C、R/T、tanφ為可變參數，根據研究區域內采樣點巖土力學參數分別指定C、tanφ的上限與下限，R/T參數由室內土工試驗及參考其他學者研究確定取值范圍，并假定它們在指定范圍內均勻分布[12-13]。

地形坡度由DEM在ArcGIS中使用坡度計算獲得。單位匯水面積的取值與每個柵格內模擬水流流向有關，流向的確定采用水文分析當中常用的D8算法，流向與柵格邊平行取邊長，與對角線平行則取柵格對角線長度。匯水面積A在數值上等于流量與柵格面積的乘積，其中流量由DEM在ArcGIS中使用自行開發的計算工具獲取。可以理解為單位匯水面積越大，徑流強度越大，對坡體沖刷能力越強，圖4為單位匯水面積示意圖，SHALSTAB模型提出者Dietrich繪制，有修改[21]。

為研究DEM精度及降雨條件對地質災害評價結果的影響，區分了兩種工況，劃分標準見表3。降雨數據來自對玉壺流域內雨量站數據的統計，2016年超強臺風“尼伯特”經過溫州帶來了短時間內的強降雨，文成縣日最大降雨量達到400 mm，而2017年期間沒有較大臺風過境，日最大降雨量為166.5 mm。圖5為反距離權重法插值得到的研究區域內降雨量分布圖。

圖4 單位匯水面積示意圖Fig.4 Schematic diagram of unit catchment area

表3 工況劃分Table 3 Working condition of division

圖5 降雨量分布圖Fig.5 Rainfall distribution map

圖6 易發性制圖柵格累計圖Fig.6 Susceptibility cartographic raster cumulative map

3 結果對比分析

圖6為各種工況下的易發性制圖，按照模型提出者的分級標準進行分類，通過統計易發性制圖(圖7)的數據得到SINMAP模型輸出值統計(表4)和柵格累計(圖6)，表4中列出了不同DEM精度與降雨量組合后模型的輸出值分類結果。

圖7 易發性制圖Fig.7 Susceptibility map

可以看出三種精度結果差異較大，主要是穩定和無條件不穩定區域的差別，降雨量較小的A工況下，5 m精度的無條件不穩定區域為12.1%，穩定區域為41.6%；12.5 m精度的無條件不穩定區域為6.5%，穩定區域為54.1%；30 m精度的無條件不穩定區域為4.5%，穩定區域為54.2%；DEM精度較高的5 m條件下識別出的不穩定區域接近12.5 m精度的兩倍，30 m精度的3倍，足見DEM精度對模型輸出值的影響非常大，主要是影響無條件不穩定區域的大小。相同降雨量不同DEM精度間的對比很明顯表現出一條規律：精度越高，模型識別出的無條件不穩定區域越多。

對比工況A下三種DEM精度輸出結果：5 m精度下無條件不穩定區域為12.1%，12.5 m精度的無條件不穩定區域為6.5%，30 m精度的無條件不穩定區域為4.5%，像元大小同樣增大2.5倍，由5 m增大到12.5 m時無條件不穩定區域減小量為5.6%，而12.5 m增大到30 m時無條件不穩定區域面積減小量僅為2%，這說明精度降低后模型識別不穩定區域能力的影響迅速減弱。

同一精度下，工況B(降雨量較大)的總體不穩定區域面積最大，這也符合推測，降雨量會影響T/R范圍進而影響輸出結果，但無條件不穩定區域相同，推測原因可能是由于該區域汛期降雨量大，日最大降雨量都在100 mm以上，地下水側向補給容易達到飽和狀態，超過界限后降雨量的增大對SINMAP模型不穩定區域的增加影響較小。相差三倍的日最大降雨量數據在5、12.5、30 m三種精度下模型穩定區域面積分別相差2%、3.6%、4.5%，無條件不穩定區域面積相同。

表4 SINMAP模型輸出統計Table 4 Statistics of SINMAP model output

4 結論

(1)DEM精度對SINMAP模型下地質災害易發性評價結果影響非常大，特別是對極端情況的影響。隨著DEM精度的提高，模型輸出的無條件不穩定區域擴大，無條件穩定區域減小，中間過渡地帶變化并不明顯。在降雨量相同的情況下，DEM精度較高的5 m條件下識別出的不穩定區域接近12.5 m精度的兩倍，30 m精度的3倍。

(2)DEM精度對SINMAP模型下地質災害易發性評價結果影響是非線性的，高精度轉為中精度，中精度轉為低精度，模型輸出無條件不穩定區域下降速率降低很快。當DEM精度降低到一定范圍后，精度的降低對模型輸出的地質災害易發性制圖影響很小。

(3)日最大降雨量超過100 mm情況下降雨量對模型輸出值的影響不大，相差3倍的日最大降雨量數據在5、12.5、30 m三種精度下模型穩定區域面積僅相差2%、3.6%、4.5%，無條件不穩定區域面積相同。