磨英飛,潘宏堅,鄧蔭萬
(廣西壯族自治區(qū)地質(zhì)環(huán)境監(jiān)測站,廣西 南寧 530201)
廣西北流市是典型的桂東南花崗巖丘陵區(qū),滑坡、崩塌地質(zhì)災(zāi)害多發(fā),開展地質(zhì)災(zāi)害易發(fā)分區(qū),在地質(zhì)災(zāi)害預(yù)警和防治規(guī)劃、國土空間規(guī)劃中具有重要應(yīng)用。目前,地質(zhì)災(zāi)害易發(fā)性評價方法有綜合指數(shù)法、證據(jù)權(quán)法、信息量法、多元統(tǒng)計分析法、邏輯回歸法等,綜合指數(shù)法、信息量法、多元統(tǒng)計分析法等在選取評價因子和權(quán)重中存在較大主觀性和人為因素,不同人員選取評價因子和賦予的權(quán)重可能會有差別,同時評價模型可靠度直接取決于地質(zhì)災(zāi)害原始數(shù)據(jù)的精度和發(fā)育數(shù)量,存在評價因子組合狀態(tài)比較多、樣本需求大、實際統(tǒng)計數(shù)量受限等缺點,具有一定的局限性。證據(jù)權(quán)法是采用貝葉斯條件概率原理、綜合各種證據(jù)層來預(yù)測某種事件發(fā)生概率的一種定量方法,在地質(zhì)災(zāi)害易發(fā)性評價中,充分考慮了地質(zhì)災(zāi)害與評價因子之間的相關(guān)關(guān)系以及評價因子之間的相互關(guān)系,評價結(jié)果較其他方法具有客觀性和可靠性。國內(nèi)學(xué)者和工作人員應(yīng)用證據(jù)權(quán)法[1-6]開展地質(zhì)災(zāi)害易發(fā)性評價,均取得了較好的效果。如胡燕等[7]采用證據(jù)權(quán)法開展巴東縣城滑坡災(zāi)害易發(fā)性評價,將巴東縣城劃分為滑坡極高易發(fā)區(qū)、高易發(fā)區(qū)、中易發(fā)區(qū)、低易發(fā)區(qū)與極低易發(fā)區(qū)5類;羅可[8]應(yīng)用證據(jù)權(quán)法進(jìn)行大埔縣崩塌滑坡地質(zhì)災(zāi)害易發(fā)性評價分區(qū),98.37%的災(zāi)害點位于中、高易發(fā)區(qū)內(nèi)。本文采用證據(jù)權(quán)法對廣西北流市南部石窩鎮(zhèn)、六靖鎮(zhèn)、清灣鎮(zhèn)[9]等地區(qū)進(jìn)行地質(zhì)災(zāi)害易發(fā)性評價和分區(qū)。
研究區(qū)位于廣西東南的玉林市北流市南部,面積397.9 km2,多年平均降雨量約2 200 mm,降雨多集中在4—9月,占年均降雨量80%。河流發(fā)育有5條小河流。
研究區(qū)屬構(gòu)造剝蝕丘陵地貌,分為高丘、低丘2種類型。
(1)高丘區(qū)山頂高程+300.0~+616.8 m,自然坡度一般30°~40°,局部45°~50°,溝谷切割深200~300 m,多呈“V”形,分布面積約192.6 km2。
(2)低丘區(qū)山頂高程+100.0~+300.0 m,自然坡度一般20°~30°,局部35°~40°,溝谷切割深50~150 m,多呈“U”形,分布面積205.1 km2。
植被以速生桉樹為主,覆蓋率較高,村莊多分布在低丘區(qū)的丘陵坡腳。
研究區(qū)地層巖性有侵入巖和變質(zhì)巖(圖1)。
侵入巖侵入時期為寒武紀(jì)—志留紀(jì)、白堊紀(jì),巖性以黑云(鉀長)二長花崗巖為主,少量黑云花崗片麻巖、片麻狀黑云二長花崗巖。
變質(zhì)巖包括寒武系淺變質(zhì)巖及晚元古代平政組花崗變粒巖。其中,寒武系淺變質(zhì)巖巖性為千枚巖、片巖、板巖,由含礫長石石英砂巖、不等粒砂巖、粉砂巖、頁巖淺變質(zhì)形成;平政組花崗變粒巖巖性為黑云變粒巖、黑云二長變粒巖、黑云斜長變粒巖及少量黑云石英片巖、含硅線堇青黑云片巖等。
花崗巖、花崗變粒巖區(qū)丘陵山坡分布第四系坡殘積粉質(zhì)黏土、砂質(zhì)黏土,厚2~3 m;下部為全風(fēng)化層砂土,以礫砂、粗砂為主,局部含黏性土、中砂、粉砂,呈散體狀,遇水易軟化崩解,厚度10~30 m。淺變質(zhì)巖區(qū)丘陵山坡第四系坡殘積黏土、粉質(zhì)黏土,厚度一般在0.5~3.0 m,局部厚5.0~8.0 m。
研究區(qū)構(gòu)造以南南西—北北東向為主,主要發(fā)育有陸川米場—石窩斷裂、大坡外—六靖斷裂、上珍斷裂、中龍斷裂、侯山斷裂、六侖斷裂、黃田斷裂7條次級斷裂。
研究區(qū)破壞地質(zhì)環(huán)境條件的人類工程活動主要有城鎮(zhèn)工程建設(shè)、居民臨坡建房、公路建設(shè)形成高切坡、高填方等。集鎮(zhèn)及村屯建設(shè)挖、填方形成高5~15 m、坡度50°~80°的高陡人工邊坡,道路建設(shè)高挖低填形成高5~20 m、坡度45°~70°的高陡邊坡,多數(shù)人工邊坡無支護(hù)、排水等防護(hù)措施,強(qiáng)降雨作用下易引發(fā)崩塌、滑坡地質(zhì)災(zāi)害。經(jīng)調(diào)查發(fā)現(xiàn),研究區(qū)有553處高陡人工邊坡。
研究區(qū)發(fā)育地質(zhì)災(zāi)害點155處。其中,崩塌30處,滑坡17處,不穩(wěn)定斜坡108處。崩塌滑坡物質(zhì)以殘坡積土體和全風(fēng)化花崗巖為主,花崗巖、花崗變粒巖區(qū)共發(fā)育地質(zhì)災(zāi)害106處,占災(zāi)害點總數(shù)的68.4%,淺變質(zhì)巖區(qū)共發(fā)育地質(zhì)災(zāi)害49處,占災(zāi)害點總數(shù)的31.6%。地質(zhì)災(zāi)害以小型為主?;麦w厚1~5 m,崩塌體厚2~3 m。地質(zhì)災(zāi)害主要發(fā)育在居民屋后人工邊坡和道路兩側(cè)的人工邊坡。
證據(jù)權(quán)法是一種條件獨立假設(shè)的前提下,基于貝葉斯條件概率的定量預(yù)測方法。證據(jù)權(quán)法最初是作為醫(yī)療診斷支持的方法,20世紀(jì)80年代末,加拿大數(shù)學(xué)地質(zhì)學(xué)家Bonham-Carter G F等[10-15]將該方法引入到礦產(chǎn)資源定量預(yù)測與評價。該方法作為一種人工智能模型,已被國內(nèi)外學(xué)者廣泛用于多元信息綜合評價,近年來又被引入到地質(zhì)災(zāi)害的評價中。其基本原理為:假設(shè)在一特定的研究區(qū)內(nèi),已知有n種二值證據(jù)圖層,1代表因子對災(zāi)害發(fā)生的證據(jù)存在,0代表不存在,通過證據(jù)權(quán)模型給出該二值化的證據(jù)因子圖層的權(quán)重,最終疊加多元圖層,實現(xiàn)地質(zhì)災(zāi)害易發(fā)性評價。證據(jù)權(quán)法的分析流程如下。
圖1 研究區(qū)地質(zhì)簡圖Fig.1 Geological sketch map of the study area
2.1.1 權(quán)重計算
(1)
式中,W+為證據(jù)因子存在區(qū)的權(quán)重值;W-為證據(jù)因子不存在區(qū)的權(quán)重值,其大小表示證據(jù)因子與地質(zhì)災(zāi)害發(fā)生和不發(fā)生的關(guān)系密切程度。證據(jù)因子和災(zāi)點正相關(guān)表示為W+>0,W-<0,負(fù)相關(guān)為W+<0,W->0,不相關(guān)時權(quán)重為0。P為不同條件發(fā)生的條件概率。
已知地質(zhì)災(zāi)害點的先驗概率為Po=D/T。證據(jù)因子權(quán)重由落入特定證據(jù)因子圖層的災(zāi)點數(shù)和全部災(zāi)點數(shù)之比與證據(jù)因子圖層面積和調(diào)查區(qū)總面積之比的比值決定。為表示證據(jù)因子對于地質(zhì)災(zāi)害發(fā)生和不發(fā)生的區(qū)分能力,可計算相對系數(shù)C=W+-W-,用來度量證據(jù)因子和地質(zhì)災(zāi)害之間的相關(guān)性大小,C值越大越有利于地質(zhì)災(zāi)害的出現(xiàn),C值越小越不利于地質(zhì)災(zāi)害的發(fā)生。
2.1.2 證據(jù)綜合
在權(quán)重值計算及分析的基礎(chǔ)上,通過證據(jù)層的優(yōu)選,選擇權(quán)重較大、與地質(zhì)災(zāi)害關(guān)系密切的證據(jù)層,剔除權(quán)重較小、與地質(zhì)災(zāi)害關(guān)系不密切的證據(jù)層;進(jìn)一步進(jìn)行證據(jù)因子相對災(zāi)點的條件獨立性檢驗,剔除地質(zhì)災(zāi)害權(quán)重相對較小而與其他證據(jù)因子相關(guān)性大的證據(jù)層。對最終篩選出的n個關(guān)于地質(zhì)災(zāi)害點條件獨立的證據(jù)因子,根據(jù)貝葉斯法則,研究區(qū)任一單元K發(fā)生地質(zhì)災(zāi)害的可能性,即對數(shù)后驗概率見式(2)。
(2)
式中,O為D的概率,O(D)=D/(T-D);D為存在地質(zhì)災(zāi)害的單元網(wǎng)格數(shù);Bi為第i個證據(jù)層;K(i)為在第i個證據(jù)因子層存在時為+,不存在時為-;Wi為第i個證據(jù)因子存在或不存在的權(quán)重。
最后計算后驗概率:
(3)
后驗概率值的大小表示易發(fā)性的高低,其值在0~1。后驗概率值越大,表示易發(fā)性越高;后驗概率值越小,表示易發(fā)性越低。
不同地域地質(zhì)環(huán)境條件和地質(zhì)災(zāi)害發(fā)育特征不同,選取的地質(zhì)災(zāi)害易發(fā)評價因子[16-19]存在差異。研究區(qū)地貌類型單一,巖組較單一,植被茂盛且類型較單一,地質(zhì)災(zāi)害主要形成于第四系殘坡積土體和全風(fēng)化花崗巖中。根據(jù)地質(zhì)災(zāi)害發(fā)育特征和形成機(jī)理,選擇地形地貌、地層巖性、斜坡結(jié)構(gòu)類型、土體厚度、構(gòu)造、人類工程活動等證據(jù)因子進(jìn)行地質(zhì)災(zāi)害易發(fā)性評價,各證據(jù)因子細(xì)分不同的證據(jù)層。
2.2.1 地形地貌
地形地貌選取坡度、坡面曲率等證據(jù)因子,將坡度劃分為 ≤10°、(10°,15°]、(15°,20°]、(20°,25°]、(25°,35°]、(35°,45°]、(45°,55°]、>55° 共8個證據(jù)層;將坡面曲率劃分為≤-7、(-7,-3]、(-3,3]、(3,7]、>7共5個證據(jù)層。
2.2.2 地層巖性
研究區(qū)地質(zhì)災(zāi)害主要發(fā)育在丘陵斜坡的第四系殘坡積粉質(zhì)黏土和全風(fēng)化花崗巖中,溝谷及河岸分布的第四系全新統(tǒng)沖洪積層地質(zhì)災(zāi)害不發(fā)育,各地層分布面積和地質(zhì)災(zāi)害發(fā)育情況如圖2所示。易崩易滑地層為晚志留世扶新單元(S3F)、六楊單元(S3Ly)侵入黑云(鉀長)二長花崗巖、寒武紀(jì)古桑單元(∈G)、平政組(Pt3p)黑云變粒巖、片巖、寒武系黃洞口組第一段(∈h1)淺變質(zhì)長石石英砂巖與頁巖,共劃分為11個證據(jù)層。
圖2 研究區(qū)各地層巖性分布面積及發(fā)育災(zāi)害數(shù)量柱狀圖Fig.2 Histogram of lithologic distribution area and development disaster of each stratum in the study area
2.2.3 斜坡結(jié)構(gòu)類型
研究區(qū)的淺變質(zhì)巖由沉積巖變質(zhì)形成,繼承了沉積巖的巖層層理,按斜坡坡向與巖層層理傾向間的夾角,劃分為順向坡(夾角≤30°)、同向斜交坡(30°,60°]、橫向坡(60°,120°]、反向斜交坡(120°,150°]、逆向坡(150°,180°]共5個證據(jù)層,侵入巖及平政組花崗變粒巖按順向坡考慮。
2.2.4 土層厚度
土層厚度考慮斜坡殘坡積層和全風(fēng)化層的厚度,劃分為≤ 3 m、(3,6 m]、(6,9 m]、(9,12 m]、>12 m共5個證據(jù)層。
2.2.5 構(gòu)造
包括與斷裂的距離和巖體節(jié)理裂隙等結(jié)構(gòu)面對地質(zhì)災(zāi)害的相關(guān)性。
與斷裂的距離劃分為[0,200 m]、(200,500 m]、(500,1 000 m]、>1 000 m共4個證據(jù)層。結(jié)構(gòu)面考慮不同傾向的結(jié)構(gòu)面對地質(zhì)災(zāi)害發(fā)育的控制或影響,劃分為[0°,30°]、(30°,60°]、(60°,90°]、(90°,120°]、(120°,150°]、(150°,180°]、(180°,210°]、(210°,240°]、(240°,270°]、(270°,300°]、(300°,330°]、(330°,360°] 共12個證據(jù)層。研究區(qū)結(jié)構(gòu)面傾向在120°~130°、170°~190°、225°~245°、265°~280°等方位區(qū)段上有比較明顯的優(yōu)勢,與120°~140°、170°~210°、225°~270°等方向的地質(zhì)災(zāi)害發(fā)育數(shù)量較多基本吻合(圖3),在0°~100°、300°~360°方位區(qū)段上地質(zhì)災(zāi)害分布零散。
圖3 研究區(qū)結(jié)構(gòu)面傾向與地質(zhì)災(zāi)害分布統(tǒng)計Fig.3 Structural inclination and distribution statistics of geological hazards in the study area
2.2.6 人類工程活動
人類工程活動考慮丘陵坡腳房屋建設(shè)強(qiáng)度、切坡強(qiáng)度,由于研究區(qū)內(nèi)道路沿線未發(fā)育地質(zhì)災(zāi)害點,同時切坡強(qiáng)度也考慮了道路的人工邊坡分布,因此不考慮交通建設(shè)強(qiáng)度指標(biāo)。
丘陵坡腳房屋建設(shè)強(qiáng)度可以采用房屋面積密度表示,選用最新的1∶1萬測繪信息和遙感影像資料,統(tǒng)計房屋面積密度,利用自然間斷點法分為4個等級,劃分為[0,23 780]、(23 780,74 312]、(74 312,157 540]、>157 540共4個證據(jù)層。
切坡強(qiáng)度采用高陡斜坡和地質(zhì)災(zāi)害的分布點密度表示,研究區(qū)調(diào)查了553處高陡斜坡和155處地質(zhì)災(zāi)害點,大部分地質(zhì)災(zāi)害點都是建房切坡引發(fā)的,高陡斜坡在強(qiáng)降雨等因素作用下可能失穩(wěn)發(fā)生崩塌滑坡,因此高陡斜坡和地質(zhì)災(zāi)害的分布點密度劃分為[0,11.41]、(11.41,28.29]、(28.29,55.67]、>55.67共4個證據(jù)層。
坡度、坡面曲率、斜坡結(jié)構(gòu)類型、土層厚度、與斷裂的距離、結(jié)構(gòu)面、丘陵坡腳房屋建設(shè)強(qiáng)度、切坡強(qiáng)度等均為連續(xù)型數(shù)據(jù),地層巖性為分類型數(shù)據(jù)。對于連續(xù)型證據(jù)因子數(shù)據(jù),采用計算累計證據(jù)權(quán)重判斷拐點的方法進(jìn)行分級。先按照較小的級距將連續(xù)性證據(jù)因子分級,然后在ArcGIS平臺下初步計算W+、W-累計權(quán)重和相對系數(shù)C,初步判斷累計權(quán)重值曲線的拐點作為分級的斷點值,將證據(jù)因子劃分為若干級,重新計算各個因子等級的權(quán)重,再繪制累計權(quán)重曲線,判斷曲線拐點和分級的合理性,并進(jìn)行多次迭代,直至分類趨于合理;分類型數(shù)據(jù)可直接分級。最終各評價因子的證據(jù)層權(quán)重及后驗概率計算結(jié)果見表1—表9。
表1 坡度證據(jù)層權(quán)重及后驗概率計算結(jié)果Tab.1 Calculation results of slope evidence layer weight and posterior probability
表2 坡面曲率證據(jù)層權(quán)重及后驗概率計算結(jié)果Tab.2 Calculation results of slope curvature evidence layer weight and a posteriori probability
證據(jù)權(quán)法的前提是各證據(jù)因子相互之間應(yīng)滿足條件獨立性假設(shè)。為確保評價指標(biāo)的相互獨立性,需對評價指標(biāo)進(jìn)行相關(guān)性分析,此次采用卡方檢驗方法進(jìn)行檢驗。首先,將各因子圖層二值化,即將相對系數(shù)C<0和C=0對應(yīng)的柵格賦值為0,將C>0的柵格單元賦值為1,在GIS軟件內(nèi)將地質(zhì)災(zāi)害點與二值模式柵格圖層疊加分析;然后,用地質(zhì)災(zāi)害點圖層匯總二值的數(shù)量,計算期望值和卡方值(二聯(lián)表計算過程簡略),各證據(jù)因子相互的卡方值見表10。
表3 地層巖性證據(jù)層權(quán)重及后驗概率計算結(jié)果Tab.3 Stratum lithology evidence layer weight and posterior probability calculation results
表4 斜坡結(jié)構(gòu)證據(jù)層權(quán)重及后驗概率計算結(jié)果Tab.4 Calculation results of evidence layer weight and posterior probability of slope structure
表5 土層厚度證據(jù)層權(quán)重及后驗概率計算結(jié)果Tab.5 Calculation results of evidence layer weight and posterior probability of soil layer thickness
表6 與斷裂的距離證據(jù)層權(quán)重及后驗概率計算結(jié)果Tab.6 Calculation results of weight and posterior probability of evidence layer of distance from fault
根據(jù)統(tǒng)計學(xué)概率,當(dāng)P(K2≥K)=0.001時,查表得卡方臨界值為10.828,即當(dāng)卡方值小于10.828時(可信度99.9%),兩因子之間不存在相關(guān)性,此時2組證據(jù)因子層可以同時參與易發(fā)性指數(shù)計算。
表7 結(jié)構(gòu)面傾向證據(jù)層權(quán)重及后驗概率計算結(jié)果Tab.7 Calculation results of evidence layer weight and posterior probability of structural plane tendency
表8 房屋建設(shè)強(qiáng)度證據(jù)層權(quán)重及后驗概率計算結(jié)果Tab.8 Calculation results of evidence layer weight and posterior probability of housing construction intensity
表9 切坡強(qiáng)度證據(jù)層權(quán)重及后驗概率計算結(jié)果Tab.9 Calculation results of evidence layer weight and posterior probability of slope cutting strength
根據(jù)表10計算結(jié)果,地形坡度、坡面曲率、地層巖性、土層厚度、切坡強(qiáng)度5個證據(jù)因子相互之間的卡方值小于臨界值10.828,滿足各證據(jù)因子相互之間獨立性假設(shè),滿足基于貝葉斯條件概率的證據(jù)權(quán)計算方法。研究區(qū)地質(zhì)災(zāi)害主要發(fā)育于第四系殘坡積土體和全風(fēng)化花崗巖、花崗變粒巖土體中,斜坡結(jié)構(gòu)類型、結(jié)構(gòu)面對地質(zhì)災(zāi)害的控制和影響僅屬個別現(xiàn)象,斷裂對地質(zhì)災(zāi)害的影響微弱。因此,選擇地形坡度、坡面曲率、地層巖性、土層厚度、切坡強(qiáng)度5個證據(jù)因子組合進(jìn)行易發(fā)性評價,也基本符合研究區(qū)的地質(zhì)災(zāi)害發(fā)育現(xiàn)狀。
研究區(qū)易發(fā)性評價最終選取地形坡度、坡面曲率、地層巖性、土層厚度、切坡強(qiáng)度等5個證據(jù)因子,采用證據(jù)權(quán)法計算各證據(jù)因子的權(quán)重和后驗概率,利用ArcGIS提取各證據(jù)層所在斜坡單元的后驗概率值,按斜坡單元對各證據(jù)因子的后驗概率值進(jìn)行疊加計算,采用自然間斷法,將研究區(qū)劃分為地質(zhì)災(zāi)害極高易發(fā)區(qū)、高易發(fā)區(qū)、中易發(fā)區(qū)、低易發(fā)區(qū)4個等級(圖4),其中極高易發(fā)區(qū)面積62.16 km2、高易發(fā)區(qū)135.22 km2、中易發(fā)區(qū)171.23 km2、低易發(fā)區(qū)29.29 km2。
表10 各證據(jù)因子相互卡方值計算結(jié)果Tab.10 Calculation results of mutual chi-square values of various evidence factors
圖4 北流市南部易發(fā)性評價分區(qū)Fig.4 Zoning map of vulnerability assessmentin the south of Beiliu City
根據(jù)易發(fā)性分區(qū)結(jié)果,極高易發(fā)區(qū)分布有69處地質(zhì)災(zāi)害點,高易發(fā)區(qū)分布有45處,中易發(fā)區(qū)分布有23處,低易發(fā)區(qū)分布有18處,極高、高易發(fā)區(qū)分布的災(zāi)害點占73.5%。通過繪制ROC曲線(圖5),計算AUC值為 0.718,說明證據(jù)權(quán)法在北流市南部易發(fā)性評價的可信度較高。
采用證據(jù)權(quán)法開展北流市南部地質(zhì)災(zāi)害易發(fā)性評價,根據(jù)北流市南部地質(zhì)災(zāi)害發(fā)育和分布特征,選擇合理的證據(jù)因子,通過檢測各證據(jù)因子獨立性,選擇地形坡度、坡面曲率、地層巖性、土層厚度、切坡強(qiáng)度5個證據(jù)因子組合作為評價指標(biāo),劃分出地質(zhì)災(zāi)害極高、高、中、低易發(fā)區(qū)。經(jīng)驗證,73.5 %的地質(zhì)災(zāi)害點落在極高、高易發(fā)區(qū)中,劃分結(jié)果具有較高的可信度。
圖5 易發(fā)性分區(qū)與地質(zhì)災(zāi)害數(shù)量ROC曲線Fig.5 Prone zoning and ROC curve of geological disaster quantity