Global inducing point variational posteriors for Bayesian neural networks and deep Gaussian processes
Authors: Sebastian W. Ober
Laurence Aitchison
Venue:  (2021)
Keywords: Bayesian neural networks
deep Gaussian processes
variational inference
inducing points

Abstract: We derive the optimal approximate posterior over the top-layer weights in a Bayesian neural network for regression, and show that it exhibits strong dependencies on the lower-layer weights. We adapt this result to develop a correlated approximate posterior over the weights at all layers in a Bayesian neural network. We extend this approach to deep Gaussian processes, unifying inference in the two model classes. Our approximate posterior uses learned "global" inducing points, which are defined only at the input layer and propagated through the network to obtain inducing inputs at subsequent layers. By contrast, standard, "local", inducing point methods from the deep Gaussian process literature optimise a separate set of inducing inputs at every layer, and thus do not model correlations across layers. Our method gives state-of-the-art performance for a variational Bayesian method, without data augmentation or tempering, on CIFAR-10 of $86.7\%$.

Under review as a conference paper at ICLR 2021
GLOBAL INDUCING POINT VARIATIONAL POSTERIORS FOR BAYESIAN NEURAL NETWORKS AND DEEP GAUS- SIAN PROCESSES
Anonymous authors Paper under double-blind review
ABSTRACT
We derive the optimal approximate posterior over the top-layer weights in a Bayesian neural network for regression, and show that it exhibits strong depen- dencies on the lower-layer weights. We adapt this result to develop a correlated approximate posterior over the weights at all layers in a Bayesian neural network. We extend this approach to deep Gaussian processes, unifying inference in the two model classes. Our approximate posterior uses learned "global" inducing points, which are defined only at the input layer and propagated through the network to ob- tain inducing inputs subsequent layers. By contrast, standard, "local", inducing point methods from the deep Gaussian process literature optimise a separate set of inducing inputs every layer, and thus do not model correlations across layers. Our method gives state-of- of-the-art performance for a variational Bayesian method, without data augmentation or tempering, on CIFAR-10 of 86.7%.
INTRODUCTION
Deep models, formed by stacking together many simple layers, give rise to extremely powerful machine learning algorithms, from deep neural networks (DNNs) to deep Gaussian processes (DGPs) (Damianou & Lawrence 2013). One approach to reason about uncertainty in these models is to use variational inference (VI) (Jordan al., 1999). VI in Bayesian neural networks (BNNs) requires the user to specify family of approximate posteriors over the weights, with the classical approach being independent Gaussian distributions over each individual weight (Hinton & Van Camp 1993 Graves 2011 Blundell et al., 2015). Later work has considered more complex approximate posteriors, for instance using a Matrix-Normal distribution as the approximate posterior for a full weight-matrix
(Louizos & Welling, 2016 Ritter et al. 2018). By contrast, DGPs use approximate posterior defined over functions the standard approach is to specify the imputs and outputs at a finite number
of "iinducing" points (Damianou & Lawrence 2013 Salimbeni & Deisenroth] 2017).
Critically, these classical BNN and DGP approaches define approximate posteriors over functions that independent across layers. An approximate posterior that factorises across layers is problematic, because what matters for a deep model the overall input-output transformation for the full model, not the input-output transformation for individual layers. This raises the question of what family of approximate posteriors should be used to capture correlations across layers. One approach for BNNs would be introduce a flexible "hypernetwork"" used generate the weights (Krueger et al. 2017 Pawlowski et al., 2017). However, this approach is likely to be suboptimal as it does not suffciently exploit the rich structure in the underlying neural network. For guidance, we consider the optimal approximate posterior over the top-layer units in deep network for regression. Remarkably, the optimal approximate posterior for the last-layer weights given the earlier weights can be obtained in closed form without choosing restrictive family of distributions. In particular, the optimal approximate posterior given by propagating the training inputs through lower layers to compute the top-layer representation, then using Bayesian linear regression map from the top-layer representation to the outputs.
Inspired by this result, we use Bayesian linear regression to define generic family of approximate posteriors for BNNs. In particular, we introduce learned "pseudo-data" at every layer, and compute the posterior over the weights by performing linear regression from the inputs (propagated from

---

Under review as a conference paper at ICLR 2021
lower layers) onto the datudo-data. We reduce the burden of working with many training inputs by summarising the posterior using small number of "inducing" points. We find that these approximate posteriors give excellent performance in the non-tempered, no-data-un ftataa-ugmeettttion regime, with performance on datasets such as CIFAR-I0 reaching 86.7%, comparable to SGMCMC (Wenzel et al. 2020). Our approach can be extended to DGPs, and explore connections the indueing point GP literature, showing that inference in the two classes of models can be unified.
2
METHODS
We consider neural networks with lower-layer weights (We)=r. We E RN--LING and top-layer weights, W L+11 E R 1111,
11,1, where the activity, F,, at layer l is given by,
Fi FI=
oF= We for E (2....,4), L}.
where o6() an elementwise nonlinearity. The outputs, Y PxN++1 depend on the top-level activity. Fz. and the output weights, WL+1. according a likelihood, P (Y W L+11 F In the following derivations, we will focus on > 1; corresponding expressions for the imput layer can be obtained by replacing ((Fo) with the inputs, X PNNo, The prior weights is independent across layers and output units (see Sec. 2.3 for the form of St).
(We) = 11 (w)),
P (ww) = =N (w3 (w 0, S.).
where w is a column of W., i.e. all the input weights to unit in layer E. To fit the parameters of the approximate posterior, Q ((W) we maximise the evidence lower bound (ELBO),
Fo (W) [log (Y, ( x log Q (( ((W)
To build intuition about how to parameterise Q ((W) we consider the optimal Q (WL+1R(Wa W2))) for any given Q ((We)/z)) We begin by simplifying the ELBO by in- corporating terms that do not depend on WL+1 into a constant, c,
C= - [log (Y, W L+1 X, (W - log Q (Wa+l (W + +c].
(4)
Rearranging these terms, we find that W 1++1 dependence can be written terms of the divergence between the approximate posterior of interest and the true posterior,
L = Bo((w)) log P (Y (Wo)
DRL (Q (W L+1 (W.)) (W Y,X, (Wo) +c]. (5)
Thus, the optimal approximate posterior is,
Q (W L+1 ( P (WL+1) Y.LHi/YXX... (WOH2) & P (L (WW+1))
(6)
and where the final proportionality comes by applying Bayes theorem and exploiting the model's conditional independencies. For regression, the likelihood is Gaussian,
P W F2) N N (yx: (Fz) wx ,AZL))
where yx is the value of single output channel for all training inputs, and A++1 is precision matrix. Thus, the posterior given in closed form by Bayesian linear regression (Rasmussen & Williams 2006).
2.1 DEFINING THE FULL APPROXIMATE POSTERIOR WITH GLOBAL INDUCING POINTS PSEUDO-DATA
We adapt the optimal scheme above to give a scalable approximate posterior over the weights at all layers. To avoid propagating all training inputs through the network, which is intractable for large datasets, we instead propagate M global inducing locations, Ug.
U U, W1.
2...., U o (U+-1) We for f = 22 ..... . L 1.
(8)
2

---

Under review as a conference paper at ICLR 2021
Next, the optimal posterior requires outputs, Y. However, no outputs available inducing locations for the output layer, let alone for intermediate layers. We thus introduce learned variational parameters to mimic the form of the optimal posterior. In particular, over weights and a "psendo-likelhhodd" .V (v1: ug, A representing noisy "pseuudo-observations" of the outputs of the linear layer at the inducing locations, u = (U wX. Substituting u into the pseudo-likelihood likelihood the approximate posterior becomes,
the product of the prior
Q (Wellee 18-1
No (AIO((U-)))) wY, A, (wi), II W
Q W
(We)e t-1
Ex (U AeVi 2ạ w
Σx
(
where v and A₂ variational parameters. Therefore, our full approximate posterior factorises as
((W) II W: (Wel
(10)
Thus, the full ELBO can be written,
L T==Egfa ((WW)/ [log P (Y|X; () + log P ((W+) - log Q ((W,)
(11)
L+1 log Q (Wol(
((w)
log (W)
(12)
The forms of the ELBO and approximate posterior suggest a sequential procedure to evaluate and subsequently optimise it: we alternate between sampling the weights using Eq. (9) and propagating the data and inducing locations using Eq. { (see Alg. [. In summary, the parameters of the
approximate posterior are the global inducing inputs, Uo. and the pseudo-data and precisions at (We (W.(W
all layers, (Ve, Ae | L+1 As each factor Q (W
} is a Gaussian, these parameters
can be optimised using standard reparameterised variational inference (Kingma & Welling 2013 Rezende et al.) 2014)) in combination with the Adam optimiser (Kingma & 2014) (Appendix A
Importantly, by placing inducing imputs on the training data (i.e. Uo X), and setting v = YA this
approximate posterior matches the optimal top-layer posterior (Eq. 60. Finally, we note that while
this posterior is conditionally Gaussian, the full posterior over all (We)
L+1 is non- non-Gaassian, and is
thus potentially more flexible than a full-covariance Gaussian over all weights.
Algorithm 1: Global inducing points for neural networks
Parameters: inducing inputs, Up. inducing outputs and precisions, (Ve, Ad) at all layers. Neural network inputs: (e.g. MNIST digits) F.
Neural network outputs: (e.g. classification logits) Fr++1
++ Do
for & {l...., + 1} do
Compute the mean and covariance over the weights at this layer
ΣW (Ne-1ST' (Ue-1)" AeO (Uc.))
M. DN 60( (U,1)) A.V.
Sample the weights and compute the ELBO
W (M,, EN = Q (we (We 1e
+ L Lllog P (We) -cl--log log N MW, EX)
Propagate the inputs and inducing points using the sampled weights,
U, = (U+11) W.
F. = ó (F(-1) We
C - C + -og P (Y L41)
3

---

Under review as a conference paper at ICLR 2021
2.2 EFFICIENT CONVOLUTIONAL BAYESIAN LINEAR REGRESSION
The previous sections were valid for a fully connected network. The extension to convolutional networks is straightforward in principle: we transform the convolution into matrix multiplica- tion by treating each patch as a separate input feature-vector, flattening the spatial and channel dimensions together into single vector. Thus, the feature-vectors have length in an__canneel x hernel_width x kernel_height and the matrix Uy contains patches_pper mage x minibatch patches. Likewise, we now have inducing outputs, v at each location in all the inducing images, so this again has length patches _ares_pe_______ image x min baiibatch. After explicitly extracting the patches, we can straightforwardly apply standard Bayesian linear regression.
However, explicitly extracting image patches is very memory intensive in DNN. If we consider standard convolution with a 3 x 3 convolutional kernel, then there pixel in the imput image, meaning factor of increase in memory consumption. Instead, we note that computing the matrices required for linear regression, & (U,-1)" Aro (U,-1) and (U,)² A V. does not require explicit extraction of image-patches. Instead, these matrices can be computed by taking the autocorrelation of the image/feature map, i.e. a convolution operation where we treat the image/feature map, as both the inputs and the weights (Appendix B for details).
3 patch centred each
2.3 PRIORS
We consider four priors in this work, which we refer to using the class names in the BNN library published with this paper. We careful to ensure that all parameters in the model have prior and approximate posterior, which is necessary to ensure that ELBOs are comparable across models.
First, we consider Gaussian prior with fixed scale, NealPrior, named because it is necessary obtain meaningful results when considering infinite networks (Neal 1996).
=
(13)
though bears strong similarities to the "He" initialisation (He et al. 2015). NealPrior is defined so to ensure that the activations retain sensible scale they propagate through the network. We compare this to the standard (0, 1) (StandardPPior)), which causes the activations to increase exponentially as they propagate through network layers (see Eq.
S N.
(14)
Next, we consider ScalePrior, which defines prior and approximate posterior over the scale,
S. -==1 II
P (se) = Gamma (sc; 2,
Q (se) = Gamma 2 + are 2 + (3))
(15)
where here we parameterise the Gamma distribution with the shape and rate parameters, and oz and 3g are non-negative learned parameters of the approximate posterior over Sg. Finally, we consider SpatiallWPriorr which allows for spatial correlations in the weights. In particular, take the covariance to the Kronecker product of an identity matrix over channel dimensions, and a Wishart-distrbuteed matrix, L over the spatial dimensions,
S ② L
P (L.) Inverse Wishart (Le) (N--1 + 1)I
,N-1 +1)
Q (Lę) InverseWishar (Lę (Ne-1 + I) N.N +1+v)
(16)
where the non-negative real number, v, and the positive definite matrix, V are learned parameters of the approximate posterior (see Appendix C).
EXTENSION TO DGPs
It is remarkable but underappreciated fact that BNNs are special cases of DGPs, with particular choice of kernel (Louizos Welling 2016 Aitchison 2019). Combining Eqs. and
1 No_ & (F S (F (F (17)
P (FF -II ( (F-.))) K (FFi))
4

---

Under review as a conference paper at ICLR 2021
local inducing
global inducing
Hamiltonian Monte Carlo
factorised
True Predicted
2
-2
2
-2
Figure 1: Predictive distributions on the toy dataset. Shaded regions represent one standard deviation.
Here, generalise our approximate posterior the DGP and link to the DGP literature. In DGP there are weights; instead we work directly with inducing outputs {U,)
P (U. I (uf O. (U
(18)
Note that here, we take the "global" inducing approach of using the inducing outputs from the previous layer, Uz-1 as the inducing inputs for the layer. this case, we need only learn the original inducing inputs, Un- This contrasts with the standard "local" inducing formulation, (as in Salimbeni & Deisenroth 2017), which learns separate inducing inputs at every layer, Ze-1. giving
P ZU = 1I ,V (u,, 0, K (Z,))).
As usual in DGPs (Salimbeni & Deisenroth 2017), the approximate posterior over Up induces an approximate posterior on F2 through the prior correlations. However, is important to remember that underneath the tractable distributions in Eqs. (77) and (18), there is infinite dimensional re-distributed function, Fe. such that Fe- Fe (-1) Standard local inducing point methods specify factorised approximate posterior over Fe by specifying the function's inducing outputs, U, F2 (Ze-1). )). at finite number of inducing input locations, Z₁-1- Importantly, the approximate posterior over a function, Fe, depends only on Ze-1. and Up. Thus, standard, local inducing. DGP approaches (c.g. Salimbeni & Deisenroth 2017), give a layervise-indeeenden approximate posterior over [Ft -+1 as they treat the inducing inputs, (Z,-1 1 as fixed, learned parameters and use inyerwise-independent approximate posterior over (Ue)
L-1 (Appendix F).
Next, need choose the approximate posterior on (U+)/=1. However, if our goal is intro- duce dependence across layers, seems inappropriate to use the standard weeewise-indepedeeee approximate posterior over (U L+1 Indeed, Appendix F we show that such a posterior implies functions in non-adjacent layers (e.g. F2 and Fe+2) are marginally independent, even with global inducing points.
To obtain more appropriate approximate posteriors, derive the optimal top-layer posterior for DGPs, which involves GP regression from activations propagated from lower layers onto the output data (Appendix E.. Inspired by the form of the optimal posterior again define approximate posterior by taking the product of the prior and a (pseudo-likelhhood
Q (U, x &h Q ( (v ( A7 Q (U, -1 (uf Σ A ET)
ΣN = -K- (U
(19)
where V and A are learned parameters, and our global inducing method, the inducing inputs, Ug-1, propagated from lower layers (Eq. T8) Importantly, setting the inducing inputs to the training data and v++ yx. the approximate posterior captures the optimal top-layer posterior for regression (Appendix Ep.. Under this approximate posterior, dependencies in U, naturally arise across all layers, and hence there dependencies between functions F3 all layers (Appendix F.

---

Under review as a conference paper at ICLR 2021
width 10
width 50
width 250
-0.25
-0.50
ELBO
factorised
-0.75
local inducing rand Iac
gl
1.00
global inducing
log ((y(x)
0
5
2
3
0
3
4
5
hidden layers
hidden layers
hidden layers
Figure 2: ELBO for different approximate posteriors as we change network depth/width on a dataset generated using linear Gaussian model. The rand gi line lies behind the global inducing line in width = 50 and width = 250.
In summary, we propose approximate posterior over inducing outputs that takes the form
L+1
Q ((U))))
II Q (U.U))))
(20)
As before, the parameters of this approximate posterior are the global inducing inputs, Up, and the pseudo-date and precisions at all layers, (Ve.Ar) 1411 Our ELBO, which we derive in takes the form
P (U. U
C = BQ(IP.UHH log (Y
log
(21)
We provide full description of method applied DGPs in Appendix D In concluding discussion of our proposed method, we note that we provide an analysis of the asymptotic complexity in Appendix L
4 RESULTS
We describe our experiments and results to assess the performance of global inducing points ('gi') against local inducing points ('li') and the fully factorised ( ('fac') approximation family. We addi- tionally consider models where we use one method to the last layer and another for the last layer, which may have computational advantages; we denote such models 'method1 -> method2'.
4.1 UNCERTAINTY IN 1D REGRESSION
We demonstrate the use of local and global inducing point methods a toy 1-D regression problem, comparing with fully factorised VI and Hamiltonian Monte Carlo (HMC; (Neal et al. 2011h) Following Hernándezz-Lbato & Adams 2015). we generate 40 input-output pairs (z, y) with the inputs v sampled i.i.d. from UQ-44-22
where ~ (O.033) We then normalised the inputs and outputs. Note that we have introduced 'gap' in the inputs, following recent work (Foong et al. 2019b] Yao et al., 2019 Foong et al.. (019a) that identifies the ability to express 'in-between' uncertainty as important quality of approximate inference algorithms. We evaluated the inference algorithms using fully-connected BNNs with hidden layers of 50 ReLU hidden units, using the NealPrior. For the inducing point methods, we used 100 inducing points per layer.
[2, 4) and the outputs generated by y -
The predictive distributions for the toy experiment can be seen in Fig. We observe that of the variational methods, the global inducing method produces predictive distributions closest to HMC, with good uncertainty in the gap. Meanwhile, factorised and local inducing fit the training data, but do not produce reasonable error bars: this demonstrates an important limitation of methods lacking correlation structure between layers.

---

Under review as a conference paper at ICLR 2021
bodion N-508, D-13
concrete N-1030, D-8
energy N-768, D
win N-1599, D-11
yacht N-508, D-6
BNN N(0, 1)
8NN Scale
ELBOs
DGP-2
D
BNN
A0 1)
Test LLs
BNN
Scale
DGP-2
DGP-5
-2.76
-2:30
-3.18
-3.07
-0.80
-0.67 -1.24
-0.84 -1.44
4.20
global inducing
Figure 3: Average test log likelihoods for BNNs the UCI datasets (in nats). Error bars represent one standard error. Shading represents different priors. connect the factorised models with the fac + gi models with thin grey line aid for casier comparison. Further to the right is better.
We provide additional experiments looking at the effect of the number of inducing points in Ap- pendix H and experiments looking at compositional uncertainty (Ustyuzhaninov et al. 2020 in both BNNs and DGPs for ID regression in Appendix
4.2
IN DEEP LINEAR NETWORKS
The lack of correlations between layers might be expected to become more problematic in deeper networks. To isolate the effect of depth on different approximate posteriors, we considered deep linear networks trained on data generated from toy linear model: input features were mapped to output feature, where the 1000 training and 100 test inputs are drawn IID from a standard Gaussian, and true outputs are drawn using weight- vector drawn IID from Gaussian with variance 1/5, and with noise variance of 0.1. We could evaluate the model evidence under the true data generating process which forms upper bound (in expectation) the model evidence and BLBO for all models.
We found that the ELBO for methods that factorise across layers - factorised and local inducing drops rapidly networks get deeper and wider (Fig. 2 This is undesirable behaviour, that wide, deep networks are necessary for good performance on difficult machine learning tasks. In contrast, we found that methods with global inducing points at the last layer decay much more slowly with depth, and perform better as networks get wider. Remarkably, global-inducing points gave good performance even with lower-layer weights drawn at random from the prior, which is not possible for any method that factorises across layers. We believe that fac - gi performed poorly at width 250 due to optimisation issues as rand >> gi performs better yet a special case of fac - gi.
know
4.3 REGRESSION BENCHMARK: UCI
We benchmark our methods on the UCI datasets in Hernández-Lobato & Adams (2015), popular benchmark regression datasets for BNNs and DGPs. Following the standard approach (Gal & Ghahramani 2015), each dataset uses 20 train-test 'splits" (except for protein with 5 splits) and the inputs and outputs are normalised to have zero mean and unit standard deviation. We focus on the five smallest datasets, App. MandM for all datasets). We consider two-layer fully-connected ReLU networks, using fully factorised and global inducing approximating families, well two- and five-layer DGPs with
expect Bayesian methods to be most relevant in small-data settings (see
7

---

Under review as a conference paper at ICLR 2021
Table 1: CIFAR-I0 classification accuracy. The first block shows our results using SpatialIWPrior, with ScalePrior in brackets. The next block shows comparable past results, from GPs and BNNs. The final block show con-comparable (sampling-basee based) methods. Dashes indicate that the figures were either not reported, are not applicable. The time is reported per epoch with ScalePrior and for MNIST, rather than CIFAR-I0 because known performance bug in the convolutions required in Sec. 2.2 with 32 x 32 (and above) images /tttps://giee cto wtt../yyotooooottooorrccccc issues 35603
test log like.
accuracy (%)
ELBO
time
factorised
-0.58 (-0.66)
80.27 (77.65)
-1.06 (-1.12)
19
local inducing
-0.62 (-0.60)
78.96 (79.46)
-0.84 (-0.88)
33
fac - gi
-0.49 (-0.56)
83.33 (81.72)
-0.91 (-0.96)
25 8
global inducing
-0.40 (-0.43)
86.70 (85.73)
-0.68 (-0.75)
65% s
Shi et (2019)
80.30%
Li et al.
2019
81.40%
Shridhar et al.
2019
73%
Wenzel et al.
2020
-0.35
88.50%
doubly-stochastic variational inference (DSVI) (Salimbeni & Deisenroth/ 2017) and global inducing. For the BNNs, we consider the standard N(00 1) prior and ScalePrior.
We display ELBOs and average test log likelihoods for the un-normalised data in Fig. where the dots and error bars represent the means and standard errors over the test splits, respectively. We observe that global inducing obtains better ELBOs than factorised and DSVI in every case, indicating that it does indeed approximate the true posterior better (since the ELBO is the marginal likelihood minus the KL to the posterior). While this the case for the ELBOs, this does not always translate to better test log likelihood due to model misspecification, as we see that occasionally DSVI outperforms global inducing by very small margin. The very poor results for factorised on ScalePrior indicate that it has difficulty learning useful prior hyperparameters for prediction, which is due to the looseness of its bound to the marginal likelihood. We provide experimental details, as well as additional results with additional architectures, priors, datasets, and RMSEs, in Appendices and M for BNNs and DGPs, respectively.
4.4 CONVOLUTIONAL BENCHMARK: CIFAR-10
For CIFAR-10, we considered a ResNet-inspired model consisting of conv2d-relu-blont htte--reluubooo-aaeeeeeee -avgool- hagpooo/-block-agggeoe
where the ResNet blocks consisted of shortcut connection in using 32 channels in all layers. In all our experiments, we
parallel with used no data augmentation and 500 inducing points. Our training scheme (see App. O ensured that our results did not reflect a 'cold posterior' (Wenzel et al. 2020). Our results are shown in Table T We achieved remarkable performance of 86.7% predictive accuracy, with global inducing points used for all layers, and with a spatial inverse Wishart prior on the weights. These results compare very favourably with comparable Bayesian approaches, i.e. those without data augmentation or posterior sharpening: past work with deep GPs obtained 80.3% (Shi et al. 2019), 9, and work using infinite-widhh neural networks define GP obtained 81.4% accuracy (Li inducing points we are approaching the accuracy of sampling based methods (Wenzel et al. 2020). which are in principle able to more closely approximate the true posterior. Furthermore, we see that global inducing performs the best in terms of ELBO (per datapoint) by wide margin, demonstrating that it gets far closer to the true posterior than the other methods. We provide additional results on uncertainty calibration in Appendix K
redd-relu-oov22d
2019). Remarkably, with only 500
RELATED WORK
Louizos & Welling (2016) attempted to pseudo-data along with matrix variate Gaussians to form an approximate posterior for BNNs; however, they restricted their analysis to BNNs, and is not clear how their method can be applied to DGPs. Their approach factorises across layers, thus
8

---

Under review as a conference paper at ICLR 2021
missing the important layerwise correlations that obtain. Moreover, they encountered an important limitation: the BNN prior implies that Ug is low-rank and it is difficult to design an approximate posterior capturing this constraint. As such, they were forced to use M < N2 inducing points, which particularly problematic the convolutional, global-inducing case where there are many patches (points) in each inducing image input.
Note that some work on BNNs reports better performance on datasets such as CIFAR-10. However, the best of our knowledge, no variational Bayesian method outperforms ours without modifying the BNN model or some form of posterior tempering (Wenzel et al. 2020), where the KL term in the relative the likelihood (Zhang et al. 2017: Bae al. 2018 Osawa et al.) 2019: Ashukha et al., 2020), which often increases the test accuracy. However, tempering clouds the Bayesian perspective, as the KI. to the posterior no longer minimised and the resulting objective is no longer lower bound on the marginal likelihood. By contrast, we use the untempered ELBO, thereby retaining the full Bayesian perspective. Dusenberry et al. (2020) report better performance on CIFAR-10 without tempering, but only perform variational inference over rank-1 perturbation to the weights, and maximise over all the other parameters. Our approach retains a full-rank parameterisation of the weight matrices.
ELBO fown-weights
Ustyuzhaninov et al. (2020) attempted to introduce dependence across layers in a deep GP by coupling inducing inputs to pseudo-outputs, which they term "inducing points inducing locations". However, as described above, the choice of approximate posterior over Ug is also critical. They used the standard approximate posterior that is independent across layers, meaning that while functions in adjacent layers were marginally dependent, the functions for non-adjacent layers were independent (Appendix F. By contrast, our approximate posteriors have marginal dependencies across U, and functions at all layers, and are capable of capturing the optimal top-layer posterior.
6 CONCLUSIONS
We derived optimal top-layer variational approximate posteriors for BNNs and deep GPs, and used them to develop generic, scalable approximate posteriors. These posteriors make use of global inducing points, which are learned only at the bottom layer and are propagated through the network. This leads to extremely flexible posteriors, which even allow the lower-layer weights to drawn from the prior. We showed that these global inducing variational posteriors lead to improved performance with better ELBOs, and state-of-the-art performance for variational BNNs on CIFAR-10.
REFERENCES
Laurence Aitchison. Why bigger is not always better: on finite and infinite neural networks. arXiv preprint arXiv:1910.08011 2019.
Arsenii Ashukha, Alexander Lyzhov, Dmitry Molchanov, and Dmitry Vetrov. Pitfalls of in-domain uncertainty estimation and ensembling in deep learning. arXiv preprint arXiv:2002.06477. 2020.
Juhan Bae, Guodong Zhang, and Roger Grosse. Eigenvalue corrected noisy natural gradient. arXiv preprint arXiv:1811.125555 2018
Maurice Stevenson Bartlett. On the theory of statistical regression. Proceedings of the Royal Society of Edinburgh, 53:260-283, 1933.
Charles Blundell, Julien Comebise, Koray Kavukcuoglu, and Daan Wierstra. Weight uncertainty in neural networks. arXiv preprint arXiv:1505.05424. 2015.
David R Burt, Carl Edward Rasmussen, and Mark van der Wilk. Convergence of sparse variational inference in gaussian processes regression. Journal of Machine Learning Research, 21:1-63, 2020.
Djalil
Chafai.
Bartlett URL
decomposition and other factorizations, Oct
2015.
http://ajaaa....aa.a.a..aoa.aaaa dee_lett-deomppsiiin and other fatereritttoiiiiiisssss
Andreas Damianou and Neil Lawrence. Deep Gaussian processes. In Artificial Intelligence and Statistics, pp. 207-215, 2013.

---

Under review as a conference paper at ICLR 2021
Michael W Dusenberry, Ghassen Jerfel, Yeming Wen, Yi-an Ma, Jasper Snoek, Katherine Heller, Balaji Lakshminarayanan, and Dustin Tran. Efficient and scalable bayesian neural nets with rank-1 factors. arXiv preprint arXivv2005.07186, 2020.
Vincent Dutordoir, Mark van der Wilk, Artem Artemev, Marcin Tomezak, and James Hensman. Trans- lation insensitivity for deep convolutional Gaussian processes. arXiv preprint arXiv: 19XXvv:902.0555 1902.05888, 2019.
Sebastian Farquhar, Lewis Smith, and Yarin Gal. Try depth instead of weight correlations: Mean-fiel field is less restrictive assumption for deeper networks. arXiv preprint arXiv:2002.03704 2020.
Andrew YK Foong, David R Burt, Yingzhen Li, and Richard E Turner. Pathologies of factorised Gaus- sian and MC dropout posteriors in Bayesian neural networks. arXiv preprint arXiv:1909.0071 2019a.
Andrew YK Foong, Yingzhen Li, José Miguel Hernández-Lobato, and Richard E Turner. 'in-between' uncertainty in Bayesian neural networks. arXiv preprint arXiv: /906. I 1537, 2019b.
Yarin Gal and Zoubin Ghahramani. Dropout Bayesian approximation: Representing model uncertainty in deep learning. arXiv preprint arXiv:150.002144 2015.
Alex Graves. Practical variational inference for neural networks. In Advances in neural information processing systems, pp. 2348-2356, 2011.
Chuan Guo, Geoff Pleiss, Yu Sun, and Kilian Q Weinberger. On calibration of modern neural networks. In Proceedings of the 34th International Conference on Machine Learning- Volume 70, pp. 1321-1330 JMLR. org, 2017.
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Delving deep into rectifiers: Surpassing human-level performance on ImageNet classification. In Proceedings of the IEEE international conference on computer vision, pp. 1026- 1034, 2015.
José Miguel Hernnnezz-Lobat and Ryan Adams. Probabilistic backpropagation for scalable learning of Bayesian neural networks. In International Conference on Machine Learning, pp. 1861-1869, 2015.
Geoffrey E Hinton and Drew Van Camp. Keeping the neural networks simple by minimizing the description length the weights. In Proceedings of the sixth annual conference on Computational learning theory, pp. 5-13, 1993.
Michael Jordan, Zoubin Ghahramani, Tommi Jaakkola, and Lawrence Saul. An introduction to variational methods for graphical models. Machine learning, 37(2): 1999.
Diederik P Kingma and Jimmy Ba. Adam: method for stochastic optimization. arXiv preprint arXiv: 1412:v:412..980.. 2014.
Diederik P Kingma and Max Welling. Auto-encoding variational bayes. arXiv:1312.61 14, 2013.
arXiv preprint
Durk Kingma, Tim Salimans, and Max Welling. Variational dropout and the local reparameterizztioo trick. In Advances in neural information processing systems, pp. 2575-2583, 2015.
Alex Krizhevsky, Geoffrey Hinton, et al. Learning multiple layers of features from tiny images. Technical report, University of Toronto, 2009.
David Krueger, Chin-Wei Huang, Riashat Islam, Ryan Turner, Alexandre Lacoste, and Aaron Courville. Bayesian hypernetworks. arXiv preprint arXiv: 2017.
Zhiyuan Li, Ruosong Wang, Dingli Yu, Simon S Du. Wei Hu, Ruslan Salakhutdinov, and Sanjeev
Arora. Enhanced convolutional neural tangent kernels. arXiv preprint 191.0009 2019.
Christos Louizos and Max Welling. Structured and efficient variational deep learning with matrix Gaussian posteriors. In International Conference Machine Learning, pp. 1708-1716, 2016.
10

---

Under review as a conference paper at ICLR 2021
Mahdi Pakdaman Naeini, Gregory Cooper, and Milos Hauskrecht. Obtaining well calibrated prob- abilities using bayesian binning. In Twenty-Ninth AAAI Conference on Artificial Intelligence, 2015.
Radford M Neal. Priors for infinite networks. In Bayesian Learning for Neural Networks, pp. 29-53. Springer, 1996.
Radford M Neal et al. MCMC using Hamiltonian dynamics. Handbook of Markov chain Monte Carlo, 2(11):2, 2011.
Kazuki Osawa, Siddharth Swaroop, Mohammad Emtiyaz E Khan, Anirudh Jain, Runa Eschenhagen, Richard E Turner, and Rio Yokota. Practical deep learning with Bayesian principles. In Advances in Neural Information Processing Systems, pp. 4289-4300 2019.
Nick Pawlowski, Andrew Brock, Matthew CH Lee, Martin Rajchl, and Ben Glocker. Implicit weight uncertainty in neural networks. arXiv preprint arXiv:1711.0129) 2017.
Judea Pearl. Probabilistic Reasoning in Intelligent Systems: Networks of Plausible Inference. Morgan Kaufmann, 1988.
Carl Edward Rasmussen and Christopher KI Williams. Gaussian Processes for Machine Learning. ISBN-13 978-0-262-18253-9, 2006.
Danilo Jimenez Rezende, Shakir Mohamed, and Daan Wierstra. Stochastic backpropagation and approximate inference in deep generative models. arXiv preprint arXiv:1401.4082, 2014.
Hippolyt Ritter, Aleksandar Botev, and David Barber. A scalable Laplace approximation for neural networks. In 6th International Conference on Learning Representations, ICLR 2018-Core Conference Track Proceedings, volume International Conference on Representation Learning. 2018.
Hugh Salimbeni and Mare Deisenroth, Doubly stochastic variational inference for deep Gaussian processes. In Advances in Neural Information Processing Systems, pp. 4588-4599, 2017.
Jiaxin Shi, Michalis K Titsias, and Andriy Mnih. Sparse orthogonal variational inference for Gaussian processes. arXiv preprint
2019.
Kumar Shridhar, Felix Laumann, and Marcus Liwicki. A comprehensive guide to Bayesian convolu- tional neural network with variational inference. arXiv preprint arXiv:1901.022 2019.
Tijmen Tieleman and Geoffrey Hinton. Lecture 6.5-rmsprop: Divide the gradient by a running average of its recent magnitude. COURSERA: Neural networks for machine learning, 4(2):26-31 31, 2012.
Brian Trippe and Richard Turner. Overpruning in variational bayesian neural networks, arXiv preprint 1rXiv:1801.02555 2018
Ivan Ustyuzhaninov, leva Kazlauskaite, Markus Kaiser, Erik Bodin, Neill DF Campbell, and Carl Hen- rik Ek. Compositional uncertainty in deep Gaussian processes. UAI, 2020.
Florian Wenzel, Kevin Roth, Bastiaan S Veeling, Jakub Swiętkowski, Linh Tran, Stephan Mandt, Jasper Snoek, Tim Salimans, Rodolphe Jenatton, and Sehastian Nowozin. How good is the Bayes posterior in deep neural networks really? arXiv preprint arXiv:2002.02405, 2020.
Jiayu Yao, Weiwei Pan, Soumya Ghosh, and Finale Doshi- Velez. Quality of uncertainty quantification for Bayesian neural network inference. arXiv preprint (rXXiv/906.096666 2019.
Guodong Zhang, Shengyang Sun, David Duvenaud, and Roger Grosse. Noisy natural gradient as variational inference. arXiv preprint arXiv:1712.023990 2017.
11

---

Under review as a conference paper at ICLR 2021
A REPARAMETERISED VARIATIONAL INFERENCE
In variational inference, the ELBO objective takes the form,
( = (Q,(w)) log P (Y[X, w) + log P (w)
(22)
Qs (w)
where w a vector containing all of the elements of the weight matrices in the full network, {We) W+t and ó (Ug, (Ve. Ar) Ar) )++ are the parameters of the approximate posterior. This objective is difficult to differentiate wrt Φ. because á parameterises the distribution over which the
expectation is taken. Following Kingma & Welling (2013) and Rezende et al. (2014), we sample from a simple, fixed distribution (e.g. IID Gaussian), and transform them to give samples from Q (w).
w(c:o)~ ~ 2 (w).
(23)
Thus, the ELBO can be written,
P (ww(:: 6))
(()) = E. log (Y]X, w(e: Ø)) + log
(24)
Qo (w (w(c:
As the distribution over which the expectation is taken is now independent of ø. we can form unbiased estimates of the gradient of [(ő) by drawing one or a few samples of
B EFFICIENT CONVOLUTIONAL. LINEAR REGRESSION
Working in one dimension for simplicity, the standard form for a convolution in deep learning is
Yiu.o Σ
Xicu+bie oe
(25)
where X the input image/feature-map is the output feature-map, W is the convolutional weights, i indexes images, and c' index channels, indexes the location within the image, and indexes the location within the convolutional patch. Later, we will swap the identity of the "patch location" and the "image location" and facilitate this, define them both to be centred on zero,
u E (- 1)/2, (W 1)/2}
8 {- 1//18-11/...... 1)/2...., (S - 1)/2}
(26)
where W is the size of an image and S is the size of a patch, such that, for example for a size kernel, 5 E (-1,0,1)
To use the expressions for the fully-connected case, we can form a new input, X' by cutting out each image patch,
X/u.cc = Xi(u+8)es
(27)
leading to
Yiu,c
Σ XX Wodtee 5c'
(28)
Note that we have used commas to group pairs indices (iu and 8c') that may be combined into single index (e.g. using a reshape operation). Indeed, combining i and into single index and combining and e' into a single index, this expression can be viewed as standard matrix multiplication,
Y X'W.
(29)
This means that we can directly apply the approximate posterior we derived for the fully-connected case in Eq. (9) to the convolutional case. To allow for this, we take
X' = = d (U,-1).
Y V
(30)
While we can explicitly compute X' by extracting image patches, this imposes a very large memory cost (a factor of for 3 x kernel, with stride 1, because there are roughly many patches as pixels in the image, and a 3 x 3 patch requires times the storage a pixel. To implement convolutional linear regression with more manageable memory cost, we instead compute the matrices required
12

---

Under review as a conference paper at ICLR 2021
for linear regression directly convolutions of the imput feature- maps. X with themselves, and as convolutions of the X with the output feature maps, Y, which we describe here.
For linear regression (Eq. 9 first need to compute,
s (UL-)) A.V,)
X (XTY))
=E
Xe Vin,e--
(31)
óe,c'
60,e²
Rewriting this in terms of X (i.e. without explicitly cutting out image patches), we obtain,
X (x"Y)
=E iM
Xi(u+tje Viaee'
(32)
This can be directly viewed as the convolution X and Y, where we treat Y as the "convolutional weights", the location within the now very large (size W) "convolutional patch", and ó as the location in the resulting output. Once we realise that the computation is a spatial convolution, is possible to fit it into standard convolution functions provided by deep-learning frameworks (albeit with some rearrangement of the tensors).
Next, need to compute,
Aeo (U.-))) = (x"x') (o (U
6c,8'e =
(33)
iu
Again, rewriting this in terms X (i.e. without explicitly cutting out image patches), we obtain,
X 'X'x''
Σ Xilut
(34)
To treat this as convolution, we first need exact translational invariance, which can be achieved by using circular boundary conditions. Note that cireular boundary conditions not typically used in neural networks for images, and we therefore only use circular boundary conditions to define the approximate posterior over weights. The variational framework does not restrict us to also using circular boundary conditions within our feedforward network, and as such, we use standard padding. With exact translational invariance, we can write this expression directly as convolution,
(x'x')
= Xieex
(35)
Se,8'd
where
- 8) E (-(S-)),.. 1)}
(36)
i.c. for a size 3 kernel, -) E (-2, 1,, where we treat Xiue as the "convolutional weights", u as the location within the "convolutional patch", and 8' 8 as the location in the resulting output "feature-map"
Finally, note that this form offers considerable benefits terms memory consumption. In particular, the output matrices are usually quite small - the number of channels is typically 32 or 64, and the number of locations within patch typically 9, giving a very manageable total size that typically smaller than 1000 1000.
C WISHART DISTRIBUTIONS WITH REAL- VALUED DEGREES OF FREEDOM
The classical description of the Wishart distribution,
Σ Wishart (I,)),
(37)
where Σ is P x P matrix, states that P 2 1 is an integer and we can generate Σ by taking the "x generated IID from standard Gaussian,
product of matrices.
E X
X N (0,1).
(38)
However, for the purposes of defining learnalle approximate posteriors, we need to be able to sample and evaluate the probability density when positive real.
13

---

Under review as a conference paper at ICLR 2021
To do this, consider the alternative, much more efficient means of sampling from Wishart distribution, using the Bartlett decomposition (Bartlett (933). The Bartlett decomposition gives the probability density for the Cholesky of a Wishart sample. In particular,
Tn
Tim
T
(39)
Tmm
(T])) = Gamma (T)) - 2)
(40)
12)
P (Tzck) = N 0T
(41)
Here, T5 is usually considered to be sampled from a x² distribution, but we have generalised this slightly using the equivalent Gamma distribution to allow for real-valued Following Chafai (0155 We need to change variables to Tjj rather than T
P = P(
(42)
Gamma (T); ø-j+1
2 :2)
2T
(43)
(T)) (V-1+1)/2-1 -
23
(44)
1/-1++1//21 ("
Th e
(45)
2 21w-j-1)/21 ("=4t1)
Thus, the probability density for T under the Bartlett sampling operation is
To 12,12
P (T) = リーj-1 kE()+I...,m)
II
e
T32/2
(46)
2
I (1
2
on-diagonals
off-diagonals
(47)
To convert this a distribution on Σ, we need the volume element for the transformation from T to
Σ,
T70
dΣ
2
m
dT,
(48)
35
J
which can be obtained directly by computing the log-determinant of the Jacobian for the transforma- tion from to E, or by taking the ratio of Eq. (46) and the usual Wishart probability density (with integral v). Thus,
P
2
Tmm j+1
(49)
j=1
Tw-m-
1 T2
II
T3./2
(50)
v-j+1
2
r (الهt)
kE(J+1....m)
22
14

---

Under review as a conference paper at ICLR 2021
Breaking this down into separate components and performing straightforward algebraic manipula- tions,
Tw-m-
(T T
(v-m-1))
(51)
P Σ
12
II
e e
(52)
kE (j+
Eij Tjhl2 Tr(E1/2
(53)
II 2112 1)/2
II
1
(54)
V2
211-1-18/2
2mm-1-1y/2
kE m
1
2(u-j+11/2
20/11//2
(55)
II 2 m2/ =2 22
(56)
where the second crooddto-laree last line was obtained by noting that j - covers the same range of integers as m j under the product. Finally, using the definition of the multivariate Gamma function,
II
V (m 1)/4 II ( - (5)
(57)
kE(j+I....m)
We thereby re-obtain the probability density for the standard Wishart distribution,
(V- 1)/2 11/22 E Tr(E)/2
P P
(58)
2ma/2T m (%)
D FULL DESCRIPTION OF OUR METHOD FOR DEEP GAUSSIAN PROCESS
Here we give the full derivation for our doubly-stochastic variational deep GPs, following Salimbeni & Deisenroth (2017). A deep Gaussian process (DGP; (Damianou & Lawrence, 2013 Salimbeni Deisenroth 2017 defines a prior over function values, F. E RT'NN, where is the layer, P is the number of input points, and N₂ is the "width" of this layer, by stacking L + 1 layers of standard Gaussian processes (we use L + layers instead of the typical L. layers to retain consistency with how define BNNs):
P PYY
L+1 IY (F.
(59)
Here, the input is Fo = X, and the output is Y (which could be continuous values for regression, or class labels for classification), and the distribution over F, PxNN, factorises into independent multivariate Gaussian distributions over each function,
P (F = I P (EX
(E.
(60)
where I is the Ath column of Fr, giving the activation of all datapoints for the Xth feature, and K a function that computes the kernel-matrix from the features in the previous layer.
(.)
To define a variational approximate posterior, we augment [Fr) with inducing points consisting of the function values (U
P (Y, (Fe, U U,) = P (Y (Fr,, U,-))
(61)
and because F, and U, are the function outputs corresponding to different inputs (Fr-1 and U²-1), they form joint multivariate Gaussian distribution, analogous to Eq. (60).
P (Fr F. U--1) = II W
(( (5))
0. ((D((()
(62)
==
15

---

Under review as a conference paper at ICLR 2021
Following Salimbeni & Deisenroth (2017) we form approximate posterior by conditioning the function values, Fe on the inducing outputs, Uz,
Q ((F2, X. U.) I Q (F& F-1, Ue-1)
= - P UF (e-
(63)
where QU,,(UZ-I) is given by Eq. i.e.
Q QU. II (u E A .
= - (U,-1) + AA)
-1
(64)
and P UF, 1 F is given by standard manipulations of Eq. 62 can be factorised in analogous fashion,
Importantly, the prior
P (F& U. U2-1) = (U.
(65)
The full ELBO can be written as
P(Y,FF x, Up) log Q ( U
[ = E. ((P
(66)
1++1 X, Un)
Substituting Eqs.
61
663
and
66
into Eq. (66). the P (Fel Ug U2-1, F (-1) terms cancel and
obtain.
PUUUU Q (U,
L C-
log P (Y|)++1))+ log
(67)
Analogous to the BNN case, we evaluate the ELBO by alternatively sampling the inducing function values Ue given U and propagating the data by sampling from P [F. Fx-1) (see Alg. 2 As in (allimenni & Deisenroth (2017), since all likelihoods we use factorise across datapoints, we only need to sample from the marginals of the latter distribution to sample the propagated function values. As in the BNN case, the parameters of the approximate posterior are the global inducing inputs, Uo. and the pseudo-data and precisions at all layers, (ee. Ark We can similarly use standard reparameterised variational inference (Kingma & Welling 2013 Rezende et al. (014] to optimise these variational parameters, as Q (U, is Gaussian. We use Adam (Kingma & Ba 2014) to optimise the parameters.
Algorithm 2: Global inducing points for deep Gaussian processes
inducing inputs, Up. inducing outputs and precisions, (Ve. A at all layers. DGP inputs: (e.g. MNIST digits) Fo
Parameters:
DGP outputs: (e.g. classification logits) FL+1
C+0
for E ......,,+. do
Compute the mean and covariance over the inducing outputs at this layer
E3 (K- (U,-1) Ac)
M. Σ4 A.
Sample the inducing oulputs and compute the ELBO
U. N' (M,, Em)= = Q
(UI)
log N (U, E Propagate the imputs using the sampled inducing outpuss,
t C log P (U,(U,-1)
F, ~ P (Fe Ue, U,-1, FE-1)
L + L+ log (Y F
E MOTIVATING THE APPROXIMATE POSTERIOR FOR DEEP GPs
Our original motivation for the approximate posterior was for the BNN case, which we then extended deep GPs. Here, show how the same approximate posterior be motivated from deep GP
16

---

Under review as a conference paper at ICLR 2021
perspective. As with the BNN case, first derive the form of the optimal approximate posterior for the last layer, in the regression case. Without inducing points, the ELBO becomes
P ( (F
C ((((
log
(68)
Q ((Fe 1++1
where we have defined generic variational posterior Q (fFe) Since we are interested in the form of Q (F (Fe))=)) we rearrange the ELBO so that all terms that do not depend on F L++ are absorbed into a constant, c:
L = E ((PO) [log P
- log Q (FL-1)+
(69)
Some straightforward rearrangements lead to a similar form to before,
C = BQ((Fe
(log [logP (FI
-x (Q (FL+i ||P (F (Ft) + tc]. (70)
from which we see that the optimal conditional posterior is given by Q (F1z_1(FF))) = Q (FL+)(Fr) P (Fz+1)) F, L). which has a closed form for regression: it is simply the standard GP posterior given by training data Y at inputs F L.. In particular, for likelihood
P PY = I W 1 AZt1),
(71)
where A 144 is the precision.
Q (FL+i = Ie Not1 N 1 A 1+1 E EL+1), Et - = (K- (FL) AL+1)
(72)
(73)
Finally, as usual in GPs (Rasmussen & Williams, 2006), the predictive distribution for test points can be obtained by conditioning using this approximate posterior.
F COMPARING OUR DEEP GP APPROXIMATE POSTERIOR TO PREVIOUS WORK
The standard approach to inference in deep GPs (e.g. Salimbeni & Deisenroth 2017 involves local inducing points and an approximate posterior over {U+)2T that is factorised across layers,
Q ((U+)2+) - NH (/,
(74)
This approximate posterior induces an approximate posterior over the underlying crinite-dimensione functions Fz at each layer, which are implicitly used to propagate the data through the network via F2 = Fe We show a graphical model summarising the standard approach in Fig. Ala
While, Salimbeni & Deisenroth (011)) point out, the function values [Fr) L++1 are correlated, the functions [ themselves are independent across layers. We note that for BNNs, this is equivalent to having a posterior over weights that factorises across layers.
One approach to introduce dependencies across layers for the functions would be to introduce the notion of global inducing points, propagating the initial Uo through the model. In fact, as we note in the Related Work section, (Ustyuzhaninov et al. (2020) independently proposed this approach to introducing dependencies, using toy problem to motivate the approach. However, they kept the form of the approximate posterior the same as the standard approach (Eq. we show the corresponding graphical model in Fig. A1b The graphical model shows that adjacent functions F, and Fe+1 share the parent node Ug, they are fact dependent. However, non-adjacent functions do not share any parent nodes, and are independent; this can be by considering the d-separation criterion [Pear] 1988) for F and Fe+1. which have parents (Ug-2, U²-1) and (U,, U respectively.
Our approach, by contrast, determines the form of the approximate posterior over U, by performing Bayesian regression using U--1 as input to that layer's 's GP, where the output data is V.. This results in posterior that depends the previous layer, Q (-=) We show the corresponding graphical model in Fig. Ale From this graphical model it is straightforward to see that our approach results in posterior over functions that correlated across all layers.
17

---

Under review as a conference paper at ICLR 2021
Up-1
U4+1
Fe-1
Ft
Fe+1
F
Fe
F E+1
(a) Graphical model illustrating the approach of Salimbeni Deisenroth (017).
The inducing inputs,
(Uc-1) 1++1 are treated as learned parameters and are therefore omitted from model.
U-11
U
U4+1
Fi-1
Ft
Fe+1
F
F2
Fi+1
(b) Graphical model illustrating the approach of Ustyuzhaninov al. 2020 the inducing outputs at the previous layer.
The inducing inputs given by
U
U²
U
F
Ft
Fe+1
Fi-1
Fr
Fi+1
(c) Graphical model illustrating approach.
Figure A1: Comparison of the graphical models for three approaches to inference in deep GPs:
Salimbeni & Deisenroth (2017). Ustyuzhaninov et al. (2020), and ours.
G PARAMETER SCALING FOR ADAM
The standard optimiser for variational BNNs is ADAM (Kingma & Ba, 2014). which we also use. Considering similar RMSprop updates for simplicity (Tieleman & Hinton 2012).
AW
(75)
E [
where the expectation over g² approximated using a moving-average of past gradients. Thus, absolute parameter changes going be of order y.. This is fine if all the parameters have roughly
18

---

Under review as a conference paper at ICLR 2021
2 points ELBO: -2.244
10 points ELBO: 2.055
20 points ELBO: 2. 2.125
40 points ELBO: 2.179
True Predicted
0
0
0
0
-2
-2
-2
-2
Figure A2: Predictive distributions on the toy dataset as the number of inducing points changes.
the same order of magnitude, but becomes serious problem if some of the parameters are very large and others are very small. For instance, if parameter is around 10-4 and 10-1 then a single ADAM step can easily double the parameter estimate, or change it from positive to negative. In contrast, a parameter around then ADAM, with 17 can make proportionally much smaller changes to this parameter, (around 0.01%). Thus, we need to ensure that all of our parameters have the same scale, especially as we mix methods, such as combining factorised and global inducing points. We thus design all our new approximate posteriors (i.e. the inducing inputs and outputs) such that the parameters have a scale around 1. The key issue is that the mean weights in factorised methods tend to be quite small - they have scale around 1/V fan-in. To resolve this issue, we store scaled weights, and we divide these stored, scaled mean parameters by the fan-in as part of the forward pass.
weights = scaled weights
(76)
Vfan-in
This scaling forces us to use larger learning rates than are typically used.
H EXPLORING THE EFFECT OF THE NUMBER OF INDUCING POINTS
In this section, we briefly consider the effect of changing the number of inducing points, M, used in global inducing. We reconsider the toy problem from Sec. 4.1 and plot predictive posteriors obtained with global inducing as the number of inducing points increases from 2 to 40 (noting that in Fig. we used 100 inducing points).
We plot the results of our experiment in Fig. A2 While two inducing points are clearly not sufficient, we observe that there remarkably very little difference between the predictive posteriors for or more inducing points. This observation is reflected in the ELBOs per datapoint (listed above each plot), which show that adding more points beyond 10 gains very little in terms closeness to the true posterior.
However, we note that this a very simple dataset: consists of only two clusters of close points with very clear trend. Therefore, would expect that for more complex datasets more inducing points would be necessary. We leave a full investigation of how many inducing points are required to obtain suitable approximate posterior, such as that found in Burt et al. (2020) for sparse GP regression, to future work.
I UNDERSTANDING COMPOSITIONAL UNCERTAINTE
In this section, we take inspiration from the experiments of Ustyuzhaninov et al. (2020) which investigate the compositional uncertainty obtained by different approximate posteriors for DGPs. They noted that methods which factorise over layers have a tendency to cause the posterior distribution for each layer to collapse (nearly) deterministic function, resulting in worse uncertainty quantification within layers and worse ELBOs. In contrast, they found that allowing the approximate posterior to have correlations between layers allows those layers to capture uncertainty, resulting in better
19

---

Under review as a conference paper at ICLR 2021
(x
F2(()
local inducing
ELBO: 4 4.921
global inducing ELBO: 7 7.483
2
0
0
0
0
1
X
Figure A3: Posterior distributions for 2-1ayer DGPs with local inducing and global inducing. The first two columns show the predictive distributions for each layer taken individually, while the last column shows the predictive distribution of the output y.
Table 2: ELBOs and variances of the intermediate functions for BNN fit to the toy data of Fig.
ELBO
V[Fi]
V[F2
V[Fy
factorised
-4.585
0.0728
0.4765
0.1926
local inducing
-5.469
0.0763
0.4473
0.0643
global inducing
2.236
0.4820
0.4877
1.0820
ELBOs and therefore a closer approximation to the true posterior. They argue that this then allows the model to better discover compositional structure in the data.
first consider a toy problem consisting of 100 datapoints generated by sampling from a two-layer DGP of width one, with sguared-exponential kernels in each layer. We then fit two two-layer DGPs to this data one using local inducing, the other using global inducing. The results of this experiment be seen in Fig. A3 which show the final fit, along with the learned posteriors
intermediate
functions Fi and F22 These results mirror those observed by Ustyuzhaninov et al. (2020) on a similar nearly deterministic posterior
experiment: local inducing, which factorises over layers, collapses over the intermediate functions, whereas global inducing provides a much broader distribution over functions for the two layers. Therefore, global inducing leads wider range of plausible functions that could explain the data via composition, which can be important in understanding the data. We observe that this behaviour directly leads better uncertainty quantification for the out-of diiitri-iititttiinn region, well as better ELBOs.
To illustrate a similar phenomenon in BNNs, we reconsider the toy problem of Sec. As is not meaningful to consider neural networks with only one hidden unit per layer, instead of plotting intermediate functions we instead look at the mean variance the functions at random input points, following roughly the experiment Dutordoir et al. (2019) in Table 1. For each layer, we consider the quantity
N;
E
v(d((z))
(77)
N₂
20

---

Under review as a conference paper at ICLR 2021
boston D-13
comorelo 030,
enorgy
kinōom
DOWer
Wane
yash!
N/0,11 1)
ELBOs
Neal
Scale
-1.44
-0.69
-1.36
-0.66
0.30
0.88
-0.48
-0.99
2.11
2.82
-0.11
-0.02
-1.13
-1.06 -1445
-1.25
-0.85
0.58
N/00.11 1)
Neal
-.80
-2.50
3.87
-3.18 0.79
-0.87
1.08
1.27
6.70
7.41
-2.44
-2.79 22882
2.81 -1.21
-0.98
-2.33
-0.86
tactorised
inducing
global inducing
gi
Figure A4: ELBOs per datapoint and average test log likelihoods for BNNs on UCI datasets.
where the expectation is random input points, which sample from standard normal. We expect that for methods which introduce correlations across layers, this quantity will be higher, as there will wider range of intermediate functions that could plausibly explain the data. We confirm this in Table 2 which indicates that global inducing leads to many more compositions of functions being considered as plausible explanations of the data. This is additionally reflected in the ELBO, which is far better for global inducing than the other, factorised methods. However, we note that the variances are far closer than we might otherwise expect for the second layer. We hypothesise that this is due to the pruning effects described in Trippe & Turner (2018), where layer has many weights that are close to the prior that are then pruned out by the following layer by having the outgoing weights collapse to zero. In fact, note that the variances in the last layer are small for the factorised methods, which supports this hypothesis. By contrast, global inducing leads to high variances across all layers.
We believe that understanding the role compositional uncertainty in variational inference for deep Bayesian models lead to important conclusions about both the models being used and the compositional structure underlying the data being modelled, and therefore an important direction for future work to consider.
J UCI RESULTS WITH BAYESIAN NEURAL NETWORKS
For this Appendix, we consider all of the UCI datasets from Hernández-Lobato & Adams (015). along with four approximation families: factorised (i.e. mean-fin -ield., local inducing, global inducing. and fac-*gi, which may offer some computational advantages to global inducing. We also considered three priors: the standard (0, prior, NealPrior, and ScalePrior. The test L.Ls and ELBOS for BNNs applied to UCI datasets are given in Fig. A4 Note that the ELBOs for the global inducing methods (both global inducing and fac- are almost always better than those for baseline methods, often by very large margin. However, as noted earlier, this does not necessarily correspond to better test log likelihoods due to model misspecification: there is not straightforward relationship between the ELBO and the predictive performance, and it is possible obtain better test log likelihoods with worse inference. We present all the results, including for the test RMSEs, in tabulated form in Appendix p
21

---

Under review as a conference paper at ICLR 2021
factorised
local inducing
fac ->g1 gi
global inducing
1,0
1.0
1.0
1.0
0.8
0.8
0.8
0.8
accurady
0.6
0.6
0.4
0.4
0.4
0.4
SpatialW Scale
0.0
0.0
0.0
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.4
0.6
0.8
1.0
contidence
confidence
contidence
confidance
Figure A5: Calibration curves for CIFAR-10
J.1 EXPERIMENTAL DETAILS
The architecture we considered for all BNN UCI experiments were fully-connected ReLU networks with 2 hidden layers of 50 hidden units each. We performed a grid search to select the learning rate and minibatch size. For the fully factorised approximation, we selected the learning rate from (3e-4, 1e-3, 3e-3, 1e-2) and the minibatch size from {32, 100, 500), optimising for 25000 gradient steps; for the other methods we selected the learning rate from {3e-3, 1e-2] and fixed the minibatch size to 10000 (as in Salimbeni & Deisenroth (2017)) optimising for 10000 gradient steps. For all methods we selected the hyperparameters that gave the best EL.B.. We trained the models using 10 samples from the approximate posterior, while using 100 for evaluation. For the inducing point methods, we used the selected batch size for the number of inducing points per layer. For all methods, we initialised the log noise variance at -3, but use the scaling trick in Appendix G to accelerate convergence, scaling by factor of 10. Note that for the fully factorised method we used the local reparameterisation trick (Kingma et al. 2015): however, for fac + gi we cannot do so because the inducing point methods require that covariances be propagated through the network correctly. For the inducing point methods, we additionally use output channel-specific precisions, AI,x which effectively allows the network to prune unnecessary neurons if that benefits the ELBO. However, only parameterise the diagonal of these precision matrices to save on computational and memory cost.
K UNCERTAINTY CALIBRATION FOR CIFAR-10
To assess how well our methods capture uncertainty, we consider calibration. Calibration is assessed by comparing the model's probabilistic assessment of its confidence with its accuracy - the pro- portion of the time that it is actually correct. For instance, gathering model predictions with some confidence (e.g. softmax probabilities in the range 0.9 to 0.95), and looking at the accuracy of these predictions, we would expect the model to be correct with probability 0.925; higher lower value would represent miscalibration.
We begin by plotting calibration curves in Fig. A5 obtained by binning the predictions in 20 equal bins and assessing the mean accuracy of the binned predictions. For well-calibrated models, we expect the line to lie on the diagonal. A line above the diagonal indicates the model is underconfident (the model is performing better than it expects), whereas a line below the diagonal indicates it is overconfident (it performing worse than it expects). While is difficult to draw strong conclusions from these plots, it appears generally that factorised is poorly calibrated for both priors, that Spatiall generally improves calibration over ScalePrior, and that local inducing with SpatialIWPrior performs very well.
To come to more quantitative conclusions, we use expected calibration error (ECE; (Naeini et al. 2015; Guo et al.) (017)), which measures the expected absolute difference between the model's confidence and accuracy. Confirming the results from the plots (Fig. A5). we find that using the more sophisticated SpatialIWPrior gave considerable improvements in calibration. While, expected, we find that our most accurate prior, Spatial[WPrior, in combination with global inducing points did very well (ECE of 0.021), the model with the best ECE is actually local inducing with PaatiallWPProo.. albeit by a very small margin. We leave investigation of exactly why this is to future work. Finally, note our final ECE value of 0.021 is considerable improvement over those for uncalibrated models
22

---

Under review as a conference paper at ICLR 2021
Table Expected calibration error for CIFAR-10
factorised
local inducing
fac - gi
global inducing
ScalePrior
0.053
0.040 0.018
0.049
0.038
SpatialIWPPior
0.045
0.036
0.021
in Guo al. 2017) (Table 1), which are in the region of 0.03-0.045 (although considerably better calibration can be achieved by post-hoc scaling of the model's confidence).
L ASYMPTOTIC COMPLEXITY
In the deep GP case, the complexity for global inducing is exactly that of standard inducing point Gaussian processes, i.e. O(M* + PM2) where M is the number of inducing points, and P can be taken to be the number of training inputs, or the size of minibatch, as appropriate. The first term, M comes from computing and sampling the posterior over Uę based on the inducing points (e.g. inverting the covariance). The second term, and PM2 comes from computing the implied distribution over Fz.
In the fully-connected BNN case, we have three terms, O(N3 + M + PN3). The first term, N where N corresponds to the width the network, arises from taking the inverse of the covariance matrix in Eq. (9). but is also the complexity e.g. for propagating the inducing points from layer the next (Eq. The second term, MN2 comes from computing that covariance Eq. 2 by taking product of input features with themselves. The third term PN2 comes from multiplying the training inputs/minibatch by the sampled inputs (Eq. I
M UCI RESULTS WITH DEEP GAUSSIAN PROCESSES
In this appendix, we again consider all of the UCI datasets from Hernández-Lobato & Adams (2015)
for DGPs with depths ranging from two five layers. We compare DSVI (Sai(maleeeen)) & Deisenroth 2017). local inducing, and global inducing. While local inducing uses the same inducing ppint architecture as /Salimbeni Deisenroth (2017), the actual implementation and parameterisation is very different and in addition they used complex architecture involving skip connections that we may not have matched exactly. As such, we do expect see differences between local inducing and Salimbeni & Deisenroth (2017).
We show our results in Fig. A6 Here, the results are not as clear-cut as the BNN case. For the smaller datasets (i.e. boston, concrete, energy, wine, yacht), global inducing generally outperforms both local inducing and DSVI, as noted the main text, especially when considering the ELBOs. We do however observe that for power, protein, and one model for kin8nm, the local approaches sometimes outperform global inducing, even for the ELBOs. We believe this is due to the fact that relatively few inducing points far fewer variational parameters than the local approaches. This may make optimisation harder in the global inducing case, especially for larger datasets where the model uncertainty does not matter as much, as the posterior concentration will be stronger. Importantly, however, our results on CIFAR-10 indicate that these issues do not arise in very large-scale, digeesimmnsional datasets, which are of most interest for future work.
used (100), in combination with the fact that global inducing has
We provide tabulated results, including for RMSEs, in Appendix P
M.1 EXPERIMENTAL DETAILS
Here, we matched the experimental setup in Salimbeni & Deisenroth (2017) closely possible. In particular, we used 100 inducing points, and full-covariance observation noise. However, our parameterisation is still somewhat different from theirs, in part because approximate posterior is defined in terms of noisy function- valuue,, while their approximate posterior was defined in terms of the function- values themselves.
23

---

Under review as a conference paper at ICLR 2021
boston
cancrete
energy
kinänm N=0192,
naval
power
protein
wine
yatht
ELBOS
0.50
-..11-0.0
-0.41
0.81
1.39
0.19
-0.16
2.41
3.50
0.02
0.05
1.07
-1.01-120
-1.17
1.16
1.83
5
-253
-2.3--3.15
3..0-0.87
-0.69 1.34
1.38
7.90
2.72
-271-0097
-0.94
local inducing
global inducing
Figure A6: ELBOs per datapoint and average test log likelihoods for DGPs on UCI datasets. The numbers indicate the depths of the models.
As the original results UCI splits, and did not provide the ELBO, we reran their code https: /tt///gttuu///// com/ ICL-ICLLIII SMLISMML//RA HoML/Doubly http//obbb-sschssttiiiaaea of epochs and noise variance to reflect the values in the paper), which gave very similar log likelihoods to those in the paper.
Salimbeni & Deisenroth (2017) used differ-
(changing the number
N MNIST 500
For MNIST, we considered a LeNet-inspired model consisting of two conv2d-relu maceeuloooooooll blocks, followed by conv2d-relu-linearr where the convolutions all have x 3 kernels with 64 channels. We trained all models using a learning rate of 33
When training very small datasets, such the first 500 training examples in MNIST, we can see variety pathologics emerge with standard methods. To help build intuition for these pathologies, introduce sanity check for the ELBO. In particular, distribution over all lower-layer parameters equal to the prior, and sets the top-layer parameters so uniform. With 10 classes, this results in average test log
could imagine model that sets the
ensure that the predictions likelihood of 2230 R log(1/10), and an ELBO (per datapoint) of approximately -2.30. We found that many combinations of the approximate posterior/prior converged to ELBOs near this baseline. Indeed, the only approximate posterior to escape this baseline for ScalePrior and SpatiallW is global inducing points. This because ScalePrior and SpatialIWPrior both offer the flexibility to shrink the prior variance, and hence shrink the weights towards zero, giving uniform predictions, and potentially zero KI. divergence. In contrast, NealPrior and StandardPrior do offer this flexibility: you always have to pay something in KI. divergence in order to give uniform predictions. We believe that this the reason that factorised performs better than expected with NealPrior, despite having an ELBO that is close to the baseline. Furthermore, it is unclear why local inducing gives very test log likelihood and performance, despite having an ELBO that is similar to factorised. For StandardPrior, all the ELBOs are far lower than the baseline, and far lower than for any other priors. Despite this, factorised and fac + gi in combination with StandardPrior appear to transiently perform better in terms of predictive accuracy than any other method. These results should sound note of caution
24

---

Under review as a conference paper at ICLR 2021
NealPrior
StandardPrior
ScalePrior
SpatiallPPrir
0.0
-2.
ELBO
5.0
factoried
local inducing
-7.5
fac
gi
global inclucing
0
test log- sggllikllhoodd
-2
1.0
test prob correct
0.8
0.6
20
40
20
40
20
40
20
40
thousand epochs
thousand epochs
thousand epochs
thousand epochs
Figure A7: The ELBO, test log likelihoods and classification accuracy with different priors and approximate posteriors reduced MNIST dataset consisting of only the first 500 training examples.
whenever try use factorised approximate posteriors with fixed prior covariances (e.g. Blundell et al. 2015 Farquhar et al., 2020). We leave a full investigation of these effects for future work.
O ADDITIONAL EXPERIMENTAL DETAILS
All the methods were implemented in PyTorch. We ran the toy experiments on CPU, with the UCI experiments being run mixture of CPU and GPU. The remaining experiments - linear, CIFAR-10. and MNIST 500 were run on various GPUs. For CIFAR-10, the most intensive our experiments, trained the models on one NVIDIA Tesla 110-PCIE-16GBB We optimised using ADAM (Kingma Ba 2014) throughout.
Factorised We initialise the posterior weight means to follow the Neal scaling (i.e. drawn from NealPrior); however, we use the scaling described in Appendix G to accelerate training. We initialise the weight variances to 1e-3/ /N N for each layer.
Inducing point methods For global inducing, we initialise the inducing inputs, Uo, and pseudo- data for the last layer, V2+1, using the first batch of data, except for the toy experiment, where we initialise using samples from N (0, 1) (since we used more inducing points than datapoints). For the remaining layers, we initialise the pseudo-data by sampling from N (0. 1). We initialise the log precision -4, except for the last layer, where we initialise
We additionally use scaling factor of 3 as described in Appendix G For local inducing. the initialisation is largely the same, except we initialise the pseudo-data for every layer by sampling from N (0, 1). We additionally sample the inducing inputs for every layer from N 10.
Toy experiment For each variational method, we optimise the ELBO over 5000 epochs, using full batches for the gradient descent. We learning rate of 1e-2. We fix the noise variance at its
25

---

Under review as a conference paper at ICLR 2021
true value, help assess the differences between each method more clearly. We use 10 samples from the variational posterior for training, using 100 for testing. For HMC, we use 10000 samples to burn in, and 10000 samples for evaluation, which subsequently thin by factor of 10. We initialise the samples from a standard normal distribution, and use 20 leapfrog steps for each sample. We hand-tune the leapfrog step sizes to be 0.0007 and 0.003 for the burn-in and sampling phases, respectively,
Deep linear network We use inducing points for the inducing point methods. We use sample from the approximate posterior for training and 10 for testing, training for 40 periods of 1000 gradient steps, using full batches for each step, with a learning rate of 1e-2.
UCI experiments
The splits that we used for the UCI datasets can be found at atttts::/gitttuu D/rpaa/ppppottttrrrataatttttttaa
CIFAR-10 The CIFAR-I0 dataset
cs wtt.::s..ootto..eee//ii////e/e
html (Krizhevsky et al. 2009)) is 10-class dataset comprising RGB, 32 x 32 images. It is divided two sets: a training set of 50,000 examples, and validation set of 10,000 examples. For the purposes of this paper, we use the validation set our test set and refer to it as such, as is commonly done in the literature. We use a batch size 500, with one sample from the approximate posterior for training and 10 for testing. For pre-processing. we normalise the data using the training dataset's mean and standard deviation. Finally, we train for 1000 epochs with a learning rate of le (see App. G for explanation of why learning rate is higher than might be expected), and we use a tempering scheme for the first 100 epochs, slowly increasing the influence of the KI. divergence to the prior by multiplying it by factor that increases from 0 to 1. In our scheme, factor in step-wise manner, meaning that for the first ten epochs it is O, then 0.1 for the next ten, 0.2 for the following ten, and so on. Importantly, we still have 900 epochs of training where the standard, untempered ELBO is used, meaning that our results reflect that ELBO. Finally, we note that we share the precisions Ae within layers instead of using separate precision for each output channel as was done in the UCI case. This saves memory and computational cost although possibly at the expense of predictive performance.
increase the
MNIST 500 The MNIST dataset (http://yand) Secaann.eecuue lecun. cww. wtt../////x/mmiiceeeeeeee is dataset of grayscale handwritten digits, each 28 28 pixels, with 10 classes. It comprises 60,000 training images and 10,000 test images. For the MNIST 500 experiments, images from the training dataset and discarded the rest. We normalised the images using the full training dataset's statistics.
trained using the first 500
Deep GPs As mentioned, we largely follow the approach of Salimbeni & Deisenroth (201) for hyperparameters. For global inducing, we initialise the inducing inputs to the first batch of training inputs, and we initialise the pseudo-outputs for the last layer to the respective training outputs. For the remaining layers, we initialise the pseudo-outputs by sampling from a standard normal distribution. We initialise the precision matrix to be diagonal with log precision zero for the output layer, and log precision -4 for the remaining layers. For local inducing, we initialise inducing inputs and pseudo data by sampling from standard normal for every layer, and initialise the precision matrices to be diagonal with log precision zero.
26

---

Under review as a conference paper at ICLR 2021
P
TABLES OF UCI RESULTS
Table 4: Average test log likelihoods nats for BNNs on UCI datasets (errors are standard error)
factorised
local inducing
global inducing
factorised - global
boston (V
-2.74 = 0.03
-2.00 + 0.04
-2.59 = 0.08
-2.60 0.02
NealPrior
-2.76 + 0.04
-2.71 + 0.04
-2.55 0.05
-2.63 0.05
ScalePrior
-3.63 + 0.03
-2.73 + 0.03
-2.50
0.04
-2.61 + 0.04
concrete -areete-Wiie 10
-3.17 + 0.02
-3.28 + 0.01
-3.08
0.01
-3.14 + 0.01
NealPrior
-3.21 0.01
-3.29 0.01
-3.14
0.01
-3.15 0.02
ScalePrior
-3.89 = 0.09
-3.12 0.01
-3.15 + 0.01
energy energy-N(00, 10
-0.76 = 0.02
-1.75 + 0.01
-0.73 + 0.03
-0.68 + 0.03
NealPrior
-0.79 = 0.02
-2.66 + 0.09
-0.72 + 0.02
00.00 - 0.02
ScalePrior
-2.55 = 0.01
-2.42 0.02
-0.33 0.02
-0.74 0.02
kin8nm (O, 1)
1.24 0.01
1.06 1 0.01
1.22 0.01
1.28 = 0.01
NealPrior
1.26 0.01
1.06 + 0.01
1.18 - 0.01
1.29 + 0.01
ScalePrior
1.23 - 0.01
1.10 0.01
1.19 0.01
1.29 + 0.00
naval (0,
7.25 2225:000
6.06 = 0.10
7.29 0.04
7.23 0.02
NealPrior
7.37 0.03
4.28 + 0.37
6.97 0.04
7.45 + 0.02
ScalePrior
6.99 0.03
2.80 + 0.00
6.63 + 0.04
7.50 + 0.02
power-N(0,))
-2.11 + 0.01
2.82 0.01
-2.00 + 0.01
-2.79 + 0.01
NealPrior
-2.81 0.01
-284 0.01
-2.82 0.01
-2.81 0.01
ScalePrior
-2.82 = 0.01
2.84 0.01
-2.82 + 0.01
-2.11 + 0.01
protein (O,,
-2.83 = 0.00
-293 0.00
-2.44 0.00
-2.81
NealPrior
-2.86 = 0.00
-2.2 0.01
-2.87 0.00
-
ScalePrior
-2.86 + 0.00
-2.11 + 0.00
-2885 + 0.00
-2.80 0.00
wine (O, 1)
-0.98 = 0.01
00999
-0.96 0.01
0.06
NealPrior
-0.99 = 0.01
0.98
-0.97 0.01
-0.96
ScalePrior
-1.22 + 0.01
0.96 0.01
0.96 0.01
yacht 10,
-1.41 0.05
22399 + 0.05
0.68 0.03
-1.12 0.02
NealPrior
-1.58 = 0.04
-1.84 0.05
0.81 0.03
-1.04 = 0.01
ScalePrior
4.12 = 0.03
2.71 + 0.22
0.79 = 0.02
-0.97 0.05
27

---

Under review as a conference paper at ICLR 2021
Table 5: Test RMSEs for BNNs on UCI datasets (errors are + standard error)
factorised
local inducing
global inducing
factorised + global
boston Coston-N(C) (O, 1)
3.60 = 0.21
3,85 = 0.26
3.13 + 0.20
3.14 = 0.20
NealPrior
3.64 0.24
3.55 + 0.23
3.14 + 0.18
3.33 + 0.21
ScalePrior
9.03 + 0.26
3.57 + 0.20
2.97 + 0.19
3.27 + 0.20
concrete (0, 1)
5.73 = 0.11
6.34 + 0.11
5.39 + 0.09
5.55 = 0.10
NealPrior
5.96 = 0.11
6.35 + 0.12
5.64 0.10
5.70 0.11
ScalePrior
12.66 + 1.04
6.48 0.12
5.56 + 0.10
5.68 = 0.09
energy NO
0.51 0.01
1.35
0.50 + 0.01
0.47 0.02
NealPrior
0.51 = 0.01
1.95
0.49 - 0.01
0.47 0.01
ScalePrior
3.02 0.05
2.67 + 0.06
0.50 + 0.01
0.49 0.01
kin8nm (0,
0.07 0.00
0.08 0.00
0.00
0.00
NealPrior
0.07 =
0.08 0.00
0.07 0.00
0.07
ScalePrior
0.07 0.00
0.08 + 0.00
0.00
0.07
naval (O, 1)
0.00 0.00
0.00 0.00
0.00 0.00
0.00
NealPrior
0.00 0.00
0.01 = 0.00
0.00
0.00 0.00
ScalePrior
0.00 - 0.00
0.01 = 0.00
0.00
0.00
power power-N (0. 1)
4.00 + 0.03
4.06 + 0.03
3.96 0.04
3.93 : 0.04
NealPrior
4.03 0.04
4.03
4.06 0.03
4.00 0.04
ScalePrior
4.06 = 0.04
4.12
4.05 0.04
4.01 0.04
protein (0, 1)
4.12
4.54 0.02
4.14
4.04
NealPrior ScalePrior
4.21 0.01
4.50
4.27 0.02
4.06
4.22 0.02
4.46
4.19
4.00 0.02
wine (O, 1)
0.65 0.01
0.63 + 0.01
0.63 + 0.01
NealPrior
0.66 0.01
0.01
0.64
0.64
0.01 0.01
ScalePrior
0.82 0.01
0.65 0.01
0.64
yacht-Ji (0, 1) 0.98 0.07
2.35 0.13
0.56 + 0.04
0.50 0.04 0.04
NealPrior
1.37 0.11
0.57 0.04
0.63 0.63
ScalePrior
14.55 + 0.59
5.75
0.56 0.04
0.63
28

---

Under review as a conference paper at ICLR 2021
Table 6: ELBOs per datapoint in nats for BNNs on UCI datasets (errors are + standard error)
factorised
local inducing
global inducing
factorised + global
boston - (0, 1)
-1.55 = 0.00
-1.54 = 0.00
-1.02 + 0,01
-1.03 + 0.00
NealPrior
-1.03 + 0.00
-0.99 0.00
-0.63 + 0.00
-0.70 + 0.00
ScalePrior
-1.54 + 0.00
-0.96 = 0.00
-0.59 + 0.00
-0.70 + 0.00
concrete (0, 1)
-1.10 = 0.00
-1.08 0.00
-0.71 + 0.00
-0.78 0.00
NealPrior
-0.88 = 0.00
-0.87 0.00
-0.99 0.00
-0.65 + 0.00
ScalePrior
-1.45 = 0.01
-0.88 0.00
-0.57 0.00
-0.63 = 0.00
energy (O, 1)
-0.13 = 0.02
-0.53 0.01
0.00
0.59 0.00
NealPrior
0.21 = 0.00
-0.33 = 0.04
0.95
0.79 = 0.01
ScalePrior
-1.12 + 0.00
-0.47 0.01
0.96+0.0
0.80 0.01
kin8nm (0,
-0.38 = 0.00
-0.43 = 0.00
0.26
-0.31 = 0.00
NealPrior
-0.35 = 0.00
-0.43 = 0.00
00.31 0.00
-0.28 0.00
ScalePrior
-0.51 + 0.00
-0.39 = 0.01
-0.29 0.00
-0.27 0.00
naval (0,
1.68 0.01
1.65 + 0.09
2.89
2.11 + 0.02
NealPrior
2.02 0.02
-0.09 0.34
2.53 0.04
2.62 0.02
ScalePrior
1.91 = 0.03
-1.42 = 0.00
2.30 0.03
2.71 = 0.02
nower-N(O,, 1)
-0.08 + 0.00
-0.06 = 0.00
0.02 0.00
-0.03 + 0.00
NealPrior
-0.05 = 0.00
-0.55 0.00
0
-0.01 0.00
ScalePrior
-0.13 = 0.00
-0.05 0.00
0.01
-0.01 0.00
protein (0.
-1.09 = 0.00
-1.14 0.00
1.00
-1.99 = 0.00
NealPrior
-1.11 = 0.00
-1.13 = 0.00
1.00
-1.99 0.00
ScalePrior
-1.13 = 0.00
-1.12 0.00
1.07 + 0.00
-1.7 0.00
wine (0, 1)
-1.48 = 0.00
-1.47 0.00
1.36+0.00
-1.36 = 0.00
NealPrior
-1.31 = 0.00
-1.30 = 0.00
-1.22
-1.23 0.00
ScalePrior
-1.46 = 0.00
-1.29 0.00
-1.22 0.00
-1.23 0.00
yacht - (0, 1)
-1.04 + 0.02
-1.30 + 0.02
0.08 + 0.01
-0.23 + 0.01
NealPrior
-0.46 = 0.02
-0.39 0.01
0.74
0.31 0.01
ScalePrior
-1.61 = 0.00
0.77 0.10
0.79
0.30 - 0.01
29

---

Under review as a conference paper at ICLR 2021
Table 7: Average test log likelihoods for of Salimbeni & Deisenroth (0117). and our implementations of local and global inducing points for deep GPs of various depths.
[dataset) {depth}
DSVI
local inducing
global inducing
boston
2
2.50 = 0.05
-2.48 + 0.05
-2.35 1: 0.05
3
-2550 0.05
-2.49 + 0.05
-2.36 + 0.05
4
-251 0.05
- 0.05
-2.35 = 0.05
5
-2.52 0.05
-2.47
-2.35 = 0.05
concrete 2
-3.14 0.01
-3.13
-3.11 0.02
3
-3.15 0.01
-3.13
-3.11 0.02
4
-3.15 = 0.01
-3.13 0.02
33.1 0.02
5
-3.11 0.01
-3.13 0.02
3.11 0.02
energy
2
-0.72 0.02
-0.71 0.03
0.69 0.03
3
-0.78 = 0.02
-0.71 0.03
0.70 0.03
-0.88 0.02
-0.71 0.03
0.70 0.03
-0.79 = 0.02 5
0.70
00.00 0.03
kin8nm 2
1.34 + 0.00
0.00
1.00
1.38 = 0.01
0.00
1.34 - 0.00
4
1.38 = 0.00
1.00
1.34 0.00
5
1.36 + 0.00
+
1.34 0.00
naval 2
6.55 + 0.22
+ 0.03
7.99 - 0.06
6.46 0.23
7.82 0.04
8.02 + 0.05
4
6.79 + 0.05
7.03
8.00 = 0.03
5
6.70 + 0.07
7.82
7.94 0.04
power 2
-2.78 = 0.01
-2.74
-2.78 = 0.01
3
-2.74 0.01
-2.72 0.01
-2.78 0.01
4
-2.74 0.01
+
-2.78 0.01
5
0.01 -2.74
-2.72 0.00
-2.78 0.01
protein
2
-2.80 0.00
-2.11 0.00
-2.82 0.00
3
-2.73 0.00
-2.78 + 0.02
-2.78 + 0.00
-2.71
0.00
-2.79 0.01
-2.78 0.00
5
-2.
2.70 + 0.00
-2.82
-2.78 = 0.00
wine
2
0.95 = 0.01
0.96
-0.96 = 0.01
3
0.95 0.01
-0.96 0.01
-0.96 + 0.01
4
0.95
0.01
-0.96 + 0.01
-0.96 0.01
5
- 0.95 + 0.01
0.01
-0.96 0.01
yacht
2
-0.41 0.03
0.27 + 0.06
-0.29 + 0.05
-0.46 0.03
-0.33 0.04
-0.25 0.06
4
-0.49 0.03
-0.00
-0.25 0.05
5
-0.49 = 0.03
-0.27 + 0.05
-023 0.04
30

---

Under review as a conference paper at ICLR 2021
Table 8: Test RMSEs for rerun of Salimbeni & Deisenroth (2017). and our implementations of local and global inducing points for deep GPs various depths.
(dataset) (depth)
DSVI
local inducing
global inducing
boston 2
2.94 + 0.17
2.87 + 0.14
2.63 + 0.14
2.94 + 0.18
2.89 + 0.14
2.68 0.14
4
2.98 = 0.18
2.88 + 0.15
2.64 0.13
2.99 0.19
2.87 387
2.63 = 0.13
concrete
5.73 0.09
5.53 0.09
5.42 0.10
5.75 = 0.09
5.14
5.46
5.74 = 0.09
5.13
5.42
5.49 = 0.09
5.52
energy
0,49 = 0.01
0
0.48 0.01
0,51 + 0.01
0.48
0.48 0.01 0.01
0.56 0.02
0.49
0
5
0.51 = 0.01
0.48
kin8nm
2
0.06 = 0.00
0.00
0.00
0.06 = 0.00
0.00
0.06
0.06 = 0.00
0.00
0.06
0.06 + 0.00
0.00
0.06:0.00
naval
2
0.00 + 0.00
0.00
0.00 0.00
0.00 0.00
0.00
0.00
4
0.00 = 0.00
0.00
0.00
0.00 + 0.00
0.00
0.00
power
3.87 = 0.04
3.74 0.04
3.87 0.04
3.75 = 0.04
3.67 0.04
3,89 0,04
3.75 + 0.04
3.68 + 0.04
3.87 0.03
3.73 + 0.04
3.67 + 0.05
3.87 + 0.04
protein
3.99 0.02
4
4.08 0.01
3.77 + 0.01
0.05
3.95 0.01
3.70 0.01
3.99
3.93
5
3.69 + 0.01
4.08
3
wine
2
0.63 = 0.01
0.63 + 0.01
0.63 + 0.01
0.63 + 0.01
0.63 +: 0.01
0.63 + 0.01
0.63 0.01
0.63 + 0.01
0.63 + 0.01
5
0.63 = 0.01
0.63 0.01
0.63 - 0.01
yacht
2
0.42 + 0.03
0.34 + 0.03
0.36 + 0.03
0.44 0.03
0.33 1 0.03
0.35 0.03
0,44 - 0.04
0.33 0.03
0.35 4.03
0.45 = 0.03
0.35 + 0.03
0.33 003+0.00
31

---

Under review as a conference paper at ICLR 2021
Table 9: ELBOs per datapoint for of Salimbeni & Deisenroth (2017) and our implementa- tions of local and global inducing points for deep GPs of various depths.
[dataset) [dephh}
DSVI
local inducing
global inducing
boston
2
-0,48 = 0.01
-0.37 + 0.01
-0.29 + 0.01
3
-0,49 0.01
-0.38 + 0.01
-0.28 1 0.01
4
-0.50 = 0.01
-0.37 = 0.01
-0.30 = 0.01
5
-0.52 = 0.01
-0.37
-0.29 + 0.01
concrete 2
-0.52 0.00
-0.42
-0.39 0.00
3
-0,53 0.00
-0.42
-0.40 = 0.01
4
-0.58 = 0.01
-0.42 0.00
-0.99 0.00
5
-0.62 0.00
-0.43 0.01
-0.88 = 0.00
energy
2
0.97 0.00
1.40
3
0,86 + 0.01
1.02
1.48 0.00
0.73 0.02
0.00
1.48 0.00
5
0.82 - 0.02
1.44
1.48 0.00
kin8nm 2
-0.18 0.00
0.06
-0.16 0.00
-0.18 = 0.00
0..6
0.16 0.00
4
-0.18 0.00
-0.16
-0.16 0.00
5
-0.19 = 0.00
0.16
0.16 0.00
naval 2
2.02 + 0.27
3.63 + 0.04
1.87 0.31
3.29
3.57 + 0.06
4
2.30 + 0.06
3,25
3.54 0.06
5
2.17 + 0.08
3.32
3.51 + 0.05
power 2
0.02 + 0.00
0.05 0.00
0.03 = 0.00
3
0.03 = 0.00
0.05 - 0.00
0.03 = 0,00
4
0.03 + 0.00
0.05 0.00
0.04 4:4 0.00
5
0.03 + 0.00
0.05 0.00
0.03+0.00
protein
2
1.06 0.00
-1.07 0.00
1.06 0.00
3
1.02 + 0.00
-1.05 0.01
-1.05 - 0.00
11.01
0.00
-1.06 0.01
-1.04 0.00
5
- -1.01 + 0.00
-1.08 + 0.00
-1.04 0.00
wine
2
-1.19 0.00
-1.18 + 0.00
-1.17 + 0.00
3
-1.19 = 0.00
-1.18 + 0.00
-1.17 + 0.00
4
-1,20 0.00
-1.18 + 0.00
-1.17 + 0.00
5
-1.20 = 0.00
-1.88 - 0.00
-1.17 + 0.00
yacht 2
1.16 + 0.01
1.78 + 0.01
1.86 + 0.01
1.09 0.01
1.80 +
1.91 0.01
4
1.07 + 0.01
1.81 0.01
5
1.07 1 0.01
1.81 0.01
0.01
32


Reviews:
[{"review_id": "xEi6lDpHqnB", "reviewer": "AnonReviewer6", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper proposes a posterior approximation for BNN that models correlations between the layers weights. The paper begins by pointing out that, for any posterior distribution approximation, the optimal conditional posterior distribution over the top-layer weights given the weights of the previous layers has a closed-form in the case of Gaussian likelihood, which, due to Bayes rule, turns out to be product of the likelihood and the top-layer prior. Based on this insight the paper proposes to model each conditional posterior distribution over intermediate layers weights given previous layers weights following the same structure, that is, as a product of a 'pseudo-likelihood' over unobserved noisy activations and the prior for that layer. In order to make inference tractable, the paper proposes the use of global inducing points as well as noisy pseudo-observations of the activations of intermediate layers which are treated as variational parameters.  The paper also describes how such procedure applies to convolutional neural networks and how it can be applied for DGPs as well.\n\nEven though the paper covers several topics, its presentation is clear and straightforward. The main issue I see with the paper is that the validation of the proposed model is done mostly in terms of accuracy and log-likelihood. As pointed out by [1], one of the main interests on posterior approximations that allow correlations between the weights is for capturing the uncertainty of the compositional structure. Based on this, I think the paper could increase its strength by describing the advantages on interpretability of the learned approximate posterior provides over a fully independent posterior approximation such as [2]. In particular, it would be useful to see how the conditional posterior distributions evolve given the previous layers weights.\n\n[1] Ivan Ustyuzhaninov, Ieva Kazlauskaite, Markus Kaiser, Erik Bodin, Neill DF Campbell, and Carl Hen-\nrik Ek. Compositional uncertainty in deep Gaussian processes. UAI, 2020.\n\n[2] Hugh Salimbeni and Marc Deisenroth. Doubly stochastic variational inference for deep Gaussian\nprocesses. In Advances in Neural Information Processing Systems, pp. 4588–4599, 2017.\n", "extra_scores": {}, "extra_text": {}}, {"review_id": "wQsHS2wQOb2", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The paper proposed a new way of doing Bayesian deep learning in which the optimal conditional posterior for the last layer weights could be reached if the inducing input $Z_0$ is chosen to be the input data $X$ and the pseudo-observation for the last layer $V^L$ is the observation $Y$. Instead of factorizing the inducing points The global inducing input $Z_0$ is propagated through the network to ensure posterior dependencies across layers. The authors also extend this idea to deep Gaussian processes so that the latent functions across layers are correlated. Experiments show better performance than previous methods on both synthetic and real datasets, without the need to anneal the weights for the KL term in ELBO.\n\nPros:\nThe idea of capturing the posterior dependencies between layers in BNN and DGP is great, and may have a good amount of significance in both areas.\n\nCons:\nThe description of the method itself is not very clear. For the BNN part, while Alg. 1 is easy to follow, how does it relates to the training process of BNN? The DGP part is even harder to understand -- what is the complete form of the variational distribution and the ELBO? How is the ELBO optimized? The paper didn't say much about these important details, and by just staring at Eq. (16) I really couldn't figure them out. Comparing to methods the paper compares with such as (Salimbeni and Deisenroth, 2017), Section 3 of this paper is not clearly written.\n\nOther questions:\n\nHow does $N_{l-1}$ in Eq. (2) disappear in Eq. (14)?\n\nDoes this method apply to GP classification?\n\nIn Sec. 4.1, 100 inducing points per layer was used on a dataset with only 40 points. What if the number of inducing points is smaller? Does the proposed method still work as expected?\n\n\n-----------------------\n\nI thank the authors for their detailed reply. The revision after the initial review solves the concerns and questions I had. Very good work. I vote for accept.", "extra_scores": {}, "extra_text": {}}, {"review_id": "E2k5OlmR8ok", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "**Summary**: \nThe paper proposes a posterior estimation for Bayesian neural networks (BNNs) and deep Gaussian processes (DGPs). The difference from the previous approaches is to use global inducing points which help to take into account correlations across layers. \n\nI find a paper well written and the results worth publishing. \n\n**Comments**:\n“We found that the ELBO for methods that factorise across layers — factorised and local inducing — drops rapidly as networks get deeper and wider (Fig. 2). This is undesirable behaviour, as we know that wide, deep networks are necessary for good performance on difficult machine learning tasks. In contrast, we found that methods with global inducing points at the last layer decay much more slowly with depth, and perform better as networks get wider.” I wonder if it can be due to the increasing tails with depth in BNNs with Gaussian prior as stated by Vladimirova et al. (2019). According to their paper, the distributions on the units across the layers become heavier-tailed but knowing that the limit is Gaussian Process (Matthews et al., 2018; Lee et al., 2018), the distributions on the units for wider neural networks are closer to Gaussians (I think the center of the distribution reminds the Gaussian more and more). And as the methods that factorise across the layers use the assumptions of Gaussian units, the error is higher for shallow NNs. So maybe it is possible to keep in mind the “heaviness” for shallow NNs induces due to the correlations between the units and to improve the results. \n\nLee et al. (2018) Deep neural networks as Gaussian processes, ICLR\n\nMatthews et al. (2018) Gaussian Process Behaviour in Wide Deep Neural Networks, ICLR\n\nVladimirova et al. (2019) Understanding priors in Bayesian neural networks at the unit level, ICML\n\n\n\n**Minor**:\n- “initalization”\n- I think you should stick either to American, either to British version of writing: “factorises”, “parameterise”, “initialization”, “initialize”\n- “work using using infinite-width neural networks”\n- “The next block show”\n- “it is difficult design” -> to\n- “downweighted” \n- “by coupling inducing inputs to inducing outputs” -> to induce \n- “the number of locations within an patch”\n- “we need to be able sample”\n- “these issues do not arises”", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: The paper proposes a variational inference method for Bayesian neural networks where the approximate posterior models the correlations between the weights at all layers, using the concept of “global” inducing points. 

Some concerns raised by the reviewers regarding how global inducing points allow us to capture uncertainty across the compositional structure and clarity of the manuscript have been addressed by the authors. However, the more general issue of interpretability has been left for future work. 

One of the main deficiencies of this paper is that it seems to ignore other scalable approaches that also provide more complex posteriors, for example, those based on stochastic gradient Hamiltonian Monte Carlo (see, e.g., https://arxiv.org/abs/1806.05490, and references therein). Overall, although there was support for this paper, it is unclear if approaches such as those presented here are really necessary. A comparison between the two methodologies maybe not only illustrative but required. 

Decision: Reject