基于成分聚类的高阶奇异谱分析及在GNSS监测序列分析中的应用

翟长治，岳　顺，李小奇

(1.武汉大学 测绘学院,湖北 武汉 430079;2.河海大学 地球科学与工程学院,江苏 南京 210098；3.上海勘测设计研究院有限公司，上海 200434)



基于成分聚类的高阶奇异谱分析及在GNSS监测序列分析中的应用

翟长治1，岳顺2，李小奇3

(1.武汉大学 测绘学院,湖北 武汉 430079;2.河海大学 地球科学与工程学院,江苏 南京 210098；3.上海勘测设计研究院有限公司，上海 200434)

摘要：高阶奇异谱分析(HSSA)相对于奇异谱分析对不同延时、嵌入维数变化有较好的鲁棒性，在信号处理、工程等领域有着广泛应用。针对高阶奇异谱分析应用中成分序列的选取凭经验过于主观的问题，在理论推导与分析的基础上，提出一种基于成分功率谱聚类的方法，并应用于实际GNSS监测序列的处理与分析，得到了具有明显以一天为周期的周期成分、趋势变化成分以及不规则变动成分，反映了不同因素造成的大桥索塔变形，为后续索塔变形的建模预报、安全评估提供了技术支撑，也说明方法的实用性。

关键词：高阶奇异谱;成分聚类;功率谱;GNSS监测序列

随着计算机技术、对地观测技术，以及各种传感器技术的发展，自动化的、实时的监测技术己经成为了可能。其中GNSS系统以精度高、全天候、实时性等优点成为这些网络的重要组成部分。这些GNSS系统在提供其常规服务的同时，也积累了大量包含区域或建筑物形变等重要信息的监测序列。但由于环境因素的复杂性以及各种随机干扰的影响，监测序列表现为各种变形与误差的叠加，因此如何有效分离与提取变形信息，为后续变形分析、预报、安全评估提供技术支持，具有重要的实用价值[1]。常用的分离与提取变形信息的方法有小波变换、现代谱估计、卡尔曼滤波等方法。近年因奇异谱分析具有稳定识别和强化周期信号的优良特性，国内外学者也尝试将其应用于变形序列的分析与变形特征提取[2]。但应用中嵌入维数、延时以及有效成分的确定存在主观性，在一定程度上限制了其应用的效果。高阶奇异谱分析(HSSA)是一种高阶统计方法，与奇异谱分析相比，具有对不同嵌入维数和延时的鲁棒性，并有盲高斯噪声性和抑制相关噪声的特性[3]。但在高阶奇异谱分析应用中同样存在重构成分的选取问题。以往研究者主要依据经验进行分类、选择，存在较强的主观性，从而会影响重构结果的客观性。针对这一问题，本文提出一种根据成分序列功率谱分析结果聚类的重构成分选择法，并对实际的GNSS观测数据进行分析处理。

1高阶奇异谱分析的原理

高阶奇异谱是一种基于高阶累积量的非线性信号的处理方法。高阶奇异谱是将原状态向量的高阶累积量进行切片处理，得到二元的高阶矩函数。然后将二元的高阶矩函数进行SVD分解，得到奇异值和左右奇异向量。高阶奇异谱重构是将原状态向量向奇异向量空间进行投影，分解出与奇异值数量相同的成分序列，然后从中选择部分成分序列进行重构，从而获得感兴趣的有用信息。令xN表示某一时段的观测数据，序列长度为N。令m1,m2,m3,m4表示序列的一、二、三、四阶矩。

(1)

式中：t1,t2,t3表示时延值;E(·)表示数学期望。四阶累积量为

Cum4(t1,t2,t3)=

(2)

由于四阶谱的计算量较大，而时延值t1,t2,t3可以任意给定，为便于进行奇异值分解处理，可以对原高阶矩函数进行切片处理，得到二元高阶矩函数。如令t2=t3，则

Cum4(t1,t2,t3)=

(3)

令Tx(i,j)=Cum4(i,i,j)，便得到高阶矩函数的切片矩阵

Tx=

对Tx进行奇异值分解，求得Tx的奇异值V和奇异向量Ek。将时间序列x1,x2,x3,…,xN中心化处理之后，排列得到共N-M+1个状态，N为序列长度，M为嵌入维数，其中第i个状态量为

(5)

矩阵中的元素与原来的时间序列的对应关系为Xji=xj+i。把式(1)看做M个变量N-M+1个样本。定义原状态向量在奇异向量空间的投影

(6)

也称为时间系数。可以由奇异向量和时间系数来重建x的成分，得到M个x的成分序列。

(7)

在实际应用中，可以根据需要从M个成分序列中选择部分序列进行重构，得到具有不同性质的重构序列[4]。

(8)

其中每个重构序列反映了原状态向量的不同性质。

2基于成分聚类的重构成分选取

应用奇异谱分析进行序列的特征提取最关键的是选择有效的成分序列进行重构。不同成分序列具有不同的特征，代表了原始序列中不同的影响，对于变形监测序列来说，就是代表不同因素引起的不同形变特征[5]。因此针对不同的应用，应选择合适的成分，选择不当，会使重构序列内容混乱，不能有效提取原序列中不同方面的影响，甚至有可能丢失有用信息，引进不必要的噪声[6-8]。针对有效成分的选取，以往的研究大都根据经验确定，具有很大的主观性。根据高阶奇异谱分析的原理，在分解时奇异值按从大到小排列，因此基于每个奇异值重构的成分序列代表原始序列中不同特征的信息，即具有不同的频率与振幅。显然可以根据成分序列的功率谱来分析成分序列，基于功率谱来选择重构成分。对于GNSS监测序列，由于环境影响因素以及随机干扰多而复杂[9-10]，在实际工作中不可能分离每种影响因素引起的变形大小，事实上环境因素也在不断变化中，表现在频率上为一个频率区间，而非一个固定的频率。故此，提出基于成分聚类的重构成分选取，即根据成分频率对成分进行聚类分析，然后根据聚类结果选取成分进行重构，达到分离提取不同特征信号的目的。

首先根据高阶奇异谱分析的原理计算得到每一个成分序列xk,1≤k≤M，采用自相关函数估计法得到各成分的主频率fi。若F看作是频率对于成分序列号的函数

(10)

由于H(i)并不是严格单调递增，i>I时，f(i)仍有可能小于H1，这时可以结合奇异值对成分序列进行取舍。考虑到高阶奇异谱下降比较迅速，当i>I时，其所对应的奇异值一般都非常小，所以可以将其舍去，即只取前i个作为聚类总体。将成分序列截取得到的聚类总体{fi|0

1)计算样品两两距离，得一距离阵记为D0，开始每个样品自成一类；

2)找出D0的非对角线最小元素，设为Dpq,则将Dp和Dq合并成一个新类，记为Dr；

3)给出计算新类与其它类的距离公式：Dkr={Dkp,Dkq}，将D(0)中第p,q行及p,q列用上面公式并成一个新行新列，新行新列对应Gr，所得到的矩阵记为D(1)。

4)对D(1)重复上述对D0的2)、3)两步得D(2)；如此下去，直到所有的元素并成一类为止。

最后根据聚类结果，把属于同一类的成分序列进行重构，得到几个不同的重构序列。

3苏通大桥GPS观测数据的高阶奇异谱分析

本文采用的实际观测数据为苏通大桥索塔顶部的GPS观测数据，GPS系统使用的是徕卡GRX1200Pro系统，采样间隔为30 s。系统共5台GPS接收机，其中4台为监测站，1台为基准站[11]。本文分析了X方向(桥轴方向)的数据，而Y,Z方向的变化性质与X方向相似[12]，可以以此类推。选取其中30 d的观测数据，为了减小数据量，采用每10 min平滑得到一个数据，最后得到分析数据序列长度为4 320，同时为了分析方便，对原数据进行了中心化处理。利用基于成分聚类的高阶奇异谱方法提取大桥索塔变形特征。

按照上述方法对数据进行高阶奇异谱分析。序列长度N=4 320，嵌入维数M=1 500。得到1 500个成分序列。对各成分序列进行功率谱分析，得到各个成分序列的主频率值。

图1显示，第50序列以前，各成分频率变化平稳，主要频率在0.05 Hz以下，而从第50序列以后，频率迅速提高，而且变化无规则。对比奇异值归一化的结果(见图2)，第50序列后的奇异值趋于0。因此本文选取前50个成分序列进行聚类分析。

对前50个成分序列进行最短距离聚类分析，得到聚类冰柱状图，如图3所示。

根据索塔变形特性和有关资料，并结合聚类分析图，把前50个成分可以分成三类较为合适，具体分类见表1。按分类结果进行重构，分别见图4～图6。

表1　最短距离分类

图1　前100个成分序列的主频率

图2　归一化奇异值图

图3　主频率聚类分析结果

图4　周期序列

从图4中可以看出，点位变化在1 d内具有明显的周期特性，称为日周期变形，日周期变形幅值为0.08 m。由于日照方向的影响，阳面的混凝土温度比阴面混凝土温度高，该温差将导致阳面的混凝土膨胀，而阴面的混凝土收缩，使塔柱顶部产生无外力作用下的扭转变形，随着时间推移，日照的方向在周期性地变化，因此这种扭转变形在一日24 h内呈周期性变化。另外，由于各天天气状况不同，温度高低有别，使各周期变形的幅值不同。

图5　趋势项序列

图6　不规则变形序列

图5显示了索塔一个月内的趋势变形，波动范围为[-0.15,0.05]m，总体具有较好的周期性。

图6为不规则变形部分，可见波动范围较小，且无固定周期性，可能是由于受两侧不对称索力、风力和桥面荷载所引起的索塔不规则变形。将三部分合成即反映了一个月内索塔的总体变形，得到图7。

图7　原序列与高阶奇异谱分析后重建序列

从图7可见，在一个月内大桥索塔在X方向的趋势项波动范围为[-0.2,0.2]m,整体上处于周期变化状态，说明索塔具有较好的周期性。原GPS监测序列表面粗糙，有奇异点，对监测分析带来了较大的影响，通过高阶奇异谱分析得到的曲线连续光滑，说明该方法具有显著的滤波效果，对原GPS序列的噪声进行了有效地去除，可以有效地提取索塔的变形特征，摒弃干扰信息。

4结论

本文针对高阶奇异谱分析应用中重构成分的选取问题，提出了一种基于成分功率谱聚类的重构成分选择方法，克服了以往评经验、主观性强的缺点。通过对苏通大桥索塔30 d GNSS实测数据处理与分析，提取除了索塔3种变形：周日变形、趋势变形与不规则变形。由于日照变化造成的日周期变化、天气变化造成的总体变形趋势和风、车辆荷载等因素造成的不规则变形。分析结果与以往大桥索塔变形分析结果基本一致，并在时域和空间域上更加细致，显示了该方法的有效性。重构成分聚类在一定程度上基于经验，仍缺乏成熟的理论支持，有待进一步深入研究。

参考文献：

[1]岳建平,方露,黎昵.变形监测理论与技术研究进展[J].测绘通报,2007(7):1-4.

[2]王解先,连丽珍.沈云中.奇异谱分析在GPS 站坐标监测序列分析中的应用[J].同济大学学报(自然科学版),2013,41(2):282-288.

[3]袁坚,肖先赐.非线性时间序列的高阶奇异谱分析[J].物理学报,1998,47(6):18-26.

[4]杨向锋,张效民,李亚安.一种基于高阶统计量的相空间重构方法[J].信号处理,2006,22(1):5-8.

[5]岳东杰.水利水电工程变形监测中GPS技术与数据处理研究[D].南京：河海大学,2006.

[6]黄大伟,戴吾蛟,文斌.相空间重构单通道ICA算法及其在变形监测数据去噪中的应用[J].测绘工程,2014,23(6):38-42.

[7]李征航,张小红,朱智勤.利用GPS进行高精度变形监测的新模型[J].测绘学报,2002,31(3):206-210.

[8]吴江淮.建筑物静态变形监测数据分析与处理[D].上海：同济大学,2007.

[9]刘成龙,张德强,黄泽纯.大跨悬索桥施工索塔变形成因分析与监测[J].西南交通大学学报,2000,35(5):501-504.

[10] 王化光,黄育龙.鄂东长江公路大桥索塔日照变形监测及应用[J].铁道建筑,2010(3):14-17.

[11] 余波,邱洪兴,王浩，等.苏通大桥结构健康监测系统设计[J].地震工程与工程振动,2009,29(4):170-177.

[12] 翟长治,岳东杰,岳顺.大型桥梁GPS监测数据的奇异谱分析[J].勘察科学技术,2014(4):20-23.

[责任编辑：刘文霞]

Higher order singular spectrum analysis based on component cluster and its applications of analysis of GNSS monitoring sequence

ZHAI Changzhi1,YUE Shun2,LI Xiaoqi3

(1.School of Geodesy and Geomatics,Wuhan University,Wuhan 430079,China;2.School of Earth Science and Engineering,Hohai University,Nanjing 210098,China;3.Shanghai Investigation,Design & Research Institnte Co.,Ltd.,Shanghai,200434,China)

Abstract：Higher Singular Spectrum Analysis (HSSA) has better robutness than singular spectrum analysis in the different delay,and embedding dimension changes.HSSA has been used in wide range of applications such as signal processing,engineering and other fields.For solving the problem of HSSA component sequences selected too subjectively,a new method is proposed based on component of the power spectral clustering,which is applied to the processing and analysis of the actual sequence of GNSS monitoring.The result has been obvious to one day for periodic components,trends and changes in composition of irregular fluctuation component,reflecting the bridge tower deformation caused by different factors.These works provide technical support for the following tower deformation modeling forecasts,safety assessment and also approve the applicability of this method.

Key words：higher-order singular spectrum;component cluster;power spectrum;GNSS monitoring sequence

中图分类号：O433.4

文献标识码：A

文章编号：1006-7949(2016)04-0046-05

作者简介：翟长治(1991-)，男，硕士研究生.

收稿日期：2015-01-15