标签：index PLR merge 序列 vec error 数据挖掘 data 分段

前言

本篇博客用于记录个人在时间序列数据挖掘中进行的 time series representation 的实践。主要采用 PLR ( piecewise linear representation ) 的方式进行时间序列的降维表达。

时间序列的降维表达对于时间序列的挖掘是非常重要的。因为时间序列通常具有非常高的维度（即数据通常很多）。但同时其又拥有很多的噪声，对于时间序列数据挖掘的下游任务（例如：分类、聚类）来说，原始的庞大时间序列数据不仅数据量大，会影响其任务的处理速度。同时，由于其噪声的存在，还容易影响其任务准确率的表现。所以，对时间序列数据进行降维表达，提取其主要趋势，减少其数据量的同时还能减少噪声数据的影响。

常用的时间序列降维表达的方式有如下类别：

直线近似：用一段段的直线来表达整个时间序列
符号表达：将每一段时间序列都编码成对应的符号
频域表达：将时域数据通过离散小波变换或离散傅里叶变换转换到频域进行表达
模型表达：假设当前的时域数据是由一个模型随机生成出来的，那么我们通过生成出来的数据反过来估计模型的参数，就能得到数据生成的过程。常见的有隐马尔可夫模型(HMMs)等
自编码器（auto-encode）

本文记录的是使用直线近似类别中的分段线性表达（PLR）的实现。

PLR（Piecewise Linear Represent）

概述来说，PLR 就是找到一段段线段来表达整个时间序列数据，线段的端点都对应于原始时间序列的点。

那么如何寻找这一段段线段呢？

分为如下三步：寻找转折点 --> 使用 bottom up 的方法融合转折点 --> 得到最终的分段点。

那么下文将从这三步描述 PLR 的实现

寻找转折点

如何定义关键转折点有非常多的标准。例如：斜率，角度，极值点等。本文参考论文[^1] 提出的使用角度的方法。作者在文章中阐述了为什么使用角度而非斜率。

因为斜率对应的函数是 \(tanx\) 其函数图像是非线性函数，在旋转同样角度（即角度差相等）的情况下，起始角度越大，斜率相差越大。

例如：\((tan60 - tan45) > (tan45 - tan30)\).

使用该方法进行转折点的衡量主要是考察局部相邻三点所构成的两条线段之间的关系。

若两条线段的角度差大于设定的阈值（某一个角度），就将三点中中间点加入转折点序列，留待之后 bottom up 的过程中使用。

这里要注意，虽然每次判定是衡量局部相邻三点所构成的两条线段的关系，但每次向前步进都是以一个点为单位，而非三个点。

实现代码：

def get_TurnPoint(data, angle_thred):
    ''' 
    :paras data:时间序列数据
    :paras angle_thred: 判定其是否为转折点的角度阈值

    :return turnPoints: 预选的转折点
    '''
    turnPoints = []
    for index in range(1, len(data)-1):
        # 这里存在一个隐性的 时间序列采样是否均匀 的问题 
        # 这里如果时间采样不均匀，记得将 atan() 中的数据除以对应时间差
        theta1 = math.atan(data[index] - data[index - 1]) * (180 / math.pi)
        theta2 = math.atan(data[index + 1] - data[index]) * (180 / math.pi)
        if theta1 * theta2 < 0:
            alpha = abs(theta1) + abs(theta2)
        else:
            alpha = abs(theta1 - theta2)
        if alpha > angle_thred:
            turnPoints.append(index)
    turnPoints.insert(0, 0)
    turnPoints.append(len(data)-1)
    return turnPoints

Bottom Up

bottom up 的关键点：

误差阈值：在融合的过程中，用于决定此次融合是否能够进行
误差衡量：融合的过程中，衡量此次融合会造成多大的误差
融合：从起点开始，不断向前汇集，删去中间点，直到分段的误差大于误差阈值。此时就以该分段的结尾作为下一个分段的起点进行下一次融合。

误差衡量的实现：

def merge_cost(data, merge_start, merge_end):
    ''' 
    本方法使用 least-squares solution to a linear matrix equation 方法进行直线拟合。
    同时获得拟合误差。
    调用 numpy.linalg.lstsq() 方法
    :param data: 时间序列数据
    :param merge_start: 合并起始点
    :param merge_end: 合并终点

    :return error: 合并产生的误差
    '''
    x = np.arange(merge_start, merge_end+1) # +1 的原因是 np.arange 为左闭右开
    y = np.array(data[merge_start:merge_end+1])
    A = np.ones((len(x), 2), np.float64)
    A[:,0] = x
    try:
        (result, residual, rank, s) = LA.lstsq(A, y, rcond=None)
        error = math.sqrt(residual[0])
    except np.linalg.LinAlgError:
        error = 0.0
    except IndexError:
        error = 0.0

    return error

代码中使用了最小二乘近似的方法计算拟合直线，同时得到误差。在此处，我想叙述一些有关最小二乘方法的有关内容：

least square

在一开始，我很疑惑如何通过求解线性方程得到直线，在思考之后有了如下答案：

一条直线的定义为 \(y=ax+b\)，此时我们将其以线性方程的方式表达：

\(\vec{y} = [\vec{x}, \vec{1}] \times [a,b]^T\) 即 \([\vec{x}, \vec{1}] \times [a,b]^T = \vec{y}\)

线性方程： \(Az = c\)，将其对应于直线方程的各个部分：

\(A\) 对应于 \([\vec{x}, \vec{1}]\)
\(z\) 对应于 \([a,b]^T\)
\(\vec{y}\) 对应于 \(c\)

那么，当我们拥有直线方程中 \(\vec{x}\) 和 \(\vec{y}\) 两组对应的数据时，想要得到直线方程，其实就是求解 \(Az=c\) 的解 \(z\)，即对应于直线中的系数 \([a,b]^T\)。关于如何使用最小二乘法求解线性方程的答案可以参考[^2]。

那么衍生一下：最小二乘法求解并不是只能做线性拟合。

当我们的矩阵 \(A\) 是以 \(A = [\vec{x}, \vec{x^2}...]\) 的方式构成时，使用最小二乘求解线性方程得到的系数 \([a_0,a_1,a_2,\cdots,a_n]^T\) 就是多项式拟合的系数。即对应如下表达式:

\[y=a_0+a_1x+a_2x^2+\cdots+a_nx^n \]

其实，从此衍生，我们可以理解到，线性组合的真正意义。

就是将空间的几个维度通过线性组合，变为该空间中的一个单独的维度。

例如上面的，我们将 \(\vec{1}, \vec{x}, \vec{x^2}, \cdots, \vec{x^n}\) 当成 n+1 个维度，其线性组合构成了 y 这一个维度。而线性组合的系数决定了新维度是该 n+1 维空间中的哪条直线。

当我们求解线性方程时，就是将已知的数据代入，找到一个最合适的将 n+1 个维度进行组合来表达 y 这一个维度的组合方式（系数）。

话归正题：

bottom up 的实现：

def bottomUp(data, turn_points, error_thred, error_estimator=merge_cost):
    ''' 
    :param data: 时间序列数据
    :param turn_points: 转折点数据
    :param error_thred: 融合过程中的误差阈值, 超过阈值后, 融合停止
    :param error_estimator: 融合过程中的误差估计方法

    :return segments
    '''
    index = 0
    segments = [turn_points[0]] # 加入起始的开始点
    len_turnPoints = len(turn_points) - 1

    while index < len_turnPoints:
        start = index
        end = index + 1
        while error_estimator(data, turn_points[start], turn_points[end]) < error_thred:
            end += 1
            if not (end < len_turnPoints): break
        segments.append(turn_points[end])
        index = end

    return segments

对于 bottom up 而言，有复杂度较低的优点，只需要一次遍历就能完成融合。即 \(O(n)\) 的时间复杂度。虽然代码中使用了两层 while 但是其只遍历了一次转折点序列。

但其缺点是其融合效果依赖于转折点寻找的过程，即其只能对现有转折点进行 merge，但如果转折点寻找过程中没有找到的点，在 bottom up 的过程中是无法找到的。

总结

分段线性表示其实能以较好且较快的方式获得压缩率较高的时间序列降维表达，对于下一步进行时间序列数据挖掘的下游任务有很大的帮助。并且其实现也非常简单。

其实我很想说，在实践的过程中，最小二乘近似的问题困扰了我一会儿。在以前从来没有想通过这个部分的内容，总是很困扰。但在实践中，我找到了\(Az=b\) 与直线方程中的各部分的对应关系。

之前老将这个 \(z\) 当成直线方程中的自变量去思考，苦苦得不到答案。但其实这个 \(z\) 的本质是需要求解的量。而直线方程中需要求解的量显然是系数，而非自变量。我们拥有的数据是 \((x,y)\) 对，要通过这个数据去求解系数。

参考

1 基于转折点和趋势段的时间序列趋势特征提取. 刘意杨等

2 numpy.linalg.lstsq 参考文档

标签：index,PLR,merge,序列,vec,error,数据挖掘,data,分段
From： https://www.cnblogs.com/StephenSpace/p/16900214.html

时间序列数据挖掘之分段线性表示（PLR）

前言