基于加權信息法的大孤山露天礦滑坡風險分析

李 旭 于慶磊 楊天鴻 董 鑫 鄧文學

（東北大學資源與土木工程學院，遼寧 沈陽 110819）

大孤山露天礦目前已進入深部開采階段，由于開采期間礦山巖體處于動態演化中，導致邊坡的穩定性受到諸多因素的影響。為此，有必要進行露天礦山的邊坡穩定性研究。基于算法預測邊坡發生失穩的高風險地區，對高風險進行重點加固防護，這樣有利于提高礦山生產的安全性和節省邊坡加固的成本。

近年來，在預測滑坡風險的案例中，許多學者在預測滑坡失穩方面做出重要貢獻。樊曉一等［1］基于三峽重點滑坡區域利用層次分析法和地形判別法，建立典型滑坡危險度判別矩陣進行危險度計算；白宇等［2］在概率論與數理統計方面與預測邊坡的極限平衡原理相結合，利用Monte-Carlo法來模擬隨機變量分布；李雪平等［3］基于歷史滑坡案例建立了Logistic回歸模型，利用ArcGIS進行影響因子提取，基于VisualC++6.0進行Logistic回歸模型分析；王萌等［4］以沐川縣的滑坡為對象，提取5個滑坡因子，采用貢獻率法計算權重；畢華興等［5］利用航片和TM數據提取并構建影響滑坡的8個因子，并建立滑坡因子專題圖，以8張滑坡因子滑坡專題圖為基礎進行抽樣調查，利用數量化理論公式進行數據處理，根據樣本得到的分值進行滑坡預測；于進慶等［6］基于GIS系統提取影響滑坡因子，并建立層次分析法的判斷矩陣，利用GIS統計獲取的滑坡因子專題圖面積得到的信息量作為層次分析法判別的重要程度進行滑坡因子計算；倪恒等［7］提出用敏感性分析思想來計算各因素對滑坡產生影響的大小，以趙樹嶺滑坡為實驗對象，考慮了5個影響因素，通過滑坡敏感性分析結合SARMA算法得出結論：水位變化＞地震荷載＞內摩擦角＞粘聚力＞巖體密度；張桂榮等［8］依托GIS平臺，采用規則網格單元劃分方法，運用信息量模型對該區斜坡穩定性進行了空間定量預測，并依信息量法的結果編制了該區的危險性預測分區圖；孫冉、李娜娜、吳麗清等［9-11］學者基于加權信息量算法，結合當地地質條件提取影響因素，并通過GIS軟件制作滑坡風險云圖。

綜上分析，眾學者已對滑坡風險預測做了諸多貢獻，但都是基于區域大、滑坡點多的情況。對于在礦山中使用這類預測方法并不多見。故在本項目中，將加權信息量法運用在礦山邊坡風險預測方面。通過對大孤山露天礦的數字高程模型（DEM）進行空間分析和地形數據提取，利用GIS的柵格圖繪制和矢量線段構建緩沖區功能提取礦山地質巖性和斷層因子。基于層次分析法給諸多因素進行重要性評價，并用GIS柵格統計進行信息量的評估，最終得出滑坡風險圖，從而達到提高礦山施工安全性，節約加固資金成本的目的。

1 大孤山鐵礦地質概況

大孤山鐵礦位于鞍山市東南12 km的千山腳下，礦區中心點坐標：東經 123°3'35.2443"、北緯41°3'5.22813929688"，占地10.6 km2，是亞洲最深的露天鐵礦。礦區主要巖性為太古界鞍山群櫻桃園組和元古界遼河群浪子山組的中深變質巖系，其次為震旦系釣魚臺石英巖及第四紀地層［12］。目前礦山（圖1）已開采至-330 m水平，開采深度420 m［13］。由于礦山進入深部開采階段，滑坡現象及巖石開裂變形情況愈加嚴重，引起大孤山鐵礦破壞的原因［14］主要有：結構面和人工開挖爆破。現階段，大孤山西北端已出現多處裂縫，且仍處在變形階段，嚴重威脅巷道運輸和礦山作業安全。

根據現場地質構造、水文地質和巖性等諸多條件，將大孤山礦坑已發生滑坡地區分為3個區域。其中I區位于西井邊坡F15斷層南部，巖性為太古代花崗巖（糜棱巖）；II區位于F15斷層以北、F14斷層以南，受2條斷層帶影響呈楔形發育；III區位于F14斷層以北，其巖性也為太古代花崗巖。邊坡破壞現狀如圖2所示。

2 滑坡評價因子的選取和層次分析法權重的確定

2.1 滑坡評價因子的選取

考慮當前礦山所能獲取到的數據情況，本文選取兩大影響層次，分別為地形地貌層次、地質層次。地形地貌層次有坡度、坡向、坡形、高程；地質層次有巖性、斷層。

由于本試驗采取加權信息量法，所以在選取因子時需利用層次分析法確定權重，并賦予到最后的信息量計算中。

2.2 層次分析法權重的確定

層次分析法是19世紀由Seaty提出的算法，該算法將多目標決策問題作為一個系統，將目標分為多個目的和準則，通過模糊的定性定量方法將其權重進行計算［15］。本算法提供3個層次，分別是目標層、決策層和方案層，目標層在本次滑坡災害風險評級里是進行滑坡風險的一個預測目標，而方案層是所選因子，決策層將所有因子分為地形地貌和地質構造2個層次。層次分析法的第二步是將所有不定性不定量數據進行重要程度比較，也稱專家打分。Seaty將打分表分為幾個重要程度，如表1所示。

地形地貌層次，參照其他學者在滑坡案例分析方面選取權重大小并予以賦值［16］；地質層次巖性和斷層的重要程度依舊參考學者在做地質地貌層次打分表進行打分［17］。層次分析法在專家打分結束后，兩兩因子按表1重要程度進行比較后便形成層次分析法判斷矩陣，如表2～表4所示。

計算得出的權重由于采用專家評審這一主觀判斷，且隨著層次分析法的判斷矩陣階數的增大，計算得出的最大特征向量比矩陣階數會大出許多，此時矩陣的不一致性會愈發的加重。為此本文對采取的打分權重進行一致性分析［18］，從而達到數據的客觀度。層次分析法在進行計算矩陣的特征向量的同時，也算出各個層次特征向量所對應的最大特征值λmax，而利用最大特征向量與矩陣階數的差值（一致性指標CI）大小可以衡量矩陣的不一致性。

式中，CI為一致性指標；n為層次分析法判斷矩陣階數。

為防止矩陣因隨機性而影響一致性指標，故結合隨機一致性指標RI表進行參考，一致性指標隨矩陣階數的不確定增大而數值發生相應的增加（表5）。

結合一致性指標和隨機引發矩陣發生偏差的指標，二者之間進行比值，便得出我們最終的檢驗系數CR。若CR＜0.1則認為該權重通過一致性分析，是符合客觀評價標準的，該數據權重可被使用；若CR＞0.1則認為不通過一致分析，則需重新取值。

水安全是水資源、水環境、水生態、水工程和供水安全五個方面的綜合效應，這五者之間并不是相互獨立的，如水資源的多寡，與水生態、供水之間有很強的關聯性，水工程、水環境與供水安全有很強的關聯性。

本文2個方案層，一個決策層進行權重賦值，其一致性比率分別為：地形地貌因子CR1為0、地質構造因子CR2為0，總決策層次CR3為0，故權重通過一致性檢驗，具有可信度。

將得到的各因子權重按照決策層所占重要程度進行歸一化處理如表6所示。

3 加權信息量法預測滑坡

3.1 信息量算法原理及滑坡像元統計

信息量法是基于信息論方法而建立起來的一種統計學方法，它不僅可以反映各種成災要素的相對敏感程度，也可反映特定成災要素中不同要素區間的致災貢獻大小（以成災要素為評價因子）［19］。其算法原理：

式中，W為研究區域內單位面積；W0為已經發生滑坡地質災害的單元面積之和；Sn為具有相同因子x1，x2，x3，…，xn組合的單元總面積；S0為具有相同因子x1，x2，x3，…，xn組合的單元中發生滑坡地質災害的總面積。

本文基于大孤山露天礦坑的數字高程模型（DEM）進行像元劃分，選取DEM柵格精度為0.1 m。為了在地圖中顯示完整，DEM數據最大高度與實際礦坑高程不符，該DEM數據是基于實際高程的500 m以上進行建立的。由于所有因子均由大孤山露天礦坑DEM數據提取，故本文以DEM數據為柵格像元劃分依據，對DEM柵格進行像元統計。在SuperMap iDesktop中具有柵格統計功能，統計得出大孤山露天礦坑總柵格單元為1 676 499個，根據現場提供的滑坡位置，進行滑坡面的繪制（圖3）。

選取圖3西北幫3個區域為滑坡位置進行統計得出，經統計得出滑坡區域柵格值為12 292個。

3.2 加權信息量因子數據的提取

基于SuperMap iDesktop自帶的GIS空間分析功能，將大孤山露天礦坑的DEM數據進行坡度、坡向、坡形、高程的提取，并對柵格進行單元值的分級。

3.2.1 地形地貌因子

坡向：坡向不同，巖體受陽光光照強度不同，土質的松散程度也受到相應的變化，坡向同樣影響著邊坡失穩。坡向提取與坡度提取方法同理，均采用GIS自帶的空間分析功能。以22.5°為間距，進行坡向等級劃分，共劃分16個等級。經統計得出90°～135°的坡向發生滑坡面積最多，占滑坡總面積34%。

坡形：通過對DEM表面曲率來獲取坡形，選取曲率大于0的值為直線型或凸形坡，小于0為凹形坡；經統計得出凸形坡與凹形坡所占滑坡面積數量相同，二者均勻分布。

高程：選取每100 m為1個等級進行劃分；共劃分6個等級，經統計得出-300～-200 m（實際高程）發生滑坡災害面積最多，占比為57.13%。

3.2.2 地質因子

斷層：距離斷層的大小同樣影響著邊坡的失穩程度，由于大孤山露天礦存在斷層數量較多，故將斷層因子考慮在內。以斷層矢量線段為中心，在該中心擴展20 m和50 m的圓形區域（緩沖區）。斷層分別構建了20 m和50 m的緩沖區，以緩沖區域進行統計得出距斷層20 m的緩沖區，發生滑坡災害的面積占比為16.31%；距斷層50 m發生滑坡災害面積占比為26.44%。

巖性：通過大孤山露天礦地質報告來進行巖性劃分，以不同巖性賦予不同柵格值進行云圖的繪制。本文巖性劃分為千枚巖、片麻狀糜棱巖、花崗斑巖、綠泥石英片巖、礦體、閃長玢巖。其中滑坡災害區域西北幫片麻狀糜棱巖占比最高，為70.33%。

由此便得到了坡度、坡向、坡形、高程、斷層、巖性各自對應的柵格圖。經信息量統計結果得出，信息量最大值（絕對值）為1.565 5，所對應因子為距斷層50 m；最小值（絕對值）為凹形坡，其信息量值為0.033 4。

3.3 加權信息量法的集成

利用的信息量計算公式（3）進行信息量統計，再結合上述層次分析法所確定的各因子的權重，待獲得上述數據后，根據式（4）可得知最終的加權信息量［20］，并統計如表7所示。

注：A1為地形地貌因子，A2為地質因子，信息量為負值時表示該因子與滑坡發生與否成負相關。

式中，Wi為該因子下層次分析法所算得的權重；I(y,xi)為該因子下的信息量的值。

4 大孤山露天礦滑坡風險云圖

根據計算得出的加權信息量的值，按照對應柵格所相應的加權信息量賦值，賦值采用GIS軟件的柵格重分級對相應區間進行相應賦值。由此就得到了相應因子的加權信息量的云圖。圖4對坡度進行信息量加權的柵格進行分級并予以賦值。

在進行重分級后的柵格，便是所求得加權信息量圖。同理可得坡向、坡形以及其他因素的加權信息量圖，如圖5所示。

將得到的所有信息量云圖進行加權計算，由于已經賦予了相應的層次分析法的權重，所以只需按照簡單的線性相加即可。這樣便得到最終的滑坡風險圖。經柵格統計得知信息量的取值范圍為-0.432 96～0.579 44之間，采用自然間斷法分為5級：低風險級、次低險級、中風險級、次高、高風險，得出大孤山露天礦山最終滑坡風險圖，如圖6所示。

5 結論

（1）本次加權信息量法采取2個決策層：地質層、地形層進行分析，通過GIS空間分析和其他工具手段提取影響災害因子，并將其運用在礦山，對礦山災害風險評價和監測預警都有著指導意義。

（2）經信息量統計，影響礦坑發生滑坡最高的因子是坡向在90°～115°因子和距斷層50 m以內的因子。

（3）通過對大孤山露天礦坑滑坡風險概率云圖分析高風險區主要在大孤山露天礦西北幫分布，這一點與現場實際考察較為相符。礦山高風險區占總礦區面積比例為7.2%。