CNTK API文档翻译(10)——使用LSTM预测时间序列数据

本篇教程展示如何用CNTK构建LSTM来进行时间序列数据的数值预测。

目标

我们使用一个连续函数的模拟数据集（本例使用正弦曲线）。对于函数y=sin(t)，我们使用符合这个函数的N个值来预测之后的M个值。

在本教程中我们将使用基于LSTM的模型。LSTM比较擅长从以往的数据中学习，因此比较适合我们的教程。

本教程分为三个部分：

生成模拟数据
构建LSTM网络模型
模型训练和评估

LSTM数据已经在很多真实的数据上得到过实践，不过在本教程我们使用简单的模拟数据，下一期我们会使用从一些物联网设备上采集的数据预测太阳能电池板每天的输出电量。

使用CNTK我们能很容易的实现模型。

import math
from matplotlib import pyplot as plt
import numpy as np
import os
import pandas as pd
import time

import cntk as C

在下面的代码中，我们通过检查在CNTK内部定义的环境变量来选择正确的设备（GPU或者CPU）来运行代码，如果不检查的话，会使用CNTK的默认策略来使用最好的设备（如果GPU可用的话就使用GPU，否则使用CPU）

# Select the right target device when this notebook is being tested:
if 'TEST_DEVICE' in os.environ:
    if os.environ['TEST_DEVICE'] == 'cpu':
        C.device.try_set_default_device(C.device.cpu())
    else:
        C.device.try_set_default_device(C.device.gpu(0))

我们设定了两种运行模式：

快速模式：isFast变量设置成True。这是我们的默认模式，在这个模式下我们会训练更少的次数，也会使用更少的数据，这个模式保证功能的正确性，但训练的结果还远远达不到可用的要求。
慢速模式：我们建议学习者在学习的时候试试将isFast变量设置成False，这会让学习者更加了解本教程的内容。

数据生成

我们需要几个辅助函数来生成模拟正弦波数据。也就是让上文说到的N和M分别是正弦波数据的有序集合。

generate_data()

在本教程中，我们会在正弦波上取N个连续的样本，输入模型，试图预测与最后一次观测值M步之后的数据。在我们的训练中，我们会生成多个这样的数据，并且预测与之相对应的值。假设一个取样包的大小是k，那么generate_data函数生成的数据如数据是：
X=[{y[1][1],y[1][2]…y[2][N]},{y[2][1],y[2][2]…y[2][N]}…{y[k][1],y[k][2]…y[k][N]}]
那么他的预测值我们设为X
L=[{y[1][N+M]},{y[2][N+M]}…{y[k][N+M]}]
split_data()

如函数名所诉，split_data会把数据分成训练数据集、验证数据集和测试数据集。

def split_data(data, val_size=0.1, test_size=0.1):
    """
    splits np.array into training, validation and test
    """
    pos_test = int(len(data) * (1 - test_size))
    pos_val = int(len(data[:pos_test]) * (1 - val_size))

    train, val, test = data[:pos_val], data[pos_val:pos_test], data[pos_test:]

    return {"train": train, "val": val, "test": test}

def generate_data(fct, x, time_steps, time_shift):
    """
    generate sequences to feed to rnn for fct(x)
    """
    data = fct(x)
    if not isinstance(data, pd.DataFrame):
        data = pd.DataFrame(dict(a = data[0:len(data) - time_shift],
                                 b = data[time_shift:]))
    rnn_x = []
    for i in range(len(data) - time_steps):
        rnn_x.append(data['a'].iloc[i: i + time_steps].as_matrix())
    rnn_x = np.array(rnn_x)

    # Reshape or rearrange the data from row to columns
    # to be compatible with the input needed by the LSTM model
    # which expects 1 float per time point in a given batch
    rnn_x = rnn_x.reshape(rnn_x.shape + (1,))

    rnn_y = data['b'].values

    # Reshape or rearrange the data from row to columns
    # to match the input shape
    rnn_y = rnn_y.reshape(rnn_y.shape + (1,))

    return split_data(rnn_x), split_data(rnn_y)



N = 5 # input: N subsequent values 
M = 5 # output: predict 1 value M steps ahead
X, Y = generate_data(np.sin, np.linspace(0, 100, 10000, dtype=np.float32), N, M)

f, a = plt.subplots(3, 1, figsize=(12, 8))
for j, ds in enumerate(["train", "val", "test"]):
    a[j].plot(Y[ds], label=ds + ' raw');
[i.legend() for i in a];

生成的数据展示出来如下：

网络模型

我们对每一个输入使用一个LSTM单元。我们有N个输入数据，每个输入数据都是正弦函数上的一个值。从LSTM单元中数出来的值是全连接网络层的输入值，从而生成一个输出值。在LSTM和全连接网络层之间我们加入一个舍弃层（DropoutLayer，关于Dropout的详情可看我的Python与人工神经网络的第七期），舍弃层会随机丢掉从LSTM中出来的百分之二十的数据，从而避免过度拟合。需要注意的是过度拟合实在训练时出现，因此舍弃层只有在训练时才有，预测的时候就不用了。

def create_model(x):
    """Create the model for time series prediction"""
    with C.layers.default_options(initial_state = 0.1):
        m = C.layers.Recurrence(C.layers.LSTM(N))(x)
        m = C.sequence.last(m)
        m = C.layers.Dropout(0.2, seed=1)(m)
        m = C.layers.Dense(1)(m)
        return m

训练网络

我们定义了一个叫next_batch的迭代器，用来生成数据包，提供给训练函数。注意，因为CNTK支持边长序列，我们的数据包必须是一个序列列表。这样也比较方便的生成小一点的数据包，也就是我在前面的文章里说过的取样包。

def next_batch(x, y, ds):
    """get the next batch to process"""

    def as_batch(data, start, count):
        part = []
        for i in range(start, start + count):
            part.append(data[i])
        return np.array(part)

    for i in range(0, len(x[ds])-BATCH_SIZE, BATCH_SIZE):
        yield as_batch(x[ds], i, BATCH_SIZE), as_batch(y[ds], i, BATCH_SIZE)

下面的代码设置了一些训练需要的参数

# Training parameters

TRAINING_STEPS = 10000
BATCH_SIZE = 100
EPOCHS = 20 if isFast else 100

关键知识

在训练之前我们需要说一些在LSTM神经网络中使用序列的关键知识，简单来说：

NTK的输入数据、输出数据和参数全部都用张量表示。每个张量有一个阶数，一个标量是一个0阶张量，一个向量是一个1阶张量，一个矩阵是一个2阶张量等等。我们经常用坐标轴来表示张量的不同维度。

每个CNTK张量有一些固定维度和一些动态的维度。张量在神经网络运行的整个生命周期在固定维度上的长度保持不变，动态维度在定义时与固定维度类似，但有一些不同：

其长度可能根据实力的不同而变化
在训练的取样确定之前，其长度通常是不确定的
他们可能是按顺序列好的

在CNTK里面，如果要运行一个重复迭代的工作，那么他的数据维度是变化的，因此在定义神经网络时我们也是不知道的。所以输入变量就只定义固定的那个维度。打个比方，如果我们的输入数据是一维的，那我们如下定义：

C.sequence.input_variable(1)

在本例中，首先我们有N个已经观测到了的数据，还有持续生成的数据，因此我们会先定义一个默认值：

x_axes = [C.Axis.default_batch_axis(), C.Axis.default_dynamic_axis()]
C.input_variable(1, dynamic_axes=x_axes)

读者应该意识到默认参数的意义，特别是循环赢球在动态维度上按语气的顺序处理时间序列数据。

# input sequences
x = C.sequence.input_variable(1)

# create the model
z = create_model(x)

# expected output (label), also the dynamic axes of the model output
# is specified as the model of the label input
l = C.input_variable(1, dynamic_axes=z.dynamic_axes, name="y")

# the learning rate
learning_rate = 0.001
lr_schedule = C.learning_rate_schedule(learning_rate, C.UnitType.minibatch)

# loss function
loss = C.squared_error(z, l)

# use squared error to determine error for now
error = C.squared_error(z, l)

# use adam optimizer
momentum_time_constant = C.momentum_as_time_constant_schedule(BATCH_SIZE / -math.log(0.9)) 
learner = C.fsadagrad(z.parameters, 
                      lr = lr_schedule, 
                      momentum = momentum_time_constant, 
                      unit_gain = True)

trainer = C.Trainer(z, (loss, error), [learner])

训练一百个周期应该能得到可以接受的结果

# train
loss_summary = []
start = time.time()
for epoch in range(0, EPOCHS):
    for x1, y1 in next_batch(X, Y, "train"):
        trainer.train_minibatch({x: x1, l: y1})
    if epoch % (EPOCHS / 10) == 0:
        training_loss = trainer.previous_minibatch_loss_average
        loss_summary.append(training_loss)
        print("epoch: {}, loss: {:.5f}".format(epoch, training_loss))

print("training took {0:.1f} sec".format(time.time() - start))

输出：

epoch: 0, loss: 0.22282
epoch: 2, loss: 0.19430
epoch: 4, loss: 0.16113
epoch: 6, loss: 0.15103
epoch: 8, loss: 0.11810
epoch: 10, loss: 0.07352
epoch: 12, loss: 0.06572
epoch: 14, loss: 0.07864
epoch: 16, loss: 0.09846
epoch: 18, loss: 0.06972
training took 12.9 sec

通常来说我们需要使用我们之前分离出来的验证数据就来做验证，不过因为输入数据的量比较小，我们可以使用所有数据来做验证。

# validate
def get_mse(X,Y,labeltxt):
    result = 0.0
    for x1, y1 in next_batch(X, Y, labeltxt):
        eval_error = trainer.test_minibatch({x : x1, l : y1})
        result += eval_error
    return result/len(X[labeltxt])

# Print the train and validation errors
for labeltxt in ["train", "val"]:
    print("mse for {}: {:.6f}".format(labeltxt, get_mse(X, Y, labeltxt)))

# Print validate and test error
labeltxt = "test"
print("mse for {}: {:.6f}".format(labeltxt, get_mse(X, Y, labeltxt)))

因为我们使用的是非常简单的正弦数据，因此理论上来说对于训练数据集、验证数据集和测试数据集他们的差值应该是差不多的，但是在真实数据中这种情况一般不会出现。当然我们也把预测的数据展示出来看看和正弦波相差多远。

# predict
f, a = plt.subplots(3, 1, figsize = (12, 8))
for j, ds in enumerate(["train", "val", "test"]):
    results = []
    for x1, y1 in next_batch(X, Y, ds):
        pred = z.eval({x: x1})
        results.extend(pred[:, 0])
    a[j].plot(Y[ds], label = ds + ' raw');
    a[j].plot(results, label = ds + ' predicted');
[i.legend() for i in a];