使用Python库valuequant处理时间序列结构变化的建模方法2

继续上期关于外生刺激复杂多变条件下时间序列数据结构变化处理的讨论。本期介绍第2种需求场景及其建模方法。

场景描述

（1）时间序列的非周期波动比较大

（2）需要识别时间序列当前趋势状态，即识别时间序列当前n阶差分正负偏向

（3）时间序列的趋势变化容易被剧烈的非周期波动掩盖

（4）对时间序列当前趋势状态的判定需要足够的时间步数据进行验证

（5）对时间序列的未来估计只考虑符合当前趋势的最新信息

（6）环境多变，数据稀疏量少，暴力型大数据算法不适用

建模方法

第一步，针对第(2)点需求，计算时间序列的n阶差分序列。因为实际应用过程中，时间序列趋势变化一般以线性变化为主，除随机变化情况外，差分序列偏向大多出现在1阶和2阶，极少出现3阶以上的差分偏向。所以，本期的建模方法根据时间序列的1阶和2阶差分序列的正负偏向分析时间序列当前趋势状态，即计算分析时间序列的1阶差分和2阶差分。

第二步，针对(1)(3)(4)点需求，本期模型相较上期提出的模型，不再按趋势明确与否进行划分，而是按分析者设置的粒子时间步MS将时间序列划分为多个粒子波段。

第三步，从1阶差分最新粒子波段开始向前回溯，逐个粒子波段按显著性水平分别进行偏向正和偏向负的单尾假设检验，判定1阶粒子波段趋势状态1阶为正（+）、负（-）或无（n）：

（a）若粒子波段至最新粒子波段的趋势状态均1阶为正（+），

        且或粒子波段的趋势状态1阶非正（+），

        则判定当前趋势状态1阶为正（+），且为当前趋势状态1阶序列；

（b）若粒子波段至最新粒子波段的趋势状态均1阶为负（-），

        且或粒子波段的趋势状态1阶非负（-），

        则判定当前趋势状态1阶为负（-），且为当前趋势状态1阶序列； 

（c） 若粒子波段至最新粒子波段的趋势状态均1阶为无（n），

        且或粒子波段的趋势状态1阶非无（n），

        则判定当前趋势状态1阶为无（n）， 且为当前趋势状态1阶序列；

（d）若，则判定当前趋势状态1阶为无（n），且为当前趋势状态1阶序列。

第四步：若当前趋势状态1阶状态为正或负，则对2阶差分进行相同的分析步骤，识别当前趋势状态2阶状态和当前趋势状态2阶序列。但不再进行3阶以后差分分析。

至第四步结束，已经解决了第(1)(2)(3)(4)点需求，并且完成时间序列当前趋势结构-即当前趋势状态和当前趋势序列的识别。建模者可结合识别结果和习惯使用的参数估计方法进行后续工作，也可按接下来的步骤进行参数估计。

第五步：满足第(5)(6)点需求，建立以下时间序列推演模型（预测模型）：

（a）当前趋势状态2阶为正（+）或负（-），计算当前趋势状态2阶序列的均值为2阶变化值G，取1阶差分最新值为1阶变化值D，取时间序列最新值为基期值A。将A、D、G作为时间序列推演模型参数，基期，则模型的未来推演序列（预测序列）的计算公式为：

（b）当前趋势状态1阶为正（+）或负（-），计算当前趋势状态1阶序列的均值为1阶变化值D，取时间序列最新值为基期值A。将A、D作为时间序列推演模型参数，基期，则模型的未来推演序列（预测序列）的计算公式为：

（c）当前趋势状态1阶为无（n），取时间序列的均值为基期值A。将A作为时间序列推演模型参数，则模型的未来推演序列（预测序列）的计算公式为：

为了方便记忆及与后续介绍的建模方法进行区分，称该方法为概率状态更新模型Probable Latest Status Model，可以直接使用Python库valuequant实现，分析者可在终端执行命令pip install valuequant安装该模块。详细使用方法请查阅BoomEvolve官网valuequant官网。

代码实现和实际应用的简单示例

构建上期曾讨论的波动幅度较大、容易被上期的严格状态更新模型忽略趋势结构变化的时间序列作为测试序列。

>>> #构建测试序列 
>>> import pandas as pd 
>>> import numpy as np 
>>> import matplotlib.pyplot as plt 
>>> series=pd.Series(20-np.arange(1,13)+np.random.normal(0,1,12),index=pd.date_range(start='2010-01-01',periods=12,freq='Q-DEC')) 
>>> aseries=pd.Series(series[-1]+np.random.normal(0,1,12),index=pd.date_range(start='2013-01-01',periods=12,freq='Q-DEC')) 
>>> series=pd.concat((series,aseries)) 
>>> aseries=pd.Series(np.arange(1,21)*1.2+series[-1]+np.random.normal(0,1,20),index=pd.date_range(start='2016-01-01',periods=20,freq='Q-DEC')) 
>>> series=pd.concat((series,aseries)) 
>>> #绘图观察测试序列 
>>> plt.plot(series) 
>>> plt.show()

​接下来使用valuequant的概率状态更新模型实现当前趋势结构-即当前趋势状态和当前趋势序列的识别以及参数估计。

>>> #调用valuequant并登录账号，登录账号可在https://www.boomevolve.com注册获得 
>>> import valuequant as vq 
>>> vq.login(name=<账户名称>,pswd=<账户密码>) 
登录成功 
>>> #构建概率状态更新模型，minstep参数是分解时间序列粒子波段的粒子时间步，significance参数是假设检验的显著性水平 
>>> model=vq.Models.probable_latest_status_model(series=series,minstep=8,significance=0.3) 
>>> #获取建模结果 
>>> model.model() {'method': 'l1', 'step': 18, 'params': {'l1_diff': 1.167286, 'value_latest': 31.54779}, 'status': {'l1': '+', 'l2': 'n'}} 
>>> #结果显示当前趋势状态status为1阶正2阶无，当前趋势状态序列为最后18步，使用1阶变化值进行未来序列推演 
>>> #可使用model.fit()命令获取样本内拟合数据，使用model.forecast(t=<需要预测的步数>)获取样本外推演数据 
>>> #绘图展示结果 
>>> plt.plot(model.data()) 
>>> plt.plot(model.fit()) 
>>> plt.plot(model.forecast(t=12)) 
>>> plt.show()

​以上是概率状态更新模型Probable Latest Status Model的简单代码实现。

注意事项：valuequant在建模过程中，会自动进行异常样本过滤。分析者也可根据建模需要，在使用前自行多加异常处理流程。

含季节项时间序列的应用示例

与上期严格状态更新模型同样的，概率状态更新模型的主要功能是识别时间序列的趋势结构变化，对含有季节项的时间序列，需要先分解时间序列季节项。该需求可以通过设置valuequant相关模型的seasonal参数满足。在上述测试时间序列的基础上，增加季节乘数，构建含季节项的测试序列，以作演示阐述。

>>> #构建含季节项的测试序列 
>>> series=series*pd.Series(np.tile([1.05,1.1,0.95,0.9],int(np.ceil(len(series)/4)))[:len(series)], index=series.index) 
>>> #绘图演示测试序列 
>>> plt.plot(series) 
>>> plt.show() 
>>> 设置seasonal=True（默认为False），分离时间序列季节项 
>>> model=vq.Models.probable_latest_status_model(series=series,minstep=8,significance=0.3,seasonal=True) 
>>> model.model() 
{'method': 'l1', 'step': 18, 'params': {'l1_diff': 1.150202, 'value_latest': 31.484396}, 'status': {'l1': '+', 'l2': 'n'}, 'seasonal': {'12': 0.901812, '3': 1.028093, '6': 1.117934, '9': 0.931923}, 'qorder': ['3', '6', '9', '12']} 
>>> #模型结果增加了关于季节乘数'seasonal'的估计 
>>> #绘图展示建模结果 
>>> plt.plot(model.data()) 
>>> plt.plot(model.fit()) 
>>> plt.plot(model.forecast(t=12)) 
>>> plt.show()

​

参数突变问题

需要注意，对于存在参数突变的情况，该方法并未提供直接识别的功能，但转折附近的2阶变化往往会使得模型将转折处附近的波段和后续波段分开（转折处附近的信息也会因此被忽略）。

>>> #构建参数突变但趋势方向不变的测试序列 
>>> series=pd.Series(np.arange(1,13)+1+np.random.normal(0,0.5,12),index=pd.date_range(start='2015-01-01',periods=12,freq='Q-DEC')) 
>>> aseries=pd.Series(np.arange(1,13)*3+series[-1]+np.random.normal(0,0.5,12),index=pd.date_range(start='2018-01-01',periods=12,freq='Q-DEC')) 
>>> series=pd.concat((series,aseries)) 
>>> # 设置anneal参数，在假设检验和均值计算中使用时间权重 
>>> model=vq.Models.probable_latest_status_model(series=series,minstep=8,significance=0.3,anneal=0.9) 
>>> model.model() 
{'method': 'l1', 'step': 10, 'params': {'l1_diff': 3.00639, 'value_latest': 49.550087}, 'status': {'l1': '+', 'l2': 'n'}}

​使用保守估计参数

概率更新状态模型与严格更新状态模型相似，由于严格状态更新模型基于1、2阶差分变化值估计，在1阶差分方向与2阶差分方向不一致时，模型的未来推演序列（预测序列）将在有限步数内出现1阶差分方向反转。如果分析者对该反转持保守态度，可以设置valuequant相关模型的参数conservative=True（默认为False），使得模型不再使用2阶变化值，而使用1阶差分变化率g，模型的未来推演序列（预测序列）的计算公式改为

具体代码实现如下：

>>> #构建存在2阶变化，且1、2阶差分方向不一致的序列 
>>> series=pd.Series(8-6*pow(0.98,np.arange(1,31))+np.random.normal(0,0.001,30),index=pd.date_range(start='2013-01-01',periods=30,freq='Q-DEC')) 
>>> #设置严格状态更新模型参数conservative=True 
>>> model=vq.Models.probable_latest_status_model(series=series,minstep=8,significance=0.3,conservative=True) 
>>> model.model() 
{'method': 'l2smooth', 'step': 18, 'params': {'l1_latest': 0.06604399999999977, 'value_latest': 4.725968, 'l2smooth_rate': 0.9783735482464369}, 'status': {'l1': '+', 'l2': '-'}, 'smooth': True} 
>>> #绘图展示结果 
>>> plt.plot(model.data()) 
>>> plt.plot(model.fit()) 
>>> plt.plot(model.forecast(t=12)) 
>>> plt.show()

​

如上左图设置conservative=True相较右图未设置conservative（默认False），使用保守估计方法可以使得未来推演序列的1阶差分绝对值趋于0而非发生方向反转。

概率状态更新模型的局限性

局限1：在处理同时存在乘数季节项和高阶差分变化的时间序列时，分离季节项乘数难免产生的误差容易导致异方差问题，进而影响分离季节项后的趋势状态识别。

局限2：容易忽视长期趋势方向。

局限3：如前文所述，容易在参数突变转折处忽略信息。

其中，概率状态更新模型的第1、2缺陷与严格状态更新模型相似，可参考上期相关内容描述，此篇不再赘述。

如分析者需要解决上述问题或需要应对与该模型不同的需求场景，请阅读本专栏后续文章其他相关建模方法的介绍。