1. 引言
中国大陆构造环境监测网络(以下简称“陆态网络”)在“十一五”期间由中国地震局、总参测绘局、中国科学院等多部门共同承担建设[1]。在山西省“陆态网络”基准站共九个站。至今已积累了十余年的时间序列,为山西地区的大地测量及地球动力学提供了丰富的基础数据。
基于eemd-spen (样本熵)联合小波阈值去噪方法是一种用于信号降噪的信号处理方法,它结合了集合经验模态分解方法(eemd)、样本熵(spen)和小波阈值处理技术。该方法已应用在各个领域,其中学者张雪应用该方法对语音进行去噪,结果表明该方法缓解分解,模态混叠现象,解决了小波分解中预先设定小波基及分解层数的数据,较好的抑制噪声[2]。集合经验模态分解方法(eemd)是把时间序列分解成m个固有模态分量(imf)和一个趋势项,通过样本熵计算得到每个imf分量样本熵值,选择合适参数做阈值,使用小波阈值去噪方法去除各imf分量小于阈值的小波系数,再重构imf分量,分到去噪的时间序列。
在本论文中,我们将基于集合经验模态分解(eemd)联合样本熵(spen)对山西太原(sxty)的gnss时间序列进行降噪分析。通过对eemd联合spen的研究和应用,我们将能够更准确地分离出gnss时间序列中的噪声和干扰,从而提高gnss数据的精度和可靠性。本研究对于山西地区的gnss数据处理和应用具有重要的意义。
2. 集合经验模态分解(eemd)去噪方法的步骤
集合经验模态分解(eemd)步骤:
1) 进行gnss站时间序列的eemd分析前,须在原始信号基础上添加白噪声,即白噪声(恒定幅值标准差并满足n次取均值结果为零) w(t)加入到带分解的时间序列中y(t)中,由此产生经过噪声处理的新时间序列[3]。公式如下:
(1)
将时间序列在叠加了白噪声的序列实施eemd解算,得到不同频域的imf构件,公式如下:
(2)
式中k为分解尺度,
为第k个独立imf分量,r称为余项。
2) 对得到imf分量进行样本熵计算,得到各分量的值。
第一步,将各分量组成一组维数为m的向量序列,
,其中
,
。这些向量代表从i点开始的m个连续的x的值。
第二步,定义向量
与
之间的距离d [xm(i)]与[xm(j)]为两者对应元素中最大差值的绝对值。即:
(3)
第三步,对于给定的xm(i),统计xm(i)和xm(j)之间距离小于等于r的j (1 ≤ j ≤ n – m, j ≠ i)的数目,并记作bi。对于1 ≤ j ≤ n − m,定义:
。 (4)
b(m)(r)为:
(5)
第四步,增加维数到m 1,计算
和
(
)距离小于等于r的个数,记为ai。
定义为:
(6)
第四步,am(r)为:
(7)
这样,bm(r)是两个序列在相似容限r下匹配m个点的概率,而am(r)是两个序列匹配m 1个点的概率。样本熵定义为:
(8)
当n为有限值时,可以用下式估计:得到各分量的样本熵值
(9)
使用小波阈值去噪,选imf分量进行n层小波分解,得到各层的小波系数,对第一层到第n层的高频小波系数采取合适的阈值进行量化,阈值为利用其去相关性去除噪声。
3) 重构imf,对小波阈值的结果进行时间序列重构,得到最终的去噪站点时间序列。
上述集合经验模态分解(eemd)理得出,分解不同频率的固有模态分量(imf),通过相似性和光滑性指标选出核实模态分量,并进行重构,得到去噪后的时间序列[4]-[7]。
3. 算例实验分析
选取山西太原(sxty)站点时间序列,数据来自中国地震局gnss数据产品服务中心平台,时间跨度为2012至2024。数据在解算过程中已进行物理改正模型(海潮、轨道采用松弛模式;地球重力场等),得到高精度时间序列(图1)。本文使用上述降噪步骤对山西太原(sxty)站n、e、u方向分别做降噪处理,来验证集合经验模态分解方法(eemd)联合样本熵(spen)对山西太原gnss时间序列降噪有效性和可行性[8]。
figure 1. time series of sxty station
图1. sxty时间序列
从山西太原(sxty)时间序列图发现n、e、u三个方向趋势变化,其中n、e方向时间序列呈线性变化趋势,而在u方向时间序列呈周期变化趋势,且存在非线性趋势[9]。
通过集合经验模态分解方法(eemd)联合样本熵(spen)进行时间序列的降噪,原始曲线中人为叠加了白噪声,使用集合经验模态分解方法(eemd)对时间序列进行分解,得依次提取出各个不同频带的固有模态函数(imf)以及剩余成分,具体分解结果如图2展示。
从图2看到固有模态分量imf1为高频信号,随着固有模态分量(imf)序号增加,频率逐渐降低。由图可得到不同频率固有模态分量(imf)关于降噪指标的结果。
采用样本熵法对不同成分做了分析,其结果显示在图3,即第一至第十内模函数(imf1至imf10)的样本熵值大致分布于0.26至1.26的区间内。
验证得到最优样本熵值为0.7,该值作为小波去噪的阈值。通过验证得到最优值为0.7,通过小波阈值法降噪,并重构时间序列,通过信噪比和能量百分比确认最优值,并得到去噪后的时间序列。
(a) n方向
(b) e方向 (c) u方向
figure 2. mode decomposition results of time series components of taiyuan (sxty) station in shanxi province
图2. 山西太原(sxty)时间序列各分量固有模态分解结果图
(a) n方向
(b) e方向 (c) u方向
figure 3. calculation results of sample entropy of shanxi taiyuan (syty)
图3. 山西太原(syty)各分量样本熵计算结果
(a) n方向
(b) e方向
(c) u方向
figure 4. comparison of time series before and after noise reduction of each component of shanxi taiyuan (sxty) time series
图4. 山西太原(sxty)时间序列各分量降噪前后时间序列对比图
图4为山西太原(sxty)各分量降噪前后时间序列对比图,为了进一步集合经验模态分解方法(eemd)联合样本熵(spen)对gnss基准站时间序列降噪的有效性,计算降噪前后序列信噪比的变化和相关系数2个参数,作为评价降噪效果指标[10]。其中信噪比rsn的计算公式如下:
(10)
table 1. the statistical results of signal-to-noise ratio and correlation coefficient of each component
表1. 各分量的信噪比、相关系数的统计结果
站名 |
信噪比rsn/db |
相关系数 |
n |
e |
u |
n |
e |
u |
sxlf |
31.188 |
31.0021 |
35.3098 |
99.9 |
99.9 |
89.6 |
计算结果见表1,n、e、u三份量降噪后的信噪比均较高,最高达35.3098 db;降噪前后序列相关性最高值达99.9%,说明本文方法整体上降噪效果明显。
4. 结论
结论部分只介绍了本文的研究内容,需要把通过本文研究得出的主要结论进行归纳和总结。
本文采用集合经验模态分解方法(eemd)联合样本熵(spen)对山西太原gnss原始位置时间序列作降噪处理,有效分离了原始数据中噪声,得到以下结论:
1) 集合经验模态分解方法(eemd)对时间序列进行分解,得到imf1~imf10的样本值分布在0.26至1.26区间,小波去噪的阈值介于之间,通过验证最优值为0.7。以此值进行小波阈值的去噪。
2) 把小波阈值设为0.7,进行时间序列的降噪,结果表明各分量信噪比有很大提升,n、e方向的前后相关系数高,降噪效果明显。
目前,本文只对山西太原站(sxty),进行时间序列的降噪,下一步可扩展到山西其他站点时间序列数据,通过更多的站点数据进一步验证该方法的有效性。
基金项目
山西省地震局科研项目(sbk-2405)资助。
notes
*通讯作者。