陳 曦，高雅萍，涂 銳

(1.成都理工大學(xué)地球科學(xué)學(xué)院，四川 成都 610059；2.中國科學(xué)院國家授時中心，陜西 西安 710600)

滑坡監(jiān)測得到的地表位移數(shù)據(jù)經(jīng)過科學(xué)處理和有效預(yù)測，可以得到滑坡體未來一段時間內(nèi)的位移趨勢變化，這對滑坡災(zāi)害防災(zāi)減災(zāi)方案的制訂有重要的數(shù)據(jù)決策支持作用，對于滑坡災(zāi)害的預(yù)報預(yù)警具有非常重要的意義。

目前，運用較多的預(yù)測模型有統(tǒng)計機器學(xué)習(xí)模型、灰色模型、回歸模型、時間序列模型等[1,4]，其中神經(jīng)網(wǎng)絡(luò)模型、灰色模型、自回歸模型等在滑坡位移預(yù)測中應(yīng)用較多[2]；數(shù)據(jù)序列的分解方法有小波分解、EMD、EEMD、VMD等[1,3,5-6]；許多學(xué)者將數(shù)據(jù)分解方法和預(yù)測模型進行組合預(yù)測，預(yù)測精度相對提高。宋麗偉[5]對位移序列進行經(jīng)驗?zāi)B(tài)分解和模態(tài)重構(gòu)，利用LSTM模型組合預(yù)測，為“階梯狀”滑坡位移的預(yù)測提供了一種可行的思路。鄢濤等[6]利用極限學(xué)習(xí)機(ELM)對大壩變形的EEMD不同分解分量引入影響因子分別預(yù)測并疊加，預(yù)測精度有所提高。但這些預(yù)測模型都需要較多的累積觀測數(shù)據(jù)進行樣本學(xué)習(xí)，在實際監(jiān)測過程中當(dāng)影響因子監(jiān)測數(shù)據(jù)積累不足或者影響因子數(shù)據(jù)誤差較大與位移相關(guān)性較少時，位移預(yù)測就會比較困難，預(yù)測精度較低，需要應(yīng)用對單位移序列擬合預(yù)測效果較好的模型，使滑坡位移預(yù)測結(jié)果有較高的精度。吳棟等[7]對比研究了LSTM 和SVM對VMD分量位移的組合預(yù)測發(fā)現(xiàn)LSTM更能有效地預(yù)測滑坡位移變化的本質(zhì)規(guī)律。簡文彬等[8]、楊大明[9]利用門限自回歸模型進行了地基沉降預(yù)測。針對這樣的問題，選擇了應(yīng)用較廣泛的經(jīng)驗?zāi)B(tài)分解(Empirical Mode Decomposition,EMD)分解方法和對非線性非穩(wěn)態(tài)序列預(yù)測效果較好的門限自回歸模型(Threshold Auto Regressive,TAR)[10-12]，以組合的形式預(yù)測滑坡位移變化。門限自回歸模型[13]是在自回歸 (AR) 模型的基礎(chǔ)上增加了門限區(qū)間約束條件,是對一類非線性的時間序列進行局部線性逼近的非線性時間序列模型。

本文以三峽庫區(qū)秭歸縣白水河滑坡為例，首先利用經(jīng)驗?zāi)B(tài)分解方法將GPS滑坡地表監(jiān)測水平累計位移分解為不同頻率位移分量，然后利用門限自回歸模型對各位移分解序列進行獨立預(yù)測，最后模態(tài)疊加得到最終預(yù)測位移序列。選取MAE、MRE、RMSE和擬合優(yōu)度4種精度指標(biāo)來準(zhǔn)確評價組合模型的精度，通過與BP(Back Propagation)神經(jīng)網(wǎng)絡(luò)模型、長短時間記憶網(wǎng)絡(luò)模型的單一或組合模型對比驗證其準(zhǔn)確性，為滑坡的位移預(yù)測提供一種新的方法。

1 主要方法

1.1 經(jīng)驗?zāi)B(tài)分解

經(jīng)驗?zāi)B(tài)分解(Empirical Mode Decomposition,簡稱EMD)以數(shù)據(jù)序列的時變特征為基礎(chǔ)，可將序列中存在的各種具有一定周期性規(guī)律的不同變化頻率的分量分離單獨存在。

對于滑坡位移數(shù)值序列，由于各種因素的影響，實測數(shù)據(jù)會存在各種模式的波動變化，也就是由若干個不同頻率的數(shù)值波形序列所混合構(gòu)成的，是不同變化頻率數(shù)值序列的復(fù)合非線性數(shù)據(jù)，EMD分解的目的就是為了分解出其中的位移本征模函數(shù)，將復(fù)雜位移序列分化為多個近似單一頻率的位移序列形式。

EMD分解方法：在原位移序列極值點基礎(chǔ)上用三次樣條曲線擬合形成上下包絡(luò)線，計算上下包絡(luò)均值作為第一個位移序列分離分量，將原位移序列減去該分量，得到1個新的位移序列。新序列若還存在局部極值，說明還不是1個無波動趨勢序列，需要繼續(xù)進行擬合分解，使用上述方法得到第一個分量后，用原始序列減去該分量，作為新的數(shù)據(jù)序列，可以得到所有剩余分量，以此類推，直至完成經(jīng)驗?zāi)B(tài)分解得到趨勢變化序列。

1.2 門限自回歸模型

門限自回歸(Threshold Auto Regressive,TAR)模型對于存在較規(guī)律波動變化的時間序列預(yù)測效果較好。由于門限的分割，該模型對數(shù)據(jù)序列可以分段處理，詳細描述序列的變化情況。

自回歸(Auto Regressive,TAR)模型是分析數(shù)據(jù)內(nèi)部相關(guān)關(guān)系的方法，假設(shè)存在觀測時序{xt}，則定義p階自回歸模型，簡稱AR(p)：

xt=φ0+φ1xt-1+…+φpxt-p+εt

(1)

式中，φp不為0，{φp}為自回歸系數(shù)序列；εt為零均值白噪聲序列。

φ0+φ1+…+φp<1對于平穩(wěn)的自回歸模型，其自相關(guān)系數(shù)滿足如下遞推關(guān)系：

ρk=φ1ρk-1+φ2ρk-2+…+φpρk-p

(2)

(3)

AR模型的參數(shù)估計采用最小二乘估計，估計自回歸系數(shù)序列{φp}，建立各參數(shù)矩陣,Y為AR模型預(yù)測矩陣，X為AR模型樣本矩陣，φ為AR模型自回歸系數(shù)矩陣，ε為零均值白噪聲序列：

(4)

線性方程組為：

Y=Xφ+ε

(5)

自回歸系數(shù)最小二乘估計為：

φ=(XTX)-1XT(Y-ε)

(6)

誤差方差的最小二乘估計為：

(7)

AR模型的定階采用貝葉斯信息準(zhǔn)則：

(8)

其中,N表示樣本序列長度,σ2表示樣本殘差的方差，使得BIC達到最小值的p即為該準(zhǔn)則下的最優(yōu) AR 模型的階數(shù)。

TAR模型是分段的AR模型，通過在觀測時序{xi}的取值范圍內(nèi)引入k-1個門限值ri(i=1,2,…,n-1)，將該整個時間序列分成n個門限區(qū)間，可用r0,rn分別表示上界和下界，并根據(jù)延遲步數(shù)d按{xi-d}值的大小分配到不同的門限區(qū)間內(nèi)，再對區(qū)間內(nèi)的xi采用不同階數(shù)的自回歸模型(AR模型)，從而形成時間序列的非線性分段動態(tài)詳細描述，其模型形式見式(9)：

(9)

(10)

2 EMD-TAR組合預(yù)測模型

首先利用經(jīng)驗?zāi)B(tài)分解方法對GPS滑坡監(jiān)測地表位移數(shù)據(jù)進行序列分解，得到不同變化頻率的位移序列，對各變形分量獨立采用門限自回歸模型進行預(yù)測，避免各分量誤差相互影響，最后將所有形變分量模態(tài)疊加得到最終滑坡地表位移數(shù)值預(yù)測序列。其優(yōu)勢是利用門限自回歸模型的對復(fù)雜非線性數(shù)據(jù)序列預(yù)測效果較好的特點，結(jié)合EMD分解算法的只用依據(jù)序列自身時間特征的特點，極大地提高了滑坡形變監(jiān)測單位移序列的精度。該文組合預(yù)測模型數(shù)據(jù)處理步驟見圖1。

圖1 組合模型預(yù)測分析步驟

所采用的滑坡位移預(yù)測模型的精度評價指標(biāo)包括平均絕對誤差(Mean Absolute Error，MAE)、均方誤差(Root Mean Square Error，RMSE)、平均相對誤差(Mean Relative Error，MRE)、擬合優(yōu)度R2,計算式為:

(11)

(12)

(13)

(14)

3 工程應(yīng)用及模型驗證

以三峽庫區(qū)秭歸縣白水河滑坡GPS地表位移監(jiān)測數(shù)據(jù)為例，滑坡體呈階梯狀向長江展布，為堆積層滑坡，坡體屬順向坡?；麦w共布設(shè)11個GPS監(jiān)測點，選擇滑坡體中下部監(jiān)測點ZG118和滑坡體北東部監(jiān)測點DX-02對滑坡體不同位置形變位移進行預(yù)測描述，數(shù)據(jù)來源于國家冰川凍土沙漠科學(xué)數(shù)據(jù)中心 (http://www.ncdc.ac.cn)。

ZG118監(jiān)測點為從2006年12月至2012年12月，共73期GPS地表位移觀測數(shù)據(jù)。為提高單序列位移的預(yù)測精度，首先利用EMD對GPS監(jiān)測點水平位移序列進行模態(tài)分解，迭代次數(shù)設(shè)置為200次，得到5個不同變化頻率的IMF分量和一個殘差分量，殘差分量作為粗差剔除。各IMF分量見圖2。

a)IMF1

將73期滑坡位移數(shù)據(jù)中的1—50期作為訓(xùn)練樣本，51—73期數(shù)據(jù)作為預(yù)測檢驗值。建立門限自回歸預(yù)測模型，門限區(qū)間個數(shù)取2，最大門限延遲量設(shè)置為5，自回歸最大階數(shù)為5，默認最小AIC值為1E+10,便于后續(xù)優(yōu)化AIC值獲得最優(yōu)多維參數(shù)。按樣本的30%～70%分位區(qū)間以1%進度搜索最優(yōu)門限值。每組分量預(yù)測步長為一步，預(yù)測后數(shù)據(jù)更新，再進行下一步預(yù)測，實現(xiàn)模型的動態(tài)更新和預(yù)測。

為驗證所建立組合預(yù)測模型的精度，分別建立EMD-TAR模型、EMD-LSTM模型、EMD-BPNN模型、TAR模型、LSTM模型對該滑坡位移序列進行預(yù)測，組合模型EMD-TAR模型、EMD-LSTM模型、EMD-BPNN模型[14-15]的預(yù)測方法為分別利用TAR模型、LSTM模型、BPNN模型預(yù)測每個IMF分量，最后時序合成得到實際滑坡預(yù)測位移序列，各模型的預(yù)測結(jié)果見圖3。

圖3 各模型預(yù)測結(jié)果曲線序列

由圖3可看出，單一TAR模型的預(yù)測效果比LSTM和BPNN更好，LSTM和BPNN單一模型會隨著時間推移預(yù)測誤差變大。相對于單一TAR模型和其他組合模型，組合EMD-TAR模型與原始位移序列擬合度較高，擬合優(yōu)度最高，較為符合原始位移的變化趨勢。為了準(zhǔn)確評價模型的整體預(yù)測精度，對比計算幾種預(yù)測模型的精度評價指標(biāo)，結(jié)果表明，所建立的EMD-TAR預(yù)測模型在MAE、RMSE、MRE和擬合優(yōu)度上均大幅提高。各模型精度評價值見表1。

表1 各模型精度評價

為更好地驗證所建立預(yù)測模型的精度和可靠性，選取滑坡體另一監(jiān)測點進行預(yù)測實驗，DX-02監(jiān)測點為從2006年12月至2010年8月，共45期GPS地表位移觀測數(shù)據(jù)，1—32期作為訓(xùn)練樣本，33—45期數(shù)據(jù)作為預(yù)測檢驗值。EMD模態(tài)分解迭代次數(shù)設(shè)置為2 000次，得到5個不同變化頻率的IMF分量和一個殘差分量，殘差分量作為粗差剔除。各IMF分量見圖4。所建立門限自回歸預(yù)測模型中門限區(qū)間個數(shù)設(shè)置為2，最大門限延遲量設(shè)置為5，自回歸最大階數(shù)為5，默認最小AIC值為1E+10。考慮到分量頻率的變化因素，其中IMF1按樣本的30%～70%分位區(qū)間以1%進度搜索最優(yōu)門限值，IMF2-4按樣本的30%～50%分位區(qū)間以1%進度搜索最優(yōu)門限值。每組分量預(yù)測步長為一步，預(yù)測后數(shù)據(jù)更新，再進行下一步預(yù)測，實現(xiàn)模型的動態(tài)更新和預(yù)測。

a)IMF1

同樣，以此監(jiān)測數(shù)據(jù)為例，分別建立EMD-TAR模型、EMD-LSTM模型、EMD-BPNN模型、TAR模型、LSTM模型對該滑坡位移序列進行預(yù)測，組合模型EMD-TAR模型、EMD-LSTM模型、EMD-BPNN模型的預(yù)測方法為分別利用TAR模型、LSTM模型、BPNN模型預(yù)測每個IMF分量，最后時序合成得到實際滑坡預(yù)測位移序列，各模型的預(yù)測結(jié)果見圖5。

圖5 各模型預(yù)測結(jié)果曲線序列

由圖5可看出，組合EMD-TAR模型與原始位移序列擬合程度最高，最符合原始位移的變化趨勢。參考分析2個算例的各模型預(yù)測結(jié)果曲線圖，樣本數(shù)據(jù)量和預(yù)測時間的變化對組合EMD-TAR模型影響較小，模型預(yù)測較穩(wěn)定。在算例一中EMD-TAR、EMD-LSTM模型預(yù)測結(jié)果誤差較小而在算例二中EMD-TAR、LSTM模型預(yù)測結(jié)果誤差較小，樣本數(shù)據(jù)量變化對LSTM及其組合模型的預(yù)測精度產(chǎn)生了影響，算例二較算例一數(shù)據(jù)量減少，EMD分解后的不同頻率分量變少，算例二滑坡位移的整體變化趨勢較算例一小，所以在算例二中，對未分解的滑坡位移進行預(yù)測的LSTM模型的預(yù)測效果相比EMD-LSTM模型更好，綜合分析得到EMD-TAR在預(yù)測精度和穩(wěn)定性上更優(yōu)。

為了準(zhǔn)確評價模型的整體預(yù)測精度，計算幾種預(yù)測模型的精度評價指標(biāo)，結(jié)果表明，該文所建立的EMD-TAR預(yù)測模型在MAE、RMSE、MRE和擬合優(yōu)度上均大幅提高。各模型精度評價值見表2。

表2 各模型精度評價

4 結(jié)論

基于EMD算法和TAR模型構(gòu)建了EMD-TAR滑坡位移預(yù)測模型，針對的是在實際監(jiān)測過程中當(dāng)影響因子監(jiān)測數(shù)據(jù)積累不足或者影響因子數(shù)據(jù)誤差較大與位移相關(guān)性較少的情況,利用EMD算法將非穩(wěn)態(tài)的GPS監(jiān)測滑坡位移序列分解為在時間尺度上具有不同變化頻率的多個分量，對AR模型進行門限改進得到TAR模型，結(jié)合TAR模型對非線性波動數(shù)據(jù)序列預(yù)測效果較好的優(yōu)勢，對滑坡位移序列進行分解、預(yù)測、再疊加得到精度較高的預(yù)測結(jié)果。經(jīng)過白水河滑坡實例數(shù)據(jù)的試驗驗證，該組合模型相比于單一TAR模型、LSTM和BPNN模型及其組合模型預(yù)測精度更高，為滑坡位移的預(yù)測提供了一種新方法。但該模型在滑坡環(huán)境影響下監(jiān)測數(shù)據(jù)積累不足或影響因子關(guān)聯(lián)較小時對位移的預(yù)測具有一定的局限性，后續(xù)可與其他顧及影響因子的模型進行組合預(yù)測研究。