基于Prophet-RF模型的GNSS高程坐標(biāo)時(shí)間序列預(yù)測(cè)分析

李 威 魯鐵定 賀小星 錢文龍

1 東華理工大學(xué)測(cè)繪工程學(xué)院，南昌市廣蘭大道418號(hào)，330013 2 華東交通大學(xué)土木建筑學(xué)院，南昌市雙港東大街808號(hào)，330013 3 中鐵第一勘察設(shè)計(jì)院集團(tuán)有限公司軌道交通工程信息化國(guó)家重點(diǎn)實(shí)驗(yàn)室,西安市西影路2號(hào),710043

近年來，不斷累積的IGS基準(zhǔn)站坐標(biāo)數(shù)據(jù)為大地測(cè)量學(xué)和地球動(dòng)力學(xué)研究提供了豐富的數(shù)據(jù)基礎(chǔ)，也為GNSS坐標(biāo)時(shí)間序列的研究提供了堅(jiān)實(shí)的數(shù)據(jù)支撐[1]。GNSS坐標(biāo)時(shí)間序列的精確預(yù)測(cè)對(duì)建筑物的變形監(jiān)測(cè)、氣候氣象的預(yù)測(cè)及區(qū)域內(nèi)地殼形變等研究具有重要意義。有研究指出，GNSS坐標(biāo)時(shí)間序列在3個(gè)方向都有明顯的周期性運(yùn)動(dòng)，其中高程方向的季節(jié)性變化最為顯著[2-3]。另外，由于高程方向的噪聲分量往往大于水平方向，且噪聲模型的組成并不單一，使得建立高精度GNSS高程坐標(biāo)時(shí)間序列預(yù)測(cè)模型的難度較大[4]。

諸多學(xué)者在GNSS坐標(biāo)時(shí)間序列的擬合與預(yù)測(cè)領(lǐng)域進(jìn)行了廣泛研究[5-8]，促進(jìn)了GNSS坐標(biāo)時(shí)間序列的研究進(jìn)展，但針對(duì)高程坐標(biāo)長(zhǎng)期預(yù)測(cè)模型的研究仍存在一些不足：1)現(xiàn)階段時(shí)間序列模型多存在建模復(fù)雜、穩(wěn)定性差、無(wú)法動(dòng)態(tài)調(diào)整模型參數(shù)等缺陷；2)模型對(duì)坐標(biāo)數(shù)據(jù)的要求較高，建模精度容易受粗差影響，且未顧及數(shù)據(jù)缺失等不利影響；3)構(gòu)建的模型缺乏長(zhǎng)期預(yù)測(cè)能力，在預(yù)測(cè)時(shí)需要滾動(dòng)預(yù)測(cè)才能保持模型的精度，增加了模型的操作難度與運(yùn)行時(shí)間。

針對(duì)這些問題，本文提出將Prophet模型應(yīng)用于GNSS高程坐標(biāo)時(shí)間序列中，該模型對(duì)數(shù)據(jù)有良好的適應(yīng)能力，在長(zhǎng)期預(yù)測(cè)中有較高的可靠性與精度，且無(wú)需太多關(guān)于時(shí)間序列的先驗(yàn)知識(shí)。本文首先檢驗(yàn)了Prophet模型在高程坐標(biāo)時(shí)間序列中的適用性；再針對(duì)Prophet模型對(duì)非線性部分預(yù)測(cè)能力較弱的缺陷，引入隨機(jī)森林(RF)模型，得到Prophet-RF組合預(yù)測(cè)模型；最后根據(jù)2種模型的特性設(shè)計(jì)多種組合方案,并利用BJFS站的實(shí)測(cè)高程坐標(biāo)數(shù)據(jù)檢驗(yàn)不同組合方案的精度與適用性。

1 基本原理與評(píng)價(jià)指標(biāo)

1.1 Prophet模型

Prophet模型由Taylor等[9]于2017年提出并發(fā)布了相關(guān)的開源軟件包，是一種可以分析時(shí)間序列的新模型，包括處理時(shí)間序列數(shù)據(jù)中的異常值和缺失值及對(duì)時(shí)間序列進(jìn)行長(zhǎng)短期預(yù)測(cè)并研究其變化規(guī)律。該模型在電離層異常探測(cè)[10]、市場(chǎng)銷量研究[11]、空氣質(zhì)量分析[12-13]等領(lǐng)域都有成功案例，但在GNSS坐標(biāo)時(shí)間序列方面幾乎沒有涉及。

Prophet模型在預(yù)測(cè)時(shí)的循環(huán)分析共分為4個(gè)部分：建立模型、預(yù)測(cè)評(píng)估、表現(xiàn)問題、可視化檢測(cè)預(yù)測(cè)結(jié)果。使用Prophet模型分析數(shù)據(jù)時(shí)將4個(gè)部分構(gòu)成一個(gè)循環(huán)體，根據(jù)不同專業(yè)及數(shù)據(jù)的具體要求對(duì)模型進(jìn)行人工調(diào)整，并結(jié)合模型的自動(dòng)預(yù)測(cè)功能，使Prophet模型較傳統(tǒng)的時(shí)間序列模型有更好的靈活性與適用性[9]。

Prophet模型是一個(gè)廣義加法模型，并且使用基于貝葉斯的曲線擬合方法對(duì)時(shí)間序列進(jìn)行平滑和預(yù)測(cè)[9-10]。模型由趨勢(shì)項(xiàng)、周期項(xiàng)及假日項(xiàng)組成，其加法模型的基本分解形式為：

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

(1)

其中，g(t)為趨勢(shì)項(xiàng)，表示時(shí)間序列中非周期部分的變化函數(shù)，包含不同程度的假設(shè)與可調(diào)節(jié)光滑度的參數(shù)[11]；s(t)為周期項(xiàng)，在Prophet模型中使用傅里葉級(jí)數(shù)對(duì)其進(jìn)行構(gòu)造，常以周、月、季度或年作為單位；h(t)為假日項(xiàng)，表示因假期或特殊原因造成的不規(guī)則影響，在GNSS坐標(biāo)時(shí)間序列中通常不存在此類情況，故在分解數(shù)據(jù)時(shí)可不必對(duì)假日項(xiàng)進(jìn)行建模；ε(t)為殘差項(xiàng)，服從正態(tài)分布，表示模型中未預(yù)測(cè)到的隨機(jī)趨勢(shì)。具體描述見文獻(xiàn)[9]。

1.2 隨機(jī)森林模型

隨機(jī)森林(random forest，RF)模型是一種具備集成思想的機(jī)器學(xué)習(xí)模型，是由一定數(shù)量的弱學(xué)習(xí)器決策回歸樹(classification and regression tree，CART)結(jié)合bagging集成學(xué)習(xí)理論和隨機(jī)特征子空間構(gòu)成的強(qiáng)學(xué)習(xí)器[14-15]。RF模型的隨機(jī)性表現(xiàn)在隨機(jī)選擇樣本及隨機(jī)選擇特征，在選擇樣本時(shí)采用的是bootstrap重抽樣方法，其特點(diǎn)是既能保證每個(gè)樣本在每次抽取時(shí)都有相同的概率被抽中，又能在不降低樣本訓(xùn)練規(guī)模的同時(shí)留出驗(yàn)證集[16]；而在特征選擇上，CART使用基尼系數(shù)進(jìn)行選擇，基尼系數(shù)越小表示特征越好。構(gòu)建RF模型的步驟如圖1所示。

圖1 隨機(jī)森林算法流程Fig.1 Random forest algorithm flow

RF模型的構(gòu)建流程為：

1)采用bootstrap重抽樣方法從原始數(shù)據(jù)集N中抽取ntree個(gè)訓(xùn)練樣本集，每個(gè)訓(xùn)練樣本集的大小約為原始數(shù)據(jù)集的三分之二。

2)對(duì)每個(gè)訓(xùn)練樣本集分別建立CART，并從全部M個(gè)特征中隨機(jī)抽取mtry個(gè)特征作為CART的特征變量，再根據(jù)基尼系數(shù)最小原則選擇最佳方案對(duì)節(jié)點(diǎn)進(jìn)行分裂。

3)每棵決策樹在生長(zhǎng)期間都按照步驟2)進(jìn)行，并保持mtry的大小不變；且每棵決策樹都要進(jìn)行完全的生長(zhǎng)，不對(duì)其進(jìn)行任何剪枝操作。

4)由前3步生成的ntree棵決策樹組成一個(gè)隨機(jī)森林。在分類問題中，RF模型對(duì)檢驗(yàn)集進(jìn)行投票并根據(jù)少數(shù)服從多數(shù)的原則進(jìn)行分類預(yù)測(cè)；而在回歸問題中，將每棵樹的結(jié)果等權(quán)相加，求得的均值作為最終預(yù)測(cè)值。

1.3 模型評(píng)價(jià)指標(biāo)

本文使用均方根誤差(RMSE)與平均絕對(duì)誤差(MAE)作為各個(gè)模型擬合與預(yù)測(cè)精度的評(píng)判標(biāo)準(zhǔn)。RMSE與MAE的定義為：

(2)

(3)

式中，g為預(yù)測(cè)值，h為原始值。RMSE和MAE的值越小，表示模型的精度越高；反之，則表示模型的精度越低。

通過引入皮爾森相關(guān)系數(shù)(Pearson correlation coefficient)和斯皮爾曼相關(guān)系數(shù)(Spearman correlation coefficient)來檢驗(yàn)各模型預(yù)測(cè)數(shù)據(jù)與原始數(shù)據(jù)的相關(guān)性。通常皮爾森相關(guān)系數(shù)對(duì)數(shù)據(jù)的要求較高，需要2個(gè)變量為線性關(guān)系且符合正態(tài)分布，而斯皮爾曼相關(guān)系數(shù)對(duì)數(shù)據(jù)沒有過多限制。根據(jù)表1可對(duì)數(shù)據(jù)的相關(guān)性強(qiáng)弱進(jìn)行判斷[17]，相關(guān)系數(shù)的絕對(duì)值越接近1，表示其相關(guān)程度越強(qiáng)。

表1 相關(guān)性強(qiáng)弱判定表

2 Prophet-RF組合模型的構(gòu)建

選用BJFS站2008～2013年的高程坐標(biāo)數(shù)據(jù)作為實(shí)驗(yàn)數(shù)據(jù)，數(shù)據(jù)來源于中國(guó)地震局GNSS數(shù)據(jù)產(chǎn)品服務(wù)平臺(tái)(http://www.cgps.ac.cn/),采樣間隔為1/365.25 a，最大采樣頻率為365.25 Hz。本文所有實(shí)驗(yàn)方案都保持相同的樣本集劃分：2008-01-01～2012-12-31的每日數(shù)據(jù)作為訓(xùn)練樣本集，2013-01-01～12-31的數(shù)據(jù)作為測(cè)試樣本集。圖2為由Prophet模型得到的預(yù)測(cè)結(jié)果，可以看出，Prophet模型能夠較好地?cái)M合原始數(shù)據(jù)中的周期性變化，且具有較好的抗粗差能力。

通過分析Prophet模型與RF模型的特性構(gòu)建不同的組合方案：方案1是基于Prophet模型自帶的分解功能，將數(shù)據(jù)分解為周期部分與非線性部分；方案2是假設(shè)在數(shù)據(jù)噪聲量級(jí)較大的情況下，使用Prophet模型對(duì)原始數(shù)據(jù)進(jìn)行擬合；方案3是假設(shè)在數(shù)據(jù)噪聲量級(jí)較小的情況下，使用RF模型對(duì)原始數(shù)據(jù)進(jìn)行擬合。

圖2 Prophet模型對(duì)BJFS站高程坐標(biāo)的預(yù)測(cè)Fig.2 Prediction of elevation coordinate of BJFS station by Prophet model

2.1 Prophet-RF組合模型方案1

Prophet模型對(duì)線性趨勢(shì)有較好的擬合與預(yù)測(cè)效果，對(duì)于非線性部分的預(yù)測(cè)效果較差。根據(jù)這一特性，使用Prophet模型自帶的分解功能將高程坐標(biāo)時(shí)間序列分解為趨勢(shì)項(xiàng)與周期項(xiàng)，其他部分看作殘差項(xiàng)。圖3中Prophet模型將BJFS站高程坐標(biāo)時(shí)間序列分解為3個(gè)主要成分：1)由邏輯回歸函數(shù)或分段線性函數(shù)擬合的趨勢(shì)項(xiàng)；2)基于傅立葉級(jí)數(shù)構(gòu)建的以a為周期的季節(jié)項(xiàng)；3)由虛擬變量表示的以周為周期的季節(jié)項(xiàng)[13]。假日項(xiàng)可單獨(dú)添加，當(dāng)數(shù)據(jù)中不存在此類情況時(shí)，采用默認(rèn)分解即可。圖4為Prophet-RF組合模型方案1的流程。

圖3 Prophet模型分解坐標(biāo)數(shù)據(jù)Fig.3 Prophet model decomposition coordinate data

圖4 Prophet-RF模型方案1流程Fig.4 Flow of scheme 1 of Prophet-RF model

Prophet-RF模型方案1的實(shí)現(xiàn)流程如下：

1)使用Prophet模型自帶的分解功能，將時(shí)間跨度為n天的坐標(biāo)時(shí)間序列Sn分解為趨勢(shì)項(xiàng)trendn、以周為周期的季節(jié)項(xiàng)weeklyn、以a為周期的季節(jié)項(xiàng)yearlyn及殘差項(xiàng)errorn。

2)趨勢(shì)項(xiàng)、以周為周期的季節(jié)項(xiàng)及以a為周期的季節(jié)項(xiàng)分別使用Prophet模型建立預(yù)測(cè)模型，殘差項(xiàng)則使用RF模型建立預(yù)測(cè)模型。將每個(gè)分量模型的前t天作為訓(xùn)練集，后k天作為測(cè)試集，且t+k=n。

3)將各分量模型得到的擬合結(jié)果與預(yù)測(cè)結(jié)果等權(quán)相加，得到Prophet-RF組合模型方案1的擬合值與預(yù)測(cè)值：

S擬合=trendt+weeklyt+yearlyt+errort

S預(yù)測(cè)=trendk+weeklyk+yearlyk+errork

(4)

2.2 Prophet-RF組合模型方案2

RF模型在某些噪聲較大的分類或回歸問題中可能會(huì)出現(xiàn)過擬合的情況，且數(shù)據(jù)中的噪聲量級(jí)往往難以判斷，故方案2在不明確數(shù)據(jù)噪聲量級(jí)的情況下，假設(shè)原始數(shù)據(jù)中噪聲量級(jí)較大。首先選用Prophet模型擬合坐標(biāo)原始數(shù)據(jù)，再對(duì)擬合數(shù)據(jù)建立RF預(yù)測(cè)模型，避免了在數(shù)據(jù)噪聲量級(jí)較大時(shí)RF模型造成數(shù)據(jù)過擬合的缺點(diǎn)。圖5為Prophet-RF組合模型方案2的流程。

圖5 Prophet-RF模型方案2流程Fig.5 Flow of scheme 2 of Prophet-RF model

Prophet-RF模型方案2的實(shí)現(xiàn)流程如下：

1)將時(shí)間跨度為n天的原始坐標(biāo)數(shù)據(jù)劃分為訓(xùn)練樣本集和測(cè)試樣本集，前t天作為訓(xùn)練樣本集，后k天作為測(cè)試樣本集，且t+k=n。

2)使用Prophet模型擬合訓(xùn)練樣本集，并將得到的擬合數(shù)據(jù)看作新的特征變量。

3)利用RF模型對(duì)Prophet模型產(chǎn)生的特征變量建立預(yù)測(cè)模型，得到最終的擬合數(shù)據(jù)與預(yù)測(cè)數(shù)據(jù)。

2.3 Prophet-RF組合模型方案3

RF模型在數(shù)據(jù)噪聲含量較少的情況下能夠保持較好的獨(dú)立性，即不易過度擬合數(shù)據(jù)，故方案3在假設(shè)原始數(shù)據(jù)噪聲量級(jí)較小的情況下，使用RF模型擬合原始數(shù)據(jù)，避免了Prophet模型對(duì)非線性部分預(yù)測(cè)能力較弱的缺陷，從而能夠充分發(fā)揮Prophet模型的預(yù)測(cè)優(yōu)勢(shì)。圖6為Prophet-RF組合模型方案3的流程。

圖6 Prophet-RF模型方案3流程Fig.6 Flow of scheme 3 of Prophet-RF model

Prophet-RF模型方案3實(shí)現(xiàn)流程如下：

1)將時(shí)間跨度為n天的坐標(biāo)原始數(shù)據(jù)劃分為訓(xùn)練樣本集和測(cè)試樣本集，前t天作為訓(xùn)練樣本集，后k天作為測(cè)試樣本集，且t+k=n。

2)使用RF模型擬合訓(xùn)練樣本集，并將得到的擬合數(shù)據(jù)看作新的特征變量。

3)利用Prophet模型對(duì)RF模型產(chǎn)生的特征變量建立預(yù)測(cè)模型，并得到最終的擬合數(shù)據(jù)與預(yù)測(cè)數(shù)據(jù)。

2.4 方案對(duì)比與分析

2.4.1 精度對(duì)比

本文選用BJFS站高程坐標(biāo)時(shí)間序列作為實(shí)驗(yàn)數(shù)據(jù)，將Prophet模型及Prophet-RF組合模型方案1、方案2、方案3得到的擬合值與預(yù)測(cè)值分別與原始數(shù)據(jù)進(jìn)行對(duì)比，同時(shí)引入經(jīng)典ARMA模型作為參照，得到各模型相應(yīng)的RMSE值與MAE值，結(jié)果見表2。

表2 不同模型的性能對(duì)比

由表2可知，ARMA模型的擬合精度最高，但預(yù)測(cè)精度最低，這是由于經(jīng)典ARMA模型容易出現(xiàn)過擬合現(xiàn)象, 且在長(zhǎng)期預(yù)測(cè)中缺乏穩(wěn)定性。而單一的Prophet模型也能得到較高的擬合精度，但預(yù)測(cè)精度卻不如組合模型，表明模型的預(yù)測(cè)能力與擬合效果沒有直接關(guān)系。Prophet-RF組合模型方案1、方案2、方案3預(yù)測(cè)數(shù)據(jù)的RMSE值較單一的Prophet模型分別提升0.7%、14.7%、16.1%，MAE值分別提升7.8%、23.4%、24.2%，其中方案1的預(yù)測(cè)精度提升較小，而方案2及方案3的預(yù)測(cè)精度都有較大提升。

為檢驗(yàn)各模型得到的擬合值、預(yù)測(cè)值與原始數(shù)據(jù)的關(guān)系，引入皮爾森相關(guān)系數(shù)和斯皮爾曼相關(guān)系數(shù)，由于ARMA模型在該數(shù)據(jù)中的適用性不強(qiáng)，故沒有分析其相關(guān)性。由表3可知，各模型的擬合值與原始數(shù)據(jù)具有強(qiáng)相關(guān)性，而預(yù)測(cè)值與原始數(shù)據(jù)具有中等強(qiáng)度相關(guān)性。通過對(duì)比分析表2和3的指標(biāo)認(rèn)為，組合模型方案3在BJFS站高程方向的預(yù)測(cè)精度最高、適用性最強(qiáng)。圖7為Prophet模型與Prophet-RF組合模型方案3預(yù)測(cè)誤差的對(duì)比，由預(yù)測(cè)值減去原始數(shù)據(jù)得到。

表3 不同模型的相關(guān)性對(duì)比

圖7 模型預(yù)測(cè)誤差對(duì)比Fig.7 Model prediction error comparison

2.4.2 組合模型分析

在方案1中，通過分解原始數(shù)據(jù)的組成成分，并根據(jù)Prophet模型和RF模型的特點(diǎn)，針對(duì)不同組成成分選用合適的預(yù)測(cè)模型，但該方案較單一的Prophet模型性能幾乎沒有提升。這可能是由于實(shí)測(cè)數(shù)據(jù)中往往存在由于各種誤差，如板塊運(yùn)動(dòng)、電離層延遲、觀測(cè)墩熱脹冷縮等造成的非線性變化，還有季節(jié)性海洋、陸地及大氣變化造成的地表水文負(fù)載等，這些都會(huì)對(duì)GNSS觀測(cè)數(shù)據(jù)產(chǎn)生負(fù)面影響[1]。所以方案1在分解數(shù)據(jù)時(shí)受到誤差的影響，無(wú)法較好地分離出各組成成分，從而導(dǎo)致模型精度提升較小。

方案2與方案3可看作一組對(duì)比實(shí)驗(yàn)，分別在數(shù)據(jù)噪聲量級(jí)較大和較小的假設(shè)下進(jìn)行。這2組方案較單一的Prophet模型預(yù)測(cè)精度都有較大提升，但方案3的預(yù)測(cè)精度更高，表明原始數(shù)據(jù)中的噪聲量級(jí)相對(duì)較小。方案2由于使用Prophet模型擬合原始數(shù)據(jù)時(shí)對(duì)非線性部分的擬合效果較差，再使用RF模型進(jìn)行預(yù)測(cè)時(shí)精度可能會(huì)受到Prophet模型擬合效果的影響。而方案3使用RF模型擬合原始數(shù)據(jù)，避免了過擬合并對(duì)非線性部分進(jìn)行了修正，再使用Prophet模型進(jìn)行預(yù)測(cè)時(shí)能夠充分發(fā)揮模型優(yōu)勢(shì)，從而提高預(yù)測(cè)精度。

3 結(jié) 語(yǔ)

高精度的預(yù)測(cè)模型有助于研究GNSS坐標(biāo)時(shí)間序列的長(zhǎng)期變化與規(guī)律。本文將Prophet模型應(yīng)用于GNSS坐標(biāo)數(shù)據(jù)中，該模型能夠分解坐標(biāo)時(shí)間序列中的組成成分，且具有較強(qiáng)的可解釋性；在Prophet模型的基礎(chǔ)上引入RF模型，以解決Prophet模型對(duì)隨機(jī)部分預(yù)測(cè)能力較弱的缺點(diǎn)。實(shí)驗(yàn)結(jié)果表明，Prophet-RF組合模型的預(yù)測(cè)精度相較于單一的Prophet模型，RMSE與MAE分別提升了16.1%與24.2%，并且對(duì)GNSS坐標(biāo)時(shí)間序列有更好的長(zhǎng)期預(yù)測(cè)能力。

由于數(shù)據(jù)往往受到噪聲及其他未知成分的影響，本文在進(jìn)行組合模型擬合時(shí)充分考慮了數(shù)據(jù)中可能存在的情況，并結(jié)合各模型的優(yōu)缺點(diǎn)設(shè)計(jì)了多種組合方案，以保證實(shí)驗(yàn)結(jié)果的可靠性，最終得到BJFS站高程方向的最優(yōu)組合方案。

本文對(duì)單個(gè)高程坐標(biāo)時(shí)間序列進(jìn)行了模型算法上的結(jié)合與研究，在今后的工作中應(yīng)當(dāng)結(jié)合更多的模型或引入更新的算法，以進(jìn)一步提高預(yù)測(cè)精度。