李 婷
(浙江省國土勘測規劃有限公司 浙江杭州 310030)
作為重要的空間數據基礎設施,連續運行參考站(CORS)為社會各行業提供穩定的位置信息服務。目前,國內外許多國家和地區已經開展或建成了CORS系統,建立全天候、全覆蓋、高精度和實時動態定位的衛星導航系統是國際大地測量的一個發展趨勢。對CORS站進行觀測時,環境因素的影響使得其觀測數據在高程方向存在較大噪聲。因此,通過剔除觀測數據中的噪聲,獲取CORS站高程方向真實值具有重要的意義。
對于CORS站高程時間序列的非線性、非平穩信號,通常利用事先給定的周期函數模型擬合趨勢項和周期項去除非平穩性,但是,這種方法得到的分析結果不準確,甚至會偏離序列真實的運動情況。近些年來,整體經驗模態分解(Ensemble Empirical Mode Decomposition,EEMD)方法被廣泛應用到時間序列分析中,經EEMD法分解后的信號可以得到若干含有特征時間尺度的本征模態函數(Intrinsic Mode Function,IMF)分量,從而對信號進行多分辨率分析[1]。吉長東等[2]將EEMD方法應用到CORS站時間序列分析中,有效地改善了信號分解過程中產生的模態分解問題,但還會存在整體平均不能完全消除噪聲的影響。本文引入一種新的噪聲自適應整體經驗模態分解(Empirical Mode Decomposition with Adaptive Noise,EEMDAN)方法,與EEMD方法不同的是,EEMDAN方法在經驗模態分解處理時添加白噪聲分解的模態分量,可以改進CORS站高程時間序列信號分解的模態混疊。最后利用某市高程CORS站點數據,驗證EEMDAN方法處理CORS站高程時間序列的有效性。
EEMD可以對待處理信號的極大值特征尺度進行信號分解,通過不斷篩選,將待分解信號分解為頻率各不相同的IMF,根據頻率的不同來反映周期長度,體現了對信號多分辨率分析的濾波過程。EEMD方法的濾波過程是通過對信號進行EEMD分解產生IMF分量,根據IMF分量特征,構造時空濾波器。
分解后的IMF可以反映信號不同范圍的特征尺度,通過特征可以對信號進行濾波。以往是通過傅里葉變換在頻域上實現濾波,其缺點是會在頻域產生諧波信號,很難將信號的頻譜有效分離。
對監測體進行變形監測采集的GNSS數據,噪聲與有用信號的時頻特性是有區別的。GNSS時間序列中通常存在由于多路徑影響產生的頻率低于結構特性的低頻噪聲,也會存在頻率高于結構特性的高頻噪聲。結合GNSS的多路徑特點,對GNSS監測數據經EEMD分解產生的不同頻率特征相關性進行分析,從而消除多路徑誤差對觀測信號的影響。
EEMDAN是對EEMD方法的改進。與EEMD分解所有 IMF分量進行整體平均不同的是,EEMDAN分解得到第一階分量后即進行整體平均,并將第一階整體平均IMF剔除后得到除第一階整體平均IMF的剩余分量,再對剩余分量進行整體平均,得到第二階、第三階直至最后一階IMF分量[3]。EEMDAN濾波過程的主要步驟如下[4]:
第一步:將白噪聲添加到信號中做經驗模態分解(Empirical Mode Decomposition,EMD),對分解的IMF分量進行平均化處理,獲取第一階IMF分量。
第二步:剔除第一階分量后,得到剩余信號,再分解后得到含有噪聲的模態分量集;在剩余信號中加入該模態分量集,組合成分解信號,繼續進行EMD分解,得到第二階分量。
第三步:剔除第二階分量后,重復第二步,直至剔除最后一階模態分量后得到剩余信號。
EEMDAN利用白噪聲頻譜分布均勻這一特性,對CORS站點的時間序列進行投影,以達到消除由分解帶來的模態混疊。EEMDAN可以將上一次分解后的剩余信號作為下一次分解信號,減小了重構誤差和整體平均噪聲在沒有消除的情況下對信號分解造成的影響。
對CORS站點時間序列進行分解,可以通過兩個指標評價其分解效果,一是分解后的IMF分量是否存在模態混疊;二是參考正交性指標(Orthogona-lity Index,IO),分解精度隨著IO的減小而提高。IO可表示為
(1)
將CORS站時間序列分解為一系列不同的IMF分量,若將Tj和CMSEj定義為第j個IMF分量的平均周期與能量密度,則可以通過參考分量的能量密度與平均周期的乘積得到信號與噪聲的臨界位置(const值),通過判別臨界位置來剔除包含噪聲的分量,完成對信號的重構。通過計算重構后的相關系數R和均方根誤差(RMSE)來量化重構效果。相關系數R越大,重構效果就越好;RMSE值越小,重構信號的精度與質量越高[5-7]。
R和RMSE可分別表示為
(2)
(3)
本文選取某市CORS站BTRD站點2016年1月1日—2019年12月31日的GNSS高程觀測數據為試驗數據,對高程觀測數據進行均值化處理,圖1所示為高程觀測數據時間序列。
圖1 BTRD站點2016—2019年高程方向上時間序列
為了便于對比,將原始信號標準差與噪聲的比值設置為10,經EEMD和EEMDAN兩種方法分解得到部分低頻IMF分量如圖2所示(實線表示的是EEMD分解結果,虛線表示的是EEMDAN分解結果)。
圖2 BTRD站點高程時間序列分解得到的部分IMF分量
本文使用IMF分量的方差貢獻率獲取CORS站點時間序列運動的周期形式。IMF分量在整個時間序列中的比重可以用該分量的方差貢獻率反映,時間序列中的主要周期運動貢獻項由方差貢獻率高的IMF分量組成。BRTD站點高程時間序列經EEMD和EEMDAN兩種方法分解的低頻分量的方差貢獻率如表1所示。
表1 不同IMF分量的方差貢獻率 單位:%Tab.1 Variance Contribution Rate in Different IMF Components 方法c4c5c6c7c8c9EEMD21.428.64.11.25.70.8EEMDAN24.826.33.93.10.91.4
由表1可知,BRTD站點高程時間序列經EEMD分解后,IMFc4與IMFc5方差貢獻率所占比重最高,分別為21.4%與28.6%,為時間序列中的主要周期貢獻項。BRTD站點高程時間序列經EEMDAN分解后,同樣是IMFc4與IMFc5的方差貢獻率占比最高,分別為24.8%與26.3%。
對于BRTD站點高程時間序列來說,可以通過IMFc4與IMFc5兩個分量表現出周期形式。為了評價EEMD和EEMDAN兩種方法的分解精度,計算方差貢獻率較高的IMF分量之間的正交指數,結果如表2所示。
由表2可知,取200次整體平均的分解精度比取100次整體平均的分解精度低,也符合文獻[8]中的推薦值。通過計算EEMD和EEMDAN兩種方法分解后的IMF分量的IO值,可知EEMDAN方法的分解精度比EEMD方法的分解精度提高了53.7%。
表2 方差貢獻率較高的IMF分量之間的正交指數Tab.2 Orthogonal Index of IMF Components in Higher Variance Contribution Rate方法N=100N=200EEMD0.054 10.068 7EEMDAN0.035 20.041 6
BRTD站點高程時間序列經EEMD和EEMDAN兩種分解方法分解后得到不同頻率的IMF分量,其中IMFc1~IMFc3是高頻分量。為了有效剔除噪聲并重構得到有用的信號,通過IMF分量的能量密度與平均周期的乘積指標來識別和剔除噪聲,其計算結果見表3。
表3 BRTD站點分解分量的能量密度與平均周期乘積Tab.3 Product of Energy Density and Average Period of Decomposition Component of BRTD Station方 法指 標c1c2c3c4c5c6c7EEMDT/d3.4067.25814.32631.43360.221177.274384.214CMSE/(cm2/d)0.3570.2140.1550.1640.0980.2860.179T×CMSE/cm21.2161.5532.2205.1555.90250.70068.774EEMDANT/d3.3577.16615.82328.16344.581159.548181.532CMSE/(cm2/d)0.2440.1950.0970.1470.1260.3270.218T×CMSE/cm20.8191.3971.5354.1405.61752.17239.574
陳仁祥等[9]根據CMSE值是否發生明顯變化來判斷信號與噪聲的臨界點,但是,這種判別方式在進行噪聲剔除時,會將雙月與近似月的分量剔除。故本文引入const值作為噪聲與信號的判別標準[10-12]。
對于BRTD站高程時間序列,經EEMDAN處理后的前3個分量的const均值為1.250,前3個分量const值與均值的最大偏差不超過34.5%,而IMFc4的const值與均值的偏差超過2倍;經EEMD處理后的前3個分量的const均值為1.663,前3個分量const值與均值的最大偏差不超過36.8%,最小偏差為7.1%,而IMFc4的const值與均值的偏差超過2倍。因此,EEMD與EEMDAN分解得到的結果相同。
const值對噪聲與信號的區別判別指標主要是針對信號中高斯白噪聲,對有色噪聲的處理效果并不理想。通過噪聲識別對提取的信號進行重構,其結果如圖3所示,重構信號的均方根誤差與相關系數見表4。
表4 信號重構后的均方根誤差與相關系數Tab.4 RMSE and Correlation Coefficient after Signal Reconstruction方法RMSE/mm相關系數/%EEMD1.57596.7EEMDAN1.29698.8
由圖3與表4可知,經EEMD與EEMDAN兩種方法重構后的時間序列均能有效提取信號中的有用信息,但是經EEMDAN處理后的信號均方根誤差比經EEMD處理后的均方根誤差減少21.5%,相關系數增加2.1%,表明經EEMDAN剔除噪聲,重構信號的結果要比經EEMD方法得到的結果更好。
圖3 兩種方法重構后時間序列
這是因為EEMDAN處理CORS站高程時間序列時,較好地利用了白噪聲頻譜均勻分布的特性,使信號能夠在遍布整個時頻空間且分布一致的白噪聲背景上,按照合適的尺度進行投影,以此來削弱分解產生的模態混疊。
以EEMD為基礎的EENDAN方法改進了對整體平均的處理,通過計算分量能量密度與平均周期的乘積、均方根誤差、相關系數和正交性指標的值,驗證了EEMDAN方法的分解精度與重構后信號的效果都要比EEMD方法好。
基于信號特有的運動特征,EEMDAN方法可以將原始信號的有用信息提取出來,為CORS站點時間序列的處理提供有益的參考。但是本文對于如何剔除信號中有色噪聲沒有做進一步的研究。