張擁軍

（鄭州市市政工程勘測(cè)設(shè)計(jì)研究院，河南 鄭州 450046）

0 引言

針對(duì)海量的GPS 坐標(biāo)時(shí)間序列進(jìn)行分析是目前研究的熱點(diǎn)之一，國(guó)內(nèi)外學(xué)者開(kāi)展了大量的研究。但GPS坐標(biāo)時(shí)間序列中不可避免地會(huì)產(chǎn)生噪聲，而噪聲可能會(huì)與GPS 坐標(biāo)時(shí)間序列中的信號(hào)產(chǎn)生混疊，淹沒(méi)有用的信號(hào)，因此將GPS 坐標(biāo)時(shí)間序列中的噪聲進(jìn)行剔除具有重要的意義。為了剔除GPS 坐標(biāo)時(shí)間序列中的噪聲，不同學(xué)者采用了多種不同的模型進(jìn)行了研究。魯鐵定等人采用改進(jìn)的經(jīng)驗(yàn)?zāi)Ｊ椒纸鈱?duì)GPS 高程時(shí)間序列中的噪聲進(jìn)行剔除研究；黃聲亨等人采用小波方法對(duì)坐標(biāo)時(shí)間序列數(shù)據(jù)中的噪聲進(jìn)行剔除；馬俊等人采用小波包系數(shù)信息熵法去除GPS 坐標(biāo)時(shí)間序列中的有色噪聲；焦雄風(fēng)等人采用抗差卡爾曼濾波（Kalman 濾波）算法對(duì)GPS 監(jiān)測(cè)數(shù)據(jù)去噪。上述方法各有優(yōu)、缺點(diǎn)及適應(yīng)的數(shù)據(jù)，該文采用奇異譜分析方法（Singular Spectrum Analysis，SSA）對(duì)GPS 坐標(biāo)時(shí)間序列噪聲進(jìn)行剔除。SSA 方法是主成分分析方法的一種特殊場(chǎng)景，是對(duì)一維的時(shí)間序列進(jìn)行分析的主成分分析方法。SSA 方法不受正弦波假定的約束，無(wú)須先驗(yàn)信息，具有穩(wěn)定識(shí)別和強(qiáng)化周期信號(hào)的優(yōu)點(diǎn)，能夠較好地從時(shí)間序列中提取信號(hào)，目前已被廣泛應(yīng)用于大地測(cè)量領(lǐng)域。

1 SSA 方法與其他方法原理及精度比較指標(biāo)

1.1 SSA 方法原理

設(shè)GPS 坐標(biāo)時(shí)間序列為,,,… , x；x為GPS 時(shí)間序列中的觀測(cè)值；為時(shí)間序列長(zhǎng)度，且＞2；選用合適的窗口長(zhǎng)度，為一個(gè)實(shí)數(shù)，且＜/2，從而構(gòu)造出滯后矩陣，其矩陣形式如公式（1）所示。

式中：x為,,,… , x中的第個(gè)實(shí)數(shù)值，其他類同。根據(jù)公式（1）可知，滯后矩陣副對(duì)角線值相等，然后計(jì)算滯后矩陣，構(gòu)造自協(xié)方差矩陣，自協(xié)方差矩陣=，其矩陣如公式（2）所示。

式中：c 為實(shí)數(shù)值，其他類同。

由于自協(xié)方差矩陣等于滯后矩陣與其轉(zhuǎn)置矩陣相乘，即=，可知為×的矩陣，因此 c為實(shí)數(shù)值，其他類同。然后利用奇異值分解分別求出自協(xié)方差矩陣的特征值 λ與對(duì)應(yīng)的特征向量E。再計(jì)算滯后矩陣在特征向量的投影 a，如公式（3）所示。

式中：為窗口長(zhǎng)度，是一個(gè)正整數(shù)；為1～中的數(shù)，為滯后矩陣中的數(shù)； λ為自協(xié)方差矩陣的特征值，通過(guò)奇異值分解得到；E為自協(xié)方差矩陣特征值對(duì)應(yīng)的特征向量；E為特征向量E中的一項(xiàng)；a 為滯后矩陣在特征向量上的投影。

奇異譜分析的重建則是由時(shí)間主成分和時(shí)間主分量來(lái)重建，重建后的值為ix，具體形式如公式（4）所示。

1.2 小波方法

小波方法的基本原理是利用多個(gè)小波基去逼近信號(hào)或可積函數(shù)。對(duì)任意信號(hào)或可積函數(shù)（），它的小波變換為該信號(hào)與小波函數(shù), ()的內(nèi)積，如公式（6）所示。

因此通過(guò)選擇合適的小波函數(shù)，就能獲取信號(hào)頻域與時(shí)域的局部信息。小波去噪主要分為3 個(gè)步驟：1） 選擇合適的小波函數(shù)，確定其分解層次，并分別計(jì)算其分解系數(shù)。2）對(duì)每個(gè)分解層次選擇恰當(dāng)?shù)拈撝?，處理其高頻系數(shù)，即剔除高頻部分的噪聲成分。3） 將處理后的高頻系數(shù)與低頻系數(shù)進(jìn)行重構(gòu)，重構(gòu)后的信號(hào)就是去噪后信號(hào)。

1.3 經(jīng)驗(yàn)?zāi)Ｊ椒纸夥?/h3>

經(jīng)驗(yàn)?zāi)Ｊ椒纸猓‥mpirical Mode Decomposition，EMD）是黃鍔提出的自適應(yīng)局部時(shí)頻分析方法。經(jīng)驗(yàn)?zāi)Ｊ椒纸馐菍⑷我庖痪S時(shí)間序列分解為不同尺度的固有模態(tài)函數(shù)，固有模態(tài)函數(shù)必須同時(shí)滿足以下2 個(gè)條件：1） 局部極值點(diǎn)和過(guò)零點(diǎn)的數(shù)目必須相等或最多相差1 個(gè)。2） 在任意時(shí)刻，局部最大值與最小值的包絡(luò)平均為0。

經(jīng)驗(yàn)?zāi)Ｊ椒纸膺^(guò)程如下：首先，確定時(shí)間序列的極大值點(diǎn)與極小值點(diǎn)，并用3 次樣條插值法擬合極大值點(diǎn)與極小值點(diǎn)，分別得到其上、下包絡(luò)線。其次，對(duì)上、下包絡(luò)線進(jìn)行均值處理，得到均值曲線。再次，用時(shí)間序列減去其均值曲線，若滿足給定閾值要求，則得到本征函數(shù)。若不滿足則重復(fù)上述過(guò)程。最后，用時(shí)間序列減去所有本征函數(shù)，當(dāng)殘差序列不超過(guò)2 個(gè)極值點(diǎn)時(shí)，經(jīng)驗(yàn)?zāi)Ｊ椒纸馔瓿伞?/p>

1.4 評(píng)價(jià)指標(biāo)

為更加精確地評(píng)價(jià)奇異譜分析去噪的有效性與優(yōu)越性，該文采用了3 個(gè)評(píng)價(jià)指標(biāo)對(duì)上述3 種方法，即SSA 方法、小波方法和經(jīng)驗(yàn)?zāi)Ｊ椒纸夥ㄟM(jìn)行定量的比較分析。評(píng)價(jià)指標(biāo)分別為噪聲序列的均方根值；去噪后信號(hào)殘差序列的均方根值；去噪后的信號(hào)與模擬信號(hào)的相關(guān)系數(shù)。

2 模擬試驗(yàn)分析

為驗(yàn)證SSA 去噪效果的有效性與優(yōu)越性，模擬1 組具有不同周期的正弦波時(shí)間序列數(shù)據(jù)，然后加入噪聲，其數(shù)據(jù)模擬模型y如公式（7）所示。設(shè)y由u與e構(gòu)成，其中u由周期分別為1200 s、500 s 和200 s 的3 種正弦波信號(hào)構(gòu)成，即u=sin（2π/1200）+sin（2π/500）+sin（2π/200），e為服從正態(tài)分布（0，0.5）的噪聲。

設(shè)模擬數(shù)據(jù)的數(shù)量為3 000 個(gè)，時(shí)間間隔單位為s。為驗(yàn)證SSA 方法的穩(wěn)健性，該文又分別將小波方法（db8 小波基）、經(jīng)驗(yàn)?zāi)Ｊ椒纸夥ǎ‥MD）與SSA 方法進(jìn)行了比較，其結(jié)果如圖1 所示。

圖1中各行子圖從上至下的含義如下：第一行為模擬數(shù)據(jù)；第二行為加入（0，0.52）噪聲后的模擬數(shù)據(jù)；第三行是采用3 種方法去噪后的數(shù)據(jù)；第四行是采用3 種方法去噪后數(shù)據(jù)與模擬數(shù)據(jù)的差值；第五行是3 種方法去噪后模擬數(shù)據(jù)與加噪后模擬數(shù)據(jù)的差值。

圖1 3 種濾波方法的模擬試驗(yàn)結(jié)果

通過(guò)對(duì)比可以看出，3 種方法均能有效地去除模擬數(shù)據(jù)中的噪聲，其形態(tài)整體相似，說(shuō)明SSA 方法具有良好的去噪能力。為了說(shuō)明SSA 方法的優(yōu)越性，該文又分別將值、值與相關(guān)系數(shù)這3 種指標(biāo)進(jìn)行比較，比較結(jié)果見(jiàn)表1。

表1 3 種方法去噪效果對(duì)比

從表1 中噪聲部分的值與相關(guān)系數(shù)可以看出，其值大小基本相當(dāng)，說(shuō)明去除的噪聲部分相同，模擬數(shù)據(jù)與去噪后的數(shù)據(jù)強(qiáng)相關(guān)，再次說(shuō)明3 種方法都有效地剔除了噪聲。但3 種方法的值有較大差異，其中SSA 方法模擬數(shù)據(jù)與去噪后的數(shù)據(jù)相差相對(duì)較小，而小波方法與EMD 方法則相差較大，說(shuō)明后2 種方法在去噪過(guò)程中剔除了模擬的正弦波信號(hào)。因此通過(guò)綜合對(duì)比可知，采用SSA 方法去噪更有優(yōu)越性。

3 試?yán)治?/h2>

該試驗(yàn)數(shù)據(jù)采用GPS 基準(zhǔn)站ZIMM 站的GPS 坐標(biāo)時(shí)間序列，該站有較長(zhǎng)時(shí)間跨度的時(shí)間坐標(biāo)序列，該文選用跨度為1999—2011 年共12 a 的數(shù)據(jù)，ZIMM 站坐標(biāo)時(shí)間序列包括趨勢(shì)項(xiàng)、周期信號(hào)以及噪聲，較適于驗(yàn)證奇異譜分析方法去噪的效果。為了方便比較，該文對(duì)ZIMM 站坐標(biāo)時(shí)間序列進(jìn)行了預(yù)處理。首先，去除ZIMM 站坐標(biāo)時(shí)間序列的趨勢(shì)項(xiàng)及周期項(xiàng)。其次，采用3 倍中誤差方法剔除ZIMM 站坐標(biāo)時(shí)間序列中的中誤差。再次，對(duì)剔除粗差后的坐標(biāo)時(shí)間序列進(jìn)行線性插值，得到處理后的坐標(biāo)時(shí)間序列。最后，采用奇異譜分析方法（SSA 方法）進(jìn)行去噪，其結(jié)果如圖2 所示。

根據(jù)圖2 可知，奇異譜分析方法能夠有效地去除噪聲，并能有效地保留GPS 坐標(biāo)時(shí)間序列中的趨勢(shì)項(xiàng)及周期項(xiàng)，且圖2（b）的噪聲序列比較平穩(wěn)，說(shuō)明奇異譜分析能夠有效地剔除噪聲。為了進(jìn)一步說(shuō)明奇異譜分析方法去噪的效果，該文又對(duì)其進(jìn)行功率譜分析，其結(jié)果如圖3 所示。

圖2 GPS 坐標(biāo)時(shí)間序列去噪前后的對(duì)比及殘差序列

根據(jù)圖3 可知，在低頻域，GPS 坐標(biāo)時(shí)間序列基本保持不變，其存在年周期信號(hào)與半年周期信號(hào)，且信號(hào)的波峰較為一致，說(shuō)明奇異譜分析方法較好地保留了其趨勢(shì)項(xiàng)與周期項(xiàng)。而在高頻域，去噪后的坐標(biāo)時(shí)間序列其功率譜密度明顯低于未去噪的坐標(biāo)時(shí)間序列，該現(xiàn)象說(shuō)明奇異譜分析方法剔除的是GPS 坐標(biāo)時(shí)間序列中的噪聲。

圖3 GPS 坐標(biāo)時(shí)間序列去噪前后的功率譜對(duì)比

4 結(jié)論

奇異譜分析是一種高效的時(shí)間序列分析方法，該方法是對(duì)一維時(shí)間序列進(jìn)行分析的方法，能夠較好地從坐標(biāo)時(shí)間序列中提取信息。該文通過(guò)模擬數(shù)據(jù)分析，驗(yàn)證了和小波方法與經(jīng)驗(yàn)?zāi)Ｊ椒纸夥ㄏ啾?，該方法去噪效果更好，能夠在最大程度地保留坐?biāo)時(shí)間序列中有用的信號(hào)。該文還通過(guò)GPS 坐標(biāo)時(shí)間序列實(shí)測(cè)數(shù)據(jù)分析，表明該方法能夠有效地剔除GPS 坐標(biāo)時(shí)間序列中的噪聲，且去噪后的坐標(biāo)時(shí)間序列更加穩(wěn)定和平滑。