本文將主要講述如何使用BLiTZ（PyTorch貝葉斯深度學習庫）來建立貝葉斯LSTM模型，以及如何在其上使用序列數據進行訓練與推理。
在本文中，我們將解釋貝葉斯長期短期記憶模型（LSTM）是如何工作的，然後通過一個Kaggle數據集進行股票置信區間的預測。
貝葉斯LSTM層
眾所周知，LSTM結構旨在解決使用標準的循環神經網絡（RNN）處理長序列數據時發生的信息消失問題。
在數學上，LSTM結構的描述如下：
我們知道，貝葉斯神經網絡的核心思想是，相比設定一個確定的權重，我們可以通過一個概率密度分布來對權重進行採樣，然後優化分布參數。
利用這一點，就有可能衡量我們所做的預測的置信度和不確定性，這些數據與預測本身一樣，都是非常有用的數據。
從數學上講，我們只需要在上面的方程中增加一些額外的步驟，也即權值和偏置的採樣，這發生在前向傳播之前。
這表示在第i次在模型第N層上權重的採樣。
這表示在第i次在模型第N層上偏置的採樣。
當然，我們的可訓練參數是和，用來表示不同的權重分布。 BLiTZ具有內置的BayesianLSTM層，可以為您完成所有這些艱苦的工作，因此您只需要關注您的網絡結構設計與網絡的訓練/測試。
現在我們看一個例子。
第一步，先導入庫
除了導入深度學習中最常用的庫外，我們還需要從blitz.modules中導入BayesianLSTM，並從blitz.utils導入variational_estimator，後者是一個用於變量訓練與複雜度計算的裝飾器。
我們還要導入collections.deque來執行時間序列數據的預處理。
import pandas as pdimport numpy as npimport torchimport torch.nn as nnimport torch.optim as optimimport torch.nn.functional as Ffrom blitz.modules import BayesianLSTMfrom blitz.utils import variational_estimatorfrom sklearn.model_selection import train_test_splitfrom sklearn.preprocessing import StandardScalerimport matplotlib.pyplot as plt%matplotlib inlinefrom collections import deque
數據預處理
現在，我們將創建並預處理數據集以將其輸入到網絡。 我們將從Kaggle數據集中導入Amazon股票定價，獲取其「收盤價」數據並將其標準化。
我們的數據集將由標準化股票價格的時間戳組成，並且具有一個形如(batch_size，sequence_length，observation_length)的shape。
下面我們導入數據並對其預處理：
#importing the datasetamazon="data/AMZN_2006-01-01_to_2018-01-01.csv"ibm="data/IBM_2006-01-01_to_2018-01-01.csv"df = pd.read_csv(ibm)#scaling and selecting dataclose_prices = df["Close"]scaler = StandardScaler()close_prices_arr = np.array(close_prices).reshape(-1, 1)close_prices = scaler.fit_transform(close_prices_arr)close_prices_unscaled = df["Close"]
我們還必須創建一個函數來按照時間戳轉換我們的股價歷史記錄。 為此，我們將使用最大長度等於我們正在使用的時間戳大小的雙端隊列，我們將每個數據點添加到雙端隊列，然後將其副本附加到主時間戳列表：
def create_timestamps_ds(series, timestep_size=window_size): time_stamps = [] labels = [] aux_deque = deque(maxlen=timestep_size) #starting the timestep deque for i in range(timestep_size): aux_deque.append(0) #feed the timestamps list for i in range(len(series)-1): aux_deque.append(series[i]) time_stamps.append(list(aux_deque)) #feed the labels lsit for i in range(len(series)-1): labels.append(series[i + 1]) assert len(time_stamps) == len(labels), "Something went wrong" #torch-tensoring it features = torch.tensor(time_stamps[timestep_size:]).float() labels = torch.tensor(labels[timestep_size:]).float() return features, labels
創建神經網絡類
我們的網絡類接收variantal_estimator裝飾器，該裝飾器可簡化對貝葉斯神經網絡損失的採樣。我們的網絡具有一個貝葉斯LSTM層，參數設置為in_features = 1以及out_features = 10，後跟一個nn.Linear(10, 1)，該層輸出股票的標準化價格。
@variational_estimatorclass NN(nn.Module):def __init__(self): super(NN, self).__init__() self.lstm_1 = BayesianLSTM(1, 10) self.linear = nn.Linear(10, 1) def forward(self, x): x_, _ = self.lstm_1(x) #gathering only the latent end-of-sequence for the linear layer x_ = x_[:, -1, :] x_ = self.linear(x_) return x_
如您所見，該網絡可以正常工作，唯一的不同點是BayesianLSTM層和variantal_estimator裝飾器，但其行為與一般的Torch對象相同。
完成後，我們可以創建我們的神經網絡對象，分割數據集並進入訓練循環：
創建對象
我們現在可以創建損失函數、神經網絡、優化器和dataloader。請注意，我們不是隨機分割數據集，因為我們將使用最後一批時間戳來計算模型。由於我們的數據集很小，我們不會對訓練集創建dataloader。
Xs, ys = create_timestamps_ds(close_prices)X_train, X_test, y_train, y_test = train_test_split(Xs,ys, test_size=.25, random_state=42, shuffle=False)ds = torch.utils.data.TensorDataset(X_train, y_train)dataloader_train = torch.utils.data.DataLoader(ds, batch_size=8, shuffle=True)net = NN()criterion = nn.MSELoss()optimizer = optim.Adam(net.parameters(), lr=0.001)
我們將使用MSE損失函數和學習率為0.001的Adam優化器
訓練循環
對於訓練循環，我們將使用添加了variational_estimator的sample_elbo方法。 它對X個樣本的損失進行平均，並幫助我們輕鬆地用蒙特卡洛估計來計算損失。
為了使網絡正常工作，網絡forward方法的輸出必須與傳入損失函數對象的標籤的形狀一致。
iteration = 0for epoch in range(10):for i, (datapoints, labels) in enumerate(dataloader_train): optimizer.zero_grad() loss = net.sample_elbo(inputs=datapoints, labels=labels, criterion=criterion, sample_nbr=3) loss.backward() optimizer.step() iteration += 1 if iteration%250==0: preds_test = net(X_test)[:,0].unsqueeze(1) loss_test = criterion(preds_test, y_test) print("Iteration: {} Val-loss: {:.4f}".format(str(iteration), loss_test))
評估模型並計算置信區間
我們將首先創建一個具有要繪製的真實數據的dataframe：
original = close_prices_unscaled[1:][window_size:]df_pred = pd.DataFrame(original)df_pred["Date"] = df.Datedf["Date"] = pd.to_datetime(df_pred["Date"])df_pred = df_pred.reset_index()
要預測置信區間，我們必須創建一個函數來預測同一數據X次，然後收集其均值和標準差。 同時，在查詢真實數據之前，我們必須設置將嘗試預測的窗口大小。
讓我們看一下預測函數的代碼：
def pred_stock_future(X_test,future_length, sample_nbr=10): #sorry for that, window_size is a global variable, and so are X_train and Xs global window_size global X_train global Xs global scaler #creating auxiliar variables for future prediction preds_test = [] test_begin = X_test[0:1, :, :] test_deque = deque(test_begin[0,:,0].tolist(), maxlen=window_size) idx_pred = np.arange(len(X_train), len(Xs)) #predict it and append to list for i in range(len(X_test)): #print(i) as_net_input = torch.tensor(test_deque).unsqueeze(0).unsqueeze(2) pred = [net(as_net_input).cpu().item() for i in range(sample_nbr)] test_deque.append(torch.tensor(pred).mean().cpu().item()) preds_test.append(pred) if i % future_length == 0: #our inptus become the i index of our X_test #That tweak just helps us with shape issues test_begin = X_test[i:i+1, :, :] test_deque = deque(test_begin[0,:,0].tolist(), maxlen=window_size) #preds_test = np.array(preds_test).reshape(-1, 1) #preds_test_unscaled = scaler.inverse_transform(preds_test) return idx_pred, preds_test
我們要將置信區間保存下來，確定我們置信區間的寬度。
def get_confidence_intervals(preds_test, ci_multiplier):global scaler preds_test = torch.tensor(preds_test) pred_mean = preds_test.mean(1) pred_std = preds_test.std(1).detach().cpu().numpy() pred_std = torch.tensor((pred_std)) upper_bound = pred_mean + (pred_std * ci_multiplier) lower_bound = pred_mean - (pred_std * ci_multiplier) #gather unscaled confidence intervals pred_mean_final = pred_mean.unsqueeze(1).detach().cpu().numpy() pred_mean_unscaled = scaler.inverse_transform(pred_mean_final) upper_bound_unscaled = upper_bound.unsqueeze(1).detach().cpu().numpy() upper_bound_unscaled = scaler.inverse_transform(upper_bound_unscaled) lower_bound_unscaled = lower_bound.unsqueeze(1).detach().cpu().numpy() lower_bound_unscaled = scaler.inverse_transform(lower_bound_unscaled) return pred_mean_unscaled, upper_bound_unscaled, lower_bound_unscaled
由於我們使用的樣本數量很少，因此用一個很高的標準差對其進行了補償。 我們的網絡將嘗試預測7天，然後將參考數據：
future_length=7sample_nbr=4ci_multiplier=10idx_pred, preds_test = pred_stock_future(X_test, future_length, sample_nbr)pred_mean_unscaled, upper_bound_unscaled, lower_bound_unscaled = get_confidence_intervals(preds_test,ci_multiplier)
我們可以通過查看實際值是否低於上限並高於下限來檢查置信區間。 設置好參數後，您應該擁有95％的置信區間，如下所示：
y = np.array(df.Close[-750:]).reshape(-1, 1)under_upper = upper_bound_unscaled > yover_lower = lower_bound_unscaled < ytotal = (under_upper == over_lower)print("{} our predictions are in our confidence interval".format(np.mean(total)))
檢查輸出圖形
現在，我們將把預測結果繪製為可視化圖形來檢查我們的網絡是否運行的很順利，我們將在置信區間內繪製真實值與預測值。
params = {"ytick.color" : "w","xtick.color" : "w", "axes.labelcolor" : "w", "axes.edgecolor" : "w"}plt.rcParams.update(params)plt.title("IBM Stock prices", color="white")plt.plot(df_pred.index, df_pred.Close, color='black', label="Real")plt.plot(idx_pred, pred_mean_unscaled, label="Prediction for {} days, than consult".format(future_length), color="red")plt.fill_between(x=idx_pred, y1=upper_bound_unscaled[:,0], y2=lower_bound_unscaled[:,0], facecolor='green', label="Confidence interval", alpha=0.5)plt.legend()
最後，我們放大一下著重看看預測部分。
params = {"ytick.color" : "w","xtick.color" : "w", "axes.labelcolor" : "w", "axes.edgecolor" : "w"}plt.rcParams.update(params)plt.title("IBM Stock prices", color="white")plt.fill_between(x=idx_pred, y1=upper_bound_unscaled[:,0], y2=lower_bound_unscaled[:,0], facecolor='green', label="Confidence interval", alpha=0.75)plt.plot(idx_pred, df_pred.Close[-len(pred_mean_unscaled):], label="Real", alpha=1, color='black', linewidth=0.5)plt.plot(idx_pred, pred_mean_unscaled, label="Prediction for {} days, than consult".format(future_length), color="red", alpha=0.5)plt.legend()
總結
我們看到BLiTZ內置的貝葉斯LSTM使得貝葉斯深度學習的所有功能都變得非常簡單，並且可以順利地在時間序列上進行迭代。 我們還看到，貝葉斯LSTM已與Torch很好地集成在一起，並且易於使用，你可以在任何工作或研究中使用它。
我們還可以非常準確地預測IBM股票價格的置信區間，而且這比一般的點估計可能要有用的多。
作者:Piero Esposito
deephub翻譯組：Alenander Zhao