作者|Christophe Pere
编译|VK
起源|Towards Datas Science

介绍

长期以来，我据说工夫序列问题只能用统计办法(AR[1]，AM[2]，ARMA[3]，ARIMA[4])。这些技术通常被数学家应用，他们试图不断改进这些技术来束缚安稳和非安稳的工夫序列。

几个月前，我的一个敌人(数学家、统计学传授、非安稳工夫序列专家)提出让我钻研如何验证和改良重建恒星光照曲线的技术。事实上，开普勒卫星[11]和其余许多卫星一样，无奈间断测量左近恒星的光通量强度。开普勒卫星在2009年至2016年间致力于寻找太阳系外的行星，称为太阳系外行星或系外行星。

正如你们所了解的，咱们将比咱们的行星地球走得更远一点，并利用机器学习进入河汉之旅。天体物理学始终是我的挚爱。

这个notebook能够在Github上找到：https://github.com/Christophe...。

RNN，LSTM，GRU，双向，CNN-x

那么咱们将在哪个模型上进行这项钻研？咱们将应用循环神经网络(RNN[5])，LSTM[6]、GRU[7]、Stacked LSTM、Stacked GRU、双向LSTM[8]、双向GRU以及CNN-LSTM[9]。

对于那些热衷于树的人，你能够在这里找到一篇对于XGBoost和工夫序列的文章，作者是jasonbrownley。github上提供了一个对于工夫序列的很好的存储库：https://github.com/Jenniferz2...

对于那些不相熟RNN家族的人，把它们看作是具备记忆效应和忘记能力的学习办法。双向来自体系结构，它是指两个RNN，它将在一个方向(从左到右)和另一个方向(从右到左)“读取”数据，以便可能最好地示意长期依赖关系。

数据

如前文所述，这些数据对应于几颗恒星的通量测量值。实际上，在每一个工夫增量(小时)，卫星都会测量来自左近恒星的通量。这个通量，或者说是光强度，随工夫而变动。这有几个起因，卫星的正确挪动、旋转、视角等都会有所不同。因而，测量到的光子数会发生变化，恒星是一个熔化的物质球(氢和氦聚变)，它有本人的静止，因而光子的发射取决于它的静止。这对应于光强度的稳定。

然而，也可能有行星，系外行星，它们烦扰恒星，甚至从恒星之间穿过卫星的眼帘(凌日办法[12])。这条通道遮住了恒星，卫星接管到的光子较少，因为它们被后面通过的行星挡住了(一个具体的例子是月球引起的日食)。

通量测量的汇合被称为光曲线。光曲线是什么样子的？以下是一些示例：

不同恒星之间的通量十分不同。有的乐音很大，有的则很稳固。通量依然出现异样。在光照曲线中能够看到孔或短少测量。咱们的指标是看是否有可能在没有测量的状况下预测光曲线的行为。

数据缩减

为了可能应用模型中的数据，有必要进行数据简化。这里将介绍两种办法，挪动平均法和窗口法。

挪动平均线：

挪动均匀包含取X个间断点并计算它们的平均值。这种办法能够缩小变异性，打消噪声。这也缩小了点的数量，这是一种下采样办法。

上面的函数容许咱们从点列表中计算挪动平均值，办法计算点的平均值和标准差的数字。

def moving_mean(time, flux, lag=5):    '''    该函数通过设定平均值，使数据去噪，缩小数据量。    @param time: (list) 工夫值列表    @param flux: (list) 浮点列表->恒星通量    @param lag: (int) 平均值个数，默认值5    @return X: (list) 工夫调整    @return y: (list) 通量按平均值从新标定    @return y_std: (list) 标准差列表    '''    # 让咱们做一些简略的代码    # 空列表    X = []    y = []    y_std = []    j = 0 # 增量    for i in range(int(len(flux)/lag)):        X.append(np.mean(time[(i+j):(i+j+lag)]))        y.append(np.mean(flux[(i+j):(i+j+lag)]))        y_std.append(np.std(flux[(i+j):(i+j+lag)]))        j+= lag            return X, y, y_stdn

能够看到函数在输出中承受3个参数。工夫和通量是工夫序列的x和y。lag 是管制计算工夫和通量平均值以及通量标准差时所思考的数据个数。

当初，咱们能够看看如何应用这个函数以及通过转换失去的后果。

# #导入所需的包matplotlib inlineimport scipyimport pandas as pdimport numpy as npimport matplotlib.pyplot as pltimport sklearnimport tensorflow as tf# 让咱们看看进度条from tqdm import tqdmtqdm().pandas()

当初咱们须要导入数据。文件kep_lightcurves.csv蕴含着数据。每颗恒星有4列，原始磁通量(“…orig”)，从新缩放的通量是原始磁通量减去均匀磁通量(“…rscl’”)、差值(“…diff”)和残差(“…_res”)。总共52列。

# 20个数据点x, y, y_err  = moving_mean(df.index,df["001724719_rscl"], 20)

df.index示意工夫序列的工夫

df[" 001724719_rscl "] 从新缩放的通量(" 001724719 ")

lag=20是计算平均值和std的数据点的个数

后面3条光照曲线的后果：

窗口办法

第二种办法是窗口办法，它是如何工作的？

你须要取很多点，在前一个例子中是20，而后计算平均值(与后面的办法没有区别)，这个点是新工夫序列的开始，它在地位20(偏移19个点)。然而，窗口不是挪动到下20个点，而是挪动一个点，用之前的20个点计算平均值，而后通过向后挪动一步，以此类推。

这不是一种下采样办法，而是一种清理办法，因为其成果是平滑数据点。

让咱们看看代码：

def mean_sliding_windows(time, flux, lag=5):    '''    该函数通过设定平均值，使数据去噪。    @param time: (list) 工夫值列表    @param flux: (list) 浮点列表->恒星通量    @param lag: (int) 平均值个数，默认值5    @return X: (list) 工夫调整    @return y: (list) 通量按平均值从新标定    @return y_std: (list) 标准差列表    '''    # 让咱们做一些简略的代码    # 空列表    X = []    y = []    y_std = []    j = 0 # 增量    for i in range(int(len(flux)-lag)):                _flux = flux[i:(i+lag)]        _time = time[i:(i+lag)]        X.append(np.mean(_time))        y.append(np.mean(_flux))        y_std.append(np.std(_flux))                         j+= 1 # 咱们只挪动一步            return X, y, y_std

你能够很容易地这样应用它：

# 应用20个点x, y, y_err  = mean_sliding_windows(df.index,df["001724719_rscl"], 40)

df.index示意工夫序列的工夫

df[" 001724719_rscl "] 从新缩放的通量(" 001724719 ")

lag=40是计算平均值和std的数据点的个数

当初，看看后果：

嗯，还不错。将lag设置为40容许“预测”或在小孔中扩大新的工夫序列。然而，如果你认真看，你会发现在红线的开始和完结局部有一个一致。能够改良函数以防止这些伪影。

在接下来的钻研中，咱们将应用挪动平均法取得的工夫序列。

将x轴从值更改为日期：

如果须要日期，能够更改轴。开普勒工作开始于2009年3月7日，完结于2017年。Pandas有一个叫做pd.data_range()的函数。此函数容许你从一直递增的列表中创立日期。

df.index = pd.date_range(‘2009–03–07’, periods=len(df.index), freq=’h’)

这行代码将创立一个频率为小时的新索引。打印后果如下所示。

$ df.indexDatetimeIndex(['2009-03-07 00:00:00', '2009-03-07 01:00:00',               '2009-03-07 02:00:00', '2009-03-07 03:00:00',               '2009-03-07 04:00:00', '2009-03-07 05:00:00',               '2009-03-07 06:00:00', '2009-03-07 07:00:00',               '2009-03-07 08:00:00', '2009-03-07 09:00:00',               ...               '2017-04-29 17:00:00', '2017-04-29 18:00:00',               '2017-04-29 19:00:00', '2017-04-29 20:00:00',               '2017-04-29 21:00:00', '2017-04-29 22:00:00',               '2017-04-29 23:00:00', '2017-04-30 00:00:00',               '2017-04-30 01:00:00', '2017-04-30 02:00:00'],              dtype='datetime64[ns]', length=71427, freq='H')

当初，对于原始工夫序列，你有了一个很好的工夫刻度。

生成数据集

因而，既然曾经创立了数据简化函数，咱们能够将它们组合到另一个函数中(如下所示)，该函数将思考初始数据集和数据集中的恒星名称(这部分能够在函数中实现)。

def reduced_data(df,stars):    '''    Function to automatically reduced a dataset     @param df: (pandas dataframe) 蕴含所有数据的dataframe    @param stars: (list) 蕴含咱们想要简化数据的每个恒星的名称的列表    @return df_mean: 蕴含由缩小平均值的数据的dataframe    @return df_slide: 蕴含通过滑动窗口办法缩小的数据    '''    df_mean = pd.DataFrame()    df_slide = pd.DataFrame()    for i in tqdm(stars):                x , y, y_std = moving_average(df.index, df[i+"_rscl"], lag=25)        df_mean[i+"_rscl_x"] = x        df_mean[i+"_rscl_y"] = y        df_mean[i+"_rscl_y_std"] = y_std                x , y, y_std = mean_sliding_windows(df.index, df[i+"_rscl"], lag=40)        df_slide[i+"_rscl_x"]= x        df_slide[i+"_rscl_y"]= y        df_slide[i+"_rscl_y_std"]= y_std        return df_mean, df_slide

要生成新的数据帧，请执行以下操作：

stars = df.columnsstars = list(set([i.split("_")[0] for i in stars]))print(f"The number of stars available is: {len(stars)}")> The number of stars available is: 13

咱们有13颗恒星，有4种数据类型，对应52列。

df_mean, df_slide = reduced_data(df,stars)

很好，在这一点上，你有两个新的数据集，其中蕴含挪动均匀和窗口办法缩小的数据。

办法

筹备数据：

为了应用机器学习算法来预测工夫序列，必须相应地筹备数据。数据不能仅仅设置在(x，y)个数据点。数据必须采纳序列[x1，x2，x3，…，xn]和预测值y的模式。

上面的函数演示如何设置数据集：

def create_dataset(values, look_back=1):    '''    函数筹备一列(x, y)数据指向用于工夫序列学习的数据    @param values: (list) 值列表    @param look_back: (int) x列表的值[x1, x2, x3，…默认值1    @return _x: x工夫序列的值    @return _y: y工夫序列的值    '''    # 空列表    _x, _y = [], []    for i in range(len(values)-look_back-1):        a = values[i:(i+look_back)]              _x.append(a)                        # 汇合x        _y.append(values[i + look_back]) # 汇合y    return np.array(_x), np.array(_y)

开始之前有两件重要的事。

1.须要从新缩放数据

当数据在[0,1]范畴内时，深度学习算法对工夫序列的预测成果更好。为此，scikit learn提供了MinMaxScaler()函数。你能够配置feature_range参数，但默认值为(0，1)。并革除nan值的数据(如果不删除nan值，则损失函数将输入nan)。

# 缩放数据num = 2 # 抉择数据集中的第三颗星values = df_model[stars[num]+"_rscl_y"].values # 提取值scaler = MinMaxScaler(feature_range=(0, 1)) # 创立MinMaxScaler的实例dataset = scaler.fit_transform(values[~np.isnan(values)].reshape(-1, 1)) # 数据将革除nan值，从新缩放并扭转形态

2.须要将数据转换为x list和y

当初，咱们将应用create_values()函数为模型生成数据。然而，以前，我更喜爱通过以下形式保留原始数据：

df_model = df_mean.save()
# 分成训练和测试集setstrain_size = int(len(dataset) * 0.8)   # 生成80%的训练数据train = dataset[:train_size] # 设置训练数据test  = dataset[train_size:] # 设置测试数据#重塑为X=t和Y=t+1look_back = 20trainX, trainY = create_dataset(train, look_back)testX, testY = create_dataset(test, look_back)# 将输出重塑为[示例、工夫点、特色]trainX = np.reshape(trainX, (trainX.shape[0], trainX.shape[1], 1))testX = np.reshape(testX, (testX.shape[0], testX.shape[1], 1))

看看后果吧

trainX[0]> array([[0.7414906],       [0.76628096],       [0.79901113],       [0.62779976],       [0.64012722],       [0.64934765],       [0.68549234],       [0.64054092],       [0.68075644],       [0.73782449],       [0.68319294],       [0.64330245],       [0.61339268],       [0.62758265],       [0.61779702],       [0.69994317],       [0.64737128],       [0.64122564],       [0.62016833],       [0.47867125]]) # x数据的第一个有20个值trainY[0]> array([0.46174275]) # 对应的y值

度量

咱们用什么指标来预测工夫序列？咱们能够应用均匀绝对误差和均方误差。它们由函数给出：

def metrics_time_series(y_true, y_pred):    '''    从sklearn.metrics计算MAE和MSE度量    @param y_true: (list) 实在值列表    @param y_pred: (list) 预测值列表    @return mae, mse: (float), (float) mae和mse的度量值    '''    mae = round(mean_absolute_error(y_true, y_pred), 2)    mse = round(mean_squared_error(y_true, y_pred), 2)    print(f"The mean absolute error is: {mae}")    print(f"The mean squared error is: {mse}")    return mae, mse

首先须要导入函数：

from sklearn.metrics import mean_absolute_error, mean_squared_error

RNNs：

你能够用几行代码轻松地用Keras实现RNN家族。在这里你能够应用这个性能来配置你的RNN。你须要首先从Keras导入不同的模型，如：

# 导入一些包import tensorflow as tffrom keras.layers import SimpleRNN, LSTM, GRU, Bidirectional, Conv1D, MaxPooling1D, Dropout

当初，咱们有从Keras导入的模型。上面的函数能够生成一个简略的模型(SimpleRNN，LSTM，GRU)。或者，两个模型(雷同的)能够重叠，或者用于双向或两个双向模型的堆栈中。你还能够增加带有MaxPooling1D和dropout的CNN局部(Conv1D)。

def time_series_deep_learning(x_train, y_train, x_test, y_test, model_dl=LSTM ,  unit=4, look_back=20, cnn=False, bidirection=False, stacked=False):    '''    生成不同组合的RNN模型。能够是简略的、重叠的或双向的。模型也能够与CNN局部一起应用。x是(样本、工夫步长、特色)的训练数据    @param x_train: (matrix) 训练数据，维度为 (samples, time steps, features)    @param y_train: (list) 预测    @param x_test: (matrix) 测试数据，维度为 (samples, time steps, features)    @param y_test: (list) 预测    @param model_dl: (model) RNN类型    @param unit: (int) RNN单元数量    @param look_back: (int) x列表中值的数量，配置RNN的形态    @param cnn: (bool) 增加cnn局部模型，默认为false    @param bidirection: (bool) 为RNN增加双向模型，默认为false    @param stacked: (bool) 重叠的两层RNN模型，默认为假    @return train_predict: (list) x_train的预测值    @return train_y: (list) 实在y值    @return test_predict: (list) x_test的预测值    @return test_y: (list) 实在y值    @return (dataframe) 蕴含模型和度量的名称    '''    #配置提前进行的回调，以防止过拟合    es = tf.keras.callbacks.EarlyStopping(        monitor='loss',  patience=5, verbose=0, mode='auto',    )        # 序列模型的实例    model = Sequential()        if cnn: # 如果cnn局部是须要的        print("CNN")        model.add(Conv1D(128, 5, activation='relu'))        model.add(MaxPooling1D(pool_size=4))        model.add(Dropout(0.2))        if not bidirection and not stacked: # 如果须要简略的模型        print("Simple Model")        name = model_dl.__name__        model.add(model_dl(unit,  input_shape=(look_back, 1)))    elif not bidirection: # 测试是否须要双向模型        print("Stacked Model")        name = "Stacked_"+model_dl.__name__        model.add(model_dl(unit,  input_shape=(look_back, 1), return_sequences=True))        model.add(model_dl(unit,  input_shape=(look_back, 1)))    elif not stacked: # 测试是否须要重叠模型        print("Bidirectional Model")        name = "Bi_"+model_dl.__name__        model.add(Bidirectional(model_dl(unit,  input_shape=(look_back, 1))))    else: # 测试是否须要双向和重叠模型        print("Stacked Bidirectional Model")        name = "Stacked_Bi_"+model_dl.__name__        model.add(Bidirectional(model_dl(unit,  input_shape=(look_back, 1), return_sequences=True)))        model.add(Bidirectional(model_dl(unit,  input_shape=(look_back, 1))))            if cnn: # 更新名称与cnn局部        name = "CNN_"+name        # 增加单层浓密层和激活函数线性来预测间断值    model.add(Dense(1))    model.compile(loss='mean_squared_error', optimizer='adam') # MSE loss能够被'mean_absolute_error'代替    model.fit(trainX, trainY, epochs=1000, batch_size=100, callbacks=[es], verbose=0)        # 做出预测    train_predict = model.predict(x_train)    test_predict = model.predict(x_test)        # 反预测    train_predict = scaler.inverse_transform(train_predict)    train_y = scaler.inverse_transform(y_train)    test_predict = scaler.inverse_transform(test_predict)    test_y = scaler.inverse_transform(y_test)        # 计算度量    print("Train")    mae_train, mse_train = metrics_time_series( train_y, train_predict)    print("Test")    mae_test, mse_test = metrics_time_series( test_y, test_predict)        return train_predict, train_y, test_predict, test_y, pd.DataFrame([name, mae_train, mse_train, mae_test, mse_test], index=["Name", "mae_train", "mse_train", "mae_test", "mse_test"]).T

此函数计算训练局部和测试局部的度量，并以数据帧的模式返回后果。举五个例子。

LSTM

# 训练模型并计算度量> x_train_predict_lstm, y_train_lstm,x_test_predict_lstm, y_test_lstm, res= time_series_deep_learning(train_x, train_y, test_x, test_y, model_dl=LSTM ,  unit=12, look_back=20)# 画出预测的后果> plotting_predictions(dataset, look_back, x_train_predict_lstm,  x_test_predict_lstm)# 将每个模型的指标保留在数据框df_results中> df_results = df_results.append(res)

GRU

# 训练模型并计算度量> x_train_predict_lstm, y_train_lstm,x_test_predict_lstm, y_test_lstm, res= time_series_deep_learning(train_x, train_y, test_x, test_y, model_dl=GRU,  unit=12, look_back=20)

重叠LSTM：

# 训练模型并计算度量> x_train_predict_lstm, y_train_lstm,x_test_predict_lstm, y_test_lstm, res= time_series_deep_learning(train_x, train_y, test_x, test_y, model_dl=LSTM ,  unit=12, look_back=20, stacked=True)

双向LSTM：

# 训练模型并计算度量> x_train_predict_lstm, y_train_lstm,x_test_predict_lstm, y_test_lstm, res= time_series_deep_learning(train_x, train_y, test_x, test_y, model_dl=LSTM ,  unit=12, look_back=20, bidirection=True)

CNN-LSTM:

# 训练模型并计算度量> x_train_predict_lstm, y_train_lstm,x_test_predict_lstm, y_test_lstm, res= time_series_deep_learning(train_x, train_y, test_x, test_y, model_dl=LSTM ,  unit=12, look_back=20, cnn=True)

后果

思考到这些数据，后果相当不错。咱们能够看出，深度学习RNN能够很好地再现数据的准确性。下图显示了LSTM模型的预测后果。

表1：不同RNN模型的后果，显示了MAE和MSE指标

  Name    | MAE Train | MSE Train | MAE Test | MSE Test--------------------------------------------------------------------            GRU |   4.24    |   34.11   |   4.15   |   31.47            LSTM |   4.26    |   34.54   |   4.16   |   31.64       Stack_GRU |   4.19    |   33.89   |   4.17   |   32.01      SimpleRNN |   4.21    |   34.07   |   4.18   |   32.41           LSTM |   4.28    |   35.1    |   4.21   |   31.9         Bi_GRU |   4.21    |   34.34   |   4.22   |   32.54  Stack_Bi_LSTM |   4.45    |   36.83   |   4.24   |   32.22        Bi_LSTM |   4.31    |   35.37   |   4.27   |   32.4Stack_SimpleRNN |   4.4     |   35.62   |   4.27   |   33.94      SimpleRNN |   4.44    |   35.94   |   4.31   |   34.37      Stack_LSTM |   4.51    |   36.78   |   4.4    |   34.28 Stacked_Bi_GRU |   4.56    |   37.32   |   4.45   |   35.34       CNN_LSTM |   5.01    |   45.85   |   4.55   |   36.29        CNN_GRU |   5.05    |   46.25   |   4.66   |   37.17   CNN_Stack_GRU |   5.07    |   45.92   |   4.7    |   38.64

表1显示了RNN系列训练集和测试集的均匀绝对误差(MAE)和均方误差(MSE)。GRU在测试集上显示了最好的后果，MAE为4.15，MSE为31.47。

探讨

后果很好，并且重现了不同恒星的光照曲线(见notebook)。然而，稳定并不是齐全重现的，峰值的强度也不雷同，通量也有轻微的偏移。能够通过留神机制进行校对。另一种办法是调整模型、层数(堆栈)、单元数(单元)、不同RNN算法的组合、新的损失函数或激活函数等。

论断

本文展现了将所谓的人工智能办法与工夫序列相结合的可能性。记忆算法(RNN、LSTM、GRU)的弱小性能使得准确再现事件的偶发稳定成为可能。在咱们的例子中，恒星通量体现出相当强烈和显著的稳定，这些办法曾经可能捕捉到。

这项钻研表明，工夫序列不再是统计办法，如ARIMA[4]模型。

参考援用

[1] Autoregressive model, Wikipedia
[2] Moving-average model, Wikipedia
[3] Peter Whittle, 1950. Hypothesis testing in time series analysis. Thesis
[4] Alberto Luceño & Daniel Peña, 2008. Autoregressive Integrated Moving Average (ARIMA) Modeling. Wiley Online Library. https://doi.org/10.1002/97804...
[5] Rumelhart, David E. et al., 1986. Learning representations by back-propagating errors. Nature. 323 (6088): 533–536. 1986Natur.323..533R.
[6] Hochreiter, Sepp & Schmidhuber, Jürgen, 1997. Long Short-Term Memory. Neural Computation. 9 (8): 1735–1780. doi:10.1162/neco.1997.9.8.1735
[7] Cho, KyungHyun et al., 2014. Empirical Evaluation of Gated Recurrent Neural Networks on Sequence Modeling. arXiv:1412.3555
[8] M. Schuster & K.K. Paliwal, 1997. Bidirectional recurrent neural networks. IEEE Transactions on Signal Processing, Volume: 45 , Issue: 11, pp. 2673–2681. DOI: 10.1109/78.650093
[9] Tara N. Sainath et al., 2014. CONVOLUTIONAL, LONG SHORT-TERM MEMORY,FULLY CONNECTED DEEP NEURAL NETWORKS. https://static.googleusercont...
[10] Ashish Vaswani et al., 2017. Attention is all you need. https://arxiv.org/abs/1706.03762
[11] Kepler mission, Nasa

原文链接：https://towardsdatascience.co...

欢送关注磐创AI博客站：
http://panchuang.net/

sklearn机器学习中文官网文档：
http://sklearn123.com/

欢送关注磐创博客资源汇总站：
http://docs.panchuang.net/