時(shí)序InSAR滑坡形變監(jiān)測(cè)與預(yù)測(cè)的N-BEATS深度學(xué)習(xí)法
——以新鋪滑坡為例

郭澳慶，胡 俊，鄭萬基，桂 容，杜志貴，朱 武，賀樂和

1. 中南大學(xué)地球科學(xué)與信息物理學(xué)院，湖南 長(zhǎng)沙 410083; 2. 長(zhǎng)沙天儀空間科技研究院有限公司，湖南 長(zhǎng)沙 410205; 3. 長(zhǎng)安大學(xué)地質(zhì)工程與測(cè)繪學(xué)院，陜西 西安 710054

滑坡是指在斜坡上的巖體土體，當(dāng)受流水沖刷、降雨沖刷、地下水運(yùn)動(dòng)及地震等因素的影響時(shí)，在重力作用下失穩(wěn)而沿斜坡向下滑動(dòng)的自然現(xiàn)象[1]。山體滑坡通常發(fā)生突然，破壞力巨大，經(jīng)常會(huì)造成重大的生命安全事故和財(cái)產(chǎn)損失。因此，滑坡形變的全面監(jiān)測(cè)和中短期實(shí)時(shí)預(yù)測(cè)技術(shù)的進(jìn)步對(duì)防災(zāi)救災(zāi)來說更具有突出的現(xiàn)實(shí)意義。

滑坡形變監(jiān)測(cè)和預(yù)測(cè)方面的研究在目前來說仍是一大熱點(diǎn)。合成孔徑雷達(dá)干涉測(cè)量技術(shù)(interferometric synthetic aperture radar,InSAR)[2-5]以其非接觸式的測(cè)量功能獲得衛(wèi)星目所能及之處的面狀形變數(shù)據(jù)，相較于利用全球?qū)Ш叫l(wèi)星系統(tǒng)數(shù)據(jù)(Global Navigation Satellite System,GNSS)和傳統(tǒng)地面監(jiān)測(cè)數(shù)據(jù)進(jìn)行滑坡預(yù)測(cè)而言，InSAR技術(shù)能夠提供無接觸式低成本的面域測(cè)量，既降低了GNSS接收站和傳統(tǒng)地面監(jiān)測(cè)儀器的安裝維護(hù)成本，又能夠獲得滑坡形變區(qū)更加全面的面域形變狀態(tài)，成為了滑坡形變的主流監(jiān)測(cè)手段之一[6-9]。

由于滑坡形變序列受到機(jī)理復(fù)雜和種類繁多的變量因子的影響[10]，其波動(dòng)呈現(xiàn)出非線性和不確定性而導(dǎo)致滑坡預(yù)測(cè)問題較難解決，在基于InSAR技術(shù)的滑坡形變預(yù)測(cè)方面的研究更是屈指可數(shù)。因此，尋找到更加適合于InSAR形變數(shù)據(jù)分析且能夠準(zhǔn)確有效完成滑坡形變預(yù)測(cè)任務(wù)的方法和技術(shù)是解決問題的方向。

目前，針對(duì)InSAR滑坡形變數(shù)據(jù)分析的方法可以被分為兩類：模型驅(qū)動(dòng)型和數(shù)據(jù)驅(qū)動(dòng)型[11]。模型驅(qū)動(dòng)型的代表有概率積分法、卡爾曼濾波、灰色系統(tǒng)理論模型、逆速法線性逼近、非線性最小二乘及有限元方法模擬等，依靠形變數(shù)據(jù)及豐富的巖土信息能夠獲得高精度且可靠的地面形變預(yù)測(cè)結(jié)果。但實(shí)際上，這類數(shù)值模擬方法受測(cè)量誤差、環(huán)境噪聲及建模假設(shè)等因素的干擾效果明顯，表現(xiàn)為模型抗差性能低，遷移應(yīng)用能力弱，各類巖土信息的獲取也消耗巨大的人力物力資源，大多數(shù)中小型工程難以承擔(dān)設(shè)備維護(hù)及災(zāi)害預(yù)警的成本。數(shù)據(jù)驅(qū)動(dòng)是當(dāng)前很熱門的數(shù)據(jù)分析方法，它利用數(shù)據(jù)本身的信息和特定任務(wù)的實(shí)際需要能夠創(chuàng)造出新的計(jì)算模型；在各種實(shí)際情況下使模型貼合特定的數(shù)據(jù)，用客觀和清晰的數(shù)據(jù)將分析過程化繁為簡(jiǎn)，因而無須受限于現(xiàn)有模型；相較于模型驅(qū)動(dòng)更具有從統(tǒng)計(jì)觀點(diǎn)來解決實(shí)際問題的特色。數(shù)據(jù)驅(qū)動(dòng)型分析方法隨著機(jī)器學(xué)習(xí)的出現(xiàn)而逐漸被研究者們接受和青睞，深度學(xué)習(xí)則是一種既能智能模擬人類大腦又具有強(qiáng)大計(jì)算和數(shù)值分析能力的工具，例如長(zhǎng)短期記憶(long short-term memory，LSTM)網(wǎng)絡(luò)模型[12]在大范圍地表沉降形變預(yù)測(cè)的應(yīng)用中大放異彩。神經(jīng)網(wǎng)絡(luò)底層擴(kuò)展分析(neural basis expansion analysis，N-BEATS)網(wǎng)絡(luò)模型[13]是一種深度神經(jīng)網(wǎng)絡(luò)，是數(shù)據(jù)驅(qū)動(dòng)模型進(jìn)行時(shí)間序列數(shù)據(jù)分析的代表。它在不同領(lǐng)域的大型時(shí)間序列數(shù)據(jù)分析中都展現(xiàn)出了先進(jìn)的預(yù)測(cè)能力，如M3、M4和TOURISM競(jìng)賽數(shù)據(jù)集。常用時(shí)序形變預(yù)測(cè)方法見表1。

表1 地質(zhì)災(zāi)害時(shí)序數(shù)據(jù)預(yù)測(cè)方法優(yōu)劣一覽

本文利用時(shí)序InSAR獲取形變監(jiān)測(cè)結(jié)果，結(jié)合N-BEATS神經(jīng)網(wǎng)絡(luò)開展滑坡形變預(yù)測(cè)工作，發(fā)揮N-BEATS模型強(qiáng)大的時(shí)間序列數(shù)據(jù)自動(dòng)分析，以及高精度、高可靠性和具有一定抗差能力的預(yù)測(cè)優(yōu)勢(shì)，將研究目標(biāo)定位于中短期滑坡小數(shù)據(jù)集，應(yīng)用于三峽庫(kù)區(qū)的Ⅰ級(jí)滑坡—新鋪滑坡，把這種有效的深度學(xué)習(xí)時(shí)間序列預(yù)測(cè)模型應(yīng)用于解決基于時(shí)序InSAR形變觀測(cè)的中短期滑坡預(yù)測(cè)難題。

1 基于N-BEATS網(wǎng)絡(luò)和時(shí)序InSAR的滑坡形變預(yù)測(cè)方法

常規(guī)測(cè)量方法在理論建設(shè)階段都有著較為完美的效果，但在實(shí)際的工程實(shí)踐中卻難以克服一些特殊的環(huán)境因素，例如滑坡發(fā)生區(qū)域的茂密植被會(huì)遮擋導(dǎo)航定位信號(hào)，導(dǎo)致相位信息斷鏈?zhǔn)笹NSS技術(shù)時(shí)而失效。地下地面常規(guī)測(cè)量方法一次建設(shè)投入成本巨大，二次建設(shè)及維護(hù)成本在滑坡發(fā)生之時(shí)被毀壞后更是不知凡幾。InSAR技術(shù)能夠高精度、少干擾獲得面狀形變數(shù)據(jù)，更適合大范圍滑坡區(qū)域的全面監(jiān)測(cè)；而具有強(qiáng)大學(xué)習(xí)能力、有效探索時(shí)序變化規(guī)律的深度學(xué)習(xí)模型N-BEATS更適合復(fù)雜非線性滑坡形變預(yù)測(cè)。

本文基于深度學(xué)習(xí)和時(shí)序InSAR的滑坡形變預(yù)測(cè)任務(wù)主要分為3個(gè)步驟：滑坡形變監(jiān)測(cè)獲得歷史數(shù)據(jù)、組織數(shù)據(jù)結(jié)構(gòu)和時(shí)間序列數(shù)據(jù)分析向前預(yù)測(cè)。首先，在SAR監(jiān)測(cè)獲得時(shí)序形變影像后經(jīng)過時(shí)序TCP-InSAR技術(shù)[28]進(jìn)行處理，該方法能夠克服在InSAR數(shù)據(jù)分析過程中，由于相位解纏處理不當(dāng)而帶來的誤差影響的問題，得到滑坡區(qū)的形變速率。然后，對(duì)時(shí)序滑坡形變影像數(shù)據(jù)進(jìn)行選擇和重新組織，找出滑坡危險(xiǎn)區(qū)特征點(diǎn)，并將特征點(diǎn)數(shù)據(jù)重新組織成為一維時(shí)間序列，并對(duì)一維時(shí)間序列數(shù)據(jù)進(jìn)行時(shí)序分解，分析其變化特征，確定預(yù)測(cè)任務(wù)數(shù)據(jù)集的輸入輸出窗口設(shè)置策略，以防損失數(shù)據(jù)特征。最后，將時(shí)間序列按照比例劃分為訓(xùn)練集和測(cè)試集，訓(xùn)練N-BEATS網(wǎng)絡(luò)模型并測(cè)試模型訓(xùn)練結(jié)果，完成時(shí)序預(yù)測(cè)。本文的主要工作是預(yù)測(cè)任務(wù)數(shù)據(jù)集設(shè)置策略的探究及N-BEATS時(shí)序預(yù)測(cè)模型建立、訓(xùn)練和預(yù)測(cè)兩個(gè)部分。

1.1 預(yù)測(cè)任務(wù)數(shù)據(jù)集設(shè)置策略

目前，基于InSAR技術(shù)的滑坡形變預(yù)測(cè)還停留在單點(diǎn)時(shí)序預(yù)測(cè)階段，對(duì)于深度學(xué)習(xí)來說，需要人為進(jìn)行優(yōu)化控制的主要是超參數(shù)的選擇，除此之外還有數(shù)據(jù)組織形式的設(shè)置策略。由于硬件設(shè)備的計(jì)算性能有限，因此一次訓(xùn)練的輸入數(shù)據(jù)不宜過大，也就是無法將全部歷史數(shù)據(jù)作為單次輸入進(jìn)行操作，那么這就需要人工根據(jù)經(jīng)驗(yàn)選取部分具有充分代表性的數(shù)據(jù)作為輸入的訓(xùn)練數(shù)據(jù)。研究者們通過豐富的試驗(yàn)和經(jīng)驗(yàn)發(fā)現(xiàn)，不同的數(shù)據(jù)組織形式對(duì)數(shù)據(jù)預(yù)測(cè)效果和訓(xùn)練能力都有很大的影響。

在輸入窗口大小的設(shè)置中，針對(duì)具有周期性特征的時(shí)間序列，其輸入窗口的長(zhǎng)度需要達(dá)到能夠體現(xiàn)真實(shí)周期的數(shù)據(jù)長(zhǎng)度，以免丟失周期信息；針對(duì)具有突發(fā)性變化的時(shí)間序列，其輸入窗口的長(zhǎng)度需要充分包括突發(fā)性狀況的整個(gè)時(shí)序變化趨勢(shì)；針對(duì)包括其他非線性特征的時(shí)間序列，也需要獲取相應(yīng)的具有充分代表性的數(shù)據(jù)長(zhǎng)度。在輸出窗口大小的設(shè)置中，需要在滿足深度學(xué)習(xí)預(yù)測(cè)結(jié)果正確的前提下同時(shí)兼顧真實(shí)預(yù)測(cè)時(shí)間頻率及時(shí)間延展需求。

對(duì)上述數(shù)據(jù)進(jìn)行規(guī)律性分析的方法在時(shí)間序列分析領(lǐng)域有很多，本文主要采取的策略是對(duì)數(shù)據(jù)進(jìn)行簡(jiǎn)單時(shí)序分解，基于Loess的季節(jié)趨勢(shì)分解方法(seasonal-trend decomposition procedure based on Loess,STL)是一種能夠?qū)⒋嬖诳山忉屝蕴卣鞯臅r(shí)間序列分解成為趨勢(shì)因子、季節(jié)因子及隨機(jī)因子，適用于本文所針對(duì)的大型緩慢滑坡形變應(yīng)用。

1.2 N-BEATS網(wǎng)絡(luò)模型架構(gòu)

N-BEATS網(wǎng)絡(luò)是由若干感知機(jī)單元組成，依據(jù)一定策略搭建成的深度神經(jīng)網(wǎng)絡(luò)。該網(wǎng)絡(luò)結(jié)構(gòu)具有兩個(gè)作用：①模擬，對(duì)原始時(shí)間序列數(shù)據(jù)的特征進(jìn)行提取并模擬；②預(yù)測(cè)，從時(shí)間維度延續(xù)模擬出的序列數(shù)據(jù)特征，從而完成預(yù)測(cè)。

N-BEATS網(wǎng)絡(luò)模型的搭建策略是從時(shí)間序列分解中獲得的靈感，因此該模型的優(yōu)勢(shì)不僅在于能夠獲得準(zhǔn)確性高的預(yù)測(cè)結(jié)果，其預(yù)測(cè)結(jié)果也具有類似于STL方法的對(duì)趨勢(shì)性、季節(jié)性等時(shí)序變化的可解釋能力。其框架結(jié)構(gòu)如圖1所示，網(wǎng)絡(luò)一次訓(xùn)練的輸入即為設(shè)定窗口大小的歷史數(shù)據(jù)序列，使序列數(shù)據(jù)依次經(jīng)過M個(gè)堆棧來分別識(shí)別數(shù)據(jù)的不同特征，每個(gè)堆棧中存在若干塊區(qū)，塊區(qū)之間通過殘差模塊進(jìn)行銜接，每個(gè)塊區(qū)中又含有幾個(gè)計(jì)算神經(jīng)元。

1.2.1 N-BEATS網(wǎng)絡(luò)架構(gòu)的工作思路

(1) 區(qū)塊層學(xué)習(xí)時(shí)間序列演化模式。首先，使原始時(shí)間序列數(shù)據(jù)經(jīng)過4層全連接層；然后，經(jīng)過線性投影層后進(jìn)行線性變換分別生成向前預(yù)測(cè)系數(shù)和向后預(yù)測(cè)系數(shù)；最后，使向前預(yù)測(cè)系數(shù)和向后預(yù)測(cè)系數(shù)分別通過偏置層輸出向前預(yù)測(cè)和向后預(yù)測(cè)結(jié)果，如圖1(a)所示。向前預(yù)測(cè)結(jié)果是沿時(shí)間序列發(fā)展方向的預(yù)測(cè)未來序列走向的結(jié)果。向后預(yù)測(cè)結(jié)果是恢復(fù)出的歷史時(shí)間序列，用于信號(hào)分離。因此，塊區(qū)中輸入為原始序列數(shù)據(jù)；輸出為向前預(yù)測(cè)和向后預(yù)測(cè)結(jié)果。

(2) 堆棧層分解剩余時(shí)間序列特征。每一個(gè)塊區(qū)的輸出都是其的向前預(yù)測(cè)和向后預(yù)測(cè)結(jié)果，如圖1(b)所示。向前預(yù)測(cè)結(jié)果被整合起來成為該堆棧的輸出結(jié)果，向后預(yù)測(cè)結(jié)果通過殘差的形式存在于塊區(qū)之間，并將塊區(qū)銜接起來，更有利于梯度的反向傳播也使得下一個(gè)塊區(qū)的預(yù)測(cè)工作更加簡(jiǎn)便。因此，本堆棧中的輸入為各塊區(qū)向前預(yù)測(cè)和向后預(yù)測(cè)結(jié)果；輸出為這部分中經(jīng)過整合得到的向前預(yù)測(cè)結(jié)果。

(3) 整合時(shí)間序列特征預(yù)測(cè)分量。各部分堆棧的輸出最終被整合了起來，疊加各個(gè)堆棧獲得的特征分量就能夠得到總體預(yù)測(cè)結(jié)果，即預(yù)測(cè)數(shù)據(jù)的序列，如圖1(c)所示。

圖1 N-BEATS網(wǎng)絡(luò)模型架構(gòu)Fig.1 N-BEATS network model architecture

1.2.2 N-BEATS網(wǎng)絡(luò)架構(gòu)優(yōu)勢(shì)

該網(wǎng)絡(luò)架構(gòu)能夠?qū)崿F(xiàn)對(duì)趨勢(shì)性、季節(jié)性等時(shí)序變化的解釋，主要核心在于對(duì)不同堆棧中的神經(jīng)單元建議不同的線性層和偏置層，架構(gòu)中主要的優(yōu)勢(shì)即可解釋性，它體現(xiàn)在趨勢(shì)項(xiàng)和周期項(xiàng)兩堆棧中。

(1) 趨勢(shì)模型。趨勢(shì)模型的特征通常是呈線性單調(diào)遞增或單調(diào)遞減，并且其變化應(yīng)是較為緩慢的。因此在趨勢(shì)項(xiàng)堆棧中(c)部分的線性層和偏置層被設(shè)置為

y=T·x

(1)

式中，T=[1,t1,t2,…,tp]，p為一個(gè)小值，用于控制線性單調(diào)遞增或遞減行為的速度，使趨勢(shì)緩慢變化；t為時(shí)間向量；x為該塊區(qū)中線性層和偏置層的輸入；y為該塊區(qū)中線性層和偏置層的輸出。

(2) 周期模型。周期模型通常被認(rèn)為是循環(huán)往復(fù)的，有漲落規(guī)律的變化。因此在周期項(xiàng)堆棧中(c)部分的線性層和偏置層被設(shè)置為傅里葉級(jí)數(shù)形式

y=S·x

(2)

(3)

式中，t為時(shí)間序列向量；H為向前預(yù)測(cè)的窗口長(zhǎng)度；x為該塊區(qū)中線性層和偏置層的輸入；y為該塊區(qū)中線性層和偏置層的輸出。

2 研究區(qū)域與數(shù)據(jù)源介紹

2.1 研究區(qū)域介紹

新鋪滑坡區(qū)是位于三峽庫(kù)區(qū)的安坪鄉(xiāng)新鋪村境內(nèi)的Ⅰ級(jí)地質(zhì)災(zāi)害區(qū)，研究區(qū)位置示意如圖2所示。圖2(a)、(b)顯示了新鋪滑坡研究區(qū)在重慶市及其與長(zhǎng)江的相對(duì)位置，圖2(c)顯示了滑坡體的具體分布情況?；乱装l(fā)區(qū)地處長(zhǎng)江右岸。新鋪滑坡由于其多期活動(dòng)，整個(gè)斜坡的坡面形成了多級(jí)臺(tái)地地貌，滑坡從上至下共有3處，分別是大坪滑坡、上二臺(tái)滑坡和下二臺(tái)滑坡。

圖2 新鋪滑坡研究區(qū)Fig.2 Xinpu landslide research area

大坪滑坡區(qū)域的巖體呈現(xiàn)出裂隙發(fā)育，傾角為20°～28°，該區(qū)域的變形破壞特征并不明顯，為順層基巖滑坡破壞模式。

上二臺(tái)滑坡區(qū)域整體土體的厚度變化較大，在降雨入滲的影響下，土體會(huì)沿著滑面產(chǎn)生緩慢的變形破壞，形成推移式淺層土體滑坡。

下二臺(tái)滑坡區(qū)域由于前緣的長(zhǎng)江河谷創(chuàng)造了較好的臨空條件，但受到地形、降水、長(zhǎng)江水位波動(dòng)及掏蝕等的作用，表面松散物質(zhì)沿坡面發(fā)生牽引式滑移破壞，最下層至江邊地段受長(zhǎng)江河水漲落和沖刷作用產(chǎn)生蠕滑變形。

該區(qū)域滑坡評(píng)估結(jié)果顯示變形破壞模式為順層基巖滑坡破壞模式、推移式淺層土體滑坡破壞模式和牽引式淺層土體滑坡破壞模式，非常適合采用InSAR技術(shù)進(jìn)行滑坡變形的長(zhǎng)期監(jiān)測(cè)。因此本文試驗(yàn)以重慶市奉節(jié)縣安坪鄉(xiāng)新鋪村的新鋪滑坡作為研究對(duì)象。

2.2 數(shù)據(jù)源介紹

本文試驗(yàn)針對(duì)新鋪滑坡形變的監(jiān)測(cè)預(yù)報(bào)收集的原始合成孔徑雷達(dá)(synthetic aperture radar,SAR)影像取自哨兵一號(hào)(Sentinel-1)衛(wèi)星星座。該星座由1A和1B兩顆星組成，每單個(gè)衛(wèi)星12 d可覆蓋地球一次，其上的C波段合成孔徑雷達(dá)能夠全天時(shí)全天候地提供對(duì)地觀測(cè)影像。

Sentinel-1A地球監(jiān)測(cè)衛(wèi)星于2014年發(fā)射并試運(yùn)行，Sentinel-1B地球監(jiān)測(cè)衛(wèi)星于2016年發(fā)射。目前，Sentinel-1數(shù)據(jù)和影像已經(jīng)能夠提供一系列的運(yùn)營(yíng)服務(wù)，例如海洋測(cè)繪、地面運(yùn)動(dòng)風(fēng)險(xiǎn)監(jiān)測(cè)、森林制圖，以及地籍管理和測(cè)繪等眾多應(yīng)用。

本文試驗(yàn)獲取了由Sentinel-1A提供的時(shí)間序列SAR數(shù)據(jù)，獲取了由2016年8月1日—2020年10月15日時(shí)間間隔為12 d的124幅SAR影像，數(shù)據(jù)相關(guān)信息見表2。

表2 Sentinel-1A時(shí)序SAR數(shù)據(jù)信息

3 新鋪滑坡形變預(yù)測(cè)試驗(yàn)

3.1 N-BEATS網(wǎng)絡(luò)模型設(shè)置

N-BEATS模型預(yù)測(cè)數(shù)據(jù)的結(jié)構(gòu)由窗口大小決定，根據(jù)窗口大小來決定由向后的k個(gè)歷史時(shí)刻作為輸入序列，向前預(yù)測(cè)m個(gè)未來時(shí)刻的輸出結(jié)果，本文試驗(yàn)需要預(yù)先根據(jù)時(shí)序數(shù)據(jù)特征設(shè)置窗口大小，即輸入數(shù)據(jù)窗口為x1,x2,…,xk，標(biāo)簽值為xk+1,xk+2,…,xk+m，以此類推至最后一個(gè)窗口的輸入數(shù)據(jù)為xn-9,xn-8,…,xn，標(biāo)簽值為xn+1,xn+2,…,xn+m，其中n為總數(shù)據(jù)長(zhǎng)度。數(shù)據(jù)需要滿足等間隔要求，若存在缺失數(shù)據(jù)的情況則需要通過插值方法補(bǔ)全。

試驗(yàn)設(shè)置3個(gè)堆棧，分別處理季節(jié)性、趨勢(shì)性及其他形變信號(hào)，每堆棧設(shè)置塊區(qū)數(shù)為2，每個(gè)塊區(qū)中設(shè)置4層全連接層，其中隱藏層單元均為64個(gè)，隨全連接層之后的線性層和偏置層參數(shù)維度設(shè)置分別為[2,2,3](季節(jié)堆棧模型、趨勢(shì)堆棧模型和其他形變堆棧模型)。N-BEATS滑坡預(yù)測(cè)網(wǎng)絡(luò)模型結(jié)構(gòu)如圖3所示。模型訓(xùn)練過程采用小批量為16的mini-batch梯度下降，訓(xùn)練所有樣本共400次。

圖3 N-BEATS滑坡預(yù)測(cè)網(wǎng)絡(luò)模型結(jié)構(gòu)Fig.3 N-BEATS landslide prediction network model structure

損失函數(shù)用于估量模型的預(yù)測(cè)值與標(biāo)簽值的差異程度，深度學(xué)習(xí)的目標(biāo)是利用優(yōu)化器來最小化損失函數(shù)的計(jì)算值，使得模型預(yù)測(cè)值與標(biāo)簽值的差異縮小至最小，從而訓(xùn)練得到最好的模型，本文中N-BEATS網(wǎng)絡(luò)模型所使用到的損失函數(shù)是平均絕對(duì)誤差損失函數(shù)(mean absolute error loss，MAELoss)，平均絕對(duì)誤差描述的是模型預(yù)測(cè)值與樣本標(biāo)簽值之間距離的平均值。優(yōu)化器為Adam，學(xué)習(xí)率設(shè)置為0.002。Adam是通過計(jì)算梯度的一階矩估計(jì)和二階矩估計(jì)來為不同參數(shù)設(shè)計(jì)獨(dú)立的自適應(yīng)學(xué)習(xí)率的一種優(yōu)化算法。在實(shí)踐中相較于其他的隨機(jī)優(yōu)化算法具有更好的性能，以上超參數(shù)設(shè)置見表3，本文中使用的網(wǎng)絡(luò)超參數(shù)應(yīng)用了N-BEATS網(wǎng)絡(luò)模型的推薦參數(shù)進(jìn)行經(jīng)驗(yàn)性多次實(shí)驗(yàn)后微調(diào)獲得。

表3 N-BEATS網(wǎng)絡(luò)模型設(shè)置一覽表

3.2 預(yù)測(cè)任務(wù)數(shù)據(jù)集設(shè)置

為探究預(yù)測(cè)數(shù)據(jù)集輸入輸出窗口設(shè)置對(duì)預(yù)測(cè)效果的影響，本文采取對(duì)數(shù)據(jù)進(jìn)行簡(jiǎn)單時(shí)序分解的策略，利用STL算法將研究目標(biāo)時(shí)間序列分解成為趨勢(shì)因子、季節(jié)因子及隨機(jī)因子，分解后P點(diǎn)的序列特征如圖4所示，圖4(a)、(b)分別為P點(diǎn)時(shí)間序列分解特征及P點(diǎn)點(diǎn)位示意?？梢娫摶聲r(shí)間序列數(shù)據(jù)成分大致為84%的趨勢(shì)項(xiàng)、8%的周期項(xiàng)及8%的不確定因子。數(shù)據(jù)表現(xiàn)出持續(xù)下降的走勢(shì)(遠(yuǎn)離衛(wèi)星方向)，并存在以年為周期的季節(jié)性形變，除此之外還在4年間有過兩次突發(fā)性形變發(fā)生，但該突發(fā)性影響所帶來的形變并不是很大，未超過1 cm。

圖4 形變時(shí)間序列分解特征Fig.4 Deformation time series decomposition feature diagram

上述時(shí)序分解成分中，趨勢(shì)項(xiàng)成分在整個(gè)時(shí)間序列中的特性表現(xiàn)為穩(wěn)定不變，無須在輸入輸出窗口大小的設(shè)置中考慮該成分；不確定成分呈現(xiàn)無規(guī)律發(fā)展的變化形式，無法在設(shè)置窗口時(shí)照顧到該成分的變化特征；周期性成分在整個(gè)時(shí)間序列中表現(xiàn)出明顯的局部性特征，因此該成分是輸入輸出窗口大小設(shè)置的重要依據(jù)，若要求網(wǎng)絡(luò)模型學(xué)習(xí)到一個(gè)整周期的形變趨勢(shì)，則周期性成分也是限制使用者需要提供訓(xùn)練樣本最低數(shù)量的重要依據(jù)，若對(duì)形變趨勢(shì)的周期性形變捕捉能力不做要求則無須滿足此數(shù)量條件。以年為周期的周期性成分在一年的時(shí)間序列中的長(zhǎng)度約為30，考慮到總樣本長(zhǎng)度為124，為保證訓(xùn)練數(shù)據(jù)的樣本集豐富的同時(shí)兼顧序列的周期性特征，分別將輸入窗口設(shè)置為[4,6,…,18,20]，輸出窗口設(shè)置為[1,2,…,5,6]，對(duì)這54種組合進(jìn)行訓(xùn)練并測(cè)試得到預(yù)測(cè)結(jié)果如圖5所示，并對(duì)預(yù)測(cè)結(jié)果計(jì)算均方根誤差見表4，可見隨著輸入和輸出窗口的增大，預(yù)測(cè)結(jié)果的準(zhǔn)確度在逐漸降低。輸入窗口越大，表明在模型訓(xùn)練時(shí)被考慮到的周期性成分更加完備，但由于同時(shí)縮小了數(shù)據(jù)集而使得模型的學(xué)習(xí)數(shù)據(jù)減少，訓(xùn)練效果反而降低；輸出窗口越大，表明期待預(yù)測(cè)到未來未知時(shí)間序列數(shù)據(jù)的長(zhǎng)度越長(zhǎng)，即更加“貪心”，則也相應(yīng)地會(huì)付出預(yù)測(cè)精度下降的代價(jià)。

圖5 54種輸入輸出窗口組合的預(yù)測(cè)結(jié)果Fig.5 Prediction results of 54 combinations of input and output windows

對(duì)比之下，本文選擇了輸入輸出窗口為[10,2]的數(shù)據(jù)集組織模式，該模式能夠在滿足預(yù)測(cè)中短期滑坡形變時(shí)長(zhǎng)的條件下達(dá)到最高精度要求。

3.3 滑坡預(yù)測(cè)結(jié)果

經(jīng)過TCP-InSAR技術(shù)處理后獲得形變區(qū)位移圖，選取影像中整體相干性大于0.8的5148個(gè)樣本點(diǎn)，將其展開并重組為時(shí)間序列數(shù)據(jù)。將數(shù)據(jù)集按照2∶1比例劃分為訓(xùn)練集數(shù)據(jù)和測(cè)試集數(shù)據(jù)，為保證預(yù)測(cè)效果表現(xiàn)力，未將測(cè)試數(shù)據(jù)集打亂順序。根據(jù)模型運(yùn)行的實(shí)際情況發(fā)現(xiàn)窗口大小為[5,1]，設(shè)置參數(shù)維度為[2,8,3]，迭代次數(shù)在300時(shí)已經(jīng)能夠穩(wěn)定地呈現(xiàn)出較為理想的預(yù)測(cè)結(jié)果。

A、B、C3點(diǎn)位置在圖6中體現(xiàn)，其形變趨勢(shì)如圖7所示，藍(lán)色曲線為數(shù)據(jù)樣本點(diǎn)的歷史觀測(cè)值，橙紅色曲線為樣本點(diǎn)未來趨勢(shì)的真實(shí)標(biāo)簽值，黃色曲線為N-BEATS網(wǎng)絡(luò)模型的預(yù)測(cè)結(jié)果?？梢姡瑪?shù)據(jù)預(yù)測(cè)結(jié)果遵循歷史數(shù)據(jù)的走向，說明了預(yù)測(cè)的正確性；放大預(yù)測(cè)結(jié)果與真實(shí)標(biāo)簽值的細(xì)節(jié)，可以看出預(yù)測(cè)結(jié)果與真實(shí)標(biāo)簽高度吻合，說明了該方法用于新鋪滑坡區(qū)的長(zhǎng)期時(shí)間序列預(yù)測(cè)的有效性。

圖6 新鋪滑坡區(qū)沉降區(qū)特征點(diǎn)位置Fig.6 Location of feature points in settlement area of Xinpu landslide area

圖7 A、B、C 3點(diǎn)形變趨勢(shì)Fig.7 Deformation trend diagram at A,B and C points

選取均方根誤差(RMSE)、平均絕對(duì)誤差(MAE)及標(biāo)準(zhǔn)差(SD)3種時(shí)間序列評(píng)價(jià)指標(biāo)對(duì)預(yù)測(cè)結(jié)果和真實(shí)標(biāo)簽的差異進(jìn)行評(píng)估，計(jì)算結(jié)果見表5?？梢姡琋-BEATS網(wǎng)絡(luò)模型預(yù)測(cè)形變趨勢(shì)結(jié)果的各項(xiàng)誤差指標(biāo)均滿足形變預(yù)測(cè)的要求，進(jìn)一步從真實(shí)數(shù)據(jù)試驗(yàn)的角度說明了該方法的可靠性。

表5 預(yù)測(cè)效果評(píng)價(jià)

現(xiàn)將時(shí)序預(yù)測(cè)得到的部分結(jié)果陳列如圖8所示。圖8列出了2019年8月22日、2019年12月8日、2020年3月25日及2020年6月17日的預(yù)測(cè)形變結(jié)果、真實(shí)觀測(cè)的形變結(jié)果及二者之間的殘差圖?？梢园l(fā)現(xiàn)，在直觀的目視判別中二者幾乎相差無幾，雖然在局部出現(xiàn)了殘差高值，但在整體范圍內(nèi)仍然表現(xiàn)出的是較為均一且不超過±1 cm的預(yù)測(cè)偏差，說明預(yù)測(cè)效果較好。

圖8 部分時(shí)間序列形變預(yù)測(cè)結(jié)果對(duì)比及其殘差圖Fig.8 Comparison and residual diagram of partial time series deformation prediction results

3.4 方法對(duì)比

圖9為STL、支持向量回歸(support vector regression,SVR)、無損卡爾曼濾波(unscented Kalman filter,UKF)、一次指數(shù)平滑法(single exponential smoothing,SES)及N-BEATS 5種時(shí)序預(yù)測(cè)方法在不同噪聲水平數(shù)據(jù)下的應(yīng)用效果對(duì)比，該節(jié)所選取的數(shù)據(jù)仍為3.1節(jié)中提到的P點(diǎn)時(shí)間序列數(shù)據(jù)?？梢姳疚膽?yīng)用的N-BEATS網(wǎng)絡(luò)模型在原始時(shí)序數(shù)據(jù)及添加毫米級(jí)和厘米級(jí)誤差序列數(shù)據(jù)的應(yīng)用效果都優(yōu)于其他方法，達(dá)到毫米級(jí)的預(yù)測(cè)精度。在為原始時(shí)序數(shù)據(jù)添加毫米級(jí)噪聲后，N-BEATS模型能夠抵抗噪聲干擾，完成高精度預(yù)測(cè)如圖9(b)所示，然而對(duì)于添加厘米級(jí)噪聲數(shù)據(jù)的抗差性能降低如圖9(c)所示，可見該模型具有一定的抵抗數(shù)據(jù)噪聲的能力，這一特點(diǎn)表明本文提出的方法對(duì)于輸入形變序列的精度要求較低，即使序列中存在厘米級(jí)誤差，仍能夠獲得相對(duì)更優(yōu)的預(yù)測(cè)成果從而擴(kuò)展N-BEATS神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)的應(yīng)用廣度。

圖9 N-BEATS與4種時(shí)序預(yù)測(cè)方法在不同噪聲水平數(shù)據(jù)下的應(yīng)用效果對(duì)比Fig.9 Comparison of N-BEATS and four time series prediction methods in different noise level data

3.5 預(yù)測(cè)區(qū)間估計(jì)

由于點(diǎn)估計(jì)存在誤差，僅對(duì)滑坡形變做出點(diǎn)估計(jì)仍不充分，因此為衡量N-BEATS神經(jīng)網(wǎng)絡(luò)模型預(yù)測(cè)滑坡形變數(shù)據(jù)的可靠度和精度水平，本文對(duì)N-BEATS網(wǎng)絡(luò)模型的預(yù)測(cè)結(jié)果進(jìn)行了置信區(qū)間估計(jì)。假設(shè)N-BEATS網(wǎng)絡(luò)模型對(duì)InSAR數(shù)據(jù)的滑坡形變預(yù)測(cè)結(jié)果服從正態(tài)分布X～N(μ,σ2)，給定可信度為95%，利用N-BEATS網(wǎng)絡(luò)模型的300個(gè)對(duì)P點(diǎn)時(shí)序數(shù)據(jù)預(yù)測(cè)結(jié)果求取均值μ的置信區(qū)間估計(jì)如圖10所示。對(duì)置信上限和置信下限分別計(jì)算均方根誤差得到該模型預(yù)測(cè)的精度水平，可見對(duì)滑坡形變做300次預(yù)測(cè)的數(shù)據(jù)以95%的可靠性保證了預(yù)測(cè)誤差在±10 mm之內(nèi)。

圖10 N-BEATS神經(jīng)網(wǎng)絡(luò)置信區(qū)間估計(jì)Fig.10 N-BEATS Neural network confidence interval estimation

4 結(jié)論與展望

本文拓展了一種利用N-BEATS深度學(xué)習(xí)網(wǎng)絡(luò)模型預(yù)測(cè)可解釋性滑坡時(shí)間序列形變的方法，試驗(yàn)結(jié)果展示出該方法在滑坡預(yù)測(cè)應(yīng)用中高可靠性、高精度及具有一定抗差能力的優(yōu)勢(shì)。三峽庫(kù)區(qū)新鋪滑坡真實(shí)試驗(yàn)結(jié)果表明以下4點(diǎn)。

(1) 針對(duì)本文所選取的時(shí)間間隔為12 d一景的Sentinel-1數(shù)據(jù)，選取單次預(yù)測(cè)輸入輸出窗口大小為[10,2]的數(shù)據(jù)結(jié)構(gòu)能夠在滿足中短期滑坡預(yù)測(cè)的前提下取得較好的預(yù)測(cè)效果。

(2) 新鋪滑坡應(yīng)用N-BEATS網(wǎng)絡(luò)模型進(jìn)行預(yù)測(cè)，以均方根誤差為1.1 mm的高預(yù)測(cè)精度完成了滑坡預(yù)測(cè)工作，證實(shí)了N-BEATS網(wǎng)絡(luò)模型在時(shí)間序列預(yù)測(cè)任務(wù)上的有效性、可靠性及其工程實(shí)踐性。

(3) 對(duì)原始滑坡形變時(shí)序數(shù)據(jù)添加毫米級(jí)和厘米級(jí)誤差驗(yàn)證N-BEATS網(wǎng)絡(luò)模型預(yù)測(cè)的抗差性能，與其余4種傳統(tǒng)時(shí)間序列預(yù)測(cè)方法相比仍能夠保持最高精度水平的預(yù)測(cè)，特別是針對(duì)存在毫米級(jí)誤差的時(shí)間序列數(shù)據(jù)該方法體現(xiàn)出了良好的降低偏差預(yù)測(cè)功能。

(4) 利用N-BEATS網(wǎng)絡(luò)模型對(duì)滑坡形變做300次預(yù)測(cè)，分析了該方法預(yù)測(cè)結(jié)果數(shù)據(jù)能夠以95%的可靠性保證預(yù)測(cè)誤差在±10 mm之內(nèi)。

雖然利用N-BEATS網(wǎng)絡(luò)模型進(jìn)行新鋪滑坡的形變預(yù)測(cè)取得了良好的效果，但在滑坡智能預(yù)測(cè)方面仍存在著許多有待解決的問題，例如預(yù)測(cè)時(shí)間跨度取決于訓(xùn)練數(shù)據(jù)的時(shí)間跨度而無法滿足超越InSAR時(shí)間分辨率的短期預(yù)報(bào)需求；目前的時(shí)序預(yù)測(cè)模型仍停留在逐點(diǎn)預(yù)測(cè)的階段，但不可否認(rèn)的是逐點(diǎn)預(yù)測(cè)非常耗時(shí)且忽略了空間相關(guān)關(guān)系這一信息；單變量預(yù)測(cè)的突發(fā)性滑坡預(yù)測(cè)結(jié)果尚存在時(shí)間延遲問題等，在未來越來越多的SAR衛(wèi)星(座)成功發(fā)射之后，有可能實(shí)現(xiàn)利用InSAR形變數(shù)據(jù)預(yù)測(cè)突發(fā)性形變的愿景。此外，InSAR數(shù)據(jù)與其他高時(shí)間分辨率數(shù)據(jù)融合也有助于實(shí)現(xiàn)形變序列數(shù)據(jù)的突發(fā)性預(yù)測(cè)。