貝葉斯神經網絡

這篇文章是對貝葉斯神經網絡的一個初步認識，可以理解下。

[隨着人工智能技術的普及，用機器學習預測市場價格波動的方法最近層出不窮。本文中，Alex Honchar 介紹了利用概率編程和 Pyro 進行價格預測的方法，相較於常規神經網絡，新方法對於數據的依賴程度更小，結果更准確。在實驗中，作者選擇了最近流行的虛擬貨幣「以太幣」作為實例進行價格預測。]

去年我曾發表過幾篇有關使用神經網絡進行金融價格預測的教程，我認為其中有一部分結果至少還挺有意思，並且值得在實際交易中加以應用。如果你閱讀過這些文章，你一定注意到一個現象：當你試圖將一些機器學習模型應用於「隨機」數據並希望從中找到隱藏規律的時候，訓練過程往往會產生嚴重的過擬合。我們曾使用不同的正則化技術和附加數據應對這個問題，但是這不僅很費時，還有種盲目搜索的感覺。

今天，我想介紹一個略微有些不同的方法對同樣的算法進行擬合。使用概率的觀點看待這個問題能夠讓我們從數據本身學習正則化、估計預測結果的確定性、使用更少的數據進行訓練，還能在模型中引入額外的概率依賴關系。我不會過多深入貝葉斯模型或變分原理的數學、技術細節，而是會給出一些概述，也更多地將討論集中在應用場景當中。文中所用的代碼可以在以下鏈接中找到：https://github.com/Rachnog/Deep-Trading/tree/master/bayesian

與此同時，我也推薦大家查閱我此前發布的基於神經網絡的財務預測教程：

1. 簡單時間序列預測（錯誤糾正完畢）

2. 正確一維時間序列預測+回測

3. 多元時間序列預測

4. 波動預測和自定義損失

5. 多任務和多模式學習

6. 超參數優化

7. 用神經網絡進行經典策略強化

為了更深入地了解概率規划、貝葉斯模型以及它們的應用，我推薦你在以下資源網站中查看：

• 模式識別和機器學習

• 黑客貝葉斯方法

• 下面即將提到的庫文件

另外，你還可能會用到下列 Python 庫：

• PyMC3 (https://github.com/pymc-devs/pymc3)

• Edward (http://edwardlib.org/)

• Pyro (http://pyro.ai/)

概率編程

這個「概率」指的是什么？我們為什么稱其為「編程」呢？首先，讓我們回憶一下我們所謂「正常的」神經網絡指的是什么、以及我們能從中得到什么。神經網絡有着以矩陣形式表達的參數（權重），而其輸出通常是一些標量或者向量（例如在分類問題的情況下）。當我們用諸如 SGD 的方法訓練這個模型后，這些矩陣會獲得固定值。與此同時，對於同一個輸入樣本，輸出向量應該相同，就是這樣！但是，如果我們將所有的參數和輸出視為相互依賴的分布，會發生什么？神經網絡的權重將與輸出一樣，是一個來自網絡並取決於參數的樣本——如果是這樣，它能為我們帶來什么？

讓我們從基礎講起。如果我們認為網絡是一個取決於其他分布的數集，這首先就構成了聯合概率分布 p(y, z|x)，其中有着輸出 y 和一些模型 z 的「內部」隱變量，它們都取決於輸入 x（這與常規的神經網絡完全相同）。我們感興趣的是找到這樣神經網絡的分布，這樣一來就可以對 y ~ p(y|x) 進行采樣，並獲得一個形式為分布的輸出，該分布中抽取的樣本的期望通常是輸出，和標准差（對不確定性的估計）——尾部越大，則輸出置信度越小。

這種設定可能不是很明確，但我們只需要記住：現在開始，模型中所有的參數、輸入及輸出都是分布，並且在訓練時對這些分布進行擬合，以便在實際應用中獲得更高的准確率。我們也需要注意自己設定的參數分布的形狀（例如，所有的初識權重 w 服從正態分布 Normal（0,1），之后我們將學習正確的均值和方差）。初始分布即所謂的先驗知識，在訓練集上訓練過的分布即為后驗知識。我們使用后者進行抽樣並得出結果。

圖源：http://www.indiana.edu/~kruschke/BMLR/

模型要擬合到什么程度才有用？通用結構被稱為變分推理（variational inference）。無需細想，我們可以假設，我們希望找到一個可以得到最大對數似然函數 p_w（z | x）的模型，其中 w 是模型的參數（分布參數），z 是我們的隱變量（隱藏層的神經元輸出，從參數 w 的分布采樣得到），x 是輸入數據樣本。這就是我們的模型了。我們在 Pyro 中引入了一個實例來介紹這個模型，該簡單實例包含所有隱變量 q_（z）的一些分布，其中 ф 被稱為變分參數。這種分布必須近似於訓練最好的模型參數的「實際」分布。

訓練目標是使得 [log(p_w(z|x))—log(q_ф(z))] 的期望值相對於有指導的輸入數據和樣本最小化。在這里我們不探討訓練的細節，因為這里面的知識量太大了，此處就先當它是一個可以優化的黑箱吧。

對了，為什么需要編程呢？因為我們通常將這種概率模型（如神經網絡）定義為變量相互關聯的有向圖，這樣我們就可以直接顯示變量間的依賴關系：

圖源：http://kentonmurray.com/

而且，概率編程語言起初就被用於定義此類模型並在模型上做推理。

為什么選擇概率編程？

不同於在模型中使用 dropout 或 L1 正則化，你可以把它當作你數據中的隱變量。考慮到所有的權重其實是分布，你可以從中抽樣 N 次得到輸出的分布，通過計算該分布的標准差，你就知道能模型有多靠譜。作為成果，我們可以只用少量的數據來訓練這些模型，而且我們可以靈活地在變量之間添加不同的依賴關系。

概率編程的不足

我還沒有太多關於貝葉斯建模的經驗，但是我從 Pyro 和 PyMC3 中了解到，這類模型的訓練過程十分漫長且很難定義正確的先驗分布。而且，處理從分布中抽取的樣本會導致誤解和歧義。

數據准備

我已經從 http://bitinfocharts.com/ 抓取了每日 Ethereum（以太坊）的價格數據。其中包括典型的 OHLCV（高開低走），另外還有關於 Ethereum 的每日推特量。我們將使用七日的價格、開盤及推特量數據來預測次日的價格變動情況。

價格、推特數、大盤變化

上圖是一些數據樣本——藍線對應價格變化，黃線對應推特數變化，綠色對應大盤變化。它們之間存在某種正相關（0.1—0.2）。因此我們希望能利用好這些數據中的模式對模型進行訓練。

貝葉斯線性回歸

首先，我想驗證簡單線性分類器在任務中的表現結果（並且我想直接使用 Pyro tutorial——http://pyro.ai/examples/bayesian_regression.html——的結果）。我們按照以下操作在 PyTorch 上定義我們的模型（詳情參閱官方指南：http://pyro.ai/examples/bayesian_regression.html）。

1. class RegressionModel(nn.Module):

2.    def __init__(self, p):

3.        super(RegressionModel, self).__init__()

4.        self.linear = nn.Linear(p, 1)

5. def forward(self, x):

6.        # x * w + b

7.        return self.linear(x)

以上是我們以前用過的簡單確定性模型，下面是用 Pyro 定義的概率模型：

1. def model(data):

2.    # Create unit normal priors over the parameters

3.    mu = Variable(torch.zeros(1, p)).type_as(data)

4.    sigma = Variable(torch.ones(1, p)).type_as(data)

5.    bias_mu = Variable(torch.zeros(1)).type_as(data)

6.    bias_sigma = Variable(torch.ones(1)).type_as(data)

7.    w_prior, b_prior = Normal(mu, sigma), Normal(bias_mu, bias_sigma)

8.    priors = {'linear.weight': w_prior, 'linear.bias': b_prior}

9.    lifted_module = pyro.random_module("module", regression_model, priors)

10.    lifted_reg_model = lifted_module()

11. with pyro.iarange("map", N, subsample=data):

12.        x_data = data[:, :-1]

13.        y_data = data[:, -1]

14.        # run the regressor forward conditioned on inputs

15.        prediction_mean = lifted_reg_model(x_data).squeeze()

16.        pyro.sample("obs",

17.                    Normal(prediction_mean, Variable(torch.ones(data.size(0))).type_as(data)),

18.                    obs=y_data.squeeze())

從上面的代碼可知，參數 W 和 b 均定義為一般線性回歸模型分布，兩者都服從正態分布 Normal（0,1）。我們稱之為先驗，創建 Pyro 的隨機函數（在我們的例子中是 PyTorch 中的 RegressionModel），為它添加先驗 ({『linear.weight』: w_prior, 『linear.bias』: b_prior})，並根據輸入數據 x 從這個模型 p(y|x) 中抽樣。

這個模型的 guide 部分可能像下面這樣：

1. def guide(data):

2.    w_mu = Variable(torch.randn(1, p).type_as(data.data), requires_grad=True)

3.    w_log_sig = Variable(0.1 * torch.ones(1, p).type_as(data.data), requires_grad=True)

4.    b_mu = Variable(torch.randn(1).type_as(data.data), requires_grad=True)

5.    b_log_sig = Variable(0.1 * torch.ones(1).type_as(data.data), requires_grad=True)

6.    mw_param = pyro.param("guide_mean_weight", w_mu)

7.    sw_param = softplus(pyro.param("guide_log_sigma_weight", w_log_sig))

8.    mb_param = pyro.param("guide_mean_bias", b_mu)

9.    sb_param = softplus(pyro.param("guide_log_sigma_bias", b_log_sig))

10.    w_dist = Normal(mw_param, sw_param)

11.    b_dist = Normal(mb_param, sb_param)

12.    dists = {'linear.weight': w_dist, 'linear.bias': b_dist}

13.    lifted_module = pyro.random_module("module", regression_model, dists)

14.    return lifted_module()

我們定義了想要「訓練」的分布的可變分布。如你所見，我們為 W 和 b 定義了相同的分布，目的是讓它們更接近實際情況（據我們假設）。這個例子中，我讓分布圖更窄一些（服從正態分布 Normal(0, 0.1)）

然后，我們用這種方式對模型進行訓練：

1.    for j in range(3000):

2.    epoch_loss = 0.0

3.    perm = torch.randperm(N)

4.    # shuffle data

5.    data = data[perm]

6.    # get indices of each batch

7.    all_batches = get_batch_indices(N, 64)

8.    for ix, batch_start in enumerate(all_batches[:-1]):

9.        batch_end = all_batches[ix + 1]

10.        batch_data = data[batch_start: batch_end]

11.        epoch_loss += svi.step(batch_data)

在模型擬合后，我們想從中抽樣出 y。我們循環 100 次並計算每一步的預測值的均值和標准差（標准差越高，預測置信度就越低）。

1.   preds = []

2. for i in range(100):

3.    sampled_reg_model = guide(X_test)

4.    pred = sampled_reg_model(X_test).data.numpy().flatten()

5.    preds.append(pred)

現在有很多經典的經濟預測度量方法，例如 MSE、MAE 或 MAPE，它們都可能會讓人困惑——錯誤率低並不意味着你的模型表現得好，驗證它在測試集上的表現也十分重要，而這就是我們做的工作。

使用貝葉斯模型進行為期 30 天的預測

從圖中我們可以看到，預測效果並不夠好。但是預測圖中最后的幾個跳變的形狀很不錯，這給了我們一線希望。繼續加油！

常規神經網絡

在這個非常簡單的模型進行實驗后，我們想要嘗試一些更有趣的神經網絡。首先讓我們利用 25 個帶有線性激活的神經元的單隱層網絡訓練一個簡單 MLP：

1. def get_model(input_size):

2.    main_input = Input(shape=(input_size, ), name='main_input')

3.    x = Dense(25, activation='linear')(main_input)

4.    output = Dense(1, activation = "linear", name = "out")(x)

5.    final_model = Model(inputs=[main_input], outputs=[output])

6.    final_model.compile(optimizer='adam',  loss='mse')

7.    return final_model

訓練 100 個 epoch：

1. model = get_model(len(X_train[0]))

2. history = model.fit(X_train, Y_train,

3.              epochs = 100,

4.              batch_size = 64,

5.              verbose=1,

6.              validation_data=(X_test, Y_test),

7.              callbacks=[reduce_lr, checkpointer],

8.              shuffle=True)

其結果如下：

使用 Keras 神經網絡進行為期 30 天的預測

我覺得這比簡單的貝葉斯回歸效果更差，此外這個模型不能得到確定性的估計，更重要的是，這個模型甚至沒有正則化。

貝葉斯神經網絡

現在我們用 PyTorch 來定義上文在 Keras 上訓練的模型：

1. class Net(torch.nn.Module):

2.    def __init__(self, n_feature, n_hidden):

3.        super(Net, self).__init__()

4.        self.hidden = torch.nn.Linear(n_feature, n_hidden)   # hidden layer

5.        self.predict = torch.nn.Linear(n_hidden, 1)   # output layer

6. def forward(self, x):

7.        x = self.hidden(x)

8.        x = self.predict(x)

9.        return x

相比於貝葉斯回歸模型，我們現在有兩個參數集（從輸入層到隱藏層的參數和隱藏層到輸出層的參數），所以我們需要對分布和先驗知識稍加改動，以適應我們的模型：

1. priors = {'hidden.weight': w_prior,

2.              'hidden.bias': b_prior,

3.              'predict.weight': w_prior2,

4.              'predict.bias': b_prior2}

以及 guide 部分：

1. dists = {'hidden.weight': w_dist,

2.              'hidden.bias': b_dist,

3.              'predict.weight': w_dist2,

4.              'predict.bias': b_dist2}

請不要忘記為模型中的每一個分布起一個不同的名字，因為模型中不應存在任何歧義和重復。更多代碼細節請參見源代碼：https://github.com/Rachnog/Deep-Trading/tree/master/bayesian

訓練之后，讓我們看看最后的結果：

使用 Pyro 神經網絡進行為期 30 天的預測

它看起來比之前的結果都好得多！

比起常規貝葉斯模型，考慮到貝葉斯模型所中習得的權重特征或正則化，我還希望看到權重的數據。我按照以下方法查看 Pyro 模型的參數：

1. for name in pyro.get_param_store().get_all_param_names():

2.    print name, pyro.param(name).data.numpy()

這是我在 Keras 模型中所寫的代碼：

1. import tensorflow as tf

2. sess = tf.Session()

3. with sess.as_default():

4.    tf.global_variables_initializer().run()

5. dense_weights, out_weights = None, None

6. with sess.as_default():

7.    for layer in model.layers:

8.        if len(layer.weights) > 0:

9.            weights = layer.get_weights()

10.            if 'dense' in layer.name:

11.                dense_weights = layer.weights[0].eval()

12.            if 'out' in layer.name:

13.                out_weights = layer.weights[0].eval()

例如，Keras 模型最后一層的權重的均值和標准差分別為 -0.0025901748 和 0.30395043，Pyro 模型對應值為 0.0005974418 和 0.0005974418。數字小了很多，但效果真的不錯！其實這就是 L2 或 Dropout 這種正則化算法要做的——把參數逼近到零，而我們可以用變分推理來實現它！隱藏層的權重變化更有趣。我們將一些權重向量繪制成圖，藍線是 Keras 模型的權重，橙線是 Pyro 模型的權重：

輸入層與隱藏層之間的部分權重

真正有意思的不止是權重的均值與標准差變得小，還有一點是權重變得稀疏，所以基本上在訓練中完成了第一個權重集的稀疏表示，以及第二個權重集的 L2 正則化，多么神奇！別忘了自己跑跑代碼感受一下：https://github.com/Rachnog/Deep-Trading/tree/master/bayesian

小結

我們在文中使用了新穎的方法對神經網絡進行訓練。不同於順序更新靜態權重，我們是更新的是權重的分布。因此，我們可能獲得有趣又有用的結果。我想強調的是，貝葉斯方法讓我們在調整神經網絡時不需要手動添加正則化，了解模型的不確定性，並盡可能使用更少的數據來獲得更好的結果。感謝閱讀：) 

原文鏈接：https://medium.com/@alexrachnog/financial-forecasting-with-probabilistic-programming-and-pyro-db68ab1a1dba