时间序列建模的步骤（多元时间序列数据的预测和建模）

简介

时间是决定企业是上升还是下降的最关键因素。 这就是为什么我们看到商店和电子商务平台的销售与节日保持一致。 这些企业分析了多年的支出数据，以了解打开大门并看到消费者支出增加的最佳时机。

但作为数据科学家，您如何进行此分析？ 别担心，你不需要建立时间机器！ 时间序列建模是一种强大的技术，可作为理解和预测趋势和模式的门户。

但即使是时间序列模型也有不同的方面。 我们在网上看到的大多数例子都涉及单变量时间序列。 不幸的是，现实世界的用例并不像那样。 有多个变量在起作用，同时处理所有这些变量是数据科学家将获得其价值的地方。

在本文中，我们将了解多元变量时间序列是什么，以及如何处理它。 我们还将进行案例研究并在Python中实现它，以便您对该主题有实际的了解。

一、单变量时间序列与多变量时间序列

本文假设您熟悉单变量时间序列，其属性以及用于预测的各种技术。 由于本文将重点关注多变量时间序列，我建议您阅读以下文章，作为单变量时间序列的一个很好的介绍：

• 创建时间序列预测的综合指南
• 使用Auto Arima构建高性能时间序列模型

但在进入多变量时间序列的详细信息之前，我会给你快速回顾一下单变量时间序列。 让我们逐一看一下它们来理解它们之间的区别。

1.1、单变量时间序列

顾名思义，单变量时间序列是一个具有单个时间相关变量的系列。

例如，查看下面的样本数据集，其中包含过去2年的温度值（每小时）。 这里，温度是因变量（取决于时间）。

如果要求我们预测接下来几天的温度，我们将查看过去的值并尝试测量和提取一种模式。 我们会注意到早上和晚上的温度较低，而下午达到峰值。 此外，如果您有过去几年的数据，您会发现它在11月到1月期间更冷，而在4月到6月则相对较热。

这些观察将有助于我们预测未来的价值观。 您是否注意到我们只使用了一个变量（过去2年的温度）？ 因此，这称为单变量时间序列分析/预测。

1.2、多变量时间序列（MTS）

多变量时间序列具有多个时间相关变量。 每个变量不仅取决于其过去的值，而且还依赖于其他变量。 此依赖关系用于预测未来值。 听起来很复杂？ 让我解释。

考虑上面的例子。 现在假设我们的数据集包括出汗百分比、露点(dew point，在一定压力下出 现第一个液滴时的温度)、风速、云量百分比等，以及过去两年的温度值。 在这种情况下，需要考虑多个变量来最佳地预测温度。 像这样的系列将属于多变量时间序列的范畴。 以下是对此的说明：

现在我们了解了多变量时间序列的样子，让我们了解如何使用它来构建预测。

二、处理多变量时间序列 - VAR

在本节中，我将向您介绍多变量时间序列预测中最常用的方法之一 ——即向量自回归模型（Vector Auto Regression，VAR）。

在VAR模型中，每个变量都是过去自身值和所有其他变量的过去值的线性函数。 为了更好地解释这一点，我将使用一个简单的可视化示例：

我们有两个变量，y1和y2。 我们需要在时间t预测这两个变量的值，从过去n值的给定数据开始。 为简单起见，我认为滞后值为1。

为了计算y1（t），我们将使用y1和y2的过去值。 类似地，为了计算y2（t），将使用y1和y2的过去值。 下面是一种表示这种关系的简单数学方法：

其中，a​1​​和a​2​​是常数项，w​∗​​是系数矩阵中的系数，e​∗​​是误差项。

这些方程类似于AR过程的方程。 由于AR过程用于单变量时间序列数据，因此未来值仅是它们自己的过去值的线性组合。 考虑AR（1）过程：

y(t)=a w×y(t−1) e

在这个例子中，我们只有一个变量y，一个常数项a，一个误差项e和一个系数w。为了适应VAR的每个等式中的多个变量项，我们将使用向量。 我们可以用以下形式写出等式（1）和（2）：

上述形式扩展到多维之后：

上述等式中的ϵt表示多变量白噪声。 对于多变量时间序列，ϵt应该是满足以下条件的连续随机向量：

• E(εt) = 0，即误差期望为0
• E(εt1,εt2‘) = σ12，εt和εt‘的期望是时间序列的标准差

三、为什么需要VAR？

回想一下我们之前看到的温度预测的例子。 可以对其进行论证，将其视为多个单变量系列。 我们可以使用像AR这样的简单单变量预测方法来解决它。 由于目标是预测温度，我们可以简单地删除其他变量（温度除外）并在剩余的单变量系列上拟合模型。

另一个简单的想法是使用我们已知的技术单独预测每个系列的值。 这将使工作变得非常简单！ 那你为什么要学习另一种预测技术呢？ 这个话题不够复杂吗？

从上面的等式（1）和（2）可以清楚地看出，每个变量都使用每个变量的过去值来进行预测。 与AR不同，VAR能够理解和使用几个变量之间的关系。 这对于描述数据的动态行为非常有用，并且还可以提供更好的预测结果。 此外，实现VAR就像使用任何其他单变量技术一样简单（您将在上一节中看到）。

四、多元时间序列的平稳性

我们从研究单变量概念得知，固定时间序列通常会给我们一组更好的预测。 如果您不熟悉平稳性的概念，请首先阅读本文：处理非平稳时间序列的温度介绍

从数学上讲，如果时间序列的数据的均值和方差没有系统变化，且严格消除了周期变化，那么这个时间序列就是平稳的，只有平稳的时间序列数据的预测才是有可能的。非平稳时间序列数据目前在理论上有缺陷，预测很难。

五、Python的实现

预测大气质量，数据来源https://archive.ics.uci.edu/ml/datasets/Air Quality#：

#import required packages import pandas as pd import matplotlib.pyplot as plt %matplotlib inline #read the data df = pd.read_csv("AirQualityUCI.csv", parse_dates=[['Date', 'Time']]) #check the dtypes df.dtypes

输出结果：

ate_Time object CO(GT) int64 PT08.S1(CO) int64 NMHC(GT) int64 C6H6(GT) int64 PT08.S2(NMHC) int64 NOx(GT) int64 PT08.S3(NOx) int64 NO2(GT) int64 PT08.S4(NO2) int64 PT08.S5(O3) int64 T int64 RH int64 AH int64 dtype: object

Date_Time列的数据类型是object，我们需要将其更改为datetime。 此外，为了准备数据，我们需要索引具有日期时间。 请遵循以下命令：

df['Date_Time'] = pd.to_datetime(df.Date_Time , format = '%d/%m/%Y %H.%M.%S') data = df.drop(['Date_Time'], axis=1) data.index = df.Date_Time

下一步是处理缺失的值。 由于数据中的缺失值被替换为值-200，因此我们必须使用更好的数字来估算缺失值。 考虑一下 - 如果缺少当前的露点值，我们可以安全地假设它将接近前一个小时的值。 有道理，对吗？ 在这里，我将用前一个值来估算-200。

您可以选择使用之前几个值的平均值替换该值，或者选择前一天同时的值（您可以在下面的评论部分中分享您对缺失值进行估算的想法）。

#missing value treatment cols = data.columns for j in cols: for i in range(0,len(data)): if data[j][i] == -200: data[j][i] = data[j][i-1] #checking stationarity from statsmodels.tsa.vector_ar.vecm import coint_johansen #since the test works for only 12 variables, I have randomly dropped #in the next iteration, I would drop another and check the eigenvalues johan_test_temp = data.drop([ 'CO(GT)'], axis=1) coint_johansen(johan_test_temp,-1,1).eig

结果：

array([ 0.17806667, 0.1552133 , 0.1274826 , 0.12277888, 0.09554265, 0.08383711, 0.07246919, 0.06337852, 0.04051374, 0.02652395, 0.01467492, 0.00051835])

我们现在可以继续创建验证集以适应模型，并测试模型的性能：

#creating the train and validation set train = data[:int(0.8*(len(data)))] valid = data[int(0.8*(len(data))):] #fit the model from statsmodels.tsa.vector_ar.var_model import VAR model = VAR(endog=train) model_fit = model.fit() # make prediction on validation prediction = model_fit.forecast(model_fit.y, steps=len(valid))

预测采用数组的形式，其中每个列表代表行的预测。 我们将把它转换成更易于呈现的格式。

#converting predictions to dataframe pred = pd.DataFrame(index=range(0,len(prediction)),columns=[cols]) for j in range(0,13): for i in range(0, len(prediction)): pred.iloc[i][j] = prediction[i][j] #check rmse for i in cols: print('rmse value for', i, 'is : ', sqrt(mean_squared_error(pred[i], valid[i])))

输出:

rmse value for CO(GT) is : 1.4200393103392812 rmse value for PT08.S1(CO) is : 303.3909208229375 rmse value for NMHC(GT) is : 204.0662895081472 rmse value for C6H6(GT) is : 28.153391799471244 rmse value for PT08.S2(NMHC) is : 6.538063846286176 rmse value for NOx(GT) is : 265.04913993413805 rmse value for PT08.S3(NOx) is : 250.7673347152554 rmse value for NO2(GT) is : 238.92642219826683 rmse value for PT08.S4(NO2) is : 247.50612831072633 rmse value for PT08.S5(O3) is : 392.3129907890131 rmse value for T is : 383.1344361254454 rmse value for RH is : 506.5847387424092 rmse value for AH is : 8.139735443605728

预测代码：

#make final predictions model = VAR(endog=data) model_fit = model.fit() yhat = model_fit.forecast(model_fit.y, steps=1) print(yhat)

,