冉佳諾，邊家文，劉文平

(1.中國地質(zhì)大學(xué) 數(shù)學(xué)與物理學(xué)院，湖北 武漢 430074；2.湖北經(jīng)濟學(xué)院 信息管理與統(tǒng)計學(xué)院，湖北 武漢 430205)

1 引 言

全球?qū)Ш叫l(wèi)星系統(tǒng)(Global Navigation Satellite System, GNSS)坐標(biāo)時間序列在大地測量以及地球動力學(xué)領(lǐng)域的深入應(yīng)用, 使其對GNSS基準(zhǔn)站數(shù)據(jù)的測量精度有了更高要求。GNSS坐標(biāo)時間序列通常包含線性趨勢信號及非線性信號。線性趨勢信號反映了接收站點由構(gòu)造應(yīng)力導(dǎo)致的構(gòu)造運動，非線性信號反映了由環(huán)境載荷導(dǎo)致的季節(jié)性形變[1,2]。傳統(tǒng)的諧波模型僅用線性趨勢信號、非線性的常數(shù)振幅季節(jié)性信號以及噪聲項來線性模擬真實站點數(shù)據(jù)特點，但通常實際基準(zhǔn)站點會對氣候及環(huán)境負(fù)載等因素呈現(xiàn)非規(guī)則響應(yīng)，致使坐標(biāo)時間序列中的季節(jié)性信號振幅隨時間發(fā)生變化。同時，GNSS坐標(biāo)時間序列中還蘊含了大量非平穩(wěn)非線性噪聲，在中國境內(nèi)主要表現(xiàn)為白噪聲加閃爍噪聲[3]，噪聲信息與干凈的GNSS坐標(biāo)時間序列的相互交疊掩蓋了時間序列本身的運動規(guī)律，導(dǎo)致在測量GNSS坐標(biāo)時間序列時產(chǎn)生較大誤差。

針對GNSS坐標(biāo)時間序列中廣泛存在的噪聲干擾問題，目前已提出了多種解決算法。 Wu等人[4]使用多尺度小波分解(Wavelet Decomposition, WD)解構(gòu)原始信號，能夠從帶噪信號中提取出線性趨勢項以及年、半年季節(jié)項等有用信息，有效分離信號與噪聲，但這種方法由于需要事先確定小波基函數(shù)以及小波分解層數(shù)，導(dǎo)致WD在信號重構(gòu)時會無法避免地吸收部分噪聲；Chen等人[5]使用奇異譜分析(Singular Spectrum Anlysis,SSA)對GNSS坐標(biāo)時間序列重構(gòu)，并從速度不確定度、譜密度等方面分析SSA的重構(gòu)能力；Xu等人[6]發(fā)現(xiàn)GNSS坐標(biāo)時間序列中的季節(jié)性變化由環(huán)境載荷引起, 使用SSA過濾后的季節(jié)性信號仍存在有色噪聲驅(qū)動的噪聲信息。Klos等人[7]指出應(yīng)用SSA對GNSS坐標(biāo)時間序列降噪時，人為選擇的延滯窗口長度對重構(gòu)結(jié)果有較為明顯的影響。

隨著深度學(xué)習(xí)在各類研究任務(wù)中取得了出色成果，在信號處理領(lǐng)域相應(yīng)的研究也逐漸火熱，深度學(xué)習(xí)算法通過從樣本中學(xué)習(xí)較為簡易的特征，并將這些特征通過網(wǎng)絡(luò)逐層組合為更加復(fù)雜的特征，進而解決那些繁雜難以處理的問題。循環(huán)神經(jīng)網(wǎng)絡(luò)(Recurrent Neural Network, RNN)是一種常用于時序數(shù)據(jù)建模的神經(jīng)網(wǎng)絡(luò)，得益于網(wǎng)絡(luò)內(nèi)部神經(jīng)元循環(huán)連接的結(jié)構(gòu)，RNN能夠充分利用時間序列的歷史前后信息來對當(dāng)前時序結(jié)果進行預(yù)測。陶濤等人[8]通過LSTM(Long Short—term Memory)神經(jīng)網(wǎng)絡(luò)并利用自編碼器實現(xiàn)了GNSS坐標(biāo)時間序列中的異常點檢測,通過自編碼器對數(shù)據(jù)進行降維，以小特征空間捕獲原始時間序列中突出的特征，并利用這些特征進行時間序列重構(gòu)，結(jié)果表明異常值通常是那些不具代表性的特征，因此通過計算得到的原始時間序列和重構(gòu)時間序列的偏差可以作為異常值檢測的指標(biāo)[8]。Tung Kieu等[9]將LSTM和CNN(Convolutional Neural Network)結(jié)合，設(shè)計了LSTM-CNN的模型,對原始時間序列進行重建，獲得了精度較高的結(jié)果。Shen等人[10]在壓制引力波信號噪聲時使用了訓(xùn)練好的LSTM模型，能夠有效去除信號中的平穩(wěn)和非平穩(wěn)噪聲，比主成分分析以及字典學(xué)習(xí)等方法性能更加優(yōu)越; 仇悅[11]應(yīng)用循環(huán)神經(jīng)網(wǎng)絡(luò)對心電信號進行消噪，獲得了比傳統(tǒng)濾波方法更好的降噪結(jié)果。因此，本文提出基于Bi-LSTM(Bidirectional Long Short term Memory)的GNSS坐標(biāo)時間序列的降噪方法，利用Bi-LSTM強大的時序數(shù)據(jù)建模能力，對具有時變季節(jié)性信號的GNSS坐標(biāo)時間序列進行建模及重構(gòu)。對比分析Bi-LSTM在仿真坐標(biāo)時間序列數(shù)據(jù)和實際GNSS站點數(shù)據(jù)的噪聲抑制效果，并通過殘差分析方法驗證算法的有效性。

2 基于Bi-LSTM的GNSS坐標(biāo)時間序列重構(gòu)

2.1 GNSS坐標(biāo)時間序列模型

GNSS坐標(biāo)時間序列單站點單分量通常以如下方式建模[7]:

(1)

其中，x(t)、y(t) 、ε(t) 分別表示在t刻下的觀測時間序列、干凈信號、噪聲；x0是初始位置；v0是初始速度；aj(t)和bj(t)(j=1,2)代表年和半年周期振幅，分別由年和半年周期振幅的均值加上隨機變量構(gòu)成，用于刻畫信號中隨時間改變的季節(jié)項；f1和f2分別為年和半年季節(jié)信號頻率。由于模型(1)中季節(jié)信號的振幅是時變的，使得已有的重構(gòu)方法的重構(gòu)精度存在較大不足[12]。 因此，本文將使用基于循環(huán)神經(jīng)網(wǎng)絡(luò)的Bi-LSTM網(wǎng)絡(luò)，學(xué)習(xí)GNSS坐標(biāo)時間序列的時序特征，擬合季節(jié)項振幅中的時序變化，提升時變振幅季節(jié)項的GNSS坐標(biāo)時間序列的重構(gòu)精度。

2.2 循環(huán)神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)

循環(huán)神經(jīng)網(wǎng)絡(luò)(Recurrent Neural Network, RNN)是一種用于時序數(shù)據(jù)建模的神經(jīng)網(wǎng)絡(luò)[13]，得益于網(wǎng)絡(luò)內(nèi)部神經(jīng)元循環(huán)連接的結(jié)構(gòu)，RNN不僅能將當(dāng)前隱藏層狀態(tài)信息傳遞到輸出層進行結(jié)果輸出，還能將信息傳遞到下一時刻的隱藏層，使其能夠充分利用時間序列的歷史前后信息來輔助預(yù)測當(dāng)前時間步結(jié)果，其結(jié)構(gòu)如圖1所示。

圖1 RNN基本網(wǎng)絡(luò)結(jié)構(gòu)Fig.1 The structure of RNN

RNN的輸入x輸出y間的關(guān)系由如下數(shù)學(xué)表達式所示：

ht=f(Wh[xt,ht-1]+bh)

(2)

yt=g(Wyht+by)

(3)

其中，xt表示在時刻t的輸入時間序列；ht為隱藏狀態(tài)向量；yt為輸出時間序列；Wh和Wy為權(quán)值矩陣，其中Wh由Whx和Whh組成，Wy=Wyh構(gòu)成；bh和by為偏置向量；f和g為非線性激活函數(shù)。由式(2)可知，某一時刻ht的狀態(tài)由當(dāng)前輸入xt與t-1時刻隱藏層信息ht-1共同決定?；谶@樣的網(wǎng)絡(luò)結(jié)構(gòu)，RNN能夠充分利用歷史信息精確地學(xué)習(xí)輸入與輸出間的映射關(guān)系，在處理具有時序關(guān)系的數(shù)據(jù)時，比其他神經(jīng)網(wǎng)絡(luò)更具有先天優(yōu)勢。

2.3 Bi—LSTM網(wǎng)絡(luò)模型

2.3.1 單一LSTM記憶單元

傳統(tǒng)RNN網(wǎng)絡(luò)由于結(jié)構(gòu)存在固有缺陷，在參數(shù)更新時會存在梯度消失以及梯度爆炸的問題，導(dǎo)致長距離的歷史信息丟失，進一步造成網(wǎng)絡(luò)極難收斂，無法訓(xùn)練出理想的模型[14]。LSTM作為一種改進的循環(huán)神經(jīng)網(wǎng)絡(luò)，在原有網(wǎng)絡(luò)結(jié)構(gòu)的基礎(chǔ)上加入了細(xì)胞狀態(tài) (Cell State)的結(jié)構(gòu)來控制全局信息的傳輸，并通過遺忘門，輸入門，輸出門三種門控單元控制細(xì)胞狀態(tài)信息值的更新。LSTM在極大程度上緩解了傳統(tǒng)RNN模型存在的長期依賴問題，減少了長距離歷史信息的丟失，輸出的預(yù)測結(jié)果更準(zhǔn)確。LSTM的結(jié)構(gòu)示意圖如圖2所示。

圖2 單個LSTM記憶單元結(jié)構(gòu) Fig.2 The structure of a single LSTM memory cell

圖3 Bi-LSTM重構(gòu)模型Fig.3 Bi-LSTM reconstruction model

LSTM的狀態(tài)關(guān)系表達式如下[15]:

運動式治理的第二種動員技術(shù)是會議和文件。盡管中央八項規(guī)定對于會議和文件都作了簡化要求，但這只是原則性規(guī)定，針對的是一般性的會議和文件。對于地方政府而言，面對特殊情況、政治任務(wù)，會議和文件不但不能減少，反而要增加。因為只有這樣，才能彰顯地方政府對該項治理運動的高度重視和扎實推進。

其中，σ和tanh分別為Sigmoid激活函數(shù)和雙曲正切激活函數(shù)；ft、it、ot分別為t時刻遺忘門、輸入門、輸出門的輸出信息；ct表示t時刻更新后的細(xì)胞狀態(tài)信息；b為偏置；W為權(quán)值矩陣。

本文采用Bi-LSTM對GNSS坐標(biāo)時間序列進行建模[16,17]，以t時刻為例，其輸出結(jié)果yt由前后時刻隱藏層狀態(tài)ht-1和ht+1以及當(dāng)前輸入時間步xt共同決定，相比單向LSTM只利用了歷史信息，Bi-LSTM更全面地考慮了時序的上下文關(guān)系，比單向LSTM結(jié)果更精確。

2.3.2 基于Bi-LSTM的信號重構(gòu)

基于Bi-LSTM的GNSS坐標(biāo)時間序列降噪網(wǎng)絡(luò)模型結(jié)構(gòu)如圖3所示，該網(wǎng)絡(luò)主要由輸入層、隱藏層、輸出層三個模塊組成。當(dāng)網(wǎng)絡(luò)處于訓(xùn)練階段時，將2.1節(jié)中所述的仿真GNSS坐標(biāo)時間序列作為輸入層的訓(xùn)練數(shù)據(jù)輸入。隱藏層由256個LSTM單元組成，用于學(xué)習(xí)輸入數(shù)據(jù)到輸出結(jié)果間的非線性映射，保存時序數(shù)據(jù)中含有的有用信息；隱藏層數(shù)、隱藏層節(jié)點數(shù)等超參數(shù)可根據(jù)實際訓(xùn)練結(jié)果合理地進行增加或減少。

(10)

所得損失結(jié)果經(jīng)由自適應(yīng)動量估計優(yōu)化器(Adaptive Moment Estimation，ADAM)[18]進行梯度優(yōu)化調(diào)整，并將這些信息輸入到隱藏層進行權(quán)值更新，再將更新后的模型用于下一步迭代訓(xùn)練，以降低損失函數(shù)值，直至結(jié)果收斂?；贐i-LSTM的GNSS坐標(biāo)時間序列重構(gòu)及評估的具體流程如下：

2)使用步驟1)生成的數(shù)據(jù)集訓(xùn)練神經(jīng)網(wǎng)絡(luò)模型，根據(jù)實驗結(jié)果調(diào)整網(wǎng)絡(luò)參數(shù)。

3)通過步驟2)訓(xùn)練好的網(wǎng)絡(luò)模型對觀測GNSS坐標(biāo)時間序列重構(gòu)，得到重構(gòu)后的坐標(biāo)時間序列。

4)計算觀測坐標(biāo)時間序列與重構(gòu)的坐標(biāo)時間序列信號間的殘差，并計算重構(gòu)時間序列與真實坐標(biāo)時間序列的標(biāo)準(zhǔn)差。

5)對步驟4)中得到的殘差數(shù)據(jù)進行殘差分析，評估指標(biāo)結(jié)果。

3 實驗結(jié)果與分析

本節(jié)設(shè)置了仿真GNSS時序數(shù)據(jù)和實際站點數(shù)據(jù)兩種實驗，分別從兩種策略評估Bi-LSTM模型在理想情況和實際情況下的GNSS坐標(biāo)時間序列的降噪性能。同時，由于Bi-LSTM模型訓(xùn)練需要一定規(guī)模的無噪、連續(xù)的數(shù)據(jù)樣本作為標(biāo)簽數(shù)據(jù)，而現(xiàn)有公用的數(shù)據(jù)庫中符合這類要求的實際數(shù)據(jù)較少，因此有必要通過仿真數(shù)據(jù)訓(xùn)練網(wǎng)絡(luò)擴充數(shù)據(jù)集，參考第2.1節(jié)中所提到的模型及理論依據(jù)，本節(jié)設(shè)計了相應(yīng)的合成仿真GNSS數(shù)據(jù)用于模型訓(xùn)練。模型訓(xùn)練完畢后，將訓(xùn)練好的Bi-LSTM網(wǎng)絡(luò)模型應(yīng)用于實際GNSS站點數(shù)據(jù)重構(gòu)，定性地分析該深度學(xué)習(xí)方法在GNSS坐標(biāo)時間序列處理任務(wù)的測試結(jié)果。

3.1 仿真實驗

3.1.1 GNSS仿真數(shù)據(jù)集

He等人[19]分析了全球671個IGS站點，90%以上GNSS坐標(biāo)時間序列Up方向的年周期項振幅在0～8 mm之間，半年周期項振幅在0～2 mm波動；而North 和 East方向的年周期項振幅范圍在0～2 mm之間，其半年周期項振幅小于1 mm；時間序列中包含的噪聲類型為白噪聲及閃爍噪聲, 閃爍噪聲振幅大小在7～20 mm之間。為了估計GNSS坐標(biāo)時間序列年和半年信號的不同時刻的季節(jié)變化程度，Klos等人將位于亞洲基準(zhǔn)站點的長期GNSS坐標(biāo)時間序列以三年長度進行分段，并計算得出連續(xù)段之間Up方向上的標(biāo)準(zhǔn)差均值為0.8 mm，North、East分量上相應(yīng)標(biāo)準(zhǔn)差均值為0.4 mm[7]。

基于此，通過GNSS坐標(biāo)時間序列模型模擬產(chǎn)生了50 000條仿真數(shù)據(jù)，每條仿真數(shù)據(jù)中的季節(jié)項信號的年周期振幅在0～8 mm之間，半年周期振幅在0～2 mm之間；同時在年周期信號振幅中加上服從均值為0、標(biāo)準(zhǔn)差為1的隨機變量，半年周期信號振幅中加上均值為0、標(biāo)準(zhǔn)差為0.5的隨機變量，以模擬GNSS坐標(biāo)時間序列的時序變化；并在仿真信號中隨機添加振幅在0～20 mm/year-k/4的閃爍噪聲，以便能夠自適應(yīng)地完成不同振幅的GNSS坐標(biāo)時間序列降噪任務(wù)。

3.1.2 仿真實驗結(jié)果分析

本文采用以下方式定量分析LSTM的重構(gòu)結(jié)果：1) 重構(gòu)后的時間序列數(shù)據(jù)與原始干凈數(shù)據(jù)的標(biāo)準(zhǔn)差(Misfit)，即算法重構(gòu)后的信號與原始信號的差值計算標(biāo)準(zhǔn)差。Misfit越小，說明降噪效果越好，與原信號擬合程度越高。2)殘差譜指數(shù)k及殘差振幅A，帶噪時間序列與算法降噪后的時間序列差值即為殘差。若殘差譜指數(shù)及殘差振幅越接近所加入的噪聲的譜指數(shù)及振幅，則提取出的噪聲越徹底[20,21]。3)由殘差振幅及殘差譜指數(shù)計算獲得的速度不確定度[22]。

圖4 各算法在閃爍噪聲振幅為10 mm/year1/4的季節(jié)信號重構(gòu)Fig.4 Seasonal signal reconstruction of each algorithm at noise amplitude of 10 mm/ year1/4

圖4為SSA, WD, Bi-LSTM三種方法在噪聲振幅為10 mm/year1/4去趨勢后的仿真GNSS坐標(biāo)時間序列季節(jié)性信號降噪結(jié)果圖。由圖4可知，WD吸收了部分有色噪聲，重構(gòu)精度明顯受到噪聲信息干擾，重構(gòu)出的信號受噪聲信息產(chǎn)生更為明顯的波動；Bi-LSTM能夠準(zhǔn)確地重構(gòu)出原時間序列。從表1的量化分析結(jié)果得知，基于Bi-LSTM方法的Misfit值最小，說明該算法降噪后的信號與原時間序列擬合程度最好；并且其殘差譜指數(shù)k和殘差振幅A等均更接近真實值，且基于k和A所得的速度不確定度也更符合實際結(jié)果。同時，功率譜密度圖(Power Spectral Density, PSD) 從能量的角度驗證了前面的結(jié)果，如圖5所示，其中橫坐標(biāo)代表頻率，縱坐標(biāo)表示功率譜密度的幅值大??；黃色線條代表了原觀測時間序列功率譜密度；天藍色、綠色和紅色線條分別表示W(wǎng)D, SSA, Bi-LSTM三種算法從觀測數(shù)據(jù)中分離出的噪聲殘差功率譜密度；紫色線條表示噪聲的功率譜密度。從功率譜密度圖上來看，與其他方法相比,Bi-LSTM在各個頻率上的功率分量大小與噪聲的功率分量大小更加接近，說明在噪聲振幅為10 mm/year1/4時，相比SSA和WD所提取的噪聲殘差, Bi-LSTM所提取的殘差與“真實”的噪聲更接近，即Bi-LSTM能夠更徹底地分離干凈時間序列和噪聲信息?；诖耍瑥姆抡鎸嶒灲Y(jié)果得知，Bi-LSTM能從涵蓋不同振幅的大規(guī)模數(shù)據(jù)中學(xué)習(xí)到時間序列的季節(jié)性變化，如圖4所示，Bi-LSTM能夠自適應(yīng)地捕捉不同時刻的季節(jié)變化特征，更精確地估計時變振幅的季節(jié)信號。

表1 閃爍噪聲振幅 10 mm/year1/4的各項指標(biāo)值

圖6 BJFS站點North方向數(shù)據(jù)各方法信號重構(gòu)Fig.6 Reconstruction results of each method of North direction at BJFS station

3.2 實際GNSS坐標(biāo)時間序列分析

本文使用北京房山站點(BJFS)2004年2月1日至2020年7月6日 的North、East、Up三個方向的GNSS數(shù)據(jù)來驗證基于Bi-LSTM的深度學(xué)習(xí)方法對于GNSS坐標(biāo)時間序列的降噪效果，該GNSS數(shù)據(jù)可由IGS全球數(shù)據(jù)中心SOPAC(http://sopac-ftp.ucsd.edu/pub/timeseries /)下載。預(yù)先對數(shù)據(jù)進行預(yù)處理剔出異常值，并對三個分量上的原觀測時間序列缺失部分進行線性插值補全，再使用Bi-LSTM、SSA、WD三種方法對北京房山站3個分量上的數(shù)據(jù)進行降噪處理。由于無噪的實際GNSS坐標(biāo)時間序列未知，因此在分析各類估計方法時不計算Misfit值。本文將Bi-LSTM處理實際站點數(shù)據(jù)的降噪結(jié)果與另外兩種廣泛應(yīng)用于GNSS坐標(biāo)時間序列重構(gòu)方法(SSA、WD)的重構(gòu)結(jié)果計算相關(guān)性，以證明Bi-LSTM處理實際數(shù)據(jù)的有效性[23]。房山站點Up方向各方法降噪后的GNSS數(shù)據(jù)結(jié)果相關(guān)系數(shù)見表2。圖6～圖8為各方法對北京房山站North、East、Up三個分量降噪后的GNSS坐標(biāo)時間序列，其中 North、East方向由線性趨勢項占主導(dǎo)地位，其季節(jié)項信號所帶來季節(jié)性變化并不明顯；而在Up方向上能夠看見明顯的趨勢變化及季節(jié)性變化。圖9為去趨勢后的Up方向各方法信號提取結(jié)果。由圖9可知，WD容易受到噪聲的干擾，Bi-LSTM對噪聲更魯棒，且能夠良好地擬合時變的季節(jié)性信號。

表2 房山站點Up方向各方法降噪后的GNSS數(shù)據(jù)結(jié)果相關(guān)系數(shù)

圖7 BJFS站點East方向數(shù)據(jù)各方法信號重構(gòu)Fig.7 Reconstruction results of each method of East direction at BJFS station

圖8 BJFS站點Up方向數(shù)據(jù)各方法信號重構(gòu)Fig.8 Reconstruction results of each method of Up direction at BJFS station

圖9 去趨勢后的Up方向各方法信號提取Fig.9 The detrended signal extracted by each method in Up direction

4 結(jié) 論

本文探索性地引入了深度學(xué)習(xí)領(lǐng)域中廣泛用于時序數(shù)據(jù)處理的Bi-LSTM循環(huán)神經(jīng)網(wǎng)絡(luò)，對GNSS坐標(biāo)時間序列普遍含有的非平穩(wěn)噪聲進行降噪。仿真和實際GNSS站點數(shù)據(jù)結(jié)果表明，本文提出的基于Bi-LSTM的降噪方法相比傳統(tǒng)方法能夠更有效地擬合GNSS坐標(biāo)時間序列中所蘊含的時變特征，從而在處理帶噪時變振幅信號時獲得更優(yōu)秀的降噪結(jié)果。同時，該Bi-LSTM模型由帶有不同時變振幅及不同噪聲振幅的GNSS坐標(biāo)時間序列訓(xùn)練而成，因此能夠自適應(yīng)地應(yīng)對實際站點中帶有各種噪聲振幅的時間序列，即使在噪聲水平較高時也能夠取得良好的降噪結(jié)果。因此基于Bi-LSTM循環(huán)神經(jīng)網(wǎng)絡(luò)的GNSS坐標(biāo)時間序列重構(gòu)滿足了后續(xù)大地測量高測量精度的要求。