机器学习（三）——时间序列法

笔者是一个痴迷于挖掘数据中的价值的学习人，希望在平日的工作学习中，挖掘数据的价值，找寻数据的秘密，笔者认为，数据的价值不仅仅只体现在企业中，个人也可以体会到数据的魅力，用技术力量探索行为密码，让大数据助跑每一个人，欢迎直筒们关注我的公众号，大家一起讨论数据中的那些有趣的事情。

我的公众号为：livandata

简而言之：
对某一个或者一组变量 x(t) 进行观察测量，将在一系列时刻t1,t2,⋯,tn 所得到的离散数字组成的序列集合，称之为时间序列。

例如: 某股票A从2015年6月1日到2016年6月1日之间各个交易日的收盘价，可以构成一个时间序列；某地每天的最高气温可以构成一个时间序列。

一些特征:
趋势：是时间序列在长时期内呈现出来的持续向上或持续向下的变动。

季节变动：是时间序列在一年内重复出现的周期性波动。它是诸如气候条件、生产条件、节假日或人们的风俗习惯等各种因素影响的结果。

循环波动：是时间序列呈现出得非固定长度的周期性变动。循环波动的周期可能会持续一段时间，但与趋势不同，它不是朝着单一方向的持续变动，而是涨落相同的交替波动。

不规则波动：是时间序列中除去趋势、季节变动和周期波动之后的随机波动。不规则波动通常总是夹杂在时间序列中，致使时间序列产生一种波浪形或震荡式的变动。只含有随机波动的序列也称为平稳序列。

严平稳：

如果对所有的时刻 t，任意正整数 k 和任意 k 个正整数：

的联合分布与：

的联合分布相同，我们称时间序列{rt} 是严平稳的。

也就是：

的联合分布在时间的平移变换下保持不变，这是个很强的条件。而我们经常假定的是平稳性的一个较弱的方式。

弱平稳：

若时间序列 {rt} 满足下面两个条件:

则时间序列 {rt} 为弱平稳的。即该序列的均值，rt与rt-l 的协方差不随时间而改变，l为任意整数。

在金融数据中，通常我们所说的平稳序列，是弱平稳的。

差分：

差分（这里为前向），就是求时间序列{rt}在 t时刻的值 rt 与 t-1 时刻的值r(t-1) 的差不妨记做 dt,则我们得到了一个新序列{dt}，为一阶差分，对新序列{dt}再做同样的操作，则为二阶差分。

通常非平稳序列可以经过d次差分，处理成弱平稳或者近似弱平稳时间序列。回头看上图，我们发现二阶差分得到的序列比一阶差分效果更好。

自相关函数 （AutocorrelationFunction, ACF）：

相关系数度量了两个向量的线性相关性，而在平稳时间序列 {rt} 中，我们有时候很想知道，rt 与它的过去值 rt-i 的线性相关性。 这时候我们把相关系数的概念推广到自相关系数。

rt 与 rt-l 的相关系数称为 rt 的间隔为 l 的自相关系数，通常记为 ρl。具体的：

这里用到了弱平稳序列的性质:

则函数： 称为rt 的样本自相关函数(ACF)。

当自相关函数中所有的值都为0时，我们认为该序列是完全不相关的；因此，我们经常需要检验多个自相关系数是否为0。

白噪声序列：

随机变量X（t）（t=1，2，3……），如果是由一个不相关的随机变量的序列构成的，即对于所有S不等于T，随机变量Xt和Xs的协方差为零，则称其为纯随机过程。

对于一个纯随机过程来说，若其期望和方差均为常数，则称之为白噪声过程。白噪声过程的样本实称成为白噪声序列，简称白噪声。之所以称为白噪声，是因为他和白光的特性类似，白光的光谱在各个频率上有相同的强度，白噪声的谱密度在各个频率上的值相同。

线性时间序列：

时间序列{rt}，如果能写成：

则我们称{rt} 为线性序列。其中at称为在t时刻的新息(innovation)或扰动(shock)。

很多时间序列具有线性性，即是线性时间序列，相应的有很多线性时间序列模型，例如接下来要介绍的AR、MA、ARMA，都是线性模型，但并不是所有的金融时间序列都是线性的。

对于弱平稳序列，我们利用白噪声的性质很容易得到rt的均值和方差：

因为 Var(rt)一定小于正无穷，因此必须是收敛序列，因此满足：

即：随着i的增大，远处的扰动 at-i 对 rt 的影响会逐渐消失。

一、自回归(AR)模型：

在第一章中，我们计算了上证指数部分数据段的ACF，看表可知间隔为1时：

根据这点我们可以建立下面的模型：

其中{at}是白噪声序列，这个模型与简单线性回归模型有相同的形式，这个模型也叫做一阶自回归(AR)模型，简称AR(1)模型

从AR(1)很容易推广到AR(p)模型：

AR(p)模型的特征根及平稳性检验。

我们先假定序列是弱平稳的，则有;

因为{at}是白噪声序列，因此有：

所以有：

根据平稳性的性质，又有E(rt) =E(rt-1)=...= u，从而：

对式2.1，假定分母不为0， 我们将下面的方程称为特征方程:

该方程所有解的倒数称为该模型的特征根，如果所有的特征根的模都小于1，则该AR(p)序列是平稳的。

下面我们就用该方法，检验上证指数日收益率序列的平稳性。

我们可以看看模型有多少阶。

可以看出，自动生成的AR模型是17阶的。关于价次的讨论在下节内容，我们画出模型的特征根，来检验平稳性。

可以看出，所有特征根都在单位圆内，则序列为平稳的！
二、AR(p)模型的定阶：

一般有两种方法来决定p：

第一种：利用偏相关函数(PartialAuto Correlation Function,PACF)

第二种：利用信息准则函数

偏相关函数判断p：

对于偏相关函数的介绍，这里不详细展开，重点介绍一个性质：

AR(p)序列的样本偏相关函数是 p 步截尾的。

所谓截尾，就是快速收敛应该是快速的降到几乎为0或者在置信区间以内。

具体我们看下面的例子，还是以之前上证指数日收益率序列。

我们看出，按照截尾来看，模型阶次p在110+，但是之前调用的自动生成AR模型，阶数为17,这一点我也有些不解。

当然，我们很少会用这么高的阶次。。

信息准则--- AIC、BIC、HQ ：

现在有以上这么多可供选择的模型，我们通常采用AIC法则。我们知道：增加自由参数的数目提高了拟合的优良性，AIC鼓励数据拟合的优良性但是尽量避免出现过度拟合(Overfitting)的情况。所以优先考虑的模型应是AIC值最小的那一个。赤池信息准则的方法是寻找可以最好地解释数据但包含最少自由参数的模型。不仅仅包括AIC准则，目前选择模型常用如下准则：

AIC = -2 ln(L) + 2 k 中文名字：赤池信息量 akaikeinformation criterion

BIC = -2 ln(L) + ln(n)*k中文名字：贝叶斯信息量 bayesian information criterion

HQ = -2 ln(L) +ln(ln(n))*k hannan-quinn criterion

下面我们来测试下3种准则下确定的p，仍然用上证指数日收益率序列。为了减少计算量，我们只计算间隔前10个看看效果。

可以看出，3个准则在第一点均取到最小值，也就是说，p的最佳取值应该在1，我们只计算了前10个，结果未必正确。

当然，利用上面的方法逐个计算是很耗时间的，实际上，有函数可以直接按照准则计算出合适的order，这个是针对 ARMA模型的，我们后续再讨论。

模型的检验：

根据式2.0，如果模型是充分的，其残差序列应该是白噪声，根据我们第一章里介绍的混成检验，可以用来检验残差与白噪声的接近程度。

我们先求出残差序列：

然后我们检查它是不是接近白噪声序列。

观察p-value可知，该序列可以认为没有相关性，近似得可以认为残差序列接近白噪声。

拟合优度及预测：

我们使用下面的统计量来衡量拟合优度：

但是，对于一个给定的数据集，R2是用参数个数的非降函数，为了克服该缺点，推荐使用调整后的R2:

它的值在0-1之间，越接近1，拟合效果越好。

下面我们计算之前对上证指数日收益率的AR模型的拟合优度。

可以看出，模型的拟合程度并不好，当然，这并不重要，也许是这个序列并不适合用AR模型拟合。

三、滑动平均(MA)模型：

MA(q)模型的形式：

c0为一个常数项。这里的at,是AR模型t时刻的扰动或者说新息，则可以发现，该模型，使用了过去q个时期的随机干扰或预测误差来线性表达当前的预测值。
MA模型的性质：

平稳性：

MA模型总是弱平稳的，因为他们是白噪声序列（残差序列）的有限线性组合。因此，根据弱平稳的性质可以得出两个结论：

自相关函数：

对q阶的MA模型，其自相关函数ACF总是q步截尾的。因此MA(q)序列只与其前q个延迟值线性相关，从而它是一个“有限记忆”的模型。

这一点可以用来确定模型的阶次，后面会介绍。

可逆性：

当满足可逆条件的时候，MA(q)模型可以改写为AR(p)模型。这里不进行推导，给出1阶和2阶MA的可逆性条件。

1阶：

2阶：

MA的阶次判定：

我们通常利用上面介绍的第二条性质：MA(q)模型的ACF函数q步截尾来判断模型阶次。示例如下：

使用上证指数的日涨跌数据（2013年1月至2014年8月）来进行分析，先取数据（点击查看原文代码）：

可以看出，序列看上去是弱平稳的。下面我们画出序列的ACF：

我们发现ACF函数在43处截尾，之后的acf函数均在置信区间内，我们判定该序列MA模型阶次为43阶。

建模和预测：

由于sm.tsa中没有单独的MA模块，我们利用ARMA模块，只要将其中AR的阶p设为0即可。

函数sm.tsa.ARMA中的输入参数中的order(p,q)，代表了AR和MA的阶次。模型阶次增高，计算量急剧增长，因此这里就建立10阶的模型作为示例，如果按上一节的判断阶次来建模，计算时间过长。

我们用最后10个数据作为out-sample的样本，用来对比预测值。

我们先来看看拟合效果,计算：

得结果0.0278706962641

可以看出，score远小于1，拟合效果不好。

然后我们用建立的模型进行预测最后10个数据

可以看出，建立的模型效果很差，预测值明显小了1到2个数量级！就算只看涨跌方向，正确率也不足50%。所以该模型不适用于原数据。

当然如果真的这么简单能预测指数日涨跌才奇怪~

关于MA的内容只做了简单介绍，下面主要介绍ARMA模型

四、ARMA模型：

在有些应用中，我们需要高阶的AR或MA模型才能充分地描述数据的动态结构，这样问题会变得很繁琐。为了克服这个困难，提出了自回归滑动平均(ARMA)模型。

基本思想是把AR和MA模型结合在一起，使所使用的参数个数保持很小。

模型的形式为：

其中，{at}为白噪声序列，p和q都是非负整数。AR和MA模型都是ARMA(p,q)的特殊形式。

利用向后推移算子B，上述模型可写成:

（后移算子B，即上一时刻）

这时候我们求rt的期望，得到：

和上期我们的AR模型一毛一样。因此有着相同的特征方程：

该方程所有解的倒数称为该模型的特征根，如果所有的特征根的模都小于1，则该ARMA模型是平稳的。

有一点很关键：ARMA模型的应用对象应该为平稳序列！ 我们下面的步骤都是建立在假设原序列平稳的条件下的~

4.1 识别ARMA模型阶次

4.1.1 PACF、ACF 判断模型阶次

我们通过观察PACF和ACF截尾，分别判断p、q的值。(限定滞后阶数50)

可以看出，模型的阶次应该为(27,27)。然而，这么高的阶次建模的计算量是巨大的。

为什么不再限制滞后阶数小一些？如果lags设置为25，20或者更小时，阶数为(0,0)，显然不是我们想要的结果。

综合来看，由于计算量太大，在这里就不使用(27,27)建模了。采用另外一种方法确定阶数。

4.2 信息准则定阶

关于信息准则，上一期有过一些介绍：

目前选择模型常用如下准则： （其中L为似然函数，k为参数数量，n为观察数）

AIC = -2 ln(L) + 2 k 中文名字：赤池信息量 akaikeinformation criterion

BIC = -2 ln(L) +ln(n)*k 中文名字：贝叶斯信息量 bayesian information criterion

HQ = -2 ln(L) +ln(ln(n))*k hannan-quinn criterion

我们常用的是AIC准则，AIC鼓励数据拟合的优良性但是尽量避免出现过度拟合(Overfitting)的情况。所以优先考虑的模型应是AIC值最小的那一个模型。

下面，我们分别应用以上3种法则，为我们的模型定阶，数据仍然是上证指数日涨跌幅序列：

为了控制计算量，我们限制AR最大阶不超过6，MA最大阶不超过4。 但是这样带来的坏处是可能为局部最优。

可以看出，AIC求解的模型阶次为(3,3)。我们这里就以AIC准则为准~至于到底哪种准则更好，小伙伴们可以分别建模进行对比~

4.2 模型的建立及预测

我们使用上一节AIC准则求解的模型阶次(3,3)来建立ARMA模型，源数据为上证指数日涨跌幅数据，最后10个数据用于预测。

同样的，先来看看拟合效果：

得结果：0.0495081497069

对比之前建立的AR、MA模型，可以发现拟合精度上有所提升，但仍然是不够看的级别。

接着来看预测效果：

 

如图，可以看出，虽然还是准头很低，不过相比之前的MA模型，只看涨跌的话，胜率为55.6%，效果还是好了不少。

五、 ARIMA模型：

到目前为止，我们研究的序列都集中在平稳序列。，即ARMA模型研究的对象为平稳序列。如果序列是非平稳的，就可以考虑使用ARIMA模型。

ARIMA比ARMA仅多了个"I"，代表着其比ARMA多一层内涵：也就是差分

一个非平稳序列经过d次差分后，可以转化为平稳时间序列。d 具体的取值，我们得分被对差分1次后的序列进行平稳性检验，若果是非平稳的，则继续差分。直到d次后检验为平稳序列。

5.1 单位根检验：

ADF是一种常用的单位根检验方法，它的原假设为序列具有单位根，即非平稳，对于一个平稳的时序数据，就需要在给定的置信水平上显著，拒绝原假设。

下面给出示例，我们先看上证综指的日指数序列：

看图形，这里显然是非平稳的。接着我们进行ADF单位根检验。

可以看出，p-value为0.1704489，大于显著性水平。原假设：序列具有单位根即非平稳。不能被拒绝。 因此上证指数日指数序列为非平稳的。

我们将序列进行1次差分后再次检验！

从图形来看，序列近似平稳序列，我们来进行ADF检验：p-value: 2.31245750144e-30

可以看出，p-value非常接近于0，拒绝原假设，因此，该序列为平稳的。

可见，经过1次差分后的序列是平稳的，对于原序列，d的取值为1即可。

5.2 ARIMA(p,d,q)模型阶次确定

上面一小节我们确定了差分次数d，接下来，我们就可以将差分后的序列建立ARMA模型。

首先，我们还是尝试PACF和ACF来判断p、q。

可以看出，模型的阶次为(27,27)，还是太高了。建模计算了太大。我们再看看AIC准则根据AIC准则，差分后的序列的ARMA模型阶次为(2,2)。因此，我们要建立的ARIMA模型阶次(p,d,q) = (2,1,2)

5.3 ARIMA模型建立及预测

根据上一结确定的模型阶次，我们对差分后序列建立ARMA(2,2)模型：

我们先看差分序列的ARMA拟合值。

得0.0397489997089

再看对差分序列的预测情况：

可以看出，差分序列ARMA模型的拟合效果和预测结果并不好，预测值非常小，这代表什么？这代表对于新的值，这里认为它很接近上一时刻的值。

这个影响可能来自模型阶次，看来模型阶次还是得尝试更高阶的。这里就不建模了，（计算时间太长），大家有兴趣可以试试高阶的模型。

然后我们将预测值还原（即在上一时刻指数值的基础上加上差分差值的预估）：

我们发现，预测序列像是比原序列延后了一段时间。 这充分表现了对最近数据的依赖。这是因为对差分序列进行预测时过于“保守”。显然，这种滞后性，在我们预测金融数据时会陷入被动。

当然，这也说明该模型并不能很好得预测指数，上面也说过，模型阶次的提高可能会提高预测的效果。