基于EWT-Prophet方法的地表沉降預測

周子琪 周世健 陶 蕊

1 東華理工大學測繪工程學院，南昌市廣蘭大道418號，330013 2 南昌航空大學校長辦公室，南昌市豐和南大道696號，330063

大規模的地下礦產開采活動易引發地面塌陷、滑坡、泥石流等地質災害，破壞地質結構，嚴重威脅礦區周邊居民的生命和財產安全[1]，因此對礦區地表進行長期的有效監測和精準預測尤為重要[2]。目前，已有學者對地表沉降預測進行研究，并取得一定的進展[3-4]。傳統的灰色類預測模型雖然使用廣泛，但適用性較差；而ARMA、ARIMA等經典模型的預測精度不穩定，需滾動預測才能保持預測精度；神經網絡預測模型雖有較高的預測精度，但存在預測過程不穩定、參數設置復雜、自適應能力較差等問題。針對上述問題，本文提出將Prophet模型引入地表沉降預測研究。Prophet模型是一種自適應分解預測、擬合插值的新模型，對原始數據有較好的自適應性，使用靈活且無需對缺失項進行插值，也無需任何先驗條件，擬合過程速度較快，有良好的預測精度[5]。本文在深入研究Prophet模型原理后，引入經驗小波變換(EWT)方法，基于EWT和Prophet模型構建一種分解-預測的自適應預測新方法，以達到提升預測精度的目的[6]。

1 研究區概況及數據來源

本文以德興礦區為研究對象，研究區位于江西省上饒市德興境內，屬于江南丘陵地區，地勢東南部偏高、西北部偏低，地貌整體起伏較大。實驗數據選取覆蓋德興礦區的30景C波段升軌Sentinel-1A影像數據，時間跨度為2018-12-29～2019-12-24，時間間隔為12 d，具體參數見表1。采用由美國航空航天局(NASA)提供的30 m分辨率SRTM1數據作為外部參考DEM，以消除地形誤差對實驗的影響。

表1 Sentinel-1A數據參數

2 SBAS沉降監測及組合預測模型建立

2.1 SBAS-InSAR技術原理

SBAS-InSAR是以傳統D-InSAR技術為基礎，用來獲取工作區地表形變時間序列的技術。該技術選取N+1幅SAR影像，通過設置合適的時間基線和空間基線，生成N個小基線集差分干涉對，去除地形相位后生成差分干涉圖并進行相位解纏，將常規D-InSAR監測的觀測結果用作單個觀測值，采用奇異值分解(SVD)方法將每組小基線數據集進行連接，以解決在時域上采樣過于稀疏的問題。同時結合穩定散射體的干涉相位信息，以獲取更高的空間分辨率，最后根據最小二乘準則計算高精度的沉降形變時間序列。

2.2 Prophet模型原理

Prophet預測模型[7]是一種高效分析時序信號的新模型，可以處理時序信號分析中產生的丟失值、異常值，對不同跨度樣本的時序信號進行長短期預測，具備較好的預測精度、良好的自適應能力和一定的抗差性和魯棒性。Prophet模型已經在GNSS高程信號預測[5]、電力需求預測[8]、空氣指數預測[9]和TEC電離層預測[10]等領域得到廣泛應用，并取得較好的預測效果，但幾乎沒有涉及地表沉降預測領域。Prophet模型的詳細原理見文獻[7]。

Prophet采用廣義加法模型來進行擬合平滑和函數預測，將地表沉降時序數據自適應分解為：

y(t)=g(t)+s(t)+εt

(1)

式中，y(t)為原始地表沉降時序數據；g(t)為沉降時序數據趨勢項，表示時間序列非線性增長(非周期項)部分的變化函數；s(t)為高程方向的震蕩周期項；εt為服從正態分布的殘差項，可表示未預測到的隨機噪聲或趨勢。SBAS-InSAR監測數據復雜，趨勢項采用邏輯回歸函數，表示為：

(2)

式中，k為地表沉降趨勢增長率；m為位移量參數；c為趨勢上限值，隨著時間t的增加，g(t)趨于c。在地殼板塊運動中，s(t)主要包含季節性的周期震蕩和微小的趨勢變化，其擬合函數通過沉降時序數據的傅里葉級數進行構造：

(3)

2.3 EWT方法原理

經驗小波變換(EWT)是一種計算量小且具備較強魯棒性的新變換方法，其核心思想是依據原始信號中的頻譜特征，進行自適應分割后構建基于傅里葉變換的正交小波濾波器組，提取出具有緊支撐傅里葉頻譜的不同AM-FM(調幅-調頻)成分。具體過程為：

1)對信號進行傅里葉變換，求得支撐區間(0，π)的傅里葉頻譜F(ω)；

2)對于自適應分割傅里葉頻譜F(ω)，依據香農定理將頻譜分解為N個頻帶與N-1個分界頻率，圖1為傅里葉頻譜分割示意圖；

圖1 傅里葉頻譜分割Fig.1 Fourier spectrum segmentation

3)根據各分界頻率構建經驗小波φ(ω)，確定經驗尺度函數和經驗小波函數；

4)對F(ω)*φ(ω)反傅里葉變換，得到不同的模態分量。

(4)

2.4 EWT-Prophet模型的沉降預測方法

圖2為EWT-Prophet組合預測模型算法的流程，圖3為EWT分解示意圖。可以看出，EWT將原始沉降數據自適應分解為1個趨勢項信號和3個周期項信號，其中趨勢項信號可有效反映原始沉降數據的趨勢性，其余周期項信號可反映原始沉降數據的微小震蕩變化。

圖2 EWT-Prophet流程Fig.2 Flow chart of EWT-Prophet

圖3 EWT分解Fig.3 EWT decomposition

本文組合模型先對沉降時序信號進行EWT分解，將原始時序信號分解為經驗尺度分量F0和N個經驗小波分量FN，再對各分量進行Prophet預測，將得到的所有預測分量疊加重構為最終的預測信號。該方法綜合了EWT方法自適應能力強、理論性強和Prophet模型擬合預測效果好、預測速度快等特點。EWT-Prophet組合預測模型對地表沉降時序數據進行預測的主要步驟為：

1)使用EWT對沉降數據進行分解：

(5)

2)對所有有效分量進行Prophet預測，將所有預測分量重構為最終的預測數據，并與單一的Prophet預測方法進行比較，驗證組合模型的可靠性。

3 實驗分析

3.1 礦區SBAS-InSAR監測及分析

實驗選取覆蓋研究區的30景升軌Sentinel-1A數據，使用GAMMA軟件進行SBAS-InSAR處理，選取2019-06-15的SAR影像作為公共主影像，再將輔影像配準到主影像上(配準精度達0.001像元)，對干涉對進行多視處理(多視系數設置為4∶1)。實驗的時間基線和空間基線閾值分別為60 d和200 m，最終構成112對干涉對用于時序分析，圖4為干涉對基線分布。根據干涉對時空基線分布圖對影像進行干涉處理，利用外部DEM數據去除地形相位的影響，生成差分干涉圖，并對其進行濾波增強處理，采用最小費用流方法進行相位解纏，得到研究區相位的真實值。通過查驗每組產生的相干系數圖、濾波后的干涉圖和解纏圖，移除干涉質量較差的干涉對，在具有高相干性的像素點上建立模型。通過奇異值分解方法反演估算形變速率并去除殘余地形，計算相位殘差并進行殘差分離，最終得到研究區地表沉降速率及地表沉降分布情況。

圖4 干涉對時空基線分布Fig.4 Spatial-temporal baseline of interferometric pairs distribution

實驗結果表明，研究區有3處明顯沉降，本文將這3處沉降區域分別標記為A、B、C，并在每個區域中分別選取3個沉降特征點作為研究對象(分別標記為#1、#2、#3)，各特征點位置如圖5所示。圖6為各區域中特征點的沉降時序，可以發現，A區域沉降最為嚴重，最大年平均沉降速率為490 mm/a，其中A#3位于沉降區中心，最大累積沉降量達440 mm，A#1和A#2在沉降區邊緣，沉降量分別為184 mm和131 mm；B區域最大年平均沉降速率達390 mm/a，位于沉降中心的特征點B#3累積沉降量為358 mm，沉降區邊緣的B#1和B#2累積沉降量分別為120 mm和118 mm；C區域最大年平均沉降速率可達233 mm/a，C#1、C#2和C#3的累積沉降量分別為148 mm、98 mm和117 mm。

圖5 A、B、C區域沉降速率Fig.5 Settlement rates in regions A, B and C

圖6 A、B、C區域特征點時間序列形變Fig.6 Time series deformation of feature points in areas A, B and C

馬濤等[12]對采用SBAS-InSAR技術獲取的時序沉降數據進行分析，并引入水準監測數據進行對比。結果表明，監測沉降數據的平均誤差為1.62 mm，符合《地面沉降干涉雷達數據處理技術規范》中10 mm的精度要求，表明SBAS-InSAR技術監測沉降數據的精度具有可靠性。本文將監測數據作為真實觀測數據進行回溯性預測研究。

3.2 地表沉降預測及結果分析

以SBAS-InSAR技術獲得的特征點沉降時序數據作為實際值，選取前27期(時間間隔為12 d)監測數據作為樣本數據進行預測研究，預測后3期的沉降情況，并將預測結果與實際監測結果進行對比分析。為證明本文提出的組合預測模型在礦區沉降預測中的可行性，分別利用Prophet模型和EWT-Prophet組合模型進行實驗，并對實驗結果進行對比分析。

表2和表3為Prophet模型和EWT-Prophet組合模型的預測結果和相對誤差統計。可以看出，組合模型的平均相對誤差分別為1.55%、2.89%、1.57%，平均相對誤差相較于單一模型分別提升48.68%、53.24%、68.21%。由此可見，引入單一Prophet模型進行對比不足以驗證本文組合方法的有效性，因此引入經典ARMA預測模型進行對比分析。

表2 Prophet模型預測值

表3 EWT-Prophet組合模型預測值

表4(單位mm)為不同方法預測結果的統計，表5(單位mm)為預測結果殘差絕對值統計。可以看出，在殘差絕對值指標中，EWT-Prophet組合模型預測結果的殘差絕對值最小，3期平均值分別為2.516 7 mm、5.065 6 mm、2.877 8 mm，表明本文方法優于單一Prophet模型和ARMA模型。

表4 不同方法預測結果

表5 預測結果殘差絕對值

為分析本文預測方法的有效性，利用均方根誤差(RMSE)、平均百分比誤差(MAPE)等指標對預測結果進行對比分析(表6)。由表可知，EWT-Prophet組合模型的精度整體優于單一Prophet模型和ARMA模型，平均均方根誤差和平均百分比誤差分別為3.97 mm和2.00%，相較于Prophet模型和ARMA模型，RMSE分別提升51.53%、59.03%，MAPE分別提升57.81%、64.85%。實驗結果表明，本文構建的組合模型具有更好的預測精度，且具有有效性和適用性。

表6 預測精度比較

4 結 語

本文以江西德興礦區為研究區域，首先采用覆蓋研究區的30景影像進行SBAS地表沉降監測，獲取該地區30期時序沉降變化與累積沉降值。然后以SBAS監測值為研究時序數據樣本，使用前27期監測值作為樣本數據進行回溯性預測研究，利用EWT-Prophet組合模型與單一Prophet模型和ARMA模型進行預測研究。結果表明，組合模型結果整體優于單一模型，本文方法預測精度更佳，相較于Prophet模型和ARMA模型，RMSE分別提升了51.53%、59.03%，MAPE分別提升了57.81%、64.85%，說明組合模型具有適用性和有效性。

本文將EWT與Prophet預測模型相結合并引入SBAS-InSAR沉降研究，給類似的沉降時序信號數據分析提供了一種自適應強、預測效果好的組合方法，但此類分解-預測的組合方法也給預測過程帶來額外的工作量。在研究過程中，礦區的沉降變化與諸多因素有關，如何建立一種高效、自適應強，并能及時反映突變情況的預測方法還有待進一步研究。