基于機器學習的長沙市滑坡災害快速風險評價①

王 璨， 肖 浩， 肖 婷， 方亞其， 劉磊磊

（1.湖南省地質災害調查監(jiān)測所，湖南 長沙 410004； 2.湖南省地質災害監(jiān)測預警與應急救援工程技術研究中心，湖南 長沙 410004； 3.有色金屬成礦預測與地質環(huán)境監(jiān)測教育部重點實驗室，湖南 長沙 410083； 4.湖南省有色資源與地質災害探查湖南省重點實驗室，湖南 長沙 410083； 5.中南大學地球科學與信息物理學院，湖南 長沙 410083）

我國地質環(huán)境脆弱、孕災條件復雜，是世界上滑坡災害最嚴重的國家之一[1］。 滑坡災害風險評價是滑坡災害風險管理的核心，是認識滑坡災害災情、制定防災政策、實施防治措施的重要依據(jù)。 國內外學者對區(qū)域滑坡災害風險評價的各環(huán)節(jié)進行了大量探索性研究。 在危險性評價中，常用的評價方法為定性分析法、確定性分析法和概率統(tǒng)計法。 其中，支持向量機、人工神經網絡、隨機森林等統(tǒng)計機器學習算法[2］因其優(yōu)良的評價效果得到廣泛應用。 目前，評價機器學習模型性能的常用指標有模型精度和接受者操作特性曲線（ROC）[3］，這些指標可以有效評價模型在樣本集數(shù)據(jù)建模預測中的優(yōu)劣，但不能反映模型在非樣本區(qū)域的評價性能。 對承災體易損性評價的研究，目前處于從定性分析到定量計算的發(fā)展階段。 在區(qū)域性易損性評價工作中，承災體系統(tǒng)相當復雜，不同研究區(qū)承災體與災害的作用機理不盡相同，導致定量評價模型并不具備很好的普適性。 風險包括災害體本身及其造成的后果，常通過精細化承災體分析計算定量風險，或采用啟發(fā)式矩陣的經驗方法對風險進行分級[4］。

本文采用機器學習模型進行長沙市滑坡災害危險性評價，基于頻率比法檢驗模型在非樣本集上的性能，并對評價結果較差的模型進行矯正；不考慮承災體的脆弱性差異，選取代表性指標因子構成易損性評價體系，實現(xiàn)易損性的快速評價；引入風險分區(qū)矩陣完成長沙市滑坡災害風險區(qū)劃，以期為長沙市的防災減災規(guī)劃提供科學依據(jù)和理論基礎。

1 研究區(qū)概況

長沙市位于湖南省東部偏北，湘江下游和長瀏盆地西緣，其地理位置介于東經111°53′～114°15′、北緯27°51′～28°40′。 該地屬亞熱帶季風濕潤氣候區(qū)，其特點是春冬多雨、夏秋多晴、嚴冬期短、暑熱期長。 域內各個地質歷史時期的地層均有出露，花崗巖體廣布，以第四系地層分布最多。 據(jù)湖南省國土資源廳提供的資料顯示，長沙市滑坡災害點共808 處（見圖1），按照地質災害災情分級標準，以小型滑坡為主（802 處），中型滑坡層占極少數(shù)（6 處）。

圖1 研究區(qū)滑坡災害點分布圖

2 研究方法

2.1 危險性評價方法

2.1.1 隨機森林模型

隨機森林模型（Random Forests， RF）是利用多棵決策樹作為分類器對樣本進行訓練與預測的分類模型[5］，算法模型如圖2 所示。 該模型方法結合了基于訓練樣本操作的裝袋算法和基于特征集操作的隨機子空間方法，將多個決策樹結合在一起，有放回地隨機選取樣本并將部分特征作為輸出[6］，其預測結果由每棵樹的結果投票產生，取投票數(shù)最多的類或取其平均值作為結果，具有較好的準確率和穩(wěn)定性。

圖2 隨機森林模型算法示意

2.1.2 極限梯度提升

極限梯度提升（XGBoost）使用梯度增強框架，是一種基于決策樹的集成方法[7］。 該方法的核心原理是逐棵建立分類或回歸樹，然后利用之前樹的殘差來訓練后續(xù)的模型。 它將之前已訓練樹的結果進行集成從而得到更好的結果。 XGBoost 預測值可按下式計算：

式中＾y（t）為最終樹群的預測值；＾y（t－1）為之前樹群的預測值；xi為樣本i對應的特征值；ft（xi）為新生成樹的預測值；t為基礎樹模型的總數(shù)。

2.2 基于FR 的危險性評價模型校正

常用的ROC 曲線和AUC（Area Under Curve，曲線下面積值）值能夠很好地反映評價模型在樣本集數(shù)據(jù)上的評價性能，但不能反映模型在非樣本集數(shù)據(jù)上的性能表現(xiàn)。 本文采用滑坡發(fā)生頻率比（Frequency ratio，F(xiàn)R）與滑坡空間概率之間的線性理論對模型預測結果進行校正[8］，以補償ROC 曲線和AUC 指標在評價非樣本集數(shù)據(jù)上的不足。

FR 表征了各子分類對于滑坡發(fā)生的重要程度，可由下式計算得到：

子分類內的滑坡數(shù)與子分類面積之比表示子分類內滑坡發(fā)生的空間概率，而研究區(qū)內滑坡數(shù)與研究區(qū)面積之比為一常數(shù)，那么滑坡發(fā)生空間概率與FR的意義就相一致了[8］。 因此可以使用FR值對危險性評價模型進行校正：當FR與危險性值成線性關系時，該危險性區(qū)劃圖能更合理地表達滑坡發(fā)生的空間概率。FR＞1 時，說明子分類有利于滑坡的發(fā)生；FR＜1 時，說明該子分類不利于滑坡發(fā)生[9］。

2.3 易損性評價方法

采用層次分析法綜合分析各指標因素，建立易損性評價模型。 層次分析法能夠將復雜的定性問題簡明化，合理計算出各指標的綜合權重系數(shù)，應用步驟如下[10］：

1） 將要素分解為目的層、準則層和方案層，構造遞階層次結構模型；

2） 根據(jù)影響因素的結構隸屬關系，引用數(shù)字1 ～9及其倒數(shù)作為標度構造判斷矩陣；

3） 進行一致性檢驗并確定各因素的總權重值，帶入下式計算得到易損性指數(shù)[11］：

式中Vi為評價單元的地質災害綜合易損性指數(shù)；i為評價單元；ωj為第i個評價單元第j個評價指標的權重值；yj為第i個評價單元第j個評價指標標準化后的取值。

3 基于機器學習的滑坡危險性評價

3.1 評價因子選取

從地形地貌、地質條件、自然環(huán)境和人類工程活動要素中提取15 個評價因子研究滑坡的形成機制。 數(shù)據(jù)分別來源于數(shù)字高程模型、長沙市1 ∶50000 地質圖、長沙市1 ∶50000 地形圖及土地利用數(shù)據(jù)。 評價區(qū)域屬于大區(qū)域，評價單元選用25 m×25 m 的柵格。

在模型計算前，通過皮爾遜相關系數(shù)對指標因子進行相關性分析，得到15 個評價因子間的皮爾遜系數(shù)見圖3，其中平面曲率、剖面曲率、地形起伏度和水流強度指數(shù)與部分因子之間相關性系數(shù)的絕對值大于0.5，表明這些因子之間相關程度較高[12］。 以剔除平面曲率等因子后剩余的11 個因子組成危險性評價指標體系，并將各因素指標按照相應標準劃分出多個二級狀態(tài)。 其中，離散型數(shù)據(jù)依據(jù)野外調查論證制定劃分標準，連續(xù)性數(shù)據(jù)采用自然間斷點法劃分等級，相應因子分布如圖4 所示。

圖3 皮爾遜系數(shù)矩陣

圖4 滑坡影響因子分布圖

3.2 模型建立

危險性評價樣本集由正樣本和負樣本構成。 正樣本數(shù)據(jù)為滑坡編錄數(shù)據(jù)中的808 個滑坡點，在滑坡點100 m 緩沖區(qū)外隨機采樣生成808 個負樣本數(shù)據(jù)，并將樣本集按7 ∶3切割為訓練集和測試集。 在Python 中分別構建RF 和XGBoost 模型，導入樣本集數(shù)據(jù)進行訓練，將訓練好的模型應用至全區(qū)實現(xiàn)危險性評價。

ROC 曲線常用來評價易發(fā)性模型精度，AUC 與評價精度正相關[13］。 一般來說，ROC 曲線越接近（0，1）點，AUC 值越接近1，模型的預測準確率越高。 采用ROC 曲線對RF 和XGBoost 模型預測結果進行評價檢驗，2 個機器學習模型的AUC 值分別為0.77 和0.76（見圖5），精度相近。

圖5 ROC 曲線

3.3 基于FR 模型校正

采用自然斷點法將危險性值進行分區(qū)，得到RF和XGBoost 模型下的滑坡危險性分區(qū)如圖6 所示。 為使分布圖能更合理表達滑坡發(fā)生的空間概率，基于FR值對其進行檢驗及校正。

圖6 不同模型生成的危險性分布圖

將危險性值等間隔劃分為50 個區(qū)間，統(tǒng)計各個區(qū)間的FR 值，圖7 為RF 和XGBoost 模型中危險性值與FR 值的關系。 可知，2 個模型中危險性值與FR 值均不成線性關系，且2 個模型都在低危險區(qū)低估了滑坡發(fā)生的概率，而在高危險區(qū)高估了滑坡發(fā)生的空間概率，因此需要對危險性分布圖進行校正。 其中，RF 模型比XGBoost 模型下的危險性預測值分布更加極端，滑坡幾乎只在高危險區(qū)中出現(xiàn)。

圖7 各危險性分布圖的頻率比分布

2 個模型FR 值與危險性值的擬合關系表達式分別為：

將原始的危險性值代入式（4）～（5）中，并將歸一化后的結果作為校正后的危險性值，從而將FR 值與危險性值之間的非線性關系轉化為線性關系，結果如圖8 和圖9 所示。 可知，RF 模型校正后的結果幾乎分布于低危險性值區(qū)，頻率分布過于極端。 XGBoost模型校正后危險性值與FR 值呈現(xiàn)良好的線性關系，且危險性值分布較均勻，表明XGBoost 模型在非樣本集上的預測性能更加優(yōu)良。 以矯正后的XGBoost 模型結果作為最終的危險性區(qū)劃結果（見圖9（b）），研究區(qū)內低、較低、中、較高和高危險區(qū)面積分別占比36.3%、28.5%、22.4%、11.5%和1.3%，低危險區(qū)到高危險區(qū)面積占比依次降低。

圖8 校正后危險性分布圖的頻率比分布

圖9 校正后模型生成的危險性分布圖

4 基于層次分析法的承災體易損性評價

4.1 評價因子選取及權重確定

城市地質災害以人口、建筑物和道路為重要的承災體，選取這3 個評價因子密度值來表征各因子的易損性程度，并構建遞階層次模型以確定各因子之間的從屬關系，相互對比之后的構造判斷矩陣如表1所示。 計算結果顯示，最大特征值λmax＝3，一致性結果CR＜0.1，滿足一致性檢驗要求[14］，將λmax對應的特征向量進行歸一化得到各指標因子權重值如表1所示。

表1 易損性評價因子判斷矩陣及權重值

4.2 易損性評價結果

基于ArcGIS 平臺，根據(jù)式（3）計算每個評價單元的易損性指數(shù)，按自然間斷法將研究區(qū)由低到高劃分為5 個等級，易損性分級與區(qū)劃如圖10 所示。 據(jù)統(tǒng)計，低、較低、中、較高和高易損區(qū)面積分別占長沙市面積的49.9%、26.6%、16.5%、4%和3%；較高-高易損區(qū)主要為城鎮(zhèn)及人口聚集區(qū)、風景區(qū)和交通干線。

圖10 承災體易損性區(qū)劃圖

5 滑坡災害風險評價

5.1 風險分區(qū)矩陣

在地質災害的風險評估中，風險是因變量，風險等級由滑坡災害危險性等級與承災體易損性等級共同決定。 采用如圖11 所示的風險矩陣對長沙市風險進行定性評估。 其中不同的數(shù)字代表不同的程度等級，1～5分別表示低、較低、中、較高和高。

圖11 風險分區(qū)矩陣

5.2 風險評價結果

風險評價結果如圖12 所示。 低、較低、中、較高和高風險區(qū)分別占研究區(qū)面積的41.4%、32.9%、21.1%、4.2%和0.4%。 可知，隨著風險等級增加，相應區(qū)劃面積遞階減小，具備較高的評價精度。

圖12 滑坡風險區(qū)劃圖

較高-高風險區(qū)主要分布在溝谷、城鎮(zhèn)和交通干線等區(qū)域。 對比長沙市行政區(qū)劃圖，可知芙蓉區(qū)、天心區(qū)東部等地處于較高-高風險區(qū)，這些區(qū)域是全市城市化程度最高、人類工程活動最強烈的地帶，聯(lián)合自然因素的影響產生了大量的潛在不穩(wěn)定斜坡。 同時，雨花區(qū)、開福區(qū)等區(qū)域的危險性等級不高，但這些區(qū)域承災體有較高的易損性等級，故相應區(qū)域也呈現(xiàn)出較高的風險態(tài)勢。

6 結論

采用RF 和XGBoost 模型對長沙市滑坡危險性進行評價，基于FR 值對結果進行檢驗和校正，以提高模型在非樣本集數(shù)據(jù)上的合理性和準確性，并采用層次分析法實現(xiàn)承災體易損性快速評價，集成危險性和易損性評價結果進而生成滑坡風險評價結果。 主要研究結論如下：

1） 經FR 值檢驗，RF 和XGBoost 模型都在低危險區(qū)低估了滑坡空間概率，而在高危險區(qū)高估了滑坡空間概率。 以頻率比法校正后，XGBoost 模型比RF 模型評價結果分布更合理。 在RF 模型（AUC ＝0.77）與XGBoost 模型（AUC ＝0.76）性能相近的情況下，以在非樣本集表現(xiàn)更好的XGBoost 模型生成危險性評價結果。

2） 危險性評價結果表明，研究區(qū)內低、較低、中、較高和高危險區(qū)面積分別占比36.3%、28.5%、22.4%、11.5%和1.3%，面積占比依次降低，精度較高。

3） 選取人口密度、道路密度和建筑物密度構成易損性快速評價指標體系，采用層次分析法得到人口密度對易損性評價的影響最大。

4） 運用數(shù)值分級方法將研究區(qū)劃分為低風險區(qū)（41.4%）、較低風險區(qū)（32.9%）、中風險區(qū)（21.1%）、較高風險區(qū)（4.2%）和高風險區(qū)（0.4%）。 長沙市較高、高風險多分布在溝谷、城鎮(zhèn)和交通干線等地，與實際情況一致，模型預測精度較高。