引言
前期介绍了一些常用的突变检测方法,如MK突变检测法、Bernaola Galvan分割算法:Manner-Kendall(M-K)—突变检验;非平稳时间序列突变检测–BG算法。今天介绍一种非参数突变检测方法–Pettitt突变检测方法。下图是本文实现的Pettitt突变检测方法检测得到的突变点以及满足的置信水平,第一栏为原始序列数据及突变点,数据是两段均匀采样,第一段的均值0.9,第二段的均值为0.6,从图中可以看出明显的突变位置,第2栏为统计量Ut及显著性水平,可以看出该突变满足显著性水平p=0.001。
去看原文
原理
实现
在这里,采用Python实现上述原理,并通过测试用例数据进行验证,具体如下:
# -*- coding: utf-8 -*-
# @Author: 武辛
# @Email: geo_data_analysis@163.com
# @Note: 如有疑问,可加微信"wxid-3ccc"
# @All Rights Reserved!
import numpy as np
import pandas as pd
import sys, math
import matplotlib.pyplot as plt
def main():
alpha = 0.05
data = pd.read_csv("data_random.csv")
Year = list(data["Year"])
X = list(data["X"])
Ut, max_idx, Kt = pettitt(X)
print("The %dth sample is change point. Kt=%.2f, Year:%d, p=%.5f" % (max_idx + 1, Kt, Year[max_idx], CalP(Kt, len(X))))
# print(Ut, max_idx, Kt)
if CalP(Kt, len(X)) < alpha:
print("The %dth sample is change point. Year:%d, p=%.9f" % (max_idx + 1, Year[max_idx], CalP(Kt, len(X))))
else:
print("Does not satisfy the significance level! P(Kt) = %.5f" % CalP(Kt, len(X)))
plot(X, Ut, Kt, max_idx)
if __name__ == '__main__':
main()
结果
第一栏为原始序列数据及突变点。第2栏为统计量Ut及显著性水平,可以看出该突变满足显著性水平p=0.001。
原文有源码,更多内容,请关注地学分析与算法。