安 北，蔣亞楠,2，曾啟菲

(1.成都理工大學地球科學學院，四川 成都 610059;2.成都理工大學地質災害防治與地質環(huán)境保護國家重點實驗室，四川 成都 610059)

滑坡地質災害在中國各個地區(qū)時有發(fā)生，具有高頻性、突發(fā)性的特點，每年給國家和人民帶來了不可估量的財產(chǎn)損失[1]，因此，積極開展滑坡預測預報研究工作對預防滑坡的產(chǎn)生，減少財產(chǎn)損失具有重要的實踐意義。

目前，在滑坡位移的預測預報工作中，大多學者應用較多的預測模型有灰色模型(GM)[2]、人工神經(jīng)網(wǎng)絡模型(ANN)[3]，也有學者對GPS監(jiān)測的滑坡位移時間序列采用不同的分解方法得到多個周期項、趨勢項、隨機項、噪聲項序列，并根據(jù)不同序列的特點采用相應的機器學習算法或非線性擬合方法進行訓練并預測，最后對應累加得到未來多期的滑坡位移預測值。劉藝梁等[4]利用經(jīng)驗模態(tài)分解(EMD)對非平穩(wěn)的滑坡位移信號進行分解并結合BP神經(jīng)網(wǎng)絡進行預測，鄧冬梅等[5]進一步提出了基于時間序列集合經(jīng)驗模態(tài)分解(EEMD)與機器學習算法結合的位移預測方法。機器學習預測模型中，最為典型的就是支持向量機(SVM)，同時粒子群算法(PSO)[5-6]、蟻群算法(ACO)[7]、差分進化算法(DE)[8]等也被引入SVM模型中以優(yōu)化參數(shù)。此外李瀟等[9]引入最小二乘法以轉化SVM中二次規(guī)劃問題并結合小波變換，證實了最小二乘支持向量機(LSSVM)預測滑坡變形的可靠性，Zhu等[10]在LSSVM基礎上進一步考慮了降雨對滑坡的影響，取得了較高的預測精度。

以上的科研工作中均取得了一定的成果，但也存在一些不足：GM無學習過程，對非線性滑坡預測精度不高，ANN易出現(xiàn)局部收斂并過于依賴樣本的情況。在累計滑坡位移時間序列的分解上，EMD分解方法得到的IMF分量個數(shù)未知，具有不可控性，若分解數(shù)量過少則考慮不夠周全，分解數(shù)量過多則會增加預測模型的工作量，同樣，受滑坡位移序列極其復雜的非線性特點，小波變換分解出的高頻子序列往往過多，加大了預測的工作量。SVM模型中，核函數(shù)參數(shù)δ及誤差懲罰因子C存在不確定性，需要進行人工設定。

本實驗中，采用變分模態(tài)分解(VMD)對累計滑坡位移進行分解，相比于EMD方法，VMD分解出的低頻分量更容易表達出信號波動趨勢，且可以確定合理的模態(tài)分解個數(shù)。其自適應性特點可根據(jù)信號確定序列的模態(tài)分解個數(shù)，匹配各模態(tài)最優(yōu)中心頻率和有限帶寬，從而實現(xiàn)各IMF分量的分離。對于低頻分量，采用二階指數(shù)平滑(SES)進行擬合，并逐期預測。而高頻分量，通過灰色關聯(lián)分析法對三峽地區(qū)滑坡的誘發(fā)因素進行分析，確定了模型的輸入影響因子包含降雨量、庫水位，與對應時間段的滑坡變形量作為BP神經(jīng)網(wǎng)絡的訓練集進行學習并預測，最后疊加得到預測累計位移值。

1 模型建立

1.1 模型概述

通過灰色關聯(lián)分析法分別將滑坡月累計降雨量，月平均庫水位高度與該月滑坡變形量進行關聯(lián)分析，證實了滑坡除受自身地質作用，如重力勢能，巖層土質條件等影響外，區(qū)域降雨對滑坡的蠕動變形也具有較強影響作用，此外，針對于三峽地區(qū)一類涉水型滑坡，庫水位高度的變化量也對滑坡變形存在一定的影響。通過對三峽地區(qū)某堆積層滑坡上安置的GPS監(jiān)測器所測定的累計滑坡變形量數(shù)據(jù)進行變分模態(tài)分解，分解為多個IMF分量，對于低頻分量，因為其逐漸穩(wěn)定遞增，可認為是滑坡變形的趨勢分量，采用二階指數(shù)平滑進行擬合，并對未來多期變化進行預測。對于高頻分量，其存在一定的周期性變化規(guī)律，通過灰色關聯(lián)分析對隨機樣本的實驗，得出滑坡變形受每年周期性的降雨與庫水位變化等外界因素的影響，通過相應的影響因子作為模型的輸入量，在BP神經(jīng)網(wǎng)絡模型中進行學習和訓練，最終得到短期的預測值。

1.2 變分模態(tài)分解原理

變分模態(tài)分解(VMD)最早由Dragomiretskiy[11]等人提出，是信號領域一種全新的自適應處理方法，該方法通過迭代搜尋變分模態(tài)的最優(yōu)解，最終將信號分解與重構，得到若干IMF分量以及降噪信號。VMD對采樣和噪聲具有較強的魯棒性，能夠減弱模態(tài)混沌現(xiàn)象。該方法目前正逐步應用于地震信號去噪、圖像處理等領域。

VMD將輸入信號f分解為離散的模態(tài)uk，這些子信號再輸入時具有特殊的稀疏性。因此每種模態(tài)可以選擇其頻譜域中的有限帶寬作為稀疏性。假設每個子信號在中心頻率ωk周圍大部分是緊湊的，并隨分解而確定，則為了評估信號的有限帶寬，進行以下調整。

a)對于每個模態(tài)，通過希爾伯特變換得到相關的分析信號，以獲得單邊頻譜。

b)對于每個模態(tài)解析得到的分析信號，通過與各自估計的中心頻率的指數(shù)混合，將模式的頻譜移至“基帶”。

c)接著通過解調信號的H1高斯平滑度估算有限帶寬，產(chǎn)生約束變分問題,見式(1):

(1)

d)式(1)中，K為分解模態(tài)個數(shù)，uk、ωk分別為分解后第k個模態(tài)分量和中心頻率，δ(t)為狄拉克函數(shù)，*為卷積運算,見式(2):

L({uk},{ωk} ,λ)=

(2)

e)式(2)中，α為二次懲罰因子，用于降低高斯噪聲的干擾，從而使uk與ωk得到優(yōu)化，將優(yōu)化的解引入ADMM算法中，并在適當?shù)那闆r下直接在Parseval / Plancherel傅里葉變換中進行優(yōu)化，來增強拉格朗日函數(shù)的鞍點，最終得到完整的變分模態(tài)算法，見式(3)—(5)：

(3)

(4)

(5)

f)式(5)中，γ為噪聲容忍度，最后通過迭代滿足式時，結束迭代，輸出各個IMF分量，見式(6)：

(6)

1.3 二次指數(shù)平滑

二次指數(shù)平滑是基于一次指數(shù)平滑上進行的修正，去解決時序信號出現(xiàn)直線趨勢時產(chǎn)生的滯后偏差等問題，利用該偏差規(guī)律建立趨勢預測模型，計算如下。

(7)

相較于傳統(tǒng)的多項式擬合方法，二次指數(shù)平滑計算更簡單，適應性也更強。

1.4 BP神經(jīng)網(wǎng)絡原理

BP神經(jīng)網(wǎng)絡是一種按照誤差逆向傳播訓練的多層前饋神經(jīng)網(wǎng)絡，主要學習過程是信號的正向傳播以及誤差的反向傳播[12-13]。其網(wǎng)絡基本結構包含輸入層、隱含層、輸出層，見圖1。

圖1 BP神經(jīng)網(wǎng)絡原理示意

圖中wmq表示輸入層單元到隱含層的權重，θq表示隱含層中第q個節(jié)點的閾值，若q足夠多則能夠以任意精度去逼近期望值，q=1…i…q，wqL表示從隱含層到輸出層權重。

正向傳播，將信號x1…xj…xm作為 BP神經(jīng)網(wǎng)絡的輸入值；從輸入層經(jīng)過各層權重wmq計算到達隱含層，在隱含層的激勵函數(shù)φ(x)中進行處理得到hq，見式(8)。

netφq=w1q×x1+w2qx2+…+wmq×xm，

hq=φ(netφq+θq)

(8)

信號再從隱含層經(jīng)過各層權重計算，將計算值代入激勵函數(shù)φ(x)得到輸出值yL：

netφL=w1L×h1+w2L×h2…+wqL×hq,

yL=φ(netφL+aL)

(9)

若輸出值與理想值偏差較大，則計算出輸出值與輸入信號對應的期望值之間的誤差值E：

(10)

為加快網(wǎng)絡學習速度并提高數(shù)據(jù)精度，本實驗的數(shù)據(jù)均進行了無量綱化。由此將月累計降雨量，滑坡變形量，月平均庫水位變化量歸一化處理為x1…xj…xm，將x1…xj…xm作為信號，通過訓練擬合效果選擇合適的激勵函數(shù)φ(x)與φ(x)，得到預測的輸出層yk。

2 工程實例

2.1 灰色關聯(lián)分析

以三峽地區(qū)白水河滑坡為例，該滑坡位于長江南岸，滑坡體處于長江寬河谷地段，南高北低，呈階梯狀向長江展布。前緣抵長江，東西兩側以基巖山脊為界，總體坡度約30°，其后緣高程達410 m，以巖土分界處為界。其南北向長度約600 m，東西向寬度約710 m，滑體平均厚度約30 m，體積約1.26×107m3，為堆積層滑坡[14]。取該滑坡區(qū)域GPS監(jiān)測器72期歷史數(shù)據(jù)進行實驗，每月1期，通過灰色關聯(lián)分析法研究影響滑坡形變的影響因素。為保證分析樣本的隨機性和可靠性，隨機選取12期滑坡的變形量與該期對應月累計降雨量，月平均庫水位變化量，使用灰色關聯(lián)分析法，比較各因素之間的發(fā)展趨勢的相似程度。具體步驟如下。

a)確定影響系統(tǒng)行為特征的比較數(shù)列，其中影響目標序列的因素ai包括月平均庫水位變化量，月累計降雨量，累計形變量，記作式(11)：

在我院檔案管理工作中應用優(yōu)秀計算機信息技術,實現(xiàn)信息全面交互,創(chuàng)建辦公自動化大環(huán)境,提高自身工作效率。

(11)

其中,ai(m)為第i個對象的第m個影響因素。將庫水位和變形量的數(shù)據(jù)以及降雨量和變形量的數(shù)據(jù)代入式(11)中，建立兩組比較數(shù)列。

b)分析兩組數(shù)據(jù)，建立相應的參考數(shù)列。確定反應系統(tǒng)特征的參考數(shù)列，以各指標的最優(yōu)值構成參考數(shù)據(jù)列，見式(12)：

a0=(a0(1),…,a0(m))

(12)

c)數(shù)據(jù)初始化，因系統(tǒng)中各因素的物理意義不同，進行灰色關聯(lián)度分析時，本實驗采用初值化法進行無量綱化，見式(13)：

(13)

d)計算參考數(shù)列對應元素的差的絕對值，見式(14)：

Δ=|t0(k)-ti(k)|，i=1,2…,n;k=1,2…,m

(14)

e)計算關聯(lián)系數(shù)，見式(15)

(15)

式(15)中，miniminkΔ表示查找每個對象，對應的影響因素的最小值；ρ為分辨系數(shù)，取值范圍為0<ρ<1，其是為了減小由于Δ最大值的數(shù)值太大而造成的失真，提高關聯(lián)系數(shù)之間的差異顯著性[15]，本實驗取ρ為0.5。

f)計算關聯(lián)度，分別計算各個參考序列對應元素的關聯(lián)系數(shù)的均值，見式(16)：

(16)

最終由式(16)得到每期滑坡變形量分別與當期月降雨量、月平均庫水位變化量關聯(lián)度，關聯(lián)結果見表1、2。

表1 滑坡變形量與庫水位關聯(lián)度

表2 滑坡變形量與降雨量關聯(lián)度

若因素趨勢曲線越相近，則其關聯(lián)程度越強，關聯(lián)度越強趨近于1，反之越弱越趨近于0。由12期的隨機試驗樣本可見，滑坡變形量與庫水位變化量的最大關聯(lián)度為0.987，最小關聯(lián)度為0.700，平均關聯(lián)度為0.877。與累計月降雨量的最大關聯(lián)度為0.995，最小關聯(lián)度為0.667，平均關聯(lián)度為0.899?？梢?，降雨和庫水位的變化對滑坡的形變具有較強的關聯(lián)作用，且?guī)焖蛔兓瘜滦巫兊挠绊懘笥诮涤陮滦巫兊挠绊?。此外，降雨對滑坡的形變存在滯后效應，因此在使用BP神經(jīng)網(wǎng)絡對高頻IMF分量學習訓練時，需要加入兩月累計降雨量作為新的影響因子。

2.2 滑坡位移分解與預測

a)根據(jù)試驗步驟，首先進行滑坡變形的位移分解，通過變分模態(tài)分解法分解72期白水河滑坡累計位移變形數(shù)據(jù)，前60期作為訓練學習樣本，后12期作為預測對比樣本。為避免產(chǎn)生多余的虛假分量，且保證分解出的IMF分量均有物理意義，通過試驗最后確定分解模態(tài)數(shù)為K=3，懲罰參數(shù)C=0.5和上升步長τ=0.2。因而得到3個模態(tài)分量，結果見圖 2。

a)原始位移

b)IMF1隨機分量

c)IMF2周期分量

d)IMF3趨勢分量

圖2中， IMF1、IMF2分量為高頻部分，代表了滑坡位移的隨機噪聲和周期部分，IMF3分量平穩(wěn)遞增，為滑坡位移的趨勢部分。

b)對IMF1、IMF2分量通過BP神經(jīng)網(wǎng)絡進行學習并預測。根據(jù)2.1節(jié)中灰色關聯(lián)分析，選用7個影響因子作為BP神經(jīng)網(wǎng)絡的輸入值，分別為：前60期對應的月累計降雨量，兩月累計降雨量，月平均庫水位變化量以及每期滑坡變形量的前四期數(shù)據(jù)，因而輸入節(jié)點數(shù)為7，輸出節(jié)點數(shù)為1。隱含層傳遞函數(shù)設置為purelin，學習速率為0.01，通過實驗設定隱含層數(shù)為3，學習結果見圖3。通過擬合學習以后，根據(jù)對應影響因素得到后12期IMF1、IMF2分量的預測值。

c)對IMF3趨勢分量的擬合由1.3節(jié)中二次指數(shù)平滑進行擬合，通過試驗，選取β=0.15，可取得較好的擬合效果。根據(jù)式(7)最終得到后12期趨勢分量預測值。

d)最終，將各自IMF分量的預測值對應相加得到滑坡位移預測值，并與真實值進行對比分析。

a)IMFI

b)IMF2

2.3 試驗驗證

為證實本模型的預測精度，采用均方根誤差(RMSE)，平均相對百分誤差(MAPE)，相關系數(shù)(R2)作為定量評估標準，見式(17)—(19)：

(17)

(18)

(19)

表3 試驗數(shù)據(jù)對比

最終在本模型中，總位移均方根誤差(RMSE)為3.14 cm，平均相對百分誤差(MAPE)為0.78%，總位移相關系數(shù)(R2)為0.986，可見本模型預測精度較高，在滑坡位移預測方面具有一定的應用價值。

3 結論

本實驗針對三峽地區(qū)“階梯狀”位移滑坡，提出了一種新的擬合預測模型：VMD-SES-BP模型。在模型中，將滑坡位移時間序列作為數(shù)字信號，通過變分模態(tài)分解將該信號分解為高頻序列IMF1、 IMF2分量以及低頻序列IMF3分量。通過灰色關聯(lián)分析法，證實三峽地區(qū)的降雨和庫水位的變化是滑坡演化過程中的主要觸發(fā)因素。將歷史監(jiān)測的60期月累計降雨量、兩月累計降雨量、月平均庫水位變化量和前4期位移變形值共同作為BP神經(jīng)網(wǎng)絡的輸入層，通過訓練學習后得到后12期IMF1、IMF2分量預測值。同時采用二階指數(shù)平滑對IMF3分量進行擬合并預測，最終將各分量預測值累加得到后12期滑坡位移預測值。通過精度分析，RMSE為3.14 cm，MAPE為0.78%，R2為0.986，證實了該模型針對三峽地區(qū)“階梯狀”滑坡位移預測上具有一定的應用價值。

此外，從實驗的結果和趨勢分析，此類“凸躍”位移的變化大部分集中在每年的7、8、9月份，在此期間，三峽地區(qū)的降雨量較大，增大了滑坡的重力勢能，加劇了滑坡形變的產(chǎn)生。同時，受水流滲透進巖體的速率影響，三峽水庫水位的下降或上漲使得地下水位與庫水位形成正負落差，同樣影響著滑坡形變量。因此，未來對三峽滑坡形變的預測預報中，應著重分析降雨和庫水位的變化。