Time Series 简明教程
Time Series - Introduction
时间序列是某段时间内一系列的观察结果。单变量时间序列由一个变量在一段时间内按周期性时间实例获取的值组成,而多变量时间序列由多个变量在一段时间内在同一周期性时间实例获取的值组成。我们所有人每天都会遇到的最简单的时间序列示例是全天或全周或全月或全年的温度变化。
时间数据分析能够让我们深入了解一个变量如何随时间变化,或者它如何依赖其他变量值的改变。变量与其之前的值和/或其他变量之间的这种关系可以分析时间序列预测并在人工智能中得到众多应用。
Time Series - Programming Languages
Time Series - Python Libraries
由于其易于编写和易于理解的代码结构以及各种开源库,Python 在从事机器学习的人员中享有很高的知名度。我们在接下来章节中将使用的几个这样的开源库已在下面介绍。
NumPy
NumPy(Numerical Python)是一个用于科学计算的库。它处理一个 N 维数组对象,并提供基本的数学功能,如大小、形状、平均值、标准差、最小值、最大值以及更复杂的一些函数,如线性代数函数和傅里叶变换。随着我们在此教程中不断前进,您将了解更多相关内容。
Pandas
此库提供诸如系列、数据框和面板等高效且易于使用的数据结构。它提升了 Python 的功能,使其从单纯的数据收集和准备转变为数据分析。Pandas 和 NumPy 这两个库大大简化了针对从小到非常大的数据集的任何操作。想要进一步了解这些函数,请查看本教程。
SciPy
SciPy(Science Python)是一个用于科学和技术计算的库。它提供函数优化、信号和图像处理、积分、插值和线性代数功能。在执行机器学习时,此库非常有用。我们将在此教程中逐步讨论这些功能。
Scikit Learn
此库是一个 SciPy 工具包,广泛用于统计建模、机器学习和深度学习,因为它包含各种可定制的回归、分类和聚类模型。它可以很好地与 Numpy、Pandas 和其他库配合使用,从而更易于使用。
Time Series - Data Processing and Visualization
时间序列是在等距时间间隔内编制索引的一系列观测。因此,在任何时间序列中都应该保持顺序和连续性。
我们将使用的该数据集是一个多变量时间序列,它具有一个受严重污染的意大利城市空气质量的约一年的时均数据。可以从以下提供的链接下载该数据集 − https://archive.ics.uci.edu/ml/datasets/air+quality 。
必须确保 −
-
该时间序列是等距的,并且
-
其中没有冗余值或间隙。
如果时间序列不连续,我们可以对其上采样或下采样。
Showing df.head()
[122] 中:
import pandas
[123] 中:
df = pandas.read_csv("AirQualityUCI.csv", sep = ";", decimal = ",")
df = df.iloc[ : , 0:14]
[124] 中:
len(df)
Out[124]:
9471
In [125]:
df.head()
Out[125]:
对于时间序列的预处理,我们要确保数据集中没有 NaN(NULL) 值;如果有,我们可以将它们用 0 或平均值,或者前一个或后一个值替换。替换比丢弃是首选,以便保持时间序列的连续性。但是,在我们的数据集中,最后几个值似乎是 NULL,因此丢弃不会影响连续性。
Dropping NaN(Not-a-Number)
In [126]:
df.isna().sum()
Out[126]:
Date 114
Time 114
CO(GT) 114
PT08.S1(CO) 114
NMHC(GT) 114
C6H6(GT) 114
PT08.S2(NMHC) 114
NOx(GT) 114
PT08.S3(NOx) 114
NO2(GT) 114
PT08.S4(NO2) 114
PT08.S5(O3) 114
T 114
RH 114
dtype: int64
In [127]:
df = df[df['Date'].notnull()]
In [128]:
df.isna().sum()
Out[128]:
Date 0
Time 0
CO(GT) 0
PT08.S1(CO) 0
NMHC(GT) 0
C6H6(GT) 0
PT08.S2(NMHC) 0
NOx(GT) 0
PT08.S3(NOx) 0
NO2(GT) 0
PT08.S4(NO2) 0
PT08.S5(O3) 0
T 0
RH 0
dtype: int64
时间序列通常被描绘成时序折线图。为此,现在我们将日期和时间列结合起来,并将其从字符串转换为 datetime 对象。这可以使用 datetime 库完成。
Converting to datetime object
In [129]:
df['DateTime'] = (df.Date) + ' ' + (df.Time)
print (type(df.DateTime[0]))
<class 'str'>
In [130]:
import datetime
df.DateTime = df.DateTime.apply(lambda x: datetime.datetime.strptime(x, '%d/%m/%Y %H.%M.%S'))
print (type(df.DateTime[0]))
<class 'pandas._libs.tslibs.timestamps.Timestamp'>
让我们看看一些变量,例如随着时间的变化温度如何变化。
Showing plots
In [131]:
df.index = df.DateTime
In [132]:
import matplotlib.pyplot as plt
plt.plot(df['T'])
Out[132]:
[<matplotlib.lines.Line2D at 0x1eaad67f780>]
In [208]:
plt.plot(df['C6H6(GT)'])
Out[208]:
[<matplotlib.lines.Line2D at 0x1eaaeedff28>]
箱形图是另一种有用的图表,允许您将有关数据集的大量信息浓缩到单个图表中。它显示一个或多个变量的平均值、25% 和 75% 的四分位数以及异常值。当异常值数量较少并且与平均值相差很远时,我们可以通过将它们设置为平均值或 75% 四分位数来消除这些异常值。
Showing Boxplots
In [134]:
plt.boxplot(df[['T','C6H6(GT)']].values)
Out[134]:
{'whiskers': [<matplotlib.lines.Line2D at 0x1eaac16de80>,
<matplotlib.lines.Line2D at 0x1eaac16d908>,
<matplotlib.lines.Line2D at 0x1eaac177a58>,
<matplotlib.lines.Line2D at 0x1eaac177cf8>],
'caps': [<matplotlib.lines.Line2D at 0x1eaac16d2b0>,
<matplotlib.lines.Line2D at 0x1eaac16d588>,
<matplotlib.lines.Line2D at 0x1eaac1a69e8>,
<matplotlib.lines.Line2D at 0x1eaac1a64a8>],
'boxes': [<matplotlib.lines.Line2D at 0x1eaac16dc50>,
<matplotlib.lines.Line2D at 0x1eaac1779b0>],
'medians': [<matplotlib.lines.Line2D at 0x1eaac16d4a8>,
<matplotlib.lines.Line2D at 0x1eaac1a6c50>],
'fliers': [<matplotlib.lines.Line2D at 0x1eaac177dd8>,
<matplotlib.lines.Line2D at 0x1eaac1a6c18>],'means': []
}
Time Series - Modeling
Introduction
一个时间序列有如下 4 个组成部分:
-
Level − 它是一个平均值,序列围绕它变化。
-
Trend − 它是一个变量随着时间增长或减小的行为。
-
Seasonality − 它是一个时间序列的周期性行为。
-
Noise − 这是由于环境因素而添加到观测中的误差。
Time Series - Parameter Calibration
Introduction
任何统计或机器学习模型都有一些参数,这些参数极大地影响对数据的建模方式。例如,ARIMA 具有 p、d、q 值。这些参数将被决定,使得实际值和建模值之间的误差最小。参数校准被称为模型拟合中最关键和最耗时的任务。因此,为我们选择最优参数非常重要。
Methods for Calibration of Parameters
有各种方式校准参数。本节将详细讨论其中一些。
Hit-and-try
一种常见的校准模型的方法是手动校准,你在其中首先可视化时间序列,直观地尝试一些参数值并反复更改它们,直到达到足够好的拟合。它要求对我们尝试的模型有一个很好的理解。对于 ARIMA 模型,手工校准是借助自相关图来进行“p”参数、偏自相关图来进行“q”参数和 ADF 测试来确认时间序列的平稳性和设定“d”参数。我们将在接下来的章节中详细讨论所有这些。
Time Series - Naïve Methods
Introduction
朴素方法,例如假设时间“t”处的预测值为时间“t-1”处变量的实际值或序列的滚动平均值,用于权衡统计模型和机器学习模型的表现如何,并强调它们的需要。
让我们尝试在本教程中在时间序列数据的一个特性中使用这些模型。
首先,我们将查看我们数据的“温度”特性的平均值及其周围的偏差。了解最大值和最小值也很有用。我们可以在此使用 numpy 库的功能。
Showing statistics
In [135]:
import numpy
print (
'Mean: ',numpy.mean(df['T']), ';
Standard Deviation: ',numpy.std(df['T']),';
\nMaximum Temperature: ',max(df['T']),';
Minimum Temperature: ',min(df['T'])
)
我们有跨越等间隔时间线的所有 9357 项观察统计数据,我们可借此了解数据。
现在我们将尝试第一个朴素方法,将当前时间的预测值设置为前一时间点的实际值,并计算均方根误差 (RMSE) 来量化此方法的性能。
Showing 1st naïve method
In [136]:
df['T']
df['T_t-1'] = df['T'].shift(1)
[137] 中:
df_naive = df[['T','T_t-1']][1:]
[138] 中:
from sklearn import metrics
from math import sqrt
true = df_naive['T']
prediction = df_naive['T_t-1']
error = sqrt(metrics.mean_squared_error(true,prediction))
print ('RMSE for Naive Method 1: ', error)
质朴方法 1 的 RMSE:12.901140576492974
让我们看看下一个朴素方法,其中当前时间的预测值等同于 предшествующих时期的平均值。我们还将计算此方法的 RMSE。
Showing 2nd naïve method
[139] 中:
df['T_rm'] = df['T'].rolling(3).mean().shift(1)
df_naive = df[['T','T_rm']].dropna()
[140] 中:
true = df_naive['T']
prediction = df_naive['T_rm']
error = sqrt(metrics.mean_squared_error(true,prediction))
print ('RMSE for Naive Method 2: ', error)
RMSE for Naive Method 2: 14.957633272839242
在此处,您还可以尝试各种前面时间段(也称为“滞后”)的数量,您想要考虑这些数量,此处保留为 3。在该数据中,您可以看到随着滞后数的增加,误差也会增加。如果滞后保持为 1,它将成为与之前使用的质朴方法相同。
Points to Note
-
您可以编写一个非常简单的函数来计算均方根误差。在此处,我们已使用软件包 “sklearn” 中的均方误差函数,然后对其进行平方根计算。
-
在 pandas 中,df[“column_name”] 还可以写成 df.column_name,但是对于此数据集,df.T 与 df[“T”] 的工作方式不相同,因为 df.T 是用于转置数据框的函数。因此,仅使用 df[“T”] 或在使用其他语法之前考虑重命名该列。
Time Series - Auto Regression
对于平稳时间序列,自回归模型将时间“t”处的变量值视为其之前“p”时间步的值的线性函数。数学上可以写成以下形式:
y_ {t} = \:C+\:\phi_{1}y_{t-1}\:+\:\phi_{2}Y_{t-2}…\phi_{p}y_{t-p}+\epsilon_{t}
我无法使用 Gemini 翻译任何内容。
其中,“p”是自回归趋势参数
\epsilon_ {t} 是白噪声,并且
y_ {t-1},y_ {t-2} \:\: …y_ {t-p} 表示先前的时期变量的值。
可以使用多种方法校准 p 的值。找到“p”的适当值的一种方法是绘制自相关图。
Note - 在对数据执行任何分析之前,我们应该以 8:2 的可用总数据集比率将数据分割为训练和测试,因为测试数据只能找出我们模型的准确性,并且假设在作出预测之前我们无法获得该数据。对于时间序列,数据点的序列非常重要,因此在分割数据时应记住不要丢失顺序。
自相关图或相关图显示变量与其自身在先前的时步关系。它使用 Pearson 相关并且显示 95% 置信区间内的相关。让我们看看我们数据的“温度”变量是怎样的。
Time Series - Moving Average
对于平稳时间序列,移动平均模型将时间“t”处变量的值视为前“q”时间步长残差误差的线性函数。残差误差是通过将时间“t”处的值与前面值的移动平均值进行比较来计算的。
在数学上可以写成 −
y_{t} = c\:+\:\epsilon_{t}\:+\:\theta_{1}\:\epsilon_{t-1}\:+\:\theta_{2}\:\epsilon_{t-2}\:+\:…+:\theta_{q}\:\epsilon_{t-q}\:
其中“q”是移动平均趋势参数
\epsilon_ {t} 是白噪声,并且
$\epsilon_{t-1}, \epsilon_{t-2}…\epsilon_{t-q}$ 是前一时间段的误差项。
“q”的值可以使用多种方法进行校准。找到“q”的恰当值的一种方法是绘制偏自相关图。
与显示直接和间接相关性的自相关图不同,偏自相关图显示变量与其自身在之前时间步长的关系,同时消除了间接相关性,让我们看看它对我们数据的“temperature”变量有何影响。
Time Series - ARIMA
我们已经了解到,对于平稳时间序列,时间“t”处的变量是先前观测或残差误差的线性函数。因此,现在是时候将这两者结合起来,建立自回归移动平均 (ARMA) 模型了。
然而,有时时间序列不是平稳的,即序列的统计特性(如均值、方差)随时间变化。而我们迄今为止学习过的统计模型假设时间序列是平稳的,因此,我们可以包括差分时间序列的预处理步骤,使其平稳。现在,对于我们正在处理的时间序列是否是平稳的,我们必须找出答案。
用于查找时间序列平稳性的各种方法正在寻找时间序列图中的季节性或趋势,检查不同时间段的均值和方差差异、增强型迪基-福勒 (ADF) 检验、KPSS 检验、赫斯特指数等。
让我们使用 ADF 检验来确定数据集中的“温度”变量是否是平稳的时间序列。
In [74]:
from statsmodels.tsa.stattools import adfuller
result = adfuller(train)
print('ADF Statistic: %f' % result[0])
print('p-value: %f' % result[1])
print('Critical Values:')
for key, value In result[4].items()
print('\t%s: %.3f' % (key, value))
ADF 统计量:-10.406056
p 值:0.000000
临界值:
1%:-3.431
5%:-2.862
10%:-2.567
现在我们已经运行了 ADF 检验,让我们解释一下结果。首先,我们将 ADF 统计量与临界值进行比较,较低的临界值告诉我们该序列很可能是非平稳的。接下来,我们查看 p 值。大于 0.05 的 p 值也表明时间序列是非平稳的。
或者,p 值小于或等于 0.05,或 ADF 统计量小于临界值表明时间序列是平稳的。
因此,我们正在处理的时间序列已经平稳。在平稳时间序列的情况下,我们将“d”参数设置为 0。
我们还可以使用赫斯特指数确认时间序列的平稳性。
In [75]:
import hurst
H, c,data = hurst.compute_Hc(train)
print("H = {:.4f}, c = {:.4f}".format(H,c))
H = 0.1660,c = 5.0740
H<0.5 的值表示反持续性行为,H>0.5 表示持续性行为或趋势序列。H=0.5 表示随机游走/布朗运动。H<0.5 的值,确认我们的序列是平稳的。
对于非平稳时间序列,我们将“d”参数设置为 1。此外,自回归趋势参数“p”和移动平均趋势参数“q”的值是根据平稳时间序列计算的,即通过在对时间序列求差后绘制 ACP 和 PACP 来计算的。
ARIMA 模型的特点在于 3 个参数 (p,d,q),现在我们已经了解了它,因此让我们对我们的时间序列建模并预测温度的未来值。
In [156]:
from statsmodels.tsa.arima_model import ARIMA
model = ARIMA(train.values, order=(5, 0, 2))
model_fit = model.fit(disp=False)
In [157]:
predictions = model_fit.predict(len(test))
test_ = pandas.DataFrame(test)
test_['predictions'] = predictions[0:1871]
In [158]:
plt.plot(df['T'])
plt.plot(test_.predictions)
plt.show()
In [167]:
error = sqrt(metrics.mean_squared_error(test.values,predictions[0:1871]))
print ('Test RMSE for ARIMA: ', error)
ARIMA 的测试 RMSE:43.21252940234892
Time Series - Variations of ARIMA
在上一章中,我们已经看到了 ARIMA 模型如何工作,以及它不能处理季节性数据或多元时间序列的局限性,因此引入了新的模型来包含这些特征。
这里提供了这些新模型的概览 −
Vector Auto Regression Moving Average (VARMA)
它是 VAR 和 VMA 的组合以及多元固定时间序列 ARMA 模型的广义版本。它的特征是“p”和“q”参数。与 ARMA 一样,通过将“q”参数设置为 0 来充当 AR 模型并且通过将“p”参数设置为 0 来充当 MA 模型,VARMA 也能够通过将“q”参数设置为 0 来充当 VAR 模型并且通过将“p”参数设置为 0 来充当 VMA 模型。
[209] 中:
df_multi = df[['T', 'C6H6(GT)']]
split = len(df) - int(0.2*len(df))
train_multi, test_multi = df_multi[0:split], df_multi[split:]
[211] 中:
from statsmodels.tsa.statespace.varmax import VARMAX
model = VARMAX(train_multi, order = (2,1))
model_fit = model.fit()
c:\users\naveksha\appdata\local\programs\python\python37\lib\site-packages\statsmodels\tsa\statespace\varmax.py:152:
EstimationWarning: Estimation of VARMA(p,q) models is not generically robust,
due especially to identification issues.
EstimationWarning)
c:\users\naveksha\appdata\local\programs\python\python37\lib\site-packages\statsmodels\tsa\base\tsa_model.py:171:
ValueWarning: No frequency information was provided, so inferred frequency H will be used.
% freq, ValueWarning)
c:\users\naveksha\appdata\local\programs\python\python37\lib\site-packages\statsmodels\base\model.py:508:
ConvergenceWarning: Maximum Likelihood optimization failed to converge. Check mle_retvals
"Check mle_retvals", ConvergenceWarning)
[213] 中:
predictions_multi = model_fit.forecast( steps=len(test_multi))
c:\users\naveksha\appdata\local\programs\python\python37\lib\site-packages\statsmodels\tsa\base\tsa_model.py:320:
FutureWarning: Creating a DatetimeIndex by passing range endpoints is deprecated. Use `pandas.date_range` instead.
freq = base_index.freq)
c:\users\naveksha\appdata\local\programs\python\python37\lib\site-packages\statsmodels\tsa\statespace\varmax.py:152:
EstimationWarning: Estimation of VARMA(p,q) models is not generically robust, due especially to identification issues.
EstimationWarning)
[231] 中:
plt.plot(train_multi['T'])
plt.plot(test_multi['T'])
plt.plot(predictions_multi.iloc[:,0:1], '--')
plt.show()
plt.plot(train_multi['C6H6(GT)'])
plt.plot(test_multi['C6H6(GT)'])
plt.plot(predictions_multi.iloc[:,1:2], '--')
plt.show()
以上代码显示了如何使用 VARMA 模型对多元时间序列进行建模,虽然此模型可能不最适合我们的数据。
Seasonal Auto Regressive Integrated Moving Average (SARIMA)
这是 ARIMA 模型的扩展,用于处理季节性数据。它将数据划分为季节性和非季节性部分,并以类似的方式对它们进行建模。它的特征是 7 个参数,对于非季节部分(p、d、q)参数与 ARIMA 模型相同,对于季节部分(P、D、Q、m)参数,其中“m”是季节周期的数量,并且 P、D、Q 与 ARIMA 模型的参数类似。这些参数可以使用网格搜索或遗传算法进行校准。
SARIMA with Exogenous Variables (SARIMAX)
这是 SARIMA 模型的扩展,可用于包含外生变量,这有助于我们对我们感兴趣的变量进行建模。
在将变量作为外生变量使用之前,对变量进行相关分析可能很有用。
[251]:
from scipy.stats.stats import pearsonr
x = train_multi['T'].values
y = train_multi['C6H6(GT)'].values
corr , p = pearsonr(x,y)
print ('Corelation Coefficient =', corr,'\nP-Value =',p)
Corelation Coefficient = 0.9701173437269858
P-Value = 0.0
皮尔森相关性显示 2 个变量之间的线性关系,要解释结果,我们首先看 p 值,如果 p 值小于 0.05,则系数的值具有显着性,否则系数的值不具有显着性。对于有显着性的 p 值,相关系数的正值表示正相关,负值表示负相关。
因此,对于我们的数据,“温度”和“C6H6”似乎具有高度正相关性。因此,我们将
[297]:
from statsmodels.tsa.statespace.sarimax import SARIMAX
model = SARIMAX(x, exog = y, order = (2, 0, 2), seasonal_order = (2, 0, 1, 1), enforce_stationarity=False, enforce_invertibility = False)
model_fit = model.fit(disp = False)
c:\users\naveksha\appdata\local\programs\python\python37\lib\site-packages\statsmodels\base\model.py:508:
ConvergenceWarning: Maximum Likelihood optimization failed to converge. Check mle_retvals
"Check mle_retvals", ConvergenceWarning)
[298]:
y_ = test_multi['C6H6(GT)'].values
predicted = model_fit.predict(exog=y_)
test_multi_ = pandas.DataFrame(test)
test_multi_['predictions'] = predicted[0:1871]
[299]:
plt.plot(train_multi['T'])
plt.plot(test_multi_['T'])
plt.plot(test_multi_.predictions, '--')
Out[299]:
[<matplotlib.lines.Line2D at 0x1eab0191c18>]
与单变量 ARIMA 建模相比,这里的预测似乎变化更大。
不用说,SARIMAX 可通过仅将相应参数设置为非零值来用作 ARX、MAX、ARMAX 或 ARIMAX 模型。
Time Series - Exponential Smoothing
在本章中,我们将讨论时间序列指数平滑涉及的技术。
Simple Exponential Smoothing
指数平滑是一种通过在一段时间内对数据分配指数递减的权重来平滑单变量时间序列的技术。
在数学上,给定时间 t 时的值 y_(t+1|t) 时,时间“t+1”时的变量值定义为 −
y_{t+1|t}\:=\:\alpha y_{t}\:+\:\alpha\lgroup1 -\alpha\rgroup y_{t-1}\:+\alpha\lgroup1-\alpha\rgroup^{2}\:y_{t-2}\:+\:…+y_{1}
其中,0≤α≤1 是平滑参数,并且
$y_{1},….,y_{t}$ 是网络流量在时间点 1、2、3、… 、t 的前值。
这是建模没有明显趋势或季节性的时间序列的一种简单方法。但指数平滑也可用于具有趋势和季节性的时间序列。
Triple Exponential Smoothing
三重指数平滑 (TES) 或霍尔特冬季方法应用指数平滑三次 - 水平平滑 $l_{t}$、趋势平滑 $b_{t}$ 和季节性平滑 $S_{t}$,其中 $\alpha$, $\beta^{*}$ 和 $\gamma$ 作为平滑参数,“m”为季节性频率,即一年中的季节数。
根据季节性分量的性质,TES 有两种类别−
-
Holt-Winter’s Additive Method − 当季节性本质上是加性的。
-
Holt-Winter’s Multiplicative Method − 当季节性本质上是乘性的。
对于非季节性时间序列,我们只有趋势平滑和水平平滑,这称为霍尔特线性趋势法。
让我们尝试对我们的数据应用三重指数平滑。
输入[316]:
from statsmodels.tsa.holtwinters import ExponentialSmoothing
model = ExponentialSmoothing(train.values, trend= )
model_fit = model.fit()
输入[322]:
predictions_ = model_fit.predict(len(test))
输入[325]:
plt.plot(test.values)
plt.plot(predictions_[1:1871])
输出[325]:
[<matplotlib.lines.Line2D at 0x1eab00f1cf8>]
在这里,我们用训练集训练了一次模型,然后我们继续做出预测。一种更现实的方法是在一个或多个时间步之后重新训练模型。当我们从训练数据“直到时间‘t’”获得时间“t+1”的预测时,下一次时间“t+2”的预测可以使用训练数据“直到时间‘t+1’”来做出,因为那时将知道时间“t+1”的实际值。这种为一个或多个未来步做出预测然后重新训练模型的方法称为滚动预测或向前验证。
Time Series - Walk Forward Validation
在时间序列建模中,随着时间的推移,预测会变得越来越不准确,因此根据实际数据重新训练模型是一种更为现实的方法,因为它可以用于进一步预测。由于统计模型的训练并不耗时,因此按步进验证是获得最准确结果的最优解决方案。
我们对数据应用一步按步进验证,并将其与我们之前获得的结果进行比较。
[333] 中:
prediction = []
data = train.values
for t In test.values:
model = (ExponentialSmoothing(data).fit())
y = model.predict()
prediction.append(y[0])
data = numpy.append(data, t)
[335] 中:
test_ = pandas.DataFrame(test)
test_['predictionswf'] = prediction
[341] 中:
plt.plot(test_['T'])
plt.plot(test_.predictionswf, '--')
plt.show()
[340] 中:
error = sqrt(metrics.mean_squared_error(test.values,prediction))
print ('Test RMSE for Triple Exponential Smoothing with Walk-Forward Validation: ', error)
Test RMSE for Triple Exponential Smoothing with Walk-Forward Validation: 11.787532205759442
我们可以看到,现在我们的模型执行得明显更好。事实上,趋势被跟踪得如此紧密,以至于在图中预测与实际值重叠。你也可以尝试对 ARIMA 模型应用按步进验证。
Time Series - Prophet Model
2017年,Facebook开源了Prophet模型,该模型能够对具有日、周、年等多强季度的时序进行建模,以及对趋势建模。它具有直观的参数,非专业级的数据科学家可以对其进行调整,以获得更好的预测。它的核心是一个加法回归模型,可以检测时序的变化点。
Prophet将时序分解为趋势分量$g_{t}$,季节分量$S_{t}$和节假日分量$h_{t}$。
$y_{t}=g_{t}s_{t}+h_{t}\epsilon_{t}$
其中,$\epsilon_{t}$ 是误差项。
类似的时间序列预测包(例如因果影响和异常检测)分别由谷歌和推特在 R 中引入。
Time Series - LSTM Model
现在,我们已熟悉时间序列的统计建模,机器学习目前盛行,所以也务必要熟悉一些机器学习模型。我们从时间序列领域最流行的模型——长短期记忆模型开始。
LSTM 是一类循环神经网络。因此,在跳到 LSTM 之前,务必要了解神经网络和循环神经网络。
Recurrent Neural Networks
这是一类专门用于处理临时数据的类神经网络。RNN 的神经元具有单元状态/记忆,输入根据此内部状态进行处理,这是通过神经网络中的循环来完成的。RNN 中的“tanh”层有周期性模块,可让它们保留信息。然而,保留的时间不会很长,因此我们才需要 LSTM 模型。
LSTM
这是一种特殊类别的循环神经网络,能够学习数据中的长期依赖性。这是因为模型的周期性模块将四个层层进行交互。
上图以黄色框表示四个神经网络层,以绿色圆圈表示逐点运算符,以黄色圆圈表示输入,以蓝色圆圈表示单元状态。一个 LSTM 模块具有一个单元状态和三个门,这给了它有选择地学习、取消学习,或者保留每个单元信息的能力。LSTM 中的单元状态帮助信息在单元之间流动,并通过只允许少量线性交互而保持不变。每个单元具有输入、输出,以及可以将信息添加到单元状态、或者从单元状态中移除信息的遗忘门。遗忘门使用一个 sigmoid 函数来确定前一个单元状态中的哪些信息应该被忘记。输入门使用“sigmoid”和“tanh”的逐点乘法运算来控制信息流向当前单元状态。最后,输出门决定哪些信息应传递到下一个隐藏状态
现在我们已经了解了 LSTM 模型的内部工作原理,让我们来实施它。要了解 LSTM 的实施,我们从一个简单的示例开始——一条直线。我们来看看 LSTM 能否学习直线关系并预测它。
首先,我们创建描述一条直线的数据集。
In [402]:
x = numpy.arange (1,500,1)
y = 0.4 * x + 30
plt.plot(x,y)
Out[402]:
[<matplotlib.lines.Line2D at 0x1eab9d3ee10>]
In [403]:
trainx, testx = x[0:int(0.8*(len(x)))], x[int(0.8*(len(x))):]
trainy, testy = y[0:int(0.8*(len(y)))], y[int(0.8*(len(y))):]
train = numpy.array(list(zip(trainx,trainy)))
test = numpy.array(list(zip(trainx,trainy)))
现在数据已经创建并分成训练和测试。让我们将时间序列数据转换为监督学习数据,根据回看期的值,这实质上是用于预测时间 t 值的滞后数。
因此像这样的时间序列 −
time variable_x
t1 x1
t2 x2
: :
: :
T xT
在回溯周期为 1 时,转化为 −
x1 x2
x2 x3
: :
: :
xT-1 xT
In [404]:
def create_dataset(n_X, look_back):
dataX, dataY = [], []
for i in range(len(n_X)-look_back):
a = n_X[i:(i+look_back), ]
dataX.append(a)
dataY.append(n_X[i + look_back, ])
return numpy.array(dataX), numpy.array(dataY)
In [405]:
look_back = 1
trainx,trainy = create_dataset(train, look_back)
testx,testy = create_dataset(test, look_back)
trainx = numpy.reshape(trainx, (trainx.shape[0], 1, 2))
testx = numpy.reshape(testx, (testx.shape[0], 1, 2))
现在,我们将训练模型。
少量训练数据显示给网络,当所有训练数据分批显示给模型并且计算误差时称为一次 epoch。Epochs 将继续运行直到误差减少为止。
In []:
from keras.models import Sequential
from keras.layers import LSTM, Dense
model = Sequential()
model.add(LSTM(256, return_sequences = True, input_shape = (trainx.shape[1], 2)))
model.add(LSTM(128,input_shape = (trainx.shape[1], 2)))
model.add(Dense(2))
model.compile(loss = 'mean_squared_error', optimizer = 'adam')
model.fit(trainx, trainy, epochs = 2000, batch_size = 10, verbose = 2, shuffle = False)
model.save_weights('LSTMBasic1.h5')
In [407]:
model.load_weights('LSTMBasic1.h5')
predict = model.predict(testx)
现在,让我们看看我们的预测值是什么。
In [408]:
plt.plot(testx.reshape(398,2)[:,0:1], testx.reshape(398,2)[:,1:2])
plt.plot(predict[:,0:1], predict[:,1:2])
Out[408]:
[<matplotlib.lines.Line2D at 0x1eac792f048>]
现在,我们应尝试以类似方式对正弦波或余弦波建模。您可以运行下面提供的代码并使用模型参数进行操作,以查看结果如何变化。
In [409]:
x = numpy.arange (1,500,1)
y = numpy.sin(x)
plt.plot(x,y)
Out[409]:
[<matplotlib.lines.Line2D at 0x1eac7a0b3c8>]
In [410]:
trainx, testx = x[0:int(0.8*(len(x)))], x[int(0.8*(len(x))):]
trainy, testy = y[0:int(0.8*(len(y)))], y[int(0.8*(len(y))):]
train = numpy.array(list(zip(trainx,trainy)))
test = numpy.array(list(zip(trainx,trainy)))
In [411]:
look_back = 1
trainx,trainy = create_dataset(train, look_back)
testx,testy = create_dataset(test, look_back)
trainx = numpy.reshape(trainx, (trainx.shape[0], 1, 2))
testx = numpy.reshape(testx, (testx.shape[0], 1, 2))
In []:
model = Sequential()
model.add(LSTM(512, return_sequences = True, input_shape = (trainx.shape[1], 2)))
model.add(LSTM(256,input_shape = (trainx.shape[1], 2)))
model.add(Dense(2))
model.compile(loss = 'mean_squared_error', optimizer = 'adam')
model.fit(trainx, trainy, epochs = 2000, batch_size = 10, verbose = 2, shuffle = False)
model.save_weights('LSTMBasic2.h5')
In [413]:
model.load_weights('LSTMBasic2.h5')
predict = model.predict(testx)
In [415]:
plt.plot(trainx.reshape(398,2)[:,0:1], trainx.reshape(398,2)[:,1:2])
plt.plot(predict[:,0:1], predict[:,1:2])
Out [415]:
[<matplotlib.lines.Line2D at 0x1eac7a1f550>]
现在,您可以继续处理任意数据集了。
Time Series - Error Metrics
对模型的性能进行量化以便将其用作反馈和比较对我们来说很重要。在本教程中,我们使用了最常见的错误指标之一:均方根误差。还有其他各种可用错误指标。本章将简要讨论它们。
Mean Square Error
它是预测值与真值之间差值的平方平均值。Sklearn 以函数形式提供它。它的单位与真值和预测值的平方相同,并且始终为正。
MSE = \frac{1}{n} displaystyle\sum\limits_{t=1}^n \lgroup y' {t}\:-y {t}\rgroup^{2}
其中,$y'_{t}$ 是预测值,
$y_{t}$ 是实际值,
n 是测试集中值的总数。
从方程式中可以清楚地看出,对于较大误差或异常值,MSE 具有更大的惩罚性。
Root Mean Square Error
它是均方误差的平方根。它也总是为正,并且在数据范围内。
RMSE = \sqrt{\frac{1}{n} displaystyle\sum\limits_{t=1}^n \lgroup y' {t}-y {t}\rgroup ^2}
其中,$y'_{t}$ 是预测值,
$y_{t}$ 表示实际值,而
n 是测试集中值的总数。
它具有单位功率,因此与 MSE 相比更具可解释性。RMSE 对较大误差的惩罚也更大。我们在教程中使用了 RMSE 度量。
Mean Absolute Error
这是预测值和真值之间绝对差值的平均值。它的单位与预测值和真值相同,且始终为正。
MAE = \frac{1}{n}\displaystyle\sum\limits_{t=1}^{t=n} | y'{t}-y_{t}\lvert
其中,$y'_{t}$ 是预测值,
$y_{t}$ 表示实际值,而
n 是测试集中值的总数。
Time Series - Applications
我们在本教程中讨论了时间序列分析,它使我们理解了时间序列模型首先根据现有观测识别趋势和季节性,然后基于该趋势和季节性预测一个值。这种分析在各个领域比如以下领域中很有用 −
-
Financial Analysis − 包括销售预测、库存分析、股票市场分析、价格估算。
-
Weather Analysis − 包括温度估算、气候变化、季节性变化识别、天气预测。
-
Network Data Analysis − 包括网络使用率预测、异常或入侵检测、预测性维护。
-
Healthcare Analysis − 包括人口普查预测、保险福利预测、病人监测。
Time Series - Further Scope
机器学习处理各种问题。实际上,几乎所有领域都具有借助机器学习实现自动化或改进的范围。正在对此进行大量工作的此类一些问题如下。
Non-Time Series Data
它是不随时间变化的数据,大多数 ML 问题都是非时间序列数据。为简单起见,我们将其进一步分类为 -
-
Numerical Data − 计算机不同于人类,只理解数字,因此对于机器学习来说,最终所有类型的数据都会转换为数字数据,例如图像数据转换为 (r,b,g) 值,字符转换为 ASCII 代码,单词转换为数字,语音数据转换为包含数字数据的文件 mfcc。
-
Image Data − 计算机视觉彻底改变了计算机的世界,它在医学、卫星成像等领域有广泛应用。
-
Text Data − 自然语言处理 (NLP) 用于文本分类、同义句检测和语言摘要。这正是使 Google 和 Facebook变得智能的原因。
-
Speech Data − 语音处理涉及语音识别和情感理解。它在赋予计算机类人特质方面发挥着至关重要的作用。