Machine Learning Mastery 深度学习时间序列教程（十三）
用于时间序列预测的有状态 LSTM 在线学习的不稳定性
原文：
machinelearningmastery.com/instability-online-learning-stateful-lstm-time-series-forecasting/
一些神经网络配置可能导致模型不稳定。
这可能会使用描述性统计量难以表征和比较同一问题的其他模型配置。
看似不稳定模型的一个很好的例子是使用有状态长短期记忆（LSTM）模型的在线学习（批量大小为 1）。
在本教程中，您将了解如何使用在线学习标准时间序列预测问题来探索有状态 LSTM 拟合的结果。
完成本教程后，您将了解：
- 如何设计一个强大的测试工具来评估 LSTM 模型的时间序列预测问题。
- 如何分析结果的总体，包括摘要统计，传播和结果分布。
- 如何分析增加实验重复次数的影响。
让我们开始吧。
用于时间序列预测的有状态 LSTM 在线学习的不稳定性 Magnus Brath 的照片，保留一些权利。
模型不稳定
当您在同一数据上多次训练相同的网络时，您可能会得到截然不同的结果。
这是因为神经网络是随机初始化的，并且它们如何适合训练数据的优化性质可导致网络内的不同最终权重。在给定相同输入数据的情况下，这些不同的网络可以反过来导致变化的预测。
因此，重复在神经网络上重复任何实验以找到平均预期表现是很重要的。
有关神经网络等机器学习算法的随机性的更多信息，请参阅帖子：
神经网络中的批量大小定义了在暴露于训练数据集的情况下网络中权重的更新频率。
批量大小为 1 意味着在每行训练数据之后更新网络权重。这称为在线学习。结果是一个可以快速学习的网络，但配置可能非常不稳定。
在本教程中，我们将探讨用于时间序列预测的有状态 LSTM 配置的在线学习的不稳定性。
我们将通过查看 LSTM 配置在标准时间序列预测问题上的平均表现来探索这一问题，该问题涉及可变数量的实验重复。
也就是说，我们将多次在相同数据上重新训练相同的模型配置，并在保留数据集上查看模型的表现，并查看模型的不稳定性。
教程概述
本教程分为 6 个部分。他们是：
- 洗发水销售数据集
- 实验测试线束
- 代码和收集结果
- 结果基本情况
- 重复与测试 RMSE
- 审查结果
环境
本教程假定您已安装 Python SciPy 环境。您可以在此示例中使用 Python 2 或 3。
本教程假设您安装了 TensorFlow 或 Theano 后端的 Keras v2.0 或更高版本。
本教程还假设您安装了 scikit-learn，Pandas，NumPy 和 Matplotlib。
接下来，让我们看看标准时间序列预测问题，我们可以将其用作此实验的上下文。
如果您在设置 Python 环境时需要帮助，请参阅以下帖子：
洗发水销售数据集
该数据集描述了 3 年期间每月洗发水的销售数量。
单位是销售计数，有 36 个观察。原始数据集归功于 Makridakis，Wheelwright 和 Hyndman（1998）。
下面的示例加载并创建已加载数据集的图。
# load and plot dataset
from pandas import read_csv
from pandas import datetime
from matplotlib import pyplot
# load dataset
def parser(x):
return datetime.strptime('190'+x, '%Y-%m')
series = read_csv('shampoo-sales.csv', header=0, parse_dates=[0], index_col=0, squeeze=True, date_parser=parser)
# summarize first few rows
print(series.head())
# line plot
series.plot()
pyplot.show()
运行该示例将数据集作为 Pandas Series 加载并打印前 5 行。
Month
1901-01-01 266.0
1901-02-01 145.9
1901-03-01 183.1
1901-04-01 119.3
1901-05-01 180.3
Name: Sales, dtype: float64
然后创建该系列的线图，显示明显的增加趋势。
洗发水销售数据集的线图
接下来，我们将了解实验中使用的 LSTM 配置和测试工具。
实验测试线束
本节介绍本教程中使用的测试工具。
数据拆分
我们将 Shampoo Sales 数据集分为两部分：训练和测试集。
前两年的数据将用于训练数据集，剩余的一年数据将用于测试集。
将使用训练数据集开发模型，并对测试数据集做出预测。
测试数据集的持久性预测（朴素预测）实现了每月洗发水销售 136.761 的错误。这在测试集上提供了较低的可接受表现限制。
模型评估
将使用滚动预测场景，也称为前进模型验证。
测试数据集的每个时间步骤将一次一个地走。将使用模型对时间步长做出预测，然后将获取测试集的实际预期值，并使其可用于下一时间步的预测模型。
这模仿了一个真实世界的场景，每个月都会有新的洗发水销售观察结果，并用于下个月的预测。
这将通过训练和测试数据集的结构进行模拟。
将收集关于测试数据集的所有预测，并计算错误分数以总结模型的技能。将使用均方根误差（RMSE），因为它会对大错误进行处罚，并产生与预测数据相同的分数，即每月洗发水销售额。
数据准备
在我们将 LSTM 模型拟合到数据集之前，我们必须转换数据。
在拟合模型和做出预测之前，对数据集执行以下三个数据变换。
- 转换时间序列数据，使其静止。具体而言，滞后= 1 差分以消除数据中的增加趋势。
- 将时间序列转换为监督学习问题。具体而言，将数据组织成输入和输出模式，其中前一时间步的观察被用作预测当前时间步的观察的输入
- 将观察结果转换为具有特定比例。具体而言，要将数据重缩放为-1 到 1 之间的值，以满足 LSTM 模型的默认双曲正切激活函数。
这些变换在预测时反转，在计算和误差分数之前将它们恢复到原始比例。
LSTM 模型
我们将使用基础状态 LSTM 模型，其中 1 个神经元适合 1000 个时期。
批量大小为 1 是必需的，因为我们将使用前向验证并对最后 12 个月的测试数据进行一步预测。
批量大小为 1 意味着该模型将使用在线训练（而不是批量训练或小批量训练）。因此，预计模型拟合将具有一些变化。
理想情况下，将使用更多的训练时期（例如 1500），但这被截断为 1000 以保持运行时间合理。
使用有效的 ADAM 优化算法和均方误差损失函数来拟合模型。
实验运行
每个实验场景将运行 100 次，并且测试集上的 RMSE 得分将从每次运行结束时记录。
所有测试 RMSE 分数都写入文件以供以后分析。
让我们深入研究实验。
代码和收集结果
完整的代码清单如下。
在现代硬件上运行可能需要几个小时。
from pandas import DataFrame
from pandas import Series
from pandas import concat
from pandas import read_csv
from pandas import datetime
from sklearn.metrics import mean_squared_error
from sklearn.preprocessing import MinMaxScaler
from keras.models import Sequential
from keras.layers import Dense
from keras.layers import LSTM
from math import sqrt
import matplotlib
import numpy
from numpy import concatenate
# date-time parsing function for loading the dataset
def parser(x):
return datetime.strptime('190'+x, '%Y-%m')
# frame a sequence as a supervised learning problem
def timeseries_to_supervised(data, lag=1):
df = DataFrame(data)
columns = [df.shift(i) for i in range(1, lag+1)]
columns.append(df)
df = concat(columns, axis=1)
return df
# create a differenced series
def difference(dataset, interval=1):
diff = list()
for i in range(interval, len(dataset)):
value = dataset[i] - dataset[i - interval]
diff.append(value)
return Series(diff)
# invert differenced value
def inverse_difference(history, yhat, interval=1):
return yhat + history[-interval]
# scale train and test data to [-1, 1]
def scale(train, test):
# fit scaler
scaler = MinMaxScaler(feature_range=(-1, 1))
scaler = scaler.fit(train)
# transform train
train = train.reshape(train.shape[0], train.shape[1])
train_scaled = scaler.transform(train)
# transform test
test = test.reshape(test.shape[0], test.shape[1])
test_scaled = scaler.transform(test)
return scaler, train_scaled, test_scaled
# inverse scaling for a forecasted value
def invert_scale(scaler, X, yhat):
new_row = [x for x in X] + [yhat]
array = numpy.array(new_row)
array = array.reshape(1, len(array))
inverted = scaler.inverse_transform(array)
return inverted[0, -1]
# fit an LSTM network to training data
def fit_lstm(train, batch_size, nb_epoch, neurons):
X, y = train[:, 0:-1], train[:, -1]
X = X.reshape(X.shape[0], 1, X.shape[1])
model = Sequential()
model.add(LSTM(neurons, batch_input_shape=(batch_size, X.shape[1], X.shape[2]), stateful=True))
model.add(Dense(1))
model.compile(loss='mean_squared_error', optimizer='adam')
for i in range(nb_epoch):
model.fit(X, y, epochs=1, batch_size=batch_size, verbose=0, shuffle=False)
model.reset_states()
return model
# make a one-step forecast
def forecast_lstm(model, batch_size, X):
X = X.reshape(1, 1, len(X))
yhat = model.predict(X, batch_size=batch_size)
return yhat[0,0]
# run a repeated experiment
def experiment(repeats, series):
# transform data to be stationary
raw_values = series.values
diff_values = difference(raw_values, 1)
# transform data to be supervised learning
supervised = timeseries_to_supervised(diff_values, 1)
supervised_values = supervised.values[1:,:]
# split data into train and test-sets
train, test = supervised_values[0:-12, :], supervised_values[-12:, :]
# transform the scale of the data
scaler, train_scaled, test_scaled = scale(train, test)
# run experiment
error_scores = list()
for r in range(repeats):
# fit the base model
lstm_model = fit_lstm(train_scaled, 1, 1000, 1)
# forecast test dataset
predictions = list()
for i in range(len(test_scaled)):
# predict
X, y = test_scaled[i, 0:-1], test_scaled[i, -1]
yhat = forecast_lstm(lstm_model, 1, X)
# invert scaling
yhat = invert_scale(scaler, X, yhat)
# invert differencing
yhat = inverse_difference(raw_values, yhat, len(test_scaled)+1-i)
# store forecast
predictions.append(yhat)
# report performance
rmse = sqrt(mean_squared_error(raw_values[-12:], predictions))
print('%d) Test RMSE: %.3f' % (r+1, rmse))
error_scores.append(rmse)
return error_scores
# execute the experiment
def run():
# load dataset
series = read_csv('shampoo-sales.csv', header=0, parse_dates=[0], index_col=0, squeeze=True, date_parser=parser)
# experiment
repeats = 100
results = DataFrame()
# run experiment
results['results'] = experiment(repeats, series)
# summarize results
print(results.describe())
# save results
results.to_csv('experiment_stateful.csv', index=False)
# entry point
run()
运行实验会在测试数据集中保存拟合模型的 RMSE 分数。
结果保存到文件“experiment_stateful.csv”。
下面提供了截断的结果列表。
自己重新运行实验可能会给出不同的结果，因为我们没有为随机数生成器播种。
...
116.39769471284067
105.0459745537738
93.57827109861229
128.973001927212
97.02915084460737
198.56877142225886
113.09568645243242
97.84127724751188
124.60413895331735
111.62139008607713
结果基本情况
我们可以从计算 100 个测试 RMSE 分数的整个人口的一些基本统计量开始。
通常，我们希望机器学习结果具有高斯分布。这允许我们报告模型的平均值和标准偏差，并在对看不见的数据做出预测时指示模型的置信区间。
下面的代码段会加载结果文件并计算一些描述性统计量。
from pandas import DataFrame
from pandas import read_csv
from numpy import mean
from numpy import std
from matplotlib import pyplot
# load results file
results = read_csv('experiment_stateful.csv', header=0)
# descriptive stats
print(results.describe())
# box and whisker plot
results.boxplot()
pyplot.show()
运行该示例将从结果中打印描述性统计量。
我们可以看到，平均而言，该配置实现了约 107 个月洗发水销售的 RMSE，标准偏差约为 17。
我们还可以看到观察到的最佳测试 RMSE 大约是 90 个销售，而更差的是不到 200 个，这是相当多的分数。
results
count 100.000000
mean 107.051146
std 17.694512
min 90.407323
25% 96.630800
50% 102.603908
75% 111.199574
max 198.568771
为了更好地了解数据的传播，还创建了一个盒子和胡须图。
该图显示了中位数（绿线），中间 50％的数据（框）和异常值（点）。我们可以看到数据差异很大，导致 RMSE 评分不佳。
在洗发水销售数据集上的 100 个测试 RMSE 分数的盒子和晶须图
还创建原始结果值的直方图。
该图表示倾斜甚至是指数分布，其质量在 RMSE 为 100 左右，长尾导向 RMSE 为 200。
结果的分布显然不是高斯分布。这是不幸的，因为平均值和标准偏差不能直接用于估计模型的置信区间（例如 95％置信度为平均值附近标准偏差的 2 倍）。
偏态分布还强调，中位数（第 50 百分位数）将是更好的集中趋势，而不是使用这些结果的均值。对于异常结果，中位数应该比平均值更稳健。
洗发水销售数据集中测试 RMSE 评分的直方图
重复与测试 RMSE
我们可以开始研究实验的摘要统计量如何随着重复次数从 1 增加到 100 而改变。
我们可以累积测试 RMSE 分数并计算描述性统计量。例如，来自一次重复的得分，来自第一次和第二次重复的得分，来自前三次重复的得分，等等至 100 次重复。
我们可以回顾一下中心趋势如何随着线图的重复次数的增加而变化。我们将看看均值和中位数。
一般来说，我们预计随着实验重复次数的增加，分布将越来越好地与基础分布相匹配，包括中心趋势，如均值。
完整的代码清单如下。
from pandas import DataFrame
from pandas import read_csv
from numpy import median
from numpy import mean
from matplotlib import pyplot
import numpy
# load results file
results = read_csv('experiment_stateful.csv', header=0)
values = results.values
# collect cumulative stats
medians, means = list(), list()
for i in range(1,len(values)+1):
data = values[0:i, 0]
mean_rmse, median_rmse = mean(data), median(data)
means.append(mean_rmse)
medians.append(median_rmse)
print(i, mean_rmse, median_rmse)
# line plot of cumulative values
line1, = pyplot.plot(medians, label='Median RMSE')
line2, = pyplot.plot(means, label='Mean RMSE')
pyplot.legend(handles=[line1, line2])
pyplot.show()
随着重复次数的增加，打印分布的累积大小，平均值和中值。截断的输出如下所示。
...
90 105.759546832 101.477640071
91 105.876449555 102.384620485
92 105.867422653 102.458057114
93 105.735281239 102.384620485
94 105.982491033 102.458057114
95 105.888245347 102.384620485
96 106.853667494 102.458057114
97 106.918018205 102.531493742
98 106.825398399 102.458057114
99 107.004981637 102.531493742
100 107.051145721 102.603907965
还创建了线图，显示随着重复次数的增加，平均值和中位数如何变化。
结果表明，正如预期的那样，平均值受异常值结果的影响大于中位数。
我们可以看到中位数在 99-100 附近看起来相当稳定。在情节结束时跳跃到 102，表明在后来的重复中出现了一系列较差的 RMSE 分数。
平均值和中值测试 RMSE 与重复次数的线图
审查结果
我们根据标准时间序列预测问题对有状态 LSTM 的 100 次重复进行了一些有用的观察。
特别：
- 我们观察到结果的分布不是高斯分布。它可能是倾斜的高斯分布或具有长尾和异常值的指数分布。
- 我们观察到结果的分布不随着重复从 1 增加到 100 而稳定。
观察结果表明了一些重要的特性：
- LSTM 和问题的在线学习的选择导致相对不稳定的模型。
- 所选择的重复次数（100）可能不足以表征模型的行为。
这是一个有用的发现，因为从 100 次或更少的实验重复中对模型做出有力的结论是错误的。
在描述您自己的机器学习结果时，这是一个需要考虑的重要注意事项。
这表明该实验有一些扩展，例如：
- 探索重复次数对更稳定模型的影响，例如使用批量或小批量学习的模型。
- 将重复次数增加到数千或更多，以试图通过在线学习解释模型的一般不稳定性。
摘要
在本教程中，您了解了如何使用在线学习分析 LSTM 模型的实验结果。
你了解到：
- 如何设计一个强大的测试工具来评估 LSTM 模型的时间序列预测问题。
- 如何分析实验结果，包括汇总统计。
- 如何分析增加实验重复次数的影响以及如何识别不稳定模型。
你有任何问题吗？ 在下面的评论中提出您的问题，我会尽力回答。
用于罕见事件时间序列预测的 LSTM 模型架构
原文：
machinelearningmastery.com/lstm-model-architecture-for-rare-event-time-series-forecasting/
使用 LSTM 直接进行时间序列预测几乎没有成功。
这是令人惊讶的，因为已知神经网络能够学习复杂的非线性关系，并且 LSTM 可能是能够直接支持多变量序列预测问题的最成功的循环神经网络类型。
最近在 Uber AI Labs 上进行的一项研究表明，LSTM 的自动特征学习功能及其处理输入序列的能力如何在端到端模型中得到利用，可用于驱动需求预测适用于公众假期等罕见事件。
在本文中，您将发现一种为时间序列预测开发可扩展的端到端 LSTM 模型的方法。
阅读这篇文章后，你会知道：
- 跨多个站点的多变量，多步骤预测的挑战，在这种情况下是城市。
- 用于时间序列预测的 LSTM 模型架构，包括单独的自编码器和预测子模型。
- 所提出的 LSTM 架构在罕见事件中的技能需求预测以及在不相关的预测问题上重用训练模型的能力。
让我们开始吧。
概观
在这篇文章中，我们将通过 Nikolay Laptev 等人回顾 2017 年题为“Uber 神经网络的时间序列极端事件预测”的论文。在 ICML 2017 时间序列研讨会上发表。
这篇文章分为四个部分;他们是：
- 动机
- 数据集
- 模型
- 发现
动机
该工作的目标是为多步骤时间序列预测开发端到端预测模型，该模型可以处理多变量输入（例如，多输入时间序列）。
该模型的目的是预测优步驾驶员对乘车共享的需求，特别是预测具有挑战性日子的需求，例如假期，其中经典模型的不确定性很高。
通常，这种类型的假期需求预测属于称为极端事件预测的研究领域。
极端事件预测已成为估算乘车共享和其他应用的峰值电力需求，交通拥堵严重性和激增定价的热门话题。事实上，有一个称为极值理论（EVT）的统计分支直接处理这一挑战。
- 优步神经网络的时间序列极端事件预测，2017。
描述了两种现有方法：
- 经典预测方法：每个时间序列开发模型时，可能根据需要拟合。
- 两步法：经典模型与机器学习模型结合使用。
这些现有模型的难度激发了对单个端到端模型的需求。
此外，还需要一个可以跨区域推广的模型，特别是针对每个城市收集的数据。这意味着在一些或所有城市训练的模型可用数据并用于在一些或所有城市做出预测。
我们可以将此概括为一个模型的一般需求，该模型支持多变量输入，进行多步预测，并在多个站点（在这种情况下为城市）中进行概括。
数据集
该模型适用于 Uber 数据集，该数据集包括美国顶级城市五年的匿名乘车共享数据。
在人口方面，美国各大城市完成旅行的五年历史记录用于提供美国所有主要假期的预测。
- 优步神经网络的时间序列极端事件预测，2017。
每个预测的输入包括每个骑行的信息，以及天气，城市和假日变量。
为了避免缺乏数据，我们使用其他功能，包括天气信息（例如降水，风速，温度）和城市级信息（例如，当前旅行，当前用户，当地假期）。
- 优步神经网络的时间序列极端事件预测，2017。
下面的图表提供了一年六个变量的样本。
模型 的缩放多变量输入来自“优步神经网络的时间序列极端事件预测”。
通过将历史数据拆分为输入和输出变量的滑动窗口来创建训练数据集。
本文未指定实验中使用的回顾和预测范围的具体大小。
时间序列建模的滑动窗口方法 取自“优步神经网络的时间序列极端事件预测”。
通过对每批样品的观察值进行标准化来缩放时间序列数据，并且每个输入序列被去除趋势，但是没有去季节化。
神经网络对未缩放的数据很敏感，因此我们将每个小批量标准化。此外，我们发现，与去调味相反，降低数据的趋势可以产生更好的结果。
- 优步神经网络的时间序列极端事件预测，2017。
模型
LSTM，例如 Vanilla LSTMs 在问题上进行了评估并表现出相对较差的表现。
这并不奇怪，因为它反映了其他地方的发现。
我们最初的 LSTM 实现相对于最先进的方法没有表现出优越的表现。
- 优步神经网络的时间序列极端事件预测，2017。
使用了更精细的架构，包括两个 LSTM 模型：
- 特征提取器：用于将输入序列提取到特征向量的模型，该特征向量可以用作做出预测的输入。
- Forecaster ：使用提取的特征和其他输入做出预测的模型。
开发了 LSTM 自编码器模型用作特征提取模型，并使用 Stacked LSTM 作为预测模型。
我们发现香草 LSTM 模型的表现比我们的基线差。因此，我们提出了一种新架构，它利用自编码器进行特征提取，与基线相比实现了卓越的表现。
- 优步神经网络的时间序列极端事件预测，2017。
在做出预测时，首先将时间序列数据提供给自编码器，自编码器被压缩为平均和连接的多个特征向量。然后将特征向量作为输入提供给预测模型以做出预测。
...该模型首先通过自动特征提取对网络进行填充，这对于在大规模特殊事件期间捕获复杂的时间序列动态至关重要。 [...]然后通过集合技术（例如，平均或其他方法）聚合特征向量。然后将最终向量与新输入连接并馈送到 LSTM 预测器做出预测。
- 优步神经网络的时间序列极端事件预测，2017。
目前尚不清楚在做出预测时究竟是什么给自编码器提供了什么，尽管我们可能猜测这是一个多变量时间序列，用于预测在预测时间间隔之前观察的城市。
作为自编码器输入的多变量时间序列将导致可以连接的多个编码向量（每个系列一个）。目前尚不清楚平均在这一点上可能采取什么角色，尽管我们可能猜测它是执行自编码过程的多个模型的平均值。
特征提取模型和预测模型概述 取自“优步神经网络的时间序列极端事件预测”。
作者评论说，可以将自编码器作为预测模型的一部分，并对此进行评估，但单独的模型可以提高表现。
但是，拥有一个单独的自编码器模块可以在我们的经验中产生更好的结果。
- 优步神经网络的时间序列极端事件预测，2017。
在展示纸张时使用的幻灯片中提供了所开发模型的更多细节。
自编码器的输入是 512 LSTM 单位，自编码器中的瓶颈用于创建 32 或 64 LSTM 单位的编码特征向量。
用于特征提取的 LSTM 自编码器的详细信息 取自“优步神经网络的时间序列极端事件预测”。
使用'_ 新输入 _'将编码的特征向量提供给预测模型，尽管未指定此新输入是什么;我们可以猜测这是一个时间序列，也许是预测区间之前的观测预测的城市的多变量时间序列。或者，从这个系列中提取的特征论文中的博客文章暗示（尽管我对这篇文章和幻灯片与此相矛盾时持怀疑态度）。
该模型接受了大量数据的训练，这是栈式 LSTM 或一般 LSTM 的一般要求。
所描述的生产神经网络模型在数千个时间序列上进行训练，每个时间序列具有数千个数据点。
- 优步神经网络的时间序列极端事件预测，2017。
在进行新的预测时，不会对该模型进行再训练。
还使用引导程序实现了估算预测不确定性的有趣方法。
它分别使用自编码器和预测模型分别估计模型不确定性和预测不确定性。输入被提供给给定模型并且使用了激活的丢失（如幻灯片中所评论的）。该过程重复 100 次，模型和预测误差项用于预测不确定性的估计。
预测不确定性估计概述 取自“优步神经网络的时间序列极端事件预测”。
这种预测不确定性的方法可能更好地描述于 2017 年论文“优步时间序列的深度和自信预测”。
发现
对该模型进行了评估，特别关注美国城市对美国假期的需求预测。
没有具体说明模型评估的具体情况。
新的广义 LSTM 预测模型被发现优于优步使用的现有模型，如果我们假设现有模型得到了很好的调整，这可能会令人印象深刻。
结果显示，与包含单变量时间序列和机器学习模型的当前专有方法相比，预测准确率提高了 2％-18％。
- 优步神经网络的时间序列极端事件预测，2017。
然后将在 Uber 数据集上训练的模型直接应用于由[约 1,500 个月的单变量时间序列预测数据集]组成的 M3-竞赛数据集的子集。
这是一种转移学习，一种非常理想的目标，允许跨问题域重用深度学习模型。
令人惊讶的是，该模型表现良好，与表现最佳的方法相比并不是很好，但比许多复杂模型更好。结果表明，可能通过微调（例如在其他转移学习案例研究中完成），该模型可以重复使用并且技巧娴熟。
LSTM 模型在优步数据上的表现和对 M3 数据集的评估 取自“优步神经网络的时间序列极端事件预测”。
重要的是，作者提出，深度 LSTM 模型对时间序列预测的最有益应用可能是：
- 有大量的时间序列。
- 每个系列都有大量的观察结果。
- 时间序列之间存在很强的相关性。
根据我们的经验，选择时间序列的神经网络模型有三个标准：（a）时间序列的数量（b）时间序列的长度和（c）时间序列之间的相关性。如果（a），（b）和（c）高，则神经网络可能是正确的选择，否则经典的时间序列方法可能效果最好。
- 优步神经网络的时间序列极端事件预测，2017。
通过本文介绍中使用的幻灯片很好地总结了这一点。
应用 LSTM 进行时间序列预测的经验教训 取自“优步神经网络的时间序列极端事件预测”幻灯片。
进一步阅读
如果您希望深入了解，本节将提供有关该主题的更多资源。
- 优步神经网络的时间序列极端事件预测，2017。
- 优步工程极端事件预测与循环神经网络，2017 年。
- 优步神经网络的时间序列建模，Slides，2017。
- 时间序列极端事件预测案例研究，幻灯片 2018。
- 时间序列研讨会，ICML 2017
- 优步时间序列的深度和自信预测，2017。
摘要
在这篇文章中，您发现了一个可扩展的端到端 LSTM 模型，用于时间序列预测。
具体来说，你学到了：
- 跨多个站点的多变量，多步骤预测的挑战，在这种情况下是城市。
- 用于时间序列预测的 LSTM 模型架构，包括单独的自编码器和预测子模型。
- 所提出的 LSTM 架构在罕见事件中的技能需求预测以及在不相关的预测问题上重用训练模型的能力。
你有任何问题吗？ 在下面的评论中提出您的问题，我会尽力回答。
用于时间序列预测的 4 种通用机器学习数据变换
原文：
machinelearningmastery.com/machine-learning-data-transforms-for-time-series-forecasting/
时间序列数据通常需要在使用机器学习算法建模之前进行一些准备。
例如，差分运算可用于从序列中去除趋势和季节结构，以简化预测问题。一些算法（例如神经网络）更喜欢在建模之前对数据进行标准化和/或标准化。
应用于该系列的任何变换操作也需要在预测上应用类似的逆变换。这是必需的，以便得到的计算表现度量与输出变量的比例相同，并且可以与经典预测方法进行比较。
在这篇文章中，您将了解如何在机器学习中对时间序列数据执行和反转四种常见数据转换。
阅读这篇文章后，你会知道：
- 如何在 Python 中转换和反转四种方法的变换。
- 在训练和测试数据集上使用变换时的重要注意事项。
- 在数据集上需要多个操作时建议的转换顺序。
让我们开始吧。
用于时间序列预测的 4 种通用机器学习数据变换 照片由 Wolfgang Staudt 拍摄，保留一些权利。
概观
本教程分为三个部分;他们是：
- 时间序列数据的变换
- 模型评估的考虑因素
- 数据转换顺序
时间序列数据的变换
给定单变量时间序列数据集，在使用机器学习方法进行建模和预测时，有四种变换很流行。
他们是：
- 电力转换
- 差分变换
- 标准化
- 正常化
让我们依次快速浏览一下以及如何在 Python 中执行这些转换。
我们还将审查如何反转变换操作，因为当我们想要以原始比例评估预测时，这是必需的，以便可以直接比较表现度量。
您是否希望在时间序列数据上使用其他变换来进行机器学习方法的建模？ 请在下面的评论中告诉我。
电力转换
功率变换从数据分布中移除偏移以使分布更正常（高斯分布）。
在时间序列数据集上，这可以消除随时间变化的方差。
流行的例子是对数变换（正值）或广义版本，例如 Box-Cox 变换（正值）或 Yeo-Johnson 变换（正值和负值）。
例如，我们可以使用 SciPy 库中的 boxcox（）函数在 Python 中实现 Box-Cox 变换。
默认情况下，该方法将以数字方式优化变换的 lambda 值并返回最佳值。
from scipy.stats import boxcox
# define data
data = ...
# box-cox transform
result, lmbda = boxcox(data)
变换可以反转，但需要一个名为invert_boxcox()的下面列出的自定义函数，它接受一个变换值和用于执行变换的 lambda 值。
from math import log
from math import exp
# invert a boxcox transform for one value
def invert_boxcox(value, lam):
# log case
if lam == 0:
return exp(value)
# all other cases
return exp(log(lam * value + 1) / lam)
下面列出了将功率变换应用于数据集并反转变换的完整示例。
# example of power transform and inversion
from math import log
from math import exp
from scipy.stats import boxcox
# invert a boxcox transform for one value
def invert_boxcox(value, lam):
# log case
if lam == 0:
return exp(value)
# all other cases
return exp(log(lam * value + 1) / lam)
# define dataset
data = [x for x in range(1, 10)]
print(data)
# power transform
transformed, lmbda = boxcox(data)
print(transformed, lmbda)
# invert transform
inverted = [invert_boxcox(x, lmbda) for x in transformed]
print(inverted)
运行该示例将在转换变换后打印原始数据集，幂变换的结果以及原始值（或接近它）。
[1, 2, 3, 4, 5, 6, 7, 8, 9]
[0\. 0.89887536 1.67448353 2.37952145 3.03633818 3.65711928
4.2494518 4.81847233 5.36786648] 0.7200338588580095
[1.0, 2.0, 2.9999999999999996, 3.999999999999999, 5.000000000000001, 6.000000000000001, 6.999999999999999, 7.999999999999998, 8.999999999999998]
差分变换
差分变换是从时间序列中去除系统结构的简单方法。
例如，可以通过从系列中的每个值中减去先前的值来消除趋势。这称为一阶差分。可以重复该过程（例如差异系列）以消除二阶趋势，等等。
通过从前一季节中减去观察值，可以以类似的方式去除季节性结构。 12 个步骤之前的月度数据与年度季节性结构。
可以使用下面列出的名为difference()的自定义函数计算系列中的单个差异值。该函数采用时间序列和差值计算的间隔，例如， 1 表示趋势差异，12 表示季节性差异。
# difference dataset
def difference(data, interval):
return [data[i] - data[i - interval] for i in range(interval, len(data))]
同样，可以使用自定义函数反转此操作，该函数将原始值添加回名为invert_difference()的差值，该值采用原始序列和间隔。
# invert difference
def invert_difference(orig_data, diff_data, interval):
return [diff_data[i-interval] + orig_data[i-interval] for i in range(interval, len(orig_data))]
我们可以在下面演示这个功能。
# example of a difference transform
# difference dataset
def difference(data, interval):
return [data[i] - data[i - interval] for i in range(interval, len(data))]
# invert difference
def invert_difference(orig_data, diff_data, interval):
return [diff_data[i-interval] + orig_data[i-interval] for i in range(interval, len(orig_data))]
# define dataset
data = [x for x in range(1, 10)]
print(data)
# difference transform
transformed = difference(data, 1)
print(transformed)
# invert difference
inverted = invert_difference(data, transformed, 1)
print(inverted)
运行该示例将打印原始数据集，差分变换的结果以及转换后的原始值。
注意，变换后序列中的第一个“间隔”值将丢失。这是因为它们在“间隔”之前的时间步长没有值，因此无法区分。
[1, 2, 3, 4, 5, 6, 7, 8, 9]
[1, 1, 1, 1, 1, 1, 1, 1]
[2, 3, 4, 5, 6, 7, 8, 9]
标准化
标准化是具有高斯分布的数据的变换。
它减去均值并将结果除以数据样本的标准差。这具有将数据转换为具有零或中心的均值的效果，其标准偏差为 1.这样得到的分布称为标准高斯分布，或标准法线，因此称为变换的名称。
我们可以使用 scikit-learn 库中的 Python 中的 StandardScaler 对象执行标准化。
此类允许通过调用fit()将变换拟合到训练数据集上，通过调用transform()应用于一个或多个数据集（例如训练和测试）并且还提供通过调用inverse_transform()来反转变换的函数。
下面应用完整的示例。
# example of standardization
from sklearn.preprocessing import StandardScaler
from numpy import array
# define dataset
data = [x for x in range(1, 10)]
data = array(data).reshape(len(data), 1)
print(data)
# fit transform
transformer = StandardScaler()
transformer.fit(data)
# difference transform
transformed = transformer.transform(data)
print(transformed)
# invert difference
inverted = transformer.inverse_transform(transformed)
print(inverted)
运行该示例将打印原始数据集，标准化变换的结果以及转换后的原始值。
请注意，期望数据作为具有多行的列提供。
[[1]
[2]
[3]
[4]
[5]
[6]
[7]
[8]
[9]]
[[-1.54919334]
[-1.161895
[-0.77459667]
[-0.38729833]
[ 0\.
[ 0.38729833]
[ 0.77459667]
[ 1.161895
[ 1.54919334]]
[[1.]
[2.]
[3.]
[4.]
[5.]
[6.]
[7.]
[8.]
[9.]]
正常化
规范化是将数据从原始范围重新缩放到 0 到 1 之间的新范围。
与标准化一样，这可以使用 scikit-learn 库中的转换对象来实现，特别是 MinMaxScaler 类。除了规范化之外，通过在对象的构造函数中指定首选范围，此类可用于将数据重新缩放到您希望的任何范围。
它可以以相同的方式用于拟合，变换和反转变换。
下面列出了一个完整的例子。
# example of normalization
from sklearn.preprocessing import MinMaxScaler
from numpy import array
# define dataset
data = [x for x in range(1, 10)]
data = array(data).reshape(len(data), 1)
print(data)
# fit transform
transformer = MinMaxScaler()
transformer.fit(data)
# difference transform
transformed = transformer.transform(data)
print(transformed)
# invert difference
inverted = transformer.inverse_transform(transformed)
print(inverted)
运行该示例将打印原始数据集，规范化转换的结果以及转换后的原始值。
[[1]
[2]
[3]
[4]
[5]
[6]
[7]
[8]
[9]]
[[0\.
[0.125]
[0.25 ]
[0.375]
[0.5
[0.625]
[0.75 ]
[0.875]
[1\. ]
[[1.]
[2.]
[3.]
[4.]
[5.]
[6.]
[7.]
[8.]
[9.]]
模型评估的考虑因素
我们已经提到了能够反转模型预测变换的重要性，以便计算与其他方法直接相当的模型表现统计量。
另外，另一个问题是数据泄漏问题。
上述三个数据转换来自提供的数据集的估计系数，然后用于转换数据。特别：
- Power Transform ：lambda 参数。
- 标准化：平均值和标准差统计量。
- 标准化：最小值和最大值。
必须仅在训练数据集上估计这些系数。
估计完成后，可以在评估模型之前使用系数对训练和测试数据集应用变换。
如果在分割成训练集和测试集之前使用整个数据集估计系数，则从测试集到训练数据集的信息泄漏很小。这可能导致对乐观偏见的模型技能的估计。
因此，您可能希望使用领域知识增强系数的估计值，例如将来所有时间的预期最小值/最大值。
通常，差分不会遇到相同的问题。在大多数情况下，例如一步预测，可以使用滞后观察来执行差异计算。如果不是，则可以在任何需要的地方使用滞后预测作为差异计算中真实观察的代理。
数据转换顺序
您可能希望尝试在建模之前将多个数据转换应用于时间序列。
这很常见，例如应用幂变换以消除增加的方差，应用季节差异来消除季节性，并应用一步差分来移除趋势。
应用转换操作的顺序很重要。
直觉上，我们可以思考变换如何相互作用。
- 应该在差分之前执行功率变换。
- 应在一步差分之前进行季节性差异。
- 标准化是线性的，应在任何非线性变换和差分后对样本进行标准化。
- 标准化是线性操作，但它应该是为保持首选标度而执行的最终变换。
因此，建议的数据转换顺序如下：
- 电力转换。
- 季节性差异。
- 趋势差异。
- 标准化。
- 正常化。
显然，您只能使用特定数据集所需的变换。
重要的是，当变换操作被反转时，必须反转逆变换操作的顺序。具体而言，必须按以下顺序执行逆操作：
- 正常化。
- 标准化。
- 趋势差异。
- 季节性差异。
- 电力转换。
进一步阅读
如果您希望深入了解，本节将提供有关该主题的更多资源。
帖子
- 如何使用 Python 进行时间序列预测数据的电源转换
- 如何使用 Python 中的差分变换删除趋势和季节性
- 如何区分时间序列数据集与 Python
- 如何在 Python 中标准化和标准化时间序列数据
蜜蜂
- scipy.stats.boxcox API
- sklearn.preprocessing.MinMaxScaler API
- sklearn.preprocessing.StandardScaler API
用品
摘要
在这篇文章中，您了解了如何在机器学习中对时间序列数据执行和反转四种常见数据转换。
具体来说，你学到了：
- 如何在 Python 中转换和反转四种方法的变换。
- 在训练和测试数据集上使用变换时的重要注意事项。
- 在数据集上需要多个操作时建议的转换顺序。
你有任何问题吗？ 在下面的评论中提出您的问题，我会尽力回答。
Python 中使用长短期记忆网络的多步时间序列预测
长期短期记忆网络或 LSTM 是一种可以学习和预测长序列的循环神经网络。
除了学习长序列之外，LSTM 的一个好处是它们可以学习进行一次性多步预测，这对于时间序列预测可能是有用的。
LSTM 的一个难点是它们配置起来很棘手，需要大量准备才能以正确的格式获取数据进行学习。
在本教程中，您将了解如何使用 Keras 在 Python 中开发用于多步骤时间序列预测的 LSTM。
完成本教程后，您将了解：
- 如何为多步时间序列预测准备数据。
- 如何开发 LSTM 模型进行多步时间序列预测。
- 如何评估多步时间序列预测。
让我们开始吧。
Python 中长期短期记忆网络的多步时间序列预测 Tom Babich 的照片，保留一些权利。
教程概述
本教程分为 4 个部分;他们是：
- 洗发水销售数据集
- 数据准备和模型评估
- 持久性模型
- 多步骤 LSTM
环境
本教程假定您已安装 Python SciPy 环境。您可以在此示例中使用 Python 2 或 3。
本教程假设您安装了 TensorFlow 或 Theano 后端的 Keras v2.0 或更高版本。
本教程还假设您安装了 scikit-learn，Pandas，NumPy 和 Matplotlib。
如果您在设置 Python 环境时需要帮助，请参阅以下帖子：
接下来，让我们看看标准时间序列预测问题，我们可以将其用作此实验的上下文。
洗发水销售数据集
该数据集描述了 3 年期间每月洗发水的销售数量。
单位是销售计数，有 36 个观察。原始数据集归功于 Makridakis，Wheelwright 和 Hyndman（1998）。
下面的示例加载并创建已加载数据集的图。
# load and plot dataset
from pandas import read_csv
from pandas import datetime
from matplotlib import pyplot
# load dataset
def parser(x):
return datetime.strptime('190'+x, '%Y-%m')
series = read_csv('shampoo-sales.csv', header=0, parse_dates=[0], index_col=0, squeeze=True, date_parser=parser)
# summarize first few rows
print(series.head())
# line plot
series.plot()
pyplot.show()
运行该示例将数据集作为 Pandas Series 加载并打印前 5 行。
Month
1901-01-01 266.0
1901-02-01 145.9
1901-03-01 183.1
1901-04-01 119.3
1901-05-01 180.3
Name: Sales, dtype: float64
然后创建该系列的线图，显示明显的增加趋势。
洗发水销售数据集的线图
接下来，我们将看一下实验中使用的模型配置和测试工具。
数据准备和模型评估
本节介绍本教程中使用的数据准备和模型评估
数据拆分
我们将 Shampoo Sales 数据集分为两部分：训练和测试集。
前两年的数据将用于训练数据集，剩余的一年数据将用于测试集。
将使用训练数据集开发模型，并对测试数据集做出预测。
作为参考，过去 12 个月的观察如下：
"3-01",339.7
"3-02",440.4
"3-03",315.9
"3-04",439.3
"3-05",401.3
"3-06",437.4
"3-07",575.5
"3-08",407.6
"3-09",682.0
"3-10",475.3
"3-11",581.3
"3-12",646.9
多步预测
我们将设计一个多步预测。
对于数据集最后 12 个月的特定月份，我们将需要进行 3 个月的预测。
这是历史观察（t-1，t-2，... t-n）预测 t，t + 1 和 t + 2。
具体而言，从第 2 年 12月开始，我们必须预测 1 月，2 月和 3 月。从 1 月开始，我们必须预测 2 月，3 月和 4 月。一直到 10 月，11 月，12 月预测从 9 月到 3 年。
需要总共 10 个 3 个月的预测，如下：
Dec, Jan, Feb, Mar
Jan, Feb, Mar, Apr
Feb, Mar, Apr, May
Mar, Apr, May, Jun
Apr, May, Jun, Jul
May, Jun, Jul, Aug
Jun, Jul, Aug, Sep
Jul, Aug, Sep, Oct
Aug, Sep, Oct, Nov
Sep, Oct, Nov, Dec
模型评估
将使用滚动预测场景，也称为前进模型验证。
测试数据集的每个时间步骤将一次一个地走。将使用模型对时间步长做出预测，然后将从测试集中获取下个月的实际预期值，并使其可用于下一时间步的预测模型。
这模仿了一个真实世界的场景，每个月都会有新的洗发水销售观察结果，并用于下个月的预测。
这将通过训练和测试数据集的结构进行模拟。
将收集关于测试数据集的所有预测并计算错误分数以总结每个预测时间步骤的模型技能。将使用均方根误差（RMSE），因为它会对大错误进行处罚，并产生与预测数据相同的分数，即每月洗发水销售额。
持久性模型
时间序列预测的良好基线是持久性模型。
这是一个预测模型，其中最后一个观察结果是持续的。由于它的简单性，它通常被称为朴素的预测。
您可以在帖子中了解有关时间序列预测的持久性模型的更多信息：
准备数据
第一步是将数据从一系列转换为监督学习问题。
这是从数字列表到输入和输出模式列表。我们可以使用一个名为series_to_supervised()的预先准备的函数来实现这一点。
有关此功能的更多信息，请参阅帖子：
该功能如下所列。
# convert time series into supervised learning problem
def series_to_supervised(data, n_in=1, n_out=1, dropnan=True):
n_vars = 1 if type(data) is list else data.shape[1]
df = DataFrame(data)
cols, names = list(), list()
# input sequence (t-n, ... t-1)
for i in range(n_in, 0, -1):
cols.append(df.shift(i))
names += [('var%d(t-%d)' % (j+1, i)) for j in range(n_vars)]
# forecast sequence (t, t+1, ... t+n)
for i in range(0, n_out):
cols.append(df.shift(-i))
if i == 0:
names += [('var%d(t)' % (j+1)) for j in range(n_vars)]
else:
names += [('var%d(t+%d)' % (j+1, i)) for j in range(n_vars)]
# put it all together
agg = concat(cols, axis=1)
agg.columns = names
# drop rows with NaN values
if dropnan:
agg.dropna(inplace=True)
return agg
可以通过将加载的系列值传入n_in值 1 和 n_out 值 3 来调用该函数;例如：
supervised = series_to_supervised(raw_values, 1, 3)
接下来，我们可以将监督学习数据集分成训练和测试集。
我们知道，在这种形式中，最后 10 行包含最后一年的数据。这些行包含测试集，其余数据构成训练数据集。
我们可以将所有这些放在一个新函数中，该函数接受加载的系列和一些参数，并返回准备建模的训练和测试集。
# transform series into train and test sets for supervised learning
def prepare_data(series, n_test, n_lag, n_seq):
# extract raw values
raw_values = series.values
raw_values = raw_values.reshape(len(raw_values), 1)
# transform into supervised learning problem X, y
supervised = series_to_supervised(raw_values, n_lag, n_seq)
supervised_values = supervised.values
# split into train and test sets
train, test = supervised_values[0:-n_test], supervised_values[-n_test:]
return train, test
我们可以使用 Shampoo 数据集对此进行测试。下面列出了完整的示例。
from pandas import DataFrame
from pandas import concat
from pandas import read_csv
from pandas import datetime
# date-time parsing function for loading the dataset
def parser(x):
return datetime.strptime('190'+x, '%Y-%m')
# convert time series into supervised learning problem
def series_to_supervised(data, n_in=1, n_out=1, dropnan=True):
n_vars = 1 if type(data) is list else data.shape[1]
df = DataFrame(data)
cols, names = list(), list()
# input sequence (t-n, ... t-1)
for i in range(n_in, 0, -1):
cols.append(df.shift(i))
names += [('var%d(t-%d)' % (j+1, i)) for j in range(n_vars)]
# forecast sequence (t, t+1, ... t+n)
for i in range(0, n_out):
cols.append(df.shift(-i))
if i == 0:
names += [('var%d(t)' % (j+1)) for j in range(n_vars)]
else:
names += [('var%d(t+%d)' % (j+1, i)) for j in range(n_vars)]
# put it all together
agg = concat(cols, axis=1)
agg.columns = names
# drop rows with NaN values
if dropnan:
agg.dropna(inplace=True)
return agg
# transform series into train and test sets for supervised learning
def prepare_data(series, n_test, n_lag, n_seq):
# extract raw values
raw_values = series.values
raw_values = raw_values.reshape(len(raw_values), 1)
# transform into supervised learning problem X, y
supervised = series_to_supervised(raw_values, n_lag, n_seq)
supervised_values = supervised.values
# split into train and test sets
train, test = supervised_values[0:-n_test], supervised_values[-n_test:]
return train, test
# load dataset
series = read_csv('shampoo-sales.csv', header=0, parse_dates=[0], index_col=0, squeeze=True, date_parser=parser)
# configure
n_lag = 1
n_seq = 3
n_test = 10
# prepare data
train, test = prepare_data(series, n_test, n_lag, n_seq)
print(test)
print('Train: %s, Test: %s' % (train.shape, test.shape))
首先运行该示例将打印整个测试数据集，即最后 10 行。还打印了训练测试数据集的形状和大小。
[[ 342.3 339.7 440.4 315.9]
[ 339.7 440.4 315.9 439.3]
[ 440.4 315.9 439.3 401.3]
[ 315.9 439.3 401.3 437.4]
[ 439.3 401.3 437.4 575.5]
[ 401.3 437.4 575.5 407.6]
[ 437.4 575.5 407.6 682\. ]
[ 575.5 407.6 682\. 475.3]
[ 407.6 682\. 475.3 581.3]
[ 682\. 475.3 581.3 646.9]]
Train: (23, 4), Test: (10, 4)
我们可以看到测试数据集第一行的单个输入值（第一列）与第二年 12 月的洗发水销售中的观察结果相符：
"2-12",342.3
我们还可以看到每行包含 4 列用于 1 个输入，3 个输出值用于每个观察。
做出预测
下一步是进行持久性预测。
我们可以在名为persistence()的函数中轻松实现持久性预测，该函数将最后一次观察和预测步骤的数量保持不变。此函数返回包含预测的数组。
# make a persistence forecast
def persistence(last_ob, n_seq):
return [last_ob for i in range(n_seq)]
然后，我们可以在测试数据集中的每个时间步骤调用此函数，从第 2 年的 12 月到第 3 年的 9 月。
下面是一个函数 make_forecasts（），它执行此操作并将数据集的训练，测试和配置作为参数，并返回预测列表。
# evaluate the persistence model
def make_forecasts(train, test, n_lag, n_seq):
forecasts = list()
for i in range(len(test)):
X, y = test[i, 0:n_lag], test[i, n_lag:]
# make forecast
forecast = persistence(X[-1], n_seq)
# store the forecast
forecasts.append(forecast)
return forecasts
我们可以调用这个函数如下：
forecasts = make_forecasts(train, test, 1, 3)
评估预测
最后一步是评估预测。
我们可以通过计算多步预测的每个时间步长的 RMSE 来做到这一点，在这种情况下给出 3 个 RMSE 分数。下面的函数evaluate_forecasts()计算并打印每个预测时间步的 RMSE。
# evaluate the RMSE for each forecast time step
def evaluate_forecasts(test, forecasts, n_lag, n_seq):
for i in range(n_seq):
actual = test[:,(n_lag+i)]
predicted = [forecast[i] for forecast in forecasts]
rmse = sqrt(mean_squared_error(actual, predicted))
print('t+%d RMSE: %f' % ((i+1), rmse))
我们可以这样称呼它：
evaluate_forecasts(test, forecasts, 1, 3)
在原始数据集的上下文中绘制预测图也有助于了解 RMSE 分数如何与上下文中的问题相关联。
我们可以首先绘制整个 Shampoo 数据集，然后将每个预测绘制为红线。下面的函数plot_forecasts()将创建并显示此图。
# plot the forecasts in the context of the original dataset
def plot_forecasts(series, forecasts, n_test):
# plot the entire dataset in blue
pyplot.plot(series.values)
# plot the forecasts in red
for i in range(len(forecasts)):
off_s = len(series) - n_test + i
off_e = off_s + len(forecasts[i])
xaxis = [x for x in range(off_s, off_e)]
pyplot.plot(xaxis, forecasts[i], color='red')
# show the plot
pyplot.show()
我们可以按如下方式调用该函数。请注意，12 个月内在测试集上保留的观察数为 12，而上述使用的 10 个监督学习输入/输出模式则为 10。
# plot forecasts
plot_forecasts(series, forecasts, 12)
我们可以通过将持久预测与原始数据集中的实际持久值相关联来使绘图更好。
这将需要将最后观察到的值添加到预测的前面。以下是具有此改进的plot_forecasts()功能的更新版本。
# plot the forecasts in the context of the original dataset
def plot_forecasts(series, forecasts, n_test):
# plot the entire dataset in blue
pyplot.plot(series.values)
# plot the forecasts in red
for i in range(len(forecasts)):
off_s = len(series) - 12 + i - 1
off_e = off_s + len(forecasts[i]) + 1
xaxis = [x for x in range(off_s, off_e)]
yaxis = [series.values[off_s]] + forecasts[i]
pyplot.plot(xaxis, yaxis, color='red')
# show the plot
pyplot.show()
完整的例子
我们可以将所有这些部分组合在一起。
下面列出了多步持久性预测的完整代码示例。
from pandas import DataFrame
from pandas import concat
from pandas import read_csv
from pandas import datetime
from sklearn.metrics import mean_squared_error
from math import sqrt
from matplotlib import pyplot
# date-time parsing function for loading the dataset
def parser(x):
return datetime.strptime('190'+x, '%Y-%m')
# convert time series into supervised learning problem
def series_to_supervised(data, n_in=1, n_out=1, dropnan=True):
n_vars = 1 if type(data) is list else data.shape[1]
df = DataFrame(data)
cols, names = list(), list()
# input sequence (t-n, ... t-1)
for i in range(n_in, 0, -1):
cols.append(df.shift(i))
names += [('var%d(t-%d)' % (j+1, i)) for j in range(n_vars)]
# forecast sequence (t, t+1, ... t+n)
for i in range(0, n_out):
cols.append(df.shift(-i))
if i == 0:
names += [('var%d(t)' % (j+1)) for j in range(n_vars)]
else:
names += [('var%d(t+%d)' % (j+1, i)) for j in range(n_vars)]
# put it all together
agg = concat(cols, axis=1)
agg.columns = names
# drop rows with NaN values
if dropnan:
agg.dropna(inplace=True)
return agg
# transform series into train and test sets for supervised learning
def prepare_data(series, n_test, n_lag, n_seq):
# extract raw values
raw_values = series.values
raw_values = raw_values.reshape(len(raw_values), 1)
# transform into supervised learning problem X, y
supervised = series_to_supervised(raw_values, n_lag, n_seq)
supervised_values = supervised.values
# split into train and test sets
train, test = supervised_values[0:-n_test], supervised_values[-n_test:]
return train, test
# make a persistence forecast
def persistence(last_ob, n_seq):
return [last_ob for i in range(n_seq)]
# evaluate the persistence model
def make_forecasts(train, test, n_lag, n_seq):
forecasts = list()
for i in range(len(test)):
X, y = test[i, 0:n_lag], test[i, n_lag:]
# make forecast
forecast = persistence(X[-1], n_seq)
# store the forecast
forecasts.append(forecast)
return forecasts
# evaluate the RMSE for each forecast time step
def evaluate_forecasts(test, forecasts, n_lag, n_seq):
for i in range(n_seq):
actual = test[:,(n_lag+i)]
predicted = [forecast[i] for forecast in forecasts]
rmse = sqrt(mean_squared_error(actual, predicted))
print('t+%d RMSE: %f' % ((i+1), rmse))
# plot the forecasts in the context of the original dataset
def plot_forecasts(series, forecasts, n_test):
# plot the entire dataset in blue
pyplot.plot(series.values)
# plot the forecasts in red
for i in range(len(forecasts)):
off_s = len(series) - n_test + i - 1
off_e = off_s + len(forecasts[i]) + 1
xaxis = [x for x in range(off_s, off_e)]
yaxis = [series.values[off_s]] + forecasts[i]
pyplot.plot(xaxis, yaxis, color='red')
# show the plot
pyplot.show()
# load dataset
series = read_csv('shampoo-sales.csv', header=0, parse_dates=[0], index_col=0, squeeze=True, date_parser=parser)
# configure
n_lag = 1
n_seq = 3
n_test = 10
# prepare data
train, test = prepare_data(series, n_test, n_lag, n_seq)
# make forecasts
forecasts = make_forecasts(train, test, n_lag, n_seq)
# evaluate forecasts
evaluate_forecasts(test, forecasts, n_lag, n_seq)
# plot forecasts
plot_forecasts(series, forecasts, n_test+2)
首先运行该示例为每个预测的时间步骤打印 RMSE。
这为我们提供了每个时间步的表现基线，我们希望 LSTM 能够表现出色。
t+1 RMSE: 144.535304
t+2 RMSE: 86.479905
t+3 RMSE: 121.149168
还创建了原始时间序列与多步持久性预测的关系图。这些行连接到每个预测的相应输入值。
此上下文显示了持久性预测实际上是多么幼稚。
具有多步持续性预测的洗发水销售数据集线图
多步 LSTM 网络
在本节中，我们将使用持久性示例作为起点，并查看将 LSTM 拟合到训练数据所需的更改，并对测试数据集进行多步预测。
准备数据
在我们使用它来训练 LSTM 之前，必须准备好数据。
具体而言，还需要进行两项更改：
- 固定。数据显示必须通过差分消除的增加趋势。
- 比例。必须将数据的比例缩小到-1 到 1 之间的值，即 LSTM 单元的激活功能。
我们可以引入一个函数使数据静止称为 difference（）。这会将一系列值转换为一系列差异，这是一种更简单的表示方式。
# create a differenced series
def difference(dataset, interval=1):
diff = list()
for i in range(interval, len(dataset)):
value = dataset[i] - dataset[i - interval]
diff.append(value)
return Series(diff)
我们可以使用 sklearn 库中的MinMaxScaler来缩放数据。
将这些放在一起，我们可以更新prepare_data()函数以首先区分数据并重缩放它，然后执行转换为监督学习问题并训练测试集，就像我们之前使用持久性示例一样。
除了训练和测试数据集之外，该函数现在返回一个缩放器。
# transform series into train and test sets for supervised learning
def prepare_data(series, n_test, n_lag, n_seq):
# extract raw values
raw_values = series.values
# transform data to be stationary
diff_series = difference(raw_values, 1)
diff_values = diff_series.values
diff_values = diff_values.reshape(len(diff_values), 1)
# rescale values to -1, 1
scaler = MinMaxScaler(feature_range=(-1, 1))
scaled_values = scaler.fit_transform(diff_values)
scaled_values = scaled_values.reshape(len(scaled_values), 1)
# transform into supervised learning problem X, y
supervised = series_to_supervised(scaled_values, n_lag, n_seq)
supervised_values = supervised.values
# split into train and test sets
train, test = supervised_values[0:-n_test], supervised_values[-n_test:]
return scaler, train, test
我们可以调用这个函数如下：
# prepare data
scaler, train, test = prepare_data(series, n_test, n_lag, n_seq)
适合 LSTM 网络
接下来，我们需要将 LSTM 网络模型与训练数据相匹配。
这首先要求将训练数据集从 2D 数组[_ 样本，特征 ]转换为 3D 数组[ 样本，时间步长，特征 _]。我们将时间步长固定为 1，因此这种变化很简单。
接下来，我们需要设计一个 LSTM 网络。我们将使用一个简单的结构，其中 1 个隐藏层具有 1 个 LSTM 单元，然后是具有线性激活和 3 个输出值的输出层。网络将使用均方误差丢失函数和有效的 ADAM 优化算法。
LSTM 是有状态的;这意味着我们必须在每个训练时代结束时手动重置网络状态。该网络将适合 1500 个时代。
必须使用相同的批量大小进行训练和预测，并且我们需要在测试数据集的每个时间步做出预测。这意味着必须使用批量大小为 1。批量大小为 1 也称为在线学习，因为在每个训练模式之后的训练期间将更新网络权重（与小批量或批量更新相反）。
我们可以将所有这些放在一个名为fit_lstm()的函数中。该函数采用了许多可用于稍后调整网络的关键参数，并且该函数返回适合 LSTM 模型以备预测。
# fit an LSTM network to training data
def fit_lstm(train, n_lag, n_seq, n_batch, nb_epoch, n_neurons):
# reshape training into [samples, timesteps, features]
X, y = train[:, 0:n_lag], train[:, n_lag:]
X = X.reshape(X.shape[0], 1, X.shape[1])
# design network
model = Sequential()
model.add(LSTM(n_neurons, batch_input_shape=(n_batch, X.shape[1], X.shape[2]), stateful=True))
model.add(Dense(y.shape[1]))
model.compile(loss='mean_squared_error', optimizer='adam')
# fit network
for i in range(nb_epoch):
model.fit(X, y, epochs=1, batch_size=n_batch, verbose=0, shuffle=False)
model.reset_states()
return model
该函数可以调用如下：
# fit model
model = fit_lstm(train, 1, 3, 1, 1500, 1)
网络配置没有调整;如果你愿意，尝试不同的参数。
在下面的评论中报告您的发现。我很想看看你能得到什么。
进行 LSTM 预测
下一步是使用适合的 LSTM 网络做出预测。
通过调用 model.predict（），可以使用拟合 LSTM 网络进行单个预测。同样，必须将数据格式化为具有[_ 样本，时间步长，特征 _]格式的 3D 数组。
我们可以将它包装成一个名为forecast_lstm()的函数。
# make one forecast with an LSTM,
def forecast_lstm(model, X, n_batch):
# reshape input pattern to [samples, timesteps, features]
X = X.reshape(1, 1, len(X))
# make forecast
forecast = model.predict(X, batch_size=n_batch)
# convert to array
return [x for x in forecast[0, :]]
我们可以从make_forecasts()函数调用此函数并更新它以接受模型作为参数。更新版本如下所示。
# evaluate the persistence model
def make_forecasts(model, n_batch, train, test, n_lag, n_seq):
forecasts = list()
for i in range(len(test)):
X, y = test[i, 0:n_lag], test[i, n_lag:]
# make forecast
forecast = forecast_lstm(model, X, n_batch)
# store the forecast
forecasts.append(forecast)
return forecasts
可以按如下方式调用make_forecasts()函数的更新版本：
# make forecasts
forecasts = make_forecasts(model, 1, train, test, 1, 3)
反转变换
在做出预测之后，我们需要反转变换以将值返回到原始比例。
这是必要的，以便我们可以计算与其他模型相当的错误分数和图，例如上面的持久性预测。
我们可以使用提供inverse_transform()函数的MinMaxScaler对象直接反转预测的比例。
我们可以通过将最后一个观察值（前几个月的洗发水销售额）添加到第一个预测值来反转差异，然后将值传播到预测值之下。
这有点儿繁琐;我们可以在函数名inverse_difference()中包含行为，它将预测之前的最后观察值和预测作为参数，并返回反向预测。
# invert differenced forecast
def inverse_difference(last_ob, forecast):
# invert first forecast
inverted = list()
inverted.append(forecast[0] + last_ob)
# propagate difference forecast using inverted first value
for i in range(1, len(forecast)):
inverted.append(forecast[i] + inverted[i-1])
return inverted
将它们放在一起，我们可以创建一个inverse_transform()函数，该函数可以处理每个预测，首先反转比例，然后反转差异，将预测恢复到原始比例。
# inverse data transform on forecasts
def inverse_transform(series, forecasts, scaler, n_test):
inverted = list()
for i in range(len(forecasts)):
# create array from forecast
forecast = array(forecasts[i])
forecast = forecast.reshape(1, len(forecast))
# invert scaling
inv_scale = scaler.inverse_transform(forecast)
inv_scale = inv_scale[0, :]
# invert differencing
index = len(series) - n_test + i - 1
last_ob = series.values[index]
inv_diff = inverse_difference(last_ob, inv_scale)
# store
inverted.append(inv_diff)
return inverted
我们可以使用以下预测来调用此函数：
# inverse transform forecasts and test
forecasts = inverse_transform(series, forecasts, scaler, n_test+2)
我们还可以反转输出部分测试数据集上的变换，以便我们可以正确计算 RMSE 分数，如下所示：
actual = [row[n_lag:] for row in test]
actual = inverse_transform(series, actual, scaler, n_test+2)
我们还可以简化 RMSE 分数的计算，以期望测试数据仅包含输出值，如下所示：
def evaluate_forecasts(test, forecasts, n_lag, n_seq):
for i in range(n_seq):
actual = [row[i] for row in test]
predicted = [forecast[i] for forecast in forecasts]
rmse = sqrt(mean_squared_error(actual, predicted))
print('t+%d RMSE: %f' % ((i+1), rmse))
完整的例子
我们可以将所有这些部分组合在一起，使 LSTM 网络适应多步骤时间序列预测问题。
完整的代码清单如下。
from pandas import DataFrame
from pandas import Series
from pandas import concat
from pandas import read_csv
from pandas import datetime
from sklearn.metrics import mean_squared_error
from sklearn.preprocessing import MinMaxScaler
from keras.models import Sequential
from keras.layers import Dense
from keras.layers import LSTM
from math import sqrt
from matplotlib import pyplot
from numpy import array
# date-time parsing function for loading the dataset
def parser(x):
return datetime.strptime('190'+x, '%Y-%m')
# convert time series into supervised learning problem
def series_to_supervised(data, n_in=1, n_out=1, dropnan=True):
n_vars = 1 if type(data) is list else data.shape[1]
df = DataFrame(data)
cols, names = list(), list()
# input sequence (t-n, ... t-1)
for i in range(n_in, 0, -1):
cols.append(df.shift(i))
names += [('var%d(t-%d)' % (j+1, i)) for j in range(n_vars)]
# forecast sequence (t, t+1, ... t+n)
for i in range(0, n_out):
cols.append(df.shift(-i))
if i == 0:
names += [('var%d(t)' % (j+1)) for j in range(n_vars)]
else:
names += [('var%d(t+%d)' % (j+1, i)) for j in range(n_vars)]
# put it all together
agg = concat(cols, axis=1)
agg.columns = names
# drop rows with NaN values
if dropnan:
agg.dropna(inplace=True)
return agg
# create a differenced series
def difference(dataset, interval=1):
diff = list()
for i in range(interval, len(dataset)):
value = dataset[i] - dataset[i - interval]
diff.append(value)
return Series(diff)
# transform series into train and test sets for supervised learning
def prepare_data(series, n_test, n_lag, n_seq):
# extract raw values
raw_values = series.values
# transform data to be stationary
diff_series = difference(raw_values, 1)
diff_values = diff_series.values
diff_values = diff_values.reshape(len(diff_values), 1)
# rescale values to -1, 1
scaler = MinMaxScaler(feature_range=(-1, 1))
scaled_values = scaler.fit_transform(diff_values)
scaled_values = scaled_values.reshape(len(scaled_values), 1)
# transform into supervised learning problem X, y
supervised = series_to_supervised(scaled_values, n_lag, n_seq)
supervised_values = supervised.values
# split into train and test sets
train, test = supervised_values[0:-n_test], supervised_values[-n_test:]
return scaler, train, test
# fit an LSTM network to training data
def fit_lstm(train, n_lag, n_seq, n_batch, nb_epoch, n_neurons):
# reshape training into [samples, timesteps, features]
X, y = train[:, 0:n_lag], train[:, n_lag:]
X = X.reshape(X.shape[0], 1, X.shape[1])
# design network
model = Sequential()
model.add(LSTM(n_neurons, batch_input_shape=(n_batch, X.shape[1], X.shape[2]), stateful=True))
model.add(Dense(y.shape[1]))
model.compile(loss='mean_squared_error', optimizer='adam')
# fit network
for i in range(nb_epoch):
model.fit(X, y, epochs=1, batch_size=n_batch, verbose=0, shuffle=False)
model.reset_states()
return model
# make one forecast with an LSTM,
def forecast_lstm(model, X, n_batch):
# reshape input pattern to [samples, timesteps, features]
X = X.reshape(1, 1, len(X))
# make forecast
forecast = model.predict(X, batch_size=n_batch)
# convert to array
return [x for x in forecast[0, :]]
# evaluate the persistence model
def make_forecasts(model, n_batch, train, test, n_lag, n_seq):
forecasts = list()
for i in range(len(test)):
X, y = test[i, 0:n_lag], test[i, n_lag:]
# make forecast
forecast = forecast_lstm(model, X, n_batch)
# store the forecast
forecasts.append(forecast)
return forecasts
# invert differenced forecast
def inverse_difference(last_ob, forecast):
# invert first forecast
inverted = list()
inverted.append(forecast[0] + last_ob)
# propagate difference forecast using inverted first value
for i in range(1, len(forecast)):
inverted.append(forecast[i] + inverted[i-1])
return inverted
# inverse data transform on forecasts
def inverse_transform(series, forecasts, scaler, n_test):
inverted = list()
for i in range(len(forecasts)):
# create array from forecast
forecast = array(forecasts[i])
forecast = forecast.reshape(1, len(forecast))
# invert scaling
inv_scale = scaler.inverse_transform(forecast)
inv_scale = inv_scale[0, :]
# invert differencing
index = len(series) - n_test + i - 1
last_ob = series.values[index]
inv_diff = inverse_difference(last_ob, inv_scale)
# store
inverted.append(inv_diff)
return inverted
# evaluate the RMSE for each forecast time step
def evaluate_forecasts(test, forecasts, n_lag, n_seq):
for i in range(n_seq):
actual = [row[i] for row in test]
predicted = [forecast[i] for forecast in forecasts]
rmse = sqrt(mean_squared_error(actual, predicted))
print('t+%d RMSE: %f' % ((i+1), rmse))
# plot the forecasts in the context of the original dataset
def plot_forecasts(series, forecasts, n_test):
# plot the entire dataset in blue
pyplot.plot(series.values)
# plot the forecasts in red
for i in range(len(forecasts)):
off_s = len(series) - n_test + i - 1
off_e = off_s + len(forecasts[i]) + 1
xaxis = [x for x in range(off_s, off_e)]
yaxis = [series.values[off_s]] + forecasts[i]
pyplot.plot(xaxis, yaxis, color='red')
# show the plot
pyplot.show()
# load dataset
series = read_csv('shampoo-sales.csv', header=0, parse_dates=[0], index_col=0, squeeze=True, date_parser=parser)
# configure
n_lag = 1
n_seq = 3
n_test = 10
n_epochs = 1500
n_batch = 1
n_neurons = 1
# prepare data
scaler, train, test = prepare_data(series, n_test, n_lag, n_seq)
# fit model
model = fit_lstm(train, n_lag, n_seq, n_batch, n_epochs, n_neurons)
# make forecasts
forecasts = make_forecasts(model, n_batch, train, test, n_lag, n_seq)
# inverse transform forecasts and test
forecasts = inverse_transform(series, forecasts, scaler, n_test+2)
actual = [row[n_lag:] for row in test]
actual = inverse_transform(series, actual, scaler, n_test+2)
# evaluate forecasts
evaluate_forecasts(actual, forecasts, n_lag, n_seq)
# plot forecasts
plot_forecasts(series, forecasts, n_test+2)
首先运行该示例为每个预测的时间步骤打印 RMSE。
我们可以看到，每个预测时间步长的分数比持久性预测更好，在某些情况下要好得多。
这表明配置的 LSTM 确实掌握了问题的技巧。
值得注意的是，RMSE 并没有像预期的那样随着预测范围的长度而变得越来越差。这是因为 t + 2 似乎比 t + 1 更容易预测。这可能是因为向下滴答比系列中记录的向上滴答更容易预测（这可以通过对结果进行更深入的分析来确认）。
t+1 RMSE: 95.973221
t+2 RMSE: 78.872348
t+3 RMSE: 105.613951
还创建了系列（蓝色）和预测（红色）的线图。
该图显示，尽管模型的技能更好，但有些预测并不是很好，而且还有很大的改进空间。
洗发水销售数据集的线图与多步骤 LSTM 预测
扩展
如果您希望超越本教程，可以考虑一些扩展。
- 更新 LSTM 。更改示例以在新数据可用时重新安装或更新 LSTM。一个 10 秒的训练时期应该足以重新训练一个新的观察。
- 调整 LSTM 。网格搜索教程中使用的一些 LSTM 参数，例如时期数，神经元数和层数，以查看是否可以进一步提升表现。
- Seq2Seq 。使用 LSTM 的编解码器范例来预测每个序列，看看它是否提供任何好处。
- 时间范围。尝试预测不同的时间范围，并了解网络的行为在不同的交付周期中如何变化。
你尝试过这些扩展吗？ 在评论中分享您的结果;我很想听听它。
摘要
在本教程中，您了解了如何为多步时间序列预测开发 LSTM 网络。
具体来说，你学到了：
- 如何开发多步时间序列预测的持久性模型。
- 如何开发 LSTM 网络进行多步时间序列预测。
- 如何评估和绘制多步时间序列预测的结果。
您对 LSTM 的多步骤时间序列预测有任何疑问吗？ 在下面的评论中提出您的问题，我会尽力回答。
用于家庭用电机器学习的多步时间序列预测
鉴于智能电表的兴起以及太阳能电池板等发电技术的广泛采用，可提供大量的用电数据。
该数据代表了多变量时间序列的功率相关变量，而这些变量又可用于建模甚至预测未来的电力消耗。
机器学习算法预测单个值，不能直接用于多步预测。可以使用机器学习算法进行多步预测的两种策略是递归和直接方法。
在本教程中，您将了解如何使用机器学习算法开发递归和直接多步预测模型。
完成本教程后，您将了解：
- 如何开发一个评估线性，非线性和集成机器学习算法的框架，用于多步时间序列预测。
- 如何使用递归多步时间序列预测策略评估机器学习算法。
- 如何使用直接的每日和每个引导时间多步时间序列预测策略来评估机器学习算法。
让我们开始吧。
用于家庭用电的机器学习模型的多步骤时间序列预测 照片由 Sean McMenemy ，保留一些权利。
教程概述
本教程分为五个部分;他们是：
- 问题描述
- 加载并准备数据集
- 模型评估
- 递归多步预测
- 直接多步预测
问题描述
'家庭用电量'数据集是一个多变量时间序列数据集，描述了四年内单个家庭的用电量。
该数据是在 2006 年 12 月至 2010 年 11 月之间收集的，并且每分钟收集家庭内的能耗观察结果。
它是一个多变量系列，由七个变量组成（除日期和时间外）;他们是：
- global_active_power ：家庭消耗的总有功功率（千瓦）。
- global_reactive_power ：家庭消耗的总无功功率（千瓦）。
- 电压：平均电压（伏特）。
- global_intensity ：平均电流强度（安培）。
- sub_metering_1 ：厨房的有功电能（瓦特小时的有功电能）。
- sub_metering_2 ：用于洗衣的有功能量（瓦特小时的有功电能）。
- sub_metering_3 ：气候控制系统的有功电能（瓦特小时的有功电能）。
有功和无功电能参考交流电的技术细节。
可以通过从总活动能量中减去三个定义的子计量变量的总和来创建第四个子计量变量，如下所示：
sub_metering_remainder = (global_active_power * 1000 / 60) - (sub_metering_1 + sub_metering_2 + sub_metering_3)
加载并准备数据集
数据集可以从 UCI 机器学习库下载为单个 20 兆字节的.zip 文件：
下载数据集并将其解压缩到当前工作目录中。您现在将拥有大约 127 兆字节的文件“household_power_consumption.txt”并包含所有观察结果。
我们可以使用read_csv()函数来加载数据，并将前两列合并到一个日期时间列中，我们可以将其用作索引。
# load all data
dataset = read_csv('household_power_consumption.txt', sep=';', header=0, low_memory=False, infer_datetime_format=True, parse_dates={'datetime':[0,1]}, index_col=['datetime'])
接下来，我们可以用'_ 标记所有缺失值？_ '具有NaN值的字符，这是一个浮点数。
这将允许我们将数据作为一个浮点值数组而不是混合类型（效率较低）。
# mark all missing values
dataset.replace('?', nan, inplace=True)
# make dataset numeric
dataset = dataset.astype('float32')
我们还需要填写缺失值，因为它们已被标记。
一种非常简单的方法是从前一天的同一时间复制观察。我们可以在一个名为fill_missing()的函数中实现它，该函数将从 24 小时前获取数据的 NumPy 数组并复制值。
# fill missing values with a value at the same time one day ago
def fill_missing(values):
one_day = 60 * 24
for row in range(values.shape[0]):
for col in range(values.shape[1]):
if isnan(values[row, col]):
values[row, col] = values[row - one_day, col]
我们可以将此函数直接应用于 DataFrame 中的数据。
# fill missing
fill_missing(dataset.values)
现在，我们可以使用上一节中的计算创建一个包含剩余子计量的新列。
# add a column for for the remainder of sub metering
values = dataset.values
dataset['sub_metering_4'] = (values[:,0] * 1000 / 60) - (values[:,4] + values[:,5] + values[:,6])
我们现在可以将清理后的数据集版本保存到新文件中;在这种情况下，我们只需将文件扩展名更改为.csv，并将数据集保存为“household_power_consumption.csv”。
# save updated dataset
dataset.to_csv('household_power_consumption.csv')
将所有这些结合在一起，下面列出了加载，清理和保存数据集的完整示例。
# load and clean-up data
from numpy import nan
from numpy import isnan
from pandas import read_csv
from pandas import to_numeric
# fill missing values with a value at the same time one day ago
def fill_missing(values):
one_day = 60 * 24
for row in range(values.shape[0]):
for col in range(values.shape[1]):
if isnan(values[row, col]):
values[row, col] = values[row - one_day, col]
# load all data
dataset = read_csv('household_power_consumption.txt', sep=';', header=0, low_memory=False, infer_datetime_format=True, parse_dates={'datetime':[0,1]}, index_col=['datetime'])
# mark all missing values
dataset.replace('?', nan, inplace=True)
# make dataset numeric
dataset = dataset.astype('float32')
# fill missing
fill_missing(dataset.values)
# add a column for for the remainder of sub metering
values = dataset.values
dataset['sub_metering_4'] = (values[:,0] * 1000 / 60) - (values[:,4] + values[:,5] + values[:,6])
# save updated dataset
dataset.to_csv('household_power_consumption.csv')
运行该示例将创建新文件'household_power_consumption.csv'，我们可以将其用作建模项目的起点。
模型评估
在本节中，我们将考虑如何开发和评估家庭电力数据集的预测模型。
本节分为四个部分;他们是：
- 问题框架
- 评估指标
- 训练和测试集
- 前瞻性验证
问题框架
有许多方法可以利用和探索家庭用电量数据集。
在本教程中，我们将使用这些数据来探索一个非常具体的问题;那是：
鉴于最近的耗电量，未来一周的预期耗电量是多少？
这要求预测模型预测未来七天每天的总有功功率。
从技术上讲，考虑到多个预测步骤，这个问题的框架被称为多步骤时间序列预测问题。利用多个输入变量的模型可以称为多变量多步时间序列预测模型。
这种类型的模型在规划支出方面可能有助于家庭。在供应方面，它也可能有助于规划特定家庭的电力需求。
数据集的这种框架还表明，将每分钟功耗的观察结果下采样到每日总数是有用的。这不是必需的，但考虑到我们对每天的总功率感兴趣，这是有道理的。
我们可以使用 pandas DataFrame 上的 resample（）函数轻松实现这一点。使用参数'D'调用此函数允许按日期时间索引的加载数据按天分组（查看所有偏移别名）。然后，我们可以计算每天所有观测值的总和，并为八个变量中的每一个创建每日耗电量数据的新数据集。
下面列出了完整的示例。
# resample minute data to total for each day
from pandas import read_csv
# load the new file
dataset = read_csv('household_power_consumption.csv', header=0, infer_datetime_format=True, parse_dates=['datetime'], index_col=['datetime'])
# resample data to daily
daily_groups = dataset.resample('D')
daily_data = daily_groups.sum()
# summarize
print(daily_data.shape)
print(daily_data.head())
# save
daily_data.to_csv('household_power_consumption_days.csv')
运行该示例将创建一个新的每日总功耗数据集，并将结果保存到名为“household_power_consumption_days.csv”的单独文件中。
我们可以将其用作数据集，用于拟合和评估所选问题框架的预测模型。
评估指标
预测将包含七个值，一个用于一周中的每一天。
多步预测问题通常分别评估每个预测时间步长。这有助于以下几个原因：
- 在特定提前期评论技能（例如+1 天 vs +3 天）。
- 在不同的交付时间基于他们的技能对比模型（例如，在+1 天的模型和在日期+5 的模型良好的模型）。
总功率的单位是千瓦，并且具有也在相同单位的误差度量将是有用的。均方根误差（RMSE）和平均绝对误差（MAE）都符合这个要求，尽管 RMSE 更常用，将在本教程中采用。与 MAE 不同，RMSE 更能预测预测误差。
此问题的表现指标是从第 1 天到第 7 天的每个提前期的 RMSE。
作为捷径，使用单个分数总结模型的表现以帮助模型选择可能是有用的。
可以使用的一个可能的分数是所有预测天数的 RMSE。
下面的函数evaluate_forecasts()将实现此行为并基于多个七天预测返回模型的表现。
# evaluate one or more weekly forecasts against expected values
def evaluate_forecasts(actual, predicted):
scores = list()
# calculate an RMSE score for each day
for i in range(actual.shape[1]):
# calculate mse
mse = mean_squared_error(actual[:, i], predicted[:, i])
# calculate rmse
rmse = sqrt(mse)
# store
scores.append(rmse)
# calculate overall RMSE
s = 0
for row in range(actual.shape[0]):
for col in range(actual.shape[1]):
s += (actual[row, col] - predicted[row, col])**2
score = sqrt(s / (actual.shape[0] * actual.shape[1]))
return score, scores
运行该函数将首先返回整个 RMSE，无论白天，然后每天返回一系列 RMSE 分数。
训练和测试集
我们将使用前三年的数据来训练预测模型和评估模型的最后一年。
给定数据集中的数据将分为标准周。这些是从周日开始到周六结束的周。
这是使用所选模型框架的现实且有用的方法，其中可以预测未来一周的功耗。它也有助于建模，其中模型可用于预测特定日期（例如星期三）或整个序列。
我们将数据拆分为标准周，从测试数据集向后工作。
数据的最后一年是 2010 年，2010 年的第一个星期日是 1 月 3 日。数据于 2010 年 11 月中旬结束，数据中最接近的最后一个星期六是 11 月 20 日。这给出了 46 周的测试数据。
下面提供了测试数据集的每日数据的第一行和最后一行以供确认。
2010-01-03,2083.4539999999984,191.61000000000055,350992.12000000034,8703.600000000033,3842.0,4920.0,10074.0,15888.233355799992
...
2010-11-20,2197.006000000004,153.76800000000028,346475.9999999998,9320.20000000002,4367.0,2947.0,11433.0,17869.76663959999
每日数据从 2006 年底开始。
数据集中的第一个星期日是 12 月 17 日，这是第二行数据。
将数据组织到标准周内为训练预测模型提供了 159 个完整的标准周。
2006-12-17,3390.46,226.0059999999994,345725.32000000024,14398.59999999998,2033.0,4187.0,13341.0,36946.66673200004
...
2010-01-02,1309.2679999999998,199.54600000000016,352332.8399999997,5489.7999999999865,801.0,298.0,6425.0,14297.133406600002
下面的函数split_dataset()将每日数据拆分为训练集和测试集，并将每个数据组织成标准周。
使用特定行偏移来使用数据集的知识来分割数据。然后使用 NumPy split（）函数将分割数据集组织成每周数据。
# split a univariate dataset into train/test sets
def split_dataset(data):
# split into standard weeks
train, test = data[1:-328], data[-328:-6]
# restructure into windows of weekly data
train = array(split(train, len(train)/7))
test = array(split(test, len(test)/7))
return train, test
我们可以通过加载每日数据集并打印训练和测试集的第一行和最后一行数据来测试此功能，以确认它们符合上述预期。
完整的代码示例如下所示。
# split into standard weeks
from numpy import split
from numpy import array
from pandas import read_csv
# split a univariate dataset into train/test sets
def split_dataset(data):
# split into standard weeks
train, test = data[1:-328], data[-328:-6]
# restructure into windows of weekly data
train = array(split(train, len(train)/7))
test = array(split(test, len(test)/7))
return train, test
# load the new file
dataset = read_csv('household_power_consumption_days.csv', header=0, infer_datetime_format=True, parse_dates=['datetime'], index_col=['datetime'])
train, test = split_dataset(dataset.values)
# validate train data
print(train.shape)
print(train[0, 0, 0], train[-1, -1, 0])
# validate test
print(test.shape)
print(test[0, 0, 0], test[-1, -1, 0])
运行该示例表明，训练数据集确实有 159 周的数据，而测试数据集有 46 周。
我们可以看到，第一行和最后一行的训练和测试数据集的总有效功率与我们定义为每组标准周界限的特定日期的数据相匹配。
(159, 7, 8)
3390.46 1309.2679999999998
(46, 7, 8)
2083.4539999999984 2197.006000000004
前瞻性验证
将使用称为前向验证的方案来评估模型。
这是需要模型进行一周预测的地方，然后该模型的实际数据可用于模型，以便它可以用作在随后一周做出预测的基础。这对于如何在实践中使用模型以及对模型有益，使其能够利用最佳可用数据都是现实的。
我们可以通过分离输入数据和输出/预测数据来证明这一点。
Input, Predict
[Week1] Week2
[Week1 + Week2] Week3
[Week1 + Week2 + Week3] Week4
...
评估此数据集上的预测模型的前瞻性验证方法在下面的函数中提供，名为 evaluate_model（）。
提供 scikit-learn 模型对象作为函数的参数，以及训练和测试数据集。提供了另一个参数n_input，用于定义模型将用作输入以做出预测的先前观察的数量。
关于 scikit-learn 模型如何适合并做出预测的细节将在后面的章节中介绍。
然后使用先前定义的evaluate_forecasts()函数，针对测试数据集评估模型所做的预测。
# evaluate a single model
def evaluate_model(model, train, test, n_input):
# history is a list of weekly data
history = [x for x in train]
# walk-forward validation over each week
predictions = list()
for i in range(len(test)):
# predict the week
yhat_sequence = ...
# store the predictions
predictions.append(yhat_sequence)
# get real observation and add to history for predicting the next week
history.append(test[i, :])
predictions = array(predictions)
# evaluate predictions days for each week
score, scores = evaluate_forecasts(test[:, :, 0], predictions)
return score, scores
一旦我们对模型进行评估，我们就可以总结表现。
下面的函数名为 summarize_scores（），将模型的表现显示为单行，以便与其他模型进行比较。
# summarize scores
def summarize_scores(name, score, scores):
s_scores = ', '.join(['%.1f' % s for s in scores])
print('%s: [%.3f] %s' % (name, score, s_scores))
我们现在已经开始评估数据集上的预测模型的所有元素。
递归多步预测
大多数预测性建模算法将采用一些观测值作为输入并预测单个输出值。
因此，它们不能直接用于进行多步骤时间序列预测。
这适用于大多数线性，非线性和整体机器学习算法。
可以使用机器学习算法进行多步时间序列预测的一种方法是递归地使用它们。
这涉及对一个时间步做出预测，做出预测，并将其作为输入馈送到模型中，以便预测随后的时间步。重复该过程直到预测到所需数量的步骤。
例如：
X = [x1, x2, x3]
y1 = model.predict(X)
X = [x2, x3, y1]
y2 = model.predict(X)
X = [x3, y1, y2]
y3 = model.predict(X)
...
在本节中，我们将开发一个测试工具，用于拟合和评估 scikit-learn 中提供的机器学习算法，使用递归模型进行多步预测。
第一步是将窗口格式的准备好的训练数据转换为单个单变量系列。
下面的to_series()功能会将每周多变量数据列表转换为每日消耗的单变量单变量系列。
# convert windows of weekly multivariate data into a series of total power
def to_series(data):
# extract just the total power from each week
series = [week[:, 0] for week in data]
# flatten into a single series
series = array(series).flatten()
return series
接下来，需要将每日电力的顺序转换成适合于监督学习问题的输入和输出。
预测将是前几天消耗的总功率的一些函数。我们可以选择用作输入的前几天数，例如一周或两周。总会有一个输出：第二天消耗的总功率。
该模型将适合先前时间步骤的真实观察结果。我们需要迭代消耗的每日功率序列并将其分成输入和输出。这称为滑动窗口数据表示。
下面的to_supervised()函数实现了这种行为。
它将每周数据列表作为输入，以及用作创建的每个样本的输入的前几天的数量。
第一步是将历史记录转换为单个数据系列。然后枚举该系列，每个时间步创建一个输入和输出对。这个问题的框架将允许模型学习根据前几天的观察结果预测一周中的任何一天。该函数返回输入（X）和输出（y），以便训练模型。
# convert history into inputs and outputs
def to_supervised(history, n_input):
# convert history to a univariate series
data = to_series(history)
X, y = list(), list()
ix_start = 0
# step over the entire history one time step at a time
for i in range(len(data)):
# define the end of the input sequence
ix_end = ix_start + n_input
# ensure we have enough data for this instance
if ix_end < len(data):
X.append(data[ix_start:ix_end])
y.append(data[ix_end])
# move along one time step
ix_start += 1
return array(X), array(y)
scikit-learn 库允许将模型用作管道的一部分。这允许在拟合模型之前自动应用数据变换。更重要的是，变换以正确的方式准备，在那里准备或适合训练数据并应用于测试数据。这可以在评估模型时防止数据泄漏。
在通过在训练数据集上拟合每个模型之前创建管道来评估模型时，我们可以使用此功能。在使用模型之前，我们将标准化和标准化数据。
下面的make_pipeline()函数实现了这种行为，返回一个可以像模型一样使用的 Pipeline，例如它可以适合它，它可以做出预测。
每列执行标准化和规范化操作。在to_supervised()功能中，我们基本上将一列数据（总功率）分成多个列，例如，七天七天的输入观察。这意味着输入数据中的七列中的每一列将具有用于标准化的不同均值和标准偏差以及用于归一化的不同最小值和最大值。
鉴于我们使用了滑动窗口，几乎所有值都会出现在每列中，因此，这可能不是问题。但重要的是要注意，在将数据拆分为输入和输出之前将数据缩放为单个列会更加严格。
# create a feature preparation pipeline for a model
def make_pipeline(model):
steps = list()
# standardization
steps.append(('standardize', StandardScaler()))
# normalization
steps.append(('normalize', MinMaxScaler()))
# the model
steps.append(('model', model))
# create pipeline
pipeline = Pipeline(steps=steps)
return pipeline
我们可以将这些元素组合成一个名为sklearn_predict()的函数，如下所示。
该函数采用 scikit-learn 模型对象，训练数据，称为历史记录以及用作输入的指定数量的前几天。它将训练数据转换为输入和输出，将模型包装在管道中，使其适合，并使用它做出预测。
# fit a model and make a forecast
def sklearn_predict(model, history, n_input):
# prepare data
train_x, train_y = to_supervised(history, n_input)
# make pipeline
pipeline = make_pipeline(model)
# fit the model
pipeline.fit(train_x, train_y)
# predict the week, recursively
yhat_sequence = forecast(pipeline, train_x[-1, :], n_input)
return yhat_sequence
模型将使用训练数据集中的最后一行作为输入以做出预测。
forecast()函数将使用该模型进行递归多步预测。
递归预测涉及迭代多步预测所需的七天中的每一天。
将模型的输入数据作为input_data列表的最后几个观察值。此列表附有训练数据最后一行的所有观察结果，当我们使用模型做出预测时，它们会添加到此列表的末尾。因此，我们可以从该列表中获取最后的n_input观测值，以实现提供先前输出作为输入的效果。
该模型用于对准备好的输入数据做出预测，并将输出添加到我们将返回的实际输出序列的列表和输入数据列表中，我们将从中输出观察值作为模型的输入。下一次迭代。
# make a recursive multi-step forecast
def forecast(model, input_x, n_input):
yhat_sequence = list()
input_data = [x for x in input_x]
for j in range(7):
# prepare the input data
X = array(input_data[-n_input:]).reshape(1, n_input)
# make a one-step forecast
yhat = model.predict(X)[0]
# add to the result
yhat_sequence.append(yhat)
# add the prediction to the input
input_data.append(yhat)
return yhat_sequence
我们现在拥有使用递归多步预测策略来拟合和评估 scikit-learn 模型的所有元素。
我们可以更新上一节中定义的evaluate_model()函数来调用sklearn_predict()函数。更新的功能如下所示。
# evaluate a single model
def evaluate_model(model, train, test, n_input):
# history is a list of weekly data
history = [x for x in train]
# walk-forward validation over each week
predictions = list()
for i in range(len(test)):
# predict the week
yhat_sequence = sklearn_predict(model, history, n_input)
# store the predictions
predictions.append(yhat_sequence)
# get real observation and add to history for predicting the next week
history.append(test[i, :])
predictions = array(predictions)
# evaluate predictions days for each week
score, scores = evaluate_forecasts(test[:, :, 0], predictions)
return score, scores
一个重要的最终函数是 get_models（），它定义了一个 scikit-learn 模型对象的字典，映射到我们可以用于报告的简写名称。
我们将从评估一套线性算法开始。我们期望这些表现类似于自回归模型（例如，如果使用七天的输入，则为 AR（7））。
具有十个线性模型的get_models()函数定义如下。
这是一个抽查，我们对各种算法的一般表现感兴趣，而不是优化任何给定的算法。
# prepare a list of ml models
def get_models(models=dict()):
# linear models
models['lr'] = LinearRegression()
models['lasso'] = Lasso()
models['ridge'] = Ridge()
models['en'] = ElasticNet()
models['huber'] = HuberRegressor()
models['lars'] = Lars()
models['llars'] = LassoLars()
models['pa'] = PassiveAggressiveRegressor(max_iter=1000, tol=1e-3)
models['ranscac'] = RANSACRegressor()
models['sgd'] = SGDRegressor(max_iter=1000, tol=1e-3)
print('Defined %d models' % len(models))
return models
最后，我们可以将所有这些结合在一起。
首先，加载数据集并将其拆分为训练集和测试集。
# load the new file
dataset = read_csv('household_power_consumption_days.csv', header=0, infer_datetime_format=True, parse_dates=['datetime'], index_col=['datetime'])
# split into train and test
train, test = split_dataset(dataset.values)
然后，我们可以准备模型字典并定义观察的前几天的数量，以用作模型的输入。
# prepare the models to evaluate
models = get_models()
n_input = 7
然后枚举字典中的模型，评估每个模型，总结其分数，并将结果添加到线图中。
下面列出了完整的示例。
# recursive multi-step forecast with linear algorithms
from math import sqrt
from numpy import split
from numpy import array
from pandas import read_csv
from sklearn.metrics import mean_squared_error
from matplotlib import pyplot
from sklearn.preprocessing import StandardScaler
from sklearn.preprocessing import MinMaxScaler
from sklearn.pipeline import Pipeline
from sklearn.linear_model import LinearRegression
from sklearn.linear_model import Lasso
from sklearn.linear_model import Ridge
from sklearn.linear_model import ElasticNet
from sklearn.linear_model import HuberRegressor
from sklearn.linear_model import Lars
from sklearn.linear_model import LassoLars
from sklearn.linear_model import PassiveAggressiveRegressor
from sklearn.linear_model import RANSACRegressor
from sklearn.linear_model import SGDRegressor
# split a univariate dataset into train/test sets
def split_dataset(data):
# split into standard weeks
train, test = data[1:-328], data[-328:-6]
# restructure into windows of weekly data
train = array(split(train, len(train)/7))
test = array(split(test, len(test)/7))
return train, test
# evaluate one or more weekly forecasts against expected values
def evaluate_forecasts(actual, predicted):
scores = list()
# calculate an RMSE score for each day
for i in range(actual.shape[1]):
# calculate mse
mse = mean_squared_error(actual[:, i], predicted[:, i])
# calculate rmse
rmse = sqrt(mse)
# store
scores.append(rmse)
# calculate overall RMSE
s = 0
for row in range(actual.shape[0]):
for col in range(actual.shape[1]):
s += (actual[row, col] - predicted[row, col])**2
score = sqrt(s / (actual.shape[0] * actual.shape[1]))
return score, scores
# summarize scores
def summarize_scores(name, score, scores):
s_scores = ', '.join(['%.1f' % s for s in scores])
print('%s: [%.3f] %s' % (name, score, s_scores))
# prepare a list of ml models
def get_models(models=dict()):
# linear models
models['lr'] = LinearRegression()
models['lasso'] = Lasso()
models['ridge'] = Ridge()
models['en'] = ElasticNet()
models['huber'] = HuberRegressor()
models['lars'] = Lars()
models['llars'] = LassoLars()
models['pa'] = PassiveAggressiveRegressor(max_iter=1000, tol=1e-3)
models['ranscac'] = RANSACRegressor()
models['sgd'] = SGDRegressor(max_iter=1000, tol=1e-3)
print('Defined %d models' % len(models))
return models
# create a feature preparation pipeline for a model
def make_pipeline(model):
steps = list()
# standardization
steps.append(('standardize', StandardScaler()))
# normalization
steps.append(('normalize', MinMaxScaler()))
# the model
steps.append(('model', model))
# create pipeline
pipeline = Pipeline(steps=steps)
return pipeline
# make a recursive multi-step forecast
def forecast(model, input_x, n_input):
yhat_sequence = list()
input_data = [x for x in input_x]
for j in range(7):
# prepare the input data
X = array(input_data[-n_input:]).reshape(1, n_input)
# make a one-step forecast
yhat = model.predict(X)[0]
# add to the result
yhat_sequence.append(yhat)
# add the prediction to the input
input_data.append(yhat)
return yhat_sequence
# convert windows of weekly multivariate data into a series of total power
def to_series(data):
# extract just the total power from each week
series = [week[:, 0] for week in data]
# flatten into a single series
series = array(series).flatten()
return series
# convert history into inputs and outputs
def to_supervised(history, n_input):
# convert history to a univariate series
data = to_series(history)
X, y = list(), list()
ix_start = 0
# step over the entire history one time step at a time
for i in range(len(data)):
# define the end of the input sequence
ix_end = ix_start + n_input
# ensure we have enough data for this instance
if ix_end < len(data):
X.append(data[ix_start:ix_end])
y.append(data[ix_end])
# move along one time step
ix_start += 1
return array(X), array(y)
# fit a model and make a forecast
def sklearn_predict(model, history, n_input):
# prepare data
train_x, train_y = to_supervised(history, n_input)
# make pipeline
pipeline = make_pipeline(model)
# fit the model
pipeline.fit(train_x, train_y)
# predict the week, recursively
yhat_sequence = forecast(pipeline, train_x[-1, :], n_input)
return yhat_sequence
# evaluate a single model
def evaluate_model(model, train, test, n_input):
# history is a list of weekly data
history = [x for x in train]
# walk-forward validation over each week
predictions = list()
for i in range(len(test)):
# predict the week
yhat_sequence = sklearn_predict(model, history, n_input)
# store the predictions
predictions.append(yhat_sequence)
# get real observation and add to history for predicting the next week
history.append(test[i, :])
predictions = array(predictions)
# evaluate predictions days for each week
score, scores = evaluate_forecasts(test[:, :, 0], predictions)
return score, scores
# load the new file
dataset = read_csv('household_power_consumption_days.csv', header=0, infer_datetime_format=True, parse_dates=['datetime'], index_col=['datetime'])
# split into train and test
train, test = split_dataset(dataset.values)
# prepare the models to evaluate
models = get_models()
n_input = 7
# evaluate each model
days = ['sun', 'mon', 'tue', 'wed', 'thr', 'fri', 'sat']
for name, model in models.items():
# evaluate and get scores
score, scores = evaluate_model(model, train, test, n_input)
# summarize scores
summarize_scores(name, score, scores)
# plot scores
pyplot.plot(days, scores, marker='o', label=name)
# show plot
pyplot.legend()
pyplot.show()
运行该示例将评估十个线性算法并总结结果。
评估每个算法并使用单行摘要报告表现，包括总体 RMSE 以及每次步骤 RMSE。
我们可以看到，大多数评估模型表现良好，整周误差低于 400 千瓦，随机随机梯度下降（SGD）回归量表现最佳，总体 RMSE 约为 383。
Defined 10 models
lr: [388.388] 411.0, 389.1, 338.0, 370.8, 408.5, 308.3, 471.1
lasso: [386.838] 403.6, 388.9, 337.3, 371.1, 406.1, 307.6, 471.6
ridge: [387.659] 407.9, 388.6, 337.5, 371.2, 407.0, 307.7, 471.7
en: [469.337] 452.2, 451.9, 435.8, 485.7, 460.4, 405.8, 575.1
huber: [392.465] 412.1, 388.0, 337.9, 377.3, 405.6, 306.9, 492.5
lars: [388.388] 411.0, 389.1, 338.0, 370.8, 408.5, 308.3, 471.1
llars: [388.406] 396.1, 387.8, 339.3, 377.8, 402.9, 310.3, 481.9
pa: [399.402] 410.0, 391.7, 342.2, 389.7, 409.8, 315.9, 508.4
ranscac: [439.945] 454.0, 424.0, 369.5, 421.5, 457.5, 409.7, 526.9
sgd: [383.177] 400.3, 386.0, 333.0, 368.9, 401.5, 303.9, 466.9
还创建了 10 个分类器中每个分类器的每日 RMSE 的线图。
我们可以看到，除了两种方法之外，其他所有方法都会聚集在一起，在七天预测中表现同样出色。
线性算法的递归多步预测线图
通过调整一些表现更好的算法的超参数可以获得更好的结果。此外，更新示例以测试一套非线性和集成算法可能会很有趣。
一个有趣的实验可能是评估一个或几个表现更好的算法的表现，前一天或多或少作为输入。
直接多步预测
多步预测的递归策略的替代方案是对每个要预测的日子使用不同的模型。
这称为直接多步预测策略。
因为我们有兴趣预测七天，所以需要准备七个不同的模型，每个模型专门用于预测不同的一天。
训练这种模型有两种方法：
- 预测日。可以准备模型来预测标准周的特定日期，例如星期一。
- 预测提前期。可以准备模型来预测特定的提前期，例如第 1 天
预测一天将更具体，但意味着每个模型可以使用更少的训练数据。预测提前期会使用更多的训练数据，但需要模型在一周的不同日期进行推广。
我们将在本节中探讨这两种方法。
直接日方法
首先，我们必须更新to_supervised()函数以准备数据，例如用作输入的前一周观察数据以及用作输出的下周特定日期的观察结果。
下面列出了实现此行为的更新的to_supervised()函数。它需要一个参数output_ix来定义下一周的日[0,6]以用作输出。
# convert history into inputs and outputs
def to_supervised(history, output_ix):
X, y = list(), list()
# step over the entire history one time step at a time
for i in range(len(history)-1):
X.append(history[i][:,0])
y.append(history[i + 1][output_ix,0])
return array(X), array(y)
此功能可以调用七次，对于所需的七种型号中的每一种都可以调用一次。
接下来，我们可以更新sklearn_predict()函数，为一周预测中的每一天创建一个新数据集和一个新模型。
函数的主体大部分没有变化，只是在输出序列中每天在循环中使用它，其中“i”的索引被传递给 to_supervised 的调用（ ） 为了准备一个特定的数据集来训练模型来预测那一天。
该函数不再采用n_input参数，因为我们已将输入修复为前一周的七天。
# fit a model and make a forecast
def sklearn_predict(model, history):
yhat_sequence = list()
# fit a model for each forecast day
for i in range(7):
# prepare data
train_x, train_y = to_supervised(history, i)
# make pipeline
pipeline = make_pipeline(model)
# fit the model
pipeline.fit(train_x, train_y)
# forecast
x_input = array(train_x[-1, :]).reshape(1,7)
yhat = pipeline.predict(x_input)[0]
# store
yhat_sequence.append(yhat)
return yhat_sequence
下面列出了完整的示例。
# direct multi-step forecast by day
from math import sqrt
from numpy import split
from numpy import array
from pandas import read_csv
from sklearn.metrics import mean_squared_error
from matplotlib import pyplot
from sklearn.preprocessing import StandardScaler
from sklearn.preprocessing import MinMaxScaler
from sklearn.pipeline import Pipeline
from sklearn.linear_model import LinearRegression
from sklearn.linear_model import Lasso
from sklearn.linear_model import Ridge
from sklearn.linear_model import ElasticNet
from sklearn.linear_model import HuberRegressor
from sklearn.linear_model import Lars
from sklearn.linear_model import LassoLars
from sklearn.linear_model import PassiveAggressiveRegressor
from sklearn.linear_model import RANSACRegressor
from sklearn.linear_model import SGDRegressor
# split a univariate dataset into train/test sets
def split_dataset(data):
# split into standard weeks
train, test = data[1:-328], data[-328:-6]
# restructure into windows of weekly data
train = array(split(train, len(train)/7))
test = array(split(test, len(test)/7))
return train, test
# evaluate one or more weekly forecasts against expected values
def evaluate_forecasts(actual, predicted):
scores = list()
# calculate an RMSE score for each day
for i in range(actual.shape[1]):
# calculate mse
mse = mean_squared_error(actual[:, i], predicted[:, i])
# calculate rmse
rmse = sqrt(mse)
# store
scores.append(rmse)
# calculate overall RMSE
s = 0
for row in range(actual.shape[0]):
for col in range(actual.shape[1]):
s += (actual[row, col] - predicted[row, col])**2
score = sqrt(s / (actual.shape[0] * actual.shape[1]))
return score, scores
# summarize scores
def summarize_scores(name, score, scores):
s_scores = ', '.join(['%.1f' % s for s in scores])
print('%s: [%.3f] %s' % (name, score, s_scores))
# prepare a list of ml models
def get_models(models=dict()):
# linear models
models['lr'] = LinearRegression()
models['lasso'] = Lasso()
models['ridge'] = Ridge()
models['en'] = ElasticNet()
models['huber'] = HuberRegressor()
models['lars'] = Lars()
models['llars'] = LassoLars()
models['pa'] = PassiveAggressiveRegressor(max_iter=1000, tol=1e-3)
models['ranscac'] = RANSACRegressor()
models['sgd'] = SGDRegressor(max_iter=1000, tol=1e-3)
print('Defined %d models' % len(models))
return models
# create a feature preparation pipeline for a model
def make_pipeline(model):
steps = list()
# standardization
steps.append(('standardize', StandardScaler()))
# normalization
steps.append(('normalize', MinMaxScaler()))
# the model
steps.append(('model', model))
# create pipeline
pipeline = Pipeline(steps=steps)
return pipeline
# convert history into inputs and outputs
def to_supervised(history, output_ix):
X, y = list(), list()
# step over the entire history one time step at a time
for i in range(len(history)-1):
X.append(history[i][:,0])
y.append(history[i + 1][output_ix,0])
return array(X), array(y)
# fit a model and make a forecast
def sklearn_predict(model, history):
yhat_sequence = list()
# fit a model for each forecast day
for i in range(7):
# prepare data
train_x, train_y = to_supervised(history, i)
# make pipeline
pipeline = make_pipeline(model)
# fit the model
pipeline.fit(train_x, train_y)
# forecast
x_input = array(train_x[-1, :]).reshape(1,7)
yhat = pipeline.predict(x_input)[0]
# store
yhat_sequence.append(yhat)
return yhat_sequence
# evaluate a single model
def evaluate_model(model, train, test):
# history is a list of weekly data
history = [x for x in train]
# walk-forward validation over each week
predictions = list()
for i in range(len(test)):
# predict the week
yhat_sequence = sklearn_predict(model, history)
# store the predictions
predictions.append(yhat_sequence)
# get real observation and add to history for predicting the next week
history.append(test[i, :])
predictions = array(predictions)
# evaluate predictions days for each week
score, scores = evaluate_forecasts(test[:, :, 0], predictions)
return score, scores
# load the new file
dataset = read_csv('household_power_consumption_days.csv', header=0, infer_datetime_format=True, parse_dates=['datetime'], index_col=['datetime'])
# split into train and test
train, test = split_dataset(dataset.values)
# prepare the models to evaluate
models = get_models()
# evaluate each model
days = ['sun', 'mon', 'tue', 'wed', 'thr', 'fri', 'sat']
for name, model in models.items():
# evaluate and get scores
score, scores = evaluate_model(model, train, test)
# summarize scores
summarize_scores(name, score, scores)
# plot scores
pyplot.plot(days, scores, marker='o', label=name)
# show plot
pyplot.legend()
pyplot.show()
首先运行该示例总结了每个模型的表现。
我们可以看到表现比这个问题的递归模型略差。
Defined 10 models
lr: [410.927] 463.8, 381.4, 351.9, 430.7, 387.8, 350.4, 488.8
lasso: [408.440] 458.4, 378.5, 352.9, 429.5, 388.0, 348.0, 483.5
ridge: [403.875] 447.1, 377.9, 347.5, 427.4, 384.1, 343.4, 479.7
en: [454.263] 471.8, 433.8, 415.8, 477.4, 434.4, 373.8, 551.8
huber: [409.500] 466.8, 380.2, 359.8, 432.4, 387.0, 351.3, 470.9
lars: [410.927] 463.8, 381.4, 351.9, 430.7, 387.8, 350.4, 488.8
llars: [406.490] 453.0, 378.8, 357.3, 428.1, 388.0, 345.0, 476.9
pa: [402.476] 428.4, 380.9, 356.5, 426.7, 390.4, 348.6, 471.4
ranscac: [497.225] 456.1, 423.0, 445.9, 547.6, 521.9, 451.5, 607.2
sgd: [403.526] 441.4, 378.2, 354.5, 423.9, 382.4, 345.8, 480.3
还创建了每个模型的每日 RMSE 得分的线图，显示了与递归模型所见的类似的模型分组。
线性算法的直接每日多步预测线图
直接提前期方法
直接提前期方法是相同的，除了to_supervised()使用更多的训练数据集。
该函数与递归模型示例中定义的函数相同，只是它需要额外的output_ix参数来定义下一周中用作输出的日期。
下面列出了直接每个引导时间策略的更新to_supervised()函数。
与每日策略不同，此版本的功能支持可变大小的输入（不仅仅是七天），如果您愿意，可以进行实验。
# convert history into inputs and outputs
def to_supervised(history, n_input, output_ix):
# convert history to a univariate series
data = to_series(history)
X, y = list(), list()
ix_start = 0
# step over the entire history one time step at a time
for i in range(len(data)):
# define the end of the input sequence
ix_end = ix_start + n_input
ix_output = ix_end + output_ix
# ensure we have enough data for this instance
if ix_output < len(data):
X.append(data[ix_start:ix_end])
y.append(data[ix_output])
# move along one time step
ix_start += 1
return array(X), array(y)
下面列出了完整的示例。
# direct multi-step forecast by lead time
from math import sqrt
from numpy import split
from numpy import array
from pandas import read_csv
from sklearn.metrics import mean_squared_error
from matplotlib import pyplot
from sklearn.preprocessing import StandardScaler
from sklearn.preprocessing import MinMaxScaler
from sklearn.pipeline import Pipeline
from sklearn.linear_model import LinearRegression
from sklearn.linear_model import Lasso
from sklearn.linear_model import Ridge
from sklearn.linear_model import ElasticNet
from sklearn.linear_model import HuberRegressor
from sklearn.linear_model import Lars
from sklearn.linear_model import LassoLars
from sklearn.linear_model import PassiveAggressiveRegressor
from sklearn.linear_model import RANSACRegressor
from sklearn.linear_model import SGDRegressor
# split a univariate dataset into train/test sets
def split_dataset(data):
# split into standard weeks
train, test = data[1:-328], data[-328:-6]
# restructure into windows of weekly data
train = array(split(train, len(train)/7))
test = array(split(test, len(test)/7))
return train, test
# evaluate one or more weekly forecasts against expected values
def evaluate_forecasts(actual, predicted):
scores = list()
# calculate an RMSE score for each day
for i in range(actual.shape[1]):
# calculate mse
mse = mean_squared_error(actual[:, i], predicted[:, i])
# calculate rmse
rmse = sqrt(mse)
# store
scores.append(rmse)
# calculate overall RMSE
s = 0
for row in range(actual.shape[0]):
for col in range(actual.shape[1]):
s += (actual[row, col] - predicted[row, col])**2
score = sqrt(s / (actual.shape[0] * actual.shape[1]))
return score, scores
# summarize scores
def summarize_scores(name, score, scores):
s_scores = ', '.join(['%.1f' % s for s in scores])
print('%s: [%.3f] %s' % (name, score, s_scores))
# prepare a list of ml models
def get_models(models=dict()):
# linear models
models['lr'] = LinearRegression()
models['lasso'] = Lasso()
models['ridge'] = Ridge()
models['en'] = ElasticNet()
models['huber'] = HuberRegressor()
models['lars'] = Lars()
models['llars'] = LassoLars()
models['pa'] = PassiveAggressiveRegressor(max_iter=1000, tol=1e-3)
models['ranscac'] = RANSACRegressor()
models['sgd'] = SGDRegressor(max_iter=1000, tol=1e-3)
print('Defined %d models' % len(models))
return models
# create a feature preparation pipeline for a model
def make_pipeline(model):
steps = list()
# standardization
steps.append(('standardize', StandardScaler()))
# normalization
steps.append(('normalize', MinMaxScaler()))
# the model
steps.append(('model', model))
# create pipeline
pipeline = Pipeline(steps=steps)
return pipeline
# # convert windows of weekly multivariate data into a series of total power
def to_series(data):
# extract just the total power from each week
series = [week[:, 0] for week in data]
# flatten into a single series
series = array(series).flatten()
return series
# convert history into inputs and outputs
def to_supervised(history, n_input, output_ix):
# convert history to a univariate series
data = to_series(history)
X, y = list(), list()
ix_start = 0
# step over the entire history one time step at a time
for i in range(len(data)):
# define the end of the input sequence
ix_end = ix_start + n_input
ix_output = ix_end + output_ix
# ensure we have enough data for this instance
if ix_output < len(data):
X.append(data[ix_start:ix_end])
y.append(data[ix_output])
# move along one time step
ix_start += 1
return array(X), array(y)
# fit a model and make a forecast
def sklearn_predict(model, history, n_input):
yhat_sequence = list()
# fit a model for each forecast day
for i in range(7):
# prepare data
train_x, train_y = to_supervised(history, n_input, i)
# make pipeline
pipeline = make_pipeline(model)
# fit the model
pipeline.fit(train_x, train_y)
# forecast
x_input = array(train_x[-1, :]).reshape(1,n_input)
yhat = pipeline.predict(x_input)[0]
# store
yhat_sequence.append(yhat)
return yhat_sequence
# evaluate a single model
def evaluate_model(model, train, test, n_input):
# history is a list of weekly data
history = [x for x in train]
# walk-forward validation over each week
predictions = list()
for i in range(len(test)):
# predict the week
yhat_sequence = sklearn_predict(model, history, n_input)
# store the predictions
predictions.append(yhat_sequence)
# get real observation and add to history for predicting the next week
history.append(test[i, :])
predictions = array(predictions)
# evaluate predictions days for each week
score, scores = evaluate_forecasts(test[:, :, 0], predictions)
return score, scores
# load the new file
dataset = read_csv('household_power_consumption_days.csv', header=0, infer_datetime_format=True, parse_dates=['datetime'], index_col=['datetime'])
# split into train and test
train, test = split_dataset(dataset.values)
# prepare the models to evaluate
models = get_models()
n_input = 7
# evaluate each model
days = ['sun', 'mon', 'tue', 'wed', 'thr', 'fri', 'sat']
for name, model in models.items():
# evaluate and get scores
score, scores = evaluate_model(model, train, test, n_input)
# summarize scores
summarize_scores(name, score, scores)
# plot scores
pyplot.plot(days, scores, marker='o', label=name)
# show plot
pyplot.legend()
pyplot.show()
运行该示例总结了每个评估的线性模型的总体和每日 RMSE。
我们可以看到，通常每个引导时间方法产生的表现优于每日版本。这可能是因为该方法使更多的训练数据可用于模型。
Defined 10 models
lr: [394.983] 411.0, 400.7, 340.2, 382.9, 385.1, 362.8, 469.4
lasso: [391.767] 403.6, 394.4, 336.1, 382.7, 384.2, 360.4, 468.1
ridge: [393.444] 407.9, 397.8, 338.9, 383.2, 383.2, 360.4, 469.6
en: [461.986] 452.2, 448.3, 430.3, 480.4, 448.9, 396.0, 560.6
huber: [394.287] 412.1, 394.0, 333.4, 384.1, 383.1, 364.3, 474.4
lars: [394.983] 411.0, 400.7, 340.2, 382.9, 385.1, 362.8, 469.4
llars: [390.075] 396.1, 390.1, 334.3, 384.4, 385.2, 355.6, 470.9
pa: [389.340] 409.7, 380.6, 328.3, 388.6, 370.1, 351.8, 478.4
ranscac: [439.298] 387.2, 462.4, 394.4, 427.7, 412.9, 447.9, 526.8
sgd: [390.184] 396.7, 386.7, 337.6, 391.4, 374.0, 357.1, 473.5
再次创建了每日 RMSE 分数的线图。
线性算法的直接每个引导时间多步预测的线图
探索混合使用每日和每时间步骤来模拟问题可能会很有趣。
如果增加用作每个引导时间的输入的前几天的数量是否改善了表现，例如，也可能是有趣的。使用两周的数据而不是一周。
扩展
本节列出了一些扩展您可能希望探索的教程的想法。
- 调谐模型。选择一个表现良好的模型并调整模型超参数以进一步提高表现。
- 调整数据准备。在拟合每个模型之前，所有数据都被标准化并标准化;探索这些方法是否必要以及数据缩放方法的更多或不同组合是否可以带来更好的表现。
- 探索输入大小。输入大小限于先前观察的七天;探索越来越少的观察日作为输入及其对模型表现的影响。
- 非线性算法。探索一套非线性和集成机器学习算法，看看它们是否能提升表现，例如 SVM 和随机森林。
- 多变量直接模型。开发直接模型，利用前一周的所有输入变量，而不仅仅是每日消耗的总功率。这将需要将八个变量的七天的二维数组展平为一维向量。
如果你探索任何这些扩展，我很想知道。
进一步阅读
如果您希望深入了解，本节将提供有关该主题的更多资源。
API
- pandas.read_csv API
- pandas.DataFrame.resample API
- 重采样偏移别名
- sklearn.metrics.mean_squared_error API
- numpy.split API
用品
摘要
在本教程中，您了解了如何使用机器学习算法开发递归和直接多步预测模型。
具体来说，你学到了：
- 如何开发一个评估线性，非线性和集成机器学习算法的框架，用于多步时间序列预测。
- 如何使用递归多步时间序列预测策略评估机器学习算法。
- 如何使用直接的每日和每个引导时间多步时间序列预测策略来评估机器学习算法。
你有任何问题吗？ 在下面的评论中提出您的问题，我会尽力回答。
