基于MOSUM的多重滤波变点检测研究

杨 超，胡 尧，b，李 扬

(贵州大学 a.数学与统计学院；b.贵州省公共大数据重点实验室，贵州 贵阳 550025)

一、引 言

时间序列中结构性突变的统计检验和估计在各种应用中都有广泛涉及，如计量经济学、移动通信、机器学习、海洋工程或神经生理学数据分析等方面[1]。在变点分析方面，近年来，学者们也做了大量的探究[2-5]。

本研究是在一般分布假设下进行，即允许其它参数变化的情况下，时间序列中存在未知个数的多时间尺度的均值变点模型。针对均值变点的研究，Horvth等提出了非参数检测方法[6-9]，进一步邓春霞研究了存在均值和方差双重变点的非参数检测方法[10]，这些方法虽然避免了参数方法中的参数估计难题，但存在缺点是需要对模型有较强的分布假设。针对多时间尺度上的变点研究，Frick等在检测方法方面做了相关的探究[11-13]，其中一些方法只需相对较弱的分布假设。然而当目标函数为检测均值变点时，其方法需要规定其它参数变化具有鲁棒性，这使得与实际并不相符。Pein等提出了一种利用似然比统计量对高斯序列中存在多时间尺度变点估计的检测方法，其方法检测均值变点的同时允许方差变化[14]。然而一般情况下，在均值对方差变化具有鲁棒性的情况下，方法检测结果并不理想。在Arlot等的研究中，提出了运用删一交叉验证和分割的检测方法，这些方法对不规则模型的参数变化具有较好的鲁棒性[15-16]。然而，在现有的文献中，还没有学者研究一种结合这三种特性的方法，即对数据有较弱的分布假设，对多时间尺度的变点检测，对模型外的其它参数变化具有较强鲁棒性。

基于此，论文考虑这样一种方法(Multiple Filter Test，MFT)，可检测均值不变的零假设或存在均值变点的备择假设。论文将为此过程设计变点检测框架，即将Messer 等人中的方法变换为基于i.i.d.随机变量的分段序列模型[17]，这种变换需要构造一个基于MOSUM型统计量的变点模型，其中也包含参数估计过程[9]。在均值不变的零假设下，并对极限过程进行推导，得出拒绝域的计算方法。而在存在均值变点的备择假设下，基于零假设下得出的拒绝域值，应用多个MOSUM过程，给出其多时间尺度的变点检测规则，最终给出变点检测算法。

针对以上探究，论文将围绕以下三个方面进行。首先，对模型提出假设条件，推导出MFT方法，同时也给出参数估计过程。此外，由模型可同时应用多个窗宽的MOSUM过程，在拒绝零假设的情况下，结合不同窗宽的MOSUM过程，提出一个估计多个时间尺度变点的算法。其算法思想为结合小窗宽的优势，可更精确地检测较短时间范围内发生的变点，结合大窗宽可对发生较小变化的变点进行更准确的识别。其次，在无均值变化的零假设情况下，探究存在其它分布参数变化对MFT方法是否具有鲁棒性。

最后，为验证该方法的有效性和实效性，论文模拟仿真MFT方法对均值变化的灵敏性以及它对其它参数变化的鲁棒性，并最后通过实证分析验证方法的实效性。

二、模型及假设

设X=(Xi)i=1，2…为i.i.d.的随机变量序列，无均值变点的零假设为μ:=E[X1]。给定具有分段均值平稳组合的模型M，假设观测值时刻点为1，2，…，3，T，T∈/{0，1}。这里假定存在变点集合C⊂{2，3，…，T}，其元素c1

X1，[1]，…，Xc1-1，[1]，Xc1，[2]，…，Xc2-1，[2]，…，Xck-1，[k+1]，XT，[k+1]

即在每一个变点cj后，都将伴随一个新均值μ[j+1]的过程X[j+1]，具体过程见图1。

从图1中可看出，图1(a)为无均值变点的随机变量时间序列，Xi～N(0，1)，T=2 000；图1(b)为存在有三个均值变点时间序列X∈M，其变点c1=500(此处均值变化较小)，而变点c2=1 200和c3=1 290(此处均值变化迅速且变化较大)，Xi～N(μ，1)，均值为μ[1]=0，μ[2]=0.3，μ[1]=2.2，μ[4]=1.4，T=2 000。当序列均值恒定时，则需要检验无变点的零假设，即使得C=Φ(如图1(a))。在备择假设下，则要估计集合C(如图1(b)中的三个变点)。注意，为了简单起见，假设随机变量是i.i.d.的，后续提出的检测方法可对这个假设条件放宽[1]。为此，基于模型构造一种均值变点检测方法。

(a)均值无变化序列

(b)均值有变化序列图1 零假设和备择假设下的时间序列图

三、检测方法探究

(一)多重过滤测试(MFT)

首先考虑无均值变点情形，下面构造零假设C=Φ的MFT。这是一个渐近方法，在模型M中，让时间和变点在参数n中呈线性增长。即在以下渐近情况下，从参数T，c1，c2，…，ck和窗宽h下引入参数nT，nc1，…，nck和nh，n∈。这里，让X∈M，MFT是基于多个MOSUM统计量，它们比较了相邻窗口观测值的经验平均值。首先考虑一个窗宽h∈{1，2，…，T/2}的情形，时间范围t取[h，T-h]，其中·为下取整函数。故对于所有的t∈[h，T-h]，构建统计量：

(1)

(2)

(3)

注意，这里L是一个高斯过程，均值为零，方差为1，因此在零假设下D(n)通常在0附近波动。由命题1的收敛性可推导出两个较好的结果。首先，由于极大值运算符的连续性，确保了极大值的收敛性，即有：

(4)

MFT的进一步解释：式(5)～式(8)中再次描述了MFT。首先，在零假设下，基于观测值X，随机游走S(n)弱收敛为标准布朗运动W。给定窗宽h∈H，过程D(n)收敛到L。收敛关键的一点，在经验(左)侧，所有函数都是基于单个过程X的，而在极限(右)侧，所有函数都是基于布朗运动W求值，又由式(5)～(8)中所有映射都是连续的，从而保证了收敛性。因此，对于拒绝域值Q的确定，可以间接模拟布朗运动W得到，并求出每次模拟所有函数{(|Lh，t|)t∈[h，T-h]|h∈H}的最大值。最后将不同窗宽大小下的D(n)的最大值作为检验统计量，从而避免了原假设检验中的多重检验。

(5)

(6)

(7)

(8)

(a)单窗宽

在实践中，当n=1，但T很大，这种方法很适用，虽然是基于渐近结果，但只需要选择足够大的最小窗宽即可。例如，对于正态分布随机变量，对于h(取约30)值，可以认为逐点分布足够接近正态分布，也类似于t分布。

(二)多重滤波算法(MFA)

这里考虑存在均值变点的备择假设。在图3中，图3(a)和图3(b)展示了数据分布及均值分布，图3(c)和图3(d)展示了过程D(1)在三种不同的窗宽大小下的值分布。由于全局最大值M小于拒绝阈值Q(图3(c))，从而MFT接受零假设，故表明无均值变点。然而，如果零假设被拒绝(图3(d))，则目标是估计C。为此，本文采用了一种启发式算法，称为多过滤器算法(MFA)，最初提出该方法思想的是Messer等人对时间序列中的变点分析[17]。针对图3(b)序列中的变点，MFA检测方法具体如下。

(a)无均值变点序列

(b)有均值变点序列

(c)无均值变点序列变点检测

(d)有均值变点序列变点检测图3 MFA算法变点检测结果图

图3说明了MFT方法同时使用多个MOSUM过程的优势所在：较小的窗宽对快速变化的均值变点较灵敏，而较大的窗宽在均值变化较小时较灵敏。这里可看出，MFA并不代表统计经验，而是一种算法，可以在拒绝零假设后用于估计变点。同时，该算法利用多个窗宽的实质是变点估计的影响只是局部的，即给定的变点过程D(n)的影响只在h-邻域内，这是因为过程D(n)构成只依赖于2h，与实际数据的分布无关。

四、参数变化鲁棒性评估

X1，[1]，…，Xnc-1，[1]，Xnc，[2]，…，XnT，[2]

(9)

(10)

命题2 给定X∈M，对于C={nc}，c∈[h，T-h]和μ[1]=μ[2](无均值变化)。在(D[h，T-h]，dSK)上，当n→，有：

(a)L值和值分布

(b)滞后项分布图4 L值和值对比图

五、模拟与实证分析

(一)模拟

在本节中，为验证MFT方法对其它参数变化的鲁棒性以及对均值变点的灵敏性性能，模拟分为两部分，第一部分对均值变点检测的灵敏性和对均值外的参数变化的鲁棒性模拟探究，第二部分对其方法的优势性进行模拟分析。

首先，对模型均值外其它参数变化的鲁棒性分析，这里模拟设计存在不同时间尺度上的三个均值变点(c1，c3，c4)与一个方差变点(c2)的随机变量序列，见图5(a)。其中均值变点c1的变化很小，而均值变点(c3，c4)彼此距离较近且变化相对较大。由于分布假设具有一般性，选用正态分布和伽玛分布的随机变量作为测试数据，应用MFA算法进行统计检验变点，取显著性水平为5%。两分布序列1 000次模拟中检测到的变点数量直方图分布，见图5(b)和图5 (c)。

(a)序列均值和方差分布

(b)正态分布序列变点检测结果

(c)伽玛分布序列变点检测结果图5 变点检测频率直方图

表1 MFT方法变点检测结果

表2 Merge MOSUM(bandwidth)方法变点检测结果

表3 WBS方法变点检测结果

最后，通过模拟1实验分析，表明了MFT检测方法对均值外的其它参数具有较强的鲁棒性；针对模拟2，从表1与表2和表3的检测结果对比可以看出，在前两组数据检测结果中，MFT检测方法的正确率都高达90%以上，相比较Eichinger等人论文中的方法在准确度和精度方面的优势很明显，在第三组数据的检测结果中，MFT方法较Eichinger等人论文中的方法在准确度方面结果略差，相对精度方面MFT方法较优。总体而言，在方法的对比上MFT占据一定优势，也说明了方法对于均值变点检测具有良好的灵敏性。

(二)实证分析

为验证MFT方法的实用性，将使用该方法与WBS方法和“Merge MOSUM(bandwidth)”方法作实例对比分析。数据选取美国IBM(国际商业机器公司，International Business Machines Corporation)于1962年1月29日至2019年8月9日的周数据(数据来源：https：//finance.yahoo.com/)，其数据包含开盘价、最高价、最低价、收盘价、成交量和复权收盘价共6种数量。为分析IBM公司股票的变动情况，选取股票收盘价作为研究对象，共3 003条数据，数据分布见图6 (a)。

(a)wbs方法变点检测结果

(b)MFT方法变点检测结果

(c)Merge MOSUM方法变点检测结果图6 股票收盘价变点检测对比图

图6展示了IBM公司股票收盘价在三种方法下的变点检测结果，这里MFT方法和“Merge MOSUM”方法基于数据量，取相同窗宽H=250，280，300。由检测结果统计得到：wbs方法共检测到4个变点，见图6(a)。时间点对应于1982年12月27日、1998年7月27日、2004年4月16日、2010年10月11日；MFT方法共检测到5个变点，时间点对应于1973年11月19日、1983年1月3日、1991年4月8日、1998年7月27日、2010年10月11日，见图6(b)；“Merge MOSUM”方法共检测到6个变点，时间点对应于1973年11月19日、1983年1月3日、1991年4月1日、1998年7月20日、2002年11月18日、2010年10月18日，见图6 (c)。为便于分析，后续将以上三种方法分别记为法1、法2、法3。

根据图6结果，除了法3中检测到的第五个变点(2002年11月18日)外，法2(5个变点)和法3(6个变点)的检测结果极为相似，相似的5个变点在法2和法3中基本在同一时期被检测到，进一步看出这5个变点中存在3个变点也同时被法1检测到，由此分析IBM公司股票数据在这5个时期点可能存在波动。经查阅资料，变点1：可能由于20世纪60年代后期，美国政府巨额的财政赤字，通货膨胀、石油价格上涨、国际收支成为美国经济的重要问题，最终导致1971年战后的布雷顿森林体系的崩溃、美元贬值和1973年的第一次石油危机。变点2：可能由于当前世界经济危机从1979年开始，经历了1980年和1981年连续两次的下降之后，政府出台了相关的应对政策，直到1983年第一季度才开始出现经济回升的迹象，但是并不显著。变点3：可能由于美国经济在经历了1990年的危机和1991年苏联解体后，于1991年3月达到谷底，同年4月开始出现回升，此后持续了“自1854年有统计以来时间最长的一次”景气扩张。变点4：可能由于自1994年4月到1998年7月美国经济持续增长了84个月，失业率(4.5%)为20世纪70年代末的最低水平，而之后由于美国受亚洲金融危机和国内消费不振的影响(互联网泡沫事件)，导致经济增长有所放慢。变点5：可能由于2009年第二季度以来，美国金融市场状况有所改善，实体经济也出现积极变化，第二季度经济降幅明显收窄，第三季度出现了连续四个季度下滑的首次增长，美国经济开始步入复苏阶段。针对法3中出现的变点5(2002年11月18日)和法1中出现的变点3(2004年4月16日)，这两个变点有且只能被其中一种方法检测到，而另外两种方法均未检测到，由此推测两变点对应的时期点可能并没有发生变化。经查阅资料显示，在这两个时期美国经济保持相对稳定，并无较大的金融事件发生。

基于以上三种方法对比分析得出，在MFT方法检测的变点集中，对应时期点的经济状态都发生了变化，而在wbs方法和“Merge MOSUM”方法的检测结果中，出现了未识别变点和误检的情况。相比较三种方法而言，MFT方法更具实效性。

六、结束语

在时间序列中存在多个时间尺度的均值变点情况下，本文提出的MFT检测方法有三个显著的优势。首先，它允许模型具有一般分布假设，只需假设随机变量为i.i.d.的分段序列，也可对同分布或独立性条件进行放宽。其次，MFT使用了一个MOSUM型统计量和渐近情形，其中MOSUM过程弱收敛到一个标准布朗运动的函数，这种方法可同时应用多个MOSUM过程，改进了单一MOSUM过程对不同时间尺度上的变点检测的缺点。最后，MFT方法对模型外参数的变化具有较强的鲁棒性，如方差或高阶矩。

MFT方法建模简单且检测效果优良，对管理部门决策者制定相关战略及防范措施具有一定的实际意义。