基于蟻群算法優(yōu)化極限學(xué)習(xí)機模型的滑坡位移預(yù)測

曹 博， 汪 帥， 宋丹青， 杜 涵， 劉光偉， 周志偉

(1.遼寧工程技術(shù)大學(xué) 礦業(yè)學(xué)院， 遼寧 阜新 123000； 2.清華大學(xué) 水沙科學(xué)與水利水電工程國家重點試驗室，北京 100084； 3.神華寶日希勒能源有限公司露天煤礦， 內(nèi)蒙古 呼倫貝爾 021000)

1 研究背景

我國是世界上滑坡災(zāi)害發(fā)育最為嚴(yán)重的國家之一[1-2]，它的頻繁發(fā)生和快速運動，每年都會造成大量的人員傷亡和財產(chǎn)損失[3]。2018年的《全國地質(zhì)災(zāi)害通報》表明，2017年全年我國共發(fā)生地質(zhì)災(zāi)害2 966起，其中滑坡災(zāi)害1 631起，占比為55.0%，因此，如何有效防止滑坡災(zāi)害成為了該研究領(lǐng)域的熱點問題。

近年來，滑坡位移預(yù)測已逐漸成為滑坡早期預(yù)警的有效手段，它允許研究人員通過室內(nèi)計算來分析滑坡歷史變形的主要特征，并對未來一段時間內(nèi)滑坡的位移趨勢進行預(yù)測[4]。由于該方法所需成本較低并且精度較高，因此得到了快速發(fā)展。目前常用的滑坡位移定量預(yù)測方法主要包括數(shù)理統(tǒng)計模型和非線性理論[5]，2009年Li等[6]提出在滑坡預(yù)測過程中使用尖點突變模型，王珣等[7]在2017年建立切線角模型進行了臨滑預(yù)報判據(jù)研究，但是這些方法通常需要對滑坡的演化狀態(tài)進行較為深入準(zhǔn)確的分析，不確定性較大。與這些模型相比，機器學(xué)習(xí)模型不需要分析復(fù)雜的滑坡機理且對于以滑坡為代表的非線性動態(tài)系統(tǒng)具有較好的擬合效果，因此應(yīng)用更加廣泛，最具代表性的模型包括2019年柳青等[8]進行滑坡位移預(yù)測預(yù)報所采用的神經(jīng)網(wǎng)絡(luò)模型、Zhou等[9]和李仕波等[10]先后于2016及2019年預(yù)測位移時采用的支持向量機模型、2015年周超等[11]預(yù)測滑坡位移時引入的混沌理論。極限學(xué)習(xí)機模型(extreme learning machine, ELM)作為近年來應(yīng)用廣泛的智能算法，其在滑坡位移預(yù)測方面的適用性得到了不少研究的驗證，如Cao等[12]利用ELM模型對白家包滑坡的位移進行了預(yù)測，同時與支持向量機模型進行了預(yù)測效果比較；Huang等[13]對滑坡位移利用混沌理論和小波變換進行分解后，再基于ELM模型進行了非平穩(wěn)時間序列的預(yù)測。在ELM模型的計算過程中，由于有個別參數(shù)是隨機選取的初始值，因此為加快收斂速度并提高準(zhǔn)確性，需要加入一定的算法進行選參優(yōu)化，如粒子群算法[10]、遺傳算法[5]、蟻群算法[14]等常見方法。然而，對于多個算法優(yōu)化效果的對比，前人的研究較為欠缺。

綜上所述，為驗證ELM模型在我國三峽庫區(qū)滑坡位移預(yù)測中的適用性，本文選擇湖北秭歸八字門滑坡為研究對象，在深入分析其變形歷史和影響因素后，完成了基于蟻群算法優(yōu)化ELM模型的位移預(yù)測過程，并通過對比驗證了不同優(yōu)化算法的優(yōu)化效果，以期為后續(xù)研究提供科學(xué)依據(jù)。

2 蟻群算法優(yōu)化極限學(xué)習(xí)機模型

2.1 基于移動平均法的時間序列分解

基于時間序列分析原理，對于通過監(jiān)測或遙感等手段得到的滑坡位移數(shù)據(jù)，可以將其看成是一種非線性時間序列，并且有如下關(guān)系[15]：

yt=St+Ct+εt

(1)

式中：yt為滑坡累計位移，mm；St、Ct和εt分別為滑體的趨勢項、周期項和隨機項位移，mm。其中，趨勢項位移代表了滑體在長時間段內(nèi)的變形趨勢，通常為單調(diào)遞增；周期項位移表示滑坡位移在周期性因素影響下的變形；隨機項位移為滑坡在隨機因素影響下產(chǎn)生的位移，由于影響因素較難量化，因此以往研究通常將該項忽略。

對于滑坡位移預(yù)測的研究，通常要利用一定的方法將滑坡的累計位移分解成不同的分量，通過對子序列進行分析預(yù)測，疊加所有子序列的預(yù)測值來得到總的滑坡位移預(yù)測值。分解方法主要包括兩種，即數(shù)理統(tǒng)計法和智能方法，后者又包括小波分解[11]、經(jīng)驗?zāi)B(tài)分解[16]、集合經(jīng)驗?zāi)B(tài)分解[17]等方法，但是這類方法總體來說較為復(fù)雜。針對滑坡位移時序的特點，本文選擇移動平均法進行總位移的分解，其計算公式如下：

Mi={M1,M2,…,Mt,…,Mn}

(2)

(3)

2.2 極限學(xué)習(xí)機

2006年Huang等[18]首先提出了極限學(xué)習(xí)機，其本質(zhì)是一種單層神經(jīng)網(wǎng)絡(luò)，但相比較而言，ELM模型的結(jié)構(gòu)更為簡單，且運算速度更快[19]，并且其中的層間連接權(quán)值及隱含層偏差不需要事先設(shè)置。在ELM算法中，隱含層神經(jīng)元閾值及輸入層與隱含層的連接權(quán)值由函數(shù)生成，訓(xùn)練過程中無需人為干預(yù)，只需進行隱含層神經(jīng)元個數(shù)設(shè)置即可獲得該條件下的最優(yōu)解。

若設(shè)模型中輸入層與隱含層之間連接權(quán)值為w，同時隱含層與輸出層間的連接權(quán)值為β，給定訓(xùn)練樣本總數(shù)為Q，ELM神經(jīng)網(wǎng)絡(luò)的輸出為T，激活函數(shù)為g(x)時，隱含層輸出矩陣H可由下式求解：

Hβ=T′

(4)

(5)

式中：H+為輸出矩陣H的廣義逆矩陣。

極限學(xué)習(xí)機模型在理論上可以獲得最小的訓(xùn)練誤差，利用ELM模型進行滑坡位移預(yù)測模型構(gòu)建的過程中，所需考慮的參數(shù)僅為隱含層神經(jīng)元個數(shù)，利用ELM模型進行位移預(yù)測主要包含以下步驟：

(1)對隱含層相關(guān)參數(shù)隨機進行初始化，確定隱含層神經(jīng)元個數(shù)，并隨機生成輸入層與隱含層間的連接權(quán)值及隱含層神經(jīng)元偏差的初始值；

(2)合理選擇一個無限可微函數(shù)作為隱含層激活函數(shù)，以便求解隱含層輸出矩陣H；

2.3 蟻群算法優(yōu)化極限學(xué)習(xí)機

Parpinelli等[20]于2002年首次提出了蟻群算法，蟻群算法類似于自然界中螞蟻的覓食過程，螞蟻群體會通過不同路徑搜尋食物，對于正確的路徑，螞蟻的數(shù)量和密度會越來越大，而對于其他路徑則相反。最后螞蟻能夠搜尋到一條抵達(dá)食物源頭的最佳路徑，即問題的最優(yōu)解[20-21]。蟻群算法最初用于旅行商難題上，隨著科研人員的不斷探索，應(yīng)用范疇幾乎涉及各個領(lǐng)域。為更好更快地找到最優(yōu)的連接權(quán)值與神經(jīng)元閾值，本文擬采用蟻群優(yōu)化算法對極限學(xué)習(xí)機模型進行改進，以達(dá)到更優(yōu)的計算效果。

使用蟻群算法優(yōu)化極限學(xué)習(xí)機模型的流程如下所示：

(1)初始化參數(shù)，計算螞蟻從一個點轉(zhuǎn)移到另一個點的概率，概率計算公式如下：

(6)

(2)通過不停迭代，尋求一個最佳的覓食路徑，確定蟻群算法全局最優(yōu)路徑，即極限學(xué)習(xí)機的權(quán)值和閾值；

(3)基于蟻群算法最優(yōu)路徑的計算，合理選擇隱含層激活函數(shù)，求解隱含層輸出矩陣H；

(4)對輸出層權(quán)值矩陣進行求解，以此為基礎(chǔ)進行滑坡位移預(yù)測。

3 研究實例

3.1 八字門滑坡概況

八字門滑坡位于湖北省宜昌市秭歸縣歸州鎮(zhèn)香溪河右岸，距離三峽大壩38 km?；赂叱谭植荚?40～270 m之間，滑體平均長度和平均寬度分別為350和400 m，厚度約為30 m，體積達(dá)400×104m3，屬于大型堆積層滑坡。根據(jù)野外鉆孔揭露的情況顯示，滑坡包括兩層滑帶，上部為次級滑面，厚度約為3 m，下部的主滑帶連續(xù)分布，中部的厚度最大?；禄鶐r為侏羅系香溪組(J1x)的砂巖和泥質(zhì)粉砂巖，呈強分化特性，局部夾雜有炭質(zhì)泥巖。八字門滑坡平面及剖面圖如圖1所示。

3.2 滑坡變形特征及影響因素分析

自從八字門滑坡發(fā)生較為明顯的變形開始，地質(zhì)災(zāi)害管理部門就對其變形進行了持續(xù)的監(jiān)測，主要包括深部位移和地表GPS位移。在滑坡的監(jiān)測剖面Ⅰ-Ⅰ′上(圖1)，安裝了3個地表監(jiān)測點(ZG110、ZG111、ZG112)和1個深部位移監(jiān)測點(QTZK3)，監(jiān)測點多年監(jiān)測結(jié)果以及相應(yīng)的庫水位和月降雨量變化見圖2，2005年3月與2004年9月相比滑坡深部監(jiān)測點的相對位移如3所示。

3個地表監(jiān)測點的高程分別約為165、191和215 m，可以代表滑坡前、中、后部的變形特點和趨勢。圖2中滑坡累計位移曲線顯示，更靠近滑坡后緣的ZG111監(jiān)測點的位移量最大，總位移超過了1 700 mm，因為該點所在位置是滑坡的陡坎處，在每年的雨季和庫水位下降期(5-9月)，該點都會發(fā)生較為明顯的階躍式變形增長，如從2009年的5月開始，在5個月的時間內(nèi)，該點的累計位移增長超過了200 mm，增長率達(dá)到了24.9%。但是位于滑坡中部的ZG110監(jiān)測點和靠近滑坡前緣的ZG112監(jiān)測點的總位移卻較小，分別約為1 200 mm和360 mm，且該兩個監(jiān)測點的位移增長趨勢較為平穩(wěn)。上述結(jié)果說明八字門滑坡的運動模式為“前進式”，即最開始由滑坡上部發(fā)生蠕變并向下擠壓，隨著位移的逐漸增加，變形開始向滑坡前緣發(fā)展，最終形成整體破壞。

圖1 香溪河八字門滑坡平面及剖面圖

圖2 2004-2012年八字門滑坡各監(jiān)測點地表位移監(jiān)測曲線以及相應(yīng)庫水位和月降雨量變化

圖3 八字門滑坡2005年3月相對于2004年9月深部監(jiān)測點的位移

選取變形量最大的ZG111監(jiān)測點對滑坡影響因素和變形間的響應(yīng)關(guān)系進行分析，在三峽水庫蓄水前，滑坡前緣被庫水淹沒，加之地下水上升帶來的侵蝕作用，滑動面前部的抗剪強度逐漸降低，引起了滑坡前部的緩慢蠕變，并在前緣逐漸形成了少量張拉裂縫，因此在這段時間內(nèi)滑坡存在漸進式變形的特征。2003年三峽水庫開始蓄水，在八字門滑坡上出現(xiàn)了較多的變形跡象；同年8月，滑坡體上道路兩側(cè)出現(xiàn)了兩個較大的裂縫，長度分別達(dá)到了30和45 m；2004年的5月至8月，庫水位短暫下降，在道路上又發(fā)現(xiàn)了多處裂縫，如圖4所示。

圖4 八字門滑坡體上公路變形實拍圖

由于水庫的周期蓄水，滑坡的累計位移也隨之發(fā)生變化。在2007年2-4月，水庫水位由156 m下降至152 m，但累計位移幾乎不變，但是在該年的5-9月，由于總降水量達(dá)到近600 mm，該點的累計位移迅速增加，變化量超過了120 mm，同時伴隨有多處地表裂縫的出現(xiàn)。這是因為短時間內(nèi)的強降水導(dǎo)致地下水位的抬升，而該滑坡主要由堆積物構(gòu)成，高孔隙率和良好的滲透性使得地表水的入滲更加容易，導(dǎo)致了異常的孔隙水壓力和水動力壓力誘發(fā)的滑坡變形。綜上所述，滑坡的變形主要受到了庫水位波動和降雨入滲的影響。

3.3 滑坡位移分解及預(yù)測

3.3.1 基于移動平均法的累計位移分解 由于ZG111測點處的累計位移最大，而且可以反映滑坡的階躍式變形特點，因此選取該點進行位移預(yù)測。總計獲取102組數(shù)據(jù)，將前90組數(shù)據(jù)作為訓(xùn)練數(shù)據(jù)集，后12組數(shù)據(jù)作為測試數(shù)據(jù)集，首先根據(jù)公式(1)、(2)將2005年7月至2013年1月的全部數(shù)據(jù)分為趨勢項位移和周期項位移，如圖5所示。由圖5可看出，提取出的趨勢項位移較為平滑，基本不存在突變現(xiàn)象，說明移動平均法可以很好地消除滑坡階躍式曲線的影響，從而反映出真正的滑坡長期變形特征。另一方面，滑坡的周期項位移也具有較明顯的周期特征，而且初步估計該曲線的周期近似為12個月，也符合實際情況。

圖5 2005年7月至2013年1月八字門滑坡ZG111測點累計位移分解結(jié)果

3.3.2 滑坡位移分量預(yù)測 趨勢項位移反映滑坡長時間以來的變形趨勢，在目前的研究中該位移多使用多項式函數(shù)進行預(yù)測。對于上節(jié)中分解得到的趨勢項數(shù)據(jù)，利用最小二乘法進行預(yù)測。選取圖5中前79組數(shù)據(jù)先進行擬合，采用二次多項式作為擬合函數(shù)，擬合之后再將公式外推獲得后面測試數(shù)據(jù)集的預(yù)測結(jié)果，利用本案例的訓(xùn)練數(shù)據(jù)最后獲得的預(yù)測公式如下：

y=0.0955t2+6.0632t+324.59 (R2=0.998)

(6)

趨勢項位移預(yù)測結(jié)果公式(6)的擬合優(yōu)度達(dá)到了0.998，表明該公式的預(yù)測效果很好，能夠有效反映滑坡變形隨時間變化的整體趨勢。

對于周期項位移數(shù)據(jù)，它通常與滑坡的周期性影響因素有關(guān)，本文選取降雨和庫水位變化作為滑坡周期性的影響因素，再參考類似研究，還需要考慮滑坡變形對兩種誘發(fā)因素的響應(yīng)滯后時間。因此，最終選擇上月降雨量、當(dāng)月降雨量、上月庫水位平均高程和當(dāng)月庫水位平均高程作為影響因素。這些影響因素被作為輸入數(shù)據(jù)集，而當(dāng)月周期項位移被作為輸出數(shù)據(jù)集，然后基于自編代碼構(gòu)建蟻群算法優(yōu)化的極限學(xué)習(xí)機模型，進行周期項位移的預(yù)測。以2012年2月至2013年1月為預(yù)測期，得出的滑坡趨勢項位移和周期項位移預(yù)測結(jié)果與實際值比較如圖6所示。由圖6可看出，周期項位移預(yù)測結(jié)果與實際較為吻合，表明預(yù)測效果也較為理想。

圖6 2012年2月至2013年1月滑坡趨勢項和周期項位移預(yù)測結(jié)果與實際值比較

4 討 論

將3.3節(jié)得到的趨勢項和周期項兩項分量相加，得到預(yù)測期蟻群算法優(yōu)化的極限學(xué)習(xí)機模型預(yù)測的最終累計位移結(jié)果，如圖7所示，同時，將單一的極限學(xué)習(xí)機模型預(yù)測結(jié)果以及采用兩種常用的優(yōu)化算法(粒子群算法和遺傳算法)對極限學(xué)習(xí)機模型進行優(yōu)化后的預(yù)測結(jié)果一并繪入圖7，并與位移實際值進行比較，以分析本文所采用的蟻群算法優(yōu)化的極限學(xué)習(xí)機模型的預(yù)測效果。上述各模型預(yù)測效果主要參數(shù)對比見表1。

圖7 各模型預(yù)測結(jié)果

分析圖7和表1可知：(1)在所有4種模型當(dāng)中，單一極限學(xué)習(xí)機模型的計算用時最短，但是其計算精度最差，平均誤差達(dá)到了23.5 mm，擬合優(yōu)度為0.973，此2項指標(biāo)均低于其他3種預(yù)測模型。然而，考慮到滑坡的累計位移已經(jīng)超過了1 000 mm，因此單一極限學(xué)習(xí)機模型的預(yù)測精度也是可以接受的，這說明該模型在階躍式滑坡位移曲線方面的預(yù)測是較為可靠的，也驗證了極限學(xué)習(xí)機模型對于本文應(yīng)用實例的適用性。(2)對于優(yōu)化后的極限學(xué)習(xí)機模型，蟻群算法優(yōu)化的極限學(xué)習(xí)機計算用時為98 s，位移預(yù)測結(jié)果平均誤差為10.1 mm，預(yù)測結(jié)果擬合優(yōu)度為0.998，各指標(biāo)均優(yōu)于其他兩種算法優(yōu)化的極限學(xué)習(xí)機模型。因此，對于本文所示的案例而言，蟻群算法具有最佳的參數(shù)選取與優(yōu)化效果，蟻群算法優(yōu)化的極限學(xué)習(xí)機模型是進行三峽庫區(qū)八字門滑坡累計位移預(yù)測的最優(yōu)模型，可在類似滑坡體的位移預(yù)測中進一步研究與推廣。

表1 各模型預(yù)測效果主要指標(biāo)對比

5 結(jié)論與展望

本文探討了蟻群算法優(yōu)化的極限學(xué)習(xí)機模型在三峽庫區(qū)堆積層滑坡位移預(yù)測中的應(yīng)用，并以秭歸縣八字門滑坡為例進行了預(yù)測效果分析，主要得到以下幾點結(jié)論：

(1)極限學(xué)習(xí)機模型能夠較好地進行三峽庫區(qū)八字門滑坡的累計位移預(yù)測，單一極限學(xué)習(xí)機模型預(yù)測的平均誤差為23.5 mm，擬合優(yōu)度可達(dá)0.973，是一種可靠的位移預(yù)測方法。

(2)蟻群算法優(yōu)化的極限學(xué)習(xí)機模型比單一的極限學(xué)習(xí)機模型的預(yù)測精度更高，說明該算法的優(yōu)化性能較好；同時與其他2種傳統(tǒng)的優(yōu)化算法(粒子群算法、遺傳算法)相比，蟻群算法在計算用時和計算精度上更優(yōu)，其預(yù)測結(jié)果與實際結(jié)果間的平均誤差為10.1 mm，擬合優(yōu)度達(dá)0.998。該模型可在類似滑坡的位移預(yù)測研究中進行推廣。

(3)在蟻群算法優(yōu)化參數(shù)選取的過程中，部分算法參數(shù)仍需進行進一步討論，算法參數(shù)值與計算結(jié)果的用時和精度間的敏感性分析將會是未來的一個研究方向。