Year: 2021

Abstract:
Ensembles of neural networks are often used to estimate epistemic uncertainty in high-dimensional problems because of their scalability and ease of use. These methods, however, are expensive to sample from as each sample requires a new neural network to be trained from scratch. We propose a new method, Generative Posterior Networks (GPNs), a generative model that, given a prior distribution over functions, approximates the posterior distribution directly by regularizing the network towards samples from the prior. This allows our method to quickly sample from the posterior and construct confidence bounds. We prove theoretically that our method indeed approximates the Bayesian posterior and show empirically that it improves epistemic uncertainty estimation over competing methods.

heading: 1 INTRODUCTION
text: Deep learning systems are becoming ubiquitous in many industries, specifically in applications with large datasets and where prediction mistakes are inexpensive. However, for many other industries, such as autonomous manufacturing, where datasets are smaller and mistakes can be catastrophic, the adoption of deep learning techniques has been slow. One of the critical challenges to applying deep learning models in these safety-critical domains is the lack of uncertainty estimation. If a model makes a prediction in a safety critical environment, we would like that model to be highly confident that prediction is correct. Uncertainty prediction can also provide significant benefit in data-poor domains by guiding exploration towards areas of uncertainty. Most deep learning models, however, are not able to estimate epistemic uncertainty – uncertainty deriving from lack of data samples. In this work, we introduce a new method for modeling epistemic uncertainty in neural networks.
We formulate this problem as a Bayesian Inference problem: given a prior distribution over functions and a set of training data, we want to model the posterior function distribution. There are many existing approaches to constructing such a posterior, such as Gaussian Processes (GPs) or Variational Inference (VI) methods, but these methods can require significant engineering to work in practice and can be expensive in high-dimensional settings. For this reason, many practitioners turn to using
ensembles of neural networks to approximate epistemic uncertainty for their ease of use in highdimensional settings. Pearce et al. (2020) illustrated that, given the correct regularization, these ensemble members approximate samples from the posterior distribution of functions. However, using this technique, every new sample of the posterior requires training a new neural network from scratch, which can be prohibitively expensive for large problems.
To address this challenge, we introduce Generative Posterior Networks (GPNs), a generative neural network model that directly approximates the posterior distribution by regularizing the output of the network towards samples from the prior distribution. By learning a low-dimensional latent representation of the posterior, our method can quickly sample from the posterior and construct confidence intervals. We prove that our method approximates the Bayesian posterior over functions and show empirically that our method can improves uncertainty estimation over competing methods.

heading: 2 RELATED WORK
text: There are two categories of uncertainty in the context of probabilistic modeling: aleatoric uncertainty and epistemic uncertainty. Aleatoric uncertainty refers to uncertainty inherent in a system, even when the true parameters of the system are known, while epistemic uncertainty refers to modeling uncertainty that can be reduced through collecting more data (Der Kiureghian & Ditlevsen, 2009; Gal, 2016). Deep learning models have shown impressive performance in quantifying aleatoric uncertainty when given enough data. On the other hand, despite many methods seeking to estimate epistemic uncertainty due to its well-motivated applications in Bayesian Optimization (Snoek et al., 2012; Springenberg et al., 2016), Reinforcement Learning (Curi* et al., 2020), and adversarial robustness (Stutz et al., 2020), quantifying it remains a challenging problem for deep learning models.
Gaussian Processes (GPs) remain one of the most effective methods for epistemic uncertainty estimation in low-dimensional, low-data settings. While scalability can be improved with interpolation methods (Hensman et al., 2013; Titsias, 2009; Wilson & Nickisch, 2015) and kernel learning methods (Wilson et al., 2016), GPs tend to perform worse than neural network-based methods on high-dimensional problems.
Another way to formulate epistemic uncertainty estimation is as a Bayesian Inference problem. Specifically, given a prior distribution over functions and some labeled data, the goal is to construct a posterior distribution in the Bayesian sense. Because this posterior distribution is usually intractable to compute exactly, approximate sampling methods are used instead. Markov-Chain Monte Carlo (MCMC) (Tanner & Wong, 1987) and Variational Inference (VI) (Blei et al., 2017) are some of the most common methods for posterior sampling. These methods, however, tend to struggle in high dimensional problems. The computational complexity of MCMC scales poorly with the dimension of the data and VI methods tend to require restrictive function families with strong independence assumptions. Bayesian Neural Networks (BNNs) (Blundell et al., 2015) are another way to perform Bayesian Inference a neural networks by modeling each parameter in the network with an independent Gaussian distribution. While BNNs perform well in low dimensional settings, their performance tends to degrade in high-dimensional settings that require large networks.
Gal & Ghahramani (2016) proved that dropout can also be used as a form of Bayesian Inference. Specifically, if a network is trained using dropout and L2 regularization, posterior samples can be collected by applying dropout at test-time. This method has surprisingly good out-of-distribution detection and can be applied to any network. However, it can be hard to tune the dropout parameter to have good uncertainty prediction without significantly degrading performance.
Neural network ensembling is another commonly used epistemic uncertainty estimation technique, which is closely related to our proposed method in this work. This technique involves training a number of completely separate neural networks on the same data to form an ensemble. Because of the stochasticity of initialization and the training process, each member of the ensemble should be a slightly different function. Specifically, the ensemble members should produce similar outputs in data-dense regions of the state-space and dissimilar outputs in data-sparse regions. Pearce et al. (2020) proved that, by using the correct regularization, each ensemble member becomes an approximate sample from the function posterior. Ensembling methods are particularly attractive to practitioners because they require very little fine-tuning and work surprisingly well in practice. However, these methods have one main drawback, namely that each new sample from the posterior
requires training a new network from scratch. Our method, on the other hand, seeks to construct a generative posterior model using similar regularization techniques to Pearce et al. (2020), allowing for quick sampling from the posterior.

heading: 3 PROBLEM STATEMENT AND BACKGROUND
text: We will focus on the problem of Bayesian Inference on the parameters of θ of some function f(·;θ). We assume we have some known prior over parameters Pprior(θ) ∼ N (µprior,Σprior) as well as access to a noisy dataset (xobs,yobs) = {(x1obs, y1obs) . . . , (xNobs, yNobs)}, where yiobs = f(xiobs;θ) + , ∼ N (0, σ ) for some unknown parameter θ. We can then define the data-likelihood as
Plike(yobs |θ,xobs) = ∏ i N (yiobs | f(xiobs;θ), σ ). (1)
The goal is to approximate samples from the posterior, Ppost(θ |yobs,xobs) ∝ Plike(yobs |θ,xobs)Pprior(θ). (2)
For ease of notation, we will drop the xobs from the likelihood.
For our analysis, it is convenient if the posterior distribution is Gaussian. To show this, we will use the fact that the product of two Gaussian PDFs is proportional to a Gaussian PDF. We make a distinction here between the data-likelihood, a function of the observed data yobs, and the parameterlikelihood, a function of the parameters θ, defined as follows:
Pparam-like(θ; yobs) = Plike(yobs |θ). (3) The parameter-likelihood PDF is equivalent to the data-likelihood PDF for all values of yobs and θ, but is notably a function of the parameters θ. For a more detailed explanation of the distinction between the parameter and data-likelihoods, see Definition 1 from Pearce et al. (2020).
When the function f is linear, it can be shown that the parameter-likelihood is also Gaussian, Pparam-like(θ; yobs) ∝ N (µlike,Σlike). (4)
While this is not true in general, Pearce et al. (2020) assume that the parameter likelihood is Gaussian so that the posterior also becomes Gaussian, Ppost(θ |yobs) ∝ Pparam-like(θ; yobs)Pprior(θ) ∝ N (µlike,Σlike)N (µprior,Σprior) = N (µpost,Σpost), (5) where the posterior covariance and mean are given by
Σpost = ( Σ−1like + Σ −1 prior )−1 , µpost = ΣpostΣ −1 likeµlike + ΣpostΣ −1 priorµprior. (6)
Note, however, that after the reformulation of the problem in Section 4.1, we no longer require this assumption for our theoretical analysis.

heading: 3.1 RANDOMIZED MAP SAMPLING (RMS)
text: Since the values µlike and Σlike are intractible to compute in practice, the challenge of computing the posterior distribution remains. However, Pearce et al. (2020) showed how we can use a method called Randomized MAP sampling (RMS) to approximately sample from the posterior distribution. With RMS we first assume access to a mechanism to compute the maximum a posteriori (MAP) solution for our posterior. Then, by randomly shifting the prior distribution we use in our MAP solver, we can estimate samples from the posterior.
Specifically, RMS samples “anchor points” θanc from some distribution θanc ∼ N (µanc,Σanc). Each anchor point defines a shifted prior distribution Panc(θanc) = N (θanc,Σprior). Now we define the MAP function which finds the MAP solution using this shifted prior:
θMAP(θanc) = arg max θ
Plike(yobs |θ)Panc(θ). (7)
By sampling different anchor points, we can then construct a probability distribution over MAP solutions P (θMAP(θanc)), which can be shown to be Gaussian. Using moment matching, Pearce et al. (2020) proved, in Theorem 1, that if we choose µanc = µprior and Σanc = Σprior + ΣpriorΣ−1likeΣprior, then P (θMAP(θanc)) = Ppost(θ |yobs). Using the approximation Σanc ≈ Σprior, a sample from posterior, θpost, can be collected by sampling θanc ∼ N (µanc,Σanc) and optimizing the following:
θpost = arg max θ logPparam-like(θ; yobs) + logPanc(θ). (8)

heading: 4 GENERATIVE POSTERIOR NETWORKS
text: While RMS can be used to construct samples from the posterior, every new sample requires solving a new optimization problem. Instead, we propose learning the MAP function itself.
Let g be a neural network parameterized by some vector, φ, that takes as input a sample from the anchor distribution, θanc ∼ Panc(θanc), and outputs θMAP(θanc). If we assume that g has enough expressive power to represent the MAP function, then g(θanc;φ) = θMAP(θanc) if:
φ = arg max φ Eθanc∼Panc(θanc) logPlike(yobs | g(θanc;φ)) + logPanc(g(θanc;φ)). (9)
In practice, however, this is a very hard optimization problem for two main reasons. (1) Optimizing in parameter space tends to be challenging for gradient-descent optimization. (2) The space of all parameters θanc is not only very large, but there are also many parameters θanc that map to the same outputs for all inputs. Moreover, we would expect the posterior parameters to be highly correlated and, therefore, need less expressive power to be represented.
To address these two problems, we will first change the problem definition slightly to find the posterior in output space as opposed to parameter space. We prove in the next section that, by changing our loss function, we can still use RMS to approximate the posterior in output space. Next, we show practically how we construct a low-dimensional representation of the anchor distribution.

heading: 4.1 THEORETICAL RESULTS
text: As mentioned above, optimizing a generative model to output parameters can be challenging. Instead, we reformulate the Bayesian Inference problem slightly and consider the posterior in output space. Concretely, consider some discretization of the input space, xsample = {x1sample, . . . , xMsample}. Let Ŷi = f(xisample;θ) be a set of sample points and Ŷ = {Ŷ1, . . . , ŶM} be a transformation of the random variable θ, where each element is the output of the function f parameterized by θ evaluated at points xsample. The prior and likelihood for this transformed random variable Ŷ can be expressed in terms of the prior and likelihood of θ:
P (Ŷ) = ∫ θ 1Ŷi=f(xisample;θ)∀i Pprior(θ)dθ, (10)
P (yobs | Ŷ) ∝ ∫ θ 1Ŷi=f(xisample;θ)∀i Plike(yobs |θ)Pprior(θ)dθ. (11)
Note that samples from the prior of Ŷ can be easily obtained by sampling θ ∼ Plike(θ) and computing Ŷi = f(xisample;θ) for all i. The goal in this reformulated problem is to find the posterior P (Ŷ |yobs) ∝ P (yobs | Ŷ)P (Ŷ). As before, we estimate the posterior using RMS and define a new MAP function in the output space:
ŶMAP(Ŷanc) = arg max Ŷ
P (yobs | Ŷ)Panc(Ŷanc) (12)
for some anchor distribution Ŷanc ∼ Panc(Ŷanc) = N (µanc,Σanc). We now need to show that, with the correct choice of µanc and Σanc, P (ŶMAP(Ŷanc)) = P (Ŷ |yobs).
To simplify the analysis, we assume that the output prior is normally distributed P (Ŷ) = N (µŶ,ΣŶ). Since the parameters are normally distributed, this approximation becomes exact when either f is linear or the final layer is infinitely wide.
Next, we need to show that the likelihood function P (yobs | Ŷ) is roughly normally distributed with respect to the vector Ŷ. Recall that each observed data point, yiobs, is equal to a noisy output of f for some parameter θ, yiobs = f(x i obs;θ) + . In the case where xobs ⊂ xsample, yiobs is evaluated on the same points as Ŷ, so the likelihood of P (yobs |θ) becomes equivalent to P (yobs | Ŷ). By the same reasoning, if we choose xsample to contain all the points in the domain X , P (yobs |θ) becomes equivalent to P (yobs | Ŷ).
Lemma 1. Let f be some neural network parameterized by θ. Let xobs,yobs correspond to vectors of observed labeled data. Let xsample be some discretization over the input space of Md points such that for any x ∈ X there exists an x′ ∈ xsample where ‖x − x′‖2 ≤ lM for some constant l. Define Ŷ := {f(x;θ) : ∀x ∈ xsample} be a vector of evaluations of f on all points in xsample for a given parameter vector θ.
Then, if we take the limit as M →∞, P (yobs | Ŷ) = P (yobs |θ).
Proof. See Appendix A.1.
It turns out that, since the observation noise is additive, the likelihood in parameter space P (yobs |θ) is also Gaussian as a function of the output vector Ŷ. Thus, if we take |xsample| → ∞ and apply Lemma 1, we get that the output likelihood P (yobs | Ŷ) is Gaussian as a function of Ŷ.
Now, since we have that both the prior and likelihood of Ŷ can be expressed as Gaussian functions of Ŷ, we also have that the posterior is Gaussian. As before, we just need to show that the distribution of MAP solutions P (ŶMAP(Ŷanc)) is also Gaussian, then apply moment matching to find the right values of µanc and Σanc such that P (ŶMAP(Ŷanc))) = P (Ŷ |yobs). Theorem 1. Let f be some neural network parameterized by θ. Let xsample be some vector of inputs sampled uniformly from XM . Define Ŷ := [f(x1;θ), . . . , f(xM ;θ)] as a transformation of the random variable θ. Assume the prior of Ŷ follows a multivariate normal distribution, P (Ŷ) = N (µŶ,ΣŶ). We assume access to some function which takes as input the center of a shifted distribution and outputs the MAP estimate of Ŷ, ŶMAP(Ŷanc).
If we choose the distribution over Ŷanc to be P (Ŷanc) = N (µanc,Σanc), where µanc = µŶ and Σanc = ΣŶ + ΣŶΣ −1 likeΣŶ, then limM→∞ P (ŶMAP(Ŷanc)) = P (Ŷ |yobs).
Proof. See Appendix A.2.

heading: 4.2 PRACTICAL IMPLEMENTATION
text: In practice, as was done in Pearce et al. (2020), we will instead use the approximation Σanc ≈ ΣŶ to make it tractable to sample from the anchor distribution. Pearce et al. (2020) argue that this approximation, in general, causes RMS to over-estimate the posterior variance.
With this assumption, the anchor distribution Panc(Ŷanc) becomes approximately equal to the prior distribution P (Ŷ). Recall that, to sample from the output prior, we simply sample a set of parameters θ ∼ Pprior(θ) and evaluate the function f at all evaluation points xsample. Now, if we want to learn the MAP function in the output space, we can construct a neural network g, parameterized by some vector φ, that takes a sample from the output anchor distribution, Ŷanc ∼ N (µanc,Σanc), along with the evaluation points xsample, and outputs θMAP(θanc). If we assume that g has enough expressive power to represent the MAP function, then g(Ŷanc;φ) = ŶMAP(Ŷanc), if:
φ = arg max φ Eθanc∼Pprior(θanc) logPθ(yobs | g(xsample,θanc;φ)) + logPanc(g(xsample,θanc;φ)) (13)
This can be re-written as:
φ = arg min φ Eθanc∼Pprior(θanc) N∑ i=1 ‖yiobs − g(xiobs,θanc;φ)‖22 + σ2 δTΣ−1Ŷ δ, (14)
where δj = g(x j sample,θanc;φ)− f(x j sample;θanc). See Appendix A.3 for a step-by-step derivation.
In practice, of course, we cannot compute the regularization term, 1N σ 2 δ TΣ−1 Ŷ δ as |xsample| → ∞. Instead, we approximate this term using a finite number of sample points from the input space X . Because of this, for high-dimensional problems, the off-diagonal terms of ΣŶ are near zero almost
always. However, we found that, in practice, assuming that the samples of Ŷ are independent did not significantly degrade performance. Using this approximation, our optimization problem becomes:
φ = arg min φ Eθanc∼Pprior(θanc),xsample∼Unif(X ) N∑ i=1 ‖yiobs − g(xiobs,θanc;φ)‖22 + β‖δ‖22, (15)
where δ = g(xsample,θanc;φ)− f(xsample;θanc) and β is a hyper-parameter.
In order for g to have enough expressive power to represent the full MAP function, it may require the size of the parameters φ to be much larger than θ. However, after applying the low-dimensional approximation of θanc described below, in our experiments, we only used roughly double the number of parameters compared to each member of the ensemble method. And since the ensemble we used has 10 members, the total number of parameters needed for our method is significantly smaller.

heading: 4.3 LOW-DIMENSIONAL EMBEDDING
text: As mentioned above, optimizing the generative function g over the the space of all anchor parameters θanc did not yield great results because the space was too large. Instead, we sought to learn a lowdimensional representation of the anchor parameters. There are two key reasons why we would expect to be able to decrease the representational power of the anchor distribution and maintain the same of fidelity in the posterior estimation. (1) Because of the nature of neural networks, there are many settings of parameters θ that would result in the same output vector Ŷ. Because we are focusing on the posterior of Ŷ, we need less representational power to model Ŷ. (2) And, maybe more importantly, because the posterior parameters are highly correlated, we would expect to need less expressive power to represent the posterior distribution than the prior distribution.
Thus, we would like to construct a simpler estimate of the prior space using a low-dimensional embedding vector, z. That is, we want our generative model to take as input the embedding vector z ∼ N (0, I), instead of θanc, in order to estimate the MAP function. To do this, we need to learn a mapping from anchor parameters θanc to embedding vectors z. For our experiments, we used a 1-1 embedding scheme where we sample k parameters from the true prior θ1, . . . ,θk ∼ Pprior(θ) and k independent samples from our embedding z1, . . . , zk. We then jointly optimize over φ and z1, . . . , zk as follows:
arg min φ,z1,...,zk Ej∈[1,...,k] N∑ i=1 ‖yiobs − g(xiobs, zj + ;φ)‖22 + β‖δ‖22 + Lreg(z1, . . . , zk), (16)
where δj = g(x j sample, zj + ;φ) − f(x j sample;θj), ∼ N (0, I) is a noise injection vector that improves the smoothness of interpolations between embedding vectors, and Lreg(z1, . . . , zk) is a regularizer to keep z1, . . . , zk roughly normally distributed, which allows us to easily sample from the embedding space.. For our experiments we use the KL divergence between z and the normal distribution, Lreg(z1, . . . , zk) = DKL(N (z̄, sz),N (0, 1)). Figure 1 illustrates how we can sample from this embedding space to construct posterior functions.

heading: 4.4 CLASSIFICATION
text: Our method uses the assumption that the prior over Ŷ is approximately normally distributed. While this approximation is exact, assuming the prior parameters are normally distributed, only when either the network is linear or the final layer is infinitely wide and does not have an activation function, for most deep networks, this is a good approximation. However, if we add a soft-max to the final output, as is usually done for classification tasks, this assumption is violated. To get around this, for classification tasks, we add the anchor loss to the pre-softmax outputs.

heading: 5 EXPERIMENTS
text: The goal of our experiments is to illustrate the ability of our method to accurately model epistemic uncertainty for out of distribution data while retaining high performance on in distribution data. We will test our method on three problems: (1) a small-scale regression task, (2) Bayesian Optimization (BO), and (3) Classification with out of distribution examples.
We compare our method against 3 competing Bayesian methods: Gaussian Processes (GPs), Bayesian Dropout (Gal & Ghahramani, 2016), and anchor-regularized neural-network ensembles (Pearce et al., 2020) with 10 ensemble members. For our BO experiments we use exact GP with a static kernel since the problems are sufficiently low-dimensional, but for our classification experiments we need to use approximate GPs. Specifically, we use a grid-interpolated GP with Deep Kernel Learning (Wilson et al., 2016) implemented with the GPyTorch library (Gardner et al., 2018). For the anchor-regularized ensembles, since every member requires a unique set of parameters (and anchor points) the number of ensemble members was limited by the memory capacity of our GPUs.

heading: 5.1 SMALL SCALE REGRESSION
text: Our first experiment is a small 1-dimensional regression problem where we can easily compute the ground truth. We provide each method with the same 6 observations. For this experiment, we compute 100 samples from the approximate posterior for each method (except for the Ensemble where we only use 10 samples). For ground truth, we use the Metropolis-Hastings MCMC algorithm. Figure 2 shows the results.

heading: 5.2 BAYESIAN OPTIMIZATION
text: The goal of Bayesian Optimization is to find the maximum of a function in as few samples as possible. This is commonly done by first constructing a posterior over possible functions given the samples we have collected so far, then, through various methods, choosing sample points that might maximize the true function. By using different methods for predicting the posterior in BO, we can evaluate how well these methods capture epistemic uncertainty.
For all methods we use Thompson Sampling and implemented the experiments using the BOTorch library by Balandat et al. (2020). We tested each method on five standard benchmark BO functions: the Branin, Hartmann-6D, Shekel, Michalewicz-10D, and Powell functions. We normalize the output of each function between a minimum value (found by random sampling) and the optimal value. To make sure that each method had a large enough sample variance, we chose hyper-parameters such that 100 prior samples from each method have a range of roughly [−1, 1]. To choose the best regularization hyper-parameter, we ran each method with 5 different regularization parameters, then chose the hyper-parameters with the highest average performance across all BO problems.
Figure 3 shows the best normalized function value sampled so far during optimization averaged across the 5 BO problems each ran for 5 trials. Our method significantly out-performs Random sampling and even slightly outperforms GPs, one of the best BO techniques for low-dimensional problems. Even though our method trails slightly behind Bayesian Dropout, these results still illustrate the uncertainty estimation power of our method, even on low-dimensional problems.

heading: 5.3 CLASSIFICATION
text: We also look at the problem of image classification. Along with our method having high accuracy, we are also interested in how well our method (and the baseline methods) is able to capture epistemic uncertainty. To measure this we split the dataset into two sets based on class labels: the in distribution dataset and the out of distribution dataset. During training, we only provide the labels for
the in-distribution dataset. During testing, along with recording the accuracy on the in distribution dataset, we also sample our models on the out of distribution dataset. For such out of distribution examples, we expect the models to show high epistemic uncertainty and, thus, the samples for these out of distribution examples should have high variance.
Note that we still provide the full set of unlabeled training data to our agent. This obviously favors our method over the baseline methods as our method is the only method that is able to take advantage of this unlabeled data. However, from a practical perspective this is a reasonable setup as many realworld datasets have far more unlabeled data than labeled data.
One desired property of an epistemic uncertainty estimator is in constructing confidence sets. For the classification setting, the confidence interval (CI) we consider is the set of possible class predictions. We would expect that this confidence interval contains the true class with high probability. We also expect that, on in-distribution examples, the confidence interval only contains the true label with high probability. To construct these confidence intervals, we take 100 samples from each model on every data point in the test data. We define the metrics CI-Any as the proportion of data points for which at least one sample contains the correct label and CI-All as the proportion of data points for which every sample predicts the correct label.
Another desired property of an epistemic uncertainty estimator is in predicting when a data point is outside the training distribution and, thus, any classification prediction made on such a data point will likely be incorrect. For out of distribution data points, we expect high sample variance from the posterior model and low sample variance from in distribution data points. Thus, variance can be used as a metric to predict which test data points are outside the training distribution. To quantify each method’s ability to perform out of distribution detection using sample variance, we construct ROC curves for each method and measure the area under the curve (AUC).
We tested all methods on both the MNIST dataset and the CIFAR-10 dataset. For MNIST, the in distribtuion classes were [0, 1, 2, 3, 4, 5] and the out of distribution classes were [6, 7, 8, 9]. For CIFAR-10, the in distribtuion classes were [“airplane”, “automobile”, “bird”, “deer”, “dog”, “ship”] and the out-of-distribution classes were [“cat”, “frog”, “horse”, “truck”]. To make sure that each method had a large enough sample variance, we chose hyper-parameters such that prior distributions had near 100% CI-All. To choose the best regularization hyper-parameter, as in the BO experiments, we ran each method with 5 different regularization parameters, then chose the model with both a high validation accuracy and CI-All.
For all other methods, we take 100 samples for measuring CI-Any, CI-All, and sample variance.
Table 1 summarizes the performance of each model using the metrics of test accuracy, CI-Any, CI-All, and the ROC AUC for out of distribution detection. Full ROC curves for each method are plotted in Figure 4. While all methods are able to achieve high prediction accuracy and high CI-All on in distribution examples, our method is the only one to achieve high CI-Any on out of distribution examples, illustrating our method’s ability to construct valid confidence intervals. Additionally, our method significantly out performs all other methods on out of distribution detection.
Figure 5 shows two informative test images, one from the in distribution dataset and one from the out of distribution dataset, for both MNIST and CIFAR-10, along with box plots of sampled PMFs from each of the learned models. While all methods predict the correct label with high precision on the in distribution example, ours has a uniquely wide prediction distribution on the out of distribution example, illustrating high predicted epistemic uncertainty.

heading: 6 CONCLUSION
text: In this paper we introduce Generative Posterior Networks (GPNs), a method to learn a generative model of the Bayesian posterior distribution by regularizing the outputs of the network towards samples of the prior. We prove that under mild assumptions, GPNs approximate samples from the true posterior. We then show empirically that our method is not only competitive with the best uncertainty predictions techniques on small Bayesian Optimization problems, it significantly outperforms these competing methods on high-dimensional classification tasks.

heading: A PROOFS
text: A.1 PROOF OF LEMMA 1
Proof. It is a standard result that linear feed-forward neural networks have bounded Lipschitz constants. We’ll denote L(θ) as the Lipschitz constant of f for a given set of parameters θ.
Recall that this likelihood is defined as follows: P (yobs | Ŷ ∝ ∫ θ 1Ŷi=f(xisample;θ)∀i (∏ i N (yiobs|f(xiobs;θ), σ ) ) Pprior(θ)dθ
Let xsample′ ⊂ xsample be the subset of points in our discretization, xsample, that are closest to the observed points, xobs; in other words, for each element xiobs ∈ xobs, there exists an element xi ∈ xsample
′ where: x′i = min
x∈xsample ‖xiobs − x‖2
We denote the elements of xsample′ as [x′1, . . . , x ′ N ]. Let Ŷ ′ ⊂ Ŷ be the corresponding elements in Ŷ; in other words, Ŷ ′i = f(x ′ i;θ) for i ∈ [1, . . . , N ].
By definition of our discretiztaion, we know that for all i ∈ [1, . . . , N ], ‖xiobs − x′i‖ ≤ lM . Using the Lipschitz-continuity of f , for any setting of parameters θ, ‖f(xiobs;θ)− Ŷ ′i ‖ ≤ lML(θ).
If we take the limit as M →∞, f(xiobs;θ) = Ŷi and
P (yobs | Ŷ) ∝ ∫ θ 1Ŷi=f(xisample;θ)∀i (∏ i N (yiobs|Ŷ ′i , σ ) ) Pprior(θ)dθ
= (∏ i N (yiobs|Ŷ ′i , σ ) )∫ θ 1Ŷi=f(xisample;θ)∀i Pprior(θ)dθ
= (∏ i N (yiobs|Ŷ ′i , σ ) ) P (Ŷ)
= (∏ i N (yiobs|f(xiobs;θ), σ ) ) P (Ŷ)
∝ P (yobs |θ)
Now we know that P (yobs | Ŷ) ∝ P (yobs |θ). But since both P (yobs | Ŷ) and P (yobs |θ) are both probability distributions of yobs, then P (yobs | Ŷ) = P (yobs |θ).
A.2 PROOF OF THEOREM 1
Proof. We will start by showing that P (Ŷ|yobs) is normally distributed. By Lemma 1, we know that limM→∞ P (yobs | Ŷ) = P (yobs | Ŷ). Thus, in the limit as M →∞,
P (Ŷ |yobs) ∝ P (yobs | Ŷ)P (Ŷ) = P (yobs |θ)P (Ŷ)
= (∏ i N (yiobs|f(xiobs;θ), σ ) ) P (Ŷ)
= (∏ i N (f(xiobs;θ)|yiobs, σ ) ) P (Ŷ)
= N (Ŷ|µlike,Σlike)N (Ŷ|µŶ,ΣŶ)
for some mean and covariance µlike,Σlike. Since the product of two multivariate Gaussians is a Gaussian, we know that:
P (Ŷ |yobs) = N (µpost,Σpost) (17) where
Σpost = ( Σ−1like + Σ −1 Ŷ )−1 , µpost = ΣpostΣ −1 likeµlike + ΣpostΣ −1 Ŷ µŶ. (18)
Now we consider at the distribution P (ŶMAP(Ŷanc)) where Ŷanc ∼ N (µanc,Σanc). We will show that, when we set
µanc = µŶ, Σanc = ΣŶ + ΣŶΣ −1 likeΣŶ
the distribution P (ŶMAP(Ŷanc)) becomes equal to the posterior distribution N (µpost,Σpost). There are three steps needed to show equality between these distributions:
1. Show that P (ŶMAP(Ŷanc)) is normally distributed for some mean and variance, µRMSpost ,Σ RMS post .
2. Show that µRMSpost = µpost.
3. Show that ΣRMSpost = Σpost.
Using the same reasoning as above, in the limit as M →∞,
P (ŶMAP(Ŷanc)) = arg max Ŷ
P (yobs|Ŷ)Panc(Ŷ)
= arg max Ŷ
N (Ŷ|µlike,Σlike)N (Ŷ|µanc,Σanc)
Since the max of a Gaussian is the mean, then
FMAP(Ŷanc) = AŶanc + b (19)
where we define:
A = ΣpostΣ −1 Ŷ
(20)
b = ΣpostΣ −1 likeµlike (21)
We will now show that E[ŶMAP(Ŷanc)] = µpost. Because we set µanc = µŶ:
E[ŶMAP(Ŷanc)] = E[AŶanc + b] = AE[Ŷanc] + b = AµŶ + b
= ΣpostΣ −1 Ŷ µŶ + ΣpostΣ −1 likeµlike
= µpost
Finally, we will show that Var[ŶMAP(Ŷanc)] = Σpost.
Var[ŶMAP(Ŷanc)] = Var[AŶanc + b] = AVar[Ŷanc]AT
= (ΣpostΣ −1 Ŷ )(ΣŶ + ΣŶΣ −1 likeΣŶ)(ΣpostΣ −1 Ŷ )T = (Σpost + ΣpostΣ −1 likeΣŶ)(Σ −1 Ŷ Σpost) = ΣpostΣ −1 Ŷ Σpost + ΣpostΣ −1 likeΣpost = Σpost(Σ −1 Ŷ + Σ−1like)Σpost
= Σpost
So, P (ŶMAP(Ŷanc)) ∼ N (µRMSpost ,ΣRMSpost ) where µRMSpost = µpost and ΣRMSpost = Σpost.
Thus, in the limit as M →∞, P (ŶMAP(Ŷanc)) = P (Ŷ|yobs).
A.3 DERIVATION OF LOSS FUNCTION
We start with
φ = arg max φ Eθanc∼Pprior(θanc) logPlike(yobs | g(xsample,θanc;φ)) + logPanc(g(xsample,θanc;φ))
If we choose µanc = µŶ and Σanc = ΣŶ, we can simplify the logarithm of the likelihood and anchor distributions as follows:
logPlike(yobs | g(xsample,θanc;φ)) = ∑ i logN (yiobs|g(xisample,θanc;φ), σ )
= − 1 2σ2 ∑ i ‖yiobs − g(xisample,θanc;φ)‖22
logPanc(g(xsample,θanc;φ)) = logN (g(xsample,θanc;φ)|µŶ,ΣŶ)
= −1 2 δTΣ−1 Ŷ δ
where δj = g(x j sample,θanc;φ)− f(x j sample;θanc).
Putting these together, we get:
φ = arg max φ
Eθanc∼Pprior(θanc) − 1
2σ2 ∑ i ‖yiobs − g(xisample,θanc;φ)‖22 − 1 2 δTΣ−1 Ŷ δ
= arg min φ
Eθanc∼Pprior(θanc) 1
N ∑ i ‖yiobs − g(xisample,θanc;φ)‖22 + 1 N σ2 δ TΣ−1 Ŷ δ

Review summaries:
NN ensembles are useful models to estimate epistemic uncertainty in the predictions. Since their computational cost is generally high, they introduce Generative Posterior Networks (GPNs) as a generative NN model to approximate the posterior distribution of the inference process. Under not very restrictive assumptions, GPNs work as well as previous methods, if not better. This is shown through extensive experiments, especially in classification problems.

This paper built upon Pearce et al 2020 that shows regularizing parameters to prior distributions in ensemble of NNs results in the approximate Bayesian inference, aiming to improve computational efficiency of the previous method. Specifically, in the method in Pearce et al, sampling one parameter from posterior requires to find a MAP solution for some anchor point, which can be computationally demanding. Therefore, the authors propose to learn MAP function from an anchor point to MAP solution, which removes the required MAP optimization procedure.  

The paper develops an uncertainty quantification method as an extension of the Bayesian ensemble method. Instead of training multiple ensemble models through MAP optimization, the proposed method tries to learn a mapping function between the prior distribution and the posterior distribution of model parameters, which avoids the complex training of ensemble models and achieves better efficiency. Due to the high-dimensional parameter space, the authors propose to learn the mapping function in the low-dimensional output space with theoretical justifications. The experiments include image classification and Bayesian optimization, where authors compare the proposed method with Bayesian ensemble and MC-dropout in terms of uncertainty estimation.