白光順，楊雪梅，朱杰勇，張世濤，祝傳兵，康曉波，孫 濱，周琰嵩

（1.昆明理工大學(xué)國土資源工程學(xué)院，云南 昆明 650093；2.自然資源部高原山地地質(zhì)災(zāi)害預(yù)報(bào)預(yù)警與生態(tài)保護(hù)修復(fù)重點(diǎn)實(shí)驗(yàn)室，云南 昆明 650093；3.云南高正地信科技有限公司，云南 昆明 650041；4.云南省地質(zhì)環(huán)境監(jiān)測院，云南 昆明 650216；5.云南南方地勘工程總公司昆明分公司，云南 昆明 650051）

0 引言

查明與地質(zhì)災(zāi)害有關(guān)的危險(xiǎn)區(qū)域是地質(zhì)災(zāi)害管理的重要工作，也是促進(jìn)研究區(qū)人民生活和基礎(chǔ)設(shè)施發(fā)展安全的重要依據(jù)[1]，基于建模評價(jià)地質(zhì)災(zāi)害易發(fā)性是重要而且有效的途徑。

應(yīng)用經(jīng)驗(yàn)式、數(shù)值模擬和統(tǒng)計(jì)方法對地質(zhì)災(zāi)害易發(fā)性建模和評價(jià)，已經(jīng)進(jìn)行了許多研究[1-10]。其中，經(jīng)驗(yàn)式方法基于現(xiàn)場觀察和專家經(jīng)驗(yàn)判斷；數(shù)值模擬計(jì)算邊坡的穩(wěn)定性；統(tǒng)計(jì)方法部分基于實(shí)地觀察和專家的先驗(yàn)知識，部分基于對地質(zhì)災(zāi)害發(fā)生的權(quán)重或概率的統(tǒng)計(jì)計(jì)算，這類方法使用統(tǒng)計(jì)技術(shù)來評估誘發(fā)地質(zhì)災(zāi)害的各種因素的相關(guān)作用，每個(gè)因素的重要性都是根據(jù)觀察到的與地質(zhì)災(zāi)害的關(guān)系來確定的。

文中使用基于貝葉斯理論的證據(jù)權(quán)法，綜合GIS 技術(shù)評價(jià)研究區(qū)地質(zhì)災(zāi)害易發(fā)性。證據(jù)權(quán)法是一種統(tǒng)計(jì)方法，最初應(yīng)用于非空間、定量的醫(yī)學(xué)診斷，以結(jié)合臨床診斷的證據(jù)來預(yù)測疾病[11-12]。在地球科學(xué)中，該方法被廣泛應(yīng)用，如：礦產(chǎn)資源潛力評估和礦床預(yù)測[13-16]，公路路基巖溶塌陷危險(xiǎn)性評價(jià)[17]和滑坡易發(fā)性和危險(xiǎn)性[1,3,18-23]。

1 研究區(qū)和數(shù)據(jù)

1.1 研究區(qū)概況

文中選擇云南高原滇中昆明盆地低山丘陵地帶這一云南省地質(zhì)災(zāi)害防治重點(diǎn)地區(qū)的典型代表，云南省省會昆明市的主要行政區(qū)之一，昆明市五華區(qū)作為研究對象，該區(qū)地質(zhì)災(zāi)害易發(fā)性評價(jià)研究具有典型代表性，可向整個(gè)云南高原昆明盆地低山丘陵區(qū)和其他低山丘陵區(qū)推廣，具有技術(shù)方法和社會經(jīng)濟(jì)意義。研究區(qū)面積381.6 km2，地勢西北高東南低，昆明盆地內(nèi)地形開闊低緩，北部山區(qū)地形崎嶇，溝壑較發(fā)育。區(qū)域年降水量的80%以上集中在6—9月，年平均降水量608.4～887.0 mm。碳酸鹽巖分布最廣，約占全區(qū)面積的38.93%，其次為砂巖、泥巖、頁巖，約占23.11%，巖漿巖主要為玄武巖，約占16.95%，主要分布在昆明盆地和其他小盆地的松散碎石土體約占11.36%，石英砂巖類約占7.56%，還發(fā)育一些巖脈；斷裂構(gòu)造較發(fā)育，以南北向構(gòu)造為主[24-25]。

1.2 數(shù)據(jù)收集

通過地質(zhì)災(zāi)害風(fēng)險(xiǎn)普查獲得了研究區(qū)地質(zhì)災(zāi)害分布數(shù)據(jù)。根據(jù)調(diào)查分析，選擇工程地質(zhì)巖組、斷裂構(gòu)造、高程、坡度、坡向、坡面曲率、距公路距離和土地利用類型等8 類因素納入評價(jià)分析。地質(zhì)數(shù)據(jù)收集自云南省地質(zhì)局1∶20 萬昆明幅、武定幅區(qū)域地質(zhì)調(diào)查報(bào)告和圖件[24-25]，12.5 m 分辨率DEM（數(shù)字高程模型）收集自ASF，道路數(shù)據(jù)收集自O(shè)SM，土地利用類型數(shù)據(jù)收集自ESA（圖1、表1）。

表1 數(shù)據(jù)簡介Table 1 Data introduction

圖1 因素基礎(chǔ)數(shù)據(jù)圖Fig.1 Basic data charts of factors

1.3 地質(zhì)災(zāi)害發(fā)育特征

現(xiàn)狀發(fā)育地質(zhì)災(zāi)害89 處，滑坡73 處，崩塌11 處，泥石流4 條，地面沉降1 處，為小—中型，無大型，中型14 處，小型75 處，主要分布在研究區(qū)低山丘陵地貌區(qū)，盆地內(nèi)僅發(fā)育1 處（圖2）。

圖2 地質(zhì)災(zāi)害分布圖（底圖為高程和山體陰影渲染）Fig.2 Map of geological hazard distribution (The bottom was rendered by elevation and hillshade)

選擇指標(biāo)“因子面積百分比A”“地災(zāi)數(shù)百分比B”和“比率(β=B/A)”表征地質(zhì)災(zāi)害的空間分布特征、主控因素和成災(zāi)特征。β 定義了地質(zhì)災(zāi)害點(diǎn)在因素分級中相對于均勻分布的豐度，β＞1 表示相對豐度更高，β＜1則相反。β＞1 的因素分級有（圖3、表2）：高程1 800～1 850 m、1 920～1 950 m 和1 950～2 000 m，坡度15°～25°、25°～35°和＞35°，坡向北東、東、南東和北，坡面曲率-0.75～-0.28（凹形）、-0.28～-0.15（凹形）、-0.15～-0.05(凹形)和0.05～0.15（凸形），石英砂巖巖組和砂巖、泥巖、頁巖巖組，距斷層距離0～50 m、300～500 m 和1 000～2 000 m，距主要公路距離0～50 m 和50～100 m，草地和裸地/稀疏植被區(qū)域。這些因素分級內(nèi)，發(fā)育了相對于均勻分布豐度更高的地質(zhì)災(zāi)害，表征這些因素分級可能是研究區(qū)地質(zhì)災(zāi)害的主控因素。

圖3 各因素分級分區(qū)和地災(zāi)點(diǎn)數(shù)量相關(guān)性統(tǒng)計(jì)圖Fig.3 Statistical charts of correlation between the factors and the number of geological hazard points

2 方法

2.1 權(quán)重計(jì)算

把研究區(qū)柵格單元化，利用條件概率計(jì)算證據(jù)因素圖層所有單元對地質(zhì)災(zāi)害發(fā)生的貢獻(xiàn)權(quán)重[13-15,26-27]。定義D為已發(fā)生地質(zhì)災(zāi)害的單元，為未發(fā)生地質(zhì)災(zāi)害的單元，B為證據(jù)因素區(qū)內(nèi)的單元，為證據(jù)因素區(qū)外的單元。

證據(jù)因素B條件下D的條件（后驗(yàn)）概率為：

式中：O(D)——證據(jù)因素B 的先驗(yàn)概率，O(D)=

P(B|D)、P(B|)——在地質(zhì)災(zāi)害發(fā)生（D）和未發(fā)生（）時(shí)，證據(jù)因素B的條件概率，取自然對數(shù)即是證據(jù)權(quán)法中的正權(quán)重（證據(jù)因素存在區(qū)的權(quán)重值）W+。

用D和B的單元數(shù)N可表示為：

同式（3）—（6）：

根據(jù)式（11）和（12），使用ArcGIS 空間分析工具執(zhí)行權(quán)重W+和W-計(jì)算。

W+的大小表明證據(jù)因素的存在與地質(zhì)災(zāi)害發(fā)生之間存在正相關(guān)關(guān)系。W-表示負(fù)相關(guān)，即證據(jù)因素存在抑制誘發(fā)地質(zhì)災(zāi)害的作用。證據(jù)因素原始數(shù)據(jù)缺失區(qū)域的權(quán)重值取0。兩個(gè)權(quán)重之間的差異Wf=W+-W-，即綜合權(quán)重，量化證據(jù)因素和地質(zhì)災(zāi)害相關(guān)性大小。如果Wf為正，則證據(jù)因素對地質(zhì)災(zāi)害有利，如果為負(fù)，則對滑坡不利。如果Wf接近于零，則表明證據(jù)因素與地質(zhì)災(zāi)害的相關(guān)性不大。

2.2 因素分級（分類）的顯著性測試

在上述權(quán)重值計(jì)算及分析的基礎(chǔ)上，實(shí)施證據(jù)因素分類的優(yōu)選，選擇類間差異顯著的證據(jù)因素類，歸并不顯著的證據(jù)因素類。選擇近似學(xué)生化檢驗(yàn)（Student-T）統(tǒng)計(jì)值進(jìn)行顯著性測試[15,28]：

式中：σW+、σW-——分別是W+和W-的標(biāo)準(zhǔn)差；

Wf——綜合權(quán)重；

σWf——綜合權(quán)重標(biāo)準(zhǔn)差。

當(dāng)測試值的絕對值 |Student-T|為1.96 和2.326 時(shí)，置信度達(dá)97.5%、99%，文中以 |Student-T|=2作為閾值。先將證據(jù)因素劃分為若干分級（分類），計(jì)算權(quán)重和標(biāo)準(zhǔn)差、Student-T，將 |Student-T|＜2的各分類視為顯著性低并歸為一類，保留 |Student-T|≥2的因素分類，然后重新計(jì)算歸并后各分類的權(quán)重值。

2.3 證據(jù)綜合

根據(jù)貝葉斯法則，任一單元K為地質(zhì)災(zāi)害的可能性，即對數(shù)后驗(yàn)概率可表示為[13-15,26,27]：

式中：Bi——第i個(gè)證據(jù)因素層；

K(i)——Wi是第i個(gè)證據(jù)因素存在或不存在的權(quán)重，在第i個(gè)證據(jù)因素層存在時(shí)是+，不存在時(shí)是-。

最后計(jì)算后驗(yàn)概率：

后驗(yàn)概率的大小作為易發(fā)性高低的指標(biāo)，值越大表示易發(fā)性越高，值越小表示易發(fā)性越低。

3 結(jié)果與討論

3.1 證據(jù)權(quán)重計(jì)算結(jié)果和相關(guān)性分析

證據(jù)權(quán)重計(jì)算結(jié)果（表2、圖4）與1.3 節(jié)可相互印證。在地形高程方面，1 800～1 850 m、1 920～1 950 m 和1 950～2 000 m 段利于地質(zhì)災(zāi)害發(fā)生，正權(quán)重0.555 0、1.175 8 和0.643 9。＞35°和15°～25°的山體斜坡較易于地質(zhì)災(zāi)害發(fā)生，正權(quán)重0.543 6 和0.378 5。坡向因素各分級權(quán)重值均不高，表明坡向?qū)Φ刭|(zhì)災(zāi)害發(fā)生的驅(qū)動(dòng)作用可能不太顯著。坡面曲率結(jié)果顯示，-0.75～-0.28（凹形）和-0.28～-0.15（凹形）兩個(gè)凹形坡分級段較易于地質(zhì)災(zāi)害發(fā)生，正權(quán)重0.569 0 和0.757 7。工程地質(zhì)巖組各巖組分類的正權(quán)重值總體不高，但砂巖、泥巖、頁巖巖組的統(tǒng)計(jì)結(jié)果仍然表現(xiàn)出對地質(zhì)災(zāi)害發(fā)生的較有利性，其正權(quán)重0.447 4，高于排在第二位的石英砂巖巖組（正權(quán)重值為0.294 7）。距斷層距離和距主要公路距離因素統(tǒng)計(jì)結(jié)果均顯示出了較明顯的距離效應(yīng)，即距斷裂或主要公路遠(yuǎn)的地區(qū)與地質(zhì)災(zāi)害發(fā)生負(fù)相關(guān)，距斷裂0～50 m 和距主要公路0～50 m、50～100 m 易于地質(zhì)災(zāi)害發(fā)生，其正權(quán)重0.797 3、0.982 0 和0.511 1。裸地或稀疏植被地區(qū)是易于地質(zhì)災(zāi)害發(fā)生的區(qū)域，其正權(quán)重0.871 9。

圖4 因素證據(jù)權(quán)重計(jì)算結(jié)果圖Fig.4 Calculation results charts of factor evidence weights

續(xù)表2

3.2 模型ROC 曲線和AUC

采用接受者操作特性曲線（Receiver Operating Characteristic Curve，ROC）和ROC 曲線下與坐標(biāo)軸圍成的面積（Area Under Curve，AUC）[29-32]評估模型擬合精度。模型擬合精度越好則AUC越接近1，0.7～0.9 時(shí)表示較好。文中建立的證據(jù)權(quán)法模型的AUC為80.4%，擬合精度優(yōu)異（圖5）。

圖5 模型預(yù)測性能ROC 曲線圖Fig.5 ROC curve of model prediction performance

3.3 易發(fā)性評價(jià)結(jié)果

綜合自然間斷點(diǎn)分級和地質(zhì)災(zāi)害分布，圈定了高易發(fā)區(qū)、中易發(fā)區(qū)和低易發(fā)區(qū)（表3、圖6），其中高易發(fā)區(qū)188.55 km2（占研究區(qū)總面積的49.41%），中易發(fā)區(qū)152.21 km2（占研究區(qū)總面積的39.88%），89.9%和9.1%的地災(zāi)點(diǎn)落入高易發(fā)區(qū)和中易發(fā)區(qū)，顯示易發(fā)性分區(qū)符合已發(fā)地質(zhì)災(zāi)害分布，模型預(yù)測性能較好。

表3 地質(zhì)災(zāi)害易發(fā)性分區(qū)表Table 3 Form of geological hazard susceptibility zoning

結(jié)合地質(zhì)環(huán)境因素特征分析西部高易發(fā)區(qū)（圖6藍(lán)色框范圍內(nèi)、圖7）主要位于砂巖、泥巖和頁巖巖組，斷裂構(gòu)造較密集，以山谷斜坡地貌為主，坡度15°～25°和＞35°較陡峭斜坡范圍成片發(fā)育且面積較廣，主要公路建于本區(qū)山谷，裸地/稀疏植被和草地連片覆蓋范圍較大。預(yù)測圈定的高易發(fā)區(qū)的這些分布特征，與上文分析得到的地質(zhì)災(zāi)害控制因素特征吻合，預(yù)測結(jié)果符合地質(zhì)災(zāi)害空間分布特征。

圖6 地質(zhì)災(zāi)害易發(fā)性柵格圖Fig.6 Grid map of geological hazard susceptibility

圖7 典型區(qū)因素和地質(zhì)災(zāi)害分布圖Fig.7 Factors and geological hazards in typical zone

4 結(jié)論

（1）“因子面積百分比A”“地災(zāi)數(shù)百分比B”和“比率β”，以及各因素各分類地質(zhì)災(zāi)害證據(jù)權(quán)重可以定量地分析各因素與地質(zhì)災(zāi)害發(fā)生的相關(guān)性。

（2）圈定高易發(fā)區(qū)188.55 km2（占總面積的49.41%），中易發(fā)區(qū)152.21 km2（占總面積的39.88%），易發(fā)性分區(qū)圖具有較好的等級區(qū)分度。

（3）通過證據(jù)權(quán)法繪制的地質(zhì)災(zāi)害易發(fā)性圖可以有效地預(yù)測該區(qū)地質(zhì)災(zāi)害，模型擬合精度AUC=80.4%。89.9%和9.1%的地災(zāi)點(diǎn)落入高和中易發(fā)區(qū)，建模結(jié)果與實(shí)際地質(zhì)災(zāi)害發(fā)育情況吻合度高，較好地揭示了研究區(qū)地質(zhì)災(zāi)害易發(fā)性特征。

（4）證據(jù)權(quán)法在研究區(qū)這類云南高原低山丘陵區(qū)有效性高，方法理論清晰，較為成熟，由數(shù)據(jù)驅(qū)動(dòng)，參數(shù)定義明確，易于一線工程師推廣使用。同時(shí)，該方法權(quán)重的估計(jì)和模型預(yù)測性能受預(yù)測因子選擇、因子數(shù)據(jù)空間分辨率、因子分級影響較大，具體工作中宜對這些問題進(jìn)行深入研究和統(tǒng)計(jì)分析。建議通過對因子分級進(jìn)行顯著性測試實(shí)施優(yōu)選，減小對權(quán)重的高估或低估，提高模型效能。