上方点击蓝字关注?

在简单了解了LSTM原理之后，本期我将以航班延误预测为例为大家介绍一下如何利用Python编程来构建LSTM模型。

这里我们要用到一个高级的深度学习链接库——Keras,它以TensorFlow或者Theano作为后端引擎，只处理模型的建立、训练和预测等功能，至于底层操作细节，Keras会帮你转化成Theano或TensorFlow相对指令。

Windows下用Anaconda安装Tensorflow和Keras

关于Tensorflow和Keras的安装，大家可以参考以下的博客。

版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议

本文链接：https://blog.csdn.net/lzy_shmily/article/details/80803433

需要注意的两点是：

• 一定要先安装Tensorflow再安装Keras

• 一定要在刚刚创建的tensorflow的环境下安装cpu版本的TensorFlow

利用Keras建立深度学习模型

利用Keras建立深度学习模型，就像在蛋糕架上面一层一层地插蛋糕。在这里，我们选择的“蛋糕架”是Sequential模型，它是多个神经网络层的线性堆叠。那么对于Sequential更多的信息，请参考Keras的中文官方文档(https://keras.io/zh/models/sequential/)。

构建LSTM模型来预测航班延误

接下来我们以航班延误预测为例来说明如何构建LSTM模型。

在航空领域，一架飞机每天要完成多个航班的飞行任务，这就形成了一个航班序列，如果飞机在一个航班任务中发生了延误，那该延误可能会沿着航班序列进行传播。

注：图片来源于网络

当然，飞机也可能通过中途加速或者沿着直线飞行从而将延误抵消。那么我们在这里利用Python构建LSTM模型来预测一架飞机在下一个时刻的到港延误。

数据来源：美国交通部运输统计局提供的公开的航班数据

代码来源：https://www.evolutionarylearn.com/paper/python-keras-tensorflow-mts/

数据预处理

• 剔除掉航班取消和改航记录

• 将数据整理成时间序列数据，也就是将计划到港时间(年-月-日-时-分)作为数据表的索引

导入模块

import time
import math
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.dates as mdates
from datetime import datetime
from keras.models import Sequential
from keras.layers import LSTM, Dense, Dropout
from sklearn.preprocessing import MinMaxScaler
from sklearn.preprocessing import LabelEncoder
from sklearn.preprocessing import OneHotEncoder
from sklearn.metrics import r2_score
from keras import backend as K  # Keras解决OOM超内存问题

建立LSTM类

LSTM类的基本框架是初始化，数据集分割，网络创建与训练，网络评估和结果可视化。

【初始化】：初始化主要是将基本的参数和数据导入进将要创建的实例。代码如下：

    def __init__(self, dataset, hyper_params):
        self.dataset = dataset  # Initialize data sets/数据集初始化
        self.num_neur = hyper_params[0]  # Initialize number of layer and number of neurons in each layer/初始化隐层数和各层神经元个数
        self.look_back = hyper_params[1]  # Initialize length of windows/初始化窗口长度
        self.epochs = hyper_params[2]  # Initialize training times/初始化训练次数
        self.batch_size = hyper_params[3]  # Initialize batch size/初始化批数
        self.selected_feature = hyper_params[4]  # Initialize the selected features/初始化选择特征
        self.train_ratio = hyper_params[5]  # Initialize the splitted ratio of training/初始化训练集分割比例
        self.feature_num = hyper_params[6]  # Initialize the number of features/初始化特征数量
        self.dropoutrate = hyper_params[7]
        self.x_train = []  # Initialize training features of training data set/初始化训练集x部分-训练特征
        self.y_train = []  # Initialize supervisory signals of training data set/初始化训练集y部分-监督信号
        self.x_test = []  # Initialize test features of training data set/初始化测试集x部分-测试特征
        self.y_test = []  # Initialize supervisory signals of training data set/初始化测试集y部分-监督信号

【数据集分割】：数据集分割包括两部分，一部分是特征选择，一部分是分割测试集和训练集。

        def create_dataset(dataset, look_back):
            dataX, dataY = [], []
            for i in range(len(dataset) - look_back):
                a = dataset[i:(i + look_back), 0:dataset.shape[1]]
                dataX.append(a)
                dataY.append(dataset[i + look_back, 0])
            return np.array(dataX), np.array(dataY)

        selected_list = feature_selection(self.selected_feature)  # Index list of selected feature/选择特征列表索引
        train_size = int(len(self.dataset) * self.train_ratio)  # Size of training data set/训练集大小
        train_data = self.dataset[0:train_size, selected_list]  # Training data set/训练集
        test_data = self.dataset[train_size - self.look_back - 1:len(self.dataset), selected_list]  # Test data set/测试集
        self.feature_num = len(selected_list)  # Update the number of selected feature/更新特征数量

        # Data set detail/具体分割后数据集
        x_train, self.y_train = create_dataset(train_data, self.look_back)
        x_test, self.y_test = create_dataset(test_data, self.look_back)

        # Reshape input to be [samples, feature_num, features]/整理特征数据的格式
        self.x_train = np.reshape(x_train, (x_train.shape[0], self.feature_num, x_train.shape[1]))
        self.x_test = np.reshape(x_test, (x_test.shape[0], self.feature_num, x_test.shape[1]))

这段代码中的测试集划分用的是test_data = self.dataset[train_size - self.look_back - 1:len(self.dataset), selected_list]，为什么不是test_data = self.dataset[train_size :len(self.dataset), selected_list]呢？这个是为了保证测试集中有更多样本，后面一种会由于存在窗口损失掉一部分数据，减一应该是为了作图时连接训练部分和测试部分。

【网络创建】：开始在蛋糕架上插蛋糕~

    def lstm(self):
        start_cr_a_fit_net = time.time()  # Record time/记录网络创建与训练时间
        self.split_dataset()  # Split the data set/数据分割

        # Create and fit the LSTM network/创建并拟合LSTM网络
        LSTM_model = Sequential()
        for i in range(len(self.num_neur)):  # 构建多层网络
            if len(self.num_neur) == 1:
                LSTM_model.add(LSTM(self.num_neur[i], input_shape=(None, self.look_back),dropout=self.dropoutrate))
            else:
                if i 1:
                    LSTM_model.add(LSTM(self.num_neur[i], input_shape=(None, self.look_back), return_sequences=True))
                    LSTM_model.add(Dropout(self.dropoutrate))
                else:
                    LSTM_model.add(LSTM(self.num_neur[i], input_shape=(None, self.look_back)))

        LSTM_model.add(Dense(1))
        LSTM_model.summary()  # Summary the structure of neural network/网络结构总结
        LSTM_model.compile(loss='mean_squared_error', optimizer='adam')  # Compile the neural network/编译网络
        LSTM_model.fit(self.x_train, self.y_train, epochs=self.epochs, batch_size=self.batch_size
                       , verbose=2)  # Fit the LSTM network/拟合LSTM网络
        end_cr_a_fit_net = time.time() - start_cr_a_fit_net
        print('Running time of creating and fitting the LSTM network: %.2f Seconds' % (end_cr_a_fit_net))

        # LSTM prediction/LSTM进行预测
        trainPredict = LSTM_model.predict(self.x_train)  # Predict by training data set/训练集预测
        testPredict = LSTM_model.predict(self.x_test)  # Predict by test data set/测试集预测
        return trainPredict, testPredict, self.y_train, self.y_test

【网络评估】这里我们的评估指标选用的是。

    def mape(self, scaler, trainPredict, testPredict):
        # Invert predictions start / 将预测值转换为正常数值
        # Create empty table like the dataset/创建一个空的数组, 结构同dataset
        trainPredict_dataset_like = np.zeros(shape=(len(trainPredict), self.dataset.shape[1]))
        # Put the predicted values in the right field/将预测值填充进新建数组
        trainPredict_dataset_like[:, 0] = trainPredict[:, 0]
        # Inverse transform and then select the right field/数据转换
        trainPredict = scaler.inverse_transform(trainPredict_dataset_like)[:, 0]

        y_train_dataset_like = np.zeros(shape=(len(self.y_train), self.dataset.shape[1]))
        y_train_dataset_like[:, 0] = self.y_train
        self.y_train = scaler.inverse_transform(y_train_dataset_like)[:, 0]

        testPredict_dataset_like = np.zeros(shape=(len(testPredict), self.dataset.shape[1]))
        testPredict_dataset_like[:, 0] = testPredict[:, 0]
        testPredict = scaler.inverse_transform(testPredict_dataset_like)[:, 0]

        y_test_dataset_like = np.zeros(shape=(len(self.y_test), self.dataset.shape[1]))
        y_test_dataset_like[:, 0] = self.y_test
        self.y_test = scaler.inverse_transform(y_test_dataset_like)[:, 0]
        # Invert predictions end/数据转换结束

        # 计算R2值
        train_R2 = r2_score(self.y_train, trainPredict)
        test_R2 = r2_score(self.y_test, testPredict)
        #trainMAPE = np.mean(np.abs(self.y_train - trainPredict) / self.y_train)
        #testMAPE = np.mean(np.abs(self.y_test - testPredict) / self.y_test)

        print("Train R2: " + str(round(train_R2, 2)))
        print("Test R2: " + str(round(test_R2, 2)))
        return trainPredict, testPredict, train_R2, test_R2

【可视化】

    def plot(self, scaler, trainPredict, testPredict):
        # Shift training predictions for plotting/转换数据结构用于作图-训练预测结果
        sub_traindataset = [[data] for data in self.dataset[:, 0]]
        trainPredictPlot = np.empty_like(sub_traindataset)
        trainPredictPlot[:, 0] = np.nan
        trainPredictPlot[self.look_back:len(trainPredict) + self.look_back, 0] = trainPredict

        # Shift test predictions for plotting/转换数据结构用于作图-测试预测结果
        sub_testdataset = [[data] for data in self.dataset[:, 0]]
        testPredictPlot = np.empty_like(sub_testdataset)
        testPredictPlot[:] = np.nan
        testPredictPlot[len(trainPredict) + self.look_back - 1:len(self.dataset), 0] = testPredict

        # plot baseline and predictions/作图
        datasety_like = np.zeros(shape=(self.dataset.shape[0], self.dataset.shape[1]))
        datasety_like[:, 0] = self.dataset[:, 0]
        y = scaler.inverse_transform(datasety_like)[:, 0]

        #dates = pd.date_range('2010-12', periods=len(y), freq='M')
        xs = [d for d in dataframe.index]
        # 配置横坐标
        #plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%Y-%m'))
        #plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[1, 7]))

        A, = plt.plot(xs, y[0:len(y)], linewidth='2', color='cornflowerblue')  # 真实值
        B, = plt.plot(xs, trainPredictPlot, linewidth='1.5', color='lightgreen')  # LSTM训练集结果
        C, = plt.plot(xs, testPredictPlot, linewidth='4', color='lightcoral')  # LSTM测试集结果

        # plt.plot(NpredYPlot,linewidth = '3',color='k')
        #plt.axvline(xs[76], linewidth='2', color='black')  # 画直线区分训练部分与测试部分
        plt.legend((A, B, C), ('real_value', 'LSTM_train', 'LSTM_test'), loc='best')
        plt.gcf().autofmt_xdate()  # 自动旋转日期标记
        plt.xticks([])#由于时间太多，不显示横坐标轴
        #plt.xlabel('Date', family='Times New Roman', fontsize=16)  # X轴
        plt.ylabel('Delay', family='Times New Roman', fontsize=16)  # Y轴

        plt.title('LSTM', family='Times New Roman', fontsize=16)  # 添加标题

        #plt.savefig(r'C:\Users\10321\Desktop\result.png', dpi=900)  # 保存图片

        plt.show()
        del trainPredictPlot, testPredictPlot

导入参数

num_neur=[15, 10]表示这个LSTM网络有两个隐藏层，第一层的神经元个数为15，第二层的神经元个数为10；select_feature = [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]表示我选择了数据集中所有的特征，不选择的设为0就好。所选特征包括：到港延误、季度、周几、出发机场、到达机场、计划和实际运行时间、距离、在空时间、计划过站时间、到达机场每小时的计划离港流量和计划到港流量。

if __name__ == "__main__":
    # Load the dataset/导入数据集
    file = r'E:/ZhangBo/data/aircraft_sequential_sarr.csv'
    dataframe = pd.read_csv(file, index_col=0)
    dataframe = dataframe.dropna(axis=0,how='any')
    encoder = LabelEncoder()
    dataset = dataframe.values
    dataset[:, 3] = encoder.fit_transform(dataset[:, 3])
    dataset[:, 4] = encoder.fit_transform(dataset[:, 4])
    dataset = dataset.astype('float32')

    # Normalize the dataset/标准化数据集
    scaler = MinMaxScaler(feature_range=(0, 1))
    dataset = scaler.fit_transform(dataset)

    # Set hyper-parameters/设定超参数
    num_neur = [15,10]  # Number of layer and number of neurons in each layer/隐藏层数和各层神经元个数
    look_back = 15 # Length of windows/窗口长度
    epochs = 30  # Training times/训练次数
    batch_size = 10 # Batch size/批数大小
    select_feature = [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]# Selected features list/被选择特征列表
    train_ratio = 0.832  # Splitted ratio of training data set/训练集分割比例
    feature_num = dataset.shape[1]  # Feature number+y/特征数量+1,也将预测项作为特征
    dropoutrate = 0.02
    # Hyper-parameter list/超参数列表
    hyper_params = [num_neur, look_back, epochs, batch_size, select_feature, train_ratio, feature_num, dropoutrate]

    # Start an LSTM model/开始一个LSTM网络
    model = lstm(dataset, hyper_params)  # Create instance of LSTM/实例化模型
    trainPredict, testPredict, y_train, y_test = model.lstm()  # Create and fit the LSTM network/创建并拟合LSTM网络
    trainPredict, testPredict,train_R2, test_R2 = model.mape(scaler, trainPredict , testPredict)  # Evaluate network performance/评估网络效果
    model.plot(scaler, trainPredict, testPredict)  # Visualization results/可视化结果
    K.clear_session()  # 关掉内存中神经网络

程序运行结果

【网络结构】

【拟合网络、训练和测试结果】

我们创建和拟合LSTM网络所需的时间为225.81秒，训练集的值为0.55，测试集的值为0.44，可以看出，模型的预测结果并不是很好，因为这里我们并没有对超参数进行选择，而且航班延误的一个很重要的影响因素——天气的相关数据我们并没有加进去。

本期结语

到这里，LSTM网络的理论及其在交通建模中的应用的相关介绍就全部结束了。

如果你有问题和不同的看法，可以在后台留言，我们互相交流，互相学习，共同进步~

编辑：庄桢

“交通科研Lab”：分享学习点滴，期待科研交流！

如果觉得还不错

点这里！???