魯鐵定 李 禎
1 東華理工大學測繪工程學院,南昌市廣蘭大道418號,330013
近20 a來,全球IGS基準站不斷積累的時間序列數據為大地測量和地球動力學研究提供了寶貴的數據基礎[1-3]。這些數據可有效反映由地球物理效應引起的長期變化趨勢和非線性變化[4],因此對GNSS坐標時間序列進行分析有助于地殼板塊運動監測[5]、大壩或橋梁變形監測[6-7]、全球或區域坐標系維護[8-9]等領域的發展。通過分析GNSS坐標時間序列,可以預測連續時間點坐標,為判斷運動趨勢提供重要依據。研究發現,GNSS坐標時間序列中高程方向的噪聲分量通常大于水平方向,且噪聲組合模型豐富,因此構建高精度的GNSS高程時間序列預測模型具有較大困難[10-12]。
目前,已有許多學者在GNSS高程坐標時間序列預測研究中基于機器學習、深度學習等算法構建相應的預測模型,并針對GNSS高程坐標時間序列的信號特征進行高精度的預測研究工作[11]。但是相關研究仍存在一定不足:1)模型預測精度容易受到異常值以及缺失值影響,抗干擾能力較弱,導致使用模型時對數據質量具有較高的要求。2)模型單次預測精度無法保證隨著預測時間跨度的增長持續保持在較高水平上,在進行長期預測時需要進行滾動預測以保持模型精度,從而會增大模型的運行載荷。3)在相同參數設置和數據集下,部分模型的預測結果存在差異,并且無法進行動態調整。
針對這些問題,本文提出一種Prophet-XGBoost組合預測模型用于GNSS高程時間序列預測,并通過SOPAC公布的實測高程數據對本文方法的適用性和有效性進行檢驗。首先檢驗XGBoost模型和Prophet模型在GNSS時間序列中的適應度;再基于單模型在長期預測中表現出的特點,針對XGBoost模型對于趨勢項部分預測能力較弱和Prophet模型對于非線性部分預測能力較弱的問題,構建Prophet-XGBoost組合預測模型;最后通過ALGO、ALRT和BRST基準站U分量坐標時間序列數據對模型預測精度進行檢驗。
XGBoost作為一種集成學習算法,可通過構建多個CART樹對數據集進行預測,然后集成多個樹模型的預測結果,從而輸出最終預測數據。與同屬于樹模型算法的GBDT算法相比,XGBoost算法更為高效,并且在算法方面也進行過改進。XGBoost使用二階泰勒展開式逼近目標函數的泛化誤差部分,可有效簡化目標函數的計算;XGBoost還可通過在目標函數中加入正則項來降低模型預測的波動性以及改善模型過擬合現象[13]。

(1)
式中,fk表示決策樹,fk(xi)表示第k棵樹賦予第i個觀測值的分數。同時,使用函數fk時,下述正則化目標函數應被最小化:
(2)
式中,l為損失函數。為防止模型過于復雜,在模型中設置懲罰項Ω:
(3)
式中,γ可控制懲罰項中枝葉數量T,λ可控制懲罰項中枝葉重量w。設置Ω(fk)項不僅可以簡化算法生成模型,同時可以防止模型過擬合。
XGBoost算法采用迭代法最小化目標函數。模型通過增加fj項在第j次迭代后得到減小的目標函數:
(4)
式(4)可以通過泰勒展開進行簡化,并且可以推導出決策樹從給定節點分裂后的損失減少公式:
Lsplit=
(5)
式中,I為當前節點中可用觀測數據的一個子集,IL和IR分別為分裂后左右節點中可用觀測數據的一個子集。函數gi和hi可表示為:
(6)
(7)
從推導Lsplit公式中可找到任意給定節點的最佳分裂,該函數只依賴于損失函數和正則化參數γ。同時,基于式(5),XGBoost算法可以優化任何損失函數,并且提供一階和二階梯度。
Prophet模型是一種基于加性模型的時間序列預測模型,由趨勢項、周期項和假日項組成[14-15],表達式為:
y(t)=g(t)+s(t)+h(t)+εt
(8)
式中,g(t)為趨勢項;s(t)為周期項;h(t)為假日項,該項在GNSS時間序列中通常不存在,因此建模過程中不予考慮;εt為殘差項,即模型未能預測到的隨機趨勢。
在Prophet模型中,趨勢項包含邏輯回歸函數和分段線性函數兩個重要函數。根據GNSS高程時間序列的趨勢變化特點,選用邏輯回歸函數來計算趨勢項。同時,Prophet模型使用傅里葉級數表示時間序列的周期性:
(9)
式中,P為時間序列周期,P=365.25時為年周期項,P=7時為周周期項。因此,Prophet模型在進行GNSS日坐標時間序列分解時會分別輸出年周期項和周周期項。
使用平均絕對誤差(MAE)和均方根誤差(RMSE)作為模型預測精度的評價指標:
(10)
(11)

通過引入Pearson相關系數判斷預測時間序列和原始時間序列之間的相關性。Pearson相關系數將時間序列視為變量,從而計算兩個時間序列之間的相關性,其計算公式為:
(12)

為驗證Prophet-XGBoost模型的有效性,選取IGS基準站實測數據作為數據集。選用ALGO、ALRT和BRST站2010~2015年U分量日坐標時間序列作為實驗數據,該數據來源于SOPAC(scripps orbit and permanent array center)。數據基本情況如表1(單位mm)所示。

表1 數據基本情況
為檢驗XGBoost模型和Prophet模型對于GNSS坐標時間序列預測的適用性,引入經典時間序列預測模型ARIMA模型作為對比。選用ALGO站數據作為實驗數據,分別用ARIMA模型、XGBoost模型和Prophet模型對數據集進行預測,將預測值與原始數據值進行對比,得到相應的MAE值和RMSE值。
該實驗數據集中2010~2014年數據為訓練集,2015年數據為驗證集,共包含2 149個日坐標數據,缺失數據天數為42,數據完整率為98.08%,訓練中最大連續缺失數據時間跨度為7 d。
由表1和四分位距值(interquartile range,IQR)計算可知,該實驗數據集中含有少量異常值。因此,該數據訓練集中包含異常值、缺失值、連續缺失值等條件,可以更好地檢驗模型對于GNSS實測數據預測的適用性。各模型精度見表2(單位mm)。

表2 不同模型精度對比
由表2可知,相較于ARIMA模型,XGBoost模型和Prophet模型在MAE和RMSE兩個精度評價指標上均有優勢。基于3個模型的預測原理可知,在長時間跨度的GNSS坐標時間序列預測研究中,ARIMA模型精度更容易受到預測時間長度的影響,而XGBoost模型和Prophet模型的預測精度隨預測跨度變化而變化的幅度較小,且XGBoost模型的算法實現更為簡便,具有更高的預測效率。因此,XGBoost模型和Prophet模型均可應用于GNSS坐標時間序列預測,并在實測數據實驗中具有適用性。圖1為ALGO站Prophet模型預測結果,其中離散點為原始時間序列值,曲線為擬合預測值。

圖1 ALGO站Prophet模型預測結果Fig.1 Prediction results of Prophet modelof ALGO station
基于XGBoost模型和Prophet模型在GNSS坐標時間序列預測中表現出的模型特性,構建Prophet-XGBoost組合預測模型,以更好地發揮兩個模型的預測優勢,削弱模型對于GNSS時間序列部分分解項預測能力不足的影響,進而實現更高精度的GNSS坐標時間序列預測。圖2為Prophet-XGBoost模型構建流程。

圖2 Prophet-XGBoost模型構建流程Fig.2 Construction process of Prophet-XGBoost model
Prophet-XGBoost模型實現流程如下:
1)通過Prophet模型對GNSS坐標時間序列U方向數據進行分解,將時間序列分解為趨勢項、周期項和殘差項,其中周期項由年周期項和周周期項等權相加構成。圖3為ALGO站Prophet模型分解結果。

圖3 ALGO站Prophet模型分解結果Fig.3 Decomposition results of Prophet modelof ALGO station
2)通過Prophet模型對分解出的趨勢項和周期項進行擬合,然后使用XGBoost模型對擬合時間序列進行預測,得到預測值yfit。圖4為ALGO站XGBoost模型擬合時間序列預測結果。

圖4 ALGO站XGBoost模型擬合時間序列預測結果Fig.4 Prediction results of fitting time series ofXGBoost model of ALGO station
3)使用XGBoost模型對Prophet模型分解時間序列過程中遺留的殘差項進行預測,得到預測值yresidual。
4)將兩個預測結果等權相加,得到預測值ypred:
ypred=yfit+yresidual
(13)
本文選用ALGO、ALRT、BRST三個IGS站U分量GNSS坐標時間序列作為實驗數據,分別采用XGBoost模型、Prophet模型和Prophet-XGBoost模型進行預測。其中,ALGO和ALRT站緯度差異較大,ALGO和BRST站經度差異較大。預測精度對比結果如表3所示。

表3 預測精度對比
由表3可知,相較于XGBoost模型和Prophet模型,ALGO站Prophet-XGBoost模型的MAE值分別降低7.4%和11.5%,RMSE值分別降低16.6%和3.5%;ALRT站Prophet-XGBoost模型的MAE值分別降低7.1%和20.5%,RMSE值分別降低7.2%和6.5%;BRST站Prophet-XGBoost模型的MAE值分別降低10.5%和9.2%,RMSE值分別降低5.6%和3.5%。相較于Prophet模型,組合模型的RMSE值降低幅度較小,這是因為Prophet-XGBoost模型進行預測時的擬合數據是由Prophet模型擬合所得。對Prophet-XGBoost模型預測結果和原始數據值進行相關性檢驗,通過引入Pearson相關系數計算ALGO、ALRT和BRST三個站的ρ值分別為0.53、0.64和0.49,其中ALGO和BRST兩個站預測結果和原始數據呈中等相關,ALRT站預測結果和原始數據呈強相關。
以ALGO站為例,對比分析XGBoost模型、Prophet模型和Prophet-XGBoost模型預測結果誤差分布。圖5為各模型預測結果絕對誤差分布情況。

圖5 各模型預測結果絕對誤差分布Fig.5 Distribution of absolute error of predictionresults of each model
由圖5可知,Prophet-XGBoost組合模型的較高精度預測結果占比最高,即相較于Prophet和XGBoost單模型預測結果,組合模型準確率得到提升。圖6為模型預測結果誤差對比。

圖6 模型預測結果誤差對比Fig.6 Error comparison of prediction resultsof each model
由圖6可知,相較于XGBoost模型,Prophet-XGBoost組合模型的誤差變化幅度和凸點值更??;Prophet模型和Prophet-XGBoost組合模型的誤差變化趨勢較為一致,但Prophet-XGBoost組合模型誤差更小。結合圖4結果可知,Prophet模型對于時間序列的趨勢變化較為敏感。
本文基于XGBoost模型和Prophet模型提出一種GNSS高程時間序列預測組合模型Prophet-XGBoost。在數據處理階段,利用Prophet模型對時間序列進行分解,將時間序列分解成3個主要組成部分,有利于XGBoost模型在預測時提取時間序列中不同組成部分的特征,減弱預測模型提取特征時受到其他組成部分的影響。得益于此,XGBoost模型可以進行更高精度的GNSS坐標時間序列預測工作,也可彌補兩個模型在單模型預測時存在的問題。此外,在計算預測結果時,也可依據XGBoost模型提供的特征分數進行定權相加。
XGBoost模型和Prophet模型在實驗過程中也存在需要改進之處:XGBoost模型對于GNSS坐標時間序列的趨勢變化不夠敏感以及Prophet模型對于GNSS坐標時間序列非線性部分預測能力較弱。同時考慮到XGBoost模型的特征項構建具有很大靈活性,研究人員可以使用非GNSS坐標時間序列數據輔助進行特征構造,如水文載荷對GNSS坐標時間序列U方向振幅具有較為顯著的影響[17]。因此,在以后研究中可以將水文載荷等因素納入特征構造條件中,從而通過更多的特征子集進行模型構建和預測?;谠撎匦?,Prophet-XGBoost模型可以作為通用預測模型,將其推廣到其他研究領域,如海平面時間序列預測、天氣預報和發病率預測,可極大增強該模型的適用性。