白 宇，孫郡慶，劉明明，李衛(wèi)鵬

(國能新疆托克遜能源有限責任公司，新疆維吾爾自治區(qū)吐魯番市，830011)

0 引言

邊坡安全穩(wěn)定問題長期困擾著露天煤礦的高效開采，是露天煤礦安全生產不可回避的關鍵問題[1-2]。近年來，隨著GPS、邊坡雷達、柔性測斜儀等新型自動化監(jiān)測設備的廣泛應用，露天煤礦邊坡變形監(jiān)測技術得到長足發(fā)展，為主動防控邊坡災害提供越來越多的現場數據支持[3-6]。相比而言，露天煤礦邊坡風險判識技術明顯滯后，大量觀測數據沒有得到充分挖掘與利用[7-9]。目前，邊坡位移觀測數據主要用于臨滑險情預警，難以在潛在風險區(qū)的早期識別方面發(fā)揮作用?；骂A測缺乏多源監(jiān)測數據融合分析手段，通常僅圍繞某個監(jiān)測點的觀測數據進行分析，邊坡狀態(tài)判識存在主觀性和盲目性。

從系統(tǒng)辨識的角度看，邊坡立體監(jiān)測網不同空間點位處觀測到的監(jiān)測數據集蘊含邊坡系統(tǒng)在外部應力作用下，由無序穩(wěn)態(tài)向有序非穩(wěn)態(tài)轉換過程中時空全域演化的信息[10-12]。為挖掘隱藏在多源監(jiān)測數據背后的邊坡系統(tǒng)演化規(guī)律，一些學者引入數據挖掘領域中的聚類分析概念，對邊坡觀測數據進行特征提取和知識挖掘。秦雨樵等[13]綜合考慮邊坡各點的位移以及點安全系數，提出基于K-means聚類算法的滑面搜索方法；雷從雨[14]將邊坡節(jié)點位移聚類分析與變點分析相融合，提出改進的特征點位移突變判據；陳波等[15]采用劃分聚類和投影聚類對邊坡位移監(jiān)測資料進行測點區(qū)域劃分和特征提取，以達到邊坡監(jiān)測數據簡約的目的。

作為按觀測時間先后順序排列的一組數據序列，邊坡位移觀測序列表現出十分明顯的非線性、非平穩(wěn)特征，不僅序列內部相鄰觀測值存在關聯(lián)性，而且序列之間存在不等長、偏移性、時間軸扭曲等結構差異。然而，現有研究工作普遍忽略了位移觀測序列結構的復雜性，或者僅對觀測序列中特定時刻的觀測值進行聚類[13-14]，或者雖然將整個觀測序列看作聚類對象，但采用歐氏距離進行同步硬性度量[15]，聚類結果的準確性和可靠性無法得到保證。為此，筆者提出一種將動態(tài)時間彎曲(DTW)距離和K-means(K均值)聚類相結合的無監(jiān)督數據挖掘方法，從多源監(jiān)測數據集中聚類“觀測域上相近、空間域上連續(xù)、變形趨勢明顯”的監(jiān)測點子集，將其所在區(qū)域圈定為潛在滑移區(qū)。實驗結果表明，上述聚類方法能夠在滑動面尚未完全形成階段，提前辨識出邊坡潛在滑移區(qū)域。

1 邊坡動力災害演化規(guī)律

露天煤礦邊坡屬于復雜非線性系統(tǒng)范疇，其動力演化遵循一般非線性系統(tǒng)演化所具有的普遍規(guī)律，即從無序到有序的演化歷程[10-11]。在邊坡災害的早期階段，各子系統(tǒng)基本保持力學平衡，整個系統(tǒng)處于最無序的平衡態(tài)；隨著邊坡高度、暴露面積、暴露時間的增加，在大氣降水、開采擾動等外部應力的作用下，局部區(qū)域巖土體單元力學平衡狀態(tài)被打破，應力重新分布，表現出短程相關和局部有序；在邊坡災害的中后期階段，隨著微損傷從小尺度到大尺度的串級發(fā)展以及主控滑動面的逐漸形成，各子系統(tǒng)趨向協(xié)同作用，表現為長程相關和全局有序，直至整個邊坡系統(tǒng)發(fā)生突變而形成滑坡。

當邊坡子系統(tǒng)發(fā)展到有序階段時，所處區(qū)域上的監(jiān)測點位移觀測序列必然表現出相似的演化行為，即有序區(qū)域具有“觀測域上相近、空間域上連續(xù)、變形趨勢明顯”的特點。因此，通過邊坡立體監(jiān)測網各測點的時空域數據進行聚類分析，可以對露天煤礦邊坡系統(tǒng)的有序程度進行定量化表征，為邊坡災害的早期辨識與晚期預警預報提供支撐。

2 位移觀測序列相似性度量方法

邊坡位移觀測序列是按觀測時間先后順序排列的數據序列，可以采用時間序列相似性度量方法進行相似性評價。時間序列相似性度量方法主要分為歐氏距離和DTW距離2種。歐式距離和DTW距離數據匹配如圖1所示。

圖1 歐式距離和DTW距離數據匹配

圖1中藍色曲線與紅色曲線代表2個時間序列，灰色曲線代表2個時間序列元素之間的相似匹配關系，曲線兩端代表2個元素。歐式距離是一對一元素相似匹配(藍色曲線上的一個元素匹配紅色曲線上的一個元素)，DTW距離是一對多元素相似匹配(藍色曲線上的一個元素有可能匹配紅色曲線上的多個元素)。

2.1 歐氏距離

歐氏距離(ED)是一種對時間序列進行同步硬性度量(“一對一”匹配)的方法。設A=(a1，a2，…，an)和B=(b1，b2，…，bn)分別代表2個具有相同長度的位移觀測序列，那么它們之間的歐氏距離定義為：

(1)

式中：d(A，B)——歐氏距離；

n——位移觀測序列的長度，即觀測時間的總次數；

ai——位移觀測序列A中第i個觀測值；

bi——位移觀測序列B中第i個觀測值。

d值越小說明序列A和B相似度越高，當d=0時，可認為A=B。

2.2 DTW距離

DTW距離是一種對時間序列進行異步相似形態(tài)度量(“一對多”匹配)的方法[16]。假定2個位移觀測序列A=(a1，a2，…，an)和B=(b1，b2，…，bn)，長度分別為n和m，則序列A和B的DTW距離計算過程如下。

2.2.1 距離矩陣

構造大小為n×m的距離矩陣：

(2)

式中：d(ai，bj)——序列A中第i個觀測值ai到序列B中第j個觀測值bj的距離。

2.2.2 彎曲路徑

定義彎曲路徑：

(3)

式中：wr——序列A中第i個點與序列B中第j個點之間的匹配關系。

并且，彎曲路徑需要滿足3個約束條件：

(1)邊界約束。彎曲路徑起于距離矩陣D的左下角，止于距離矩陣D的右上角，即w1=(1，1)，wk=(n，m)。

(2)單調性約束。對于相鄰元素wr=(i，j)和wr+1=(i′，j′)，需要滿足i≤i′和j≤j′，即彎矩路徑只能沿著時間軸單向移動。

(3)連續(xù)性約束。對于相鄰元素wr=(i，j)和wr+1=(i′，j′)，需要滿足i′≤i+1和j′≤j+1，即彎曲路徑上任意2個相鄰的元素在距離矩陣D中也相鄰。

2.2.3 DTW距離

定義序列A和B之間的DTW距離：

DTW(A，B)=min(∑(i，j)∈Wd(ai，bj))

(4)

由式(4)可知，DTW(A，B)對應的彎曲路徑為最優(yōu)路徑，在此路徑上的累積距離最小?；趧討B(tài)規(guī)劃，求解式(2)等價于構造累積距離矩陣Γ，其中的元素服從如下遞歸關系：

(5)

式中：γ(i，j)——累積距離。

當完成累積距離矩陣Γ計算后，元素γ(n，m)便為序列A和B的DTW距離，即：

DTW(A，B)=γ(n，m)

(6)

進一步，以Γ的右上角為起點，即wk=(n，m)，以累積距離最小為選取原則，逆序尋找W中各元素，直至Γ的左下角，即w1=(1，1)，可以得到完整的最優(yōu)彎曲路徑。

DTW距離不僅能實現不等長時間序列的度量，而且對時間序列的偏移、形狀縮放、振幅變化、時間點異常數據具有較強的魯棒性。因此，筆者采用DTW距離對位移觀測序列的相似性進行分析。

3 位移觀測序列相似性聚類方法

聚類分析一直是數據挖掘和機器學習研究的熱點之一，其目的是將數據集中的數據點按各自的相似度劃分為不同的簇類，并確保簇內的數據點相似程度最大、簇間數據點相似程度最小，從而判識數據集合的結構，提取需要的重要信息[17]。為深度挖掘位移觀測序列背后蘊含的邊坡系統(tǒng)演化規(guī)律，采用基于DTW距離的K-means聚類算法，對位移觀測序列進行相似性聚類分析。

3.1 算法原理

假定數據集X由n個位移觀測序列組成，Χ=(Χ1，Χ2，…，Χn)，每個序列含m個觀測值，Χi=(χi1，χi2，…，χim)。設定聚類數k，則聚類算法目標是找到一組簇類集合{S1，S2，…，Sk}，使得簇間觀測序列盡可能不同，而簇內觀測序列盡可能相似。其中，每個簇的聚類中心由cj=(cj1，cj2，…，cjm)表示，1≤j≤k。具體K-means算法框架如下。

(1)初始化階段。

步驟1：在數據集X中，隨機選取k個觀測序列作為初始聚類中心。

(2)聚簇階段。

步驟2：遍歷數據集X中所有的觀測序列，根據各觀測序列到每個聚類中心的DTW距離，將觀測序列分配到距離最近的聚類中心所代表的簇中，構建k個子集{S1，S2，…，Sk}。

步驟3：計算每個子集的位移觀測均值，作為新的聚類中心。

重復步驟2～3，直到k個聚類中心不再變化，或達到指定的迭代次數。在此過程中，每個觀測序列有且僅屬于一個簇，并且與該簇聚類中心之間的距離最小。

3.2 聚類中心計算

每個簇的聚類中心由該簇內觀測序列各個時間點觀測值的均值構成，相應計算公式如下：

(7)

式中：cjt——簇Sj聚類中心cj第t個元素；

χlt——簇Sj內觀測序列Χl第t個時間點的觀測值；

|Sj|——簇Sj中觀測序列的個數。

4 邊坡潛在滑移區(qū)智能辨識

依據露天煤礦邊坡系統(tǒng)所表現出的從無序到有序的動力演化特征，提出一種基于多源數據融合的邊坡潛在滑移區(qū)辨識方法，通過對監(jiān)測數據集的聚類分析，將“觀測域上相近、空間域上連續(xù)、變形趨勢明顯”的監(jiān)測點選定為有序測點，將有序測點所在區(qū)域圈定為潛在滑移區(qū)。具體辨識流程如下。

(8)

(2)對監(jiān)測數據集F進行數據清洗，消除數據集中的缺失值、重復值、異常值。

(3)設定聚類數k，采用基于DTW距離的K-means聚類算法，對F觀測域上的位移序列進行相似性度量及聚類分析，將F分解為k個簇，每個簇中測點位移觀測序列行為相近。計算同簇測點位移觀測序列的均值，獲得每個簇的聚類中心。

(4)對各簇中心進行排序，選出[t0，tn]區(qū)間位移變化最大的簇中心，并對該簇測點的空間連續(xù)性進行評估。如果存在連續(xù)性，則將該簇測點選定為有序測點，將該簇測點所在區(qū)域圈定為潛在滑移區(qū)?；禄瑒用嬷悄鼙孀R的流程如圖2所示。

5 實驗研究與結果分析

采用有限元數值建模方式仿真模擬典型邊坡漸進失穩(wěn)全過程，輸出邊坡不同演化階段的位移場作為實驗測試樣本，對提出的邊坡潛在滑移區(qū)辨識方法展開適用性研究。

5.1 數值建模與計算

采用文獻[18-19]中的邊坡算例：邊坡高度H=20 m，坡角β=26.57°，楊氏模量E=100 MPa，泊松比ν=0.3，容重?=20 kN/m3，黏聚力c=10 kPa，內摩擦角φ=20°。土體材料為理想彈塑性材料，強度準則采用平面應變條件下莫爾-庫侖結合DP準則(DP4準則)。邊坡尺寸及網格劃分如圖3所示，邊坡左、右兩側邊界為法向約束，底邊為雙向固定約束，共劃分四邊形單元5 686個、節(jié)點5 878個。采用有限元強度折減方法，按照c/ω，tanφ/ω方式不斷折減土體材料強度參數，模擬邊坡漸進破壞過程及其位移場演化信息，直至邊坡發(fā)生整體失穩(wěn)。其中，ω為折減系數，可看作廣義時間，其初始值ω0=1.0，步長△ω=0.002，第m步的折減系數ωm=ωm-1+△ω。

圖3 邊坡尺寸及網格劃分

折減系數分別為ω90=1.18、ω150=1.30、ω185=1.36和ω192=1.384時，邊坡等效塑性應變云圖和局部化帶擴展路徑如圖4所示。由圖4可以看出，ω90、ω150、ω185和ω192時刻分別對應于邊坡漸進演化的早期、中期、晚期和失穩(wěn)階段；當邊坡處于前3個階段時，雖然滑帶自坡腳向坡體內部逐漸擴展，但整個滑帶尚未貫通，滑移區(qū)也沒有完全形成，無法通過塑性云圖和滑帶擴展路徑提前辨識潛在滑移區(qū)；當邊坡進入失穩(wěn)階段、滑帶完全貫通后，通過塑性云圖和滑帶擴展路徑可以準確獲取滑移區(qū)范圍。

圖4 不同折減系數工況下邊坡內部等效塑性應變云圖及局部化帶擴展路徑

5.2 監(jiān)測樣本構建

將邊坡模型中的所有節(jié)點選做監(jiān)測點，以ωi0為觀測起始時間，以ωn為觀測截止時間，從數值模型計算結果中提取觀測數據，構建邊坡監(jiān)測數據集。

(9)

5.3 辨識算法測試

以觀測起始時間為分類標準，將監(jiān)測數據集細分為自邊坡漸進演化的早期、中期、晚期啟動觀測的監(jiān)測數據集，按照前文所述方法構建監(jiān)測數據集。設定聚類數k=3，開展監(jiān)測數據集觀測域聚類。聚類結果通過測點空間分布圖和簇中心曲線圖可視化展示。其中，同簇測點及其簇中心曲線使用相同顏色。

5.3.1 邊坡漸進演化早期階段啟動的觀測數據集(簡稱“早期監(jiān)測數據集”)

圖5 同簇測點空間分布與簇中心曲線(早期監(jiān)測數據集聚類結果)

5.3.2 邊坡漸進演化中期階段啟動的觀測數據集(簡稱“中期監(jiān)測數據集”)

圖6 同簇測點空間分布與簇中心曲線(中期監(jiān)測數據集聚類結果)

5.3.3 邊坡漸進演化晚期階段啟動的觀測數據集(簡稱“晚期監(jiān)測數據集”)

5.3.4 辨識結果分析

通過對早期、中期和晚期監(jiān)測數據集的聚類挖掘，可以得到如下結論：

(2)在監(jiān)測點數量足夠多的前提下，邊坡潛在滑移區(qū)辨識效果將主要取決于觀測的起始時間和觀測的累積時長。觀測起始時間越早，滑移區(qū)可辨識的時間也越早，但所需的觀測時長也越多。以ω90為觀測起始時間，于ω155時刻可辨識出潛在滑移區(qū)，需要65個時間步；以ω150為觀測起始時間，于ω175時刻可辨識出潛在滑移區(qū)，需要25個時間步；以ω185為觀測起始時間，滑移區(qū)辨識僅需要5個時間步，但辨識時間也最晚。

(3)現場位移監(jiān)測數據中的本底噪聲對辨識效果有很大影響，只有當觀測累積位移值明顯大于本底噪聲時，才有可能對潛在滑移區(qū)進行聚類辨識。因此，選用觀測精度高的監(jiān)測設備對于潛在滑移區(qū)的早期辨識十分重要。

6 結論

(1)邊坡立體監(jiān)測網不同空間點位處觀測到的監(jiān)測數據集蘊含邊坡系統(tǒng)在外部應力作用下，由無序穩(wěn)態(tài)向有序非穩(wěn)態(tài)轉換過程中時空全域演化信息。為深度挖掘隱藏在多源監(jiān)測數據背后的邊坡系統(tǒng)演化規(guī)律，提出一種將DTW距離和K-means聚類相結合的無監(jiān)督數據挖掘方法，從多源監(jiān)測數據集中聚類“觀測域上相近、空間域上連續(xù)、變形趨勢明顯”的監(jiān)測點子集，將其所在區(qū)域圈定為潛在滑移區(qū)。

(2)采用強度折減有限元數值模型仿真模擬典型邊坡漸進失穩(wěn)全過程，輸出邊坡不同演化階段的位移場作為實驗測試樣本，對潛在滑移區(qū)辨識方法展開適用性研究。實驗結果表明，該聚類方法能夠捕捉邊坡系統(tǒng)從無序到有序的漸進演化特征，可以在滑動面尚未完全形成階段，提前辨識出邊坡潛在滑移區(qū)域。

(3)需要指出的是，采用數值模型構建監(jiān)測樣本與邊坡現場觀測數據之間存在差異，并且將所有節(jié)點都作為監(jiān)測點是理想情況。因此，加快構建露天煤礦邊坡地表-深部立體監(jiān)測網，開展基于現場有限測點數據的潛在滑移區(qū)早期識別，將成為未來研究的重點。