本示例说明如何使用长短期记忆网络(LSTM)预测时间序列数据。
要预测序列未来时间步长的值,可以在序列回归LSTM网络训练序列,其中的响应是训练序列,其值偏移了一个步长。也就是说,再输入序列的每个时间步长处,LSTM网络都学会预测下一时间步长的值。
要预测将来多个时间步长的值,请使用该 predictAndUpdateState
功能一次预测一个时间步长,并在每次预测时更新网络状态。
本文使用数据集chickenpox_dataset。该示例训练了一个LSTM网络,根据前几个月的病例数来预测水痘病例数。
加载示例数据
chickenpox_dataset包含单个时间序列,时间步长对应于月份,值对应于案例数。输出是一个单元格数组,其中每个元素都是单个时间步长。将数据重塑为行向量
data = chickenpox_dataset;
data = [data{:}];
figure
plot(data)
xlabel('Month'); ylabel('Cases'); title('Monthly Cases of Chickenpox');
图1. 时间序列
划分训练集和测试集,序列的前90%数据为训练集,后10%数据为测试集。
numTimeStepsTrain = floor(0.9*numel(data));
dataTrain = data(1:numTimeStepsTrain+1);
dataTest = data(numTimeStepsTrain+1:end);
标准化数据
为了更好地拟合数据并且防止训练过程中的参数偏离,我们必须归一化训练数据使其成为零均值和单位方差的数据。同样地,我们在预测的时候,也需要对测试集数据进行归一化。
mu = mean(dataTrain);
sig = std(dataTrain);
dataTrainStandardized = (dataTrain - mu) / sig;
准备预测变量和响应
要预测序列未来时间步长的值,请将响应指定为训练序列,其值偏移一个时间步长。也就是说,在输入序列的每个时间步长处,LSTM网络都学会预测下一时间步长的值。预测变量是没有最终时间步长的训练序列。
XTrain = dataTrainStandardized(1:end-1);
YTrain = dataTrainStandardized(2:end);
定义LSTM网络架构
创建一个LSTM回归网络,指定LSTM层具有200个隐藏单元。
numFeatures = 1;
numResponses = 1;
numHiddenUnits = 200;
layers = [...
sequenceInputLayer(numFeatures) %输入层为输入数据特征个数:1
lstmLayer(numHiddenUnits) %隐藏层为神经元个数:200
fullyConnectedLayer(numResponses) %全连接层为响应个数:1
]
确认一下训练过程的参数选择,设定训练器为 ‘adam’并且训练250个epochs。为了防止梯度爆炸,设定梯度阈值为1,将初始学习率设为0.005,在训练完125个epochs后设定学习率为0.005*0.2。
options = trainingOptions('adam', ...
'MaxEpochs', 250, ...
'GradientThreshold', 1, ...
'InitialLearnRate', 0.005, ...
'LearnRateSchedule', 'piecewise', ...
'LearnRateDropPeriod', 125, ...
'LearnRateDropFactor', 0.2, ...
'Verbose', 0, ...
'Plots', 'training-progress';
)
训练LSTM网络
使用 trainNetwork
训练LSTM网络。
net = trainNetwork(XTrain, YTrain, layers, options);
预测时间步长
要预测未来多个时间步长的值,使用 predictAndUpdateState
功能一次预测一个时间步长,并在每次预测时更新网络状态。对于每个预测,将先前的预测作为函数的输入。使用与训练数据相同的参数标准化测试数据。
dataTestStandardized = (dataTest-mu) / sig;
XTest = dataTestStandardized(1:end-1);
要初始化网络状态,首先要对训练数据xTrain进行预测,使用训练响应的最后一个时间步进行第一个预测YTrain(end)。循环遍历其余预测并将先前的预测输入 predictAndUpdateState
。
如果数据量比较大,对于长序列和大型网络来说,在GPU上的预测通常要比在CPU上预测得要快。对于单个时间步长预测,使用CPU预测即可。将 ExecutionEnvironment
得选项’predictAndUpdateState’设置为’cpu’。
net = predictAndUpdateState(net, XTrain);
[net, YPred] = predictAndUpdateState(net, YTrain(end));
numTimeStepsTest = numel(XTest);
for i = 2:numTimeStepsTest
[net, YPred(:, i)] = predictAndUpdateState(net, YPred(:, i-1), 'ExecutionEnvironment', 'cpu');
end
因为上面我们对测试集数据做了标准化处理,所以这里要做去标准化
yPred = sig*YPred + mu;
计算预测数据和真实数据间得RMSE
YTest = dataTest(2:end);
rmse = sqrt(mean(YPred-YTest).^2);
rmse = single 248.5531
画出训练部分原始数据+预测数据
figure
plot(dataTrain(1:end-1));
hold on;
idx = numTimeStepsTrain : (numTimeStepsTrain+numTimesStepsTest);
plot(idx, [data(numTimesStepsTrain), YPred], '.-');
hold off;
xlabel('Month'); ylabel('Cases'); title('Forecast'); legend(['Observed', 'Forecast']);
图3. 原始数据+预测数据
比较预测数据和真实数据
figure
subplot(2, 1, 1);
plot(YTest); hold on;
plot(YPred, '.-'); hold off;
legend(['Observed', 'Forecast']);
ylabel('Cases'); title('Forecast');
subplot(2, 1, 2);
stem(YPred - YTest);
xlabel('Month'); ylabel('Error'); title('RMSE = ' + rmse);
使用观测值来更新网络状态
如果能够得到预测结果之间的真实步长时间,就可使用真实值而不必使用预测值来更新网络状态。
首先,初始化以下网络状态。当我们对一个新的时间序列做预测的时候,我们需要使用 resetState
来重置一下网络,从而防止使用原先的数据集所做的预测会影响到对新的数据所做的预测。在更新完网络状态后,通过对训练数据进行预测然后初始化网络状态。
net = resetState(net);
net = predictAndUpdateState(net, XTrain);
对每个时间步长进行预测。对于每段预测,使用原先的观测值步长来预测后来的观测值步长。设定函数predictAndUpdateState
的参数 ExecutionEnvironment
的选项为’cpu’。
YPred = [];
numTimeStepsTest = numel(XTest);
for i = 1:numTimeStepsTest
[net, YPred(:, i)] = predictAndUpdateState(net, XTest(:, i), 'ExecutionEnvironment', 'cpu');
end
使用原先的方法去标准化
YPred = sig*YPred + mu;
计算RMSE
rmse = sqrt(mean(YPred-YTest).^2);
rmse = 158.0959
将预测结果和真实结果进行比较
figure;
subplot(2, 1, 1)
plot(YTest); hold on;
plot(YPred, '.-'); hold off;
legend(['Observed', 'Predicted']);
ylabel('Cases'); title('Forecast with Updates');
subplot(2, 1, 2);
stem(YPred - YTest);
xlabel('Month'); ylabel('Error'); title('RMSE = ' + rmse);
图5. 预测和真实结果的比较
我们使用真实数据而不是使用预测数据来更新网络状态的好处就是,最后预测结果更加精确。