使用Pytorch和Pyro实现贝叶斯神经网络(Bayesian Neural Network)_pytorch贝叶斯网络分类_gang_akarui的博客

本文翻译自： https://www.noconote.work/entry/2019/01/01/160100
最近概率模型和神经网络相结合的研究变得多了起来，这次使用Uber开源的Pyro来实现一个贝叶斯神经网络。概率编程框架最近出了不少，Uber的Pyro基于Pytorch，Google的Edward基于TensorFlow，还有一些独立的像PyMC3,Stan,Pomegranate等等。
这次选择使用Pyro实现的原因是，Pyro基于Numpy实现，加上对Pytorch的支持很友好，可以自动计算梯度，动态计算，这些好处当然都是Pytorch带来的。官网链接： https://pyro.ai/，介绍到此为止，接下来开始一步一步写代码了。
Step1：导入需要的模块
import numpy as np
import torch
import torch.nn as nn
import torch.nn.functional as F
from torchvision import datasets
from torchvision import transforms
import pyro
from pyro.distributions import Normal
from pyro.distributions import Categorical
from pyro.optim import Adam
from pyro.infer import SVI
from pyro.infer import Trace_ELBO
import matplotlib.pyplot as plt
import seaborn as sns
sns.set()
Step2：准备数据集，这次使用的是MNIST 
train_loader = torch.utils.data.DataLoader(
        datasets.MNIST("/tmp/mnist", train=True, download=True,
                       transform=transforms.Compose([transforms.ToTensor(),])),
        batch_size=128, shuffle=True)
test_loader = torch.utils.data.DataLoader(
        datasets.MNIST("/tmp/mnist", train=False, transform=transforms.Compose([transforms.ToTensor(),])
Step3：定义神经网络模型 
class BNN(nn.Module):
    def __init__(self, input_size, hidden_size, output_size):
        super(BNN, self).__init__()
        self.fc1 = nn.Linear(input_size, hidden_size)
        self.out = nn.Linear(hidden_size, output_size)
    def forward(self, x):
        output = self.fc1(x)
        output = F.relu(output)
        output = self.out(output)
        return output
net = BNN(28*28, 1024, 10)
Step4：基于正态分布，初始化weights和bias 
def model(x_data, y_data):
    # define prior destributions
    fc1w_prior = Normal(loc=torch.zeros_like(net.fc1.weight),
                                      scale=torch.ones_like(net.fc1.weight))
    fc1b_prior = Normal(loc=torch.zeros_like(net.fc1.bias),
                                      scale=torch.ones_like(net.fc1.bias))
    outw_prior = Normal(loc=torch.zeros_like(net.out.weight),
                                       scale=torch.ones_like(net.out.weight))
    outb_prior = Normal(loc=torch.zeros_like(net.out.bias),
                                       scale=torch.ones_like(net.out.bias))
    priors = {
        'fc1.weight': fc1w_prior,
        'fc1.bias': fc1b_prior, 
        'out.weight': outw_prior,
        'out.bias': outb_prior}
    lifted_module = pyro.random_module("module", net, priors)
    lifted_reg_model = lifted_module()
    lhat = F.log_softmax(lifted_reg_model(x_data))
    pyro.sample("obs", Categorical(logits=lhat), obs=y_data)
Step5：为了近似后验概率分布，先定义一个函数 
def guide(x_data, y_data):
    fc1w_mu = torch.randn_like(net.fc1.weight)
    fc1w_sigma = torch.randn_like(net.fc1.weight)
    fc1w_mu_param = pyro.param("fc1w_mu", fc1w_mu)
    fc1w_sigma_param = F.softplus(pyro.param("fc1w_sigma", fc1w_sigma))
    fc1w_prior = Normal(loc=fc1w_mu_param, scale=fc1w_sigma_param)
    fc1b_mu = torch.randn_like(net.fc1.bias)
    fc1b_sigma = torch.randn_like(net.fc1.bias)
    fc1b_mu_param = pyro.param("fc1b_mu", fc1b_mu)
    fc1b_sigma_param = F.softplus(pyro.param("fc1b_sigma", fc1b_sigma))
    fc1b_prior = Normal(loc=fc1b_mu_param, scale=fc1b_sigma_param)
    outw_mu = torch.randn_like(net.out.weight)
    outw_sigma = torch.randn_like(net.out.weight)
    outw_mu_param = pyro.param("outw_mu", outw_mu)
    outw_sigma_param = F.softplus(pyro.param("outw_sigma", outw_sigma))
    outw_prior = Normal(loc=outw_mu_param, scale=outw_sigma_param).independent(1)
    outb_mu = torch.randn_like(net.out.bias)
    outb_sigma = torch.randn_like(net.out.bias)
    outb_mu_param = pyro.param("outb_mu", outb_mu)
    outb_sigma_param = F.softplus(pyro.param("outb_sigma", outb_sigma))
    outb_prior = Normal(loc=outb_mu_param, scale=outb_sigma_param)
    priors = {
        'fc1.weight': fc1w_prior,
        'fc1.bias': fc1b_prior,
        'out.weight': outw_prior,
        'out.bias': outb_prior}
    lifted_module = pyro.random_module("module", net, priors)
    return lifted_module()
Step6：定义Optimizer
 *这里多说两句，优化器还是喜闻乐见的Adam；这个SVI函数是干嘛的呢，它使用变分推理（Variational inference）来逼近后验概率分布，说到变分推理自然离不开ELBO，这里的Trace_ELBO函数是继承的父类ELBO，还有其他种类的ELBO，比如JitTrace_ELBO，TraceGraph_ELBO等等，具体选哪个，去看文档就好了。顺便推荐一个变分推理的教程，徐亦达老师的视频教程非常棒。 
optim = Adam({"lr": 0.01})
svi = SVI(model, guide, optim, loss=Trace_ELBO())
Step7：训练Bayesian Neural Network 
n_iterations = 5
loss = 0
for j in range(n_iterations):
    loss = 0
    for batch_id, data in enumerate(train_loader):
        loss += svi.step(data[0].view(-1,28*28), data[1])
    normalizer_train = len(train_loader.dataset)
    total_epoch_loss_train = loss / normalizer_train
    print("Epoch ", j, " Loss ", total_epoch_loss_train)
Step8：使用训练好的模型，预测结果
 *Bayesian Neural Network要分别针对每个类别进行采样，然后计算预测值中概率最大的类别，作为结果。这里有10个数字，所以我们设置采样的次数为10 
n_samples = 10
def predict(x):
    sampled_models = [guide(None, None) for _ in range(n_samples)]
    yhats = [model(x).data for model in sampled_models]
    mean = torch.mean(torch.stack(yhats), 0)
    return np.argmax(mean.numpy(), axis=1)
print('Prediction when network is forced to predict')
correct = 0
total = 0
for j, data in enumerate(test_loader):
    images, labels = data
    predicted = predict(images.view(-1,28*28))
    total += labels.size(0)
    correct += (predicted == labels).sum().item()
print("accuracy: %d %%" % (100 * correct / total))
*当然也可以把计算出的概率归一化一下，这样使得结果更具有解释性 
def predict_prob(x):
    sampled_models = [guide(None, None) for _ in range(n_samples)]
    yhats = [model(x).data for model in sampled_models]
    mean = torch.mean(torch.stack(yhats), 0)
    return mean
def normalize(x):
    return (x - x.min()) / (x.max() - x.min())
def plot(x, yhats):
    fig, (axL, axR) = plt.subplots(ncols=2, figsize=(8, 4))
    axL.bar(x=[i for i in range(10)], height= F.softmax(torch.Tensor(normalize(yhats.numpy()))[0]))
    axL.set_xticks([i for i in range(10)], [i for i in range(10)])
    axR.imshow(x.numpy()[0])
    plt.show()
效果如下： 
x, y = test_loader.dataset[0]
yhats = predict_prob(x.view(-1, 28*28))
print("ground truth: ", y.item())
print("predicted: ", yhats.numpy())
plot(x, yhats)
ground truth:  7
predicted:  [[  39.87391   -86.72453   112.25181   102.81404   -89.95473    52.947186
  -167.88455   402.66757   -29.799454   95.75331 ]]
                    本文翻译自：https://www.noconote.work/entry/2019/01/01/160100最近概率模型和神经网络相结合的研究变得多了起来，这次使用Uber开源的Pyro来实现一个贝叶斯神经网络。概率编程框架最近出了不少，Uber的Pyro基于Pytorch，Google的Edward基于TensorFlow，还有一些独立的像PyMC3,Stan,Pomegranate等等。...
				我们介绍了带变分推理的贝叶斯卷积神经网络，这是卷积神经网络（CNN）的一种变体，其中权重的难处理的后验概率分布是由Backprop的Bayes推断的。 我们证明我们提出的变分推断方法是如何实现的性能相当于频率论推理在几个数据集（MNIST，CIFAR10，CIFAR100），如所描述的相同结构。
贝叶斯vs频频方法中的过滤器权重分布
整个CNN的全贝叶斯视角
该存储库包含两种类型的贝叶斯lauer实现：
 BBB（Backprop的Bayes）： 基于。 该层分别对所有权重进行采样，然后将其与输入组合以从激活中计算出一个样本。
 BBB_LRT（使用本地重新参数化技巧的Backprop进行Bayes操作）： 这一层与本地重新参数伎俩结合贝叶斯通过Backprop。 这个技巧使得可以直接从激活中的分布中采样。
制作自定义贝叶斯网络？
 要创建自定义贝叶斯网络，请继承layers.m
				贝叶斯回归模型
这是尝试在Python中实现类似的库。
 它允许使用lme4语法（的子集）指定贝叶斯回归模型。 给定这样的描述和pandas数据框，该库将针对或生成模型代码和设计矩阵。
 。 显示库的笔记本已用于使模型适合数据。
这是系统可以处理的一些示例公式：
y ~ x
 人口层面的影响
y ~ 1 + x
 y ~ x1:x2
 变量之间的相互作用
y ~ 1 + x0 + (x1 | z)
y ~ 1 + x0 + (1 + x1 | z)
 y ~ 1 + x0 + (1 + x1 || z)
 组系数之间没有相关性
y ~ 1 + x0 + (1 + x1 | z1:z2)
 按多个因素分组（未测试）
 y ~ 1 + x0 + (x1 | z0) + (1 + x2 || z1)
 以上的组合
可以在各种粒度
				PyTorch中的贝叶斯LSTM实现
灵感来自：Uber的时间序列的深度和自信预测（2007） 
 由于贝叶斯神经网络提供了可量化的不确定性和置信区间，因此与等效的频频方法不同，贝叶斯神经网络引起了人们的极大兴趣。
 该存储库演示了PyTorch中的一种实现，并通过一个预测建筑能耗的真实示例总结了贝叶斯LSTM（长期短期记忆）网络的一些关键功能。 本示例中使用的Appliances能量预测数据集来自UCI机器学习存储库（  ）
 随附的笔记本直接从Google Colab共享。 结果，交互式可视化尚未转移到GitHub。
 请单击“在Colab中打开”按钮或单击此处，以在Google Colab中查看笔记本：  : //colab.research.google.com/drive/1pwMzsdRPwTO8oRVU0LwY9hs_z-ye67su
				pytorch贝叶斯网络Hyperparameters are the parameters in models that determine model architecture, learning speed and scope, and regularization.超参数是确定模型架构，学习速度和范围以及正则化的模型中的参数。
The search for optimal hyperpar...
				文章目录前言什么是贝叶斯神经网络How to train BNNBNN的损失函数
看了网上不少贝叶斯神经网络的文章，不少文章写的有点马虎，甚至一些说的不清不楚的文章，评论区许多人称赞是好文章，不禁让人怀疑他们是否真的看懂了文章。
本文将总结贝叶斯神经网络，首先，我将简单介绍一下什么是贝叶斯神经网络（BNN）；接着我将介绍BNN是怎么训练的，BNN是怎么预测的；最后，我会介绍BNN背后的运作原理。
如果您在阅读过程中发现了什么错误，请务必在评论区中指出，避免错误的观点在网上流传。
什么是贝叶斯神经网络
然后，你需要定义 ResNet18 的基本单元，即 Residual Block。Residual Block 包含两个卷积层，并且在两个卷积层之间使用批量归一化层和 ReLU 激活函数。
```python
class ResidualBlock(nn.Module):
    def __init__(self, in_channels, out_channels, stride=1, downsample=None):
        super(ResidualBlock, self).__init__()
        self.conv1 = nn.Conv2d(in_channels, out_channels, kernel_size=3, stride=stride, padding=1, bias=False)
        self.bn1 = nn.BatchNorm2d(out_channels)
        self.relu = nn.ReLU(inplace=True)
        self.conv2 = nn.Conv2d(out_channels, out_channels, kernel_size=3, stride=1, padding=1, bias=False)
        self.bn2 = nn.BatchNorm2d(out_channels)
        self.downsample = downsample
    def forward(self, x):
        residual = x
        out = self.conv1(x)
        out = self.bn1(out)
        out = self.relu(out)
        out = self.conv2(out)
        out = self.bn2(out)
        if self.downsample:
            residual = self.downsample(x)
        out += residual
        out = self.relu(out)
        return out
接下来，你需要定义 ResNet18 的主体部分，即多个 Residual Block 的堆叠。在这里，你可以使用 PyTorch 的 nn.Sequential 容器来简化模型的定义。
```python
class ResNet(nn.Module):
    def __init__(self, block, layers, num_classes=10):
        super(ResNet, self).__init__()
        self.in_channels = 64
        self.conv = nn.Conv2d(3, 64, kernel_size=7, stride=2, padding=3, bias=False)