屠水云，張鐘遠(yuǎn)，付弘流，徐世光，鄧明國(guó)，何例春，劉金宇

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

0 引言

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

1 研究方法

1.1 頻率比（FR）

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

式中：FRi——頻率比值；

li——某個(gè)評(píng)價(jià)因子i類(lèi)屬性區(qū)間發(fā)生地質(zhì)災(zāi)害的個(gè)數(shù)；

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

si——某個(gè)評(píng)價(jià)因子i類(lèi)區(qū)間的面積；

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

2 實(shí)例分析

2.1 研究區(qū)概況

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

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

2.2 評(píng)價(jià)指標(biāo)因子選取

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

2.3 評(píng)價(jià)因子相關(guān)性分析

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

表1 評(píng)價(jià)指標(biāo)因子相關(guān)性系數(shù)矩陣Table 1 Correlation coefficient matrix of evaluation index factors

2.4 評(píng)價(jià)指標(biāo)因子分析

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

表2 評(píng)價(jià)指標(biāo)因子分級(jí)、頻率比、確定性系數(shù)Table 2 Evaluation index factor classification,frequency ratio and certainty coefficient

圖2 評(píng)價(jià)指標(biāo)因子分級(jí)圖Fig.2 Grading of evaluation index factors

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

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

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

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

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

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

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

道路修建開(kāi)挖坡體改變?cè)械刭|(zhì)環(huán)境，破壞巖土體結(jié)構(gòu)[25]，以200 m 等距離提取道路緩沖區(qū)，得到6 個(gè)類(lèi)別。

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

通過(guò)對(duì)因子類(lèi)別進(jìn)行分類(lèi)后，利用式（1）對(duì)各評(píng)價(jià)因子類(lèi)別進(jìn)行頻率比計(jì)算，當(dāng)頻率比大于1 時(shí)，說(shuō)明該因子類(lèi)別對(duì)地質(zhì)災(zāi)害發(fā)育具有促進(jìn)作用，如表3所示。

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

2.5 邏輯回歸分析

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

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

3 易發(fā)性分區(qū)與評(píng)價(jià)

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

基于GIS 平臺(tái)，將評(píng)價(jià)指標(biāo)因子圖層自定義添加屬性字段，對(duì)應(yīng)輸入計(jì)算的確定性系數(shù)，利用柵格疊加得到確定性系數(shù)模型評(píng)價(jià)圖，利用自然斷點(diǎn)法將沿河縣地質(zhì)災(zāi)害易發(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所示。利用柵格計(jì)算器按照公式（3）計(jì)算得到CF-LR 模型地質(zhì)災(zāi)害發(fā)生概率圖，利用自然斷點(diǎn)法將其分為低易發(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ì)災(zāi)害易發(fā)性等級(jí)的頻率比值均從極低易發(fā)區(qū)到極高易發(fā)區(qū)明顯增大，表明有效評(píng)價(jià)了研究區(qū)地質(zhì)災(zāi)害易發(fā)性。CF 模型和CF-LR 模型計(jì)算的極高易發(fā)區(qū)頻率比值分別占總頻率比值為63.3%和71.6%。說(shuō)明CF-LR 模型比單一CF 模型評(píng)價(jià)精度更高。

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

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

3.2 精度驗(yàn)證

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

圖4 ROC 曲線Fig.4 ROC curve

4 結(jié)論

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

（2）CF 模型評(píng)價(jià)低易發(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 模型評(píng)價(jià)低易發(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ì)災(zāi)害易發(fā)性等級(jí)的頻率比值從極低易發(fā)區(qū)到極高易發(fā)區(qū)明顯增大，均有效評(píng)價(jià)了研究區(qū)地質(zhì)災(zāi)害易發(fā)性。CF 模型和CF-LR 模型計(jì)算的極高易發(fā)區(qū)頻率比值分別占總頻率比值為63.3%和71.6%。

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

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