CORS高程時間序列非線性趨勢項模式識別與估計

張恒璟，王 蕾，金澤林，蘭文琦，齊 昕

CORS高程時間序列非線性趨勢項模式識別與估計

張恒璟，王 蕾，金澤林，蘭文琦，齊 昕

（遼寧工程技術大學 測繪與地理科學學院，遼寧 阜新 123000）

針對傳統線性模型擬合連續運行參考站（CORS）高程時間序列使測站運動趨勢單一化的問題，提出一種基于CORS高程時間序列本身的非線性趨勢項定量識別與估計方法。該方法首先用均值聚類整體經驗模態分解方法，對CORS高程時間序列進行分解；其次，基于排列熵理論對非線性趨勢項定量識別，統計決策適合的趨勢項模型，并估計模型參數。實驗結果表明：基于排列熵的非線性趨勢項定量識別方法避免了人為判別誤差，準確地反映了序列自身的變化趨勢；通過比較傅里葉、正弦函數、線性多項式三種模型的擬合效果，發現傅里葉模型擬合效果最佳。

高程時間序列；趨勢項；排列熵；模型擬合；定量估計

0 引言

傳統的序列趨勢項是由固定的模型擬合得出[1]，經過擬合后能夠得出序列的周期項、常數項、線性速率等信息。針對如地震、海嘯、更換儀器等操作產生的序列影響，通常在固定模型的基礎上，加入各項改正，采用最小二乘求解擬合參數，為后續時頻分析提供基礎。

國內外研究表明，連續運行參考站（continuously operating reference stations, CORS）高程時間序列是非線性非平穩的序列[2-5]，全球范圍內各個CORS所處地理位置不同、受到地球物理因素差異性的影響，導致各個CORS的運動趨勢不完全相同，采用統一的趨勢項擬合模型使各個CORS運動趨勢單一化。且有研究指出[6-7]，趨勢項是一個基于序列自身存在的單調擬合函數，在長期時間序列趨勢項識別中，應選擇基于序列本身識別的方法。排列熵算法作為一種新穎的突變檢測方法，可以檢測動態突變和時間序列隨機性，并定量評估信號序列中包含的隨機噪聲。其優點是：①計算簡單，抗噪聲能力較強；②該算法對時間敏感，可以實現高分辨率；③算法的輸出結果直觀，突變信息更易識別。依據該方法的獨特優勢，本文引入具有自適應性的均值聚類整體經驗模態分解（-means clustering-ensemble empirical mode decomposition,-EEMD）算法，對CORS高程時間序列進行分解，結合趨勢項定義和排列熵理論，對非線性趨勢項進行定量識別與合成，解決了使CORS運動趨勢單一化問題，更好地反映CORS運動趨勢，豐富了CORS運動趨勢提取方法。

1 K-EEMD與排列熵理論

1.1 K-EEMD算法

-EEMD是文獻[8]針對整體經驗模態分解（ensemble empirical mode decomposition, EEMD）信號分解時，信噪比低及模態混疊問題，所提出的一種改進的EEMD算法。該方法首先在EEMD分解前，向原始序列中隨機添加白噪聲；其次，引入統計學均值聚類思想，對EEMD算法100次迭代計算產生的本征模態函數（intrinsic mode function, IMF）分量矩陣，以歐氏距離作為相似程度判斷依據，對IMF分量進行均值聚類分析，不斷迭代直至聚類中心不再發生改變；最后，保留各IMF分量聚類簇中最多的一類結果取平均值作為最終IMF分量。實驗結果表明，該方法獲取了CORS高程時間序列中頻率從高至低的IMF分量，減弱了EEMD算法的模態混疊現象，提高了信噪比，基本消除了EEMD算法得到的低頻IMF分量中殘留高頻噪聲的現象。

1.2 排列熵理論

排列熵（permutation entropy, PE）是一種通過計算時間序列的復雜性程度以及檢驗序列隨機性和突變的方法。該算法設計概念簡單，計算速度快，并且特別適用于非線性非平穩數據，表現出較好的魯棒性[9-11]，其理論如下：

2 非線性趨勢項定量識別與合成

利用國內外BJFS、SYDN、JOEN、HOB2站近20 a的高程時間序列進行實驗。實驗數據來自斯克里普斯軌道和常駐陣列中心（Scripps Orbit and Permanent Array, SOPAC）網站提供的原始（RAW）數據，利用四分位差法對原始數據進行粗差剔除后，采用三次樣條插值法進行插值，已進行去均值處理。如圖1所示，對四個站高程時間序列數據進行-EEMD實驗，圖1中僅給出中低頻的IMF分量。

圖1 K-EEMD分解部分IMF分量

經過-EEMD分解之后，殘差值（residual）可作為各序列非線性趨勢項的基礎項，基礎項呈現單調變化趨勢。但在分解過程中，可能存在過度分解或者長周期項與趨勢項無法識別的現象，此時可以進行趨勢項合成實驗[12]，將低頻IMF分量依次并入基礎項，如圖2所示。

此時，可以采用人工判讀對各合成項進行單調性判別，若合成趨勢項滿足單調特性，則將證明該IMF分量為趨勢項；若產生周期運動趨勢，則該IMF分量為當前時間尺度內無法識別的長周期項。

BJFS站加入IMF11、IMF10后，保持單調向上趨勢，加入IMF9分量后，產生了周期震蕩，即BJFS站的趨勢項為residual+IMF11+IMF10，IMF9分量為長周期項；SYDN站加入IMF11后，保持單調向下趨勢，加入IMF10后，無法準確判斷其是否為單調變化，加入IMF9后，產生了周期震蕩，即SYDN站的趨勢項為residual+IMF11，IMF10無法確定，IMF9為長周期項；JOEN站加入IMF11、IMF10后，保持單調向上趨勢，加入IMF9后則產生了周期震蕩，即JOEN站的趨勢項為residual+IMF11+ IMF10，IMF9為長周期項；HOB2站加入IMF11后，仍保持單調向下趨勢，加入IMF10、IMF9后，產生了周期震蕩，即HOB2站的趨勢項為residual+ IMF11，IMF10和IMF9為長周期項。這種方法缺乏理論指導，且數據量大、較小數據變化在曲線圖中無法表現，導致部分IMF分量無法準確判斷，分別計算各站合成項的排列熵，排列熵統計結果如表1所示。

表1 排列熵統計

當排列熵為0時，原始序列相空間重構后排序固定，表示序列穩定嚴格單調變化，定量地判別了合成項的單調變化趨勢。從排列熵指標不難看出：BJFS站的趨勢項為residual+IMF11+IMF10；SYDN站趨勢項為residual+IMF11+IMF10；JOEN站趨勢項為residual+IMF11+IMF10，HOB2站趨勢項為residual+IMF11。該方法與人工判別結論一致，同時能準確判別出SYDN站IMF10分量為趨勢項的一部分。

排列熵指標顯示，SYDN站余項并非嚴格單調變化，查看余項數據發現，其確實不為嚴格單調變化。這種現象可能是-EEMD分解的端點效應，加入IMF11后，趨勢項嚴格單調變化，也證明了SYDN站的趨勢項應該為residual+IMF11，而不是單純以余項為該站趨勢項。

3 非線性趨勢項定量估計

通過定量識別和合成實驗，得到了經過-EEMD分解的各個CORS非線性趨勢項，結合趨勢項定義選用三種模型：正弦函數模型、多項式模型、傅里葉模型對CORS非線性趨勢項進行擬合。

1）線性多項式模型。在經典的時間序列擬合模型中，通常認為趨勢項是一條直線，所以用線性多項式直接進行擬合，選取最小二乘進行逼近，求解線性方程的系數進而完成擬合。線性多項式的表達式為

2）正弦函數模型。通過本文-EEMD分解時間序列識別和合成趨勢項后，每個CORS站都存在單調非線性趨勢項，采用經典線性擬合方法可能會出現偏差，當三角函數的周期足夠大時，在一定時間尺度上可以保持單調趨勢，因此可以采用正弦函數對非線性趨勢項進行擬合。正弦函數模型的表達式為

3）傅里葉模型。傅里葉擬合是一種通過傅里葉技術不斷逼近數據的曲線擬合方法，可看做多個正弦函數和余弦函數疊加的過程，經過傅里葉變換，對曲線進行擬合。傅里葉模型的表達式為

采用3種擬合方法對4個CORS站高程時間序列擬合結果如圖3所示，圖3中離散點為SOPAC給出的各站運動趨勢，實線為本文識別的該測站非線性趨勢項，余下3種依次為正弦函數擬合、線性多項式擬合和傅里葉函數擬合結果。

從圖3可知，當趨勢項近似線性變化時，3種擬合模型擬合結果相近；當趨勢項變化為非線性時，傅里葉擬合模型明顯優于其他兩種模型。為驗證本文提取趨勢項的正確性，將線性擬合結果與SOPAC給出的趨勢項參考值進行比對，如表2所示。

圖3 趨勢項擬合

表2 線性擬合結果對比

對4個測站的非線性趨勢項線性擬合效果進行觀察，所有測站提取的擬合結果均在SOPAC給出的參考范圍內，說明本文提出的非線性趨勢項識別方法和識別結果是可靠的。為定量比較3種擬合方案的優劣性，采用均方根誤差、擬合優度兩種精度指標，進行擬合精度評定，其評定結果如表3所示。

1）均方根誤差。反映了擬合估計值與真實值間的偏離程度，數值越小表示偏離程度越小。均方根誤差的計算公式為

2）擬合優度。擬合優度表示擬合效果的好壞，越接近1代表該模型的擬合效果越好。擬合優度的計算公式為

表3 擬合精度統計

分析表3可知：當測站趨勢項為非線性時，傅里葉模型擬合效果最佳，正弦函數模型次之，線性多項式擬合效果最差；當測站趨勢線接近線性時，3種擬合模型均有較好的擬合效果，但傅里葉模型擬合最佳。

4 結束語

針對傳統線性模型擬合CORS站高程時間序列使測站運動趨勢單一化的問題，本文結合-EEMD算法和排列熵，提出了一種基于CORS高程時間序列本身的非線性趨勢項定量識別與估計方法。4個CORS站近20 a的高程時間序列實驗結果表明：本文提出的CORS站高程時間序列非線性趨勢項定量識別與估計方法，可以有效地識別并分離非線性趨勢項，得到基于序列自身特性的單調擬合函數。

在研究過程中發現：相同區域的測站具有相似的運動特性，不同地區的測站間運動差異較大。下一步將結合測站所處地理位置對引起測站高程非線性運動的地球物理因素進行定量分析，獲取各測站高程方向上更為準確的運動規律。

[1] 姜衛平. GNSS 基準站網數據處理方法與應用[M]. 武漢: 武漢大學出版社, 2017.

[2] 程鵬飛, 成英燕, 秘金鐘, 等. 2000國家大地坐標系建立的理論與方法[M]. 北京: 測繪出版社, 2014.

[3] 張鵬, 蔣志浩, 秘金鐘, 等. 我國 GPS 跟蹤站數據處理與時間序列特征分析[J]. 武漢大學學報(信息科學版), 2007, 32(3): 251-254.

[4] 劉宗強, 黨亞民, 楊強, 等. 顧及共模誤差的四川連續 GPS 基準站坐標時間序列噪聲分析[J]. 大地測量與地球動力學, 2018, 38(5): 510-515.

[5] 武曙光. 區域 CORS 站坐標時間序列特征分析[D]. 武漢: 武漢大學, 2017.

[6] HUANG N E, WU M L C, LONG S R, et al. A confidence limit for the empirical mode decomposition and Hilbert spectral analysis[J]. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences, 2003, 459(2037): 2317-2345.

[7] WU Z, HUANG N E. A study of the characteristics of white noise using the empirical mode decomposition method[J]. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences, 2004, 460(2046): 1597-1611.

[8] 張恒璟, 齊昕, 文漢江. 基于 K 均值聚類 EEMD 的 CORS 高程時間序列信號分析方法[J]. 大地測量與地球動力學, 2019, 39(10): 1053-1057.

[9] 鄭近德, 程軍圣, 楊宇. 改進的EEMD算法及其應用研究[J]振動與沖擊, 2013, 32(21): 21-26.

[10] 雷亞國, 孔德同, 李乃鵬, 等. 自適應總體平均經驗模式分解及其在行星齒輪箱故障檢測中的應用[J]. 機械工程學報, 2014, 50(3): 64-70.

[11] 李長偉, 雷文平, 龐博, 等. 基于EEMD的相關排列熵的滾動軸承故障診斷[J]. 組合機床與自動化加工技術, 2020(8): 1-4.

[12] BANDT C, POMPE B. Permutation entropy: a natural complexity measure for time series[J]. Physical Review Letters, 2002, 88(17): 1-4.

Pattern identification and estimation of nonlinear trend items in elevation time series of CORS

ZHANG Hengjing, WANG Lei, JIN Zelin, LAN Wenqi, QI Xin

(School of Geomatics, Liaoning Technical University, Fuxin，Liaoning 123000, China)

Aiming at the problem that the traditional linear model fits the Continuously Operating Reference Stations (CORS) elevation time series to make the station movement trend singular, this paper proposes a quantitative identification and estimation method of nonlinear trend items based on the CORS elevation time series itself. This method first uses the-means clustering global empirical mode decomposition method to decompose the CORS elevation time series. Secondly, based on the permutation entropy theory, quantitatively identify the nonlinear trend item, make statistics on the trend item model suitable for decision-making, and estimate the model parameters. The experimental results show that the quantitative identification method of nonlinear trend item based on permutation entropy avoids human discrimination errors and accurately reflects the changing trend of the sequence itself. By comparing the fitting effects of Fourier, sine function, and linear polynomial models, it is found that the Fourier model has the best fitting effect.

elevation time series; trend item; permutation entropy; model fitting; quantitative estimates

P228

A

2095-4999(2022)02-0152-07

張恒璟，王蕾，金澤林，等. CORS高程時間序列非線性趨勢項模式識別與估計[J]. 導航定位學報, 2022, 10(2): 152-158.（ZHANG Hengjing, WANG Lei, JIN Zelin, et al. Pattern identification and estimation of nonlinear trend items in elevation time series of CORS[J]. Journal of Navigation and Positioning, 2022, 10(2): 152-158.）

10.16547/j.cnki.10-1096.20220219.

2021-04-27

張恒璟（1982—），男，河南鄧州人，博士，副教授，研究方向為空間大地測量數據處理與GNSS高程非線性運動。

王蕾（1997—），男，安徽阜陽人，碩士研究生，研究方向為空間大地測量數據處理。