时间序列分析之结构突变点检测

点击蓝字 关注我们

 HAPPY NEW YEAR 

结构突变现象

在过去的几十年里,时间序列的结构性改变问题一直是统计学和计量经济学中的研究热点,这种突变可能表示状态之间发生的转换,我们将时间序列的这些部分称为时间序列的状态,或控制过程中参数不变的时间段。两个连续的不同状态通过一个变化点来区分。变点检测的目标是通过发现变化点来识别这些状态边界。变点检测在时间序列的建模和预测中非常有用,在医疗状况监测、气候变化检测、语音和图像分析以及人类活动分析等领域中也有广泛的应用。

01. 变点理论

变点理论是统计学中的一个经典分支,其基本定义是在一个序列或过程中,当某个统计特性(分布类型、分布参数)在某时间点受系统性因素而非偶然性因素影响发生变化,我们就称该时间点为结构突变点。结构突变点简称变点(Changepoint),也可以称为断点(Breakpoint),是指“模型中某个或某些量突然发生变化的点”,即在序列或过程中的某个未知时刻,某些特性出现了变化,这些特性可以是分布或分布参数,也可能是数字特征。具体定义如下:

假设存在一个数据集,每个数据观测值相互独立,如果在某一时刻,模型中的某个或某些变量突然发生了变化,即存在一个时间点,在该点之前,数据集符合一个分布,在该点之后,数据集符合另外一个分布,则该点为该数据集的变点。

变点识别即利用一定的统计指标或统计方法,对时间序列的状态进行观测,以便准确有效的估计出变点的位置。变点问题作为统计学中的一个重要课题,最早应用于工业质量控制,后延展到金融经济、计算机、气象学、流行病学等多个领域。

从监测内容来看,分为单变点研究和多变点研究。目前大多数研究集中于单变点区域,即假设研究的时序数据中至多只有一个变点。但实际情况一般都存在多个变点,由此衍生了多变点的研究问题,但在研究阶段,大部分学者都是先确定变点数目,再探究变点位置,由此研究方法大打折扣。

从监测目的来看,分为事中变点(连续抽样)和事后变点(固定样本),前者指连续观察某一随机过程,监测到变点时停止检验,不运用到未来数据,主要用于事件预警;后者从已获得的时序数列中检测过去的变点位置,主要用作历史检验。

02. 离线变点检测

离线(Offline)变点检测是指从已经获得的时序数列中检测过去的变点位置,主要用作历史检验。接下来介绍几种常见的离线变点检测算法。

2.1

 控制图方法 

传统的变点检测大多基于统计原理,包括最小二乘法、极大似然法、贝叶斯方法等。随着统计控制过程(SPC)的兴起,控制图技术迅速发展开来,并广泛应用于工业生产之中。控制图即运用典型的数理统计方法,判断产品是否偏离典型分布,甄别是否存在异动。

目前在统计过程控制中最成熟的是三大控制图:休哈特控制图(Shewhart)、累积和控制图(CUSUM)以及指数加权滑动平均控制图(EWMA)。其中最传统的休哈特控制图(mu-3sigma,mu+3sigma)在金融投资中的应用是布林带,EWMA则同指数加权均线比较接近。

CUSUM图作为工业应用检验变点的三大控制图之一,其利用假设检验和极大似然估计的相关统计原理,构建累积和统计量,不断累积观察值与基线水平的差值,放大观察数据出现的波动,从而更加迅速敏感地探测到微小的异常情况,检验出变点位置,其最大的特点是对系统性变化的敏感性,不需要积累太多的样本,因而能较好的控制风险。

由CUSUM推导过程,得到CUSUM统计量,即

其中X为金融序列,令其中k为允偏量。上下预警指标分别为若预警指标分别达到上下阈值h,则对应地产生上变点和下变点。 

选取2000个随机数,其中前1000个服从标准正态分布,后1000个服从均值偏离为0.2的,方差为0的正态分布。在序列分布已知的情况下,我们需要将分布开始变化的位置即变点估计出来。

在此我们设定允偏量k=0.1,阈值h=20,根据上述CUSUM统计量的计算得到得到该序列的上CUSUM统计量,并观察该累计值是否超过阈值。在不断累积过程中,CUSUM控制图在不到1100的位置将这微小的差距识别了出来,效果相对还比较明显,一般统计方法很难甄别如此微小的变化。

· 优点:方法简单,较容易实施,业界研究与应用比较多。· 缺点:需要调参,确定合适的阈值,太小模型太敏感,太大模型太粗糙。

2.2

 邹氏检验 

邹氏检验法(Chow test)是由邹至庄提出的,用于判断结构在预先给定的时点是否发生了变化的一种方法。邹氏转折点检验的目的:检验整个样本的各子样本中模型的系数是否相等。如果模型在不同的子样本中模型的系数不同说明该模型中存在转折点。

检验方法如下:

(1)建立多元线性回归模型:

(2)检验原假设:(3)对以上模型用OLS进行估计得回归方程:

其中,残差平方和分别为和.(4)若原假设成立,则两个回归模型可合并为一个,两组样本观测值可合并成一组样本观测值,回归模型及回归方程如下:

其中,残差平方和为.(5)则Chow统计量为:

在H0成立条件下,F统计量服从自由度为(k, n1+n2-2k)的F分布,其中k为估计参数的个数。给定显著水平α查第一自由度k,第二自由度n1+n2-2k的F分布表的临界值Fα,当F>Fα,拒绝H0,认为两个子样本所反映的经济关系显著不同,经济结构发生了变化;反之认为经济结构关系比较稳定。

· 优点:方法简单,借助数理统计知识有较好的理论基础。

· 缺点:检验时需要事先知道突变点位置且检测出的变点位置较多,不具有区分度。

2.3

 直接密度比估计 

由于前面方法依赖于预先设计的参数模型,并且在现实世界中变化点检测场景中灵活性较低,因此最近的一些研究通过直接估计概率密度的比率而引入了更灵活的非参数变化,而无需执行密度估计。这种密度比估计思想的基本原理是,知道两种密度意味着知道密度比。然而,知道比率并不一定意味着知道这两种密度,因为这种分解并不是唯一的。

因此,直接密度比估计比密度估计简单得多。根据这一想法,研究者设计了直接密度比估算方法。这些方法通过非参数高斯核模型对两个子区间和之间的密度比进行建模,如下所示:

其中,和表示和的样本概率分布,表示f-散度:

显然,两分布之间相异性度量值越高,该点越有可能是变点。根据f-散度计算方法的不同,可分为以下三类变点检测算法:

(1)KLIEP

相异性度量的一个流行选择是Kullback-Leibler(KL)散度:

该问题是一个凸优化问题,因此可以简单地获得唯一的全局最优解θ,例如通过梯度投影法。投影梯度下降在每一步沿负梯度方向移动,并投影到可行参数上,得到KL散度近似值为:

(2)uLSIF另一种直接密度比估计器是uLSIF(无约束最小二乘重要拟合),它使用Pearson(PE)散度作为相异度量,如下所示:

作为uLSIF训练标准的一部分,密度比模型拟合到平方损失下的真实密度比。PE散度的近似值如下:

(3)RuLSIF根据分母密度p(Y)的条件,密度比值可以是无界的。为了克服这个问题,考虑0≤ α<1上的α-相对PE散度用作相异性度量,称为相对uLSIF(RuLSIF)方法:利用α-相对密度比的估计量g(Y), α-相对PE散度可近似为

在随机生成的数据集上进行实证,检测到的变点结果如下:

03. 在线变点检测

在线(Online)变点检测是指连续观察某一随机过程,监测到变点时停止检验,不运用到未来数据,主要用于事件预警。在现实生活中,在线检测过程更有意义,但关于变点的在线检测算法还不是很多,下面介绍一种最成熟的变点在线检测算法。3.1 贝叶斯在线变点检测—Bayesian Online Changepoint Detection

对于变点的在线检测研究,较为常用的是贝叶斯变点检测算法。与回顾性分割不同,此算法关注的是因果预测滤波;只考虑已经观察到的数据,生成序列中下一个不可见数据的精确分布。

假设一系列的观察值可以划分为不重叠的乘积分区。分区之间的划分称为变更点。进一步假设对于每个分区ρ,其中的数据是来自某个概率分布,每个,符合i.i.d条件。我们将时间a和时间b之间的连续观测集合表示为。在变点之间的间隔的先验概率分布表示为。

根据给定的运行长度,利用边缘概率密度求和,可知:

后验概率为

因此,

预测分布仅取决于最近的数据x。

假设每一序列之前的概率值是已知的,那么的概率取决于两种可能性,即变化点发生或变化点没有发生,行程长度或。对两种可能性进行判断并计算在这两种情况下的概率,需要计算风险函数(Hazard Function):

值得注意的两点:

如果是时间尺度为λ的离散指数(几何)分布,过程是无记忆的,风险函数是常数。

是先验。

算法如下:

Matlab程序如下:

function R=bocd_cal(data, mean0, prec0, hazard, T)mean_params = [mean0];prec_params = [prec0];R = zeros(T + 1, T + 1);R(1,1) = 1;message = [1]; for t = 2: T+1 %# 2. Observe new datum.x = data(t-1);%# 3. Evaluate predictive probabilities.for ii = 1: t-1pis(ii) = normpdf (x,mean_params(ii), 1/prec_params(ii) + 1);end %# 4. Calculate growth probabilities.growth_probs = pis .* message * (1 - hazard);%# 5. Calculate changepoint probabilities.cp_prob = sum(pis .* message * hazard);%# 6. Calculate evidencenew_joint = [cp_prob, growth_probs];%# 7. Determine run length distribution.R(t, 1: t) = new_joint;evidence = sum(new_joint);R(t, 🙂 =R(t, :)/evidence;%# 8. Update sufficient statistics.offsets = 2: t+1;new_mean_params = (mean_params .* offsets + x) ./ (offsets + 1);new_prec_params = prec_params + 1;mean_params = [mean0, new_mean_params];prec_params = [prec0,new_prec_params];%# Setup message passing.message = new_joint;endend

检测结果如下:

 上图(顶部)从具有固定方差和变化均值的高斯分布生成的合成数据,变化点用红色虚线表示。(底部)使用对数色标表示每个时间点的RL后验概率,颜色越深表示概率越高,锯齿状红线是每个时间点的最大概率。

参考文献:

[1] Aminikhanghahi S,  Cook D J. A survey of methods for time series change point detection[J]. Knowledge and Information Systems, 2017, 51(2):339-367.

[2] Liu S, Yamada M, Collier N, Sugiyama M. Change-point detection in time-series data by relative density-ratio estimation[J]. Neural Networks, 2013, 43(1):72-83.

[3]

[4]#fearnhead2007line.

[5]

了解更多概率论与数理统计的相关知识,

和我们一起加入网课学习吧!选课网址:https://www.icourse163.org/learn/USTB-1003768006

 HAPPY NEW YEAR 

关注我们 大数学极客

一周一个概率知识,让数学更有趣