穆 柯,謝婉麗,劉琦琦,嚴(yán) 明,楊 惠,李嘉昊,黃 煜,朱榮森
(1.西北大學(xué) 地質(zhì)學(xué)系,陜西 西安 710069;2.西北大學(xué) 大陸動力學(xué)國家重點實驗室,陜西 西安 710069)
我國疆域幅員遼闊,山地眾多,地形地貌復(fù)雜,長期受到滑坡災(zāi)害的威脅。據(jù)統(tǒng)計,2021年全國共發(fā)生滑坡2 335起,占全年發(fā)生地質(zhì)災(zāi)害總數(shù)的48.93%,為所有地災(zāi)種類之最。
耀州區(qū)是陜西省中部連接關(guān)中與陜北的交通節(jié)點,也是重要的礦業(yè)基地,長期以來對資源的不合理開發(fā)導(dǎo)致滑坡災(zāi)害頻發(fā),阻礙著當(dāng)?shù)厣a(chǎn)建設(shè)與城市發(fā)展。針對滑坡的突發(fā)性與高危害性,科學(xué)的風(fēng)險評價與防災(zāi)管理是減輕其損失的最佳方案,而對滑坡的易發(fā)性進(jìn)行精準(zhǔn)的劃分是滑坡防治與管理工作的重要一環(huán)。
經(jīng)過國內(nèi)外學(xué)者數(shù)十年來的研究,關(guān)于滑坡易發(fā)性研究的理論與技術(shù)日新月異,目前廣泛使用的方法主要為三種類型:層次分析(AHP)法[1]、加權(quán)線性組合法[2]等定性分析方法;信息量[3]、確定性系數(shù)(CF)[4]、模糊信息[5]等統(tǒng)計預(yù)測模型方法,BP神經(jīng)網(wǎng)絡(luò)[6]、邏輯回歸(LR)[7]、支持向量機(jī)(SVM)[8-9]、隨機(jī)森林(RF)[10]等機(jī)器學(xué)習(xí)模型。實例如許沖等基于層次分析法對汶川震區(qū)的滑坡進(jìn)行了易發(fā)性評價[11];謝婉麗等引入模糊信息法在大西安地區(qū)進(jìn)行了地質(zhì)災(zāi)害易發(fā)性區(qū)劃[12-13];NHU等人采用邏輯模型樹等五種機(jī)器學(xué)習(xí)方法對伊朗庫爾德斯坦省的滑坡進(jìn)行易發(fā)性預(yù)測,并對結(jié)果進(jìn)行分析比較,最終認(rèn)為邏輯模型樹的評價精度最高[14];BEHNIA等人開發(fā)了數(shù)據(jù)驅(qū)動的隨機(jī)森林模型對育空地區(qū)的滑坡進(jìn)行了易發(fā)性評價,認(rèn)為隨機(jī)森林在滑坡易發(fā)性評價中表現(xiàn)出色[15]。
近年來,隨著滑坡易發(fā)性研究的深入與實際工作的展開,傳統(tǒng)評價方法的精度愈發(fā)難以滿足實際需求,越來越多的學(xué)者將目光轉(zhuǎn)向耦合模型,選擇合理的預(yù)測模型以提升預(yù)測精度是近幾年滑坡易發(fā)性研究的熱點領(lǐng)域。如BINH等人提出了一種徑向基函數(shù)與旋轉(zhuǎn)森林集成的耦合模型(RFRBF),并在實際案例中與其他五種機(jī)器學(xué)習(xí)模型進(jìn)行對比,認(rèn)為RFRBF的預(yù)測精度最高[16]。ZHENG等人在對日本山區(qū)的滑坡易發(fā)性評估中對比了四種集成學(xué)習(xí)算法(ML),其中SVM-boosting表現(xiàn)最佳,SVM-stacking的性能最低,表明耦合模型不一定會提升預(yù)測性能[17]。白志剛等人基于熵指數(shù)與隨機(jī)森林耦合模型對渝東北地區(qū)地質(zhì)災(zāi)害易發(fā)性評價,認(rèn)為耦合模型的加入明顯提升了預(yù)測精度[18]。上述案例表明不同預(yù)測模型的耦合不一定會提升滑坡易發(fā)性評價的精度,在實際應(yīng)用中的效果仍然存在爭議,尚未形成完整的理論體系?;谏鲜鲅芯楷F(xiàn)狀,本文嘗試將LR-RF耦合模型應(yīng)用于耀州區(qū)的滑坡易發(fā)性評價,用于檢驗與傳統(tǒng)LR模型相比,耦合模型在實際評價中的預(yù)測精度是否有所提升。
邏輯回歸(LR)模型是廣泛應(yīng)用于地質(zhì)災(zāi)害領(lǐng)域的回歸分析模型,具有假設(shè)簡單、數(shù)據(jù)要求低、輸出結(jié)果便捷等優(yōu)點,廣泛應(yīng)用于統(tǒng)計建模領(lǐng)域[19],但是需要大量而準(zhǔn)確的數(shù)據(jù)對數(shù)學(xué)模型的構(gòu)建進(jìn)行支撐,當(dāng)數(shù)據(jù)較少時,評價結(jié)果不準(zhǔn)確[20]。隨機(jī)森林(RF)模型是一種基于集成學(xué)習(xí)的算法,它通過隨機(jī)抽樣對多個決策樹進(jìn)行集成,使用投票機(jī)制進(jìn)行預(yù)測[21]。因此為了在較少量數(shù)據(jù)的前提下提高模型的預(yù)測精度,本研究將LR模型的空間屬性引入隨機(jī)森林模型中,以耀州區(qū)為研究區(qū),選取坡度等八個環(huán)境因子,建立二者的耦合模型(LR-RF)對研究區(qū)進(jìn)行滑坡易發(fā)性評價,并與傳統(tǒng)LR模型評價結(jié)果進(jìn)行對比,以探討LR-RF耦合模型的預(yù)測性能及適用性,研究成果可以為當(dāng)?shù)氐幕路乐闻c風(fēng)險管理工作提供理論支撐。
耀州區(qū)位于陜西省中部黃土高原與汾渭地塹的過渡地帶,處于108.579 0°~109.092 2°E,34.805 7°~35.220 8°N之間,面積1 617 km2,海拔236~1 732 m,氣候類型屬溫帶季風(fēng)氣候,年均降水量554.5 mm。區(qū)內(nèi)發(fā)育漆水河等五條主要河流,強(qiáng)烈的切割侵蝕作用塑造了區(qū)內(nèi)殘塬溝壑密布的地貌特征,全區(qū)大部分斜坡表面均為第四系黃土。
區(qū)內(nèi)共有滑坡災(zāi)點71處,其分布現(xiàn)狀見圖1。由圖1不難看出,研究區(qū)災(zāi)點多分布于河谷階地區(qū)與殘塬溝壑區(qū);按行政區(qū)劃則在天寶路街道最多。區(qū)內(nèi)復(fù)雜多變的地貌特征為滑坡發(fā)生提供了地形基礎(chǔ),遍布全區(qū)的濕陷性黃土等脆弱巖土體提供了物源條件,人類工程活動是重要的誘發(fā)因素[22-24]。
圖1 滑坡點分布圖(審圖號:陜S(2021)023號,底圖無修改,下同)
耀州區(qū)滑坡點數(shù)據(jù)來源于實地調(diào)查,地理環(huán)境等數(shù)據(jù)來源見表1。
表1 評價指標(biāo)數(shù)據(jù)來源
根據(jù)研究區(qū)滑坡分布特征與前人經(jīng)驗,本次易發(fā)性評價選擇30 m×30 m的柵格單元作為基本評價單元,在ArcGIS中將研究區(qū)劃分為大約181萬個柵格單元。
評價因子的選擇對評價分析至關(guān)重要,綜合分析研究區(qū)滑坡發(fā)育近況和影響要素,這里主要從地理環(huán)境、地質(zhì)環(huán)境、人類活動三個方面出發(fā):
2.2.1 地理環(huán)境因子
(1)高程。研究區(qū)高程范圍為536 m~1 732 m,高差較大,不同高程區(qū)間對滑坡的作用程度有所不同,在ArcGIS中將其分為7級(圖2a)。
(2)坡度。研究區(qū)坡度分布范圍為0°~73°,坡度跨度區(qū)間較大,在ArcGIS中將其分為7級(圖2b)。
(3)坡向。研究區(qū)的坡向為-1°~360°,以45°為間隔對其進(jìn)行等距分級(圖2c)。
(4)距河流距離。距河流距離的大小對滑坡的發(fā)生影響較大,區(qū)內(nèi)五條河流及其支流均為渭河水系,多為NW-SE向。使用ArcGIS的歐式距離分析,以500 m為間隔將其劃分為5級(圖2d)。
(5)植被指數(shù)(NDVI)。植被指數(shù)由landsat8遙感影像通過軟件進(jìn)行地理校正、大氣輻射校正后計算獲取,其取值在[-1,1],NDVI與植被發(fā)育程度呈正相關(guān)。在ArcGIS中將其分為5級(圖2e)。
2.2.2 地質(zhì)環(huán)境因子
(1)工程地質(zhì)巖組。不同地層巖性的物理力學(xué)性質(zhì)對滑坡的誘發(fā)程度有顯著差異,強(qiáng)度較小的巖層對滑坡的誘發(fā)程度較高。根據(jù)研究區(qū)的地質(zhì)資料,研究區(qū)的地層巖性可分為五類巖組(圖2f)。
(2)距構(gòu)造距離。研究區(qū)地質(zhì)構(gòu)造主要為北部山區(qū)的三條背斜。使用ArcGIS的歐氏距離分析,以500 m為間隔將其劃分為5級(圖2g)。
2.2.3 人類活動因子
(1)距道路距離。距道路距離的遠(yuǎn)近會對滑坡產(chǎn)生不同的影響,利用ArcGIS的歐氏距離算法,以200 m為間距對其進(jìn)行分級(圖2h)。
圖2 評價因子分級圖
表2 因子間相關(guān)系數(shù)表
進(jìn)行評價前需要評估選取的評價因子的相關(guān)性,以保證因子間相互獨立不互相干涉。借助皮爾森系數(shù)(PCC)對各個因子進(jìn)行獨立性檢驗[25]。假設(shè)樣本數(shù)據(jù)集(Ai,Bj)=(a1,b1),(a2,b2),…,(an,bn),則PCC計算公式如下:
(1)
當(dāng)0≤|PCC|<0.3時,表明因子間線性不相關(guān);當(dāng)|PCC|≥0.3時,表明因子間線性相關(guān)。計算結(jié)果顯示評價因子間的|PCC|均比0.3小,證明各因子之間線性不相關(guān)。
邏輯回歸(LR)模型可以在多組數(shù)據(jù)中構(gòu)建目標(biāo)數(shù)據(jù)的回歸關(guān)系,從而對目標(biāo)數(shù)據(jù)發(fā)生的概率進(jìn)行預(yù)測[26]。在滑坡易發(fā)性評價中,對滑坡發(fā)生賦值為1,滑坡未發(fā)生賦值為0,并將其作為邏輯回歸模型的目標(biāo)變量。記滑坡發(fā)生的概率為P,其計算公式為:
(2)
式中:β0,β1,β2,…,βn為邏輯回歸系數(shù);X1,X2,…,Xn為自變量;P為評估結(jié)果,即發(fā)生滑坡的概率。
隨機(jī)森林(RF)模型是一種最早由LEO Breiman提出的集成學(xué)習(xí)算法。其工作原理為:首先從建模數(shù)據(jù)中隨機(jī)抽取部分樣本,并生成對應(yīng)的樣本集和決策樹;接著借助決策樹對抽取樣本集進(jìn)行訓(xùn)練,每棵決策樹都能獨立地得出預(yù)測結(jié)論;隨后匯總?cè)繘Q策樹,組成RF算法,對新數(shù)據(jù)進(jìn)行分類和預(yù)測。模型輸出結(jié)果由全部決策樹的獨立結(jié)果投票產(chǎn)生,獲得票數(shù)最多的即為預(yù)測結(jié)果[27],其運行結(jié)構(gòu)見圖3。
在進(jìn)行LR模型構(gòu)建之前,首先需要選取樣本數(shù)據(jù)構(gòu)建正負(fù)樣本集,正樣本為研究區(qū)的滑坡點,即滑坡發(fā)生的樣本點;負(fù)樣本為非滑坡點,即滑坡不發(fā)生的樣本點。
采樣方法為:以研究區(qū)的71個滑坡點為建模正樣本,在其周圍建立1 000 m的緩沖區(qū),在緩沖區(qū)外以1:3的比例隨機(jī)生成213個點作為非滑坡點(負(fù)樣本)。71個正樣本與213個負(fù)樣本共同組成訓(xùn)練樣本,標(biāo)記滑坡點為1,非滑坡點為0。使用各因子的歸一化值作為模型指標(biāo):
(3)
(4)
式中:Sij為i因子j分級下的滑坡災(zāi)害影響面積,Sij′為該分級的總面積,Rij為求取歸一化指數(shù)Xn的中間變量,Rn為Rij的合集,Xn為歸一化指數(shù)。
將數(shù)據(jù)樣本導(dǎo)入R語言軟件中,使用glm函數(shù)構(gòu)建LR模型,運算結(jié)果見表3。
表3 LR模型詳情表
由表3可知,各因子的顯著系數(shù)值均小于0.05,滿足獨立性條件,其預(yù)測發(fā)生滑坡的概率P為:
(5)
式中:X1j,X2j,X3j,…,X8j為單元格對應(yīng)的高程、坡度等8個評價因子通過式(3)和式(4)計算的歸一化值。在ArcGIS中完成式(5)的計算,得到研究區(qū)的滑坡發(fā)生概率,將其劃分為五級,得到基于LR模型的滑坡易發(fā)性分區(qū)圖(圖4)。
圖4 基于LR模型的滑坡易發(fā)性分區(qū)圖
(1)樣本選取。與LR模型構(gòu)建類似,建立正負(fù)樣本集是RF模型構(gòu)建之前的首要步驟。因為對研究區(qū)進(jìn)行預(yù)測之前不能確定滑坡易發(fā)區(qū)的空間位置,為了避免影響模型預(yù)測精度,非滑坡點的選取需要避開滑坡易發(fā)區(qū)。在傳統(tǒng)的建模過程中,大多是通過在滑坡周圍建立緩沖區(qū),在緩沖區(qū)外選取非滑坡點作為負(fù)樣本[28]。鑒于滑坡發(fā)育的隨機(jī)性與地質(zhì)環(huán)境的復(fù)雜性,極大增加了負(fù)樣本的獲取難度。此處結(jié)合前文中獲得的LR模型易發(fā)性分區(qū),在其去除高易發(fā)區(qū)的剩余區(qū)域內(nèi)作為非滑坡點選取區(qū),以提升其為非滑坡點的概率,從而提升模型的預(yù)測精度。
根據(jù)前人經(jīng)驗,在構(gòu)建模型之前,采用1∶3的比例選取負(fù)樣本點[29]。在LR模型去除高易發(fā)區(qū)的區(qū)域選取非滑坡點213個,與滑坡點71個共同組合為模型的正負(fù)樣本集,隨后以7∶3的比例將284個正負(fù)樣本劃分為201個訓(xùn)練樣本數(shù)據(jù)與83個測試樣本數(shù)據(jù)。
(2)參數(shù)選取。隨機(jī)森林模型只有在最優(yōu)參數(shù)組合下,才能發(fā)揮最理想的性能,預(yù)測效果才能達(dá)到最佳。在模型建立過程中,對預(yù)測精度影響最大的參數(shù)有兩個:最大特征數(shù)(簡稱mtry),其含義為組成隨機(jī)森林的每棵決策樹建立過程中可以使用的最大特征數(shù)量,其取值與參與模型構(gòu)建的因子數(shù)量有關(guān),本文選取了8種因子參與模型構(gòu)建,故最大特征數(shù)取值范圍為1~8的自然數(shù);決策樹數(shù)量(簡稱ntree),即隨機(jī)森林算法中的決策樹數(shù)量[27],決策樹數(shù)量越多,模型誤差會趨于穩(wěn)定,計算量也會增加。因此,對選取最大特征數(shù)與決策樹數(shù)量的取值是建模的關(guān)鍵所在。
通過在R語言中進(jìn)行循環(huán)迭代,可得不同最大特征數(shù)下的OOB袋外誤差(圖5)。圖中OOB袋外誤差含義為隨機(jī)森林模型分類錯誤的樣本占總樣本的比例,其值越小,即表明模型預(yù)測效果越好。由圖5可知,當(dāng)最大特征數(shù)取值為4時,誤差最小,因此本次建模mtry取值為4。
圖5 特征數(shù)與袋外誤差關(guān)系圖
取最大特征數(shù)為4,在R語言中進(jìn)行迭代,獲得不同決策樹下的袋外誤差(圖6)。圖中1線是指在循環(huán)迭代中,樣本1(滑坡點)被誤分類為0(非滑坡點)的誤差,同理0線是指樣本0(非滑坡點)被誤分類為1(滑坡點)的誤差。由圖6可知,模型中三條曲線的預(yù)測誤差在決策樹數(shù)量取值為200時趨于穩(wěn)定,故在本此建模中ntree取值為200。
最終確定ntree取值為200,mtry取值為4進(jìn)行LR-RF模型的構(gòu)建。
圖6 OOB誤差迭代圖
(3)模型建立。使用上一步獲得的201個訓(xùn)練樣本數(shù)據(jù),在R語言中基于randomforest函數(shù)構(gòu)建RF模型,并將83個測試數(shù)據(jù)集代入以進(jìn)行精度測試。測試結(jié)果顯示,模型的預(yù)測準(zhǔn)確率為91.24%,kappa值為0.838 1,誤差為5.53%,表明模型預(yù)測精度較高[24]。混淆矩陣通過統(tǒng)計實際與模型預(yù)測結(jié)果的不一致性來評估模型精度,LR-RF模型的測試混淆矩陣如表4所示,在測試樣本中,62個非滑坡點預(yù)測成功60個,預(yù)測失敗2個,錯誤率為3.2%;21個滑坡點預(yù)測成功18個,預(yù)測失敗3個,成功率為85.71%。
表4 混淆矩陣
(4)模型預(yù)測。經(jīng)過上述模型參數(shù)優(yōu)化,選取mtry為4,ntree為200,建立LR-RF模型,并將研究區(qū)全部181萬個評價柵格數(shù)據(jù)導(dǎo)入訓(xùn)練好的模型中,輸出結(jié)果為模型預(yù)測研究區(qū)發(fā)生滑坡的概率。在ArcGIS中將其分為5級,繪制基于LR-RF模型的滑坡易發(fā)性分區(qū)圖(圖7)。
圖7 基于LR-RF模型的滑坡易發(fā)性分區(qū)圖
將耀州區(qū)的滑坡點與兩種模型的滑坡易發(fā)性評價結(jié)果進(jìn)行疊加,從災(zāi)點密度和ROC曲線兩個方面對LR-RF模型與LR模型的滑坡易發(fā)性評價結(jié)果進(jìn)行驗證。
滑坡密度是各滑坡易發(fā)等級滑坡數(shù)量與該等級面積的比值,可以直觀地反應(yīng)不同的滑坡密度的差異[25]。在上文中已經(jīng)獲得LR模型與LR-RF模型的易發(fā)性分區(qū)結(jié)果,將已發(fā)生的地質(zhì)災(zāi)害點與獲得的易發(fā)性分區(qū)進(jìn)行統(tǒng)計,統(tǒng)計結(jié)果見表5。從表5中可以清晰的看出,LR-RF模型的高易發(fā)區(qū)包含了42.35%的滑坡點,高于LR模型的35.21%;在滑坡密度方面,LR-RF模型的滑坡密度為0.57處/ km2,也高于LR模型的0.34處/km2。由此可以說明LR-RF模型的預(yù)測結(jié)果中易發(fā)性較高的地方滑坡點更密集,與實際災(zāi)點分布更貼合,有更高的預(yù)測成功率。
表5 滑坡易發(fā)性分區(qū)統(tǒng)計表
在研究區(qū)已發(fā)生的滑坡中隨機(jī)選擇30%作為精度測試樣本,并建立500 m的緩沖區(qū),隨后在緩沖區(qū)外隨機(jī)生成相同數(shù)量的非滑坡點,二者共同組成測試正負(fù)樣本,用于對上述兩種模型獲得的易發(fā)性評價結(jié)果精度的進(jìn)行對比分析。
ROC曲線是一種反映模型預(yù)測精度的經(jīng)典方法,其橫軸為特異性,縱軸為靈敏度,AUC值(ROC曲線下面積)的大小用于評判模型的預(yù)測精度,其值越接近1,表明測試模型預(yù)測精度越高[26]。LR-RF和LR模型的ROC曲線見圖8,LR-RF模型的AUC值為0.912 3,與LR模型的0.889 5有較大提升。這是因為LR模型在樣本數(shù)據(jù)較少時缺少有力的數(shù)據(jù)支撐,在局部預(yù)測不精確。LR-RF模型在選取非滑坡樣本時吸收了LR模型的空間屬性,同時RF模型在分類預(yù)測中可以很好地反應(yīng)樣本數(shù)據(jù)間的非線性關(guān)系,因而表現(xiàn)出較強(qiáng)的預(yù)測精度。而且LR-RF模型在滑坡密度方面也有很好地適配性,說明該模型從適應(yīng)性和準(zhǔn)確度上都很適合對研究區(qū)進(jìn)行滑坡易發(fā)性評價。
圖8 不同模型預(yù)測結(jié)果的ROC曲線
精確而及時的滑坡易發(fā)性評價在滑坡防治工作中有著至關(guān)重要的作用,也是經(jīng)久不衰的研究熱點。近年來集成學(xué)習(xí)算法在滑坡易發(fā)性評價研究中展示了豐富的發(fā)展?jié)摿Γ梢苑从郴屡c評價因子間的非線性關(guān)系。本項研究為了解決LR模型在數(shù)據(jù)源較少時的精度下降問題,在LR模型的評價基礎(chǔ)上引入集成學(xué)習(xí)算法中的RF模型進(jìn)行耦合來對研究區(qū)展開滑坡易發(fā)性評價。借助滑坡密度與ROC曲線可以看出,與傳統(tǒng)的LR模型相比,LR-RF模型可以有效減少偏差和錯誤分類,有效提升了預(yù)測準(zhǔn)確率與精度。
本文以耀州區(qū)為研究區(qū),分別建立了傳統(tǒng)邏輯回歸(LR)模型與邏輯回歸-隨機(jī)森林(LR-RF)的耦合模型對研究區(qū)進(jìn)行滑坡易發(fā)性評價,并與傳統(tǒng)LR模型評價結(jié)果進(jìn)行對比,對LR-RF耦合模型與LR模型在研究區(qū)滑坡易發(fā)性評價中的預(yù)測性能及適用性進(jìn)行了探討,主要結(jié)論如下:
(1)基于耀州區(qū)的地質(zhì)環(huán)境條件與滑坡分布特征,選取高程、坡度、坡向等八個地質(zhì)環(huán)境因子,采用皮爾森相關(guān)性系數(shù)驗證了各個評價因子間的獨立性,確立了八個評價因子構(gòu)成的研究區(qū)滑坡易發(fā)性評價體系。
(2)使用在邏輯回歸(LR)模型分區(qū)基礎(chǔ)上選取非滑坡點作為負(fù)樣本的采樣方法,在R語言中構(gòu)建耀州區(qū)LR-RF滑坡易發(fā)性評價模型,對該區(qū)進(jìn)行了滑坡易發(fā)性評價和區(qū)劃。評價結(jié)果表明,高、較高及中易發(fā)區(qū)覆蓋了研究區(qū)東南部的漆水河等河谷和中部殘塬邊緣區(qū)域,包含了73.23%的滑坡點,預(yù)測效果較好。
(3)分別從滑坡密度和ROC曲線兩項指標(biāo)對兩種模型結(jié)果精度進(jìn)行驗證:從滑坡密度來看,LR-RF模型高易發(fā)區(qū)的滑坡密度為0.57處/ km2,高于LR模型的0.34處/ km2;從ROC曲線來看,LR-RF模型的AUC值為0.912 3,也大于LR模型的0.889 5。LR-RF模型的兩項精度驗證指標(biāo)均優(yōu)于LR模型,表明耦合模型的預(yù)測結(jié)果優(yōu)于傳統(tǒng)模型。在研究區(qū)采用LR-RF模型進(jìn)行滑坡易發(fā)性評價擁有更好的評價精度和預(yù)測能力,可為研究區(qū)的防災(zāi)減災(zāi)工作提供一定理論參考和技術(shù)指導(dǎo)。