时间序列分段算法 [Time series Breakout Detection]

标签: 时间序列 算法 time | 发表时间:2015-12-15 12:01 | 作者:
出处:http://www.iteye.com

在时间序列分析中,断点检测(breakout detection)是一个很基本的问题。

通过捕捉时序数据中的断点(breakout),来发现时序数据所表示的系统在过去是否发生了某种事件(event),进而为系统诊断提供必要的数据支持。

 

为了实现对时序断点的检测,我们首先需要对时序的整体时序做拟合。

这里我们通过一条直线来拟合一段时序,如果时序的趋势发生了变化,则用多条直线来拟合整条时序数据。

如下是对一条波动规律明显的时序做拟合之后的结果。

每个红色线条的转折点,就是我们找到的断点。

 

以上数据是我们在实验环境下,为了检测算法效果而人工构造的一条时序。

那么,该算法在实际情况下表现如何?

一下是一条实际的股票价格时序数据。我们通过该算法进行断点检测,并将断点红红色线条连起来的效果:

 

 

 

算法介绍:

算法所使用的关键即使:

1. 单变量线性回归,用来拟合某一段时序

2. 动态规划算法,  用来全局最大化断点检测效果。

 

算法核心代码:

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "lsp.h"

static double loss(double * s, int n){
    int i;
    double t;
    double g0   = 0.0, g1   = 0.0; 
    double h00  = 0.0, h01  = 0.0, h10  = 0.0, h11  = 0.0;
    double hv00 = 0.0, hv01 = 0.0, hv10 = 0.0, hv11 = 0.0;
    double l0, l1;

    // grad and hessian matrix
    for (i = 0; i < n; i++){
        t    = s[i];
        g0  += t;
        g1  += t * (1.0 + i);
        h00 += 1.0;
        h01 += 1.0 + i;
        h11 += (1.0 + i) * (1.0 + i);
    }
    h10 = h01;

    // inverse of hessian
    t = h00 * h11 - h01 * h10;
    hv00 = h11 / t;
    hv01 = hv10 = -h01 / t;
    hv11 = h00 /t;

    // the theta
    l0 = hv00 * g0 + hv01 * g1;
    l1 = hv10 * g0 + hv11 * g1;

    // sqare loss
    t = 0.0;
    for (i = 0; i < n; i++){
        t += (l0 + l1 * (i + 1) - s[i]) * (l0 + l1 * (i + 1) - s[i]);
    }
    return t;
}

int * lsp(double * ts, int n, int min_size, double beta, int *ol){

    if (!ts || min_size < 2 || n < 2 * min_size || !ol){
        return NULL;
    }

    // prev breakout point
    int * prev = (int*)malloc(sizeof(int) * (n + 1));
    memset(prev, 0, sizeof(int) * (n + 1));

    // number of breakout point
    int * num  = (int*)malloc(sizeof(int) * (n + 1));
    memset(num, 0, sizeof(int) * (n + 1));

    // F scores
    double * F = (double*)malloc(sizeof(double) * (n + 1));
    memset(F, 0, sizeof(double) * (n + 1));

    // loss 
    double * lossv = (double*)malloc(sizeof(double) * (n + 1));
    memset(lossv, 0, sizeof(double) * (n + 1));

    for (int s = 2 * min_size; s < n + 1; ++s){
        for (int t = min_size; t < s - min_size + 1; ++t){
            //double ls = loss(ts + prev[t], t - prev[t]);
            double ls = lossv[t];
            double rs = loss(ts + t, s - t);
            double as = loss(ts + prev[t], s - prev[t]);
            double score = (as - ls - rs) * (t - prev[t]) * (s - t) /    \
                           ((s - prev[t]) * (s - prev[t])) - num[t] * beta;
            score += F[t];
            if (score > F[s]){
                num[s] = num[t] + 1;
                F[s] = score;
                prev[s] = t;
                lossv[s] = rs;
            }
        }
    }

    int k = num[n];
    *ol = k;
    int * re = (int*)malloc(sizeof(int) * k);
    memset(re, 0, sizeof(int) * k);
    int i = n;
    while(i > 0){
        if (prev[i])
            re[--k] = prev[i];
        i = prev[i];
    }

    free(prev);  prev  = NULL;
    free(num);   num   = NULL;
    free(F);     F     = NULL;
    free(lossv); lossv = NULL;
    return re;
}

 

算法复杂度上限为:O(n * n * n), 如果时序中的断点明显,切规律则算法复杂度接近 O(n * n)。

 

 



已有 0 人发表留言,猛击->> 这里<<-参与讨论


ITeye推荐



相关 [时间序列 算法 time] 推荐:

时间序列分段算法 [Time series Breakout Detection]

- - ITeye博客
在时间序列分析中,断点检测(breakout detection)是一个很基本的问题. 通过捕捉时序数据中的断点(breakout),来发现时序数据所表示的系统在过去是否发生了某种事件(event),进而为系统诊断提供必要的数据支持. 为了实现对时序断点的检测,我们首先需要对时序的整体时序做拟合. 这里我们通过一条直线来拟合一段时序,如果时序的趋势发生了变化,则用多条直线来拟合整条时序数据.

Long time no see,英式中文

- cow - 刁民公園
昨天「星期日檔案」探討港式英語. 節目中好幾個嘉賓包括楊鐵樑先生都說long time no see是港式英語. 楊官的英文程度不會令人懷疑,但long time no see有其背景,不能算作港式英語.

妄谈时间序列表格型大数据系统设计

- - Solrex Shuffling
一直在特定领域的分布式系统一线摸爬滚打,曾取得一些微不足道的成绩,也犯过一些相当低级的错误. 回头一看,每一个成绩和错误都是醉人的一课,让我在兴奋和懊恼的沉迷中成长. 自己是个幸运儿,作为一个 freshman 就能够有机会承担许多 old guy 才能够有的职责. 战战兢兢、如履薄冰的同时,在一线的实作和思考也让我获得了一些珍贵的经验,却直至今日才够胆量写出来一晒.

自建基于 Linux 的 Time Capsule 服务

- Class4 - Gracecode.com
2011-02-14 更改 netatalk 部分配置参数. 2011-02-10 初始化版本. Apple 的 Time Machine 提供了非常完备的数据备份方案,但每次使用外接磁盘备份很 麻烦也很容易忘记备份. 同时,个人对于 Time Capsule “高昂的”售价(虽然它的确值那些钱)颇有微词.

Google Analytics Real-Time 实时报告开测

- Ray - 谷奥——探寻谷歌的奥秘
先进的互联网越来越强调实时性,Google Analytics终于也开始了这方面的尝试. 以前我们必须等到第二天才能看到第一天的统计数据,今天随着Google Analytics Real-Time的上线,实时的统计数据报告已经不是梦了. 实时统计对于衡量社交媒体对网站的贡献率是非常重要的,在博文发出之后,你即可实时看到各个社交网络引入的流量:.

夏娃的时间:Just enjoy the time of EVE

- Eastar Lee - ACG批评
08年的《夏娃的时间》作为吉浦康裕的第三部动画作品,也是最长的一部作品,在剧情主旨和背景设定上都带有其出道作《水之语》的很多影子. 人类和机器人之间关系,可以说是SF作品中永恒的主题之一,无数电影动漫小说都对这个主题进行过探讨,要把这个老调常弹的主题说出新意并不是一件容易的事情,但是《夏娃的时间》做到了.

Popcorn Time – 开源电影播放器

- - 小众软件
Popcorn Time 是一个很火爆的开源电影播放器,无需完全下载就可以直接播放 yts 里的高清电影,并且能够加载外挂字幕文件. 由于有版权争论, Popcorn Time 以开源的方式让其不会被版权机构干掉,不过…不知道老外是否想到了天朝早已有类似服务的商业公司. Popcorn Time 已经提供了 Windows、OS X、Linux 三个平台的客户端,无广告,有搜索、分类.

『Time』一周摄影图片精选:June 25 - July 08, 2011

- 翼翔 - 『Leica中文摄影杂志』
链接: 《Times》新闻摄影图片精选     新闻摄影@Leica.org.cn. 美国佛罗里达州肯尼迪航天中心,“亚特兰蒂斯”号航天飞机将于8日发射,这是“亚特兰蒂斯”号最后一次执行太空任务. The shuttle atlantis on mission sts-135 sits on launch pad 39 after the rotating service structure has been rolled back at kennedy space center, cape canaveral, florida, usa..

『Time』一周摄影图片精选:July 30 - August 05,2011

- L - Leica中文摄影杂志
链接: 《Times》新闻摄影图片精选     新闻摄影@Leica.org.cn. 一名学生在斋月晨祷阅读“古兰经”,印度尼西亚. A student read the Quran before morning prayers in the holy month of Ramadan at the Al-Mukmin Islamic boarding school in Solo, Indonesia, Tuesday..

『Time』一周摄影图片精选:August 06 - 13, 2011

- Roger - Leica中文摄影杂志
链接: 《Times》新闻摄影图片精选     新闻摄影@Leica.org.cn. 一名摘葡萄工在倾倒整篮的白葡萄,法国Domaine de Rombeau葡萄庄园. 摄影师:Guillaume Horcajuelo. 大雨的积水导致学校被淹,两名小学生在水中搬运自己的桌椅,印度. 摄影师:Mukesh Gupta.