机器学习贝叶斯神经网络与Python

机器学习的当前趋势

机器学习目前有三大趋势：概率编程，深度学习和“大数据”。在PP内部，许多创新都在使用变分推理来进行规模扩展。在本文中，我将展示如何在PyMC3中使用变分推理来拟合简单的贝叶斯神经网络。我还将讨论如何在概率编程和深度学习之间架起桥梁，为将来的研究开辟非常有趣的途径。

概率性编程规划

概率编程允许非常灵活地创建自定义概率模型，主要关注从数据中洞察和学习。该方法本质上是贝叶斯方法，因此我们可以指定先验来告知和约束我们的模型，并以后验分布的形式得到不确定性估计。使用MCMC采样算法，我们可以从这个后验中抽取样本，非常灵活地估计这些模型。PyMC3和Stan是目前用于构建和评估这些模型的最先进工具。然而，采样的一个主要缺点是它通常非常慢，特别是对于高维模型。这就是为什么最近开发出的变分推理算法几乎和MCMC一样灵活，但速度要快得多。这些算法不是从后验中抽取样本，而是将后验的分布(如正态分布)拟合到样本上，将样本问题转化为优化问题。

不幸的是，当涉及到传统的ML问题，如分类或（非线性）回归时，概率编程常常处于第二位（在准确性和可扩展性方面），而不是像集成学习这样的算法方法(例如随机森林或梯度增强回归树) 。

深度学习

现在正处于第三次文艺复兴时期，深度学习一直成为头条新闻，它主宰了几乎所有的物体识别基准，在雅达利游戏（Atari games）中大出风头，在围棋(Go)中击败了世界冠军李世石（Lee Sedol）。从统计学角度来看，神经网络是非常好的非线性函数近似和表示学习者表现。虽然主要用于分类，但它们已经扩展到使用自编码器进行无监督学习以及各种其他有趣的方式（例如，递归网络或MDN来估计多模态分布）。为什么他们这么好用？没有人真正知道，因为统计特性仍未被完全理解。

深度学习的一大创新之处在于能够训练这些极其复杂的模型。这取决于几个支柱：

• 速度：促进GPU允许更快的处理。
• 软件：Theano和TensorFlow等框架允许灵活地创建抽象模型，然后可以对其进行优化并编译到CPU或GPU。
• 学习算法：对数据子集的训练 - 随机梯度下降 - 允许我们在大量数据上训练这些模型。dropout等技术可以避免过度拟合。
• 架构性：很多创新来自于改变输入层，如卷积神经网络，或输出层，如MDN。

弥合深度学习和概率编程

一方面，我们有概率编程，它允许我们以非常有原则和易于理解的方式构建相当小的模型，以深入了解我们的数据； 另一方面，我们有深度学习，它使用许多启发式方法来训练巨大且高度复杂的模型，这些模型在预测方面令人惊叹。变分推理的最新创新允许概率编程来扩展模型复杂性以及数据大小。因此，我们正处于能够将这两种方法结合起来以帮助解开机器学习的新创新的风口浪尖。

虽然这将允许概率编程应用于更广泛的有趣问题，但我相信，这种桥梁也为深度学习的创新带来了巨大的希望。

• 预测中的不确定性：正如我们将在下面看到的，贝叶斯神经网络告诉我们其预测的不确定性。我认为不确定性是机器学习中一个未被充分认识的概念，因为它显然对现实世界的应用很重要。但它在训练中也很有用。例如，我们可以针对最不确定的样本专门训练模型。
• 表示中的不确定性：我们还得到了我们权重的不确定性估计，这可以告诉我们网络的学习表示的稳定性。
• 先验的正则化：权重通常是L2正则化以避免过度拟合，这非常自然地成为权重系数的高斯先验。然而，我们可以想象各种其他先验，如spike-and-slab来加强稀疏性（这更像是使用L1-norm）。
• 通过知情先验迁移学习：如果我们想训练一个网络使用新的对象识别数据集，我们可以通过将知情先验集中在从其他预训练网络检索到的权重上来引导学习。
• 分层神经网络：在概率编程中一种非常强大的方法是分层建模，它允许将子群上学到的东西汇集到整个总体中。应用于神经网络，在分层数据集中，我们可以训练单个神经网络专门研究子群，同时仍然了解整体的表示。例如，想象一个经过训练的网络，可以根据汽车的图片对汽车模型进行分类。我们可以训练分层神经网络，其中训练子神经网络以仅从单个制造商分辨模型。人们的直觉是，来自某一制造商的所有汽车都有一定的相似性，因此训练专门针对品牌的个别网络是有意义的。但是，由于各个网络连接在更高层，他们仍会与其他专门的子网络分享有关对所有品牌有用的特征的信息。有趣的是，网络的不同层可以通过层级的各个级别来通知。
• 其他混合架构：我们可以更自由地构建各种神经网络。例如，可以使用贝叶斯非参数来灵活地调整隐藏层的大小和形状，以在训练期间最优地将网络架构扩展到当前的问题。目前，这需要更多的超参数优化和大量的知识。

PyMC3中的贝叶斯神经网络

生成数据

首先，让我们生成一些玩具数据 - 一个简单的二元分类问题，它不是线性可分的。

%matplotlib inline
import theano
import pymc3 as pm
import sklearn
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from warnings import filterwarnings
filterwarnings('ignore')
sns.set_style('white')
from sklearn import datasets
from sklearn.preprocessing import scale
from sklearn.model_selection import train_test_split
from sklearn.datasets import make_moons
X, Y = make_moons(noise=0.2, random_state=0, n_samples=1000)
X = scale(X)
X = X.astype(float)
Y = Y.astype(float)
X_train, X_test, Y_train, Y_test = train_test_split(X, Y, test_size=.5)
fig, ax = plt.subplots(figsize=(12, 8))
ax.scatter(X[Y==0, 0], X[Y==0, 1], label='Class 0')
ax.scatter(X[Y==1, 0], X[Y==1, 1], color='r', label='Class 1')
sns.despine(); ax.legend()
ax.set(xlabel='X', ylabel='Y', title='Toy binary classification data set');

模型规格

神经网络非常简单。基本单元是一个感知器，它只不过是逻辑回归。我们并行使用其中许多，然后将它们叠加起来以获得隐藏层。在这里，我们将使用2个隐藏层，每个隐藏层有5个神经元，足以解决这个简单问题。

def construct_nn(ann_input, ann_output):
 n_hidden = 5
 
 # Initialize random weights between each layer
 init_1 = np.random.randn(X.shape[1], n_hidden).astype(float)
 init_2 = np.random.randn(n_hidden, n_hidden).astype(float)
 init_out = np.random.randn(n_hidden).astype(float)
 
 with pm.Model() as neural_network:
 # Weights from input to hidden layer
 weights_in_1 = pm.Normal('w_in_1', 0, sd=1, 
 shape=(X.shape[1], n_hidden), 
 testval=init_1)
 
 # Weights from 1st to 2nd layer
 weights_1_2 = pm.Normal('w_1_2', 0, sd=1, 
 shape=(n_hidden, n_hidden), 
 testval=init_2)
 
 # Weights from hidden layer to output
 weights_2_out = pm.Normal('w_2_out', 0, sd=1, 
 shape=(n_hidden,), 
 testval=init_out)
 
 # Build neural-network using tanh activation function
 act_1 = pm.math.tanh(pm.math.dot(ann_input, 
 weights_in_1))
 act_2 = pm.math.tanh(pm.math.dot(act_1, 
 weights_1_2))
 act_out = pm.math.sigmoid(pm.math.dot(act_2, 
 weights_2_out))
 
 # Binary classification -> Bernoulli likelihood
 out = pm.Bernoulli('out', 
 act_out,
 observed=ann_output,
 total_size=Y_train.shape[0] # IMPORTANT for minibatches
 )
 return neural_network
# Trick: Turn inputs and outputs into shared variables. 
ann_input = theano.shared(X_train)
ann_output = theano.shared(Y_train)
neural_network = construct_nn(ann_input, ann_output)

Normal先验有助于调整权重。通常我们会在输入中添加一个常量b，但是为了保持代码的整洁，我在这里省略了它。

变分推理：缩放模型复杂度

我们现在可以运行像NUTS这样的MCMC采样器，在这种情况下工作得非常好，但正如我已经提到的，当我们将模型扩展到具有更多层的更深层架构时，这将变得非常缓慢。

相反，我们将使用ADVI变分推理算法，这是最近添加到PyMC3，并更新使用运算符变分推理（OPVI）框架。这样做速度更快，并且可以更好地扩展。

from pymc3.theanof import set_tt_rng, MRG_RandomStreams
set_tt_rng(MRG_RandomStreams(42))
%%time
with neural_network:
 inference = pm.ADVI()
 approx = pm.fit(n=50000, method=inference)

在性能方面，考虑到NUTS确实不是太好。接下来我们使它更快一些。为了让它真正运行起来，我们可能需要在GPU上运行神经网络。

由于使用样本更方便，我们可以使用sample方法快速地从变分近似中提取样本（这只是从正态分布中采样，所以与MCMC完全不同）：

trace = approx.sample(draws=5000)

绘制目标函数（ELBO）我们可以看到优化会逐渐改善拟合。

plt.plot(-inference.hist)
plt.ylabel('ELBO')
plt.xlabel('iteration');

现在我们训练了模型，让我们使用后验预测检查（PPC）来预测保留集。

# Replace arrays our NN references with the test data
ann_input.set_value(X_test)
ann_output.set_value(Y_test)
with neural_network:
 ppc = pm.sample_ppc(trace, samples=500, progressbar=False)
# Use probability of > 0.5 to assume prediction of class 1
pred = ppc['out'].mean(axis=0) > 0.5

让我们来看看我们的预测：

fig, ax = plt.subplots()
ax.scatter(X_test[pred==0, 0], X_test[pred==0, 1])
ax.scatter(X_test[pred==1, 0], X_test[pred==1, 1], color='r')
sns.despine()
ax.set(title='Predicted labels in testing set', xlabel='X', ylabel='Y');
print('Accuracy = {}%'.format((Y_test == pred).mean() * 100))

Accuracy = 95.0%

我们的神经网络已经做得很好了。

让我们看看分类器学到了什么

这样，我们在整个输入空间的网格上评估类概率预测。

grid = pm.floatX(np.mgrid[-3:3:100j,-3:3:100j])
grid_2d = grid.reshape(2, -1).T
dummy_out = np.ones(grid.shape[1], dtype=np.int8)
ann_input.set_value(grid_2d)
ann_output.set_value(dummy_out)
with neural_network:
 ppc = pm.sample_ppc(trace, samples=500, progressbar=False)

概率面

cmap = sns.diverging_palette(250, 12, s=85, l=25, as_cmap=True)
fig, ax = plt.subplots(figsize=(14, 8))
contour = ax.contourf(grid[0], grid[1], ppc['out'].mean(axis=0).reshape(100, 100), cmap=cmap)
ax.scatter(X_test[pred==0, 0], X_test[pred==0, 1])
ax.scatter(X_test[pred==1, 0], X_test[pred==1, 1], color='r')
cbar = plt.colorbar(contour, ax=ax)
_ = ax.set(xlim=(-3, 3), ylim=(-3, 3), xlabel='X', ylabel='Y');
cbar.ax.set_ylabel('Posterior predictive mean probability of class label = 0');

预测值的不确定性

到目前为止，我展示的一切都可以用非贝叶斯神经网络完成。每个类标签的后验预测均值应与最大似然预测值相同。不过，我们也可以通过后验预测的标准差来了解我们预测中的不确定性。这是看起来像：

cmap = sns.cubehelix_palette(light=1, as_cmap=True)
fig, ax = plt.subplots(figsize=(14, 8))
contour = ax.contourf(grid[0], grid[1], ppc['out'].std(axis=0).reshape(100, 100), cmap=cmap)
ax.scatter(X_test[pred==0, 0], X_test[pred==0, 1])
ax.scatter(X_test[pred==1, 0], X_test[pred==1, 1], color='r')
cbar = plt.colorbar(contour, ax=ax)
_ = ax.set(xlim=(-3, 3), ylim=(-3, 3), xlabel='X', ylabel='Y');
cbar.ax.set_ylabel('Uncertainty (posterior predictive standard deviation)');

我们可以看到非常接近决策边界的地方，我们对于预测哪个标签的不确定性是最高的。您可以想象，将预测与不确定性联系起来是许多应用（如医疗保健）的关键属性。为了进一步提高精度，我们可能希望主要根据来自高不确定区域的样本来进行模型训练。

Mini-batch ADVI

到目前为止，我们已经对所有数据进行了模型训练。显然，这不会像ImageNet那样扩展。此外，对mini-batch数据（随机梯度下降）的训练避免了局部最小值，并且可以加快收敛速度。

幸运的是，ADVI也可以在mini-batch上运行。它只需要一些设置：

minibatch_x = pm.Minibatch(X_train, batch_size=32)
minibatch_y = pm.Minibatch(Y_train, batch_size=32)
neural_network_minibatch = construct_nn(minibatch_x, minibatch_y)
with neural_network_minibatch:
 inference = pm.ADVI()
 approx = pm.fit(40000, method=inference)
plt.plot(-inference.hist)
plt.ylabel('ELBO')
plt.xlabel('iteration');

如您所见，mini-batch ADVI的运行时间要低得多。它似乎收敛的更快。

总结

希望本文演示了PyMC3: ADVI中可用的一个非常强大的新推理算法。我还认为，如前所述，弥合概率编程和深度学习之间的差距可以为这一领域的创新开辟许多新途径。