楊厚明 魯鐵定,2 孫喜文 何錦亮

1 東華理工大學(xué)測(cè)繪與空間信息工程學(xué)院,南昌市廣蘭大道 418號(hào),330013

2 東華理工大學(xué)自然資源部環(huán)鄱陽(yáng)湖區(qū)域礦山環(huán)境監(jiān)測(cè)與治理重點(diǎn)實(shí)驗(yàn)室,南昌市廣蘭大道 418號(hào),330013

受多路徑效應(yīng)、鐘差、電離層延遲等因素的影響,GNSS坐標(biāo)時(shí)間序列存在各種噪聲[1-2],呈現(xiàn)明顯的非線性變化,不能準(zhǔn)確反映測(cè)站的實(shí)際運(yùn)動(dòng)信息。常用的GNSS坐標(biāo)時(shí)序降噪方法主要有小波分解(WD)[3-4]、經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)[5]等。EMD在處理非線性、非平穩(wěn)信號(hào)方面應(yīng)用廣泛,但存在一定的端點(diǎn)效應(yīng)和模態(tài)混疊現(xiàn)象,影響降噪效果。Wu等[6]和Yeh等[7]對(duì)EMD進(jìn)行優(yōu)化和改進(jìn),分別提出集合經(jīng)驗(yàn)?zāi)B(tài)分解(EEMD)和互補(bǔ)集合經(jīng)驗(yàn)?zāi)B(tài)分解(CEEMD)。兩者都能有效抑制模態(tài)混疊現(xiàn)象,但計(jì)算效率低且過(guò)程繁瑣。Dragomiretskiy等[8]提出一種新的信號(hào)時(shí)頻分析處理方法——變分模態(tài)分解(VMD),該方法可以有效分離IMF分量、劃分信號(hào)的頻域,避免了EMD方法中模態(tài)混疊等缺陷,具有較好的魯棒性。魯鐵定等[9]將VMD應(yīng)用到變形監(jiān)測(cè)數(shù)據(jù)的降噪中,效果顯著。但VMD方法需要預(yù)先設(shè)置模態(tài)分解個(gè)數(shù)K和二次懲罰因子α,這2個(gè)參數(shù)在大多數(shù)情況下是依據(jù)經(jīng)驗(yàn)選取的,然而實(shí)測(cè)信號(hào)復(fù)雜多變,若參數(shù)設(shè)定不當(dāng)會(huì)對(duì)分解效果產(chǎn)生嚴(yán)重影響。

基于以上研究,本文先利用麻雀搜索算法(SSA)優(yōu)化VMD,然后結(jié)合WD方法提出一種GNSS坐標(biāo)時(shí)間序列降噪方法IVMD-WD,并結(jié)合仿真信號(hào)和GNSS實(shí)測(cè)數(shù)據(jù)實(shí)驗(yàn)驗(yàn)證該方法的有效性。

1 算法原理

1.1 VMD基本原理

VMD處理非線性和非平穩(wěn)信號(hào)時(shí)效果較好[8],其實(shí)質(zhì)是一個(gè)變分問(wèn)題的構(gòu)造和求解過(guò)程。VMD構(gòu)造的約束變分問(wèn)題可表示為[9]:

(1)

式中,t為時(shí)間,f為原始信號(hào),uk為模態(tài)函數(shù),ωk為各模態(tài)的實(shí)際中心頻率,e-jωkt為每個(gè)解析信號(hào)的預(yù)估中心頻率,‖·‖2為求L2范數(shù),s.t.表示約束條件。

使用二次懲罰因子α和Lagrange乘子λ(t),將約束變分問(wèn)題轉(zhuǎn)換為無(wú)約束變分問(wèn)題,進(jìn)而求得式(1)的最優(yōu)解,得到增廣Lagrange表達(dá)式為:

L({uk},{ωk},λ)=

(2)

1.2 SSA基本原理

SSA是Xue等[11]受麻雀覓食行為啟發(fā)而提出的一種新型群體智能優(yōu)化算法。相較于粒子群算法、遺傳算法等優(yōu)化算法,SSA收斂速度更快、求解精度更高、穩(wěn)定性和魯棒性更好[11-12]。麻雀種群按其職能分為發(fā)現(xiàn)者和追隨者,發(fā)現(xiàn)者一般占比為10%～20%,主要負(fù)責(zé)為種群尋找食物和提供覓食的區(qū)域和方向,剩余麻雀均為依賴(lài)發(fā)現(xiàn)者來(lái)獲取食物的追隨者。此外,部分麻雀帶有預(yù)警機(jī)制,一般占比為10%～20%,當(dāng)遇到危險(xiǎn)時(shí),會(huì)發(fā)出報(bào)警信號(hào),麻雀會(huì)進(jìn)行反捕食。詳細(xì)過(guò)程參見(jiàn)文獻(xiàn)[12]。

1.3 SSA優(yōu)化VMD

對(duì)VMD計(jì)算結(jié)果影響較大的是模態(tài)分解數(shù)K和懲罰因子α,其余參數(shù)一般設(shè)置為默認(rèn)值。當(dāng)K過(guò)小時(shí),會(huì)使信號(hào)欠分解;當(dāng)K過(guò)大時(shí),會(huì)造成信號(hào)過(guò)分解并產(chǎn)生模態(tài)混疊現(xiàn)象[10]。本文利用SSA對(duì)VMD的參數(shù)K和α進(jìn)行優(yōu)化,可快速準(zhǔn)確地獲取優(yōu)化后的參數(shù)。

選用包絡(luò)熵為SSA優(yōu)化VMD的適應(yīng)度函數(shù)。包絡(luò)熵可以較好地反映原始信號(hào)的稀疏特性和不確定性,當(dāng)信號(hào)中噪聲較多時(shí),熵值較大;反之,熵值較小。包絡(luò)熵的原理[10]為:

(3)

式中,N為信號(hào)的采樣點(diǎn)數(shù),pj為a(j)的歸一化形式,a(j)為信號(hào)x(j)經(jīng)Hilbert解調(diào)后得到的包絡(luò)信號(hào)。

SSA優(yōu)化VMD的過(guò)程見(jiàn)圖1。

圖1 SSA優(yōu)化VMD的流程

1.4 WD基本原理

WD可以避免傳統(tǒng)傅里葉變換的缺陷,具有較好的時(shí)頻局部分析和多分辨率分析性能,因此在非線性非平穩(wěn)信號(hào)研究中應(yīng)用廣泛[3,5]。

WD通過(guò)一組高通和低通濾波器將原始信號(hào)分解為低頻和高頻分量,然后將低頻分量分解。小波基函數(shù)表示為[13]:

(4)

式中,φ(t)為基小波或者母小波函數(shù),經(jīng)過(guò)尺度因子a和平移因子b變換后的φa,b(t)統(tǒng)稱(chēng)為小波函數(shù)。

小波分解過(guò)程中最重要的就是確定小波基函數(shù)和分解層數(shù),本文選用正則性較好的db4小波進(jìn)行分解,最佳層數(shù)由文獻(xiàn)[3]的方法確定。

1.5 IVMD-WD方法構(gòu)建

使用SSA優(yōu)化影響VMD分解效果的2個(gè)關(guān)鍵參數(shù)K和α。對(duì)于VMD分解后的多個(gè)IMF分量,采用多尺度排列熵(multi-scale permutation entropy,MPE)作為判斷噪聲和信號(hào)的標(biāo)準(zhǔn),具體計(jì)算步驟參考文獻(xiàn)[14]。IVMD-WD方法的降噪流程見(jiàn)圖2,具體步驟如下:

圖2 IVMD-WD方法的降噪流程

1)初始化SSA參數(shù),設(shè)置麻雀種群數(shù)量為30,最大迭代次數(shù)為10,考慮計(jì)算效率和算法精度,將K的取值范圍設(shè)為[3,10],α的取值范圍設(shè)為[100,3 500]。

2)根據(jù)得到的最優(yōu)參數(shù)組合[K,α],對(duì)原參考信號(hào)進(jìn)行VMD處理,得到K個(gè)IMF。

3)計(jì)算各IMF的MPE,設(shè)定MPE的閾值,根據(jù)閾值大小判斷出有效IMF并重構(gòu)為信號(hào),剩余IMF重構(gòu)為噪聲。

4)將步驟3)中重構(gòu)后的高頻噪聲部分利用小波分解再次降噪,利用相關(guān)系數(shù)判斷得到有效信號(hào),將步驟3)中得到的低頻信號(hào)與小波分解處理后的信號(hào)重構(gòu)為最終的降噪信號(hào)。

2 仿真信號(hào)算例分析

2.1 仿真信號(hào)降噪實(shí)驗(yàn)

進(jìn)行仿真信號(hào)實(shí)驗(yàn),并將IVMD-WD方法與EMD、WD和EEMD方法進(jìn)行對(duì)比。實(shí)驗(yàn)過(guò)程中,EMD、WD和EEMD方法均利用相關(guān)系數(shù)法分離噪聲和信號(hào)。仿真信號(hào)由3個(gè)周期項(xiàng)和高斯白噪聲組成,采樣間隔設(shè)置為1 s,采樣點(diǎn)數(shù)為1 024,并加入信噪比(signal noise ratio, SNR)為 6 dB 的高斯白噪聲,仿真信號(hào)的分量波形見(jiàn)圖3,表達(dá)式如式(5)所示:

(5)

圖3 仿真信號(hào)各分量波形

本文方法進(jìn)行降噪時(shí)需先利用SSA尋找VMD的最優(yōu)參數(shù)組合,采用包絡(luò)熵作為適應(yīng)度函數(shù),適應(yīng)度函數(shù)值在SSA尋優(yōu)過(guò)程中隨迭代次數(shù)的變化見(jiàn)圖4。由圖4可知,適應(yīng)度值在迭代到第2次時(shí)達(dá)到最小,此時(shí)對(duì)應(yīng)的[K,α]=[7,2 730]為最優(yōu)參數(shù)組合。

圖4 適應(yīng)度值收斂圖

利用SSA得到的最優(yōu)參數(shù)組合對(duì)信號(hào)進(jìn)行VMD處理,得到的7個(gè)IMF見(jiàn)圖5。由圖5可知,低頻分量主要集中在前2階模態(tài)中。

圖5 IVMD分解圖

為有效分離出低頻信號(hào)和高頻噪聲,需計(jì)算出各IMF的MPE值。計(jì)算MPE時(shí)需設(shè)置適當(dāng)?shù)膮?shù),本文參考文獻(xiàn)[15]對(duì)關(guān)鍵參數(shù)進(jìn)行設(shè)定:尺度因子s=12、嵌入維數(shù)m=6、延遲時(shí)間τ=1。計(jì)算不同尺度因子下各IMF的排列熵均值作為最終的MPE值,MPE值越接近于1,表明時(shí)間序列具有越大的隨機(jī)波動(dòng)性和不規(guī)則性。MPE閾值一般取0.6,將大于0.6的IMF視為噪聲分量,小于0.6的視為低頻信號(hào)分量[15]。VMD得到的各IMF的MPE值見(jiàn)表1。

表1 IVMD得到的各IMF的MPE值

由表1可見(jiàn),IMF1～I(xiàn)MF7的MPE值逐漸增大,表明序列的隨機(jī)波動(dòng)性越來(lái)越大,噪聲成分逐漸增多,符合圖5的波形描述。IMF1、IMF2的MPE值均小于0.6,因此,將IMF1、IMF2重構(gòu)為低頻信號(hào),IMF3～I(xiàn)MF7重構(gòu)為高頻噪聲,然后利用小波分解處理高頻噪聲部分。同時(shí),使用EMD、WD、EEMD方法對(duì)仿真信號(hào)進(jìn)行降噪處理,并與本文方法進(jìn)行對(duì)比,見(jiàn)圖6?？梢钥闯?EMD、WD、EEMD方法降噪后,信號(hào)波形與原始序列擬合效果較差,而本文方法降噪后的波形與仿真信號(hào)更加接近,曲線更為光滑,避免了EMD降噪過(guò)程中的端點(diǎn)效應(yīng)和模態(tài)混疊問(wèn)題,能夠提取更多的有效信號(hào),降噪效果更好。

圖6 4種降噪信號(hào)波形對(duì)比

2.2 仿真信號(hào)降噪效果評(píng)價(jià)

為了評(píng)價(jià)上述4種方法的降噪效果,選取均方根誤差(RMSE)、信噪比(SNR)、平均絕對(duì)誤差(MAE)和互相關(guān)系數(shù)(R)作為評(píng)價(jià)指標(biāo):

(6)

(7)

(8)

(9)

表2 4種方法的降噪效果

由表2可見(jiàn),WD降噪效果略?xún)?yōu)于EMD和EEMD,而本文方法相比于EMD、WD和EEMD方法,各降噪評(píng)價(jià)指標(biāo)均為最優(yōu),RMSE分別降低15.42%、14.35%和15.06%,MAE分別降低15.62%、13.37%和15.18%,SNR分別提高1.44 dB、1.35 dB和1.48 dB,R分別提高0.028 2、0.025 6和0.027 0。IVMD-WD方法的RMSE和MAE最小,表明其降噪后的序列偏差更小;SNR和R最大,表明其提取的有效信號(hào)最多,降噪信號(hào)波形與原始信號(hào)波形更相似,擬合度更高。仿真信號(hào)實(shí)驗(yàn)結(jié)果表明,IVMD-WD方法的降噪效果優(yōu)于EMD、WD和EEMD方法。

3 GNSS實(shí)測(cè)數(shù)據(jù)降噪分析

為進(jìn)一步驗(yàn)證本文方法的可靠性和適用性,選取SOPAC(Scrips Orbit and Permanent Array Center) 提供的美國(guó)西海岸10個(gè)GNSS基準(zhǔn)站網(wǎng)原始坐標(biāo)時(shí)間序列進(jìn)行分析,觀測(cè)時(shí)間為2000～2023年,采樣間隔為1/365.25 a。

限于篇幅,僅以GOBS站為例進(jìn)行詳細(xì)分析,圖7為該站N、E、U方向的原始時(shí)間序列與經(jīng)4種方法降噪后的信號(hào)對(duì)比結(jié)果?？梢钥闯? EMD、WD和EEMD降噪后信號(hào)波形與原始時(shí)間序列的波形存在一定差異,而IVMD-WD方法去噪后的時(shí)間序列能夠更好地?cái)M合原始時(shí)間序列,可以有效地反映局部運(yùn)動(dòng)趨勢(shì),且周期振蕩較小。圖8為利用4種方法得到的降噪信號(hào)與原始時(shí)間序列的殘差對(duì)比。由圖可知,IVMD-WD方法的平均殘差最小,表明本文方法有效濾除了高頻噪聲,降噪效果最顯著。

圖7 GOBS站點(diǎn)3個(gè)方向降噪效果對(duì)比

圖8 降噪信號(hào)與原始信號(hào)的殘差結(jié)果

圖9為10個(gè)站點(diǎn)采用上述4種方法降噪后信號(hào)的4種評(píng)價(jià)指標(biāo)值對(duì)比。由圖9(a)和9(b)可知,EMD方法降噪效果在部分站點(diǎn)上要優(yōu)于EEMD,整體來(lái)說(shuō),兩者降噪效果相差不大;WD降噪效果整體上略低于EMD,略?xún)?yōu)于EEMD;相比于EMD、WD和EEMD方法,本文方法降噪后得到的RMSE和MAE最小。由圖9(c)和9(d)可知,本文方法降噪后的SNR和R最大。綜合可知,本文IVMD-WD方法降噪后的各項(xiàng)評(píng)價(jià)指標(biāo)均為最優(yōu),表明本文方法降噪效果最顯著,能夠有效提取出更多的有用信號(hào)。

圖9 4種降噪方法的評(píng)價(jià)指標(biāo)

4 結(jié) 語(yǔ)

本文針對(duì)GNSS坐標(biāo)時(shí)間序列的非線性、非平穩(wěn)特性構(gòu)建了一種新的降噪算法IVMD-WD。該方法利用麻雀搜索算法優(yōu)化VMD參數(shù),使用多尺度排列熵作為噪聲和信號(hào)的篩選標(biāo)準(zhǔn),結(jié)合小波分解算法對(duì)VMD進(jìn)行改進(jìn)。通過(guò)仿真信號(hào)和10個(gè)GNSS基準(zhǔn)站的實(shí)測(cè)數(shù)據(jù)進(jìn)行降噪分析。結(jié)果表明,與傳統(tǒng)的EMD、WD和EEMD方法相比,本文IVMD-WD方法降噪后的各指標(biāo)均為最優(yōu),能夠更有效地剔除原始時(shí)間序列中的高頻噪聲并保留有用信號(hào),在降噪性能上具有明顯優(yōu)勢(shì)。