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

CORS高程時(shí)間序列非線性趨勢項(xiàng)模式識(shí)別與估計(jì)

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

（遼寧工程技術(shù)大學(xué) 測繪與地理科學(xué)學(xué)院，遼寧 阜新 123000）

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

高程時(shí)間序列；趨勢項(xiàng)；排列熵；模型擬合；定量估計(jì)

0 引言

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

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

1 K-EEMD與排列熵理論

1.1 K-EEMD算法

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

1.2 排列熵理論

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

2 非線性趨勢項(xiàng)定量識(shí)別與合成

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

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

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

此時(shí)，可以采用人工判讀對各合成項(xiàng)進(jìn)行單調(diào)性判別，若合成趨勢項(xiàng)滿足單調(diào)特性，則將證明該IMF分量為趨勢項(xiàng)；若產(chǎn)生周期運(yùn)動(dòng)趨勢，則該IMF分量為當(dāng)前時(shí)間尺度內(nèi)無法識(shí)別的長周期項(xiàng)。

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

表1 排列熵統(tǒng)計(jì)

當(dāng)排列熵為0時(shí)，原始序列相空間重構(gòu)后排序固定，表示序列穩(wěn)定嚴(yán)格單調(diào)變化，定量地判別了合成項(xiàng)的單調(diào)變化趨勢。從排列熵指標(biāo)不難看出：BJFS站的趨勢項(xiàng)為residual+IMF11+IMF10；SYDN站趨勢項(xiàng)為residual+IMF11+IMF10；JOEN站趨勢項(xiàng)為residual+IMF11+IMF10，HOB2站趨勢項(xiàng)為residual+IMF11。該方法與人工判別結(jié)論一致，同時(shí)能準(zhǔn)確判別出SYDN站IMF10分量為趨勢項(xiàng)的一部分。

排列熵指標(biāo)顯示，SYDN站余項(xiàng)并非嚴(yán)格單調(diào)變化，查看余項(xiàng)數(shù)據(jù)發(fā)現(xiàn)，其確實(shí)不為嚴(yán)格單調(diào)變化。這種現(xiàn)象可能是-EEMD分解的端點(diǎn)效應(yīng)，加入IMF11后，趨勢項(xiàng)嚴(yán)格單調(diào)變化，也證明了SYDN站的趨勢項(xiàng)應(yīng)該為residual+IMF11，而不是單純以余項(xiàng)為該站趨勢項(xiàng)。

3 非線性趨勢項(xiàng)定量估計(jì)

通過定量識(shí)別和合成實(shí)驗(yàn)，得到了經(jīng)過-EEMD分解的各個(gè)CORS非線性趨勢項(xiàng)，結(jié)合趨勢項(xiàng)定義選用三種模型：正弦函數(shù)模型、多項(xiàng)式模型、傅里葉模型對CORS非線性趨勢項(xiàng)進(jìn)行擬合。

1）線性多項(xiàng)式模型。在經(jīng)典的時(shí)間序列擬合模型中，通常認(rèn)為趨勢項(xiàng)是一條直線，所以用線性多項(xiàng)式直接進(jìn)行擬合，選取最小二乘進(jìn)行逼近，求解線性方程的系數(shù)進(jìn)而完成擬合。線性多項(xiàng)式的表達(dá)式為

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

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

采用3種擬合方法對4個(gè)CORS站高程時(shí)間序列擬合結(jié)果如圖3所示，圖3中離散點(diǎn)為SOPAC給出的各站運(yùn)動(dòng)趨勢，實(shí)線為本文識(shí)別的該測站非線性趨勢項(xiàng)，余下3種依次為正弦函數(shù)擬合、線性多項(xiàng)式擬合和傅里葉函數(shù)擬合結(jié)果。

從圖3可知，當(dāng)趨勢項(xiàng)近似線性變化時(shí)，3種擬合模型擬合結(jié)果相近；當(dāng)趨勢項(xiàng)變化為非線性時(shí)，傅里葉擬合模型明顯優(yōu)于其他兩種模型。為驗(yàn)證本文提取趨勢項(xiàng)的正確性，將線性擬合結(jié)果與SOPAC給出的趨勢項(xiàng)參考值進(jìn)行比對，如表2所示。

圖3 趨勢項(xiàng)擬合

表2 線性擬合結(jié)果對比

對4個(gè)測站的非線性趨勢項(xiàng)線性擬合效果進(jìn)行觀察，所有測站提取的擬合結(jié)果均在SOPAC給出的參考范圍內(nèi)，說明本文提出的非線性趨勢項(xiàng)識(shí)別方法和識(shí)別結(jié)果是可靠的。為定量比較3種擬合方案的優(yōu)劣性，采用均方根誤差、擬合優(yōu)度兩種精度指標(biāo)，進(jìn)行擬合精度評定，其評定結(jié)果如表3所示。

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

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

表3 擬合精度統(tǒng)計(jì)

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

4 結(jié)束語

針對傳統(tǒng)線性模型擬合CORS站高程時(shí)間序列使測站運(yùn)動(dòng)趨勢單一化的問題，本文結(jié)合-EEMD算法和排列熵，提出了一種基于CORS高程時(shí)間序列本身的非線性趨勢項(xiàng)定量識(shí)別與估計(jì)方法。4個(gè)CORS站近20 a的高程時(shí)間序列實(shí)驗(yàn)結(jié)果表明：本文提出的CORS站高程時(shí)間序列非線性趨勢項(xiàng)定量識(shí)別與估計(jì)方法，可以有效地識(shí)別并分離非線性趨勢項(xiàng)，得到基于序列自身特性的單調(diào)擬合函數(shù)。

在研究過程中發(fā)現(xiàn)：相同區(qū)域的測站具有相似的運(yùn)動(dòng)特性，不同地區(qū)的測站間運(yùn)動(dòng)差異較大。下一步將結(jié)合測站所處地理位置對引起測站高程非線性運(yùn)動(dòng)的地球物理因素進(jìn)行定量分析，獲取各測站高程方向上更為準(zhǔn)確的運(yùn)動(dòng)規(guī)律。

[1] 姜衛(wèi)平. GNSS 基準(zhǔn)站網(wǎng)數(shù)據(jù)處理方法與應(yīng)用[M]. 武漢: 武漢大學(xué)出版社, 2017.

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

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

[4] 劉宗強(qiáng), 黨亞民, 楊強(qiáng), 等. 顧及共模誤差的四川連續(xù) GPS 基準(zhǔn)站坐標(biāo)時(shí)間序列噪聲分析[J]. 大地測量與地球動(dòng)力學(xué), 2018, 38(5): 510-515.

[5] 武曙光. 區(qū)域 CORS 站坐標(biāo)時(shí)間序列特征分析[D]. 武漢: 武漢大學(xué), 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 高程時(shí)間序列信號(hào)分析方法[J]. 大地測量與地球動(dòng)力學(xué), 2019, 39(10): 1053-1057.

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

[10] 雷亞國, 孔德同, 李乃鵬, 等. 自適應(yīng)總體平均經(jīng)驗(yàn)?zāi)Ｊ椒纸饧捌湓谛行驱X輪箱故障檢測中的應(yīng)用[J]. 機(jī)械工程學(xué)報(bào), 2014, 50(3): 64-70.

[11] 李長偉, 雷文平, 龐博, 等. 基于EEMD的相關(guān)排列熵的滾動(dòng)軸承故障診斷[J]. 組合機(jī)床與自動(dòng)化加工技術(shù), 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高程時(shí)間序列非線性趨勢項(xiàng)模式識(shí)別與估計(jì)[J]. 導(dǎo)航定位學(xué)報(bào), 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—），男，河南鄧州人，博士，副教授，研究方向?yàn)榭臻g大地測量數(shù)據(jù)處理與GNSS高程非線性運(yùn)動(dòng)。

王蕾（1997—），男，安徽阜陽人，碩士研究生，研究方向?yàn)榭臻g大地測量數(shù)據(jù)處理。