基于CF與CF-LR模型的地質(zhì)災害易發(fā)性評價

屠水云，張鐘遠，付弘流，徐世光，鄧明國，何例春，劉金宇

（1.云南地礦工程勘察集團公司,云南 昆明 650000；2.昆明理工大學國土資源工程學院,云南 昆明 650000；3.銅仁市自然資源局,貴州 銅仁 554300）

0 引言

地質(zhì)災害易發(fā)性評價是以地質(zhì)環(huán)境條件為基礎，參考地質(zhì)災害現(xiàn)狀的靜態(tài)因素來預測一定區(qū)域內(nèi)發(fā)生地質(zhì)災害的可能性[1]。地質(zhì)災害易發(fā)性評價方法分為定性和定量兩類。定性方法主要包括專家評分[2]、層次分析[3]等。隨著數(shù)據(jù)獲取的便利、計算能力的提升以及評估模型的日趨完善，定量評價方法應用更為廣泛，定量方法主要有信息量[4]、確定性系數(shù)[5]、證據(jù)權[6]、邏輯回歸[7]、支持向量機[8]、決策樹[9]、隨機森林[10]、神經(jīng)網(wǎng)絡[11]等。其中確定性系數(shù)方法計算嚴密，可以解決多源數(shù)據(jù)類型的合并問題和影響因子內(nèi)部不同特征區(qū)間對地質(zhì)災害易發(fā)性的影響[12]，但單一的確定性系數(shù)評價法沒有考慮每個評價因素對地質(zhì)災害易發(fā)性的影響差異。邏輯回歸( Logistic Regression，LR) 可以使用簡單的線性回歸來描述自然現(xiàn)象之間的復雜非線性關系，并根據(jù)影響因素與歷史災害點之間的關系確定影響因素的權重。文章基于地理信息系統(tǒng)，將研究區(qū)劃分為柵格，選取海拔、坡度、坡向、地形曲率、歸一化植被指數(shù)（Normalized Difference Vegetation Index，NDVI）、工程地質(zhì)巖組、斷層、道路、水系這9 個孕災、誘災因素作為評價指標因子，采用頻率比法（Frequency Ratio，F(xiàn)R）、確定性系數(shù)法（Certainty Factor，CF）量化評價指標因子，基于確定性系數(shù)法進行邏輯回歸運算，計算研究區(qū)網(wǎng)格地質(zhì)災害發(fā)生的概率，得到地質(zhì)災害易發(fā)性分區(qū)圖。

1 研究方法

1.1 頻率比（FR）

頻率比是建立在假設地質(zhì)條件、孕育地質(zhì)災害的概率相似的地區(qū)。頻率比重點考慮因子類別與地質(zhì)災害發(fā)生可能性的空間相關性，定量表示環(huán)境因子各屬性區(qū)間對地質(zhì)災害發(fā)生的相對影響程度[13?15]，計算方法如式（1）。

式中：FRi——頻率比值；

li——某個評價因子i類屬性區(qū)間發(fā)生地質(zhì)災害的個數(shù)；

L——研究區(qū)內(nèi)的總數(shù)；

si——某個評價因子i類區(qū)間的面積；

S——研究區(qū)總面積。

FRi大于 1 表明該環(huán)境因子屬性區(qū)間利于地質(zhì)災害發(fā)育，值越大表示對地質(zhì)災害發(fā)育的貢獻也越大；反之，F(xiàn)Ri小于 1 表明該環(huán)境因子屬性區(qū)間不利于地質(zhì)災害發(fā)育。

1.2 確定性系數(shù)模型（CF）

確定性系數(shù)模型假設將來發(fā)生地質(zhì)災害的條件和過去發(fā)生地質(zhì)災害的條件相同。CF計算公式為:

式中：CF——地質(zhì)災害發(fā)生的確定性系數(shù)；

PPa——地質(zhì)災害在因子分類數(shù)據(jù)a中發(fā)生的條件概率，研究中通常用因子分類a中的地質(zhì)災害個數(shù)與因子分類a的面積比值表示；

PPS——地質(zhì)災害在整個研究區(qū)中發(fā)生的先驗概率，以研究區(qū)地質(zhì)災害總個數(shù)與研究區(qū)總面積比值表示。

由式（2）可知，CF的變化區(qū)間為[?1，1]。正值表示地質(zhì)災害發(fā)生的確定性大，越接近1 越易于發(fā)生地質(zhì)災害；負值表示地質(zhì)災害發(fā)生確定性小，越接近?1 越不易于發(fā)生地質(zhì)災害；值為 0 時表示條件概率和先驗概率相同，不確定是否會發(fā)生地質(zhì)災害[5]。

1.3 基于確定性系數(shù)的邏輯回歸模型（CF-LR）

邏輯回歸模型是研究二分類因變量常用的多元統(tǒng)計分析方法。自變量Xi為控制災害發(fā)生的影響因子。因變量Y屬于二分類變量，通常 0 代表地質(zhì)災害不存在，1 代表地質(zhì)災害存在。用線性回歸來描述自然現(xiàn)象之間復雜的非線性關系，揭示因變量和多個自變量之間的多元回歸關系，將每個評價因子視為自變量，能很好解決滑坡易發(fā)性評價中出現(xiàn)的二分類變量問題[16]，邏輯回歸函數(shù)如式（3）：

式中：P——地質(zhì)災害發(fā)生的概率；

Z——地質(zhì)災害發(fā)生概率的目標函數(shù)，表達為各因素自變量x1，x2，x3,···,xn的線性組合；

β1,β2,···,βn——邏輯回歸系數(shù)；

β0——常數(shù)表示在不受任何有利或不利于地質(zhì)災害發(fā)生因素影響的條件下，地質(zhì)災害發(fā)生與不發(fā)生概率之比的對數(shù)值[17]。

通過確定性系數(shù)模型計算得到各評價因子類別的值，將其結(jié)果作為邏輯回歸模型中的自變量，建立回歸方程，進行邏輯回歸運算，得到各評價因子的邏輯回歸系數(shù)，以此進行確定性系數(shù)–邏輯回歸模型（CF-LR）進行地質(zhì)災害易發(fā)性評價。

2 實例分析

2.1 研究區(qū)概況

研究區(qū)沿河土家族自治縣位于貴州省東北部，隸屬銅仁市，南北長98.28 km，東西寬53 km，行政區(qū)域總面積2 483.51km2，占貴州省總面積的1.4%，占銅仁市總面積的13.7%。沿河縣境內(nèi)有烏江及其支流洪渡河、暗溪河、白泥河、壩坨河等26 條河流，河道長548.7 km，河網(wǎng)密度0.23 km/km2。地貌輪廓明顯受地質(zhì)構造控制，全縣地貌“軸部成山，翼部成谷”。區(qū)內(nèi)出露地層從老到新有震旦系、寒武系、奧陶系、志留系、二疊系、三疊系及第四系。受烏江切割和地層、巖性、構造的影響，在內(nèi)外營力綜合作用下，形成山巒疊障、溝谷縱橫、復雜多樣的地形地貌景觀。區(qū)內(nèi)歷史地質(zhì)災害以滑坡、崩塌為主，共計130 處，滑坡、崩塌分別占全縣地質(zhì)災害的55.38%、33.85%。研究區(qū)地理位置及地質(zhì)災害分布如圖1所示。

圖1 研究區(qū)地理位置及地質(zhì)災害點分布Fig.1 Geographical location and distribution of geological hazard in the study area

2.2 評價指標因子選取

結(jié)合研究區(qū)的地質(zhì)背景、地質(zhì)災害形成條件及發(fā)育特征，初步選取海拔、坡度、坡向、地形曲率、歸一化植被指數(shù)（NDVI）、工程地質(zhì)巖組、距斷層距離、距道路距離、距水系距離9 個影響因素作為評價指標因子。數(shù)據(jù)源為沿河縣地質(zhì)災害數(shù)據(jù)庫、地理空間數(shù)據(jù)云平臺獲取研究區(qū)30 m×30 m 數(shù)字高程模型（Digital Elevation Model，DEM)、1∶50 000 的地質(zhì)圖、Google 影像地圖，利用ArcGIS 平臺通過DEM 數(shù)據(jù)提取分析得到研究區(qū)坡度、坡向、地形曲率、河流網(wǎng)評價因子圖層，通過Google 影像地圖矢量化得到道路數(shù)據(jù)，利用 landsat8影像獲得該區(qū)的歸一化植被指數(shù)（NDVI）專題圖。

2.3 評價因子相關性分析

影響地質(zhì)災害發(fā)育的因素之間存在一定的關聯(lián)，當評價因子之間存在多重共線問題時，會降低模型的預測精度，因而需對評價因素進行相關性分析。利用ArcGIS 計算相關矩陣如表1所示，相關性系數(shù)絕對值最大為0.324，說明本文選取的9 個評價指標因子之間相關性較弱，均可納入研究區(qū)評價模型[18]。

表1 評價指標因子相關性系數(shù)矩陣Table 1 Correlation coefficient matrix of evaluation index factors

2.4 評價指標因子分析

工程地質(zhì)巖組為離散型因子，根據(jù)野外地質(zhì)調(diào)查以及已有分類標準進行分類，連續(xù)型指標因子分類根據(jù)地質(zhì)災害比例進行等距離劃分，各指標因子分級如圖2所示，利用式（1）進行頻率比計算確定性系數(shù)計算，利用式（2）進行確定性系數(shù)計算，結(jié)果見表2。

表2 評價指標因子分級、頻率比、確定性系數(shù)Table 2 Evaluation index factor classification,frequency ratio and certainty coefficient

圖2 評價指標因子分級圖Fig.2 Grading of evaluation index factors

海拔高度與降雨量、植被類型、植被覆蓋等有著密切的關系，影響著人類工程活動程度，因此海拔間接影響著地質(zhì)災害的發(fā)育[19]，海拔高度209～1 408 m，將其分為6 個類別。

坡度定量描述地面的傾斜程度，它的大小對斜坡表面徑流量、斜坡表體土層剩余下滑力等都影響巨大，一定程度上影響著地質(zhì)災害發(fā)育的規(guī)模與強度[20]，研究區(qū)坡度最高達75°，以8°等間距分為5 類，大于40°為1 類，共計6 個類別。

不同坡向與巖體結(jié)構面的組合關系差異導致地質(zhì)災害發(fā)育的程度不同[21]，將研究區(qū)坡向分為9 個類別。

地形曲率是局部地形曲面在各個截面方向上形狀、凹凸變化的反映，其值為正時表明邊坡是凸面坡，為 0 時表明為平面坡，為負時表明邊坡為凹面坡[22]，由于研究區(qū)平面坡（曲率等于0）面積極小，所以用曲率為?0.2～0.2 代表近似平面坡，將其分為凹坡（

歸一化植被指數(shù)（NDVI）是遙感影像中近紅外波段（NIR）的反射值和紅光波段（R）的反射值的差與兩者之和的比值，NDVI 值的范圍為 [?1,1]，負值表示對可見光高反射，地面為江、河、湖泊等水體或有雪覆蓋，0 表示NIR 和R 近似相等，為巖石或裸地等，正值表示有植被覆蓋，數(shù)值越大表示植被覆蓋率越高[23]，研究區(qū)NDVI 在?0.02～0.54 之間，將其分為5 個類別。

巖土體是地質(zhì)災害發(fā)生的物質(zhì)來源基礎，巖石類型、堅硬程度決定巖土體的力學強度、抗風化能力和抗侵蝕能力[19]，研究區(qū)工程地質(zhì)巖組分為5 類，分別為堅硬巖組、較堅硬巖組、較軟巖組、軟巖組和軟硬相間巖組。

地質(zhì)構造影響著巖體結(jié)構及其組合特征，對山區(qū)地質(zhì)災害發(fā)育起著重要的控制作用[24]，利用ArcGIS 領域分析功能將研究區(qū)斷層以300 m 等距離提取緩沖區(qū)，得到6 個類別。

道路修建開挖坡體改變原有地質(zhì)環(huán)境，破壞巖土體結(jié)構[25]，以200 m 等距離提取道路緩沖區(qū)，得到6 個類別。

河流的侵蝕、側(cè)蝕作用影響地質(zhì)災害的發(fā)育、且河流是控制坡面侵蝕的重要原因[26]，將研究區(qū)河流200 m等距離提取緩沖區(qū)，得到6 個類別。

通過對因子類別進行分類后，利用式（1）對各評價因子類別進行頻率比計算，當頻率比大于1 時，說明該因子類別對地質(zhì)災害發(fā)育具有促進作用，如表3所示。

表3 頻率比大于1 的屬性區(qū)間Table 3 Attribute intervals with frequency ratio greater than 1

2.5 邏輯回歸分析

利用ArcGIS 以500 m 距離制作災點緩沖區(qū)，在500 m 以外提取隨機點130 個非地質(zhì)災害點，與災害訓練樣本組成訓練集共計260 個點。將9 個評價指標因子的屬性提取至訓練集樣本，導出后替換成評價因子的CF 值導入SPSS 軟件中進行邏輯回歸運算，各評價因子分類級別的CF 值作為自變量，是否發(fā)生滑坡災害作為因變量（0 表示未發(fā)生地質(zhì)災害，1 值表示已發(fā)生地質(zhì)災害），LR-CF 模型的邏輯回歸運算結(jié)果如表4所示，其計算得到的所有評價指標因子的邏輯回歸系數(shù)均為正數(shù)，表明所有評價指標因子對模型均起正向作用。在邏輯回歸計算過程中，顯著性sig≤ 0.05 則該回歸系數(shù)有效，評價指標因子具有統(tǒng)計意義[22]。

表4 邏輯回歸系數(shù)和顯著性Table 4 Logistic regression coefficient and significance

3 易發(fā)性分區(qū)與評價

3.1 CF、CF-LR 分區(qū)結(jié)果

基于GIS 平臺，將評價指標因子圖層自定義添加屬性字段，對應輸入計算的確定性系數(shù)，利用柵格疊加得到確定性系數(shù)模型評價圖，利用自然斷點法將沿河縣地質(zhì)災害易發(fā)性區(qū)劃為低易發(fā)區(qū)、中易發(fā)區(qū)、高易發(fā)區(qū)、極高易發(fā)區(qū)，其面積（頻率比）分別為361.265 km2（0.159）、784.269 km2（0.414）、895.197 km2（1.003）、442.779 km2（2.718），如圖3（a）和表5所示。利用柵格計算器按照公式（3）計算得到CF-LR 模型地質(zhì)災害發(fā)生概率圖，利用自然斷點法將其分為低易發(fā)區(qū)、中易發(fā)區(qū)、高易發(fā)區(qū)、極高易發(fā)區(qū)，其面積（頻率比）分別為671.252 km2（0.142）、467.758 km2（0.327）、927.527 km2（0.741）、507.145 km2（3.051），如圖3（b）和表5所示。CF 模型和CF-LR 模型地質(zhì)災害易發(fā)性等級的頻率比值均從極低易發(fā)區(qū)到極高易發(fā)區(qū)明顯增大，表明有效評價了研究區(qū)地質(zhì)災害易發(fā)性。CF 模型和CF-LR 模型計算的極高易發(fā)區(qū)頻率比值分別占總頻率比值為63.3%和71.6%。說明CF-LR 模型比單一CF 模型評價精度更高。

表5 地質(zhì)災害易發(fā)性評價頻率比值Table 5 Frequency ratio of geological hazard susceptibility evaluation

圖3 地質(zhì)災害易發(fā)性區(qū)劃Fig.3 Division of geological hazard susceptibility

3.2 精度驗證

本文使用ROC 曲線來表示擬合數(shù)據(jù)和實測數(shù)據(jù)之間的關系，評價成功率或有效性以AUC值來表示（圖4）。曲線中縱軸為敏感度，即實際地質(zhì)災害數(shù)量百分比累加量，橫軸為特異性，即易發(fā)性面積百分比累積量，ROC曲線下面積AUC值越大表明模型評估效果越好[27?28]。CF 模型和CF-LR 模型AUC值分別為0.722 和0.818，說明CF 和CF-LR 評價模型均能夠較為客觀準確地對沿河縣地質(zhì)災害易發(fā)性進行評價，且CF 法進行邏輯回歸后的CF-LR 模型評價精度更高。

圖4 ROC 曲線Fig.4 ROC curve

4 結(jié)論

（1）文中從選取的9 個地質(zhì)災害影響因素的各類別的頻率比值可以看出，在海拔209～800 m，坡度8°～32°，坡向朝向北、東北、東、西南，地形曲率小于?0.2，NDVI 為0.1～0.3，較軟質(zhì)巖、軟質(zhì)巖、軟硬相間質(zhì)巖，距斷層900 m、道路和河流800 m 以內(nèi)對沿河縣地質(zhì)災害發(fā)育具有促進作用。

（2）CF 模型評價低易發(fā)區(qū)、中易發(fā)區(qū)、高易發(fā)區(qū)、極高易發(fā)區(qū)，其面積（頻率比）分別為361.265 km2（0.159）、784.269 km2（0.414）、895.197 km2（1.003）、442.779 km2（2.718）；CF-LR 模型評價低易發(fā)區(qū)、中易發(fā)區(qū)、高易發(fā)區(qū)、極高易發(fā)區(qū)，其面積（頻率比）分別為671.252 km2（0.142）、467.758 km2（0.327）、927.527 km2（0.741）、507.145 km2（3.051）。CF 模型和CF-LR 模型地質(zhì)災害易發(fā)性等級的頻率比值從極低易發(fā)區(qū)到極高易發(fā)區(qū)明顯增大，均有效評價了研究區(qū)地質(zhì)災害易發(fā)性。CF 模型和CF-LR 模型計算的極高易發(fā)區(qū)頻率比值分別占總頻率比值為63.3%和71.6%。

（3）CF 模型和CF-LR 模型AUC值分別為0.722 和0.818，均能夠較為客觀準確地對沿河縣地質(zhì)災害易發(fā)性進行評價。但單一CF 法沒有考慮評價因素對地質(zhì)災害易發(fā)性的影響差異，在此基礎上，LR 法用線性回歸來表示評價因子之間復雜非線性關系，考慮了評價因子的權重，使得AUC值提高了0.096，CF-LR 模型具有更高的評價精度。

由于研究區(qū)的地質(zhì)災害研究樣本偏少，不為理想研究實驗區(qū)，將影響評價效果和精度，對地質(zhì)災害易發(fā)性評價的精度還需進一步探索。