基于邏輯回歸的四川青川縣區(qū)域滑坡災(zāi)害預(yù)警模型

方然可，劉艷輝，蘇永超，黃志全

（1.華北水利水電大學(xué)，河南 鄭州 450045；2.中國地質(zhì)環(huán)境監(jiān)測院（自然資源部地質(zhì)災(zāi)害技術(shù)指導(dǎo)中心），北京 100081；3.洛陽理工學(xué)院，河南 洛陽 471023）

四川省青川縣位于四川盆地的北部邊緣，龍門山斷裂帶的地震多發(fā)帶，轄區(qū)內(nèi)山高坡陡，地形切割強烈，地質(zhì)構(gòu)造復(fù)雜，是滑坡等地質(zhì)災(zāi)害高發(fā)區(qū)。近年來，特別是“5.12”汶川地震后，多次區(qū)域強降水作用下，滑坡災(zāi)害頻發(fā)，點多面廣的防災(zāi)形勢日益嚴(yán)峻，居民正常生產(chǎn)生活受到嚴(yán)重威脅。在區(qū)域尺度上，開展降水誘發(fā)的區(qū)域滑坡災(zāi)害預(yù)警十分關(guān)鍵。自2003年以來，中國廣泛開展各級地質(zhì)災(zāi)害氣象預(yù)警，在提高社會公眾的防災(zāi)減災(zāi)意識和有效減輕地質(zhì)災(zāi)害造成人員傷亡和財產(chǎn)損失等方面起到了積極作用[1-3]。

預(yù)警模型是成功開展地質(zhì)災(zāi)害氣象預(yù)警的關(guān)鍵，為此，大量學(xué)者開展了不懈的努力和探索。區(qū)域地質(zhì)災(zāi)害預(yù)警模型主要分為統(tǒng)計預(yù)警模型和動力預(yù)警模型兩大類，能在區(qū)域尺度得到實際應(yīng)用的主要是統(tǒng)計預(yù)警模型，特別是基于統(tǒng)計原理的臨界降水閾值模型[4-6]，在美國、中國香港、日本等滑坡早期預(yù)警系統(tǒng)中都得到了廣泛應(yīng)用。2009年以來，劉傳正等[2-3]提出了顯式統(tǒng)計預(yù)警原理，多位學(xué)者針對不同地區(qū)的特點，構(gòu)建了相應(yīng)的統(tǒng)計預(yù)警模型[7-13]，有效支撐了中國大陸各級地質(zhì)災(zāi)害氣象預(yù)警工作。但預(yù)警模型研究受制于研究區(qū)滑坡誘發(fā)機理復(fù)雜、調(diào)查監(jiān)測數(shù)據(jù)不足等限制，仍存在區(qū)域預(yù)警精度有限、精細(xì)化不足等問題。

2015年，原國土資源部和中國氣象局在四川省青川縣建立了第一個國家級地質(zhì)災(zāi)害氣象預(yù)警試驗區(qū)，針對地質(zhì)災(zāi)害預(yù)警的地質(zhì)環(huán)境專項調(diào)查監(jiān)測、降雨加密站點布設(shè)與監(jiān)測等工作逐步展開并取得了豐碩成果[14-15]。多年的數(shù)據(jù)積累和系統(tǒng)的預(yù)警專項調(diào)查，為區(qū)域滑坡災(zāi)害預(yù)警模型研究奠定了雄厚的數(shù)據(jù)基礎(chǔ)；人工智能和大數(shù)據(jù)的蓬勃發(fā)展，為預(yù)警模型發(fā)展提供了新的方法。

本文以四川省青川縣為研究區(qū)，開展了基于邏輯回歸算法的區(qū)域滑坡災(zāi)害預(yù)警模型研究并校驗，推動了人工智能在滑坡災(zāi)害預(yù)警領(lǐng)域中的應(yīng)用，有效支撐了我國正在開展的區(qū)域地質(zhì)災(zāi)害氣象預(yù)警業(yè)務(wù)，具有重要理論意義和應(yīng)用價值。

1 數(shù)據(jù)獲取與清洗

本研究區(qū)為四川省青川縣，位于川、甘、陜交界處，轄36個鄉(xiāng)鎮(zhèn)，總?cè)丝诩s22萬。研究數(shù)據(jù)主要來源于四川省青川縣1:5萬地質(zhì)災(zāi)害與地質(zhì)環(huán)境調(diào)查成果、青竹江流域地質(zhì)災(zāi)害調(diào)查成果、地質(zhì)災(zāi)害災(zāi)情直報系統(tǒng)、青川縣加密雨量監(jiān)測數(shù)據(jù)以及青川預(yù)警試驗區(qū)補充調(diào)查監(jiān)測成果等[14-16]。

1.1 滑坡和非滑坡（正負(fù)樣本）數(shù)據(jù)

據(jù)調(diào)查成果[16]，截至2015年底，青川縣共發(fā)生崩塌、滑坡、泥石流和不穩(wěn)定斜坡災(zāi)害1 672處，其中崩塌262處、滑坡643處、泥石流45處和不穩(wěn)定斜坡722處。災(zāi)害類型以不穩(wěn)定斜坡最多，其次為滑坡、崩塌，泥石流數(shù)量相對少。本文研究對象為廣義滑坡，包含了崩塌、滑坡和不穩(wěn)定斜坡。青川縣廣義滑坡占總數(shù)的97.3%，是青川縣地質(zhì)災(zāi)害的最主要災(zāi)種。

正樣本是指已經(jīng)發(fā)生滑坡的點，正樣本的采樣一般以歷史滑坡編目數(shù)據(jù)為依據(jù)進(jìn)行篩選。篩選標(biāo)準(zhǔn)為具有確定的空間地理坐標(biāo)和時間坐標(biāo)（一般精確到日）。負(fù)樣本是指沒有發(fā)生滑坡的點，無法直接獲取。本研究以正樣本為基礎(chǔ)，通過時空采樣確定負(fù)樣本[17]。空間上，在正樣本120 m（3倍的預(yù)警網(wǎng)格單元大小）緩沖區(qū)外空間隨機采樣（正負(fù)樣本采樣比例1∶2），確定負(fù)樣本的空間位置；然后在汛期范圍內(nèi)隨機時間，確定負(fù)樣本的時間屬性。最終確定2010—2018年滑坡災(zāi)害正負(fù)樣本1 826個:613個滑坡點（正樣本）和1 213個非滑坡點（負(fù)樣本）。

1.2 地質(zhì)環(huán)境數(shù)據(jù)

孕育滑坡災(zāi)害的地質(zhì)環(huán)境數(shù)據(jù)主要包括地形地貌、地層巖性、地質(zhì)構(gòu)造、溝谷水系，以及人類工程活動等。據(jù)相關(guān)分析結(jié)果[16]，地形地貌是滑坡發(fā)育的重要條件，地層巖性是滑坡發(fā)育的物質(zhì)基礎(chǔ)條件，控制了青川縣滑坡災(zāi)害的空間分布；地質(zhì)構(gòu)造特別是斷裂帶的分布對滑坡災(zāi)害發(fā)育影響明顯；溝谷水系廣泛分布，直接影響著斜坡穩(wěn)定性；人類工程活動，特別是道路切坡、房屋建筑等，對自然斜坡進(jìn)行改造，誘發(fā)或加劇了滑坡災(zāi)害的發(fā)生發(fā)展。因此，在前期研究成果基礎(chǔ)上，本文選取了坡度、地貌、地層巖性、斷裂、房屋分布等11個地質(zhì)環(huán)境因子參與模型構(gòu)建。

1.3 降雨誘發(fā)數(shù)據(jù)

降雨是研究區(qū)滑坡災(zāi)害發(fā)生的主要誘發(fā)因素。本研究降雨數(shù)據(jù)主要來源于氣象部門，收集整理了青川縣域43個加密雨量站點逐日雨量數(shù)據(jù)，通過空間插值，以公里網(wǎng)格為單元構(gòu)建雨量數(shù)據(jù)庫，保證了降雨數(shù)據(jù)的站點密度和監(jiān)測精度。選取當(dāng)日雨量，前1～15日逐日雨量等16個降雨因子參與模型構(gòu)建。

1.4 訓(xùn)練樣本集

以200 m×200 m為單元對研究區(qū)進(jìn)行網(wǎng)格剖分（約82 000個網(wǎng)格單元），分別匹配27個影響因子（包括坡度等11個地質(zhì)環(huán)境因子，當(dāng)日雨量等16個降雨因子）作為輸入特征參數(shù)，滑坡是否發(fā)生（滑坡正樣本為1，非滑坡負(fù)樣本為0）作為輸出特征參數(shù)。構(gòu)建了青川縣區(qū)域滑坡災(zāi)害預(yù)警訓(xùn)練樣本集，樣本個數(shù)1 826個（圖1），樣本集輸入特征及參數(shù)見表1。

圖1 青川縣區(qū)域滑坡災(zāi)害預(yù)警訓(xùn)練樣本集Fig.1 Training sample set of regional landslide disasters in Qingchuan County

表1 訓(xùn)練樣本輸入特征及參數(shù)Table1 Input characteristics and parameters of the training samples

1.5 數(shù)據(jù)清洗

數(shù)據(jù)清洗是指對數(shù)據(jù)進(jìn)行缺失值插補和異常值（噪音值）識別處理[18]。在機器學(xué)習(xí)模型構(gòu)建過程中，訓(xùn)練樣本集的好壞直接決定了模型準(zhǔn)確性和泛化能力。原始數(shù)據(jù)集中常存在人工錯誤、數(shù)據(jù)傳輸誤差、設(shè)備故障和地質(zhì)信息模糊等問題，嚴(yán)重影響著訓(xùn)練樣本集的質(zhì)量。因此，對原始數(shù)據(jù)集進(jìn)行數(shù)據(jù)清洗，是有效提升模型精度的必要基礎(chǔ)。

數(shù)據(jù)清洗方法要根據(jù)數(shù)據(jù)實際情況進(jìn)行選擇。例如，若空值較少可直接刪除，若空值所占比例較大，則不可直接刪除，一般可用均值代替補全。對于離群值問題，因離群值不一定是錯值，一般通過降低權(quán)重或是增加此類樣本方式進(jìn)行清洗[19-20]。

2 模型構(gòu)建與檢驗

2.1 算法與工具

邏輯回歸算法是一種常用的非線性二分類因變量回歸統(tǒng)計模型，在機器學(xué)習(xí)領(lǐng)域也得到了廣泛使用。該算法通過極大似然估計法估計參數(shù)，具有一致的漸進(jìn)正態(tài)性[21-22]。與一般線性回歸算法的區(qū)別在于邏輯回歸算法通過Sigmoid 函數(shù)可以把輸出結(jié)果約束在[0，1]之內(nèi)：由圖2所示，邏輯回歸算法輸出Y值介于0～1 之間，因此可以引入閾值的概念，從而對輸出Y進(jìn)行二分類。例如設(shè)置閾值為0.5，算法輸出大于0.5 則判定為1，反之則判定為0。

圖2 Sigmoid 函數(shù)Fig.2 Sigmoid function

滑坡災(zāi)害預(yù)警時，影響滑坡是否發(fā)生的地質(zhì)環(huán)境、降雨等因素可以作為自變量，滑坡發(fā)生或不發(fā)生可以作為分類因變量。設(shè)P為滑坡發(fā)生的概率，取值范圍為[0,1]，1-P即為滑坡不發(fā)生的概率。P/（1-P）為滑坡發(fā)生與不發(fā)生的概率比值，對其取自然對數(shù)ln[P/（1-P）]：

式中：P-滑坡發(fā)生概率；

x1，x2，···，xm-影響因子；

β0，β1，···，βm-邏輯回歸算法的回歸系數(shù)。

國內(nèi)外學(xué)者在該領(lǐng)域也開展了相應(yīng)探索，Lee S 等[23]和Ohlmacher G C 等[24]分別在韓國龍仁地區(qū)和美國堪薩斯州東北地區(qū)，選擇誘發(fā)滑坡的致災(zāi)影響因子，采用了邏輯回歸算法開展滑坡預(yù)測。李鐵鋒等[25]以長江三峽為研究區(qū)，把有效雨量與Logistic 回歸模型結(jié)合，建立了預(yù)警模型。孫德亮[26]應(yīng)用邏輯回歸算法，設(shè)置了不同閾值開展對比分析，建立了滑坡災(zāi)害易發(fā)性分區(qū)模型等。

本文選取邏輯回歸算法開展滑坡災(zāi)害預(yù)警模型訓(xùn)練。模型訓(xùn)練過程通過Python 語言實現(xiàn)，調(diào)用了sklearn 庫里的LogisticRegression模型。

2.2 模型構(gòu)建與優(yōu)化

按照4∶1的比例，將青川縣訓(xùn)練樣本集（圖1，表1）劃分為訓(xùn)練集和測試集，進(jìn)行訓(xùn)練和校驗。采用貝葉斯優(yōu)化算法、五折交叉驗證開展模型參數(shù)優(yōu)化。目前最常用的模型參數(shù)優(yōu)化方法包括傳統(tǒng)方法和超參數(shù)優(yōu)化算法。傳統(tǒng)方法，也稱為網(wǎng)格搜索法，方法優(yōu)化精度和速度成反比。為了參數(shù)優(yōu)化更高效，出現(xiàn)了超參數(shù)優(yōu)化算法—貝葉斯優(yōu)化算法[27]。貝葉斯優(yōu)化算法采用了高斯過程，通過增加樣本數(shù)量來擬合目標(biāo)函數(shù)分布，目標(biāo)函數(shù)通過交叉驗證精度來進(jìn)行優(yōu)化，每次迭代都會輸出一次超參數(shù)，在尋找最優(yōu)值的過程中優(yōu)化超參數(shù)。通過貝葉斯優(yōu)化算法，對邏輯回歸模型的正則化指數(shù)C值進(jìn)行參數(shù)優(yōu)化，最終得到模型最優(yōu)參數(shù)，C=2。

基于優(yōu)化后的Logistic模型，分別設(shè)置0.25，0.5，0.75 三種不同的閾值對模型結(jié)果進(jìn)行二分類，得到混淆矩陣，見表2。

根據(jù)表2，對比分析三種閾值分類結(jié)果的精度，可見，當(dāng)閾值取0.5時，模型總精度最高。

2.3 模型檢驗

模型檢驗從3個方面來進(jìn)行：

（1）準(zhǔn)確率（Accuracy，ACC），表達(dá)的是模型的精度。模型準(zhǔn)確率，是判斷模型預(yù)測分類正確的樣本數(shù)和總樣本數(shù)的比值。另外，還有精確率（Precision）、召回率（Recall）和F1值等指標(biāo)。

（2）ROC曲線和AUC值，表達(dá)的是模型泛化能力。ROC（Receiver Operating Characteristic）曲線是一個畫在二維平面上的曲線；AUC（Area Under Curve），即ROC曲線下的面積。通常，AUC的值介于0.5 到1.0之間，AUC值越大，說明模型表現(xiàn)越好。

表2 不同閾值下的Logistic 回歸分類結(jié)果混淆矩陣Table2 Confuse matrix of the result of the logistic regression classification under different thresholds

（3）學(xué)習(xí)曲線（Learning Curve），描述模型擬合程度，判斷模型是否存在過擬合或欠擬合。

使用測試集對Logistic 回歸模型的準(zhǔn)確率和模型泛化能力指標(biāo)進(jìn)行評估。評估結(jié)果見表3和圖3。可見，模型準(zhǔn)確率為0.943，AUC值為0.980，模型準(zhǔn)確率和泛化能力均較好。

表3 Logistic 回歸模型分類Table3 Logistic regression model classification report

圖3 邏輯回歸模型學(xué)習(xí)曲線和ROC曲線Fig.3 Learning curve and ROC curve of the logistic regression model

開展區(qū)域滑坡實際預(yù)警時，按訓(xùn)練樣本特征屬性格式，輸入研究區(qū)各預(yù)警網(wǎng)格單元27個特征屬性，調(diào)用基于邏輯回歸算法訓(xùn)練保存好的模型開展概率預(yù)測，依據(jù)模型輸出概率P劃分滑坡災(zāi)害預(yù)警等級。輸出概率的預(yù)警等級，可參考地質(zhì)災(zāi)害氣象風(fēng)險預(yù)警標(biāo)準(zhǔn)[28]中的預(yù)警等級劃分，也可根據(jù)研究區(qū)具體情況微調(diào)確定。考慮到青川縣基于邏輯回歸算法的預(yù)警模型訓(xùn)練中閾值設(shè)定為0.5，結(jié)合研究區(qū)具體情況，將地質(zhì)災(zāi)害氣象預(yù)警概率等級劃分標(biāo)準(zhǔn)調(diào)整如下，當(dāng)輸出概率P≥40%且P＜60%時，發(fā)布黃色預(yù)警；當(dāng)輸出概率P≥60%且P＜80%時，發(fā)布橙色預(yù)警；當(dāng)輸出概率P≥80%時，發(fā)布紅色預(yù)警，見表4。

表4 預(yù)警等級劃分Table4 Early warning level division

3 結(jié)論

（1）以清洗后的2010—2018年滑坡災(zāi)害為正樣本，以1∶2的比例采樣負(fù)樣本，以200 m×200 m為網(wǎng)格單元，構(gòu)建了青川縣區(qū)域滑坡災(zāi)害訓(xùn)練樣本集。訓(xùn)練樣本數(shù)量為1 826個，其中正樣本為613個，負(fù)樣本為1 213個。各訓(xùn)練樣本包含了27個輸入特征參數(shù)（地質(zhì)環(huán)境因子、降雨因子等），1個輸出特征參數(shù)（滑坡正樣本為1，非滑坡負(fù)樣本為0）。

（2）基于1 826個訓(xùn)練樣本，采用邏輯回歸算法開展滑坡災(zāi)害預(yù)警模型學(xué)習(xí)訓(xùn)練。以訓(xùn)練樣本集的80%作為訓(xùn)練集，20%作為測試集，進(jìn)行5 折交叉驗證，采用精確度（Accuracy）、ROC曲線和AUC值校驗?zāi)Ｐ蜏?zhǔn)確度和模型泛化能力。采用貝葉斯優(yōu)化算法對模型參數(shù)進(jìn)行優(yōu)化。結(jié)果顯示，設(shè)置閾值為0.5時，混淆矩陣的總精度最高，模型準(zhǔn)確率為0.943，AUC值為0.980，模型準(zhǔn)確率和模型泛化能力均較好。

（3）開展區(qū)域滑坡災(zāi)害實際預(yù)警時，可調(diào)用基于邏輯回歸算法訓(xùn)練好的模型輸出概率，根據(jù)概率分段確定預(yù)警等級，輸出概率P≥40%且P＜60%時，發(fā)布黃色預(yù)警；輸出概率P≥60%且P＜80%時，發(fā)布橙色預(yù)警；當(dāng)輸出概率P≥80%時，發(fā)布紅色預(yù)警。

今后，將在青川縣區(qū)域滑坡災(zāi)害預(yù)警中進(jìn)一步校驗?zāi)Ｐ偷臏?zhǔn)確性。