基于優(yōu)化經(jīng)驗?zāi)B(tài)分解和最小二乘支持向量機(jī)的邊坡位移預(yù)測

易智文

（萍鄉(xiāng)市山口巖水利樞紐管理中心，江西 萍鄉(xiāng)，337000）

0 前言

水庫邊坡滑坡是普遍存在的地質(zhì)災(zāi)害，嚴(yán)重威脅庫區(qū)內(nèi)人民正常的生產(chǎn)生活[1]。為降低滑坡災(zāi)害風(fēng)險，有必要開發(fā)相應(yīng)的預(yù)警預(yù)報系統(tǒng)。受地質(zhì)結(jié)構(gòu)、氣候條件等影響，邊坡位移時間序列具有顯著的非線性非穩(wěn)態(tài)特征[2]，若能表征這種非線性演化規(guī)律，則有助于實現(xiàn)位移時序的精準(zhǔn)預(yù)測。

目前，常用的邊坡位移預(yù)測方法可分為數(shù)學(xué)統(tǒng)計模型和機(jī)器學(xué)習(xí)模型。數(shù)學(xué)統(tǒng)計模型包括灰色預(yù)測模型、自回歸模型等[3]。機(jī)器學(xué)習(xí)模型近年來得到眾多學(xué)者青睞，常用的方法包括神經(jīng)網(wǎng)絡(luò)模型、支持向量機(jī)模型、長短期記憶神經(jīng)網(wǎng)絡(luò)模型、循環(huán)神經(jīng)網(wǎng)絡(luò)模型、隨機(jī)森林模型等[4]。通過這些方法，研究者們?nèi)〉昧嗽S多卓有成效的成果，但仍存在不足，如仍難以反映位移時序的非線性特征。一些學(xué)者將信號分解方法嵌入預(yù)測模型以解決上述問題，常用的方法包括小波分析、經(jīng)驗?zāi)B(tài)分解、變分模態(tài)分解等[5]，但這類方法分解結(jié)果隨機(jī)性較大，誤差較高，難以賦予分量實際物理意義。

基于此，本文引入經(jīng)軟篩分停止準(zhǔn)則優(yōu)化后的經(jīng)驗?zāi)B(tài)分解，結(jié)合最小二乘支持向量機(jī)，實現(xiàn)非線性位移時序精準(zhǔn)預(yù)測。該模型首先自適應(yīng)地將位移-時間曲線分解為若干分量，采用K-means 法將其聚類為趨勢性位移、周期性位移和隨機(jī)性位移。而后，再分別采用最小二乘法和LSSVM 模型對這3 種位移分別擬合預(yù)測。最后，累加求和3 種預(yù)測位移，即可得到累計位移預(yù)測值。通過對山口巖水庫監(jiān)測位移進(jìn)行預(yù)測，發(fā)現(xiàn)SSSC-EMD-LSSVM 模型預(yù)測效率和精度高，是一種可靠的非線性位移時序預(yù)測模型。

1 方法與原理

1.1 基于軟篩分停止準(zhǔn)則的經(jīng)驗?zāi)B(tài)分解

經(jīng)驗?zāi)B(tài)（Empirical mode decomposition，EMD）分解[6]是一種自適應(yīng)信號處理方法，可高效、精確地處理非線性非平穩(wěn)數(shù)據(jù)，故被廣泛應(yīng)用于各個領(lǐng)域。經(jīng)驗?zāi)B(tài)分解不受基函數(shù)約束，可將處理對象分解成不同尺度、相互獨立的本征模態(tài)函數(shù)（Intrinsic mode function，IMF），彌補(bǔ)了傅里葉變換和小波變換的局限性。在分解過程中，本征模態(tài)函數(shù)需嚴(yán)格滿足2 個條件：①對于整個時間序列，函數(shù)局部極值點和零穿越點數(shù)目相等或最多相差1 個；②對于任意時間序列，函數(shù)局部最大值點連成的上包絡(luò)線和最小值點連成的下包絡(luò)線均值為零。假定有一組初始時間序列x（t），EMD 分解過程如下：

（1）尋找局部極大值和極小值點，采用三次樣條函數(shù)進(jìn)行銜接，計算上包絡(luò)線emax（t）和下包絡(luò)線emin（t）的平均值m（t），計算公式如下：

（2）將初始時間序列x（t）與均值m（t）相減，得到差值h（t），見式（2）。若差值h（t）不滿足本征模態(tài)函數(shù)基本條件，則重復(fù)上述步驟。若差值滿足條件，則將h（t）作為第1 個本征模態(tài)分量IFM1，記為C1（t），將初始時間序列x（t）與該IMF 的差值作為剩余分量，記為r1（t）。

（3）將剩余分量r1（t）作為新的時間序列，對其重復(fù)步驟（1）～（2），即可得到本征模態(tài)分量IMF2，IMF3，…，IMFn，和1 個最終剩余分量rn（t）。根據(jù)上述分解過程，初始時間序列x（t）可表示為：

包絡(luò)曲線擬合參數(shù)、邊界條件參數(shù)和篩分停止標(biāo)準(zhǔn)參數(shù)直接影響分解精度以及效率。目前，針對包絡(luò)曲線和邊界條件參數(shù)已開展大量工作，少有研究報道了篩分停止準(zhǔn)則參數(shù)的優(yōu)化選取方法。在EMD 分解過程中，若篩選迭代次數(shù)太少，即“欠篩”，則可能導(dǎo)致模態(tài)分量包含過多單分量信號；若篩選迭代次數(shù)過多，即“過篩”，則使得模態(tài)分量包含互不相關(guān)的故障信號。因此，以往的研究常預(yù)先設(shè)置篩分閾值，但這一過程嚴(yán)重取決于研究者的經(jīng)驗，可能會導(dǎo)致不準(zhǔn)確的分解結(jié)果。彭丹丹等[7]提出了一種軟篩分停止準(zhǔn)則，可以自適應(yīng)優(yōu)化選取篩分參數(shù)，在一定程度上彌補(bǔ)了上述不足。

基于軟篩分停止準(zhǔn)則的經(jīng)驗?zāi)B(tài)分解考慮整體能量特性和局部沖擊特性，首先需要定義目標(biāo)函數(shù)fik描述包絡(luò)均值信號，如下式：

式中：均方根RMSik和超峭度EKik可分別由式（6）和（7）計算：

式中：n 為采樣點，總數(shù)由Ns表示；mik（t）為第i 個IMF 經(jīng)過k 次篩分后的包絡(luò)線的均值；為mik（t）的算術(shù)平均值。

根據(jù)既定目標(biāo)函數(shù)，提出啟發(fā)式機(jī)制自適應(yīng)地確定每次篩分過程的最優(yōu)（或次優(yōu)）篩分迭代次數(shù)，評價指標(biāo)定義為初始判斷系數(shù)P1、過程判斷系數(shù)P2和結(jié)果判斷系數(shù)P。軟篩分停止準(zhǔn)則流程如下：

（1）預(yù)設(shè)最大篩分迭代次數(shù)Imax，令第i 個IMF 分量經(jīng)k-1 次迭代運行后得到的平滑信號hik（t）等于ri（t）。其中，ri（t）表示初始時間序列x（t）減去第i 個分量ci（t）后的殘余項。而后，令目標(biāo)函數(shù)fik、初始判斷系數(shù)P1及迭代次數(shù)k 均設(shè)置為0。其中，P1表示為fik-1小于fik時，輸出為1，則進(jìn)入步驟（2）。

（2）令f=fik，P0=P1，k=k+1，用hik-1計算mik（t），再將hik-1（t）減去mik（t）得到hik（t）。而后，計算fk，若fik〉f，則P1=1，否則令P1=0，P=P0+P1，進(jìn)入步驟（3）。

（3）判斷是否符合篩分停止的條件：①極值點個數(shù)和零點個數(shù)相差不超過1；②fik-2〈fik-1且fik-1〈fik。若不滿足條件，則返回步驟（2）；若滿足條件，則進(jìn)入步驟（4）。

（4）對比篩分次數(shù)k 和預(yù)設(shè)篩分最大次數(shù)Imax。若k〈Imax，則將hik（t）輸出為第i 個模態(tài)分量IFMi，篩分結(jié)束；若k〈Imax，則進(jìn)入步驟（5）。

（5）判定結(jié)果判斷系數(shù)P，若P〈2，則返回步驟（2）；若P≥2，則將hik-2（t）輸出為第i 個模態(tài)分量IFMi，篩分結(jié)束。

基于啟發(fā)式機(jī)制的篩分停止準(zhǔn)則能夠自適應(yīng)確定最優(yōu)迭代次數(shù)，避免了人為干擾，解決了“欠篩”和“過篩”問題，提高了分解效率和精度。因此，本文采用經(jīng)過軟篩分停止準(zhǔn)則優(yōu)化后的經(jīng)驗?zāi)B(tài)分解（SSSC-EMD）對位移時間序列進(jìn)行分解。

1.2 最小二乘支持向量機(jī)

最小二乘支持向量機(jī)（Least squares support vector machine，LSSVM） 是一種改進(jìn)支持向量機(jī)（Support vector machine，SVM）[8]，不僅能夠解決SVM 存在的二次回歸的問題，還能提高計算效率和預(yù)測精度。LSSVM使用不同的最優(yōu)目標(biāo)函數(shù)，引入平等約束，并采用平方誤差替換原始損失函數(shù)以及控制計算數(shù)量，整體結(jié)構(gòu)具有樣本小、風(fēng)險小的特點。LSSVM 具體實現(xiàn)流程如下：

假設(shè)樣本集長度為l，樣本集T=（｛xi，y）i｝ （i=1，2，…，l），輸入量為xi∈Rn，輸出量為yi∈R。采用非線性函數(shù)φ（x）將樣本集映射至高維空間，則最優(yōu)線性回歸函數(shù)表達(dá)式為：

式中：w 為權(quán)向量；b 為偏置量。

與VSM 不同的是，LSSVM 基于結(jié)構(gòu)風(fēng)險最小原則，將誤差平方ξi作為模型的損失函數(shù)，將約束條件轉(zhuǎn)換為等式約束，則LSSVM 問題的目標(biāo)函數(shù)J 為：

式中：ξi為松弛因子；C 為正則化參數(shù)，建立拉格朗日方程優(yōu)化求解：

式中：ai為拉格朗日乘子?；诰€性方程條件，對L（w，b，ξ，a）求各階偏導(dǎo)，令各偏導(dǎo)分量為0 并代入式（10），消去w 和ξ 項，得到b 和a：

式中，I 為l 階單位矩陣，y=[y1，y2，…，yl]T，H=[1，2，…，l]T，a=[a1，a2，…，al]T。K 為核函數(shù)，滿足Mercer 條件，本次研究采用徑向基核函數(shù)（RBF）：

式中：σ 為核函數(shù)寬度。將高維特征空間中的點積運算替換為原空間中的核函數(shù)，得到LSSVM 回歸函數(shù)為：

懲罰系數(shù)ξ 和核函數(shù)寬度σ 直接影響LSSVM 回歸模型的精度。為此，本文采用交叉驗證法進(jìn)行尋優(yōu)。

1.3 邊坡位移預(yù)測流程

根據(jù)基于軟篩分停止準(zhǔn)則的經(jīng)驗?zāi)B(tài)分解和最小二乘支持向量機(jī)，對邊坡位移時間序列進(jìn)行預(yù)測，具體流程見圖1。

圖1 基于SSSC0-EMD 和LSSVM 的位移預(yù)測流程圖

（1）采用SSSC-EMD 自適應(yīng)分解位移-時間曲線。一般而言，邊坡位移由趨勢性位移、周期性位移和隨機(jī)性位移組成。殘余分量與趨勢項對應(yīng)，IMF 分量與周期性位移和隨機(jī)性位移對應(yīng)。需要指出的是，若分解出的IMF 數(shù)大于所需變量數(shù)，則采用K-means 法進(jìn)行分類，然后將同一類分量疊加求和，賦予其物理意義（即周期項或隨機(jī)項）。

（2）分別在趨勢性位移、周期性位移和隨機(jī)性位移內(nèi)劃分訓(xùn)練集和測試集。

（3）采用最小二乘法擬合趨勢性位移，輸出位移預(yù)測結(jié)果，比較預(yù)測值與實測值。

（4）采用灰色關(guān)聯(lián)分別尋找周期性位移集和隨機(jī)性位移集的輸入歷史序列數(shù)和位移預(yù)測數(shù)據(jù)點數(shù)。而后，建立LSSVM 模型，應(yīng)用交叉驗證進(jìn)行參數(shù)尋優(yōu)。最后，按照既定比例不斷迭代計算，建立效果最好的預(yù)測模型。將測試集輸入至對應(yīng)的預(yù)測模型，得到位移預(yù)測結(jié)果，比較預(yù)測值與實測值。

（5）將各預(yù)測位移分量累加求和，實現(xiàn)邊坡位移預(yù)測。

為評價模型的預(yù)測精度，引入均方根誤差（RMSE），平均絕對誤差（MAE），判定系數(shù)（R2）：

式中：yi為實際監(jiān)測位移；為位移預(yù)測；為實際監(jiān)測位移均值。

2 工程實例分析

2.1 工程背景

山口巖水利樞紐地處江西省萍鄉(xiāng)市蘆溪縣，位于贛江支流袁河上游。發(fā)電廠房位于下游河道約200 m處，廠址區(qū)自然邊坡坡度約25° ～35°，坡角約50°。運營過程中，廠房邊坡出現(xiàn)顯著蠕滑裂縫，極易誘發(fā)滑坡災(zāi)害。因此，亟待開展滑坡位移監(jiān)測、預(yù)測。

嚴(yán)格按照《混凝土壩安全監(jiān)測技術(shù)規(guī)范》（DL/T5178）采用TS15 全站儀開展邊坡位移監(jiān)測，對W14測點59 個月（2018 年5 月至2023 年3 月）的位移監(jiān)測數(shù)據(jù)進(jìn)行研究。W14 測點布置圖如圖2 所示，位移監(jiān)測數(shù)據(jù)如圖3 所示。訓(xùn)練集和測試集的劃分嚴(yán)重影響模型結(jié)果，研究表明訓(xùn)練集長度越長，測試精度越高，但相應(yīng)的計算時間更長。為協(xié)調(diào)計算精度和時間，依據(jù)前人研究經(jīng)驗[9]，將訓(xùn)練集與測試集的比例設(shè)置為8：2，即前47 個月（2018 年5 月至2022 年3 月）的數(shù)據(jù)作為訓(xùn)練集，后12 個月（2022 年4 月至2023 年3 月）的數(shù)據(jù)作為測試集。

圖2 監(jiān)測點布置圖

圖3 W14 測點累計位移

2.2 位移分解

采用SSSC-EMD 分解位移時間序列，結(jié)果如圖4所示。從圖中可知，W14 測點累計位移時間序列共分解出2 項高頻分量IMF1和IMF2，1 項低頻分量IMF3，1 項殘余分量r。由于IMF 分量數(shù)為3，故采用K-means 將其聚類為2 類，選用Sqeuclidean 作為聚類距離計算方法。聚類輪廓系數(shù)是表征聚類效果的指標(biāo)之一，聚類輪廓系數(shù)越高，分量之間相關(guān)程度越高。根據(jù)計算結(jié)果可知，若將IMF1和IMF2聚為一類，MF3單獨為一類，聚類輪廓系數(shù)為0.73，表明聚類效果較好[10]。因此，將高頻分量IMF1和IMF2累加求和，作為周期性位移；將IMF3作為隨機(jī)性位移；將殘余分量r 作為趨勢性位移。

圖4 W14 測點累計位移SSSC-EMD 分解結(jié)果

2.3 位移預(yù)測

2.3.1 趨勢性位移預(yù)測

結(jié)合圖5 可知，隨著時間的推移，殘余分量r 表現(xiàn)為單調(diào)遞增的趨勢，故可采用最小二乘法對其進(jìn)行擬合預(yù)測。結(jié)果表明，3 次多項式擬合結(jié)果精度最高，RMSE 為0.005 mm，MAE 為0.004 mm，R2達(dá)1.00。

圖5 趨勢性位移預(yù)測結(jié)果

2.3.2 周期性位移預(yù)測

建立LSSVM 模型預(yù)測周期性位移時間序列。首先，采用灰色關(guān)聯(lián)法確定歷史序列數(shù)為5，位移預(yù)測數(shù)據(jù)點數(shù)為1，再運用交叉驗證法對模型參數(shù)尋優(yōu)，C=1.26，σ=1.22，最后訓(xùn)練模型并輸出預(yù)測值。預(yù)測結(jié)果如圖6 所示，通過與監(jiān)測結(jié)果比較可知，RSME 為0.037 mm，MAE 為0.033 mm，R2為0.96。

圖6 周期性位移預(yù)測結(jié)果

2.3.3 隨機(jī)性位移預(yù)測

建立LSSVM 模型預(yù)測隨機(jī)性位移時間序列。首先，采用灰色關(guān)聯(lián)法確定歷史序列數(shù)為6，位移預(yù)測數(shù)據(jù)點數(shù)為1，再運用交叉驗證法對模型參數(shù)尋優(yōu)，C=1.05，σ=1.93，最后訓(xùn)練模型并輸出預(yù)測值。預(yù)測結(jié)果如圖7所示，通過與監(jiān)測結(jié)果比較可知，RSME 為0.018 mm，MAE 為0.017 mm，R2為0.99。

圖7 隨機(jī)性位移預(yù)測結(jié)果

2.3.4 累計位移預(yù)測

根據(jù)上述分解流程可知，累計位移為趨勢性位移、周期性位移和隨機(jī)性位移之和。因此，將趨勢性位移預(yù)測值、周期性位移預(yù)測值和隨機(jī)性位移預(yù)測值累加求和，即可得到累計位移預(yù)測值，結(jié)果如圖8 所示。通過與實測值對比可知，位移預(yù)測模型的RSME 為0.045 mm，MAE 為0.039 mm，R2為0.97。

圖8 累計位移預(yù)測結(jié)果

2.4 模型比較

本節(jié)同時也采用BP 神經(jīng)網(wǎng)絡(luò)和傳統(tǒng)最小二乘支持向量機(jī)回歸模型（LSSVM）對累計位移預(yù)測，以分析SSSC-EMD-LSSVM 模型的優(yōu)越性。為避免不同尋優(yōu)方法的影響，這2 種模型的超參數(shù)也采用交叉驗證法確定。

表1 不同模型預(yù)測精度比較結(jié)果

3 結(jié)論

（1）針對非線性位移時間序列難以準(zhǔn)確預(yù)測的難題，本文結(jié)合優(yōu)化經(jīng)驗?zāi)B(tài)分解和最小二乘支持向量機(jī)，將時間序列分層分解為多個分量，分別預(yù)測后再疊加重構(gòu)，即可得到累計位移預(yù)測值。

（2）以山口巖水庫為例開展位移時間序列預(yù)測。與傳統(tǒng)BP 神經(jīng)網(wǎng)絡(luò)和LSSVM 模型相比，本文提出的SSSC-EMD-LSSVM 模型預(yù)測精度更優(yōu)，RMSE 為0.045 mm，MAE 為0.039 mm，R2為0.97，說明這是一種性能較好的非線性位移時間序列預(yù)測模型。