小波分析在GPS时间序列分析中的运用

董 伟

(中铁第四勘察设计院集团有限公司，湖北 武汉 430063)

1 概述

GPS已被广泛运用于地壳形变的研究中。将GPS资料制成时间序列，通过分析GPS的连续观测资料便可以得到测站点与时间及空间的关系，并推估出各种构造，从而获得与地壳形变之间的关联性［1］。

目前常用的GPS时间序列分析方法对不平稳的时间序列的分析无能为力，而小波分析能有效地分析地学现象中不平稳的时间序列。根据小波的多尺度特点，可以将信号进行分解和重构，对每一层进行处理，从而获得所需要的部分。

本文在小波分析的基础上对GPS时间序列进行了去除观测噪声的研究。

2 基本理论

小波分析的基础是进行小波变换，而小波变换就是通过小波函数系去表示或逼近一个信号，可由基小波函数通过平移和伸缩构成。

基小波函数可设为ψ(t)∈L2(R)，其中L2(R)表示平方可积的实数空间，小波变换基底则可采用下式进行定义:

式中，a为尺度因子，b为平移因子。

对于给定的能量有限信号或函数 f(t)∈L2(R)，其连续小波变换可得到该信号或者函数与小波函数的内积:

式中，Wf(a，b)为小波变换系数的复共轭函数。而对于离散的数据，只需要将上式的积分形式改成求和形式即可。该变换过程在实际运用中主要是为了获得小波系数，从而通过利用这些系数来分析时间序列的时频变化特征。

本文在进行数据处理时所采用的是morlet小波函数，在实际计算过程中可以通过选择一套离散化的尺度变量来计算小波变换，选取尺度参数原则上可以通过下式确定［2］:

式中，a0为可分解的最小的尺度，δj为离散尺度之间的间距，J为决定最大尺度的因子，N为数据的个数，其中a0一般选择两倍的采样间隔，即a0=2δt，选取较小的δj可以得到更好的尺度分辨率，但计算和画图会变得非常慢，对于morlet小波，本文将选择δj=0.125。

由于GPS时间序列是一维数据，根据Mallat算法，如果已知双尺度方程中的滤波器系数，就可以快速计算出各尺度的逼近和细节。

3 实例分析

本文使用的GPS时间序列为“陆态网络”中KKN4站时间序列，利用小波分析GPS时间序列的主要工作流程见下图1。

图1 GPS时间序列小波分析流程图

由于坐标时间序列的小波分析要求原始坐标时序具有零均值的特性，因此需要对部分原始的GPS时间序列进行消除线性趋势项和常数项的处理，利用消除线性趋势项和常数项得到小波方差图和小波功率谱图，从而可以判断出数据中的周期性质。

在对KKN4站原始时间序列数据进行分析时，从图2可以看出KKN4站原始时间序列中竖直方向有着较为明显的线性趋势，在去除其线性趋势项后的时间序列变化量在零的上下波动，且分布较为均匀，已经具备零均值特性，可用于后续小波分析。

在选取KKN4站2010～2013年的数据进行小波分析时首先获得小波方差图和小波功率谱，由于小波函数通常是一个复数，所以一个实数序列经过小波变换以后可以得到的小波系数也是复数，再利用实部、虚部可以得到小波功率谱，小波方差和功率谱的具体获得过程在此不作详细介绍。

通过小波方差图分析可以确定信号中不同的尺度扰动相对强度及其存在的主要时间尺度，图3为得到的小波方差图和功率谱图。从小波方差图可以明显看出，大致在1/4，1/2、1和2处出现波峰，这与GPS时间序列中存在的季节性、半周年、周年和两周年周期相吻合，这也说明GPS时间序列中确实存在着多个不同时间尺度的周期，同时从小波功率谱图也可以明显看出其周期性。

根据图3的分析结果，确定将GPS时间序列进行4层小波分解，并确定阈值对每一层的高频部分进行去噪处理，之后再对时间序列进行重构。

图2 KKN4站原始时间序列和去趋势项时间序列

图3 KKN4站GPS时间序列小波方差和功率谱

根据小波的多尺度特点，可以将信号进行分解和重构，即用小波将信号进行分解，然后对每一层进行相关处理，最后重构信号。在对GPS的时间序列数据进行分析时，需要考虑长周期运动、同震变形以及震后变形。GPS单站、单分量坐标序列一般用下列周期性模型来表示:

式中，ti为时间，以年为单位;a为地壳位置;b为线性变换率;c、d、e、f为周年和半周年运动振幅;gj为地震造成的同震偏移;H(ti)为阶跃函数;vi为残差值，代表观测值与预测值间的差异，也就是包含GPS时间序列的噪声信号。图4为KKN4站2010～2013年间原始序列及进行4层分解后的高频、低频序列。

图4 GPS时间序列的原始序列和4层分解的高频、低频序列

本文分步采用软阈值和硬阈值两种方法进行去噪处理，然后根据信噪比大小来确定去噪的好坏。

硬阈值处理方法:

软阈值处理方法:

式中，|x|为小波变换的系数，T为预先选定的阈值，考虑到T值的选取会有不同，本文利用GPS数据解算的3倍中误差来限定，设定T=3σ，σ为中误差。

由于影响小波去噪处理效果的因素非常多，选择不一样小波基函数、不一样的阈值以及不一样的分解尺度，最终的去噪效果会有一定的差异，所以在进行降噪处理时必须以一些详细的指标来衡量去噪效果。

本文中所采用的评价方法是信噪比。信噪比是衡量信号中噪声量度的传统方法，其定义式［3］为:

式中，SNR为信噪比，Ps为原始信号功率，Pz为噪声功率。

本文利用软阈值和硬阈值两种方法进行去噪处理及信号重构后信号机信噪比对比图见图5［4－5］。

图5 GPS时间序列与硬阈值法、软阈值法去噪对比

从图5可以非常清楚地看到软阈值法去噪结果比硬阈值法去噪结果更加平滑，采用前者比采用后者进行处理的信号的信噪比小，但采用前者进行降噪处理的同时也掩盖了某些由于闪烁噪声和含随机漫步噪声而导致的点位突变，滤掉了部分变化较大的坐标值。虽然采用软阈值进行降噪处理便于分析基准站点的运动趋势，但缺少对基准站点整体稳定性分析，因此在具体运用中需要根据实际情况选择合适的软阈值或硬阈值降噪，才能获得更好的结果。

4 结语

在利用小波进行GPS时间序列分析中，小波分析可以对信号进行不同尺度的分析，而且还可以将不同特效的噪声非常有效地进行分离。另外小波分析还可以通过对GPS时间序列采用高低频分析相结合的方式进行数据处理，能有效地提取季节性、半周年、周年及两周年的周期项，再通过对信号进行降噪处理，从而可以更好地根据实际需要来分析GPS时间序列。

本文将GPS时间序列进行多尺度分析，尽量保持了边缘位置精确，同时在去噪方面做到了折中处理。分析处理结果表明，小波分析可以很好地用于GPS时间序列的分析中，但由于处理方式不同，去噪处理最终的结果也会不同，因此现场实际使用过程中需要根据具体情况选取合适的处理方式。至于如何深入地进行多尺度分解、去噪及重构，还需要进一步完善。

［1］李婧.“陆态网络”基准站坐标时间序列变化特性分析［D］.郑州:解放军信息工程大学，2013.

［2］Torrence C，Compo G P.A practical guide to wavelet analysis.［J］.Bulletin of the American Meteorological Society，1998，79(1):61－78.

［3］陈强，黄声享，王韦.小波去噪效果评价的另一指标［J］.测绘信息与工程，2008，(5):13－14.

［4］周亚，王立峰，张思慧，等.IGS连续运行参考站高程时间序列功率谱分析［J］.太赫兹科学与电子信息学报，2014，(1):103 －107.

［5］范朋飞.高精度GPS站点坐标时间序列分析与应用［D］.西安:长安大学，2013.