Markov Chain Score Ascent: A Unifying Framework of Variational Inference with Markovian Gradients
Authors: Kyurae Kim
Jisu Oh
Jacob R. Gardner
Adji Bousso Dieng
Hongseok Kim
Venue: NeurIPS 2022 Accept (2022)
Keywords: variational inference
Bayesian inference
inclusive Kullback-Leibler divergence
Markov chain gradient descent
Markov chain

Abstract: Minimizing the inclusive Kullback-Leibler (KL) divergence with stochastic gradient descent (SGD) is challenging since its gradient is defined as an integral over the posterior. Recently, multiple methods have been proposed to run SGD with biased gradient estimates obtained from a Markov chain. This paper provides the first non-asymptotic convergence analysis of these methods by establishing their mixing rate and gradient variance. To do this, we demonstrate that these methods—which we collectively refer to as Markov chain score ascent (MCSA) methods—can be cast as special cases of the Markov chain gradient descent framework. Furthermore, by leveraging this new understanding, we develop a novel MCSA scheme, parallel MCSA (pMCSA), that achieves a tighter bound on the gradient variance. We demonstrate that this improved theoretical result translates to superior empirical performance.


MARKOV CHAIN SCORE ASCENT:
A Unifying Framework of Variational Inference with Markovian Gradients
Kyurae Kim*
Jisu Oh Sogang University
Sogang University msca@hdsogang. ac. kr
jisohhsoogang ac kr
Jacob R. Gardner University of Pennsylvania jacobrgäseas. upenn. edu
Adji Bousso Dieng Princeton University
Hongseok Kim
Sogang University hongseokdsogang. ac kr
adj1pprinceton.ed
Abstract
Minimizing the inclusive Kullback-Leibler (KI) divergence with stochastic gradient descent (SGD) is challenging since its gradient is defined as an in- tegral over the posterior. Recently, multiple methods have been proposed run SGD with biased gradient estimates obtained from a Markov chain. This paper provides the first non-asymptotic convergence analysis of these methods by establishing their mixing rate and gradient variance. To do this, we demon- strate that these methods-which we collectively refer to as Markov chain score ascent (MCSA) methods-can be cast as special cases of the Markov chain gra- dient descent framework. Furthermore, by leveraging this new understanding. we develop a novel MCSA scheme, parallel MCSA (pMCSA), that achieves a tighter bound on the gradient variance. We demonstrate that this improved theoretical result translates to superior empirical performance.
1 Introduction
Bayesian inference aims to analyze the posterior distribution of an unknown latent variable z from which data is observed. By assuming model p (x z), the posterior (z) is given by Bayes' rule such that x (z) & (xx (z) (z) where (z) represents our prior belief on z. Instead of working directly with π, variational inference (VI, Blei et al. 2017) seeks variational approx- imation q A) E 2, where is variational family and À are the variational parameters, that the most similar according to discrepancy measure d (zz q (::2))
The apparent importance choosing the right discrepancy measure has led to quest spanning decade (Bamler et al., 2017; Dieng et al., 2017; Geffner & Domke, 2021b; Lorbatd et al., 2016; Li & Turner, 2016; Regli & Silva, 2018; Ruiz & Titsias, 2019; Salimans et al., 2015; Wan et al. 2020; Wang et al., 2018; Zhang et al., 2021). So far, the exclusive (or reverse, backward) Leibber (KL) divergence (:; || œ) has seen "exclusive" use, partly because it is defined an integral over q (z; which can be approximated efficiently. In contrast, the inclusive (or forward) KL defined as an integral over π as
(Z) dz Ez-m(-)
(z)
dx Il ())
(z) log
log q (z;
q 2z;
K. Kim is currently with the University of Pennsylvania. JJ. Oh currently with North Carolina State University. *Corresponding author.
36th Conference on Neural Information Processing Systems (NeurIPS 2022).

---

Since our goal is to approximate π with (=:A) but the inclusive KL involves an integral over æ, we end up facing a chicken-and-egg problem. Despite this challenge, the inclusive KL has con- sistently drawn attention due to its statistical properties, such as better uncertainty estimates due its mass covering property (MacKay, 2001; Minka, 2005; Trippe & Turner, 2017).
Recently, Naesseth et al. (2020); Ou & Song (2020) have respectively proposed Markovian score climbing (MSC) and joint stochastic approximation (JSA). These methods minimize the inclu- sive KL using stochastic gradient descent (SGD, Robbins & Monro 1951), where the gradients are estimated using Markov chain. The Markov chain kernel K,, (2, is z-iivariint (Robert & Casella, 2004) and is chosen such that it directly takes advantage of the current variational ap- proximation q (-;2,). Thus, the quality of the gradients improves over time as the KL divergence decreases. Still, the gradients are non-asymploticall biased and Markovian across adjacent il- erations, which sharply contrasts MSC and JSA from classical black-box VI (Kucukelbir et al., 2017; Ranganath et al., 2014), where the gradients are unbiased and independent. While Naes- seth et al. (2020) have shown the convergence of MSC through the work of Gu Kong (1998), this result is only asymptotic and does not provide practical insight into the performance of MSC.
In this paper, we address these theoretical gaps by casting MSC and JSA into a general frame- work we call Markov chain score ascent (MCSA), which we show is a special case of Markov chain gradient descent (MCGD, Duchi et al. 2012). This enables the application of the non- asymptotic convergence results of MCGD (Debavelaere el al., 2021; Doan et al., 2020a,b; Duchi et al., 2012; Karimi et al., 2019; Sun et al., 2018; Xiong et al., 2021). For MCGD methods, the fundamental properties affecting the convergence rate are the ergodic convergence rate (p) of the MCMC kernel and the gradient variance (G). We analyze and G of MSC and JSA, enabling their practical comparison given a fixed computational budget (N). Furthermore, based on the recent insight that the mixing rate does not affect the convergence rate of MCGD (Doan et al., b0,0a,b), we propose a novel scheme, parallel MCSA (pMCSA), which achieves lower variance by trading off the mixing rate. We verify our theoretical analysis through numerical simulations and compare MSC, and pMCSA on general Bayesian inference problems. Our experiments show that our proposed method outperforms previous MCSA approaches.
Contribution Summary
We provide the first non-asymptotic theoretical analysis of two recently proposed inclusive KL minimization methods (Section 4), MSC (Theorems 1 and 2) and JSA (Theorem 3).
To do this, we show that both methods can be viewed as what we call "Markov chain score ascent" (MCSA) methods (Section 3), which are a special case of MCGD (Proposition 1).
3 In light of this, we develop a novel MCSA method which we call parallel MCSA (pMCSA, Section 5) that achieves lower gradient variance (Theorem 4).
We demonstrate that the improved theoretical performance of pMCSA translates to supe- rior empirical performance across a variety of Bayesian inference tasks (Section 6).
2 Background
2.1 Inclusive Leibler Minimization with Stochastic Gradients
VI with SGD The goal of VI is to find the optimal variational parameters À identifying (x E Q that minimizes some discrepancy measure D (π, q (:;2).. A typical way to perform VI is to use stochastic gradient descent (SGD, Robbins & Monro 1951), provided that the opti- mization objective provides unbiased gradient estimates g (2) such that we can repeat the update
2 = 2, - Y/ Ygg( (24-1),
where V1...., .r is a stepsize schedule.
Inclusive KL Minimization with SGD For inclusive KL minimization, g should be set as
g (2) = VA ll (::A) = V2H[T, (::A] = -E2- [2;2))
where HI [z, (:) is the cross- cross-entropy between and q 2:) which shows the connection with cross- entropy methods (de Boer et al., 2005), and s (2: z) = V, log (z; 2) is known as the score gradient. Since inclusive KL minimization with SGD equivalent to ascending towards the direction of the score, Naesseth et al. (2020) coined the term score climbing. To better conform with the optimization literature, we instead call this approach score ascent as in gradient ascent.

---

Table 1: Convergence Rates of MCGD Algorithms
Algorithm
Stepsize Rule
Gradient Assumption
Rate
Reference
Mirror Descent' 7, /
F(lg(nn) 1 G² o
G² logT Ieep-ivT
Duchi et al (2012)
Corollary
7/ = 2/(1+
G² Iog
Doan et al. (2020a)
SGD-Nesterov
(l <G
P. 2LV1+1
vT
Theorem 2
SGD*
y min 24/] Il < G <G +1)
7, = y/t
G2 logT
Doan et al. (2020b)
Theorem 1,2
Notation: JJ is the o o-field formed by the iterates y, up to the rth MCGD iteration. || is the dual norm of x: ?B, is the stepsize of the momentum: 2 is the Lipschitz smoothness constant; Pμ is the strong convexity constant.
2.2 Markov Chain Gradient Descent
Overview of MCGD Markov chain gradient descent (MCGD, Duchi al. 2012; Sun et 2018) is family of algorithms that minimize function defined as f (2) J (A.r) II (dn). where y) is random noise, and II (dh) is its probability measure. MCGD repeats the steps
1 2, - Y/1 M1), n. Pay-1 (n--.)),
(1)
where PA is a Il-invariant Markov chain kernel that may depend 2, The noise of the gra- dient is Markovian and non-nonnnanneee asymptotically biased, departing from vanilla SGD. Non Non-asymon non-asymptotic convergence of this general algorithm has recently started to gather attention as by Debavelaere et al. (2021); Doan et al. (2020a,b); Duchi et al. (2012); Karimi et al. (2019); Sun et al. (2018).
Applications of MCGD MCGD encompasses an extensive range of problems, including dis- tributed optimization (Ram et al., 2009), reinforcement learning (Doan et al., 2020a; Tadić & Doucet, 2017; Xiong et al., 2021), and -ntittitonnmmimizzation (Karimi et al., 2019), to name few. This paper extends this list with inclusive KL VI through the MCSA framework.
Markov Chain Score Ascent
First, we develop Markov chain score ascent (MCSA), framework for inclusive KL minimiza- tion with MCGD. This framework will establish the connection between MSC/JSA and MCGD.
3.1 Markov Chain Score Ascent as a Special Case of Markov Chain Gradient Descent
As shown in Equation (1), the basic ingredients of MCGD are the target function f (2, 17), the gradient estimator (A, 7),, and the Markov chain kernel Pa (77, +). Obtaining MCSA from MCGD boils down designing g and P2 such that f (2) = dgu(2z ( q (+:2)) The following proposition provides sufficient conditions on and P3. to achieve this goal.
Proposition 1. Let y where II is defined as
(z (2)) z (N)) and a Markov chain kernel PA (n, be Iliivariant
(()) π (z) x (2(2) x ... x π (z))..
Then, by defining the objective function f and the gradient estimator g to be [[x) and g (A.. - Σn (z)))
f(2. - Σn logq
where [π] ] the entropy ,f MCGD results in inclusive KL minimization as
En f n) = (x ( |l q (): and En n) ] (x ll (=:
Proof. See the full proof in page 19.
This simple connection between MCGD and VI paves the way toward the non-asymptotic anal- ysis of JSA and MSC. Note that N here can be regarded as the computational budget of each MCGD iteration since the cost of (i) generating the Markov chain samples (1)) ( (N) and (ii) computing the gradient g will linearly increase with N.
In addition, the MCGD framework often assumes P to be geometrically ergodic. An exception is the analysis of Debavelaere et al. (2021) where they work with polynomially ergodic kernels.
3

---

Assumption 1 (Uniform Geometric Ergodicity). The Markov chain kernel P is geometrically ergodic as
drv 11th (n, ....) - Cp""
for some positive constant C.
3.2 Non-Asymptotic Convergence of Markov Chain Score Ascent
Non-Asymptotic Convergence Through Proposition Assumption 1 and some technical assumptions on the objective function, we can apply the existing convergence results of MCGD to MCSA. Table provides list of relevant results. Apart from properties of the objective func- tion (such as Lipschitz smoothness), the convergence rates are stated in terms of the gradient bound G, kernel mixing rate p. and the number of MCGD iterations T. We focus on G and p as they are closely related to the design choices of different MCSA algorithms.
Convergence and the Mixing Rate p Duchi et al. (2012) was the first to provide an analysis of the general MCGD setting. Their convergence rate dependent on the mixing rate through the 1/ log, term. For MCSA, this result is overly conservative since, on challenging problems, mixing can be slow such that p ~ 1. Fortunately, Doan et al. (2020a,b) have recently shown that it is possible to obtain rate independent of the mixing rate For example, in the result of Doan et al. (2020b), the influence of decreases a rate of (1/T). This observation critical since it implies that trading "slower mixing rate" for "lower gradient variance" could be profitable. We exploit this observation in our novel MCSA scheme in Section 5.
Gradient Bound G Except for Doan et al. (2020b), most results assume that the gradient is bounded for Vr, À as (2. || G. Admittedly, this condition is strong, but it is similar to the
bounded variance assumption E [l] < G used in vanilla SGD, which also known to be strong it contradicts strong convexity (Nguyen et al., 2018). Nonetheless, assuming G can have practical benefits beyond theoretical settings. For example, Geffner Domke (2020) use to compare the performance different VI gradient estimators. In similar spirit, we will obtain the gradient bound G of different MCSA algorithms and compare their theoretical performance.
4 Demystifying Prior Markov Chain Score Ascent Methods
In this section, we will show that MSC JSd both qualify as MCSA methods. Furthermore, we establish (1) the mixing rate of their implicitly defined kernel P and (ii) the upper bound on their gradient variance. This will provide insight into their practical non-asymptotic performance.
4.1 Technical Assumptions
To cast previous methods into MCSA, we need some technical assumptions.
Assumption 2 (Bounded importance weight). The importance weight ratio (z) π (z) /4 (z) is bounded by some finite constant as w* < oo for all A E A such that r 11w) < 1.
This assumption necessary to ensure Assumption 1, and can be practically ensured by using variational family with heavy tails (Domke & Sheldon, 2018) or using defensive mixture (Hes- terberg, 1995; Holden et al., 2009) as
Qdet. (z: a q (z; (1-
where 0 a < and v-) is heavy tailed distribution such that sup π (z) (z) co. Note that qdet. is only used in the Markov chain kernels and q (:;2) is still the output of the VI procedure. While these tricks help escape slowly mixing regions, this benefit quickly vanishes as A converges. Therefore, ensuring Assumption 2 seems unnecessary in practice unless we absolutely care about ergodicity. (Think of the adaptive MCMC setting for example. Brofos al. 2022; Holden et al. 2009).
Model (Variational Family) Misspecification and w* Note that w* is bounded below ex- ponentially by the inclusive KL as shown in Proposition 2. Therefore, w* will be large (i) in the initial steps of VI and (ii) under model (variational family) misspecification.
Assumption 3. (Bounded Score) The score gradient is bounded for E A and Vz E 2 such that || (2;z) 1 L some finite constant L 0.
Although this assumption strong, it enables us compare the gradient variance of MCSA methods. We empirically justify the bounds obtained using Assumption 3 in Section 6.2.

---

4.2 Markovian Score Climbing
MSC (Algorithm 4 in Appendix B) is a simple instance of MCSA where nu 2, and Py, = , is the conditional importance sampling (CIS) kernel (originally proposed by Andrieu et al. (2018)) where the proposals are generated from (c) Although MSC uses only single sample for the Markov chain, the CIS kernel internally operates N proposals. Therefore, N in MSC has a different meaning, but it still indicates the computational budget.
Theorem 1. MSC (Naesseth et al., 2020) is obtained by defining
1 (n, dny) K dz'
with 7, 2,, where K2 (z, .) is the CIS kernel with Qaet. (::²) as its proposal distribution. Then, given Assumption 2 and 3, the mixing rate and the gradient bounds are given as
drv (1) (n, I.) = =(1= (1-
N-1
and E || I F s 123
2w 2
where w" sup, (2) qdef. (z; 2).
Proof. See the full proof page 21.
Discussion Theorem shows that the gradient variance of MSC is insensitive to N. Although the mixing rate does improve With when w" is large due to model misspecification and lack of convergence (see the discussion in Section 4.1), this will be marginal. Overall, the performance of MSC cannot be improved by increasing the computational budget N.
Meanwhile, Naesseth et al. also provide a Rao-Blackwellized version MSC we denote as MSC-RB, Instead of selecting single z, by resampling over the N internal proposals, they suggest forming an weeerrance-wwightee estimate (Robert & Casella, 2004). The theoretical properties of this estimator have been concurrently analyzed by Cardoso et al. (2022). Theorem 2. (Cardoso al., 2022) The gradient variance of MSC-RB is bounded as
Pho-Blackwellization
MB Fi-1 4
N-1 dz2 (T ll q(: + o (N 3/2 + /--//--)] + ||A.|
where = S y) 2 (2w* + N 2) is the mixing rate of the Rao- Bao-Blackwellized CIS kernel, and (z q) - (dz) is the x² divergence.
Proof. See the full proof in page 22.
The variance of MSC-RB decreases as O (V/-1), which is more encouraging than vanilla MSC. However, the first term depends on the x² divergence, which is bounded below exponentially by the KL (Agapiou et al., 2017). Therefore, the variance of MSC-RB will large on challenging problems where the x² divergence is large, although linear variance reduction is possible.
4.3 Joint Stochastic Approximation
ISA (Algorithm 5 in Appendix B) was proposed for generative models where the likelihood factorizes into each datapoint. Then, subsampling can be used through random-sea scan version of the independent Netropolis-Hasting (IMH, Hastings 1970) kernel. Instead, we consider the general version of JSA with a vanilla IMH kernel since it can be used for any type of likelihood. At each MCGD step, JSA performs multiple Markov chain transitions and estimates the gradient by averaging all the intermediate states, which is closer to how traditional MCMC' is used.
Independent Metropolis-Hasteee
Similarly to MSC, the IMH kernel in JSA generates pro- posals from q 2.) To show the geometric ergodicity of the implicit kernel P. we utilize the geometric convergence rate of IMH kernels provided by Mengersen & Tweedie (1996, Theorem 2.1) and Wang (2022). The gradient variance, on the other hand, is difficult to analyze, especially the covariance between the samples. However, we show that, even if we ignore the covariance terms, the variance reduction with respect to N is severly limited in the large w/ regime. To do this, we use the exact n step marginal IMH kernel derived by Smith & Tierney (1996) as
KJ (z, dz')==) Tn ( w (z) v (z)) π (z') dz' + 2nh (w (z)) 8g (dz'),
(2)
where w (z) = 7 (z) /qaof. (z; x v y = max (x,y), and for R (v) = {z' | w (z') < {D}
n
- (0) du and 1(w) = [
π (dz'). (1- w (z')
Tn (w) =
(3)
3
R(w)
w
5

---

Theorem 3. JSA (Ou & Song, 2020) is obtained by defining
F dn)) = K 11-1)+1 (z(), dz ()) K2 (k-1)+2 (z(2, dz' (2)).... K2. (k-1))+ (z(N), dz (N)
with 7,, = (z 2 (22) zi ) Then, given Assumption and 3, the mixing rate and the gradient
variance bounds are
drv (R. (n.)). s ((,N) ANN and IF [|BB si-1| S 1?
2 N
+
3 1
+ O (4/m* +Nin) + Coor l|ell2
where μ (a:), Covv Na the covariance between the samples, w' = sup, π (z) /qact (z; 2), and C (r, > 0 is a finite constant.
En Em Cov( ( ss (aiz
ss (A: zi
(m)
F is the sum of
Proof. See the full proof page 24.
Discussion As shown in Theorem 3, JSA benefits from increasing N in terms of faster mix- ing rate. However, under lack of convergence and model misspecification (large w*), the vari- ance improvement becomes marginal. Specifically, in the large w* regime, the variance reduc- tion is limited by the constant 1/2 term. This is true even when, ideally, the covariance between the samples is ignorable such that Cov R O. In practice, however, the covariance term Cov will be positive, only increasing variance. Therefore, in the large w" regime, JSA will perform poorly, and the variance reduction by increasing N is fundamentally limited.
5 Parallel Markov Chain Score Ascent
Algorithm 1: PMCSA
Our analysis in Section 4 suggests that the statistical performance of MSC, MSC-RB, and JSA are heav- ily affected by model specification and the state of convergence through w". Furthermore, for JSA, large w" abolishes our ability to counterbalance the inefficiency by increasing the computational budget N. However, P and G do not equally impact conver- gence, recent results on MCGD suggest that gradient variance is more critical than the mixing rate (see Section 3.1). We turn to leverage this understanding to overcome the limitations of previous methods.
Input: initial samples zo initial parameter Ag number of iterations T,
2o
(N)
stepsize schedule y/
for t 1,2, . T' do
for n = 1, 2, ... N do
K KA,1 (2)
(n)
z
end
1 - S (dizz (1)
g g(2)
2 4-1 -1g (2,-1)
5.1 Parallel Markov Chain Score Ascent
end
We propose novel scheme, parallel Markov chain score ascent (pMCSA, Algorithm 1), that embraces a slower mixing rate in order to consistently achieve an (V/)) variance reduction, even on challenging problems with large w",
Algorithm Description Unlike JSA that uses N sequential Markov chain states, pMCSA op- erates parallel Markov chains. To maintain a similar per-iteration cost with JSA, it performs only single Markov chain transition for each chain. Since the chains are independent, the Metropolis-Hastiee rejections do not affect the variance of pMCSA
Theorem 4. pMCSA, our proposed scheme, is obtained by setting
R²* (1, dn') KŽ (z(), dz (1) K (z(2), dz (2)
K (z(N), dz (N))
with y = (2((1, , (N). Then, given Assumption 2 and 3, the mixing rate and the gradient variance bounds are
drv (1.* (q, II) s (N) , and E E|ller.Mccsal F- | 1
+ IIMI
where μ Exs (2; wt sup, π (z) /uall (z; 2) and C (N)+ 0 is finite constant.
Proof. See the full proof in page 27.
Discussion Unlike JSA and MSC, the variance reduction rate of pMCSA independent of w². Therefore, it should perform significantly better on challenging practical problems. If we consider the rate of Duchi et al. (2012), the combined rate constant with respect N since cancels out. In practice, however, we observe that increasing N accelerates convergence quite dramatically. Therefore, the mixing rate independent convergence rates by Doan et al. (2020a,b) appears to better reflect practical performance. This is because (i) the mixing rate p is conser- vative global bound and (ii) the mixing rate will improve naturally as MCSA converges.

---

10
Variance
10'
PMCSA
USA
MSC
MSC-RB
10
10"
10"
2k 4k 6k 8k
2k 4k 8k10k1 2k 4k 8k10k1 2k 6k 8k10k
Iteration
Iteration
Iteration
Iteration
Figure Gradient variance versus iteration and computational budget (N). pMCSA
not only achieves the least gradient variance, but its variance also scales better with N.
The colors range from light (N = 2²) to dark (N = 2²) representing the computational budgets N [23, 2 29, 2?]. The target distribution is a 50-D multivariate Gaussian with 500.
The error bands are the 80% quantiles obtained from 8 replications.
60
PMCSA (ours)
Final KL
40
JSA
MSC
20
MSC-RB
SGD
Momentum
Nesterov
ADAM
10
10
10"
10"
10
10 10
100 10
10
10"
10
10
10 10
Stepsize (y)
Stepsize (7)
Stepsize (γ)
Stepsize (y)
Figure Optimizer stepsize (y) versus final KL. pMCSA the least sensitive to opti- mizer hyperparameters and results in stable convergence. The final KL is obtained at the 10 iteration. The target distribution is a 100-D Gaussian with v = 500. The error bands are the 80% quantiles, while the solid lines are the median of 20 replications.
5.2 Computational Cost Comparison
Table 2: Computational Costs
The four schemes using the CIS and IMH kernels have different costs depend- ing on N organized in Table
Applying P
Estimating
(z,
( A)
(z:
(z,x) (z:
Cost of Sampling Proposals For the CIS kernel used by MSC, controls the number of internal proposals sampled from q (=: For JSA and pMCSA, the IMH kernel only uses single sample from 2=; but applies the kernel N times. On the other hand, pMCSA needs
# Eval. # Eval. Samples # Grad. # Grad.
ELBO
0
0
N
N
N
MSC
N
1
N
N
1
0
1
MSC-RB
N
1
N
N - 1
0
N
JSA
N N+1
1
N
0
N
pMCSA
N
2N
N
0
N
twice more evaluations of (:; How-
ever, this added cost is minimal since it is dominated by that of evaluating p (z,x).
Cost of Estimating the Score When estimating the score, MSC computes v, log (z; A) only once, while JSA and our proposed scheme compute N times. However, MSC-RB also computes the score times. Lastly, notice that MCSA methods not differentiate through the likelihood (z,x), unlike ELBO maximization, making its per-iteration cost significantly cheaper.
6 Evaluations
6.1 Experimental Setup
Implementation For the realistic experiments, we implemented MCSA methods on top of the Turing (Ge et al,, 2018) probabilistic programming framework. For the variational family, we use diagonal multivariate Gaussians (mean- field family) with the support transformation of Kucukelbir et al. (2017). We the ADAM optimizer by Kingma & Ba (2015) with stepsize of 0.01 in all experiments. The budget is set to = 10 for all experiments unless specified.
Baselines We compare (i) pMCSA (ours, Section 5), (ii) JSA (Ou & Song, 2020), (iii) MSC (Naesseth et al., 2020), (iv) MSC with with Blackwellizztiin (MSC-RB, Naes- seth al. 2020), and (v) evidence lower-bound maximization (ELBO, Kucukelbir et al. 2017; Ranganath et al. 2014) with the path derivative estimator (Roeder et al., 2017).
Available at https: //gi h/ggthtb...... httgg//ttb.coo/eeeeeee wwwp:///eoooooaa////aaaaaaaaa
7

---

Table 3: Test Log Predictive Density Bayesian Neural Network Regression
ELBO
MCSA Variants
D2.
xx Nuain
=
N
pMCSA (ours)
JSA
MSC
MSC-RB
yacht
403
277
-2.45 #0.01
-2.44
#0.01
-2.49
#0.01
-3,00
#0.05 -2.98 #0.04 -2.98 #0.04
concrete
503
927
-3.25
20.01
-3.24
20.01
-3.20
=0.01
-3.33
=0.02 +3.32
20.02 +3.33 #0.02
airfoil
353
1352
-2.53 +0.02
-2.56 #0.02
-2.27
#0.02
-2.51
#0.02 -2.53
+0.01 -2.51 #0.01
energy
503
691
-2.42 10.02
-2.40
#0.02
-1.92
#0.03
-2.38
#0.02 -2.37 #0.02 -2.36 #0.02
wine
653
12
1439
-0.96 10011
-0.96 10.01
-0.95
#0.01
-0.97
+0.01 -0.97
#0.01 -0.97 #0.01
boston
753
14
455
2.72
-2.70 =0.03
-2.69
=0.02
-2.82
=0.02 -2.80
#0.03 -2.78 #0.02
sml
1203
23
3723
-1.32 20.01
-1.25 #0.02
-1.22
#0.01
-1.72
20.01 -1.97 #0.02 -1.95 #0.02
gas
6503
129
2308
-0.06 #0.01 0.13 #0.03
-0.09
#0.02
-0.47 #0.03 -0.47 20.04 -0.50 20.03
1 D3: Dimentionality of a, Dx= Number of features, Nuaini Number of training data points. denotes the 95% bootstrap confidence intervals obtained from replications.
3 Bolded numbers don't have enough evidence to be distinguished from the best performing method under significance threshold (Friedman test with Nemenyi post-hoc test, Demsar 2006).
6.2 Simulations
Setup First, we verify our theoretical analysis on multivariate Gaussians with full-rank covari- ances sampled from Wishart distribution with degrees of freedom (values of are in the figure captions). This problem is challenging since an IMH (used by PMCSA, JSA) or CIS (used by MSC, MSC-RB) kernel with diagonal Gaussian proposal will mix slowly due to large w".
Gradient Variance We evaluate our theoretical analysis of the gradient variance. The vari- ance is estimated from 512 independent Markov chains using the parameters generated by the main MCSA procedure. The estimated variances are shown in Figure 1. We make the follow- ing observations: (i) pMCSA has the lowest variance overall, and consistently benefits from increasing N. (ii) MSC does not benefit from increasing N whatsoever. (iii) MSC-RB does not benefit much from increasing until the x² divergence between π and q (:;2) has become small. (iv) JSA does not benefit from increasing N until q (:: has sufficiently converged (when w* has become small). These results confirm our theoretical analyses in Sections 4 and
Robustness Against Optimizers Since the convergence of most sophisticated SGD optimiz- ers has yet to be established for MCGD, we empirically investigate their effectiveness. The re- sults using SGD (Bottou et al., 2018; Robbins & Monro, 1951), Momentum (Polyak, 1964), Nes- terov (Nesterov, 1983), ADAM (Kingma& & Ba, 2015), and varying stepsizes are shown in Figure 2. Clearly, pMCSA successfully converges for the broadest variety of optimizer settings. Overall, most MCSA methods seem to be the most stable with ADAM, which points out that establishing the convergence of ADAM for MCGD will be a promising direction for future works.
6.3 Bayesian Neural Network Regression
Setup For realistic experiments, we train Bayesian neural networks (BNN, Neal 1996) for regression. We use datasets from the UCI repository (Dua & Graff, 2017) with 90% ran- dom train-test splits and run for T = 5 10* iterations. We use the model and forward propagation method of Hernandez- Adams (2015) with a 50-unit hidden layer (see Appendix C.I).
airfoil
gas
pMCSA
ELBO N=1
ELBO N-10
2
4
2
2
4
11
1
&
Distribution of variational pos- terior mean of the BNN weights. pMCSA re- sults in much less pruning. The density was estimated with a Gaussian kernel and the band- width was selected with Silverman's rule.
Figure 3:
Results The results are shown in Table PMCSA achieves the best performance com- pared all other MCSA methods. Also, its overall performance is comparable to exclusive KL minimization methods (ELBO) unlike other MCSA methods. Furthermore, on airfoil and energy, pMCSA improves over ELBO by 0.29 nat and 0.48 nat. Even on gas where pMCSA did not beat ELBO, its performance is comparable, and dominates all other MCSA methods by roughly 0.4 nat. Additional experimental results, including root mean-square error (RMSE) results and plots with respect to the wall clock time, can be found in Appendix E.1
8

---

Table 4: Test Log Predictive Density Robust Gaussian Process Regression
ELBO
MCSA Variants
D2
Dx
Nuain
N =
pMCSA (ours)
JSA
MSC
MSC-RB
yacht
287
6
277
-3.63 40.02
-3.31 #0.04
-3.29 #0.05
-3.25
#0.04
-3.27 20.05
airfoil
353
6
1352
-3.14 #0.01
-2.63 #0.01
-2.83
#0.04
-2.77
#0.02
-2.73 #0.02
boston
472
13
455
-2.98
20.01
-2.96 #0.02
-3.00
#0.03
-3.00
20.03
-2.96 20.03
energy
703
8
691
-2.75
#0.01
-2.58 20.03
-2.78
#0.04
-2.70
#0.04
-2.72 #0.05
concrete
939
8
927
-3.68 #0.01
-3.49 10.01
-3.69 #0.02
-3.59
#0.04
-3.57 #0.02
wine
1454
11
1439
-1.02
#0.01
-0.94
20.02
-1.04
#0.01
-1.00
=0.02
-0.99 #0.02
gas
2440
128
2308
0.18 #0.02
-0.86 #0.02
-1.10 #0.03
-1.10 20.04 -1.06 #0.02
1 Dy: Dimentionality of a, Dx: Number of features, Nusin Number of training data points.
+ denotes the 95% bootstrap confidence intervals obtained from replications.
3 Bolded numbers don't have enough evidence to be distinguished from the best performing method under .05 significance threshold (Friedman test with Nemenyi post-hoe test, Demšar 2006).
Weight Pruning When using VI, BNNs have been known to underfit data, which Hoffman (2017); MacKay (2001); Trippe & Turner (2017) associated with "weight pruning." That is, the variational posterior converges to the zero-mean prior. Furthermore, Coker et al. (2022); Huix al. (2022) have shown that this is guaranteed to happen under certain conditions. However, these results are strictly based on exclusive KL minimization. On the other hand, Figure 3 shows that pMCSA does not suffer from weight pruning, which suggests that pruning is an artifact of using the exclusive KL. Additional plots are shown in Figure 6 (Appendix E.1).
6.4 Robust Gaussian Process Regression
Setup We train Gaussian processes (GP) with a Student-T likelihood for robust regression. We use datasets from the UCI repository (Dua & Graff, 2017) with 90% random train-test splits. We use the Matérn 5/2 covariance kernel with automatic relevance determination (Neal, 1996) (see Appendix C.2). We run all methods with = 2 10 iterations. For prediction, we use the mode (:: for the hyperparameters and marginalize the latent function over q 2=) (Rasmussen & Williams, 2006). We consider ELBO with only N since differentiating through the likelihood makes its per-iteration cost comparable to MCSA methods with N 10.
Results The results are shown in Table 4. Except for gas, pMCSA achieves better predictive densities than all other methods. This suggests that, overall, the exclusive KL may be less ef- fective in terms of uncertainty quantification for GP posteriors. While ELBO achieves the best performance on gas, among MCSA methods, pMCSA dominates. Our encouraging regression results suggest that incorporating methods such as inducing points (Snelson & Ghahramani, 2005) into MCSA may lead to an important new class of GP models. Additional experimen- tal results, including RMSE results and plots with respect to the wall clock time, can be found in Appendix E.2. Note that in Appendix E.2, MCSA methods appear worse in terms of RMSE compared to the exclusive KL since the inclusive KL puts less probability volume around the posterior mode.
7 Related Works
Inclusive KL minimization Our MCSA framework generalizes MSC (Naesseth al., 2020) and JSA Ou & Song (2020), which are inclusive KL minimization based on SGD and Markov chains. Similar to MCSA is the method of Li et al. (2017). However, the convergence of this method is not guaranteed since it uses short Markov chains, disqualifying for MCSA. Other methods based on biased gradients have been proposed by Bornschein & Rengio (2015); Le (2020), but these are specific for deep generative models. On different note, Jerfel et al. (2021) use boosting instead of SGD minimize the inclusive KI., which gradually builds a complex variational approximation from simple variational family.
Beyond the KL Divergence Discovering alternative divergences for VI has been an active research area. For example, the x² (Dieng et al., 2017), f (Wan et al., 2020; Wang et al., 2018), (Hernandez-Lobhe
et al., 2016; Li & Turner, 2016; Regli Silva, 2018), reguarlized impor- tance ratio (Bamler et al, 2017) divergences have been studied for VI. However, for gradient estimation these methods involve the importance ratio w (z) = x (z), which leads signif- icant variance and low signal-to-noise ratio under model misspecification (Bamler et al,, 2017;

---

Geffner & Domke, 2021a,c). In contrast, under stationarity, the variance of pMCSA is o²/N (a is the variance of the score over the posterior) regardless of model misspecification. Mean- while, Geffner & Domke (2021b); Ruiz & Titsias (2019); Salimans et al. (2015); Zhang et al. (2021) construct implicit divergences formed by MCMC. With the exception of Ruiz & Titsias (2019), most of these approaches aim to maximize auxiliary representations of the classic ELBO. There- fore, their property is likely to be similar to the exclusive KL.
Adaptive MCMC and MCSA As pointed out by Ou & Song (2020), using q (:: within the MCMC kernel makes MCSA structurally equivalent to adaptive MCMC. In particular, Andrieu & Thoms (2008); Brofos et al. (2022); Gabrié et al. (2022); Garthwaite et al. (2016) discuss the use of stochastic approximation in adaptive MCMC. Also, Andrieu & Moulines (2006); Brofos et al. (2022); Giordani & Kohn (2010); Habib & Barber (2019); Holden et al. (2009): Keith et al. (2008); Neklyudov et al. (2019) specifically discuss adapting the propsosal of IMH kernels, and some of them use KL divergence minimization. These methods focus showing ergodicity the samples (n, in our context) not the convergence of the variational approximation q (: A). In this work, we focused the convergence of q (::A), which could advance the adaptive MCMC side of the story.
8 Discussions
This paper presented new theoretical framework for analyzing inclusive KL divergence min- imization methods based on running SGD with Markov chains. Furthermore, we proposed pMCSA, new MCSA method that enjoys substantially low variance. We have shown that this theoretical improvement translates into better empirical performance.
Limitations Our work has three main limitations. Firstly, since our work aims to under- stand existing MCSA methods, it inherits their current limitations. For example, minibatch subsampling challenging for models with factoriatee likelihoods (Naesseth et al., 2020). Secondly, our theoretical analysis in Section 4 requires Assumption 3, which is strong, but re- quired to connect with MCGD. An important future direction would be to relax the assumptions needed by MCGD. Lastly, our MCSA framework does not include models with parameterized posteriors such as variational autoencoders.
Parameterized Posteriors On problems with parameterized posteriors, the target posterior moves around. Therefore, quickly chasing the moving posterior with fast converging MCMC kernels is as important achieving low variance. Because of this, trading bias and variance is less straightforward compared to the "static" setting we consider. Furthermore, the usage of expensive MCMC kernels could be beneficial as suggested by Zhang al. (2022).
Towards Alternative Divergences In Section 6, we have shown that minimizing the inclu- sive KL is competitive against minimizing the exclusive KL on general Bayesian inference prob- lems. Although Dhaka et al. (2021) have shown that the inclusive KL fails on high-dimensional problems, this is only the case under the presence of strong correlations. Before entirely ditching the inclusive KL, it is essential to ask, "how correlated are posteriors really in practice?" Fur- thermore, the true performance of alternative divergences often masked by the limitations of the inference procedure (Geffner& & Domke, 2021a,c). Given that PMCSA significantly advances the best-known performance of inclusive KL minimization, is possible that similar improve- ments could be extracted from other divergences. To conclude, our results motivate further development of better inference algorithms for alternative divergence measures.
Acknowledgments and Disclosure of Funding
We thank Hongseok Yang for pointing us to relevant related work, Guanyang Wang for insightful Eestrooois-Hassting algorithm, Geon Park and Kwanghee
discussions about the independent Choi for constructive comments that enriched this paper, Christian A. Naesseth for comments about Rao-Black Blackwellized MSC, and Gabriel V. Cardoso for comments about BR-SNIS. We also acknowledge the Department of Computer Science and Engineering of Sogang University for providing computational resources.
K. Kim was supported in part by the EPSRC through the Big Hypotheses grant [EP/RO18537/1] J. R. Gardner was supported by NSF award Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Science and ICT under Grant [NRF-2021R) R1 A2C10954355)
H. Kim was supported by the Basic
10

---

References
Agapiou, S., Papaspiliopouloss O., Sanz-Alonso, D., & Stuart, A. M. 2017. Importance Sampling: Intrinsic Dimension and Computational Cost. Statistical Science, 32(3), 405-431.
Andrieu, Christophe, & Moulines, Éric. 2006. On the Ergodicity Properties of Some Adaptive MCMC Algorithms. The Annals of Applied Probability, 16(3).
Andrieu, Christophe, & Thoms, Johannes. 2008. A Tutorial on Adaptive MCMC. Statistics and Computing, 18(4), 343-373.
Andrieu, Christophe, Lee, Anthony, & Vihola, Matti. 2018. Uniform Ergodicity of the Iterated Conditional SMC and Geometric Ergodicity of Particle Gibbs Samplers. Bernoulli, 24(2).
Bamler, Robert, Zhang, Cheng, Opper, Manfred, & Mandt, Stephan. 2017. Perturbative Black Box Variational Inference. In: Advances in Neural Information Processing Systems, vol. 30.
Curran Associates, Inc.
Blei, David M., Kucukelbir, Alp, & McAuliffe, Jon D. 2017. Variational Inference: A Review for Statisticians. Journal of the American Statistical Association, 112(518), 859-877.
Bornschein, Jörg, & Bengio, Yoshua. 2015 (May). Reweighted Wake-Sleep. In: Proceedings of the International Conference on Learning Representations.
Bottou, Léon, Curtis, Frank E., & Nocedal, Jorge. 2018. Optimization Methods for Large-Scale Machine Learning. SIAM Review, 60(2). 223-311.
Brofos, James, Gabrie, Marylou, Brubaker, Marcus A., Lederman. Roy R. 2022. Adaptation of the Independent Metropolis-Hasting, Sampler with Normalizing Flow Proposals. Pages 5949-5986 of: Proceedings of the International Conference on Artificial Intelligence and Statis- tics. PMLR, vol. 151. ML Research Press.
Cardoso, Gabriel, Samsonov, Sergey, Thin, Achille, Moulines, Eric, & Olsson, Jimmy. 2022. BR- SNIS: Bias Reduced Nermlizrma Importance Sampling. Advances in Neural information Processing Systems (to be presented), Sept.
Coker, Beau, Bruinsma, Wessel P., Burt, David R., Pan, Weiwei, & Doshi-Velez, Finale. 2022.
Wide Mean-Field Bayesian Neural Networks the Data. Pages 5276-5333 Proceedings of the International Conference on Artificial Intelligence and Statistics. PMLR.
de Boer, Pieter-Tjerk, Kroese, Dirk P., Mannor, Shie, & Rubinstein, Reuven Y. 2005. A Tutorial on the Cross-Entropy Method. Annals of Operations Research, 134(1), 19-67.
Debavelaere, Vianney, Durrleman, Stanley, & Allassonnière, Stéphanie, 2021. On the Conver- gence of Stochastic Approximations under Subgeometric Ergodic Markov Dynamic. Elec- tronic Journal of Statistics, 15(1).
Demšar, Janez. 2006. Statistical Comparisons of Classifiers over Multiple Data Sets. Journal of Machine Learning Research, 7(1), 1-30.
Dhaka, Akash Kumar, Catalina, Alejandro, Welandawe, Manushi, Andersen, Michael R., Hug- gins, Jonathan, & Vehtari, Aki. 2021. Challenges and Opportunities High Dimensional Variational Inference. Pages 7787-7798 of: Advances in Neural Information Processing Sys- tems, vol. 34. urran Associates, Inc.
Dieng, Adji Bousso, Tran, Dustin, Ranganath, Rajesh, Paisley, John, & Blei, David. 2017. Vari- ational Inference via Upper Bound Minimization. Pages 2729-2738 of: Advances in Neural Information Processing Systems, vol. 30. Curran Associates, Inc
Doan, Thinh T., Nguyen, Lam M., Pham, Nhan H., & Romberg, Justin. 2020a (Oct.). Convergence Rates of Accelerated Markov Gradient Descent with Applications in Reinforcement Learning. Tech. rept. arXiv:2002.02877 [math]. ArXiv.
Doan, Thinh T., Nguyen, Lam M., Pham, Nhan H., Romberg, Justin. 2020b (Apr.). Finite-Time Analysis of Stochastic Gradient Descent under Markov Randomness. Tech. rept. arXiv:2003.10973. ArXiv.
Domke, Justin, & Sheldon, Daniel R. 2018. Importance Weighting and Variational Inference. In: Advances in Neural Information Processing Systems, vol. 31. Curran Associates, Inc
Dua, Dheeru, & Graff, Casey, 2017. UCI Machine Learning Repository.
11

---

Muchi,John C., Agarwal, Alekh, Johansson, Mikael, & Jordan, Michael I. 2012. Ergodic Mirror Descent. SIAM Journal on Optimization, 22((). 1549-1578.
Gabrié, Marylou, Rotskoff, Grant & Vanden-Bijnden Eric. 2022. Adaptive Monte Carlo Aug- mented with Normalizing Flows. Proceedings of the National Academy of Sciences, 119(10),
e2109420119.
Garthwaite, P. H., Fan, Y., & Sisson, S. A. 2016. Adaptive Optimal Scaling of Seetropolis-Haasting Algorithms Using the dobbins-Monro Process. Communications in Statistics Theory and Methods, 43(17), 5098-5111.
Ge, Hong, Xu, Kai, & Ghahramani, Zoubin. 2018. Turing: A Language for Flexible Probabilis- tic Inference. Pages 1682-1690 of: Proceedings of the International Conference on Machine Learning. PMLR, vol. 84, ML Research Press.
Geffner, Tomas, & Domke, Justin. 2020. A Rule for Gradient Estimator Selection, with an Ap- plication to Variational Inference. Pages 1803-1812 of: Proceedings of the International Con- ference on Artificial Intelligence and Statistics. PMLR, vol. 108. ML Research Press.
Geffner, Tomas, & Domke, Justin. 2021a. Empirical Evaluation of Biased Methods for Alpha Divergence Minimization. In: Proceedings of the Sympusium on Advances in Approximate
Bayesian Inference.
Geffner, Tomas, & Domke, Justin. 2021b. MCMC Variational Inference via Uncorrected Hamil- tonian Annealing. Pages 639-651 of: Advances in Neural Information Processing Systems, vol. 34. Curran Associates, Inc.
Geffnnr, Tomas, & Domke, Justin. 2021c. On the Difficulty of Unbiased Alpha Divergence Min- imization. Pages 3650-3659 of: Proceedings of the International Conference on Machine Learn- ing. PMLR, vol. 139. ML Research Press.
Giordani, Paolo, & Kohn, Robert. 2010. Adaptive Independent Metropolis-Hastinge by Fast Estimation of Mixtures of Normals. Journal of Computational and Graphical Statistics, 19(2), 243-259.
Gu, Ming Gao, & Kong, Fan Hui. 1998. A Stochastic Approximation Algorithm with Markov Chain Monte-Care Method for Incomplete Data Estimation Problems. Proceedings of the National Academy of Sciences, 95(13), 7270-7274.
Habib, Raza, & Barber, David. 2019. Auxiliary Variational MCMC. In: Proceedings of the Inter- national Conference on Learning Representations.
Hastings, W. K. 1970. Monte Carlo Sampling Methods Using Markov Chains and Their Appli-
cations. Biometrika, 57(1), 97-109.
Hernandez-Lobato Jose, Li, Yingzhen, Rowland, Mark, Bui, Thang, Hernandez-Lobato. Daniel, & Turner, Richard. 2016. Black-Box Alpha Divergence Minimization. Pages 1511-1520 of: Proceedings of the International Conference Machine Learning. PMLR, vol. 48. ML Research Press.
Hernandez- Jose Miguel, & Adams, Ryan. 2015. Probabilistic Backpropagation for Scal- able Learning of Bayesian Neural Networks. Pages 1861-1869 of: Proceedings of the Interna- tional Conference on Machine Learning. PMLR, vol. 37. ML Research Press.
Hesterberg, Tim. 1995. Weighted Average Importance Sampling and Defensive Mixture Distri- butions. Technometrics, 37(2), 185-194.
Hoffman, Matthew D. 2017. Learning Deep Latent Gaussian Models with Markov Chain Monte Carlo. Pages 1510-1519 of: Proceedings of the International Conference on Machine Learning. PMLR, vol. 70. ML Research Press.
Holden, Lars, Hauge, Ragnar, & Holden, Marit. 2009. Adaptive Independent Metropolis-
Hastings. The Annals of Applied Probability, 19(1.
Huix, Tom, Majewski, Szymon, Durmus, Alain, Moulines, Eric, & Korba, Anna. 2022 (July). Variational Inference of Overparameterized Bayesian Neural Networks: Theoretical and Em- pirical Study. Tech. rept. arXiv:2207.03855, ArXiv.
Jerfel, Ghassen, Wang, Serena,
Clara, Heller, Katherine A., Ma, Yian, & Jor- dan, Michael I. 2021. Variational Refinement for Importance Sampling Using the Forward Kullback-Leible Divergence. Pages 1879-1822 of: Proceedings of the International Conference on Uncertainty in Artifical Intellience. PMLR, vol. 161. ML Research Press.
12

---

Karimi, Belhal, Miasojedow, Blazej, Moulines, Eric, & Wai, Hoi-To. 2019. Non-Asympiotic Anal- ysis of Biased Stochastic Approximation Scheme. Pages 1944-1974 of: Proceedings of the An- nual Conference on Learning Theory. PMLR, vol. 99. ML Research Press.
Keith, Jonathan M., Kroese, Dirk P., & Sofronov, George Y. 2008. Adaptive Independence Sam- plers. Statistics and Computing, 18(4), 409-420.
Kingma, Diederik P., & Ba, Jimmy. 2015. Adam: A Method for Stochastic Optimization. In: Proceedings of the International Conference on Learning Representations.
Kucukelbir, Alp. Tran, Dustin, Ranganath, Rajesh, Gelman, Andrew, & Blei, David M. 2017. Au- lomatic Differentiation Variational Inference. Journal of Machine Learning Research, 18(14),
1-45.
Le, Tuan Anh, Kosiorek, Adam R., Siddharth, N., Teh, Yee Whye, & Wood, Frank. 2020. Revis- iting Reweighted Wake-Sleep for Models with Stochastic Control Flow. Pages 1039-1049 of: Proceedings of the International Conference on Uncertainty Artifical Intelligence. PMLR, vol.
115. ML Research Press.
Li, Yingzhen, & Turner, Richard E. 2016. Rényi Divergence Variational Inference. In: Advances in Neural Information Processing Systems, vol. 29, Curran Associates, Inc.
Li, Yingzhen, Turner, Richard E., & Liu, Qiang. 2017 (May). Approximate Inference with Amor- tised MCMC. Tech. rept. arXiv:1702.0834 [cs, stat]. ArXiv.
MacKay, David J.C. 2001 (June). Local Minima, Symmetry-Breakke ational Free Energy Minimization. Technical Report.
and Model Pruning in Vari-
Mengersen, K. L.,, & Tweedie, R. L. 1996. Rates of Convergence of the Hastings and Metropolis Algorithms. The Annals of Statistics, 24(1), 101-121.
Minka, Tom. 2005 (Jan.). Divergence Measures and Message Passing. Tech, rept. MSR-TR-2005- 173. Microsoft Research.
Naesseth, Christian, Lindsten, Fredrik, & Blei, David. 2020. Markovian Score Climbing: Vari- ational Inference with KL((Plq). Pages 15499-15510 of: Advances in Neural Information Pro- cessing Systems, vol. 33. Curran Associates, Inc.
Neal, Radford M. 1996. Bayesian Learning for Neural Networks. Lecture Notes in Statistics, vol. 118. New York, NY: NY: Springer New York.
Neklyudov, Kirill, Egorov, Evgenii, Shvechikov, Pavel, & Vetrov, Dmitry, 2019 (June). Metropolis- Hastings View on Variational Inference and Adversarial Training. Tech. rept. arXiv:1810.07151 [cs, stat]. ArXiv.
Nesterov, Yurii Evgen evich. 1983. A Method of Solving a Convex Programming Problem with Convergence Rate O(1/k²). Doklady Akademii Nauk SSSR, 269(3). 543-547.
Nguyen, Lam, Nguyen, Phuong Ha, van Dijk, Marten, Richtarik, Peter, Scheinberg, Katya, & Takac, Martin. 2018. SGD and Hogwild! Convergence without the Bounded Gradients As- sumption. Pages 3750-3758 of: Proceedings of the International Conference on Machine Learn- ing. PMLR, vol. 80. ML Research Press.
Ou, Zhijian, & Song. Yunfu. 2020. Joint Stochastic Approximation and Its Application to Learn- ing Discrete Latent Variable Models. Puges 929-938 of: Proceedings of the International Con- ference on Uncertainty in Artifical Intelligence. PMLR, vol. 124. ML Research Press.
Polyak, B.T. 1964. Some Methods of Speeding up the Convergence of Iteration Methods. USSR Computational Mathematics and Mathematical Physics, 4(5), 1-17.
Ram, S. Sundhar, Nedić, A., & Veeravalli, V. V. 2009. Incremental Stochastic Subgradient Algo- rithms for Convex Optimization. SIAM Journal on Optimization, 20(2), 691-717.
Ranganath, Rajesh, Gerrish, Sean, & Blei, David. 2014. Black Box Variational Inference. Pages
814-822 of: Proceedings of the International Conference on Artificial Intelligence and Statistics. PMLR, vol. 33. ML Research Press.
Rasmussen, Carl Edward, & Williams, Christopher K. 2006. Gaussian Processes for Machine Learning. Adaptive Comput. Mach. Learn. Cambridge, Mass: MIT Press.
Regli, Bean-Baptiste & Silva, Ricardo. 2018 (May). Alpha-Bete Divergence for Variational Infer- ence. Tech. rept. arXiv:1805.01045. ArXiv.
13

---

Robbins, Herbert, & Monro, Sutton. 1951. A Stochastic Approximation Method. The Annals of Mathematical Statistics, 22(3), 400-407.
Robert, Christian & Casella, George. 2004. Monte Carlo Statistical Methods. Springer Texts in Statistics. New York, NY: Springer New York.
Roeder, Geoffrey, Wu, Yuhuai, & Duvenaud, David K. 2017. Sticking the Landing: Simple, Lower-Variance Gradient Estimators for Variational Inference. In: Advances in Neural Infor- mation Processing Systems, vol. 30. Curran Associates, Inc.
Ruiz, Francisco, & Titsias, Michalis. 2019. Contrastive Divergence for Combining Variational Inference and MCMC. Pages 5537-5545 of: Proceedings of the International Conference on Machine Learning. PMLR, vol. 97. MI Research Press.
Salimans, Tim, Kingma, Diederik, Welling, Max. 2015. Markov Chain Monte Carlo and Vari- ational Inference: Bridging the Gap. Pages 1218-1226 of: Proceedings of the International Conference on Machine Learning. PMLR, vol. 37. ML Research Press.
Smith, Richard L., & Tierney, Luke. 1996. Exact Transition Probabilities for the Independence Metropolis Sampler. Tech.. rept.
Snelson, Edward, & Ghahramani, Zoubin. 2005. Sparse Gaussian Processes Using Pseudo- Inputs. In: Advances in Neural Information Processing Systems, vol. 18. MIT Press.
Sun, Tao, Sun, Yuejiao, & Yin, Wotao. 2018. On Markov Chain Gradient Descent. In: Advances in Neural Information Processing Systems, vol. 31. Curran Associates, Inc.
Tadić, Vladislav B., & Doucet, Arnaud. 2017. Asymptotic Bias of Stochastic Gradient Search.
The Annals of Applied Probability, 27(6).
Trippe, Brian, Turner, Richard. 2017. Overpruning in Variational Bayesian Neural Networks.
Tech. rept. arXiv:1801 800.06230.... ArXiv.
Wan, Neng, Li, Dapeng, & Hovakimyan, Naira. 2020. F-Divergence Variational Inference. Pages 17370-1732 17379 of: Advances in Neural Information Processing Systems, vol. 33, Curran Asso- ciates, Inc
Wang, Dilin, Liu, Hao, & Liu, Qiang. 2018. Variational Inference with Tail-Adaptive f- Divergence. In: Advances in Neural Information Processing Systems, vol. 31. Curran Asso- ciates, Inc.
Wang.
Guanyang. 2022. Exact Convergence Rate Analysis of the Independent Metropolis-
Hastings Algorithms. Bernoulli, 28(3), 2012-2033.
Xiong, Huaqing, Xu, Tengyu, Liang, Yingbin, & Zhang, Wei. 2021. Non-Asymptotic Conver- gence of Adam-Type Reinforcement Learning Algorithms under Markovian Sampling. Pro- ceedings of the AAAI Conference on Artificial Intelligence, 35(12), 10460-10468.
Zhang, Guodong, Hsu, Kyle, Li, Jianing, Finn, Chelsea, & Grosse, Roger B. 2021. Differentiable Annealed Importance Sampling and the Perils of Gradient Noise. Pages 19398-19411 19410 of: Ad- vances in Neural Information Processing Systems, vol. 34. Curran Associates, Inc.
Zhang, Liyi, Blei, David M., & Naesseth, Christian A. 2022 (Sept.). Transport Score Climb- ing: Variational Inference Using Forward KL and Adaptive Neural Transport. Tech. rept. arXiv:2202.018 ArXiv.
14

---

Checklist
1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope? [Yes]
See the contribution summary in Section 1.
(b) Did you describe the limitations of your work? [Yes] See Section 8.
(c) Did you discuss any potential negative societal impacts of your work? [No]
(d) Have you read the ethics review guidelines and ensured that your paper conforms
to them? [Yes]
2. you are including theoretical results...
(a) Did you state the full of assumptions of all theoretical results? [Yes] The key assumptions are stated in each proof.
(b) Did you include complete proofs of all theoretical results? [Yes]
See Appendix D.
3. you ran experiments...
(a) Did you include the code, data, and instructions needed to reproduce the main experimental results (either in the supplemental material or as a URL)? [Yes] included in the supplementary material.
(b) Did you specify all the training details (e.g., data splits, hyperparameters, how they were chosen)? [Yes]
See Section 6. Additional details can be found in the code in the supplementary ma- terial.
(c) Did you report error bars (e.g., with respect to the random seed after running periments multiple times)? [Yes]
See the text in Section 6.
(d) Did you include the total amount of compute and the type of resources used (e.g.. type of GPUs, internal cluster, or cloud provider)? [Yes] See Appendix
4. If you are using existing assets (e.g., code, data, models) curating/releasing sets...
(a) If your work uses existing assets, did you cite the creators? [No]
(b) Did you mention the license of the assets? [No]
(c) Did you include any new assets either in the supplemental material or as a URL?
[N/A]
(d) Did you discuss whether and how consent was obtained from people whose data you're 'e're using/curating? [N/A]
(e) Did you discuss whether the data you are using/curating contains personally iden-
tifiable information or offensive content? [N/A]
5. If you used crowdsourcing or conducted research with human subjects...
(a) Did you include the full text of instructions given to participants and screenshots, if applicable? [N/A]
(b) Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals, if applicable? [N/A]
(c) Did you include the estimated hourly wage paid to participants and the total amount spent on participant compensation? [N/A]
15


Reviews:
[{"review_id": "ZWMOeOVdJnh", "reviewer": "Reviewer_ceLV", "summary": "The paper examines variational inference with the inclusive KL, KL(/pi, q), where required samples from /pi are obtained using MCMC.\nIn particular, the authors provide theoretical analysis of non-asymptotic convergence of two published methods, and establish mixing rates and bounds on gradient variance. \nThis is done by casting those methods as instances of Markov Chain Gradient Descent, and applying and adapting existing theoretical results.\nFurther based on their findings, the authors propose a novel scheme that trades off mixing rate with for reduced variance, and through that performs favorably.\nA number of experiments on simple real-world data confirm the theoretical findings.\n", "questions": "One question that is not really discussed in the paper is: \"Is it worth giving up the scalability of the standard variational inference for the improved performance of the discussed methods?\". I think since all papers in the MCSA space lack large-scale experiments, this really remains to be seen, but I'd like to hear the authors comment on this.\n\n289 \"We consider ELBO with only 𝑁 = 1 since differentiating through the likelihood makes its per-iteration cost comparable to MCSA methods with 𝑁 = 10.\" Should be elaborated as it is not clear where this number comes from.\n\n294 \"Our encouraging regression results suggest that incorporating methods such as inducing points (Snelson & Ghahramani, 2005) into MCSA may lead to an important new class of GP models.\" not clear what that means and what the basis of this is.", "limitations": "One problem with the results are the strong assumptions necessary (the authors do state this), however, it seems very hard to derive comparable results without any assumptions, and the experiments validate some of the findings.\nIn total, the paper appears to be honest about strengths and weaknesses of the analysis method.", "rating": 7, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "This is a very nice paper, well written, tackling a worthwhile problem, and executed with a very high level of sophistication.\n\nThe paper begins by casting two existing inclusive KL score climbing methods under the framework of Markov Chain Gradient Descent, a method to optimize an expectation of a function over an intractable distribution that works via using MCMC to draw samples from the distribution.\nThis is done via Proposition 1, which is a relatively simple re-writing exercise of the involved quantities.\nDespite being simple, this connection allows the authors to apply and adapt existing theoretical results to establish results on mixing and gradient variance (Thm1/2).\nThe authors discuss the implications, and in particular the downsides of the both analysed methods.\n\n\nOriginality:\nThe authors connect existing works to the established field of Markov Chain Gradient Descent, a method to optimize an expectation of a function over an intractable distribution that works via using MCMC to draw samples from the distribution.\nThis is done via Proposition 1, which is a relatively simple re-writing exercise of the involved quantities.\nDespite being simple, this highly original re-casting allows the authors to apply and adapt existing theoretical results to establish results on mixing and gradient variance (Thm1/2).\nOn top of this, the authors take their findings and use them to motivate an original method that overcomes shortcomings of the analysed methods.\n\n\nQuality: \nThe whole paper is carried out at a very high standard.\nAll claims are supported by proofs and their implications are discussed in the main text.\nI appreciate the discussions of the theoretical results in plain text.\nIt is particularly nice that the analysis of the existing work leads to a new algorithm that performs favorably in practice -- it is very nice to have theoretical results motivate new methodology.\nFor the new method, the authors provide a comparison of computation costs -- though it would have been nice to see this confirmed in the experiments, by e.g. plotting performance as a function of wallclock time instead of iterations.\nIt is also great that the authors compare various optimizers and step-sizes.\n\n\nClarity:\nThe paper is very well written and organized.\nThe proof contains justifications for most steps, which makes it very nice to read the derivations.\n\n\nSignificance: \nThe results are significant from a theoretical / unification perspective, as well is the proposed algorithm from a practical perspective (due to its strong performance).\nThe only real critique I have for the paper is the lack of modern large-scale experiments.\nThat is, deep generative models with complex architectures, e.g. on images, 3d scenes or scientific domains, and how MCSA perform here in general.\nThe experiments considered are pretty much all more or less solved problems with only marginal gains to be made.\nA strong large-scale results would significantly improve the impact of this paper.\n\n114 missing s in iterations\neq under 162, should say that L is score bound\nline 284, verb missing (\"use\"?)\n171 alg 4 is in appendix", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "6vBW3I2NFaK", "reviewer": "Reviewer_seMc", "summary": "This paper presents a unified framework of Markov Chain Score Ascent (MCSA) methods and a theoretical study of MCSA gradient variance for minimizing the forward KL divergence in variational inference. Based on the theoretical analysis, the paper proposes a simple modification to use parallel Markov chains, named parallel MCSA (pMCSA). Experiment results show that pMCSA achieves better empirical performance for a number of Bayesian inference problems. ", "questions": "- Are there scenarios where the mixing rate is important and JSA could attain better performance?\n- Does pMCSA also trade bias for variance during early training iterations? \n- It would be beneficial to also plot MSC-RB in Figure 1, 2 to match the computational costs between different methods. ", "limitations": "The authors adequately addressed the limitations of their work. The authors did not discuss any potential negative societal impacts of the work, but as the work is mostly theoretical it is unlikely to have immediate negative societal impacts.    ", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths:\n- The paper unifies previous works on forward KL minimization under the same framework.\n- The study on the improvements of gradient variance w.r.t. $N$ for prior methods is useful and insightful.\n- The proposed pMCSA algorithm overall achieves better empirical performance than JSA, MSC and MSC-RB. \n\nWeaknesses:\n- Analysis is largely based on the more general Markov chain gradient descent methodology, of which MCSA is only a special case. This limits the novelty and scope of the work. \n- The proposed modification is relatively simplistic and similar to other parallel MCMC algorithms. Switching the kernel to MH in MSC and increasing sample size also effectively yields the same algorithm. \n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "56U-l452I6-", "reviewer": "Reviewer_mF4s", "summary": "This paper provides a non-asymptotic convergence analysis of methods that minimize the inclusive Kullback-Leibler (KL) divergence with stochastic gradient descent (SGD). The authors prove that these Markov chain score ascent methods are special cases of a unified framework. Based on this framework, a new scheme is proposed.", "questions": "None.", "limitations": "Yes, the authors have addressed the limitations and potential negative societal impact of their work.", "rating": 6, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Unfortunately, I am somewhat of an outsider to the field of this paper (cannot check the correctness of the main part such as the non-asymptotic convergence of the proposed framework) and I feel I cannot accurately assess the importance of the proposed framework and the novel scheme. However, this paper is definitely well-organized and technically solid.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "20sKzh4nc2J", "reviewer": "Reviewer_PECh", "summary": "The authors provide non-asymptotic convergence rates of Markov Chain Score Ascent (MCSA) methods by casting them as instances of Markovian Chain Gradient Descent for which such rates have been recently established. Based on these results the authors hypothesize that the convergence of MCSA depends more heavily on the gradient variance than on the mixing rate of the Markov Chain and propose a parallel MCSA (pMCSA) scheme which leverages these insights. In contrast to JSA which performs $N$ transition on a single chain, pMCSA performs a single transition on $N$ independent Markov chains, hence lowering the variance of the gradient in exchange for a slower mixing rate.\n", "questions": "- Based on the analysis JSA and MSC are both heavily influenced by model misspecification (of the variational approximation) via $w^*$ while pMCSA is independent of $w^*$. Is this effect expected to vanish with increasing capacity of the variational distribution (with reasonable initialization) after some initial training time?\n\nRegarding the convergence rates (mentioned above):\n- I was specifically wondering about the $log T$ term in the numerator of the rate for Mirror Descent which I was not able to match to the rate stated in the original paper.\n- I was also confused by the statement that the influence of mixing rate on the convergence rate decreases as $O(1/T^2)$ even though the stated convergence rate does depend on the mixing rate at all. I might be missing something here and would be grateful for clarification.\n\nSome minor thing I noticed:\n- Equation after line 97. Product sign between $\\pi(z^{(1)}$ and $\\pi(z^{(2)}$ is missing\n- The gradient bound $G$ is mentioned on Line 114 but has not yet been introduced\n- Line 165L MSA -> MSC\n", "limitations": "The authors clearly state their assumptions and limitations.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "**Originality:**\nThe work is original in that it is, to the best of my knowledge, the first that discusses non-asymptotic convergence and provides convergence bounds for MCSA methods. The paper casts MCSA methods as MCGD to leverage recently established results on non-asymptotic convergence for MCGD. Hence the theoretical contribution lies in deriving bound for the gradient variance and mixing rate of the individual methods such that they can be used with the existing theoretical findings. I believe this work is interesting to the inference community and might help to inform the further development of MSCA algorithms and their analysis.\n\n**Clarity:**\nGenerally, due to the technical nature of the paper it was not always straightforward to me to follow the argument. I found it especially hard to follow the discussion about the convergence rates. The authors are mainly interested in the rates in terms of the gradient bound $G$, mixing rate $\\rho$, and number of iterations $T$ and hence can express the convergence rates in a simplified form in O-notation. Even though the authors gave references to the exact theorems and corollaries which establish the convergence rates in previous work it was not straightforward for me to verify the exact correspondences, especially for the Mirror descent case (please see questions below). \n\n**Significance**\nGiven the theoretical nature of the paper, I think the empirical evaluation sufficiently demonstrates the improvement of pMCSA over JSA and MSC. The first experiment seems to confirm the theoretical results but only considers two different sampling budgets. If computationally feasible, it would be nice to see an array of different sampling budgets, e.g. N=8,16,32,62,128, to further confirm the theoretical findings empirically.\n\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}]

Author rebuttal:

Decision: Accept