融合波动率的时序数据符号聚合近似方法探究

鲁法明，王 琳，包云霞，李 昂，曾庆田

(1.山东科技大学 计算机科学与工程学院，山东 青岛 266590；2.山东科技大学 数学与系统科学学院，山东 青岛 266590)

0 引言

时间序列是对某个物理量进行等时间间隔观测所得到的数值序列，时间序列分析广泛应用于工业制造、金融、气候监测和医疗等领域。例如，在工业制造领域，采集并分析生产过程中的设备状态时序数据，可以对生产线异常做出预警，提高产品质量和生产效率[1-3]；在经济金融领域，股票的价格走势、利率等信息也构成时间序列，对其进行分析可以在一定程度上预测股票收益或检测股票操纵异常[4-6]；在气候环境领域，对降雨量、河流水位、温度等时间序列进行分析可以为气候应对等问题提供依据[7-9]；在医疗领域，持续监测病人的心电(Electrocardiogram, ECG)和脑电(Electroencephalogram, EEG)活动情况便得到相应的时间序列，对其进行分析可以及时发现病情的异常情况[10-11]；在公共管理领域，对公交系统客流量和车流量时序数据进行分析可以为公交管理和路线调整提供依据[12]。

时间序列数据通常具有维度高和数据量大等特点，直接分析和挖掘原始数据的运算量较大，因此时间序列的降维表示成为业界研究热点之一。目前，比较成熟的时间序列降维表示方法有离散傅里叶变换(Discrete Fourier Transform, DFT)[13]、离散小波变换(Discrete Wavelet Transformation, DWT)[14]、分段线性近似(Piecewise Linear Approximation, PLA)[15]、分段聚合近似(Piecewise Aggregate Approximation, PAA)[16]和符号聚合近似(Symbolic Aggregate approximation, SAX)[17]等。其中，SAX广泛应用于时间序列数据的挖掘与分析，其将时间序列转换为字符串序列，在达到降维效果的同时较好地保持了时间序列的原始信息[16]。

SAX[16-17]将原始时间序列分为若干时间段，每个时间段内的观测值用一个符号聚合表示，从而快速有效地对时间序列降维。降维后的符号序列仍然可以反映时间序列的总体变化态势，且具有快速高效等优点。然而，SAX在时间序列的符号化过程中丢失了原始时序数据的许多统计和形态信息。统计信息主要指极值、方差等信息，形态信息主要包括变化趋势和波动信息等。

为减少统计信息的丢失量，LKHAGVA等[18]提出一种基于极大值、极小值和均值的扩展符号聚合近似(Extension of Symbolic Aggregate approximation, ESAX)方法，该方法同样将每个子序列的最大值和最小值加入符号序列，使原来每个子序列对应的符号变成3个。ESAX弥补了部分丢失的统计信息，强调极值信息的重要性，适用于部分对数据极值敏感的场景。

同样是针对统计信息丢失的问题，钟清流等[19]在SAX的基础上加入方差来描述时间序列的发散程度。该方法将时间序列的符号特征看做矢量，将各个子序列的方差和均值特征分别看做组成该矢量的两个分量，在一定程度上弥补了时序数据在方差方面丢失统计信息的问题，然而当某些应用场景中存在形态特征不同、均值和方差相同的时序数据时，该方法需要进一步优化。

针对SAX在符号化过程中丢失形态特征信息的问题，李海林等[20]提出SF_SAX。该方法在SAX的基础上，采用最小二乘法对子序列进行直线拟合，用直线的斜率近似值表示该子序列的形态特征，然后按文献[19]的方法将均值和斜率表示为符号矢量。SF_SAX能较好地识别时间序列在形态方面的差异，然而当压缩比较大且序列波动较大时，拟合效果需要进一步改进。

为捕捉时间序列中数据值的变化趋势信息，SUN等[21]提出一种基于趋势距离的时间序列符号化(SAX_TD)方法，该方法将每个子序列的起点值和终点值与其均值的差称为趋势因子，通过计算趋势因子之间的距离衡量时间序列趋势的变化，能较好地识别时间序列的变化趋势；与文献[21]类似，同样使用子序列的起始点来构建趋势，季海娟等[22]提出基于始末距离的时间序列符号聚合近似表示方法(SAX_SM)，该方法在计算子序列均值时去掉起始点，而将起始点加入子序列映射的字符序列中；李海林等[23]用数值导数描述子序列的趋势特征，导数为正数表示子序列呈上升趋势，为负数表示子序列呈下降趋势，为0则表示无变化。这3种方法在一定程度上弥补了符号化过程中趋势信息丢失的问题。

前述各种改进的SAX在对时序数据进行降维的同时，分别从不同角度捕捉传统SAX丢失的统计信息和形态信息，各有其适用的场合。本文着眼于弥补时间序列符号化过程中丢失的数值波动信息。具体而言，本文提出一种融合波动信息的时间序列符号聚合近似方法，定义了一种新的波动率指标来同时刻画时间序列的波动幅度和变化趋势等形态特征，用融合了波动率的符号矢量近似刻画子序列，在此基础上给出一种新的时间序列距离度量方法。然后，以该度量方法为基础，提出时间序列的相似性计算和分类方法，并在公开的加利福尼亚大学河滨分校(University of California, Riverside，UCR)时序数据集和凯斯西储大学(Case Westem Reserve University, CWRU)的轴承故障诊断数据集上对该分类学习的效果进行实验，对所提方法的适用场景进行分析。

1 SAX原理概述

SAX是LIN等[16-17]提出的基于PAA的符号化表示方法。其在时间序列数据近似服从正态分布的前提下，使用符号序列聚合表示原始时间序列，在降维的同时保留时间序列的总体变化态势，并保证符号空间相似模式之间的距离满足真实距离的下界要求。

简而言之，SAX将时间序列转换为字符序列主要经过以下3步：

步骤1将原始时间序列规格化，转换成均值为0、标准差为1的序列。

(1)

表1 字符集大小α=3,…,10时各分割点的取值情况

假设某时间序列数据经过规格化处理后的结果如图1所示，图中在128个观测点取值。按照SAX设置分段数量为8，字符集大小为3，首先采用PAA算法得到各段子序列的均值，然后通过符号聚合近似后得到符号序列abcccbaa，该符号序列对应的图形表示如图2所示。对比图1和图2可见，采用SAX进行降维后，数据维度从原始的128个观测值降到只有8个观测值，而时间序列的总体变化态势在降维前后基本保持一致。

(2)

(3)

式(3)中β取值如表1所示。

由上述步骤可见，SAX采用PAA计算出的子序列均值代替这一段子序列，容易丢失原序列部分信息(如方差、变化趋势和波动情况等)，而且数据的压缩比越大，信息丢失得越多。另外，该方法只能保留原时间序列的总体变化趋势，无法描述各段的局部形态信息。针对SAX的不足，下面给出融合波动率的时间序列符号聚合近似方法。波动率从一定程度上同时刻画时间序列的局部波动幅度、变化趋势等特征，可以弥补SAX部分丢失的信息。

2 融合波动率的时间序列符号聚合近似方法

2.1 时间序列的二元符号矢量序列近似表示

时间序列的波动情况同时隐含时间序列的数值和形态方面的特征，这些特征往往与时间序列的状态或模式变化有关，捕捉时序数据的波动信息对时间序列分析和异常检测等应用具有重要作用。观察图1所示的时间序列，可以直观地发现各个子序列局部的状态或模式变化，但SAX处理后仅保留各个子序列的均值，导致丢失了时间序列信息，从图形上看，图2相比图1丢失了很多信息。为更好地刻画时间序列的波动信息，本文在SAX的基础上加入波动率来量化时间序列的波动特征，其定义和计算公式如下：

定义1给定时间间隔0，记xi和xi+1分别为第i时刻和第i+1时刻的观测值，则称式(4)的计算结果为本时间段内的波动率。

(4)

由式(4)可见，当时间序列各个观测点的取值整体呈上升或下降趋势时，波动率的取值较大；当没有明显的上升和下降趋势，但局部波动幅度较大时，波动率的取值也较大。反之，当某个时间段内各个观测点的取值均稳定在某个恒定值附近时，波动率的取值接近零。显然，该式客观地刻画了时间序列的波动情况，后文实验结果也表明这一方法在捕捉时间序列形态信息方面的有效性。

结合上述波动率计算公式，融合波动率之后的时间序列符号化方法称为融合波动率的时间序列符号聚合近似方法(time series Symbol Aggregation approximation method for fusion VOLAtility，SAX_VOLA)，其具体过程如下：

给定长度为n的时间序列C={c1,c2,c3,…,cn}，SAX_VOLA对时间序列进行符号表示时主要有以下4步：

步骤1与SAX相同，首先对原时间序列进行规格化，将其转换成均值为0、标准差为1的序列。

例如，对于图1中的时间序列，仍然设置分段数量为8，字符集大小为3，先采用PAA算法得到各段子序列的均值，再进行符号聚合近似，得到的符号序列仍然为abcccbaa。SAX_VOLA针对每个子序列计算其波动率，并将波动率取值与对应子序列的符号组合为一个二元矢量，得到的二元符号矢量序列为

c·i+0.334·j,c·i+0.375·j,

c·i+0.365·j,b·i+0.426·j,

a·i+0.279·j,a·i+0.246·j}。

2.2 时间序列的距离度量

(5)

(6)

上述时间序列的距离度量方法在传统符号距离的基础上融入波动率距离，而波动率在一定程度上同时刻画了时间序列的波动幅度和变化趋势，因此相比传统的符号化方法，该方法通常能更加准确地度量时间序列的距离。为验证这种距离度量的准确性，后文将从某公开的数据源中寻找含有分类标签的多个数据集，基于该相似性度量指标对时间序列进行分类学习，通过分类的准确性评估本文所提符号聚合方法和时间序列距离度量方法的有效性。

3 实验评估

3.1 UCR数据集实验与分析

UCR时间序列档案库[24]是使用最为广泛的时间序列数据源之一，它含有128个时间序列数据集，具体包括ECG数据集、Trace故障数据集、Wafer传感器数据集等。其中的Trace数据集[25]是流程工业领域核电站仪器故障的综合数据集(完整的数据集一共有16类，此处UCR时间序列档案库收录了Trace数据集中的4类[26])；Wafer数据集[27]是在半导体制造领域的硅晶片加工期间，通过各种传感器收集的线上测量值所构成的数据集，分为正常和异常两类。文献[21-22]以UCR为数据源，对多种时间序列符号聚合近似表示方法在时序数据分类方面的准确性进行对比。为验证本文所提SAX_VOLA的有效性，同时便于不同方法间进行比较，采用文献[21-22]采纳的20个时序数据集进行实验评估。具体选择的UCR中的数据集如表2所示，每个数据集包括训练集和测试集两部分，不同数据集类别标签的数目从2～50不等，时间序列长度从60～637不等，这些数据集在类别数量和序列长度上的多样性可以更客观地对各种符号化方法的适用性做出评估。

表2 选择UCR的20组数据集信息表

续表2

本文借助时间序列数据分类的准确性对比不同符号聚合近似表示方法。在分类器的选择方面，鉴于K最邻近(K-Nearest Neighbor, KNN)分类器[28]原理简单，差错率较低，本文用KNN分类器作时间序列数据的分类算法。考虑到1-NN分类器直接按照训练集中与其距离最近的时间序列类别进行分类，采用KNN分类器中的1-NN分类算法对时间序列进行分类。

就时间序列距离度量而言，分类过程中分别选取原始时间序列的欧式距离，以及SAX,SAX_TD[21],SAX_SM[22]和本文所提SAX_VOLA方法中给出的距离计算公式作为距离度量的依据。其中，SAX_TD用每段子序列的起止点与均值的差构建趋势距离，进而计算出两条时间序列符号化后的距离；SAX_SM直接用子序列之间起止点的差构建趋势距离，进而求得两条时间序列符号化之后的距离。

除了前述不同的时间序列距离计算方法，符号化表示时字符集的大小α和时间序列的分段数w也会对分类的准确率产生影响。为减小这些参数对评估结果的干扰，对各种不同的距离计算方法进行多次实验，设α=3,…,10，w=2,…,n/2(n为时间序列的长度，w每次取值为前一次的2倍)，选择分类准确率最高的结果作为相应符号聚合近似表示的最终实验结果(如果不同参数值取得了相同的分类准确率，则选w较小的参数)。表3所示为不同方法在不同数据集上得到的分类准确率，其中Eucild表示欧式距离，加粗表示各数据集上所取得的最好分类准确率。

表3 各种方法在不同数据集上的分类准确率

续表3

由表3可见，本文所提SAX_VOLA及其对应的时间序列距离度量方法，在19组数据集上比直接用原始时序数据的欧式距离进行分类的效果好，在制造业相关的Trace数据集上提升效果最好，达到17%，在Yoga数据集上提升最低，仅为0.1%；在19组数据集上优于SAX，在OliveOil数据集上提升较大，高达73.3%，在Wafer数据集上提升较小，只有0.1%；在17组数据集上优于SAX_TD，1组数据集上持平，在Trace数据集上提升最高，为18%，在Wafer数据集上提升最小，为2%；在18组数据集上略优于SAX_SM方法，在Trace数据集上提升最高，为22.7%，在Wafer数据集上提升最小，仅1%；在17组数据集上优于所有方法。具体对比如图3所示，可见基于SAX_VOLA的时间序列分类效果在本文所列多数数据集上优于传统的欧氏距离、SAX、SAX_TD和SAX_SM方法。

由图3和表3可见，本文所提方法在Coffee数据集上比其他方法的分类准确率略低。该数据集中时间序列的图像如图4所示(图中横坐标表示时间序列的长度，纵坐标表示时间序列上每个点的对应值)，可见该时间序列局部波动不显著。作为对比，图5和图6给出了SAX_VOLA表现较好的Wafer和Trace数据集中的时间序列图像，其最明显的特点就是局部波动明显，或者有明显的上升或下降趋势，与前面提及的波动率能够有效捕捉子序列的上升/下降趋势信息及波动幅度信息的事实相符。一般而言，当原始时间序列有明显的上升、下降趋势或局部波动幅度较大时，本文所提SAX_VOLA的效果更好。这种现象在制造业领域生产线的监测数据上较常见，因此本文方法在该类情况下会取得更准确的分析结果。

下面分析本算法的时间性能。实验所用计算机的配置为CPU i5-4200M、8 G内存、Windows操作系统。图7所示为前述4种符号化方法在ECG200，GunPoint，OliveOil，Trace数据集上的时间成本(包括符号化降维和相似度计算)对比图，可见本文所提SAX_VOLA的时间成本略高于其他3种方法，原因是本文提出的波动率指标计算复杂度高于始末距离、趋势距离和SAX无附加指标的计算复杂度，这也是为了抽取更多序列信息而付出的代价。然而，随着分段数的增加，由于段内采样点的减少会降低波动率计算的代价，SAX_VOLA与SAX_TD的时间成本差距逐渐缩小。当然，无论本文所提SAX_VOLA还是传统的符号化方法，由于取值点大大减少，在进行时间序列相似度计算和分类等任务时，时间效率均明显优于未符号化时对原始时间序列的处理效率。

3.2 轴承故障诊断数据集实验与分析

CWRU滚动轴承数据中心的轴承故障诊断数据集[29]是世界公认的轴承诊断标准数据集之一，为进一步验证SAX_VOLA在工业领域的应用，本次实验选取CWRU轴承数据中的驱动端数据。被诊断轴承型号为深沟球轴承SKF6205，系统采样频率为12 kHz，负载为1 HP，选取两个周期。被诊断轴承的缺陷位置有滚动体损伤、外圈损伤和内圈损伤3种，损伤直径分别为0.007 inch，0.014 inch，0.02 inch，共9种损伤状态，加上正常状态共计10种状态，每种状态选取100个样本，然后随机挑选30%的数据作为测试集，字符集α和w的设置与各种方法在UCR时间序列档案库中数据集上的设置相同。表4所示为各种方法在轴承故障数据集上的最高分类准确率与取得最高分类准确率所耗费的时间，可见SAX_VOLA取得的分类准确率略优于其他方法，时间成本比略低于SAX_TD。导致这一现象的原因与3.1节相同，为了保证时间序列分类的准确性，需要在符号化过程中尽量抽取原始时间序列更多的信息，因此将耗费更多的计算时间。

表4 各种方法在轴承故障数据集上的分类准确率和耗时

4 结束语

针对制造业等领域的时间序列数据降维问题，本文提出一种新的时间序列符号聚合近似方法，通过引入波动率指标同时量化时间序列的局部波动幅度和变化趋势等信息，弥补了传统SAX在符号化过程中丢失的波动信息。实验结果表明，当时间序列有明显的上升和下降趋势，或者局部存在频繁的波动时，本文方法在时间序列分类问题上的准确率上通常优于传统方法。

然而，在时间序列的降维过程中，本文采用等长分割的办法，可能导致识别出的序列模式不准确，后续工作可以尝试根据波动率对时间序列进行不等长分割，进而更加精确地捕捉时间序列模式和状态的变化信息。另外，时间序列分析已经应用于智能制造[30-31]、业务过程管理[32-34]等诸多领域，如何在这些领域推广和应用本文方法也是下一步研究的重点。