孫 琳，任娜娜，李云安，胡樂健

(中國地質(zhì)大學(xué)(武漢)工程學(xué)院，湖北 武漢 430074)

0 引言

巖溶即喀斯特(KARST)。地下水和地表水對可溶性巖石的破壞和改造作用，這種作用及其所產(chǎn)生的地貌現(xiàn)象和水文地質(zhì)現(xiàn)象總稱為巖溶[1]。巖溶本身是一種常見的地質(zhì)現(xiàn)象，但是隨著我國經(jīng)濟(jì)的發(fā)展，基礎(chǔ)工程的建設(shè)，巖溶塌陷已經(jīng)成為我國主要地質(zhì)災(zāi)害之一。據(jù)記載，人類活動造成的巖溶塌陷存在著逐年遞增的趨勢。在巖溶地區(qū)，巖溶塌陷不明顯，但在人類工程的影響下，存在巖溶塌陷的可能性；因此進(jìn)行巖溶發(fā)育強弱預(yù)測、并對重點區(qū)域進(jìn)行巖溶塌陷的排查具有重要的現(xiàn)實意義。

地質(zhì)災(zāi)害預(yù)測是地質(zhì)災(zāi)害評估中重要的一環(huán)，一直廣受關(guān)注。國內(nèi)學(xué)者對巖溶塌陷不確定性預(yù)測也開始了一些研究。雷明堂等[2]率先將GIS技術(shù)應(yīng)用在巖溶塌陷不確定性評價中。隨后一些知識驅(qū)動模型應(yīng)用在巖溶塌陷預(yù)測中，其中權(quán)重的定量化取值最終依靠的是專家的理論和經(jīng)驗。主要方法有層次分析法(AHP)[3-5]、模糊層次分析法[6]、模糊綜合評判法[7-8]、灰色聚類法[9]等。但是基于知識驅(qū)動的模型受人為主觀性影響較大，不同情況下必將損失掉空間數(shù)據(jù)本身深層的信息。不同學(xué)者又將數(shù)據(jù)驅(qū)動模型運用到巖溶塌陷預(yù)測當(dāng)中，其中通過計算來確定各影響因素的權(quán)重。數(shù)據(jù)驅(qū)動模型主要包括人工神經(jīng)網(wǎng)絡(luò)(ANN)[10]、邏輯回歸[11]、證據(jù)權(quán)法[12]等。這些模型預(yù)測方法在滑坡、泥石流災(zāi)害區(qū)劃預(yù)測中已經(jīng)進(jìn)行大量的應(yīng)用。但是這些模型在巖溶塌陷區(qū)劃中只有少量的應(yīng)用，基本呈現(xiàn)出單點開花，個別學(xué)者應(yīng)用的現(xiàn)象，尤其是數(shù)據(jù)驅(qū)動模型迄今為止只有個別學(xué)者研究使用過，不具備說服力。因此，運用證據(jù)權(quán)模型對巖溶路基發(fā)育強弱進(jìn)行預(yù)測，依托GIS技術(shù)并且有著深厚的統(tǒng)計學(xué)理論基礎(chǔ)，這將為巖溶研究者提供一個更科學(xué)、直觀的研究背景。

本文選取丹霞互通樞紐巖溶區(qū)域，通過計算巖溶發(fā)育的先驗概率，客觀解釋了各影響因素與巖溶的關(guān)系，選取條件獨立的證據(jù)，利用證據(jù)權(quán)模型結(jié)合GIS技術(shù)確定巖溶發(fā)育強弱的后驗概率，從而得出在丹霞互通樞紐建設(shè)影響下巖溶塌陷易發(fā)分區(qū)。最后進(jìn)行預(yù)測概率的ROC精度評價分析。

1 方法介紹

證據(jù)權(quán)法(the weight of evidence，wofe)是基于貝葉斯準(zhǔn)則的一種離散的多元統(tǒng)計方法?；谪惾~斯準(zhǔn)則計算出各證據(jù)層因子的權(quán)重，有效避免了人的主觀性。證據(jù)權(quán)法最開始應(yīng)用于醫(yī)學(xué)領(lǐng)域，在20世紀(jì)80年代BONHAM-CARTER[13]和AGTERBERG等人[14]將證據(jù)權(quán)法應(yīng)用于地質(zhì)礦產(chǎn)領(lǐng)域，現(xiàn)已基本成熟。隨后又有學(xué)者將其應(yīng)用于滑坡、泥石流等地質(zhì)災(zāi)害領(lǐng)域。證據(jù)權(quán)首先進(jìn)行先驗概率的計算，就是用已發(fā)現(xiàn)的巖溶點進(jìn)行計算；然后進(jìn)行各個證據(jù)因子權(quán)重的計算；之后還要進(jìn)行證據(jù)因子的優(yōu)選；最后運用GIS將各權(quán)重因子疊加得出巖溶發(fā)育的后驗概率圖并進(jìn)行證據(jù)因子獨立性檢驗。

證據(jù)權(quán)模型原理如下。假設(shè)研究區(qū)域柵格數(shù)為N[T]，且每個柵格最多只有一個巖溶點，巖溶點的單元數(shù)為N[K]，所以任意單元發(fā)生巖溶的先驗概率公式為：

(1)

(2)

當(dāng)Wj+>0，而Wj-<0，表示該證據(jù)有利于巖溶形成；反之，不利于巖溶形成，因此這時該證據(jù)與巖溶形成無關(guān)。且Wj與證據(jù)層關(guān)系如下所示:

(4)

證據(jù)因子優(yōu)選，在計算每個證據(jù)層的權(quán)重后，需要將證據(jù)層二值化處理。但是一個證據(jù)層可能存在多級、多個權(quán)重，就需要確定最優(yōu)切值(cut of fvaule)。最優(yōu)切值可以通過對比值C(contrast)來確定。C=Wj+-Wj-表示一組空間證據(jù)與一組預(yù)測目標(biāo)之間的空間關(guān)聯(lián)度量。C值大于0為正相關(guān)，且C值越大，關(guān)聯(lián)越大；反之亦然。但是，在一些特定情況下，訓(xùn)練點較少，沒有明確的最大C值，可以采用學(xué)生化檢驗[15]。學(xué)生化檢驗(the studentizedC)stud(C)是一種顯著性檢驗。

(5)

式中：s2(Wj+)和s2(Wj-)表示為Wj+和Wj-的方差。計算方式如下：

(6)

在顯著性檢驗之后，取顯著性水平α=0.05，stud(C)>1.96，認(rèn)為95%與目標(biāo)預(yù)測相關(guān)，表明優(yōu)于其他切值。當(dāng)存在多個stud(C)大于1.96，可選取N{Bj∩K}最大者[12]。

最后可計算后驗概率。為了計算簡便，引入定義的幾率O{K}，其計算如下：

(8)

經(jīng)過推導(dǎo)，可以得出：

(10)

2 工程實例

本文利用ArcGIS Desktop 10.2軟件等功能和Excel 2013，成功地對位于廣東省韶關(guān)市仁化縣丹霞樞紐互通區(qū)域進(jìn)行了公路路基巖溶發(fā)育的后驗概率計算，并給出巖溶發(fā)育分區(qū)圖和預(yù)測精度評價。

2.1 研究區(qū)概況

擬建丹霞樞紐互通位于廣東省韶關(guān)市仁化縣，路線經(jīng)打鐵沖村與新建村。互通區(qū)設(shè)置有A、B、C、D、E、F、G、H共8個匝道，總計長度為6 258.7 m(本文只選取其中交匯的最主要一部分)?；ネ▍^(qū)以橋梁及路基形式通過，線路上部為沖洪積平原，地形平緩，地勢相對平坦，以水田、魚塘為主，有一條水溝穿過；丘陵低山，山體高大、山坡較陡，植被發(fā)育，場區(qū)有粵贛高速、水泥路及在建的韶贛鐵路與外界相通，交通條件較好。場區(qū)地貌上屬于低山丘陵單元區(qū)(圖1)。

圖1 研究區(qū)范圍及其概況Fig.1 The overview of the research area

2.2 證據(jù)因子分析

2.2.1證據(jù)因子數(shù)據(jù)準(zhǔn)備

2.2.2地層巖性

基于研究區(qū)的地層地質(zhì)圖和鉆孔數(shù)據(jù)，且?guī)r溶發(fā)育在可溶性巖石中，所以地層圖中忽略第四系粉質(zhì)黏土。地層分別是上白堊統(tǒng)南雄組(K2n)含礫砂巖和泥盆系上統(tǒng)錫礦山組(D3x)角礫狀灰?guī)r。表1中是通過計算確定證據(jù)因子和巖溶的關(guān)系，其中錫礦山組角礫狀灰?guī)r地層的stud(C)>1.96，可將該地層作為本證據(jù)圖層的二值化處理切值。

表1 地層巖性權(quán)重及對比度

2.2.3地下水

巖溶發(fā)育的一大基本條件就是可溶性水。本文根據(jù)原始數(shù)據(jù)計算出地下水位距離地面的距離，并做出地下水距地面距離分區(qū)圖。根據(jù)實際情況，地下水距離地面最小為1.25 m，最大為18.75 m。地下水距地面距離均勻劃分，可將地下水距地面的距離分為5個區(qū)域，分別是1.25～4.75 m、4.75～8.25 m、8.25～11.75 m、11.75～15.25 m、15.25～18.75 m。由表2可知，區(qū)域1.25～4.75 m可作為本證據(jù)圖層的最優(yōu)切值。且該條件符合理論經(jīng)驗認(rèn)知，地下水水位高的地方巖溶發(fā)育強烈。

表2 地下水距地面距離權(quán)重及對比度

2.2.4頂板厚度

頂板厚度是巖溶發(fā)育頂部距巖層頂面的距離。頂板厚度能夠表明巖溶在可溶性巖層中發(fā)育的位置。頂板厚度對巖溶發(fā)育有一定的影響，其中頂板越厚，水循環(huán)越差，巖溶可能就不發(fā)育。根據(jù)實際情況，將頂板厚度分為6個區(qū)域，分別是：0～3 m、3～6 m、6～9 m、9～12 m、12～18 m、18～35 m。由表3可知，證據(jù)優(yōu)選后，區(qū)域0～3 m可作為本證據(jù)圖層的最優(yōu)切值。并且符合理論經(jīng)驗認(rèn)知，在巖層分界面，巖溶發(fā)育強烈。

表3 頂板厚度權(quán)重及對比度

2.2.5斷層

斷層可以影響地下水流向，對巖溶發(fā)育有一定的影響。在研究區(qū)有一條斷層破碎帶WF1，物探成果發(fā)現(xiàn)成帶狀分布，向下延伸深度大，近北東走向，視寬度約25 m，影響寬度約60 m，構(gòu)造垂向延伸深度較大，構(gòu)造水平展布長度較大。根據(jù)距斷層距離分區(qū)，分別是：0～100 m、100～200 m、200～300 m、300～400 m、<200 m，<300 m。由表4可知，證據(jù)優(yōu)選后，區(qū)域0～100 m可作為本證據(jù)圖層的最優(yōu)切值。

2.2.6證據(jù)因子優(yōu)選

每個證據(jù)因子幾乎都有多級分區(qū)。這就需要將每個證據(jù)圖層進(jìn)行二值化處理。二值化處理之前按照方法確定最優(yōu)切值。經(jīng)計算(表1～表4)，選取的最優(yōu)切值為：地層巖性中灰?guī)r(D3x)、地下水距地面距離1.25～4.75 m、頂板厚度0～3 m、距斷層距離<100 m。

表4 斷層緩沖區(qū)權(quán)重及對比度

2.3 條件獨立性檢驗和后驗概率

2.3.1條件獨立性檢驗

(11)

式中：m——各證據(jù)圖層存在狀態(tài)的組合個數(shù)；

Pi——第i種組合狀態(tài)下的后驗概率；

N[i]——第i種組合狀態(tài)下的柵格數(shù)。

2.3.2危險區(qū)劃分

通過ArcGIS軟件生成巖溶發(fā)育的后驗概率圖，巖溶發(fā)育的地方容易產(chǎn)生巖溶塌陷，尤其是在公路路基施工設(shè)備等增加地面荷載的情況下，所以在施工過程中應(yīng)該提高警惕。在本文中巖溶發(fā)育后驗概率圖也可叫做巖溶塌陷風(fēng)險后驗概率圖。后驗概率圖在分級過程中需要確定分級閾值，通常采用CAPP(Cumulative Area—Posterior Probability)曲線的拐點確定閾值[17-18]。如圖2所示，考慮到先驗概率為0.015，則優(yōu)先選擇離先驗概率最近的拐點，閾值分別為0.007 6、0.014 2、0.022 0。則將后驗概率分四級，分別對應(yīng)穩(wěn)定、較穩(wěn)定、較危險、危險(圖3)。

圖2 后驗概率面積頻率曲線Fig.2 The cumulative area posterior probability curve

圖3 巖溶塌陷風(fēng)險后驗概率圖Fig.3 The posteriori probability map showing the risk of karst collapse

2.4 ROC精度評價

本文采用ROC(receiver operating characteristic curve)曲線[19-20]進(jìn)行巖溶預(yù)測準(zhǔn)確性評價。ROC預(yù)測曲線橫坐標(biāo)假陽性為預(yù)測巖溶區(qū)域(不包括真實巖溶區(qū)域)百分比，縱坐標(biāo)真陽性為實際巖溶區(qū)域百分比。ROC曲線下面積(AUC，area under curve)為預(yù)測概率(圖4)，預(yù)測概率為66.2%。

圖4 ROC預(yù)測曲線Fig.4 The ROC curve

在公路路基巖溶發(fā)育區(qū)，引入常用在礦產(chǎn)資源預(yù)測中的證據(jù)權(quán)模型，定量的分析計算每個證據(jù)因子的權(quán)重，并據(jù)此得出線性公路路基巖溶發(fā)育強度分區(qū)圖，由于人類工程活動易引發(fā)巖溶塌陷，公路路基巖溶發(fā)育強度分區(qū)圖也可叫做人類工程作用下公路路基巖溶塌陷發(fā)生危險分區(qū)圖。危險分區(qū)圖劃分為四個區(qū)，分別為穩(wěn)定(0.001 8～0.007 6)、較穩(wěn)定(0.007 6～0.014 2)、較危險(0.014 2～0.022 0)、危險(0.021 9～0.039 3)。其中穩(wěn)定區(qū)占27.1%、較穩(wěn)定區(qū)占35.4%、較危險區(qū)占32.9%、危險區(qū)占4.6%，可根據(jù)后驗概率圖，注意危險區(qū)施工。根據(jù)ROC的精度評價，預(yù)測成功率為66.2%(圖4)。巖溶最發(fā)育的地方基本上都在危險區(qū)，并結(jié)合預(yù)測成功率可知該分區(qū)圖基本符合實際結(jié)果。

3 結(jié)論

(1)需要指出的是，巖溶發(fā)育強度預(yù)測和巖溶塌陷預(yù)測略有不同，一個是以研究區(qū)內(nèi)已知巖溶發(fā)育點為訓(xùn)練點，一個是以研究區(qū)內(nèi)已知巖溶塌陷為訓(xùn)練點。在本研究區(qū)內(nèi)巖溶塌陷點基本不存在，所以預(yù)測巖溶發(fā)育強烈區(qū)在人類工程作用下為巖溶塌陷危險區(qū)。并且本文將巖溶發(fā)育強烈區(qū)作為公路路基施工過程中容易塌陷的區(qū)域，需注意做好防范措施。如果巖溶塌陷的訓(xùn)練點足夠多，可以直接預(yù)測巖溶塌陷危險區(qū)。影響巖溶發(fā)育的因素有很多，地表水流域，地下水變化等，但是由于數(shù)據(jù)缺少，所以無法進(jìn)行計算。并且采用的是巖溶鉆孔數(shù)據(jù)，鉆孔數(shù)據(jù)具有很強的目的性，隨機性不夠強。如有條件，可進(jìn)行巖溶塌陷危險區(qū)預(yù)測。

(2)在進(jìn)行證據(jù)因子優(yōu)選時，每個證據(jù)圖層幾乎都有多級分區(qū)，證據(jù)因子優(yōu)選后，每個證據(jù)圖層只有兩種狀態(tài)，這種情況就忽略了證據(jù)層中其他分區(qū)的狀態(tài)，這就造成了最終結(jié)果的系統(tǒng)誤差，這是證據(jù)權(quán)模型的弊端。為了解決這種弊端，成秋明等[21]提出模糊證據(jù)權(quán)模型，這樣就可以最大限度利用證據(jù)圖層中的信息。在計算證據(jù)條件獨立性檢驗時，運用的omnibus test檢驗方法是一種近似的檢驗，大致符合條件獨立性。