时间序列建模的向量自动回归(VAR)实战教程

1,249 阅读9分钟

向量自回归(VAR)是一种用于多变量时间序列分析的统计模型,特别是在变量之间有相互影响的关系的时间序列。VAR模型与单变量自回归模型不同,因为它允许对多变量时间序列数据进行分析并作出预测。VAR模型经常被用于经济学和天气预报中。

使用VAR模型的基本要求是 。

  • 至少有两个变量的时间序列。
  • 变量之间的关系。

它被认为是一个自回归模型,因为该模型所做的预测取决于过去的数值,这意味着每个观测值被建模为其滞后值的函数。

ARIMA系列与VAR模型的基本区别在于,所有ARIMA模型都用于单变量时间序列,而VAR模型则用于多变量时间序列。此外,ARIMA模型是单向模型,这意味着因变量受其过去或滞后值本身的影响,而VAR是双向模型,这意味着因变量受其过去的值或另一变量的值的影响,或受两者的影响。

关于时间序列的更多理解,请参考这些文章: - 时间序列数据分析的一般概述。

要了解更多关于时间序列建模的信息,请参考这些文章:------。

这篇文章将是另一个时间序列建模的指南,但这次将使用多变量时间序列数据。同时,我们将通过一些测试,这是了解多变量时间序列所需要的。在开始建模部分之前,让我们先了解一下模型背后的数学原理。

什么是向量自回归(VAR)?

一个典型的单变量时间序列的自回归模型(AR(p))可以表示为

图片来源

其中

  • y_t_-i表示早期的变量值。
  • A_i_是一个时间不变的_(k_×_k_)矩阵。
  • _e__t_是一个误差项。
  • c是模型的截距。

这里的p级是指,最多使用y的p个滞后期。

我们知道,VAR模型处理的是多变量的时间序列,这意味着会有两个或更多的变量相互影响。因此,VAR模型的方程随着时间序列中变量数量的增加而增加。

假设有两个时间序列变量,y1和y2,那么为了计算y1(t),VAR模型将使用两个时间序列变量的滞后期。

例如,有两个时间序列变量(y1和y2)的VAR(1)模型的方程式将是这样的。

图片来源

其中,Y{1,t-1}是yy1的第一个滞后值,Y{2,t-1}是y2的第一个滞后值。

而带有y1和y2时间序列变量的VAR(2),模型的方程会是这样的。

图片来源

在这里,我们可以清楚地了解模型的方程将如何随着变量和滞后值的增加而增加。例如,一个有3个时间序列变量的VAR(3)模型方程将看起来像。

图片来源

因此,这就是P值将如何增加模型方程的长度,而变量的数量将增加方程的高度。

在Python中实现矢量自回归(VAR)

让我们用python建立一个基本的VAR模型。

为了建立这个模型,我们可以使用python的statsmodel包,它提供了大部分的模块来进行时间序列分析,并提供了一些数据与包来进行时间序列分析的练习。

导入库

import numpy as np
import pandas as pd
import statsmodels.api as sm
from statsmodels.tsa.api import VAR
from statsmodels.tsa.base import  datetools

导入数据

mdata = sm.datasets.macrodata.load_pandas().data
mdata.head()

输出。

在这里我们可以看到我们的数据是怎样的。在这里我考虑了三个变量real gdp real cons和realinv,以便进一步建模处理。在进行进一步处理之前,还需要使用年份和季度列来制作一个日期时间值。

mdata = mdata[['year','quarter']].astype(int)
mdata = mdata[['year','quarter']].astype(str)
quarterly = mdata["year"] + "Q" + mdata["quarter"]
quarterly = datetools.dates_from_str(quarterly)
quarterly

输出。

现在我们可以使用日期时间值作为我们数据的索引。

mata = mdata[['realgdp','realcons','realinv']]
mdata.index = pandas.DatetimeIndex(quarterly)
mdata.head()

输出。

视觉化的数据。

mdata.plot()

输出。

在这里我们可以看到,realgdp和realcons都有很高的相关性,而realinv和其他变量之间也有轻微的相关性。因为我们看到的趋势,我们可以理解,所有的人都在稳步增长,直到最后一个点,但他们都有某种程度的衰退。

所以我们可以在它们的对数值中检查时间序列的模式。

data = np.log(mdata).diff().dropna()
data.plot()

输出。

在这里我们可以看到归一化的数值是如何在时间序列中被绘制出来的。realgdp和realcons之间的关系很强,因为它们遵循相同的线型,这清楚地表明,realgdp的变化将影响realcons,反之亦然。为了进行更多的分析,我们可以对时间序列进行一些测试,从统计学上确定变量之间的关系。

格兰杰因果关系测试

通过使用格兰杰因果关系检验,我们可以在建立模型之前找到变量之间的关系,因为众所周知,如果变量之间没有关系,我们可以放弃这些变量,单独进行建模。如果它们之间存在关系,我们需要在建模部分考虑该变量。

在数学中,测试提供了变量之间的p值,如果p值大于0.05,那么我们就需要接受无效假设,如果p值小于0.05,我们就需要拒绝无效假设。

Statsmodel也提供了一个模块来进行测试,所以接下来使用statsmodel,我将进行格兰杰因果关系测试。

from statsmodels.tsa.stattools import grangercausalitytests

data = mdata[["realgdp", "realcons"]].pct_change().dropna()
#Performing test on for realgdp and realcons.
gc_res = grangercausalitytests(data, 12)

输出。

这里我们可以看到,每个滞后期的P值都是零。所以现在,让我们继续进行realgdp和real inv之间的因果关系检验。

data = mdata[["realgdp", "realinv"]].pct_change().dropna()

输出。

这里我们可以看到每个滞后期的P值都高于0.05,这意味着我们需要接受无效假设。另外,我们可以说,如果我们在realcons和realinv之间进行检验,也会发生类似的事情。

data = mdata[["realcons", "realinv"]].pct_change().dropna()

输出。

在这里,我们得到的P值比以前高。这意味着realcons正在影响realinv。在图表中,我们估计realinv和其他数值之间没有明显的关系,因此我们可以理解格兰杰因果关系检验的意义。

协整检验

协整有助于找出两个或多个时间序列之间的统计联系。当两个或更多的时间序列是协整的,它们就有一个长期的、有统计学意义的关系。

我们可以使用statsmodel软件包进行这一测试。

data = mdata[["realgdp","realcons", "realinv"]].pct_change().dropna()
from statsmodels.tsa.vector_ar.vecm import coint_johansen

def cointegration_test(data, alpha=0.05): 

   """Perform Johanson's Cointegration Test and Report Summary"""
  out = coint_johansen(data,-1,5)
  d = {'0.90':0, '0.95':1, '0.99':2}
   traces = out.lr1
   cvts = out.cvt[:, d[str(1-alpha)]]
def adjust(val, length= 6): return str(val).ljust(length)

   # Summary
    print('Name   ::  Test Stat > C(95%)    =>   Signif  \n', '--'*20)
    for col, trace, cvt in zip(data.columns, traces, cvts):
        print(adjust(col), ':: ', adjust(round(trace,2), 9), ">", adjust(cvt, 8), ' =>  ' , trace > cvt)
cointegration_test(data)

输出。

在这里,我们可以看到变量对整个系统的意义。

所以,在这里我们已经看到了我们如何使用这两个测试;接下来,我们可以进一步进行建模部分。

建模

我们可以直接将预处理过的数据放入VAR模块进行建模。

var = VAR(data)

在这一步之后,程序中出现了一件事:如何选择顺序。通常要做的一件事是检查最佳适合的滞后值。我们需要比较不同的AIC(Akaike信息准则)、BIC(Bayesian信息准则)、FPE(Focused Prediction Error)和HQIC(Hannan-Quinn信息准则)。这些都是有助于选择最适合的滞后值的参数。

Statsmodel提供select_order模块来分析这些值。

x= var.select_order()
x.summary()

输出。

这里我们可以看到与AIC、BIC、FPE和HQIC相结合的最小值是以*号给出的。这里我们可以看到第三行和第一行都有这个符号。这里我选择了第三行,这意味着滞后值VAR(p)的值是3,因为它被建议使用AIC和其他参数产生最小值的组合。

results = var.fit(3)
#We can check the summary of the model by.
results.summary()

输出。

这里我们可以看到所有的系数标准误差值t检验和模型在每个滞后期的概率,直到3滞后期。然后,在最后,有一个混淆矩阵,显示了变量之间的相关性。在结果中,我们发现realgdp和realinv之间的相关性很高,在概率上的相关效应也是我们可以交叉检查的。我们还可以绘制模型,这将是了解模型性能的一个更好方法。

可视化的输入。

results.plot();

输出。

我们也可以通过模型绘制我们的预测值。

results.plot_forecast(20);

输出。

在这里,我们可以看到每个变量的模型的结果。未来20步的预测值的线条是以这样一种稳定的方式进行的,所以在这里我们可以看到我们决定的滞后期是相当令人满意的,并提供了良好的结果。

在这里,我们已经看到,在数据是多变量的情况下,我们如何进行时间序列建模。我们已经看到在这种情况下,我们如何理解一个数据中呈现的两个时间序列变量之间的关系。这种情况可以有很多例子。例如,大气温度的变化可以由湿度和季节因素引起,也可以由股票市场数据引起。我鼓励你在现实世界的场景数据集中使用这个模型,这样我们就能更深入地了解多变量时间序列分析。

参考文献:

The postHands-On Tutorial on Vector AutoRegression(VAR) For Time Series Modelingappeared first onAnalytics India Magazine.