Scalable Stochastic Gradient Riemannian Langevin Dynamics in Non-Diagonal Metrics
Authors: Hanlin Yu
Marcelo Hartmann
Bernardo Williams
Arto Klami
Venue: Accepted by TMLR ()

Abstract: Stochastic-gradient sampling methods are often used to perform Bayesian inference on neural networks. It has been observed that the methods in which notions of differential geometry are included tend to have better performances, with the Riemannian metric improving posterior exploration by accounting for the local curvature. However, the existing methods often resort to simple diagonal metrics to remain computationally efficient. This loses some of the gains. We propose two non-diagonal metrics that can be used in stochastic-gradient samplers to improve convergence and exploration but have only a minor computational overhead over diagonal metrics. We show that for fully connected neural networks (NNs) with sparsity-inducing priors and convolutional NNs with correlated priors, using these metrics can provide improvements. For some other choices the posterior is sufficiently easy also for the simpler metrics.

Published in Transactions on Machine Learning Research (08/2023)
Scalable Stochastic Gradient Riemannian Langevin Dynamics in Non-Diagonal Metrics
Hanlin Yu
hanlin.gnu@heiee
Deportment of Computer Science, University of Helsinki, Finland
Mareelo Hartmann
mareed..aartmamms@ccttttink
Department of Computer Science, University of Helsinki, Finland
Bernardo Williams
bernardo. hamaadd.wllllammmmmoreeeoccceei.
Department of Computer Science, University of Helsinki, Finland
Arto Klami
arto arto.lamm@helsinkkff
Department of Computer Science, University of Helsinki, Finland
Reviewed on OpenReview: https: // openreview. net/ forum? id- dXAuvo60GI
Abstract
Stochastic-gradiett sampling methods are often used to perform Bayesian inference on neural nelworks. has been observed that the methods which notions of differential geometry are included tend to have better performances, with the Riemannian metric improving posterior exploration by accounting for the local curvature. However, the existing methods often resort to simple diagonal metrics to remain computationally efficient. This loses some of the gains. We propose two non-diagonal metrics that can be used in stochastic gradiiettgrdeettt samplers to improve convergence and exploration but have only a minor computational overhead over diagonal metrics. We show that for fully connected neural networks (NNs) with sparsity-inducing priors and convolutional NNs with correlated priors, using these metrics can provide improvements. For some other choices the posterior is sufficiently easy also for the simpler metrics.
1 Introduction
Bayesian methods are increasingly being used for large-scale models, especially deep neural networks, with significant literature both sampling-base based algorithms (Welling & Tch, 2011; Chen et al., 2014; Zhang et al., 2020; Vono et al., 2022) and distributional approximation algorithms (Graves, 2011; Bhundell et al., 2015; Huix et al., 2022; Nazaret & Blei, 2022). Even though sampling-based methods are often considered inefficient in the general Bayesian inference literature, they are highly practical for neural networks due to the stochastic-gradieet sampling methods following closely the process of stochastic optimization (Wenzel et al., 2020). In effect, we get efficient samplers by relatively minor modification of optimization algorithms. There are two prominent familics of stochastic-gradion sampling methods, based on Langevin dynamics (Welling & Teh, 2011) and Hamiltonian dynamics (Chen al., 2014). We build on the former that has the advantage of ease of implementation and fewer dyper-parameters but our contributions focusing accounting for local curvature could likely be used also in Hamiltonian dynamics.
The samplers based on Langevin dynamics operate by numerical integration of a particular stochastic difler- ential equation (SDE) (Ma et al., 2015). Our starting point the Stochastic Gradient Riemannian Langevin Dynamics (SGRLD) method where the underlying SDE is defined on manifold characterised by a metric tensor G()) D (for a problem of parameters) and the integration is carried along the manifold. This formulation allows improved sampling by careful choice of the metric (Ma et al., 2015). Perhaps the optimal

---

Published in Transactions on Machine Learning Research (08/2023)
choice is to use metrics that accurately capture the local curvature of the target distribution and e.g. Giro- lami & Calderhead (2011) have demonstrated that using the Fisher information matrix as the metric tensor provides advantages in exploration for sufficiently small problems where it is technically feasible. However, this comes with a substantial computational challenge: Already constructing the Fisher information matrix is computationally costly and the samplers need to repeatedly compute its inverse and inverse square root.
For large-scale problems the practical methods sacrifice flexibility for computational efficiency by using simplified metrics. The early works used constant G that was further restricted to be diagonal, resulting in the basic SGLD method (Welling & Teh, 2011). More recent methods use it to be diagonal to retain computational efficiency. In a suitably chosen diagonal metric the computational overhead is negligible as element-wise operations can be used throughout the algorithm. Furthermore, the metric is typically not formed using differential geometric arguments, and is often adaptively tuned during optimization. For example, Li et al. (2016) used G(0) motivated by the RMSprop preconditioner used commonly for optimization, and Wenzel et al. (2020) constructed a diagonal G(Ø) with separate scaling for each layer of neural network. These methods remain the standard choices for inference for general cases. Some attempts of using non-diagonal metrics have been done, but all of the current solutions are structure-dependentt and only applicable for specific scenarios. Marccau-Caron & Ollivier (2017) and Nado et al. (2018) provide computationally more efficient approximations for the Fisher information metric in form of quasi-diagonal approximations and the K-FAC approximation, respectively, but they still incur a higher computational cost and require explicit construction of the approximation for a specific network structure. Recently, Lange et al. (2023) constructed a non-diagonal preconditioner that accounts for the curvature induced by the batch normalization operation and showed that it improves efficiency for networks using that operation, but they do not provide general recipe for forming the metric for all other structurcs.
G(0) but assume
In this work we provide two new SGRLD variants using non-diagonal metrics G(() that - when implemented in suitable manner are computationally efficient and applicable for general deep neural networks with arbitrary network structurcs. The resulting algorithms yield well-defined operations for different network layer lypes without specific treatments, in contrast to the previous allempts. The first method builds on the Monge metric Hartmann et al. (2022) introduced in the context of Lagrangian Monte Carlo (Lan et al., 2015). This metric is combination of a diagonal matrix and rank-one matrix that depends only the gradients, taking the form ID + V'VEVET with & being the log posterior, having differential geometric justification. Despite of full-matrix form, the metric is fast to compute and has efficient inverse. We incorporate the metric into SGRLD and show that we can implement many operations in element- wise manner just as in the case of diagonal metrics. The sampler is generally applicable as the metric docs not rely on c.g. the structure of a neural network, but is found to be quite sensitive to a Luning parameter a² controlling the metric.
The second metric builds on the Shampoo method Gupta et al. (2018) proposed for optimization. They construct full matrix preconditioner G(Ø) by using Kronecker products of layer-wise matrices that can be efficiently updated during optimization. We incorporate their preconditioner The updates are considerably more efficient than direct inversion of the full matrix but the sampler still has some computational overhead compared to diagonal metrics, in our experiments around factor of 1.3 - 2.4. We argue that such small overhead is justified by consistent good performance in our experiments.
metric tensor in SGRLD,
We evaluate the metrics in problems of varying complexity. We observe that for fully connected neural networks with Gaussian priors the metric may not be that important. For heavy-tailed priors often recom- mended for these networks (Fortuin et al., 2022), such as the horseshoe prior (Carvalho et al., 2009) used by Ghosh et al. (2019) and Popkes et al. (2019), the posterior is more challenging and the choice of the metric is more important. For convolutional neural networks, we also observe differences between the metrics.
2 Background
We denote by 0 E DR. the parameters of a model, here a neural network, and want to infer the posterior
p(0 X) = exp(- (xp(-U(()/T)
(1)
Z
2

---

Published in Transactions on Machine Learning Research (08/2023)
expressed in terms of log-potential U(0) - logg((X 0) to facilitate usage of gradient-based algorithms. Here Z is the intractable normalization constant and 7 is a temperature parameter that is sometimes used to improve predictive performance (Wenzel et al., 2020). We write the equations in this general form but use т = corresponding to standard Bayesian inference in all experiments. We use the notation & for a stochastic estimate for quantity 2, estimated from subsample of observations. In particular, Ve((()) is the estimate of the full gradient VQU(() and we refer to the estimate as the stochastic gradient.
Stochastic Gradient Riemannian Langevin Dynamics (SGRLD) (Girolami Calderhead, 2011; Patterson Teh, 2013; Ma al., 2015) builds on the Stochastic Differential Equation (SDE) (Särkkä & Solin, 2019)
d0 = - 1 VoU(() (2) kaW rT))dtt
where
D a (G(0) r2(0) = Σ k-1 00
Here W is Brownian motion, G(0) E R DxD metric tensor characterising the manifold we are operating on, and the remaining parameters correspond those of Equation 1. The SDE has differential geometric justification, as it can be seen as performing sampling on the manifold (Girolami & Calderhead, 2011). With G(()) ID. we recover the original SGLD proposed by Welling & Teh (2011). Under some mild conditions, it can be proven that we obtain approximate samples from the posterior distribution by following the trajectories induced by the SDE (Ma et al., 2015). In principle any integrator could be used, but in practice the Euler-Maruyama integrator is used throughout the SGLD literature. The SDE is discretized with step-size h₁ > 0 and the samples are obtained by iterating the update rule using stochastic gradients
00 - ((,)), (2Th,G(0,)) RR + =T(0,)
where R2 ~ ID). The total integration time becomes S Σ h2 where T is the number of iterations.
This integrator has two core computational challenges. The second and third terms require computing the metric's inverse and square root of the inverse, which in general has O(D*) complexity. The other fundamental challenge is that in the last term r(0) involves derivatives of the inverse and is challenging even if the inverse itself could be obtained efficiently. However, that term can luckily be omitted in practical algorithms. Girolami & Calderhoad (2011) showed that the term naturally disappears with the assumption that the metric tensor does not depend on the model parameters but is constant G())) = M. Moreover, Li et al. (2016) provides the following theorem that indicates we get bounded estimation error for quantities of interest even if completely ignoring the term, and consequently considering r(() explicitly long we verify have similar bound for any given metrics consider.
can derive practical samplers by not
Theorem 2.1. Denote operator norm as ||- With assumption on smoothness and boundedness as provided in the Appendir A.1, after ignoring the r(0) terms during discrelized updates, for any quantity & estimated as the empirical expectation over the posterior samples, the mean-square error of the estimate for a general SGRLD sampler is bounded as
12 ST
ET
h2)2
T h -I
1)
E ( (6- C Σ
2
+
+
+
Σ
ST
S
ST
for some constant C > 0 independent of {he}, where AV, is an operator that is defined in the Appendix A. where we also replicate the proof in our nolation.
2.1 Previous SGRLD methods
Most previous SGRLD methods build on Euler-Maruyama integration as in Equation 3, only differing in the choice of the metric. Next we briefly explain the key methods commonly used in practice from the perspective of the metric, always presenting both the metric and the inverses required for performing the integration to highlight their computational properties. Importantly, all of these use diagonal metrics. To simplify notations, we will use ğ to denote (o((()) divided by the number of data points in the training set.
3

---

Published in Transactions on Machine Learning Research (08/2023)
SGLD Welling & Teh (2011) used constant metric
G(O.) = C(O.) - -( In.
(4)
This is highly efficient, but the metric completely naive.
PSGLD Li et al. (2016) constructed G(0) inspired by the RMSprop preconditioner. We formulate it
G( = diag (u(0))
G(O,) diag (1 Ø v(0)).
G(O.) à diag 2 Ve(0)),
where (O)) is initialized as a zero vector and updated with exponential moving average (EMA) with suitable as
V(O.) - AV (1 2 (0.) = v (V(0)) +c.
The notation diag(-) denotes the operation that generates diagonal or block-diagonal matrix from an arbitrary number of scalars or matrices, and e is a small constant added mainly for numerical reason. Here the square ?, square root v- and division Ø are all performed element-wise and hence the algorithm has linear cost; we consistently use the notation where vi refers to element-wise operation and the mairix square root is denoted by G1/2.
Wenzel's SGLD Wenzel et al. (2020) constructed layer-wise metric. The implementation by Fortuin et al. (2021) follows
G(O.)) diag (10(0)))
G(O.) 1 diag (ml))(()))))(mnitt)))
G(O.) diag (((1/V((0))) In )
L
where the number of layers and D; the number of paramoters for the ith layer, and
a (0.): V mean(V((0.))),, o (0,); or/min((((0..))
where the update rule of V( is similar to above, but initialized to be vector full of ones. c has similar role to above and the operations are element-wise. The preconditioner is updated periodically, in practice after every epoch.
3 SGRLD in Non-diagonal metrics
For more efficient exploration of the posterior, for instance when the parameters are strongly correlated or the posterior involves arcas of strong curvature, we should be using full matrix G(Ø), rather than diagonal one that can only scale the dimensions. Even though the computation for arbitrary metrics is prohibitively expensive, we construct specific metrics of low cost that still capture the correlations.
Next we introduce two such metrics that lead to efficient SGRLD samplers without assuming the metric to be diagonal. The first. based on a metric derived from the perspective of differential geometry that has efficient inverses due to being a combination a diagonal matrix and a rank-one term, proposed by (Hartmann et al., 2022) for sampling algorithms with exact gradients. The latter uses efficient Kronecker products construct metric that sufficiently efficient compute for large neural networks, originally proposed in the imiiaization context (Gupta et. al., 2018; Anil et al., 2021) and generalized here for sampling.

---

Published in Transactions on Machine Learning Research (08/2023)
20
20
-10
-10
-20
-20
-5 0
-3 0
0,
Illustration of how the Monge metric captures the local curvature of the target density, here the funnel distribution (see Section 4.2). The two plots illustrate the local metric around two separate parameter settings, in form of geodesic paths (red arrows) for different initial velocities sampled from a Euclidean ball and the surfaces of the final positions (solid red line). The metric helps in reaching the narrow funnel in y-direction (left) and similar Fuclidean metric the flat areas (right).
Figure 1:
3.1 SGRLD in Monge metric
Hartmann et al. (2022) recently proposed metric named Monge metric. They embed the target probability distribution in higher-dimensiona Euclidean space M {(O), where O S, (0, (U)()), yielding first fundamental form on M, which can be seen as a natural metric induced by the graph of the function. It thus has a direct differential geometric justification. Figure illustrates the geodesics induced by the metric for the funnel distribution, demonstrating how it characterises the local curvature. See Section 4.2 for details of the example and empirical validation of how SGRLD in the Monge metric, as described next, also improves sampling from funnel.
The metric is given by G(f) = T 7 for lighter notation) G(0))1 = 1- T+a IV VUVUT and equally easy square root of that. Hartmann et al. (2022) used it in the context of Lagrangian Monte Carlo (Lan et al., 2015) and showed that the metric helps to achieve good exploration of complex target distributions, but their sampler does not scale for large problems due to requiring exact gradients and delerminants for rejection checks.
and hence has the inverse (dropping dependence
Here we use that metric in the context of SGRLD, which requires some additional development. Even though the expressions above are efficient compared to direct matrix inversion, the cost of actually forming the square metric from the direct application of -horisnnMorrison inversion would still be challenge for very large models. Next show that by careful manipulation of the equations we avoid forming the matrix explicitly and can implement SGRLD with element-wisc operations, eventually matching the computational cost of diagonal approximations. We also show that standard exponential moving average of stochastic gradients is applicable for forming the metric tensor.
To derive the efficient element-wise operations, we start by defining helper functions to simplify notations:
a
f-1(x) = 1 a² |
/ +2+ 22 ( 1)
f-4(x)= ||0||
We can then express the metric G(O), its inverse and inverse square root as
G(O.) ID + VVV4V1,
G(O.) / /-
(5)
G(O,) = To + + (1
5

---

Published in Transactions on Machine Learning Research (08/2023)
where Vlo is initialized to zero and updated based the stochastic gradient using standard EMA
v4. AV4-1 +1 A)ge-
Observe that all computations required by the integrator are of the form "( for some vector and either n -1 -1/2. Using this notation, the Ith entry of the required product can be computed using
[G(
[z fo (V4) VV))VL [ar] fa (VI)) V4),,
where (.--) is the inner product that combines the last element of Equation with Since f(-) only depends the current moving average of gradients, can implement the exact updates for batch by iterating through all the parameters once to calculate the needed quantities, followed by individual updates for cach of the parameters. In practice, all operations are still computed with standard lensor operalions.
Both the memory and computation complexity is linear in D, which is the same as using diagonal precon- ditioners like RMSprop as in Li et al. (2016). However, it is important to note that the metric does not correspond a diagonal approximation. As shown by Equation 5, the effect of the metric is addilive instead of multiplicative. For all diagonal samplers the G(0)" terms in Equation 3 take the forms where elements of the gradient multiplicd by a scaling factor, whereas here we retain the gradient itself and add scaled version of moving average to it.
For the specific case of probabilistic models, i.e. ( = - log Ty (3) 0), the Fisher information matrix characterizes the local curvature. As shown by Hartmann et al. (2022), Monge metric relates to the Fisher information matrix in expectation (with respect to Y) such that its expected value is identity matrix plus a² times the Fisher information matrix. However, this connection is somewhat thcoretical as practical use we evaluate for the observed data y rather than computing the expectation. It also important to note that despite similar expressions the outer product used for constructing the Monge metric does not correspond to the empirical estimate of Fisher information discussed by Kunstner et al. (2019). The Monge metric uses the outer product of full-data gradients (or their stochastic estimates), whereas the empirical Fisher estimate is an expectation of the outer products of per-sample gradients, as detailed Appendix A.8. The expectation of the outer product has different computational properties than the outer product of expectations, and the obscrvations Kunstner et al. (2019) make do not directly apply. However, the Monge metric retains the problematic scaling of empirical Fisher due to use of 'squared" gradients, and the scaling parameter a² is needed in part to account. for this.
We will later observe that the performance depends strongly on the choice of a² and hence the metric introduces a new hyper parameter that needs to be selected carefully. When a² = 0 the metric reduces to identity and differs increasingly more from the Euclidean one for larger values. Note that a value here does not have the same interpretation as in Lagrangian Monte Carlo, as there an implicit scaling depending on the number of data points in the training set. The validity of the sampler for all a² 2 shown by the following theorem, with the proof provided in the Appendix A.2.
Theorem 3.1. For SGRLD with the Monge metric, after ignoring the I(0) terms during discretized updates, we can bound the approrimation error as defined in Theorem 2.1 as
(3--0) ")-))= - SC
Σ
h2 5 ST
E 2
+ ST
[Σ T -=
h2)2) + +1(a97(1-))))
S
We note that the bound is provided solely for ensuring asymptotic validity of the sampler and it does not directly say anything about the empirical accuracy of the sampler for different a. There no reason believe the bound to be particularly tight and the large exponent in a is not necessarily problematic already a small can dramatically change the metric.

---

Published in Transactions on Machine Learning Research (08/2023)
3.2 SGRLD in Shampoo metric
An alternative to using a metric with efficient inverses is to use metric for which the inverses can be efficiently updated iteratively during the sampling process. Shampoo (Gupta et al., 2018; Anil et al., 2021) shows how to do this in the context of optimization, deriving a preconditioner possessing a close relationship with full-matrix Adagrad (Duchi et al., 2011). They additionally provide details for a version related to full-matrix RMSprop. We build on that version.
We denote the tensor rank of the Ith parameter as di, and the individual shapes of the Ith parameters (mr)1..... (ni) J. We can then expressed the metric
G(O,) = diag( क. H(m))
G(O.) -1 diag( [M के H(t)), G(O,)- diag( ( [m क5 H(t),
where we use i to index the dimension of the parameter, [Hi]] Elp, and
[[ A(H (1 1)
Here (1) := -ja-19j. where are all possible indexes; see Appendix A.9 for details the notation. Each [Hi] is therefore square matrix, with both the width and the height given by the ith dimension the shape of the Ith paramcter. We used a slight abuse of notation here, by using [ğe] denote the gradients with respect to the corresponding parameter but reshaping them here to match the shape of the parameter that is now tensor.
Next, we will show how this metric can be used in the context of SGRLD. In the original Shampoo, cach HH)] shi is calculated numerically after adding small constant to avoid numerical issucs. For SGRLD
also need access to [H 44, and hence we first compute [H(], 4d; mumerically, and then compute [Hj]] 24 by taking the square. The Kronecker product formulation allows efficient updates with standard tensor operations (Anil et al., 2021).
The computational cost is then C(E ( and the memory cost is C(EL Emh (m))., for some positive constants C. These costs are strictly larger than those of the diagonal metrics or the Monge metric, but in practice the metric can be updated relatively efficiently (e.g. periodically, after hundreds of steps), and hence the computation remains manageable. To further reduce computational complexity, we can divide the tensors into smaller blocks and treat them as individual tensors instead (Anil et al., 2021).
The Shampoo metric also relates to the empirical Fisher information discussed by Kunstner et al. (2019). As shown in Appendix A.8, it has a connection with the square root of the empirical Fisher matrix and it hence avoids the inverse gradient scaling problem, but lacks direct interpretation in terms of curvature. The effect of taking a square root of the metric is further discussed by Staib et al. (2020).
The following thoorem shows the asymptotic validity of the sampler, with the proof in the Appendix A.3. The bound depends on the EMA parameter À in the same way as the bound for the Monge metric and otherwise follows that of Theorem 2.1.
Theorem 3.2. For SGRLD with the Shampoo metric, after ignoring the r(0) terms during discretized updates, we can bound the approzimation error as defined in Theorem 2.1 as
(6--)) 0
h EIAK112 52
1 ET k2)2 ST
+0(2(1-1)3)
S
4 Experiments
These experiments characterise the behaviors of the proposed samplers and additionally provide information about the nature of the posteriors of some network architectures.
7

---

Published in Transactions on Machine Learning Research (08/2023)
4.1 Experimental setup
Comparison methods We evaluate the two proposed metrics against other commonly-used SGRLD algorithms applicable for arbitrary network structures explained in Section 2. We refer to these methods by their metrics, using Identity for Welling & Teh (2011), RMSprop for the pSGLD method of Li et al. (2016) and Wenzel for Wenzel et al. (2020).
Experimental details The code that can be used to reproduce all experiments can be found in https: //gi thub. wtt.///oo/////cee /SG//LLD Concerning neural network experiments, our implementations for all methods are built on top of Fortuin et al. (2021). For RMSprop, Wenzel and Shampoo we use À = 0.99 and e = 1 matching the choices of Fortuin et al. (2022), whereas for Monge we use = 0.9 based on good performance in preliminary experiments. For all methods we select constant learning rate based on performance (log-probability) on a separate validation sct, and for Monge we additionally select a² that controls the metric within the same process. For all methods we use fixed schedule for learning rate to isolate the effect of the metric. For instance, Zhang et al. (2020) showed that cyclical rates can improve sampling, but since the interaction between such schedules and the underlying metric is unknown, we do not employ them here. We use no data augmentation to isolate its effect on likelihood tempering (Kapoor et al., 2022).
We run the samplers for a total of 400 epochs using a batch size of 100. The first 1000 steps are treated as burn-in, and the actual samples used for evaluation are collected after that with a thinning interval of 100 steps. Each setup is repeated for 10 independent runs. Additional details, like the learning rates for each case, are provided in the Appendix A.4 and A.5.
Evaluation measures Except for the funnel experiment where we compare against a known ground truth, we evaluate the samplers from the perspective of predictive accuracy and quantification of predictive uncertainty, following the general observations of Wilson & Izmailov (2020). In particular, we evaluate the log-probability logp(y)XX X)) and classification accuracy computed as the ensemble average over the posterior samples. We compute both for test data not used during inference or validation, averaging the measures over ten independent runs, and also report the running times. hog-probability is a proper scoring rule (Gneiting & Raftery, 2007) that requires correctly reprosenting the predictive uncertainty and is generally used as the primary metric for evaluating quality of neural networks providing predictions with uncertainty estimates (Lasshminarayanan et al., 2017; Tran et al., 2022) and for comparing alternative inference methods (Zhang et al., 2020). However, we note that it does not directly measure how well the samplers explore the true posterior.
We also use a measure that quantifies the curvature of the posterior. Srinivas et al. (2022) proposed a notion, denoted by C(O), that measures the model's curvature with respect to the data distribution by quantifying how much the local space deviates from Euclidean. It can be efficiently computed for arbitrary neural networks and they proposed it for training models of lower curvature, but we use their measure for estimating the average curvature of the regions the samplers visit during inference as the expected value of C(O) over the posterior samples.
We use this measure for two purposes. It primarily used as an approximation of the inference complexity. It not an exact mcasure as the samplers are not guarantcod to perfectly sample from the true posterior, but still gives indication of how complex the posterior is for different architectures. We also use the average curvature of different samplers as additional indication of their differences. Lower higher average curvature does not directly indicate a sampler to be more correct since poor sampler could equally well be unable to reach areas of high curvature or to get stuck there, but differences between methods are still informative.
Finally, we also measured the expected calibration error (ECE) by Naeini et al. (2015) for the predictions. Rahaman & Thiery (2021) explains that ensemble mothods, including sampling- basod Bayesian neural net- works, are not necessarily good in terms of ECE and hence this evaluation provided as a complementary evidence in Appendix A.7, rather than as a part of the main results.
8

---

Published in Transactions on Machine Learning Research (08/2023)
Identity
RMSprop
Monge
Shampoo
Figure 2: Funnel distribution in four metrics. The better metries explore the funnel clearly better, though even the best one has difficulties reaching the very end. The right subplots show the challenging marginal of the 2D distribution (left subplot), indicating the true marginal with blue line and the samples with yellow histograms.
4.2 Funnel
The funnel distribution (O.,.... It, sp(0p)/n-1))((00) |0, where sp(x)) = log(1 e*), is known to have difficult geometry for sampling and has been used e.g. by Hartmann et al. (2022) and Neal (2003) for sampler evaluation. We consider a funnel distribution with D 2. In Figure 1 already used this distribution for visual demonstration of how the Monge metric captures the local curvature of the target distribution. In the following, we show how SGRLD algorithms with identity, RMSprop, Monge and Shampoo metrics behave in this task. When the gradient is a vector of one dimension, as is here, ignoring technicalities like periodic updates and the role of e, Shampoo metric coincides with full matrix RMSprop.
We draw 2e samples from the funnel, using well chosen hyperparameters for each method. The gradients that here are known analytically available for the samplers are corrupted by noise of standard deviation of one, to emulate a scenario where stochastic gradients are used. Figure 2 shows that the identity metric is unable to reach the area of high curvature. The better metrics clearly help, but even the Monge metric that is here the best may not be sampling from the exact correct distribution. Hartmann et al. (2022) showed that LMC can solve the problem nearly perfectly in this metric in combination with exact gradients, but SGRLD is not a perfect algorithm for the problem with stochastic gradients.
4.3 Fully-connected networks
We use fully comnected neural networks of size 784-N-N-10 on MNIST dataset (LeCun et al., 2010), where we use N E [400, 800, 1200], with sctup inspired by Li et al. (2016) and Korattikara et al. (2015). We evaluate all the samplers in problems of varying difficulty by (a) changing the prior distribution of the weights and biases and (b) changing the network size by modifying N. Tables (log-probability and accuracy) and 2 (average curvature and time) report the results for the commonly used Gaussian prior ~ (O, of o²In) and for the horseshoe prior 0 ~ AOO,, o² I2D where À follows half Cauchy distribution with mean 0 and variance 1, as implemented in Fortuin et al. (2021). Following Fortuin et al. (2022), the prior scale is adjusted based on network structure.
The gist of the results is that one of the nia-diaaonal metrics is the best terms of log- bog-pobbabiliy in all cases, the horseshoe prior is clearly better, and the non-diagonal metrics help more when the posterior is challenging. Even though the numerical differences are somewhat small, the standard deviations (computed over 10 runs) are even smaller and the differences wetween samplers are reliable. Next analyse the results in more detail from different perspectives.
Effect of prior For the Gaussian prior the posterior is easy, shown by the low average curvature for all methods. For these cases the choice of the metric is not particularly critical, but the proposed Shampoo metric is still generally the best. For the Monge metric when hidden unit size is 400 or 800 the optimal a² = 0 and hence the metric reduces to G(O) = Tp, but with the largest network it reaches an overall best performance with a² 0.75.
For the Horseshoe prior the posterior is considerably more challenging, with average curvature around 3 times higher, and see that there is clear benefit. in using this prior the ang-probabilities are clearly
9

---

Published in Transactions on Machine Learning Research (08/2023)
Gaussian prior
Horseshoe prior
Metric
logp(y|
Acc.
log logp(y| X)
Acc.
Identity
F0.1310, 0.0004
[0.9685, 0.0005]
-0.0750, 0.0004
[0.9823, 0.0003
Wenzel
00.1315, 0.0002]
[0.9686 0.0003]
-0.0768, 0.0009)
[0.9807, 0.0003)
400
RMSprop
F0.1307, 0.0005
[0.9689, 0.0004]
-0.0678, 0.0011]
[0.9798, 0.0008]
Monge
Identity (a
0) is optimal
-0.0629, 0.0009]
[0.9831, 0.0005]
Shampoo
0..128 0.0003]
[0.9689, 0.0004]
-0.0641, 0.0009
0.9820, 0.0006]
Identity
-0.1667, 0.0004]
0.9587, 0.0003
-0.0798, 0.0003
0.9818, 0.0004
Wenzel
F..1667, 0.0004
[0.9583, 0.0005]
0.0839, 0.0008
0.9789,, 0.0003
800
RMSprop
F0.1643, 0.0003]
[0.9614, 0.0003]
-0.0656, 0.0013
[0.9800, 0.0005]
Monge
Identity (a²
= 0) is optimal
-0.0612, 0.0014]
[..9844, 0.0008]
Shampoo
(001620,, 0.0004]
[0.9605, 0.0003)
-0.0603, 0.0011]
0.9884,, 0.0004]
Identity
0..193,, 0.0004
0..9544, 0.0004
-0.0809, 0.0004
0.9814, 0.0002
Wenzel
0.1193,, 0.0004
[0.9519, 0.0003]
-0.1006, 0.0009
0.9776, 0.0005
1200 RMSprop
0..1886, 0.0005]
0.95688 0.0003]
-0.0631, 0.0011
[0.9808, 0.0007]
Monge
00.1780, 0.0005]
0.95600, 0.0005]
-0.0694, 0.0005
0.9833, 0.0004]
Shampoo
-0.1854, 0.0003]
0.95559, 0.0004
-0.0564, 0.0007]
0.9834, 0.0005]
Table 1: Inference accuracy for MNIST with fully connected networks of varying size N and prior, marking the best metric for each configuration. Each entry is given as mean followed by standard deviation.
Gaussian prior
Horseshoe prior
N Metric
Curv.
Time
Curv.
Time
Identity
5.93
8.1
14.65
8.3
Wenzel
5.82
7.7
14.79
9.5
400
RMSprop
6.30
8.8
25.76
10.7
Monge
a² 0 is optimal
15.08
12.9
Shampoo
6.28
18.8
17.84
20.1
Identity
4.69
19.9
15.01
24.8
Wenzel
4.68
20.3
14.45
26.2
800
RMSprop
4.67
24.2
20.09
29.3
Monge
a² = 0 is optimal
12.09
32.1
Shampoo
4.53
46.5
12.74
53.3
Identity
4.16
43.5
15.28
56.3
Wenzel
4.10
49.1
11.41
58.2
1200
RMSprop
4.11
56.9
16.96
58.3
Monge
6.92
59.2
16.78
60.1
Shampoo
3.93
92.9
10.07
89.1
Table 2: Curvature and computational efficiency for MNTST with fully connected networks of varying size N and prior. The time is given seconds per epoch.
better than using the Gaussian prior. For this more challenging case, both of the proposed metrics provide substantial improvement over the baselines, especially in terms of log- pogobbbbbb,lity, and the optimal Monge metric always has a² > 0, with 1.25, 0.5 and 0.075 being the optimal values for the three sizes. We also see clear differences between the samplers in terms of the average curvature C(O), revealing that the pSGLD sampler using the RMSprop metric behaves in different manner than the others. For = 1200 with some learning rate the Wenzel metric had high variance over the runs, with some bad runs explaining the poor final performance.
Effect of network size In terms of the network size the results between the samplers are relatively consistent, but an important observation is that with Gaussian priors the smallest networks are the best whereas with the horseshoe the largest ones outperform the smaller ones but only if using the better
10

---

Published in Transactions on Machine Learning Research (08/2023)
0.984
00.07
0.982
log-prebability
accuracy
-0.08
0.978
-0.09
0.976
-0.10
Isenthy
Waned
RMSppop
-0.11
0.972
0.12
0.970
100
200
300
100
200
300
400
epochs
epochs
Log-pro pog-probyyility (top) and accuracy (bottom) of different samplers on MNIST with hidden layer size 400 and horseshoe prior. Shaded regions show #.96 standard deviations computed over replicates.
Figure 3:
0.984
-0.05
log proppabbbiityy
accuracy
0.982
-0.070
0.980
idetrry
-0.075
Weszsl
RMSpru
Monge
Shampon
0.978
0.0
0.5
1.0
1.5
2.0
0.0
11.5
1.0
1.5
2.0
a²
a²
Figure Log-pue bog-prbbability (left) and accuracy (right) in Monge metric with varying a² on MNIST with hidden layer size 400 and horseshoe prior, in comparison to other metrics. Shaded areas show 1.96 standard deviations computed over 10 replicates.
samplers. A practitioner using the identity metric might conclude that the best they do for this data is horseshoe prior with N 400, but in the better metrics can effectively also larger networks.
Convergence Figure 3 illustrates convergence of the samplers (averaging predictions over samples col- lected thus far) for the case of N 400 and the horseshoe prior, showing that all metrics result in similar convergence. This plot is shown as a function the iteration. The Shampoo metric takes around 1.3 2.4 times longer per iteration, but all other metrics share roughly the same cost.
The RMSprop metric that was earlier seen to have clearly different C(O) compared to others is here revealed to overfit. it quickly reaches good log- log-probability but the accuracy does not increase monotonically confirming that unusual average curvature can be sign of poor sampling.
Monge metric parameter Figure 4 shows how the Monge metric depends on its tuning parameter a for 400 with the horseshoe prior. For a 0 it matches identity but for a² > O we see clear improvement over the diagonal metric. For a > 1.25 the performance starts to degrade possibly due to numerical issucs, and the curve suggests may be possible further improve the performance if these issues could be better avoided. As noted earlier, in some other cases a² is the best choice and the metric becomes identity.
4.4 ResNet architecture
As another common network architecture we consider the ResNet (He al., 2016), known have relatively easy optimization landscape (Li et al., 2018). We consider alternative prior choices to inspect the potential effect the samplers: (a) independent Gaussian priors 0 oN) and (b) correlated normal
11

---

Published in Transactions on Machine Learning Research (08/2023)
Independent Gaussian prior
Metric
logp(y
Acc.
Curv.
Time
Identity
--0.4627, 0.0039
00.591, 0.0015
10.18
7.6
Wenzel
H0.4869, 0.0047]
[0.8533, 0.0023]
9.12
7.5
RMSprop
-0.4751, 0.0023]
(088566, 0.0014]
10.27
8.8
Shampoo
-0.4606, 0.0043]
(085596, 0.0018)
10.63
11.4
Correlated normal prior logp(ulX)
Metric
Acc.
Curv.
Time
Identity
-0.4437, 0.0040)
00.6611, 0.0025]
9.91
11.6
Wenzel
--04611, 0.0050]
[0.8614, 0.0021
8.91
11.8
RMSprop
[-0.4500, 0.0041]
(0.8642, 0.0018]
9.94
14.1
Shampoo
[-0.4365, 0.0034]
[0.8668, 0.0023]
10.39
15.4
Table 3: CIFAR10 with ResNet architecture. The best Monge metric always matches the identity metric and is not shown separately. Boldface indicates the best metric. For log-probability and accuracy, each entry is given as mean followed by standard deviation. The time is given as seconds per epoch.
prior proposed by Fortuin et al. (2022). We use the CTFAR10 data (Krizhevsky & Hinton, 2009) and Google ResNet-20 as implemented by Fortuin al. (2021).
The experimental results are shown in Table 3. Similar to the previous experiments, the Shampoo metric is the best in both cases and the computational overhead over the alternatives is manageable. For the Monge metric the optimal choice is here always to resort to identity metric with a² = 0 and hence no improvement is observed. In terms of posterior difficulty, we observe that switching to correlated priors makes the curvature slightly lower, matching the hypothesis of Fortuin et al. (2022), while also improving the log- log-probability and accuracy terms of all motries. Shampoo results in the best performances both cases, while Wenzel and RMSprop are worse than using identity metric. While the standard deviations are larger than the MNIST case, Shampoo still yields consistent improvements over identity especially when using correlated Gaussian prior, with bog-probability roughly equal to that of identity metric plus two times the standard deviation.
5 Discussion
We used relatively small data and models for evaluation in interest of manageable overall computation. SGLD-based methods have previously been applied on larger problems as well (Zhang et al., 2020), and we showed that the computational overhead of the new metrics small, with the Monge metric being essentially as fast as the diagonal ones. We see no reason why the methods could not be applied also for larger networks, and we could also take into use various other tricks for maximal scaliability, e.g. only doing posterior inference 2010; Ober & Rasmussen, 2019; Watson et al., 2021) using proper metrics only within individual layers (Gupta et al., 2018; Anil et al., 2021). Our metries could most likely also be combined with various sampler advances, like the adaptive drift by Kim et al. (2020).
over the last few layers (Lazaro-Gredilla & Figuciras-Vide
The Shampoo metric performed well in all of our experiments and is currently our practical suggestion. The Monge metric is equally fast as the diagonal ones and some cases worked extremely well, but requires careful choice of a² and for the easier posteriors a² was actually the best. This in part caused by mumerical stability issues; the method appears in principle to work well with larger values but often has issues with gradient explosion and shows degraded accuracy already before failing completely. This is seen also in Figure 4 where the performance suddenly drops after a² = 1.25. The likely core reason for these problems is the inverse scaling with the 'squared' gradient, which has previously been pointed out as problematic for methods using empirical estimates of Fisher information matrix as preconditioners (Kunstner et al., 2019). It is likely that better normalization or adaptive control of possibly clement-wise o² will resolve these issucs, but we do not have good solutions ready yet. The current experiments validate the metric has potential but further development needed for providing robust practical method.
12

---

Published in Transactions on Machine Learning Research (08/2023)
We note that in the main experiments purposefully ignored the samplers' different running times to ensure that they have equal storage cost and computational cost during evaluation, without requiring method- specific thinning intervals or other tricks that would make the interpretations of the results harder. As shown by Tables 2 and 3, the Shampoo metric takes, depending on the case, 1.3 2.4 times longer to compute than the Identity metric. In Appendix A.6 we provide additional empirical results where all samplers are restricted use the same total computation time, to provide an alternative perspoctive for sampler efficiency in terms of wall-clock time. We observe that despite some differences, the Monge and Shampoo metrics still retain the advantages over the previous methods.
Our focus was more in demonstrating that it is possible to have computationally efficient metrics that avoid making the diagonal assumption to encourage further research into this direction, instead of arguing that these specific metrics would necessarily be optimal in any specific sense. We used two metrics derived from different perspectives and Lange et al. (2023) recently provided an additional example of constructing a computationally efficient non-diagonal metric specific the batch normalization operation (Ioffe & Szegedy, 2015), remonstrating that the design space of possible choices is broad. It is likely that other metrics, possibly derived from still different perspectives, can also be made computationally more efficient, but the practical design process remains largely manual.
6 Conclusion
Sampling should be done using a metric that somehow accounts for the curvature of the postcrior distribution. Current stochestli--radienn samplers usually do this with diagonal metrics to retain computational elficiency and use adaptive metrics motivated by the optimization literature. We showed that possible to use non- diagonal metrics while retaining the same or almost the same computational cost, providing two practical metrics that can improve sampling. One metric was consistently good but has some compulational overhead. The other worked extremely well in some cases but is sensitive to a tuning parameter, and further research is needed for robust practical implementations.
Despite theoretical advantages, Riemannian samplers are not commonly used today. For instance, none of the leading probabilistic programming environments supports them. In deep learning approximate second- order optimization methods are already widely used and the interest typically in predictive performance that is easier to evaluate compared to verifying reliability of parameter estimates that can be difficult for Riemannian samplers. This suggests - somewhat counter-intuitively we might see Riemannian met entering practice faster in problems of higher dimensionality, now that scalable approaches are available.
Finally, we showed that the need for improved sampling algorithms depends on the model architecture and in particular the prior distribution. Much of the literature has considered networks with simple priors, in particular the Gaussian prior, and for those priors already simple metrics are sufficient. For more advancod priors we can see larger benefits with better metrics.
Acknowledgments
This work was supported by the Research Council of Finland, by projects 324852, 345811 and 348952 as well as the Flagship programme: Finnish Center for Artificial Intelligence FCAL.
References
Rohan Anil, Vineet Gupta, Tomer Koren, Kevin Regan, and Yoram Singer. Scalable Second Order Opti- mization for Deep Learning, 2021.
Charles Blundell, Julien Cornebise, Koray Kavukcuoglu, and Daan Wierstra. Weight Uncertainty in Neu- ral Networks. In Proceedings of the 32nd International Conference on Machine Learning - Volume 37, ICML'15, pp. 1613 1622. JMLR.rrg, 2015.
Carlos M. Carvalho, Nicholas G. Polson, and James G. Scott. Handling Sparsity via the Horseshoc. In Proceedings of the Twelth International Conference on Artificial Intelligence and Statistics, volume of Proceedings of Machine Learning Research, pp. 73 80. PMI.R, 2009.
13

---

Published in Transactions on Machine Learning Research (08/2023)
Changyou Chen, Nan Ding, and Lawrence Carin. On the Convergence of Stochastic Gradient MCMC Algorithms with High-Order Integrators. In Advances in Neural Information Processing Systems, pp. 2278 2286, 2015.
Tianqi Chen, Emily B. Fox, and Carlos Guestrin. Stochastic Gradient Hamiltonian Monte Carlo. In Pro- ceedings of the 31st International Conference on Machine Learning, Proceedings of Machine Learning Research, pp.
2014.
Joshua V. Dillon, Ian Langmore, Dustin Tran, Eugene Brevdo, Srinivas Vasudevan, Dave Moore, Brian Patton, Alex Alemi, Matt Hoffman, and Rif A. Saurous. Tensorflow distributions, 2017.
John Duchi, Flad Hazan, and Yoram Singer. Adaptive Subgradient Methods for Online Learning and Stochastic Optimization. Journal of Machine Learning Research, 12(61):2121 2159, 2011.
Vincent Fortuin, Adrià Garriga-Alonso, Mark van der Wilk, and Laurence Aitchison. BNNpriors: library for Bayesian neural network inference with different prior distributions. Software Impocts, 9:100079, 2021.
Vincent Fortuin, Adrià Garriga-Alonso, Sebastian W. Ober, Florian Wenzel, Gunnar Rätsch, Richard E. Turner, Mark van der Wilk, and Laurence Aitchison. Bayesian Neural Network Priors Revisited. In International Conference on Learning Representations, 1CLR 2022, 2022.
Soumya Ghosh, Jiayu Yao, and Finale Doshi-Velen Model Selection in Bayesian Neural Networks via Horseshoe Priors. Journal of Machine Leurning Research, 20(182):1 46, 2019.
Mark Girolami and Ben Calderhead. Riemann manifold Langevin and Hamiltonian Monte Carlo methods. Journal of the Royal Statistical Society: Series B (Statistical Methodology). 73(2):123 224-114, 2011.
Tilmann Gneiting and Adrian E. Raftery. Strictly Proper Scoring Rules, Prodiction, and Estimation. Journal of the American Statistical Association, 102(477):359 378, 2007.
Alex Graves. Practical variational inference for neural networks. In Advances in Neural Information Pro-
cessing Systems, volume 24. Curran Associates, Inc., 2011.
Vineet Gupta, Tomer Koren, and Yoram Singer. Shampoo: Preconditioned Stochastic Tensor Optimization. In International Conference on Machine Learning (ICML), volume 80 of Proceedings of Machine Learning Research, pp. 1842-1850. PMLR, 2018.
Marcelo Hartmann, Mark Girolami, and Arto Klami. Lagrangian manifold Monte Carlo on Monge patches. In Proceedings of The 25th International Conference on Artificial Intelligence and Statistics, volume 151 of Proceedings of Machine Learning Reseorch, pp. 4764-4781. PMLR, 2022.
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Deep Residual Learning for Image Recognition. In 2016 IEEE Conference on Computer Vision and Pattern Recognition, CVPR 2016, Las Vegas, NV, USA, June 27-30, 2016, pp. 770-778. TEEE Computer Society, 2016.
Tom Huix, Szymon Majewski, Alain Durmus, Eric Moulines, and Anna Korba. Variational Inference of overparameterized Bayesian Neural Networks: theoretical and empirical study, 2022.
Sergey loffe and Christian Szegedy. Batch Normalization: Accelerating Deep Network Training by Reducing Internal Covariate Shift. In Proceedings of the 32nd International Conference on Machine Learning, volume 37 of Proceedings of Machine Learning Research, pp. 448 456, Lille, France, 07 Jul 2015. PMLR.
Sanyam Kapoor, Wesley J. Maddox, Pavel Ixmailov, and Andrew Gordon Wilson. On Uncertainty, Tem- pering, and Data Augmentation in Bayesian Classification. In Advances in Neural Information Processing Systems, 2022.
Sehwan Kim, Qifan Song, and Faming Liang. Stochastic Gradient Langevin Dynamics algorithms with Adaptive Drifts. Xiu preprint ar arXio:2009.095.0 2020.
14

---

Published in Transactions on Machine Learning Research (08/2023)
Anoop Korattikara, Vivek Rathod, Kevin Murphy, and Max Welling. Bayesian Dark Knowledge. In Advances in Neural Information Processing Systems, pp. 3438-3446. MIT Press, 2015.
Alex Krizhevsky and Geoffrey Hinton. Learning multiple layers of features from tiny images. Technical Report 0, University of Toronto, Toronto, Ontario, 2009.
Frederik Kunstner, Philipp Hennig, and Lukas Balles. Limitations of the empirical Fisher approximation for natural gradient descent. In Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
Balaji Lakshminarayanan, Alexander Pritzel, and Charles Blundell. Simple and Scalable Predictive Un- certainty Estimation using Deep Ensembles. Advances in Neural Information Processing Systems, volume 30. Curran Associates, Inc., 2017.
Shiwei Lan, Vasilcios Stathopoulos, Babak Shahbaba, and Mark Girolami. Markov Chain Monte Carlo From Lagrangian Dynamics. Journal of Compulational and Gruphical Stalistics, 24(2):357 378, 2015. ISSN 10618600, 15372715.
Susanna Lange, Wei Deng, Qiang Ye, and Guang Lin. Batch Normalization Proconditioning for Stochastic Gradient Langevin Dynamics. Journal of Machine Learning, 2(1):65 82, 2023.
Vidaal. Miguel Greiila and Aníbal R. Figueiras- Vidal. Marginalized Neural Network Mixtures for Large- Scale Regression. IEEF Transactions Neural Networks, 2185):334-1135) 13551, 2010.
Yann LeCun, Corinna Cortes, and CJ Burges. MNIST handwritten digit database. ATT Labs [Online].
Available: htlp://wann http//yannneeeare www.nneeecnc.o..o com/erdll sonddb/mnis,, 2, 2010.
Chunyuan Li, Changyou Chen, David Carlson, and Lawrence Carin. Preconditioned Stochastic Gradient Langevin Dynamics for Deep Neural Networks. In AAAI Conference Artificial Intelligence, AAAI16, pp. 1788 1794, 2016.
Hao Li, Zheng Xu, Gavin Taylor, Christoph Studer, and Tom Goldstein. Visualizing the Loss Landscape of Neural Nets. In Aduances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018.
Yi-An Ma, Tiançi Chen, and Emily B. Fox. A Complete Recipe for Stochastic Gradient MCMC. In Advances in Neural Information Processing Systems, pp. 2917-2925. MIT Press, 2015.
Gaétan Marceau-Cart and Yann Ollivier. Natural Langevin Dynamics for Neural Networks. In Geometric Science of Information - Third International Conference, GSI 2017, Proceedings, volume 10589 of Lecture Notes in Computer Science, pp. 451-459. Springer, 2017.
Zachary Nado, Jasper Snoek, Roger Grosse, David Duvenaud, Bowen Xu, and James Martens. Stochastic Gradient Langevin Dynamics that exploit Neural Network Structure. In ICLR Workshop, 2018
Mahdi Pakdaman Nacini, Gregory Cooper, and Milos Hauskrecht. Obtaining Well Calibrated Probabilities Using Bayesian Binning. Proceedings of the AAAT Conference on Artificial Intelligence, 29(1), February 2015.
Achille Nazaret and David Blei. Variational Inference for Infinitoly Deep Neural Networks. In Proceedings of the 39th International Conference on Machine Learning, volume 162 of Proceedings of Machine Learning Research, pp. 16447-16461. PMLR, 17-23 Jul 2022.
Radford M Neal. Slice sampling. The annals of statistics, 31(3):705 767, 2003.
Sebastian W. Ober and Carl E. Rasmussen. Benchmarking the Neural Linear Model for Regression. Symposium on Advances in Approcimate Bayesion Inference, 2019.
15

---

Published in Transactions on Machine Learning Research (08/2023)
Yaniv Ovadia, Emily Fertig, Jie Ren, Zachary Nado, D. Sculley, Sebastian Nowozin, Joshua Dillon, Balaji Takshminarayanan, and Jasper Snoek. Can you trust your model' uncertainty? Evaluating predictive uncertainty under dataset shift. In Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
Sam Patterson and Yee Whye Teh. Stochastic Gradient Riemannian Langevin Dynamics on the Probability Simplex. In Advances in Neural Information Processing Systems, volume 26, 2013.
K. B. Petersen and M. S. Pedersen. The Matrix Cookbook, nov 2012. URL. http://www2. http:://ww22.oompute dtu. ck/pubdb/pubs pttbbbbbpuba//24-ffflll 3274-ful1
html. Version 20121115.
Anna-Lena Popkes, Hiske Overweg, Ari Ercole, Yingzhen Li, Jos'e Miguel Hernández=Lobalo, Yordan Za- ykov, and Cheng Zhang. Interpretable Outcome Prediction with Sparse Bayesian Neural Networks in Intensive Care. arXiv preprint arXiv:1905. 190.0:9905,0999, 2019.
Rahul Rahaman and Alexandre Thiery. Uncertainty Quantification and Deep Ensembles. Tn Advances in Newral Information Processing Systems, volume 34, pp. 20063 20075. Curran Associates, Inc., 2021.
Suraj Srinivas, Kyle Matoba, Himabindu Lakkaraju, and Francois Fleuret. Efficient Training of Low- Curvature Neural Networks. In Advances in Neural Information Processing Systems, 2022.
Matthew Staib, Sashank J. Reddi, Satyen Kale, Sanjiv Kumar, and Suvrit Sra. Escaping Saddle Points with Adaptive Gradient Methods, 2020.
Simo Särkkä and Arno Solin. Applied Stochastic Differential Equations. Institute of Mathematical Statistics Textbooks. Cambridge University Press, 2019.
Ba-Hien Tran, Simone Rossi, Dimitrios Milios, and Maurizio Filippone. All You Need is Good Functional Prior for Bayesian Deep Learning. Journal of Machine Learning Research, 23(74):1 2022.
Maxime Vono, Vincent Plassier, Alain Durmus, Aymeric Dieuleveut, and Eric Moulines. Qlsd: Quantised Langevin Stochastic Dynamics for Bayesian Federated Learning. In Proceedings of The 25th International Conference on Artificial Intelligence and Statistics, volume 151, pp- 6459 6500. PMTR, 28-30 Mar 2022.
Joe Watson, Jihao Andreas Lin, Pascal Klink, Joni Pajarinen, and Jan Peters. Latent Derivative Bayesian Last Layer Networks. In International Conference on Artificial Intelligence and Statistics (AISTATS), volume 130 of Proceedings of Machine Learning Research, pp. 1198-1206, 2021.
Max Welling and Yee Whye Teh. Bayesian Learning via Stochastic Gradient Langevin Dynamics. In International Conference on Machine Learning, 2011.
Florian Wenzel, Kevin Roth, Bastiaan S. Veeling, Jakub Swiatkowski, Linh Tran, Stephan Mandt, Jasper Snoek, Tim Salimans, Rodolphe Jenatton, and Sebastian Nowozin. How Good is the Bayes Posterior in Deep Neural Networks Really? In International Conference on Machine Learning (ICML). JMLR.org, 2020.
Andrew G Wilson and Pavel Izmailov. Bayesian Deep Learning and a Probabilistic Perspective of Gen- cralization. In Advances in Neural Information Processing Systems, volume 33, pp. 4697-4708. Curran Associates, Inc., 2020.
Ruqi Zhang, Chunyuan Li, Jianyi Zhang, Changyou Chen, and Andrew Gordon Wilson. Cyclical Stochastic Gradient MCMC for Bayesian Deep Learning. In International Conference on Learning Representations (ICLR), 2020.
A Appendix
This Appendix provides proofs for the theorems presented in the paper in Sections A.1-A.3, additional technical details concerning the experiments in Sections A.4-A.5, additional empirical results Scction A.6- A.7, discussions on the relationships of proposed metrics with empirical Fisher in Section A.8 and clarification the formulation of ShampooSGT.D in Section A.9.
16

---

Published in Transactions on Machine Learning Research (08/2023)
A.1 Convergence analysis
The proof is based on Chen et al. (2015) and Li el al. (2016). We impose the same mild assumption of smoothness and boundedness from their paper.
Assumption A.1. tp and its up to 3rd-order derivatives, th are bounded by a function V, i.e. | S
C₁ for k = 20, Cki Pk: > 0. Furthermore, the expectation of V on {0,} bounded: sup, E (P,)) S co, and V smooth such that supse(0.1) (s 0 +(1 s)Y) s (() + (Y)), V0.Y.p max for some C 0.
Consider test function of interest o0)). The posterior average with respect to the invariant measure p(0) related to the Stochastic Differential Equation (SDE) is defined as Φ fxx (())
Our geometric MCMC samplers give us samples (0.) We approximate o with a: 7 ET_ ((0,).
Consider the SDE for SGRLD with exact gradients:
d0 = - vV v (T) dW + +r)(d)
local generator can be written as, where a b = a b and A B tr(A (B)),
L = 1==G-(o))-Voo.uu rr())) Vo. + +G(()) Vo. Vo.
We define functional that solves the following Poisson Equation, C (v,))) = (00)) - Z.
(6)
According to the assumptions, e exists.
Consider the SDE for Ricmannian SGLD after ignoring Gamma:
2dW d0 +
local generator can be written as
- Vo. + TG(0,) Ve. Vo.
Let Li be the generator of the SDE with exact gradients. Observe that
2, +AV+ + AL.,
(7)
where AV, = (-GG(()) V(,U().)) Vo, and AL, = -Tr(O)) Ve.- We use the Fuler-Maryyama integrator, which is a first order integrator. We have E (v)(O))) = A ( O(h))
= =1 (())))))))()))
Sum over and plug in Equation we have
T
Σ B(v(O.)) = Σ (0 (((+))))+ b. Σ hr + E hhe + h
Divide both side by Se, plug in Equation and re-arrange the terms,
1 T-1
1 (E(((or) ((b)) + ((((0)))-
6 è θ
-
ST
ST
1
1
Σr h2 ST
hAV (---))
h, (e- C
ST
ST
1=1
t=1
17

---

Published in Transactions on Machine Learning Research (08/2023)
Therefore, there is positive constant C > 0 satisfying
(ā-0) 0) SC
(E(((O))) (00)))
(((((o))
(b())))++
S7
S
A2
Σ h St
T
h ST
E h
2
FIAVI1
+
Σ
=())))
+
ST
Aa
A1 is bounded according to assumptions, A2 is bounded by O(yhs) because of the added Gaussian noise. Based on these observations, after simplifications, we obtain the following theorem.
Theorem A.2. Denote operator norm as
|| ||. The MSE of the algorithm can be bounded as
h S
+ 1 EE R2)
T h
E E σ) SC
Σ
BIAVI
+ ST
+
rI(),)
S1
S
for some C >>> 0.
Ag term can be further relaxed as
D
T
h
rr((0:))
))
2
Aa S
.I
ST
As step sizes are non-zero, it follows that
he 1
0
ST
for any t. Since D is bounded by definition, if each rI, term is bounded by relatively small quantity, we can have reasonably good estimate of Ø.
A.2 Derivations for Gamma term in Monge metric
Assume first and second order gradients are bounded. Use the following to simplify notations,
Vi = xV= (1-
where Vi is the moving average after the current update and used to form the metric, V is the moving average before the current update and VI is the current stochastic gradients. Further, denote individual entries of the second order derivative for the stochastic estimate of gradients as Ligh, or in matrix form as 221. Since it is a weighted average of bounded terms, is also bounded.
The derivations follow that
2 00,- = -a
VIVIT
(G(0)
(LD a
00,
+ a² |/] 171,
1V1; + (1 A)
-1 Vigoh,
+ 20
1 a² 2112
(1+a2
D (.(0)= Σ
0
(G(() ")
00,-
1=1
(1- 1, (V9))
1L41T (1- E ((.(41) 20
= --0
1 1/12
(1 (1+ (1+121V77[132
18

---

Published in Transactions on Machine Learning Research (08/2023)
Then the T()) term in matrix form can be written
r(0) (1 1) (-a² V S
+ 2a
a VII2
(1+ a² (++1115113)2)
where
81 = Σ
1))
J=1
D
82 = VI (V
That is, the I(0) term for the Monge metric is tractable and can be expressed in matrix form.
Since first order gradients VI are bounded and VI exponential moving average of VI, follows that VI is also bounded. Since second order gradients (V
also bounded, we can conclude that V/Il3) and S2VI/(1+ a² ((]13) are bounded. Therefore, (0 = ( - A))
for all indexes i.
We have for the Monge metric
hi +
(ET h?)? S4.
E ( ②) sc
Σ BLVVP
(1aa22(1-1)))
S
ST
A.3 Derivations for Gamma term in the Shampoo metric
To simplify notation, drop dependencies on time and write
H XH + (1 ())
G( diag ( [ 24 1
Observe that every entry in every gl0 is combination of the gradients. Since the gradients are bounded according to assumptions, every entry in every g( is bounded. Therefore, every Hi is also bounded, and for any indexes x,
0 (1)) 0 A).
00
Recall the Gamma term is,
D
( = Σ
())
00,
1=1
We have
G(
II
(
2%
II 2))
i=1
where 2% n is the modulo of x on n.
Use (0) denote G(0)2d, then we have
M(0) = diag( (0 H=))
We can that each M()jja is the product of d individual terms. It follows that
G(0)-1 M(0)
19

---

Published in Transactions on Machine Learning Research (08/2023)
and
a
00, M( - -1-
since it is the of product of O(1 A) and some bounded terms.
Then each @
(G(O)=)) can equivalently be expressed as
00g
a (M
((M(0) 22 ((M )im. (M
E
0
Lmm,no ((M(0) )m @M
00p
M()) is a positive definite matrix, therefore M() z and M (0) z are also positive definite and unique.
Recall 0(x-1) = -x-- (Petersen & Pedersen, 2012), so (MM()) lik are bounded.
( Dmm
It also holds that (Petersen & Pedersen, 2012)
a(X)) ) n-1 axg (xT) (xa-1-)))
a(M(0)) Jom has closod-form expression, and is bounded. thercfore, (M)m) To summarize,
(M )5k (0) )mm OM
= 0(1 - ()),
((M(O) 1i
iM)
00:
and it trivially follows that
(,(0) 011-A)
for all indexes j. Thus, we have for the Shampoo metric
E ( ")²= sc Σ BIAVI
h ST
1
ET 1=1
11))
+
(1(7((1119)
ST
S7
A.4 Experimental details for sampling from funnel
We use SGRLD samplers with identity, RMSprop, Monge and Shampoon metrics, collecting a total of 2. samples for each.
Concerning hyperparameters, for identity metric we employ a constant step size of 0.001. For RMSprop metric, we employ constant step size of 0.0025, X 0.995 and e 0.0. For Monge metrie, we employ constant step size 0.003, a² = 0.1 and À = 0.7. For Shampoo metric, employ constant step size 0.003, X = 0.9995, le 6 and update interval of 1.
We observe that when using Monge metric, the sampler might suffer from numerical issues when reaching the bottom the funnel. When using other metrics, the samplers do not seem to reach those locations. In order to resolve the numerical issue, we resort identity metrie up to a scaling when the norm of the preconditioned gradients larger than 1000.0.
For the scatter plots, 1000 random samples from all the samples chosen and plotted.
A.5 Experimental details for neural network experiments
Following the practice of Wenzel et al. (2020), we use a bijection between learning rate and step size h { In, and use & the tuning parameter to be validated. We separate validation set to tune the
20

---

Published in Transactions on Machine Learning Research (08/2023)
Gaussian prior
Horseshoe prior
N Metric
Step size
Step size
Identity
0.05
0.25
Wenzel
0.075
0.75
400
RMSprop
0.00025
0.0005
Monge
Identity (a² = 0) is optimal
0.25
Shampoo
0.0025
0.005
Identity
0.05
0.5
Wenzel
0.05
0.5
800
RMSprop
0.0005
0.0005
Monge
Identity (a² =0) 0) is optimal
0.25
Shampoo
0.0075
0.005
Identity
0.025
0.25
Wenzel
0.05
0.25
1200 RMSprop
0.0005
0.0005
Monge 0.025
0.25
Shampoo
0.0075
0.005
Table 4: Learning rates of MNIST experiments
Independent Gaussian prior
Correlated normal prior
Metric
Step size
Step size
Identity
0.1
0.1
Wenzel
0.5
0.75
RMSprop
0.0005
0.0005
Shampoo
0.01
0.01
Table 5: Learning rates of CIFAR10 experiments
hyperparameters. The learning rates are tuned based on grid in the form of [1x 10a, 2.5 1.5 105, 7.5x 10a] for different integers æ. If in at least one run for given learning rate a sampler completely breaks down such that it cannot finish sampling, we conclude that the learning rate is not applicable. Concerning a² for Monge, we first try a² = 1.0, a² = 0.5 and a² = 0.1. If all these choices result in worse or nearly identical performances as identity and the performance generally becomes better as a² decreases, we conclude that the optimal choice is a² = 0. Otherwise, we employ grid search in terms of potentially better a² values.
For fully connected neural networks, following Ghosh et al. (2019), place horseshoe prior both the weights and biases.
We report the best learning rates of each sampler in each experimental setup in Table 4 and Table 5. For the RMSprop metric, following implementations in e.g. Fortuin et al. (2021) and TensorFlow Probability (Dillon et al., 2017), the gradients of the prior are included when calculating the metric. For the Shampoo metric, employ an update interval of 100 steps.
The rumning times of the algorithms using variant of our code base that yields minimal computational workload for training. For MNIST experiments, the code was run on single Intel Xeon @ Gold 6230 CPU @ 2.10GHz core. For CIFAR10 experiments, the code was run on a single NVIDIA @ Tesla V100-SXMP-32GB GPU with 10 Intel Xeon Gold 6230 CPU 2.10GHz cores. measure the time per step, and use that to estimate time per epoch.
estimated by separate 3 rums of the algorithm for 20 epochs each,
A.6 Experiments for equal computation time
Using non-diagonal metries typically results in small increase in running time, as observed in Tables and To ensure the slower algorithms not given unfair advantage, we here report the results for a setup where
21

---

Published in Transactions on Machine Learning Research (08/2023)
Gaussian prior
Horseshoe prior
N Metric
logp(y)
Acc.
log logp(yl X)
Acc.
Identity
F0.1311, 0.0004
[0.9684, 0.0005]
-0.0750, 0.0004]
[0.9823, 0.0003]
Wenzel
F0.1315, 0.0002]
[0.9686, 0.0003]
-0.0777, 0.0009
[0.9806, 0.0003]
400
RMSprop
H..1310, 0.0005]
[0.9689, 0.0006]
-0.0680, 0.0013]
[0.9792, 0.0007]
Monge
Identity (a
0) is optimal
-0.0650, 0.0006]
[0.9826, 0.0006]
Shampoo
00.1304, 0.0002]
[0.9684, 0.0004]
-0.0691, 0.0009]
[0.8812, 0.0008]
Identity
-..1667, 0.0004]
00.9577, 0.0003]
-0.0798, 0.0003
0.9818, 0.0004
Wenzel
0..1668, 0.0004
[0.9583, 0.0004]
-0.0822, 0.0008
[0.9777, 0.0002]
800
RMSprop
F0.1648, 0.0002]
[0.9612, 0.0003]
-0.0656, 0.0013]
[0.9801, 0.0007]
Monge
Identity (a²
= 0) is optimal
-0.0628, 0.0010]
[.9835, 0.0007]
Shampoo
[0.1635, 0.0003]
[0.9596, 0.0004]
-0.0663, 0.0009)
0.98199 0.0006]
Identity
F..1932, 0.0004]
00.9514, 0.0004
-0.0809, 0.0004
0.9814, 0.0002
Wenzel
F..1935, 0.0004]
[0.9519, 0.0003]
-0.1009, 0.0009
0.9745,, 0.0005
1200 RMSprop
(..1894, 0.0006]
0.95666 0.0003]
10.0630, 0.0011]
[0.9808, 0.0007]
Monge
00.1785, 0.0005]
[0.5560, 0.0005]
-0.0699, 0.0005
0.9832, 0.0005]
Shampoo
-0.1878, 0.0003]
00.9541, 0.0003]
-0.0599, 0.0005
0.9831, 0.0003]
Table 6: MNIST with fully connected networks of varying size N and prior for setup where each algorithm is constrained to use the same total computation time, marking the best metric for each configuration. Each entry is given as mcan followed by standard deviation.
each sampler is constrained to have equal total rumning time, matching the time it took to run the fastest sampler in the main experiments. We note that this may result in worse performance for the samplers which take longer to run than they really are, since they are effectively using fewer samples to form the ensemble, and remark additionally that the exact computation time naturally depends on implementation details and hence not all differences in terms of efficiency are caused by the algorithm details.
In previous experiments, we recorded evaluation metrics after every epoch. Since the time points that the samplers finish rumning may be different, the best learning rate may also change. We therefore further tune the learning rates of the samplers based on performances on the validation set, assuming that the running times of the samplers under different learning rales are always given by the times reconstructed from the estimated times from before. In most cases the best learning rate coincide with the learning rates as reported in Table 4. The exceptions for fully connected neural networks experiments are: Gaussian prior, N = 400, the best learning rate for Shampoo is 0.005; N 800, the best learning rate for Shampoo is 0.0025; N 1200, the best. learning rate for Shampoo is 0.0025. Horseshoe prior, N = 400, the best learning rate for RMSprop is 0.00075. The exceptions for ResNet experiments are: Gaussian prior, the best learning rate for RMSprop is 0.00075 and the best learning rate for Shampoo is 0.025. Correlated normal prior, the best learning rate for RMSprop is 0.00075 and the best learning rate for Shampoo is 0.025. The final results are again reported for the test set.
Table reports the results for the MNIST data, and in general the results are in line with the main results reported in Table 1. That is, we still have either Monge or Shampoo metric with the best log- brobbbility for all cases, for horseshoc prior also the best accuracy is obtained by the Monge metric in all cases, and the Shampoo metric for the largest network with the horseshoe prior remains among the best overall result.
Table reports the results for the CIFARIO data. The Shampoo metric that has the slowest per-iteration computation remains as the metric with the best classification accuracy for both priors, but in terms of log-probability the identity metric overtakes all compared methods. However, it is worth noting that it was close to being the best already in the main experiments.
22

---

Published in Transactions on Machine Learning Research (08/2023)
Independent Gaussian prior
Metric
log Dog((y/X
Acc.
Identity
-0.4634, 0.0039
[0.8589, 0.0013
Wenzel
-0.4869, 0.0047]
[0.8533, 0.0023]
RMSprop
-0.4828, 0.0041
[0.8562, 0.0018
Shampoo
-0.4858, 0.0047]
[0.8595, 0.0023]
Correlated normal prior
Metric
Dgg(yyX)
Acc.
Tdentity
-0.4437, 0.0040)
[0.8641, 0.0025
Wenzel
--.4624, 0.0050]
[0.8612, 0.0020
RMSprop
[-.44594, 0.0030
[0.8631, 0.0026]
Shampoo
[-0.4506, 0.0023]
[0.8697, 0.0015]
Table 7: CTFAR10 with ResNet architecture for setup where each algorithm is constrained to the same total computation time. The best Monge metric always matches the identity metric and is not shown separately. Boldface indicates the best metric. Each entry is given mean followed by standard deviation.
Gaussian prior
Horseshoe prior
N Metric
ECE
ECE
Identity
[0.0430, 0.0004
[0.0281, 0.0004]
Wenzel
(0.0438, 0.0004]
[0.0248, 0.0005
400
RMSprop
[0.0447, 0.0003]
[0.0124, 0.0009]
Monge
Identity (a²
0) is optimal
[0.0180, 0.0009]
Shampoo
0..0211, 0.0003)
00.014,, 0.0008)
Identity
0.0525, 0.0003
[0.0316, 0.0004]
Wenzel
0.05113 0.0003
[0.0273, 0.0004
800
RMSprop
00.5586, 0.0004]
[0.0106, 0.0004]
Monge
Identity (a.
0) is optimal
[0.0178, 0.0009]
Shampoo
0.05449 0.0004)
[0.0136, 0.0006
Identity
0.0570, 0.0004
(0.0318, 0.0005
Wenzel
0.0581, 0.0002]
[0.0320, 0.0005
1200 RMSprop
0.06699 0.0004
[0.0095, 0.0005
Monge
0.0544, 0.0004
00.0252, 0.0005
Shampoo
0..0639, 0.0004
[0.0114, 0.0006
Table 8: ECE for the MNIST experiments
A.7 Expected Calibration Error
This section complements the empirical results of the main paper by presenting the Expected Calibration Error (ECE) (Naeini et al., 2015) measure for the calibration the predictions, defined as
ECE = P()) lo: -eil
The predictions are divided into K bins, with P(i) being the empirical probability that predictions fall into the ith bin. o, is the true fraction of positive instances in bin e, is the mean of predicted probabilities of being positive in bin i. ECE is defined based on the intuition that the level of confidence of model should be the same as the probability of correct predictions. Principally, lower ECE is better.
The results serve as further verification that the samplers work as intended, while revealing small differences between the samplers. However, we note that there are some criticisms regarding use of ECE context of deep learning and ensembles, and we refrain from drawing strong conclusions based on these results. ECE is not proper scoring rule, and can be minimized by model making uninformative predictions (Ovadia
23

---

Published in Transactions on Machine Learning Research (08/2023)
Independent Gaussian prior
Correlated normal prior
Metric
ECE
ECE
Identity
[0.0919, 0.0032
[0.0861, 0.0020)
Wenzel
[0.1006, 0.0036
[0.0962, 0.0038
RMSprop
[0.0962,
0.0021
[0.0925, 0.0021
Shampoo
[0.0902, 0.0017
[0.0856, 0.0024
Table 9: ECE for the CIFAR10 experiments
Metric
Average confidences
Identity
0.9544
Wenzel
0.9568
RMSprop
0.9684
Monge
0.9648
Shampoo
0.9663
Table 10: Average confidences of MNIST with fully connected neural networks, hidden unit size 400, horse- shoe prior
et al., 2019), whereas Rahaman & Thiery (2021) shows that ensembling (as done in sampling- basod deep learning methods) not necessarily improving ECE.
Tables 8 and 9 report the ECE for the MNIST and CIFAR1O experiments, respectively. On MNIST with horseshoc prior that was previously identified as a more challenging sampling task, the RMSprop, Monge and Shampoo samplers result in better ECE compared to Identity and Wenzel. For MNIST with Gaussian prior the differences are generally small. For CIFAR10 experiments Identity and Shampoo have the best ECE.
We additionally ran the samplers for 3 independent runs on MNIST, hidden unit size 400, horseshoe prior, using variant of our code recording the confidences (the largest predicted probabilities of the model) and report the averages of the confidences in Table 10. These nunbers that are approximately 96% for all samplers can be compared with the accuracies that are approximately 98% as indicated in Table 1, confirming the general observation of Rahaman & Thiery (2021); all the samplers are somewhat under-confident
A.8 Discussion on the relationship of the metrics with Empirical Fisher
The Fisher information matrix characterises the local curvature of a probabilistic model, and is formed using expectation of the gradients over the possible outcomes of the model. Kunstnor et al. (2019) discusses methods that replace the expectation by conditioning on the observed data y instead, resulting in the so- called empirical Fisher
Σ Vo low.g(p(nkee I 0) Va logp(ynl T)
(8)
computed as the sum of gradients for individual data points. They argue it may not be a good choice as a curvature matrix of preconditioner, but contains useful information about the gradient noise in stochastic optimization. Here we describe how both Monge and Shampoo metrics relate to the empirical Fisher, but also have certain important differences. Both metrics use gradients of the joint density rather than the likelihood, but even if ignoring the prior the expressions are different.
The Monge metric (ignoring stochastic estimation of gradients for now) is
Io a²
Σ Vo ( [yn, smn)
Vo Vo log p(O Iy, F..)
(9)
where the latter term resembles the empirical Fisher but is different: We are now taking the outer product of the full-data gradients, rather than summing over the outer products. Within MongeSGLD compute
24

---

Published in Transactions on Machine Learning Research (08/2023)
stochastic estimate of this quantity, but again the expression only considers the observed data y instead of taking the expectation over their distribution. The core difference between the two expressions is that the second term in Equation 9 includes cross-terms between individual samples while Equation 8 does not.
The Shampoo metric has a connection with full matrix RMSprop and the full matrix RMSprrp can be interpreted as using an estimate, up to scaling, of
Ve log Wn, In Ve p(0 Tn)
That is, Shampoo is related to the empirical estimate for the Fisher matrix in Equation but takes the square root of that. As discussed by Kunstner et al. (2019), this removes possible scaling issues but makes the connection to Fisher information more vague.
A.9 Clarification the formulation of ShampooSGLD
For the convenience of notations, we using superscript i when introducing the formulation of Sham- pooSGL.D In expressions like [Ai],, is not the tensor component, but just an extra index referring to each dimension of the gradient, indexing the matrices, For each i, [H is a square matrix, with both height and width being the shape of the gradients in that dimonsion. These square matrices are constructed based on the moving average of [0=]]²) which is obtained from the gradients. We are using a slight abuse of notation here, because when describing the previous methods we used Ô. to denote the stochastic gradient concate- nated and reshaped to a D x 1 column vector, where D is the lotal nunber of parameters; but here we are using [ōr] to denote the gradients of the /th parameter the same shape of the original parameter.
25


Reviews:
[{"review_id": "ubgtMePJav", "reviewer": "Reviewer_3sW2", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "### Strengths and Weaknesses\n\nThis is an interesting submission which attempts to improve stochastic-gradient-based\nposterior samplers using non-diagonal metrics. \nNon-diagonal metrics can better capture the curvature of the log-probability,\nleading to improved exploration and higher sample efficiency in downstream \ntasks.\nAs such, I think this paper is of interest to the machine learning\ncommunity and a good fit for TMLR. \n\nThe major strength of this submission is the relative simplicity of the \nnon-diagonal metrics, which are based on well-known approximate curvature \nmatrices from the optimization literature, and \nThis means fast implementations are already available and can be leveraged\nimmediately to develop better samplers.\n\nHowever, the submission has several significant weaknesses. From a theoretical\nperspective, the approximate curvature matrices seem to worsen the sample complexity\nof estimation compared to using the identity matrix (Theorems 3.1 and 3.2).\nMoreover, the novelty/relevance of these theorems are clear to me from \nthe submission, although I am not an expert in this area.\n\nEmpirically, the Shampoo metric does lead to improved log-likelihood, but\nessentially the same accuracy as using the identity matrix on both MNIST\nand CIFAR-10.\nThus, it's not clear these improvements are meaningful in practice.\nMoreover, no calibration experiments are provided, so its impossible to \njudge the quality of the uncertainty provided by the samples.\nLog-likelihoods are not a suitable metric for judging calibration, since they \nencourage posterior collapse (the MAP estimate maximizes predictive log-likelihood).\nThus, I find the experimental evaluation to be weak at best.\n\nFinally, the utility of the Monge metric is not well justified by the experiments\nor theory.\nOn CIFAR-10 and MNIST with a Gaussian prior, the curvature matrix\ndefaulted to identity ($\\alpha^2 = 0$), so it's not clear why Monge is useful\nin general.\nThis matches Theorem 3.1, where the complexity suggests $\\alpha^2 = 0$ \nis optimal. \nFinally, is well-known that the outer-product of gradient is a poor estimate\nof the Fisher information [1], so the justification here seems very weak.\nPerhaps this is why $\\alpha^2 = 0$ is often optimal.\n\nTo summarize:\n\n**Strengths**:\n\n- The non-diagonal metrics are fast to compute, improve the log-likelihood on most tasks,\n    and well-known in optimization. \n- The paper is well-written and polished. \n    \n\n**Weaknesses**:\n\n- The theory seems relatively weak and doesn't justify using non-diagonal metrics. \n- No calibration experiments are provided for MNIST/CIFAR-10, so it's unclear\n    non-diagonal metrics lead to improved uncertainty estimation.\n- The Monge metric often reverts to the identity matrix, introduces an additional \n    tuning parameter, and isn't well-justified. \n\n### Correctness\n\nAs mentioned, this paper is outside of my area of research and makes use\nof tools I am not familiar with, so I cannot comment on correctness of the \ntheoretical results.\nI did not check the derivations in the appendix for correctness.\n\n### Writing and Presentation:\n\nThe paper is well-written and I noticed essentially no typographical issues. \nI congratulate the authors on producing such a well-polished text. \n\n\n### Additional References:\n\n\n[1] Kunstner, Frederik, Philipp Hennig, and Lukas Balles. \"Limitations of the empirical Fisher approximation for natural gradient descent.\" Advances in neural information processing systems 32 (2019).\n\n[2] Yao, Yuling, et al. \"Yes, but did it work?: Evaluating variational inference.\" International Conference on Machine Learning. PMLR, 2018.\n\n### Questions\n\n- Theorem 3.1: So the claim here is that increasing $\\alpha$ reduces $\\mathbb{E} \\|\\Delta V_t\\|^2$? \n    Otherwise, it seems like the tightest bound in 3.1 comes from setting $\\alpha = 0$ and using the diagonal metric.\n    I think a precision definition of $\\Delta V_t$ should be included in in the statement so we can see this trade-off, if it exists.\n    Currently the text only says \"where $\\Delta V_t$ is an operator\" (Theorem 2.1), which isn't helpful. \n\n- Section 3.2: I don't understand your notation. I see that $l$ indexes layers of the model, \n    and $i$ indexes the up to the rank of the tensor; $[H_{t}^i]_l$ is ith \n    tensor component, but what is its shape (is it actually a tensor or just a matrix) \n    and how do you update it with the gradient? \n    Its not clear to me how your \"slight abuse of notation\" works. \n    Also, how do you take the Kronecker product of a $d_i$ tensors? \n    That should produce a tensor whose shape is the product of the component shapes, right? \n\n- Theorem 3.1 and 3.3: What does \"ignoring the $\\Gamma(\\theta)$ terms mean? \n    Does it mean dropping the $\\Gamma$ term in Eq. (3), which you say can be done safely \n    for practical algorithms?\n    Or does it mean dropping the last term in Theorem 2.1 from the analysis?\n\n- Can you provide calibration information for fully connected networks? Log-likelihood is \n    not a good substitute for calibration information, since the MAP estimate, \n    which incorporates no uncertainty information, will maximize the log-likelihood.\n    Thus, log-likelihoods typically favor posterior collapse around the MAP and \n    don't capture uncertainty information very well.\n\n- Table 3: Monge reduces to identity, which seems to match Theorem 3.1.\n    It doesn't seem like choosing $\\alpha^2 > 0$ is supported by any strong \n    evidence---empirical or theoretical.\n    Can you comment on this?", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "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": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This submission proposes a method for sampling from approximate\nBayesian posterior distributions using on stochastic Langevin dynamics.\nThe main contribution is two new approximate curvature matrices \nwhich can be estimated online and used to improve \nsample quality and the exploration ability of the sampler. \nThe first approach, called the Monge metric, uses a diagonal plus rank-one matrix,\nwhere the rank-one component is estimated online using stochastic gradients\nof the log-probability.\nThe second, called Shampoo metric, is based on a Kronecker factored curvature\nmatrix from the optimization literature. \nThe authors characterize the sample efficiency of these two approaches\nand evaluate their practical performance of MNIST and CIFAR-10.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "I would like to see calibration experiments evaluating the quality of uncertainty for the the neural network experiments. I think such experiments are critical to judging the quality of the posterior approximations provided by the sampling methods. Using non-diagonal curvature matrices should lead to improved sampling performance and improved posterior uncertainties. Otherwise, non-diagonal metrics are not worth the additional computational complexity. See Yao et al. [2] for a basic reference on evaluating calibration for approximate inference methods.\n\n", "broader_impact_concerns": "None.", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "jVOCSyNELI", "reviewer": "Reviewer_jecC", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths\n- The paper is well-written with clear motivation and an extensive overview of related methods.\n- The experimental setup is justified. I especially appreciated the toy example in Section 4.2 that shows the strengths of incorporating non-diagonal information.\n- The authors included code and experiment details in the Appendix.\n\nWeakness\n- The improvements for incorporating non-diagonal information seem to be limited. For example, Monge always prefers $\\alpha = 0$ for CIFAR10 experiments and it seems like extensive tuning for $\\alpha$ is necessary. \n- There are some minor concerns on the experiments section, which I address below.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "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": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "The paper proposes two different approaches for using non-diagonal metrics with stochastic gradient samplers. The first approach utilizes a Monge metric, which approximates the Fisher information matrix as the identity plus a low-rank matrix. The second approach exploits Kronecker factors, as in Shampoo, to approximate the Fisher information. Both approaches efficiently support key operations like inversion and square roots. The authors empirically demonstrate that these proposed approaches can yield improvements on MNIST-MLP and CIFAR10-RESNET benchmarks in terms of log-likelihood and accuracy with reasonable computational overhead.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "- In the introduction, authors claim that there is no stochastic-gradient sampler with non-diagonal $G(\\boldsymbol{\\theta})$ for training large neural network. It would be helpful to clarify why methods such as [1] do not fall into this category in the manuscript. At least for the models the authors experiment on, I believe that [1] can be employed.\n- For Funnel problem in Section 4.2, I anticipate that Shampoo-based methods would be more beneficial in capturing the correlations. However, I could not find the Shampoo results and it would be helpful to include it.\n- The accuracies reported in Table 1 for all methods (including the baseline) seems to be lower than I would anticipate [2, 3]. Could the authors clarify why this might be the case? This also applies to the results presented in Table 2 .\n\n[1] Nado, Z., Snoek, J., Grosse, R., Duvenaud, D., Xu, B., & Martens, J. (2018). Stochastic gradient Langevin dynamics that exploit neural network structure.\n\n[2] Blundell, C., Cornebise, J., Kavukcuoglu, K., & Wierstra, D. (2015, June). Weight uncertainty in neural network. In International conference on machine learning (pp. 1613-1622). PMLR.\n\n[3] Fortuin, V., Garriga-Alonso, A., Ober, S. W., Wenzel, F., Rätsch, G., Turner, R. E., ... & Aitchison, L. (2021). Bayesian neural network priors revisited. arXiv preprint arXiv:2102.06571.\n", "broader_impact_concerns": "I don't have any concerns on the ethical implications of the work.", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "Quu9qcUImH", "reviewer": "Reviewer_fKj6", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths\n-----------\n* The paper rigorously supports its claims.\n* Some experiments are given, but there are some shortcomings (see weaknesses).\n\nWeaknesses\n---------------\n* Most of the methodology is effectively preexisting. In particular, the metrics effectively come from preexisting literature. Importantly, these metrics were previously tested on much larger scale experiments (e.g. Shampoo was used to train some large text to video diffusion models, if I remember correctly).\n* Connected to the above: the experiments are rather small scale and toy. They also really don't show that much improvement (e.g. see horseshoe accuracy numbers). Also, the additional time (e.g. from the shampoo metric) is 2-3x, which is also a rather large burden.\n* Why is the Sherman-Morrison Formula intractable here? It has been used to train DEQs (https://arxiv.org/pdf/1909.01377.pdf) for instance, which seem to have much higher dimensionality.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "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": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper applies two (existing) computationally cheap Riemannian metric methods to Riemannian langevin monte carlo sampling. There is some theory (showing that the metric matters beyond the 0th order term), and some experiments show that the proposed methods better sample from the distribution.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "N/A", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper presents two methods for performing Stochastic Gradient Riemannian Langevin Dynamics in larger scale deep networks, where the curvature is modeled using something more sophisticated than a simple diagonal matrix. The methods include a diagonal plus rank-one update and a KFAC-like Shampoo based sampler.  Drawing stochastic samples from large deep networks is certainly of interest to many within the community, for uncertainty quantification, generative modeling, Bayesian deep learning, etc.  The reviewers felt that the claims were all supported by the evidence and all of them agreed on acceptance.  There were some questions regarding treatment of previous work (e.g. I agree that KFAC-SGLD and the Shampoo-based method are rather similar).  I would ask the authors to tone down the wording regarding this comparison in their paper (as they promised in the discussion period).  Given the consensus of the reviewers, the recommendation is accept the paper.  
Decision: Accept (poster)