三種GNSS 高程時(shí)序降噪方法的效果對(duì)比分析

范小猛，胡川，張重陽，李成洪

(重慶交通大學(xué) 土木工程學(xué)院,重慶 400074)

0 引言

在全球衛(wèi)星導(dǎo)航系統(tǒng)(GNSS)坐標(biāo)時(shí)間序列中，不僅包含真實(shí)信號(hào)，還含有白噪聲(WN)、閃爍噪聲(FN)、隨機(jī)漫步噪聲(RWN)等多種有色噪聲.合理的分離出信號(hào)和噪聲，有助于獲取更準(zhǔn)確的測(cè)站運(yùn)動(dòng)參數(shù)估值，對(duì)于精化速度場(chǎng)、建立地球框架以及揭示構(gòu)造變形運(yùn)動(dòng)等有著重要的意義[1-2].

經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)、整體經(jīng)驗(yàn)?zāi)B(tài)分解(EEMD)和小波降噪是三種較常用的噪聲分離方法.HUANG等[3]在1998 年討論了EMD 在時(shí)間序列分析中的應(yīng)用問題.張雙成等[4]將基于EMD 的時(shí)間序列分析法引入GNSS 坐標(biāo)序列分析，對(duì)中國境內(nèi)及周邊9 個(gè)國際GNSS 服務(wù)(IGS)站三個(gè)坐標(biāo)分量上的序列數(shù)據(jù)進(jìn)行處理，有效地分離出線性趨勢(shì)項(xiàng)、周期項(xiàng)和噪聲項(xiàng)，證明了EMD 法降噪的有效性.隨著研究的不斷深入，研究人員發(fā)現(xiàn)EMD 存在模式混疊的現(xiàn)象[5-6]，為了改善這個(gè)問題，WU 等[7]在EMD 的基礎(chǔ)上提出了EEMD.張恒璟等[8]采用EMD 和EEMD分別對(duì)BJFS 站和URUM 站10 a 左右的高程數(shù)據(jù)進(jìn)行處理，分離序列趨勢(shì)項(xiàng)和周期項(xiàng)運(yùn)動(dòng)，發(fā)現(xiàn)EEMD可以明顯減弱模式混疊現(xiàn)象，但不能完全消除模式混疊，同時(shí)EEMD分解時(shí)所加入的WN 會(huì)導(dǎo)致其分解結(jié)果不唯一.小波分析具有多分辨率、低熵性、去相關(guān)性等特點(diǎn)，得到了GNSS 時(shí)間序列分析研究人員的重視[9-11].魏冠軍等[12]采用小波變換閾值降噪對(duì)高程分量的坐標(biāo)時(shí)序進(jìn)行處理，有效減弱了噪聲信息，證明了該方法降噪的可行性，但樣本序列數(shù)據(jù)時(shí)間長度在3 a 左右，沒有涉及到更長時(shí)間跨度時(shí)序的降噪分析.馬俊等[13]提出了小波包系數(shù)信息熵法，對(duì)KMIN 站近5 a 的三個(gè)坐標(biāo)分量上的時(shí)序進(jìn)行了降噪處理，得到了僅受WN 影響的時(shí)間序列，消除了有色噪聲的影響，但同樣沒有對(duì)更長時(shí)間跨度的坐標(biāo)時(shí)間序列進(jìn)行分析.

目前，已經(jīng)有眾多學(xué)者對(duì)EMD、EEMD 和小波降噪有了深入研究，但少有文獻(xiàn)考慮到數(shù)據(jù)長度可能會(huì)對(duì)降噪效果產(chǎn)生的影響.隨著GNSS基準(zhǔn)站的不間斷觀測(cè)，相關(guān)機(jī)構(gòu)已積累了20 多年的坐標(biāo)時(shí)間序列數(shù)據(jù)，為研究上述方法在不同時(shí)序長度下的降噪效果提供了可能.針對(duì)上述問題，本文通過對(duì)比分析EMD、EEMD 和小波降噪在不同GNSS 高程坐標(biāo)時(shí)序長度下的降噪效果，采用信噪比(SNR)、相關(guān)系數(shù)和均方根誤差(RMSE)定量評(píng)價(jià)它們的降噪性能，為相關(guān)研究提供一定參考.

1 降噪方法

1.1 經(jīng)驗(yàn)?zāi)B(tài)分解

EMD 是一種基于非線性、非平穩(wěn)數(shù)據(jù)的自適應(yīng)時(shí)頻信號(hào)分解方法，在不同尺度上分析數(shù)據(jù)，提取有用信息[3].其基本思想是通過篩選，把原始信號(hào)X(t)分解成m個(gè)不同時(shí)間尺度的本征模態(tài)分量(IMF)和殘差項(xiàng)rn，分解流程如下：

1)尋找原始信號(hào)中所有的局部極大、極小值點(diǎn)，利用三次樣條函數(shù)分別擬合上述點(diǎn)，得到關(guān)于時(shí)間的上、下包絡(luò)線，計(jì)算兩條包絡(luò)線的平均包絡(luò)線m1(t)，用原始信號(hào)X(t)減去平均包絡(luò)得到新的待分解信號(hào)h1(t)

2)判斷h1(t)是否滿足以下條件：一是待分解的信號(hào)中的極值點(diǎn)的個(gè)數(shù)與過零點(diǎn)的個(gè)數(shù)應(yīng)相等或最多相差1；二是在任意點(diǎn)上，上下包絡(luò)線的均值為0，即上下包絡(luò)線應(yīng)對(duì)稱.如果滿足上述條件，將h1(t)作為分解得到的第一個(gè)分量 I MF1(t)，然后進(jìn)行第3)步操作，否則將h1(t)作為新的原始信號(hào)重復(fù)第1)步操作.

3)計(jì)算剩余信號(hào)r1=X(t)-IMF1(t)，然后將r1作為新的待分解信號(hào).

4)判斷r1是否滿足以下條件：r1為單調(diào)變化的殘差序列，且不能繼續(xù)分解得到IMF 分量.如果滿足，EMD 分解過程結(jié)束，殘差項(xiàng)rm=r1，否則重復(fù)流程1)～3)操作.原始信號(hào)經(jīng)過m次分解可以表示為

為了獲取不含噪聲的GNSS 坐標(biāo)時(shí)間序列，需要對(duì)上式中分解出的IMF 分量進(jìn)行重構(gòu).首先，計(jì)算各分量與原始序列間的相關(guān)系數(shù)；然后識(shí)別第一個(gè)取局部極小值的相關(guān)系數(shù)對(duì)應(yīng)的分量IMFk，將其作為信號(hào)與噪聲之間的分界層；最后，將分界層與之前的IMF 分量層相加作為高頻噪聲，然后用原始信號(hào)減去重構(gòu)高頻噪聲得到降噪后信號(hào)[3].

1.2 整體經(jīng)驗(yàn)?zāi)B(tài)分解

EMD 的缺點(diǎn)之一是模式混疊：信號(hào)被EMD 分解為m個(gè)不同的IMF 分量后，不同頻率的信號(hào)會(huì)存在于同一個(gè)IMF 分量中，即一個(gè)IMF 分量包含著不同尺度的信號(hào).EEMD 通過在分解前向原始信號(hào)中添加白噪聲w(t)，生成待分解信號(hào)，改善模式混疊問題.

式中，下標(biāo)i表示第i次分解.其余分解步驟與EMD 相同，此處不再贅述.將多次分解得到的尺度相同的IMF 分量整體取平均，作為當(dāng)前尺度的IMF 分量，從而得到最終的一系列不同尺度的IMF 分量數(shù)據(jù)，按照相關(guān)系數(shù)準(zhǔn)則重構(gòu)信號(hào)即可得到降噪之后的GNSS坐標(biāo)時(shí)間序列.本文根據(jù)文獻(xiàn)[14]實(shí)驗(yàn)結(jié)果，將加入的WN 與原始信號(hào)的標(biāo)準(zhǔn)差之比設(shè)為0.1，整體平均次數(shù)設(shè)置為100 次.

1.3 小波降噪

多分辨率分析是小波分析的重要概念，它是利用特定的小波基函數(shù)的伸縮變換，將GNSS 坐標(biāo)時(shí)間序列在不同分辨率上進(jìn)行分解，得到能夠反映時(shí)間序列概貌的低頻部分以及能夠反映時(shí)間序列細(xì)節(jié)的高頻部分，進(jìn)而對(duì)不同分辨率上的信號(hào)進(jìn)行分析[12].小波閾值降噪是在多分辨率分析的基礎(chǔ)上，通過設(shè)置一定的閾值對(duì)分解之后的各個(gè)尺度上的信號(hào)進(jìn)行處理，達(dá)到分離噪聲的目的.小波閾值降噪的具體步驟如下：

1)小波分解.選取合適的小波基函數(shù)，對(duì)GNSS坐標(biāo)時(shí)間序列進(jìn)行小波分解，得到不同尺度上的小波系數(shù).

2)閾值降噪.對(duì)分解得到的各個(gè)尺度上的信號(hào)進(jìn)行非線性閾值降噪處理，剔除小波系數(shù)小于閾值的小波項(xiàng).

3)重構(gòu)信號(hào).將剔除小波項(xiàng)之后的各尺度上的信號(hào)重構(gòu)得到降噪后的GNSS 坐標(biāo)時(shí)間序列.

經(jīng)過實(shí)驗(yàn)對(duì)比，本文采用Coif3 小波對(duì)原始數(shù)據(jù)進(jìn)行分解，根據(jù)軟閾值規(guī)則對(duì)數(shù)據(jù)進(jìn)行降噪處理.

2 實(shí)驗(yàn)數(shù)據(jù)和方案

2.1 實(shí)驗(yàn)數(shù)據(jù)

由于我國境內(nèi)的IGS 站時(shí)間序列在水平方向上主要以線性為主，而在高程方向上表現(xiàn)出復(fù)雜的非線性運(yùn)動(dòng)[15-16]，因此本文以高程方向上的坐標(biāo)時(shí)序?yàn)槔M(jìn)行分析.為保證實(shí)驗(yàn)所需時(shí)序長度，選取我國境內(nèi)BJFS、KMIN、LHAS、SHAO、URUM、WUHN 共6 個(gè)IGS 站點(diǎn)所對(duì)應(yīng)的高程方向單日解坐標(biāo)時(shí)序，各站點(diǎn)完整數(shù)據(jù)的起止歷元均為1999.17—2020.39，截取其中2000—2019 年的觀測(cè)數(shù)據(jù)作為本次實(shí)驗(yàn)的原始數(shù)據(jù).

采用線性擬合去除原始數(shù)據(jù)中的趨勢(shì)項(xiàng)，并且根據(jù)3σ 準(zhǔn)則剔除孤立值.另外，為了避免由于數(shù)據(jù)缺失造成模式混疊現(xiàn)象，以及滿足小波分析要求數(shù)據(jù)具有均勻取樣和零均值的特性，利用線性插值對(duì)數(shù)據(jù)缺失的部分進(jìn)行插補(bǔ)，得到適用于降噪分析的連續(xù)數(shù)據(jù).

2.2 實(shí)驗(yàn)方案

實(shí)驗(yàn)方案共分為三種：

方案一：將每個(gè)站點(diǎn)共20 a 的坐標(biāo)時(shí)間序列數(shù)據(jù)依次分為4 段，分別為2000—2004、2005—2009、2010—2014、2015—2019 時(shí)間段，得到24 個(gè)時(shí)序長度為5 a 的樣本序列；

方案二：將每個(gè)站點(diǎn)共20 a 的坐標(biāo)時(shí)間序列數(shù)據(jù)平均分為2 段，分別為2000—2009、2010—2019，得到12 個(gè)時(shí)序長度為10 a 的樣本序列；

方案三：將每個(gè)站點(diǎn)2000—2019 共20 a 的坐標(biāo)時(shí)間序列數(shù)據(jù)分別整體作為一段，得到6 個(gè)時(shí)序長度為20 a 的樣本序列.

降噪效果通過降噪后信號(hào)與原始信號(hào)之間的SNR、相關(guān)系數(shù)以及RMSE 來判斷.其計(jì)算公式分別為：

式中：X(t)為原始序列；為降噪后的真實(shí)信號(hào)；N為采樣點(diǎn)個(gè)數(shù).SNR 體現(xiàn)的是噪聲信號(hào)在信號(hào)中所占的比重，其值越大，表示噪聲含量越低；相關(guān)系數(shù)越接近1，表示兩者的相似度越高，降噪效果越好；RMSE 表示的是兩者之間的偏差程度，其值越小代表降噪效果越好.

3 實(shí)驗(yàn)結(jié)果與分析

3.1 降噪結(jié)果

經(jīng)過計(jì)算，三種方法在不同樣本序列長度下所分解出的分量個(gè)數(shù)相同，從5 a 到20 a 分別為10、11、12 個(gè).對(duì)于EMD 和EEMD 而言，根據(jù)各IMF 分量與原樣本序列的相關(guān)系數(shù)得到的信號(hào)與噪聲間分界層不同：EMD 法得到的不同樣本序列分界層數(shù)分別為第5、4、5 層，而EEMD 法得到的不同樣本序列分界層數(shù)分別為第2、5、5 層，根據(jù)上述分界層對(duì)樣本序列進(jìn)行降噪處理.

圖1 和圖2 描述了不同IGS 站點(diǎn)高程方向上經(jīng)過降噪的坐標(biāo)時(shí)間序列.從圖中可以發(fā)現(xiàn)三種方法都可以有效的對(duì)坐標(biāo)時(shí)間序列數(shù)據(jù)進(jìn)行降噪，降噪后信號(hào)表現(xiàn)出明顯的周期性變化，曲線變化也較為平滑.其中，圖1 中BJFS 站的5 a-2 樣本序列以及KMIN 站的5 a-4 樣本序列的光滑度明顯優(yōu)于其余相同時(shí)間跨度的樣本序列，其原因可能是在這兩個(gè)樣本數(shù)據(jù)中缺失歷元數(shù)較少，數(shù)據(jù)較為完整.例如BJFS 站的5 a-2樣本序列的缺失歷元數(shù)為54，僅占樣本總數(shù)的3%，遠(yuǎn)小于其它曲線光滑度較低的樣本序列缺失歷元數(shù)，說明時(shí)間序列中缺失數(shù)據(jù)的存在會(huì)對(duì)降噪效果產(chǎn)生一定影響.

圖1 BJFS 站和KMIN 站降噪結(jié)果對(duì)比

圖2 SHAO 站和URUM 站降噪結(jié)果對(duì)比

圖2 中SHAO 和URUM 兩個(gè)站點(diǎn)，EMD 降噪后的曲線在個(gè)別時(shí)間段呈現(xiàn)出明顯的異常波動(dòng).其中SHAO 站的5 a-2 樣本序列數(shù)據(jù)在2007 年前后波動(dòng)最大，最大達(dá)到了80 mm，5 a-4 樣本序列數(shù)據(jù)中異常波動(dòng)較多，但波動(dòng)幅度較小，其余站點(diǎn)也存在相似的異常波動(dòng)情況.這種異常波動(dòng)現(xiàn)象可能是原始時(shí)間序列數(shù)據(jù)中存在長時(shí)間的連續(xù)中斷，如SHAO 站連續(xù)中斷的歷元數(shù)為411，缺值率達(dá)到了22%，即便是后續(xù)對(duì)缺失處進(jìn)行了插補(bǔ)，依然會(huì)對(duì)EMD 降噪結(jié)果造成偏差.

此外，在上述出現(xiàn)異常波動(dòng)的時(shí)間點(diǎn)，EEMD 和小波降噪得到的序列趨勢(shì)相似，振幅為3～5 mm，趨勢(shì)較平緩，說明即使數(shù)據(jù)存在長時(shí)間的連續(xù)缺失，這兩種方法也可以有效減弱插值的影響，削弱信號(hào)中存在的噪聲.

3.2 降噪結(jié)果對(duì)比分析

為了進(jìn)一步說明EMD、EEMD 和小波降噪的效果，分別計(jì)算不同樣本長度下各方法降噪結(jié)果的SNR、相關(guān)系數(shù)和RMSE，其中樣本長度為5 a 和10 a 時(shí)，將各指標(biāo)的平均值作為最終結(jié)果，如表1 所示.

由表1 可知：

表1 各站點(diǎn)評(píng)價(jià)參數(shù)統(tǒng)計(jì) m

1) EMD 受初始數(shù)據(jù)質(zhì)量的影響較大.在數(shù)據(jù)較完整的站點(diǎn)，EMD 降噪序列各評(píng)價(jià)指標(biāo)略差或持平于EEMD 和小波降噪結(jié)果；而在缺失比例大的站點(diǎn)，各評(píng)價(jià)指標(biāo)值明顯較差.

2)相比于EMD，EEMD 不僅改善了模式混疊問題，還在一定程度上提高了EMD 分解的精度，即便是在缺值數(shù)較多的SHAO 站，EEMD 降噪序列的相關(guān)系數(shù)值大部分在0.7 以上，降噪效果穩(wěn)定.

3)在5 a 和10 a 的樣本時(shí)間序列數(shù)據(jù)中，小波降噪序列的評(píng)價(jià)參數(shù)值明顯優(yōu)于EMD 和EEMD；當(dāng)時(shí)間跨度為20 a 時(shí)，EEMD 降噪序列的評(píng)價(jià)參數(shù)與小波相近.這說明在這三種方法中，小波降噪的適用性更好，而EEMD 更適用于對(duì)時(shí)間跨度大的坐標(biāo)時(shí)間序列進(jìn)行降噪分析.

4)大部分站點(diǎn)的EEMD 降噪序列和小波降噪序列隨著時(shí)間跨度的增加，三個(gè)評(píng)價(jià)參數(shù)值均有所增加，難以根據(jù)某一個(gè)參數(shù)準(zhǔn)確說明其降噪性能的強(qiáng)弱.

為了進(jìn)一步探究EMD、EEMD 和小波降噪對(duì)有色噪聲的剔除能力，分別計(jì)算三種方法降噪序列的譜指數(shù)，譜指數(shù)越接近零，在一定程度上說明有色噪聲的振幅越小，對(duì)真實(shí)信號(hào)的影響也越小[13]，結(jié)果如圖3 所示.

圖3 三種方法得到的降噪序列的譜指數(shù)

由圖3 可以得出以下結(jié)論：1)EEMD 降噪后時(shí)序的譜指數(shù)大于EMD 降噪序列的譜指數(shù)，剔除有色噪聲的能力更好，但是否因EEMD 方法所加入WN 產(chǎn)生的影響還有待進(jìn)一步研究；2)小波分析降噪在不同時(shí)間跨度的時(shí)序數(shù)據(jù)下都可以更好的抑制有色噪聲，使得降噪后的信號(hào)更加可靠.

4 結(jié)論

本文通過設(shè)置不同樣本長度的降噪實(shí)驗(yàn)，對(duì)EMD、EEMD 以及小波分析三種方法的降噪效果進(jìn)行對(duì)比，并利用SNR、RMSE 和相關(guān)系數(shù)定量分析了其降噪性能.實(shí)驗(yàn)表明：

1)EMD 易受坐標(biāo)時(shí)間序列自身質(zhì)量的影響，降噪效果相對(duì)較差；EEMD 和小波降噪在時(shí)間序列存在長時(shí)間中斷，質(zhì)量較差時(shí)，依然可以得到較好的降噪結(jié)果.

2)隨著時(shí)序樣本長度的增加，三種方法的降噪性能均有提升.其中，小波降噪的適用性最廣泛，在不同的樣本數(shù)據(jù)中，降噪性能均優(yōu)于EMD 和EEMD；EEMD 在時(shí)間跨度較長的情況下降噪性能更好，與小波降噪性能接近；EMD 降噪效果相對(duì)最差.

3)小波降噪相比于EMD 和EEMD，所得到的降噪序列中含有的有色噪聲更少，更有利于對(duì)數(shù)據(jù)進(jìn)行后續(xù)的分析和使用.