Challenges and Opportunities in High Dimensional Variational Inference
Authors: Akash Kumar Dhaka
Alejandro Catalina
Manushi Welandawe
Michael Riis Andersen
Jonathan H. Huggins
Aki Vehtari
Venue: NeurIPS 2021 Poster (2021)
Keywords: Variational Inference
BBVI
Importance Sampling
Divergences
Pre-asymptotics

Abstract: Current black-box variational inference (BBVI) methods require the user to make numerous design choices – such as the selection of variational objective and approximating family – yet there is little principled guidance on how to do so. We develop a conceptual framework and set of experimental tools to understand the effects of these choices, which we leverage to propose best practices for maximizing posterior approximation accuracy. Our approach is based on studying the pre-asymptotic tail behavior of the density ratios between the joint distribution and the variational approximation, then exploiting insights and tools from the importance sampling literature. Our framework and supporting experiments help to distinguish between the behavior of BBVI methods for approximating low-dimensional versus moderate-to-high-dimensional posteriors. In the latter case, we show that mass-covering variational objectives are difficult to optimize and do not improve accuracy, but flexible variational families can improve accuracy and the effectiveness of importance sampling – at the cost of additional optimization challenges. Therefore, for moderate-to-high-dimensional posteriors we recommend using the (mode-seeking) exclusive KL divergence since it is the easiest to optimize, and improving the variational family or using model parameter transformations to make the posterior and optimal variational approximation more similar. On the other hand, in low-dimensional settings, we show that heavy-tailed variational families and mass-covering divergences are effective and can increase the chances that the approximation can be improved by importance sampling. 

Challenges and Opportunities in High-dimensional Variational Inference
Akash Kumar Dhaka* Aalto University, Silo AI akash, dhakaßaalto..ft
Alejandro Catalina*
Aalto University
alejandro watalinaaaalllee..cmca.ammmaaa
Manushi Welandawe
Boston University
manushi.. welandawe@bu.een
Michael Riis Andersen Technical University of Denmark miri@dtu.de dk
Jonathan H. Huggins Boston University huggins@bu. edu
Aki Vehtari Aalto University aki vehtar veetari@@alto .fi
Abstract
Current black-box variational inference (BBVI) methods require the user to make numerous design choices -uch as the selection of variational objective and ap- proximating family-yet there is little principled guidance on how to do so. We develop a conceptual framework and set of experimental tools to understand the effects of these choices, which we leverage to propose best practices for max- imizing posterior approximation accuracy. Our approach based on studying the pre-asymptotic tail behavior of the density ratios between the joint distribu- tion and the variational approximation, then exploiting insights and tools from the importance sampling literature. Our framework and supporting experiments help to distinguish between the behavior of BBVI methods for approximating low-dimensional versus moderate-to- hhiit--to-iig--ddniiioaaa posteriors. In the latter case, we show that mass-covering variational objectives are difficult to optimize and do not improve accuracy, but flexible variational families improve accuracy and the effectiveness of importance sampling-al at the cost of additional optimization challenges. Therefore, for aadderatett--hgh--dmeeeee
posteriors recommend using the (mode-seeking) exclusive KL divergence since it is the easiest to optimize, and improving the variational family using model parameter transformations to make the posterior and optimal variational approximation more similar. On the other hand, in din-eninniional settings, we show that heavy-tailed variational families and mass-covering divergences are effective and can increase the chances that the approximation can be improved by importance sampling.
Introduction
A great deal of progress has been made in black-box variational inference (BBVI) methods for Bayesian posterior approximation, but the interplay between the approximating family, divergence measure, gradient estimators and stochastic optimizer is non-trivial and more for high- dimensional posteriors
10 31). While the main focus in the machine learning literature has
Equal contribution.
35th Conference on Neural Information Processing Systems (NcurIPS 2021).

---

Exclusive KL
Inclusive Target
Exclusive KL
Inclusive KL
Target
2
5
10
50
0.6
0.4
0.2
0
0.0
0.6
0.4
0.2
0.0
10
15
10
15
5
10
15
0
10
-2
0
Distance from the mode
(a)
(h)
and inclusive (mass-covering) divergences. (a) The typical 2D illustration (correlation 0.9) gives the impression that the inclusive divergence would provide better approximation. (b) For correlated Gaussian targets in dimensions D 2, 5, 10, 50, the marginal distributions of the distance from the mode for samples drawn from the approximation (red) and the target (blue). The intuition from the low-dimensional examples does not carry over to higher dimensions: although the importance ratios still bounded, even for a lower correlation level, the overlap in typical sets of the target and the approximations gets worse both for exclusive and inclusive divergences.
Figure Illustration of mean-field approximation with exclusive
been on improving predictive accuracy, the choice of method components becomes even more critical when the goal is obtain accurate summaries of the posterior itself.
In this paper, we show that, while the choice approximating family and divergence is often motivated by low-dimensional illustrations, the intuition from these examples do not necessarily carry over to higher-dimensional settings. By drawing a connection between importance sampling and the estimation of common divergences used in BBVI. for understanding the reliability of BB VI in terms of the ae-aasmmptoic behavior of the density ratio between the target and the approximate distribution. When this density ratio is heavy-tailed, unbiased estimators exhibit a large bias with high probability, in addition to high variance. Such heavy tails occur when there is mismatch between the typical sets of the approximating and target distributions. In higher dimensions, even over-dispersed distributions miss the typical set of the target [18 27. Thus, illustrated in Fig. [ the benefits of heavy- tailed approximate families and divergences favoring mass-covering diminish as dimensionality of the target distribution increases. Building on these insights, make the following main contributions:
able to develop comprehensive framework
1. We develop a conceptual and experimental framework for predicting and empirically evaluating the reliability of BBVI based on the choice of variational objective, approximating family, and target distribution. Our framework also incorporates the Pareto diagnostic [27] as a simple and practical approach for obtaining empirical and conceptual insights into the asymptott convergence rates estimators of common divergences and their gradients.
2. We validate our framework through extensive empirical study using simulated data and many commonly used real datasets with both Gaussian and non-Gaussian target distributions. We consider the exclusive and inclusive Kullback-Leeille (KI. divergences 21 24) tail-adaptive f-divergence [29], x² divergence [7], and a-divergences (12), and the resulting variational approximation for isotropic Gaussian and Student-t and normalising flow families.
3. Based on our framework and numerical results, we provide justified recommendations design choices for different scenarios, including low- to noderate-dimeessionll and high-dimensional posteriors.

---

2 Preliminaries and Background
Let p(0, Y) be joint distribution of a probabilistic model, where 0 TRD is vector of model parameters and Y is the observed data. In Bayesian analysis, the posterior p(0) = p(0 Y) = p(0, Y))p(Y) (where p = S [00, Y dd0 ) is typically the object of interest, but most posterior summaries of interest are not accessible because the normalizing integral, in general, is intractable.
Variational inference approximates the exact posterior p(0 Y) using a distribution g E from family of tractable distributions Q. The best approximation is determined by minimizing a divergence D(p q), which measures the discrepancy between p and q:
qx- arg min Dpp || q),
qxE
where E R.K is a vector parameterizing the variational family Q. Thus, the properties of the resulting approximation are determined by the choice of variational family Q as well as the choice of divergence D.
The family Q is often chosen such that quantities of interest (e.g., moments of q) can be computed efficiently. For example, q can be used to compute Monte Carlo or importance sampling estimates of the quantities of interest. Let (0) = p(0, Y) \)/((0) denote the density ratio between the joint and approximate distributions. For a function R + R. the biased self-normalized importance sampling estimator for the posterior expectation Ee [0))] is given by
w(O))
[(0) 2 w( (0.)).
where .I..... As ~ q are independent. Using importance sampling can allow for computation of more accurate posterior summaries and to go beyond the limitations of the variational family. For example, it possible to estimate the posterior covariance even when using a mean-field variational family.
Pareto Smoothed Importance Sampling. Since importance sampling estimates have very high variance, Pareto smoothed importance sampling (PSIS) can be used to substantially reduce the variance with small additional bias (27). This procedure modifies and stabilises extreme importance ratios using a generalized Pareto distribution fit to the upper tail of the distribution of the ratios.
Variational families. vector À E R for D-dimensional inputs E RR. Typical choices of include mean-field Gaussian and Student's t families 33 [4], full and low rank Gaussians [15 22 mixtures of exponential families T7 19, and normalising flows 25 We focus on the most popular mean-fice field and normalizing flow families. Mean-field families assume independence across the D dimensions: g(0) II qi (O.). where each 9: typically belongs to some exponential family or other simple class of distributions. Normalising flows [I] provide more flexible families that can capture correlation and non- linear dependencies. A normalizing flow defined via the transformation of a probability density through sequence of invertible mappings. By composing several maps, simple distribution such as a mean-field Gaussian can transformed into a more complex distribution 25
Let gx (0) be approximating family parameterised by K-dimensional
f-divergences. The most commonly used divergences are examples of f-divergences 28 For a convex function f satisfying f(1) the f-divergence is given by
( q) Eang { (p(e) (0) (0)
The exclusive Kullback-Leibler (KL. divergence corresponds to s(w) = - log(ur). ) the inclusive KL divergence corresponds fw) w log( (w), the x² divergence corresponds to f(w)) = (w 1)?.
and the general a-divergences correspond to (w - w)/ (a(a - 1)}. We also consider the adaptive f-divergence proposed by Wang et al. (29]
Loss estimation and stochastic optimization. In all the cases we consider, minimizing the f- divergence is equivalent to minimizing the loss function L || q) ==
(although, see Wan al. (28) for a different approach). Let L(A) == ( qx) denote the loss as function of the variational parameters X. The loss and its gradient G()) VxL(A) can both be approximated using, respectively, the Monte Carlo estimates
Ī(A)= 1 E. f(w(O.)) and G(A) gg,))(,),
(2)
3

---

where 01.....89 gg are independent draws from gA and RK + R is an appropriate gradient-lis, function that depends on f and w. The two most popular gradient estimators the literature are the score function and the reparameterization gradient estimator 220 30 The score function gradient corresponds to g(0) = (f)))))
estimator ((()) V log (0). It is a general-purpon
that applies to both discrete and continuous distributions q. but it is known to suffer from high variance. When this estimator is used for the mass-covering divergences such as the inclusive KL and general a-divergences with o > 1, the importance weights are usually replaced with self- normalized importance weights wo wo(0))/ Ei= w(()).. The reparameterization gradient [20 requires expressing the distribution QA as a deterministic transformation of a simpler base distribution r such that (x) ~ 9A with ~ r. This allows writing an expectation with respect to y as an expectation over the simpler distribution r. The reparameterization estimator corresponds to using g(zx) = V (z+)) (for 2, r) in place of g(0), where w implicitly depends on À as well. In the case of the adaptive f-divergence, the importance weights w(02),.... w(0s) are sorted, and the gradients corresponding to each sample are then weighed by the empirical rank. The gradient estimates can be used in a stochastic gradient optimization scheme such that
X+1 < X' + n,(()')
(3)
where 7/ is the step size. In practice, more stable adaptive stochastic gradient optimisation methods such RMSProp Adam p 13), which smooth or normalize the noisy gradients, are often used.
Numerous prior work have studied some of the challenges tied to optimizing these divergence measures under the presence of noisy gradient estimates I1 29 31. Particularly, when dealing with mass-covering divergences, the gradient estimates can become so noisy that convergence not possible in practice, as we will illustrate later on.
3
A Asssssig the Reliability of Black-box Variational Inference
3.1
Conceptual framework
How can we determine - both conceptually and experimentally what is required to obtain reliable estimates of the variational divergence and optimal variational approximation? As we have seen, the most common variational divergences and their Monte Carlo gradient estimators can be expressed in terms of the density ratio w()). Reliable black-box variational inference ultimately depends on the behavior of (0) since (1) accurate optimization requires low-variance and (nearly) unbiased gradient estimates (A), and (2) determining convergence and validating the quality of variational approxima- tions can require accurate estimates T.(2) of variational divergences [14 [5]. While asymptotically (in the number of iterations and Monte Carlo sample size S) there may be no issues with stochastic optimization or divergence estimation, in practice black-box variational inference operates in the pre-asympootic regime. Therefore, the reliability of black-box variational inference depends on the pre-asymptotie behavior of the (0), and how it interacts with the choice of variational objective and gradient estimator.
Before accounting for the effects of the objective and gradient estimator, first consider the behavior of the density ratio wf), which can also be interpreted as an importance sampling weight with (0) as the proposal distribution [cf. 16. 29, Pickands 23 proved, under commonly satisfied conditions, that for u tending to infinity, the distribution of w()) w(()) > u is well-approximated by the paremmerreettr generalized Pareto distribution GPD(u, a, k). which for l: > 0 has density p(w u, o, k) = o [11 + k(w - w)/a ///////-111/a--1
111/ where w is restricted to (u, oo). Since w0 > O, this implies its distribution is heavily skewed to the right with a power-law tail. Consider the idealized scenario of estimating the mean of w (0) ~ GPD(u, (u, ar k). We assume the mean is finite, which is equivalent to assuming A < 1 since [1/k] determines the number of finite moments. Because of the heavy right skew, most of the mass of w(f) is below its mean. Therefore, even after averaging a large number of samples, most empirical estimates Σ w(O.) will be smaller than the true mean. Figure Ba illustrates this behavior for different values k: even with million samples, the empirical mean is far below the true mean when h > 0.7. The highly variable sizes of the confidence intervals based on 10,000 replications further highlight the instability of the estimator. So, even though the empirical mean is an unbiased estimator, in the pre-asymptotic regime (before the generalized central limit theorem applicable (5),, in practice the estimates are heavily biased downward with high probability. If (0) is not generalized Pareto distribution, we can instead treat & as
-1

---

1.2
mean) /True mean
k=0.9
0.8
Estimated Estimatedmean)
k-0.8 -5.95
0.6
0.2
k=0 .9
10
10
10
10'
Sample size
Figure 2: The ratio of estimated mean and true mean for different values of k shape parameter of a generalized Pareto distribution and confidence intervals in finite sample size simulation.
the tail index k := inffe 0 (w < oo), which encodes the same tail behavior as GPD(u, o, k). Crucially, we should expect to much larger than 0 when there is a significant mismatch between the target distribution and the variational family. Since selecting variational family that can match the typical set tends to more difficult in higher dimensions, we should expect & to be larger for higher-dimensional posteriors.
We can generalize our observations about pre-asymptotic estimation bias to the estimators Σ(A) and (()). For the loss estimator, replace w(.)) with f(w(0.)), where f f(w) is polynomial in w and log w for the class of losses we consider. If the dominant term of f(w) is of order w", the tail behavior will be similar to generalized Pareto with Ka ak. Thus, L(A) will have larger pre-as aeeetmtmptoicc bias as a increases. For example, estimation of the mass-covering inclusive KL (where 1) and, more generally, mass-covering a-divergences with >> 0 will suffer from large pre-asymptotic bias. On the other hand, for the mode-seeking exclusive KI, f(w) = log(w), can expect all moments to be finite and therefore much smaller pre-asympootic bias.
Similar considerations apply the gradient estimator, with the details depending on the specific estimator used. However, when using self-normalized weights for a-divergences, can expect large pre-asymptotic bias whenever w0) has such bias since self-normalization involves estimating the mean of w0)) This bias will affect the accuracy of the solution found using stochastic optimization. Thus, the quality of the solutions found can only partially be improved by using smaller step size since smaller step sizes will only reduce the effects of large estimator variance, but not the effeets from a large bias. We provide more details on the behavior of the score function and reparameterized gradients for each of the divergences in ???? following Geffner and Domke
In summary, our framework makes two key predictions:
(P1) Estimates and gradients of mode-seeking divergences (in particular exclusive KI. divergence with log dependence wr) have lower variance and are less biased than those of mass-covering divergences (in particular a-divergences with a > O, with polynomial dependence on w).
(P2) The degree of polynomial dependence on w determines how rapidly the bias and variance will increase as approximation accuracy degrades in particular, high dimensions.
Because the adaptive f d-divergence depends directly on the (ordered) weights, we expect to behave similarly to the mass-covering divergences.
3.2 Experimental framework
In the light of potentially large non-asymptotic bias arising from the heavy right tail of w()), it is important verify the pre-asymptotic behavior of the Monte Carlo estimators used in variational inference. We follow the approach developed by Vehtari et al. [77 for importance sampling and compute empirical estimate & of the tail index & by fitting generalized Pareto distribution to the observed tail draws. In the importance sampling setting, Vehtari et 27 show that the minimal sample size to have small error with high probability scales as S ((exp k)2)). Vehtari
5

---

Exclushe KL
Chr2
Exclshe KL
Chrz
Excls KL
C122
incluske KL
Indushv KL
Inclusiv
50
3.0
25
2.0
V DAt, q)
Pareto k
0.7
-1.0
10 Dimension
20
50
10
20
40
50
Dimension
Dimension
(h)
Figure 3: Results for correlated Gaussian targets of dimension D ..... 50 using either the exclusive or inclusive KI. divergence as the variational objective. Each line in the plots corresponds fitting and evaluating the same divergence measure indicated in the legend. Each result the average of 50 independent simulations. Quantiles are computed from simulating 100, 000 draws. (a) Bias and variance of the gradients of the optimised f-divergence for one parameter 0a for increasing dimensions at the end of the optimisation for correlated Gaussian targets of dimension D) = 1,..., 20 and mean field Gaussian as variational approximation. (b) The ratio of the f-divergence estimate to the true value. (c) The & values for the variational approximations.
al. 27 also demonstrate that i provides practical pre-asymptotic convergence rale estimate even when the variance is infinite and generalized central limit theorem holds. While estimating i in general requires larger sample size than is commonly used to estimate the stochastic gradients, we can still use it to diagnose and identify the challenges with different divergences. If k > 0.7, the minimal sample size to obtain reliable Monte Carlo estimate is so large that is usually infeasible in practice. This cutoff is in agreement with our findings shown in Fig. 3a Thus, together with our conceptual framework, we have a third key prediction:
(P3) The è value can be used diagnose pre-asymptotic reliability of variational objectives. In particular, the a-divergence with a > will become unreliable when max(1, 0.7, even if is bounded (by very large constant).
3.3 Verification of Pre-asymptotic (Un)reliability
compute most of the relevant quantities such as the loss function in closed form. Specifically, we fit a mean-field Gaus- sian to a Gaussian with constant 0.5 correlation factor using the inclusive KL, exclusive KL, x2. and 1 1/2-divergences. We vary the dimensionality D from to 50, which a surrogate for the degree of mismatch between the optimal variational approximation and the target distribution. To find the optimal divergence-bar based approximation, we optimize the closed-for foosed-form expression for the divergences between two Gaussians. Hence, the complexities and uncertainty due to the stochastic optimization. Due to space limitations, we focus representative cases of the approximations from optimising the mode-seeking exclusive KI. divergence and the mass-covering inclusive KI. divergence. Results for the other divergences are included in the appendix.
We first verify our three key predictions in simple setting where
consider the best-case scenario and ignore
(P1) Mode-seeking divergences are more stable and reliable than mass-covering ones. Figure 3b shows that as the approximation- target mismatch increases with dimension, the bias in and variance of the divergence estimates increases substantially for the inclusive KL and x² but only moderately for the exclusive KL. Similarly, Fig, 2 shows that gradient bias and variance increases with dimension for inclusive KL and x2 but not exclusive KL.
(P2) Degree of polynomial dependence on w determines sensitivity to approximation-taagel mismatch. Figure 3b shows that divergence estimates resulting from optimising higher polynomials of w become more and more unstable as dimensions increases.

---

Gaussian
Planar Flow NVP Flow
Gaussian MF Student-t
Planar Flow NVP Flow
MF Stueen-t
Exclus ive KL
Exclusive KL
10
2-
10a
0.7
10a
10a
FDiv
Diiv
00
1015
Pareto k
2-
10*
(
0.7
10a
0-
10"4
Inclusive KL
Inclusive KL
1045
2-
10a
0.7
10a
10a
0 10
20
30
40
50
10
20
30
40
50
D Dimension
Dimension
(a)
(b)
Figure 4: Results for increasing dimensions of the robust regression model. (a) Pareto i values for BBVI approximations. (b) Relative error of covariance estimates for BBVI (solid lines) and after PSIS correction (dashed lines).
(P3) i diagnoses pre-asymptotic reliability. Figure 3e shows that the h values grow rapidly for the inclusive KL-based approximation, particularly for higher-degree dependence on w, which agrees with predicted behavior and large bias and variance of the inclusive KL and x². In contrast, the k values remain fairly stable for the exclusive KL-based approximation, again in agreement with predicted and observed bias and variance behavior.
Experiments
In this section, describe series of experiments study how our pre-asymptotic framework be used for assessing the reliability of black-box variational approximations for practical applications and developing best-practices. For all posterions,
mean-field Gaussian and Student-t families, planar flow (25) with 6 layers and non-volume preserving (NVP) flow [S] with 6 stacked neural networks with 2 hidden layers of 10 neurons each for both the translation and scaling operations with a standard Gaussian distribution for the latent variables. We use Stan (26 for model construction. For stochastic optimization we use RMSProp with initial step size of 10 run for either T Taxx iterations or until convergence was detected using a modified version of the algorithm by Dhaka et al. [6]. For the exclusive KL we use 10 draws for gradient estimation per iteration, while for the other divergences we use 200 draws, and a warm start at the solution of the exclusive KI. In practice, we found the optimisation for x2 divergence extremely challenging, with the solution failing to converge even for moderate dimensions D R 10. Therefore, we only include results for the KI. divergences and the adaptive f-divergence. We compare the accuracy of approximated posterior moments to ground-truth computed either analytically or using the dynamic Hamiltonian Monte Carlo algorithm in Stan (26). Specifically, consider the estimates jt and E for, respectively, the posterior mean μ and covariance matrix Σ. We also consider the mean and covariance estimates produced by PSIS and compute k. The experiments were carried on a laptop and an internal cluster with only CPU capability. The code for the experiments will be made available after acceptance using MIT license.
4.1 Haived posteriors
First, we study the toy robust. regression model previously used by Huggins al. 14 given by
3g N(O, 10),
Un In, 8 ~ (1
1),

---

NF Coussian Excl, KL
NVP Flow Excl, KL
Planar Flow Excl, KL
NF Coussion Incl. KL
NVP Flow Incl KL
Planar Flow Incl. KL
Max dmercionsiity
30
20
30
10x
23x10
$2x10
12x1e
27x10
61x10
14x10*
32v10
72x10
Láx1c
37x10
#3x1e
19x20
44x10
19x10"
Step size
Figure 5: Maximum dimensionality converged per step size for the robust regression model.
where yn E Zn E RDD are the target and predictors respectively, 3 denotes the unknown coefficients, and D is varied from 2 to 50. We generated data from the same model with covariates generated from a zero-mean Gaussian with constant correlation of 0.4. The Student's t leads to the posterior having heavy tails, making it more challenging target distribution. We use Tmax 10,000.
Mode-seeking divergences are easier to optimize. Figure 4a shows that the estimated tail index generally increases with the dimension as expected. In particular, the & values when using normalizing flows, which are more challenging to optimize. is low for D < 20 when using exclusive KL, but infinite when using either the inclusive KL or f-divergence. From Fig. 4b we can that exclusive KI. provides also more accurate and reliable posterior approximations than the inclusive KI. and adaptive /-divergence, particularly for the normalizing flows. This observation is consistent with the prediction (P3) of the proposed framework. The better performance for normalizing flows corroborates the relative ease of stochastic optimization with the exclusive KI. divergence compared to the inclusive KL or the adaptive f-divergence - despite the fact that we used 20 times as many Monte Carlo samples estimate the gradients for the inclusive KL and the f-divergence compared to the exclusive KL. To further illustrate the relative difficulty of optimizing the inclusive KL divergence, Fig. 5 shows the largest dimension for which the stochastic optimization converged as a function of the step-size. For most step-sizes, the combination of normalizing flows and the inclusive KL divergence only converged for D 2. whereas convergence is possible in higher dimensions for simpler variational families. These observations are consistent with predictions (PI)-(P2) of the proposed framework.
Adaptive f-divergence interpolates between the exclusive and inclusive KI. divergence, but is difficult optimize. In low dimensions, the adaptive f-divergence behaves somewhere between the two KL divergences as seen in ?? as it was designed to [29]. As confirmed by Fig. 4 For higher-dimensional posteriors, we expect
behave more like the exclusive KL, but it less stable due to its functional dependence on the importance weights.
Normalizing flows can be effective but are challenging to optimize. Fig. also shows that nor- malizing flows can be quite effective when used with exclusive KL to ensure stable optimization. However,
be seen in ??, when using out-of-the-box optimization with no problem-specific wning (as we have done for a fair comparison), the normalizing flows approximations can have pathological features even in low dimensions.
4.2 Realistic models and datasets
We now study how the choice of divergence and approximating family compare across diverse range of benchmark posteriors. We compare variational approximations for models and datasets from posteriorde ordb" in terms of accuracy of the estimated moments and predictive likelihood. We used an 80/20 training/test split on all datasets to compute the predictive likelihoods. We use T Tiax = 15,000.
8

---

MF Gaussian
Planar Flow NVP F Flow
MF Gaussian
Planar Flow NVP Flow
MF Student-t
MF S Student-t
Exclusive KL
10
2
10
0.)
10
0
10a
Diiv
Pareto
2
10a
0.
0-
10a
Inclus lve KL
10a
2-
109
0.7
dogs (5)
ark (7)
mesquitte (8)
eschools ncp (10)
8schools ap (10)
nes1996 (12)
diamonds (26)
radon, unp (90)
dogs (5)
ark (7)
mes quite (8)
Bschools, ncp (10)
8schools cp (10)
nes 1996 (11)
diamonds (26)
radon unp (90)
(b)
Figure 6: Results for posteriordb experiments. Dimensionality of each dataset is given in paren-
theses. (a) Pareto values for BBVI approximations, (b) Relative error of mean and covariance estimates for BBVI using exclusive KL (circles) and after PSIS correction (triangles).
Table 1: Predictive likelihood results on posteriordb datasets using Mean Field Gaussian approximation. The results denote the likelihood with the variational approximation solution obtained and after PSIS correction to the solution. Bold (underline) indicates best- best-performing method(s) (variational method(s))
Name dogs arK mesquite
HMC
Excl. KL.
Exel. KL+PSIS
Incl. KL.
Incl. KL+PSIS -70.5.14.1
-71.111.2
-71.2 113
-71.7E1..
-110+3.55
-32.4:L0.6
-34.3.10.7
.3..4
-35.2.00.8
-34.910.8
-1681=127
-2512=140
18618+186
-00
nes 19s1996
-112.911.7
-112.
-27.9911.8
2140.5559.7
499.3145.6
diamonds
22.113.11
-2.611..
1.5+12
-3196.6157.9
31.99455.7
radon
-234.4+1.1
13.3
-325.0+19.5
-3774+1.9
-3.0,5
Exclusive KL remains the most reliable for realistic posteriors. The results are summarized in Fig. 6 where the same pattern is seen: the exclusive KL is superior for higher-dimensional posteriors (e.g., > 10) or when combined with normalizing flows, while inclusive KL is better for lower-dimensional posteriors. Despite the superior performance of the exclusive KL divergence, the large values for indicate that fitting approximations based normalizing flows remains challenge in high dimensions. The performance for the adaptive f-divergence is comparable to the inclusive KL divergence. Table D shows that the exclusive KI. divergence consistently outperforms the inclusive KI. divergence terms of predictive accuracy, but can be significantly worse than HMC.
Importance sampling shows the relative errors of the first two moments for the variational approximation (dots) and after correcting the estimates using PSIS (triangles). In some cases, the PSIS correction dramatically improved the accuracy of the normalizing flows.
substantially improve accuracy. Focusing exclusive KL, Fig. 6b
Reparameterizatiion is an important tool for improving accuracy. The 8-schools model is low- dimensional (D 10), but the funnel-shaped posterior makes inference challenging for variational approximations [4 31 As has been noted previously in the literature, and is clear from Figs. 6a 6h) reparameterizing the model so that the posterior better matches the variational family can be effective way improve the accuracy of the approximation. See ?? for illustration.

---

Discussion
Our conceptual framework based on the pre-asymptotie behavior of the density ratios importance weights w along with our comprehensive experiments lead number of important takeaways for practitioners looking to obtain reasonably accurate posterior approximations using black-box variational inference:
The instability of cass-covering divergences like inclusive KL and x² means that, given currently available methodology, users are better off using the exclusive KL divergence except for easy low-dimensional posteriors. The reliance of the adaptive f f-divergence on importance weights leads to similar instability.
Importance sampling appears to almost always be beneficial for improving accuracy, even when the diagnostic is large. However, a large k does suggest the user should not expect even the PSIS-corrected estimates to be partieularly accurate.
Using normalizing flows particularly NVP flows together with exclusive KL. and PSIS provides the best and most consistent performance across posteriors of varying dimensionality and difficulty, We therefore suggest this combination as a good default choice.
Our results suggest an important direction for future work is improving the stability of optimization with normalizing flows, which still tend to have some pathological behaviors unless they are very carefully tuned since such tuning significantly detracts from the benefits of using BBVI.
Limitations
While our experiments included range common statistical model types, our findings may not generalize to all types of posteriors or to other variational families. For example, we did not explore semi-implicit methods or applications neural networks. We also did not investigate alternative divergences such as those used in weeighte autoencoders.
Acknowledgments
We thank Sihan Liu for pointing out and fixing errors in the code. We thank Ari Hartikainen and PyStan developers for help with posteriordb, Academy of Finland (grants 298742, and 313122) and Finnish Center for Artificial Intelligence for partial support of the research. We also acknowledge the computational resources provided by the Aalto Science-IT project. AKD would like to thank Tomas Geffner and Eero Siivola for useful discussions.
References
[1] Abhinav Agrawal, Daniel R. Sheldon, and Justin Domke. Advances in black-box VI: nor- malizing flows, importance weighting, and optimization. In Advances in Neural Information Processing Systems 33: Annual Conference Neural Information Processing Systems 2020, NeurlPS 2020, December 6-12, 2020, virtual, 2020.
[2] Robert Bamler, Cheng Zhang, Manfred Opper, and Stephan Mandt. Perturbative black box variational inference. In Advances in Neural Information Processing Systems, volume 30, pages 5079-5088, 2017.
[3] D. M. Blei, Alp Kucukelbir, and Jon D McAuliffe. Variational Inference: A Review for Statisticians. Journal of the American Statistical Association, 1185118)599-87,, 2017.
[4] Jörg Bornschein and Yoshua Bengio. Reweighted wake-sleep. In 3rd International Conference on Learning Representations, ICLR 2015, San Diego CA, USA, May 7-9, 2015, Conference Track Proceedings, 2015.
[5] Louis Chen and Qi-Man Shao. Normal approximation under local dependence, The Annals of Probability, 12(3):198522022 2004.
10

---

[6] Akash Kumar Dhaka, Alejandro Catalina, Michael R Andersen, Mäns Magnusson, Jonathan Huggins, and Aki Vehtari. Robust, accurate stochastic optimization for variational inference. In Advances in Neural Information Processing Systems, volume 33, pages 10961-10973, 2020.
[7] Adji Bousso Dieng, Dustin Tran, Rajesh Ranganath, John Paisley, and David Blei. Variational inference via kchi upper bound minimization. In Advances in Neural Information Processing 2017.
Systems 30, pages
[8] Laurent Dinh, Jascha Sohl-Dickstein, and Samy Bengio. Density estimation using real nvp. In International Conference on Learning Representations, 2017.
[9] John Duchi, Elad Hazan, and Yoram Singer. Adaptive Subgradient Methods for Online Learning and Stochastic Optimization. Journal of Machine Learning Research, 12:2121-2159, 2011.
[10] Tomas Geffner and Justin Domke. On the difficulty of unbiased alpha divergence minimization, 2020.
[11] Tomas Geffner and Justin Domke. Empirical evaluation of biased methods for alpha divergence minimization. In Symposium on Advances in Approximate Bayesian Inference, AABI 2020, 2020.
[12] Jose Hernandez-Lobatte Yingzhen Li, Mark Rowland, Thang Bui, Daniel and Richard Turner. Black-box alpha divergence minimization. In Proceedings The 33rd International Conference on Machine Learning, volume 48, pages 1511-1521 PMLR, 2016.
[13] G. E. Hinton and Tijmen Tieleman. Lecture 6,5 - Rmsprop: Divide the gradient by running average of its recent magnitude. In Coursera: Neural networks for machine learning. 2012.
[14] Jonathan H Huggins, Mikolaj Kasprzak, Trevor Campbell, and T. Broderick. Validated Varia- tional Inference via Practical Posterior Error Bounds. In AISTATS, 2020.
[15] Alp Kucukelbir, Rajesh Ranganath, Andrew Gelman, and D. M. Blei. Automatic Variational Inference Stan. Advances in Neural Information Processing Systems, Advances in Neural Information Processing Systems, 2015.
[16] Yingzhen Li and Richard E Turner. Rényi divergence variational inference. In Advances in Neural Information Processing Systems, volume 29, pages 1073-1081, 2016.
[17] Wu Lin, Mohammad Emtiyaz Khan, and Mark Schmidt. Fast and Simple Natural-Gradient Variational Inference with Mixture of Exponential-fanne family Approximations. In International Conference Machine Learning, 2019.
[18] David C. MacKay. Information Theory, Inference and Learning Algorithms. Cambridge University Press, 2003.
[19] Andrew C. Miller, Nicholas J. Foti, and Ryan P. Adams. Variational boosting: Iteratively refining posterior approximations. In Proceedings of the 34th International Conference on Machine Learning, volume 70 of Proceedings of Machine Learning Research, pages 2420-2429. PMLR, 2017.
[20] Shakir Mohamed, Mihaela Rosca, Michael Figurnov, and Andriy Mnih. Monte Carlo gradient estimation in machine learning. arXiv preprint arXiv: arXiv: 1006. 2019.
[21] Christian A. Naesseth, Fredrik Lindsten, and David M. Blei. Markovian score climbing: Variational inference with KL(p/lq). CoRR, abs/2003. 10374, 2020.
[22] Victor M H Ong, David Nott, and Michael S Smith. Gaussian Variational Approximation With a Factor Covariance Structure. Journal of Computational and Graphical Statistics, 27(3):
465 478, 2018. doi: 10.1080/1061 18000.001..399044
[23] James Pickands. Statistical inference using extreme order statisties. Annals of Statistics, 119-131, 1975.
[24] Dennis Prangle. Distilling importance sampling, 2021.

---

[25] Danilo Rezende and Shakir Mohamed. Variational inference with normalizing flows. In Pro- ceedings of the 32nd International Conference on Machine Learning, volume 37 of Proceedings of Machine Learning Research, pages 1530-1538. PMLR, 2015.
[26] Stan Development Team. Stan modeling language users guide and reference manual. 2.26, 2020. httph:t/t//mmommm./ stan. org
[27] Aki Vehtari, Daniel Simpson, Andrew Gelman, Yao Yuling, and Jonah Gabry. Pareto smoothed importance sampling. arXiv preprint 150i 02640,, 2019.
[28] Neng Wan, Dapeng Li, and Naira Hovakimyan. I-Divergence Variational Inference. In Advances in Neural Information Processing Systems, 2020.
[29] Dilin Wang, Hao Liu, and Qiang Liu. Variational inference with tail-adaptive f-divergence. In Advances Neural Information Processing Systems, volume 31, pages 5737-5747, 2018.
[30] Ming Xu, Matias Quiroz, Robert Kohn, and Scott A. Sisson. Variance reduction properties of the reparameterization trick. In Proceedings of Machine Learning Research, volume 89 of Proceedings of Machine Learning Research, pages 271 1-2720. PMLR, 2019.
[31] Yuling Yao, Aki Vehtari, Daniel Simpson, and Andrew Gelman. Yes, but did it work?: Evalu- ating variational inference. In Proceedings of the 35th International Conference Machine Learning, volume 80 of Proceedings of Machine Learning Research, pages 5581-5590, PMLR, 2018.
12


Reviews:
[{"review_id": "mVvHe6sFILu", "reviewer": "Reviewer_EXRc", "summary": "The paper discusses the accuracy of variational inference and the choices that need to be me made with regard to divergence and approximating family, particularly as they relate to the dimensionality of the problem. Based on an analysis of the pre-asymptotic behavior of the density ratios involved, the authors propose a framework based on Pareto smoothed importance sampling (PSIS) to help practitioners analyze the success of their approximate inference procedure. The authors provide experiments demonstrating general findings for variational inference for a range of divergences and approximate posteriors. ", "questions": "", "limitations": "", "rating": 7, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Overall, I found the paper to be quite clear and well-written. I believe the paper to be quite significant as I am not aware of much work on trying to provide a framework for understanding variational inference and the choices that need to be made, particularly in a way that could be used by a practitioner with relative ease. I believe the insights on dimensionality, mode-seeking vs mode-covering, and the differences in the quality of approximation from different approximate posteriors (e.g. MF Gaussian, Student's t, NVP), are novel, valuable and provide good guidance for practitioners and those who regularly use VI in their research.\n\nThe main issue I had with the paper was the lack of a proper review of prior work. The paper does reference some prior work in passing, but it does not really discuss the prior work (for instance [14, 30]) in the detail that I think would be required. For instance, [30] quite similarly proposes using the $\\hat{k}$ obtained from PSIS to evaluate how successful VI is - although they did not analyze different divergences from the standard mode-seeking one nor the effect of dimensionality. I think the paper would really require a proper discussion of these papers, as well as other attempts to assess the quality of VI approximations.\n\nOverall, as someone who focuses on VI in their research, I think this is quite a good paper and represents a useful and needed direction for the field. Therefore I believe the paper deserves acceptance, but I would of course be curious to hear what the other reviewers say.\n\nBelow are some more minor/tangential comments/questions:\n- The authors frequently make reference to the fact that they are interested in the \"pre-asymptotic\" regime. To me it would seem intuitive that this would depend quite significantly on the initialization of the approximating distribution, particularly as the divergences and dimensionality are changed. Do the authors have any intuition as to what the best way of initializing the approximate posterior are for these different cases?\n- For more complex models, it is often the case that there are symmetries that result in many modes that are identical when it comes to predictions, for instance permuting neurons in the same layer of a neural network would result in such modes. In the case of a unimodal approximate posterior, it would in principle therefore be possible to model one of these modes while missing the others in a way that the predictions are actually reflective of the true posterior. It would seem to me that in such a case the $\\hat{k}$ would be quite bad even though the predictions are fine - do the authors have any ideas on how this could be addressed?\n- It seems as though the labels are incorrect for Figs. 2 and 3. For instance, in line 138 Fig. 3a is referred to when I believe Fig. 2 is the correct one. \n- In Table 1, which approximate posterior is being used?\n- l. 259 - \"but it less\"", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "The authors addressed some limitations of their work in that the discussion is limited to a typical f-divergences used in the literature and does not discuss all types of posteriors seen in the literature, such as semi-implicit posteriors. While not an extensive discussion, I found this description of limitations to be adequate for the goals of the work. I am not aware of potential negative societal impacts for this work.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "8"}, "extra_text": {}}, {"review_id": "VEJQlKNRbGD", "reviewer": "Reviewer_c7n8", "summary": "This paper provides a statistical tool to analyze different variational inference methods in the f-divergence family. It showed that for high dimensional and low dimensional problems, the behavior of the variational inference methods may differ, and the analysis based on the GPD statistics can suggestion different strategies for different situations. The predictions from this theory were validated by some commonly used models (e.g., linear regression). In summary, this paper can provide useful guidelines for variational inference practitioners.", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The main novelty of this paper is to use GPD to model the density ratio, which reveals important properties of different variational inference methods in the f-divergence family, which includes the inclusive KL and exclusive KL approaches. This provides a practical tool for people to use for variational inference diagnosis. One particular conclusion from this analysis is that, for moderate-to-high-dimensional posteriors, exclusive KL divergence is a better choice than inclusive KL divergence. These conclusion can be very helpful for people to choose the right variational inference method in different scenarios.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "Here are some more detailed commons of the paper:\n\n1. The three predictions (P1, P2, P3) are reasonable and are verified by some experiments. However, more rigorous proofs or justifications are needed to generalize the conclusions to other applications scenarios. For example, the authors mentioned that application to neural networks is not covered theoretically or empirically. \n\n2. It would be better to provide a table to summarize the conclusions in the paper, for low dimensional and high dimensional problems. Properties and behaviors of different variational inference method could be included in the table, together with the recommended inference method in the f-divergence family.\n\n3. A minor technical issue, in line 103, for the normalized importance weight, there should be a 1/S term in the denominator.\n\n4. A problem in notation: D is used for divergence, but also used for dimension.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "4 hours"}, "extra_text": {}}, {"review_id": "Ms_2rnwN8F", "reviewer": "Reviewer_8K18", "summary": "The authors propose a conceptual framework and design guidelines for variational inference based on the pre-asymptotic tail behaviour of the unnormalized importance weight. This framework is used to study the behaviour of commonly variational families and f-divergences in mid-to-high dimensions. The authors observe that mass-covering divergences, while often superior theoretically and for low dimensional problems, are hard to optimize in mid-to-high dimensions compared to the mode-seeking exclusive KL-divergence. These observations are consistent with the prediction made using their framework based on estimating the tail-index of the best-fit generalized pareto distribution.", "questions": "", "limitations": "", "rating": 7, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "**Originality and Significance:**\nAs far as I can tell, the proposed framework does not contain novel contributions on a technical level, i.e. the work builds heavily on previous work on pareto smoothed importance sampling (PSIS) [1] and PSIS based diagnostics [2] (as adequately cited by the authors), and the authors use standard estimators for the different f-divergences and their gradients.\nThat said, I believe there is value, especially for new partitioners, in presenting the above in a coherent and unified framework. The presented framework makes it easy to study the behaviour of commonly used variational families and divergences and helps making educated predictions based on the order of the (estimated) tail index and f-divergence. Hence, I do believe the setup and presentation are original and valuable for the community.\n\n**Clarity:**\nThe writing of the paper is clear and easy to follow.\nOne thing that I noticed is that the authors do not explain how to actually use the PSIS correction to improve the estimator. Adding a paragraph explaining this would make the paper more self-contained and enable readers not familiar with [1, 2] to better understand the experiments.\n\n**Quality:**\nOverall, the methodology is technically sound and the observations and conclusions drawn by the authors seems correct but sometimes slightly to general. \n\nSpecifically, the statements P1, P2 made in section 3.3 are too general to be supported by the actual experiments.\n- While *P1* mention 'mode-seeking'-divergences and 'mass-covering'-divergences the experiments only evaluate the exclusive KL-divergence as an example for 'mode-seeking' and the inclusive KL-divergence and $\\chi^2$-divergence as examples for 'mass-covering' divergences. The experiments are still insightful, but the claim should be less general to be supported by the actual experiments.\nThe main text mentions additional results for the other divergences in the appendix, however, the appendix does only include additional results for the evaluation presented in Figure 3c, but not for for the evaluations presented in Figure 3a and 3b.\n- The claim made in *P2* would be better supported if the plot in the main text would include the results for the $1/2$-divergence (mentioned on Line 194) and the importance weight $w$ itself. The corresponding plot for Figure 3c can be found in the appendix (Figure B.1.). \n   - Is there a reason why these results are not reported? \n   - For each of the line plots in Figure 3, what was the variational objective used for training? The figures in the appendix are explicit about what variational objective was used. It would be helpful to mention this in the figure caption in the main text as well. \n  - For clarity it would also be helpful to include the degree of the variational objective in the legend similar to the figure in the appendix.\n\nIt is also unclear to me if the results are averaged over multiple independently trained models. Given that this paper consideres non-convex stochastic optimization, I consider this necessary to reach any conclusions. Specifically,\n- Figure 4, Figure 6 and Table 1 do not report any metric of variability and I can not find any statement in the text mentioning if the presented values are averages over multiple runs.\n- How are the confidence intervals computed in Figure 2?\n- How is the reported variance computed in Figure 3?\n\nSome clarification would be helpful here.\n\n\n**Additional Questions and Comments:**\n- Regarding the intuition gained by observing the distance from the mode (Figure 1): The overlap between the distribution over distances from the mode is already less for the inclusive KL-divergence in 2D and this difference seems to carries over to the higher dimensional examples.  Even the 2D illustration does not give the impression that the approximation based on the inclusive KL would do better w.r.t. this metric (overlap between distributions over distances). For other metrics, i.e. the variance of the resulting importance weight, I'd agree that the approximation based on the inclusive KL-divergence seems preferable. Overall, I'm not sure how Figure 1 demonstrates that 'the benefits of heavy-tailed approximate families and divergences favoring mass-covering behaviour diminishes as the dimensionality of the target distribution increases'.\nAm I missing something here?\n- Figure 1; caption: The approximation is red for the exclusive KL and green for the inclusive KL but the caption only mentions that the approximation is shown in red.\n- Line 38, 29; some words seem to be missing in that sentence. 'the benefits of heavy-tailed approximate families and divergences favoring mass-covering *behavior* diminishes as *the* dimensionality of the target distribution increases.\n- Line 205: I think this should reference Figure 3a instead of Figure 2\n\n\n*[1] Aki Vehtari, Daniel Simpson, Andrew Gelman, Yao Yuling, and Jonah Gabry. Pareto smoothed importance sampling. arXiv preprint arXiv:1507.02646, 2019.*\n\n*[2] Yuling Yao, Aki Vehtari, Daniel Simpson, and Andrew Gelman. Yes, but did it work?: Evalu- ating variational inference. In Proceedings of the 35th International Conference on Machine Learning, volume 80 of Proceedings of Machine Learning Research, pages 5581–5590. PMLR, 2018.*\n", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "adequate", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "5"}, "extra_text": {}}, {"review_id": "A3ld2tFGJo4", "reviewer": "Reviewer_x3v3", "summary": "In this paper, the authors study the pre-asymptotic behavior of the density ratios between the joint distribution and the variational approximation. The authors use generalized Pareto distribution to model this density ratio in the pre-asymptotic regime and make three major predictions for how BBVI methods should behave. They verify their predictions through extensive experimentation, and finally, provide actionable advice to BBVI practitioners (use flows with exclusive KL and PSIS as a good default.) ", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "I appreciate the authors for taking up this important direction in BBVI research. In the recent years, the number of components in BBVI methods have increased dramatically and it is important that we, as a community, analyze the existing methods. This paper takes some very concrete steps towards it.\n\nOverall, I like the paper. It is clearly written for most part. However, I have a few concerns and seek explanations.\n\n- Figure 6. I am not sure how to make sense of $\\hat k$ when using flows. Even with exclusive KL, the $\\hat k$ value reaches $\\infty$ for moderate dimensions (less than 100). Still, the relative errors seem to suggest they work better than the other methods. How reliable is $\\hat k$ when comparing different methods while keeping the variational objective same? Since $\\hat k$ value for flows seems to be $\\infty$ for moderately large models, the use of $\\hat k$ as a diagnosis tool is rather limited?  \n - Figure 4b. I don't see the Planar flow and NVP flow lines (Green lines) for Fdiv and Inclusive KL? Were they not calculated because the $\\hat k$ is $\\infty$? If so, I believe something is not consistent. There is Relative error estimates for NVP Flow even though the $\\hat k$ estimate is $\\infty$ for exclusive KL and D = 50. Maybe I missed the explanation somewhere; it will be great if the authors can point it out or provide one. \n - Table 1. I am not entirely sure the bold and the underline are supposed to be. Also, is the entire table only for NVP based method? I could not find, in text or in the table caption, the family that was used for the reported numbers. Are there any error bars for these numbers? \n - Figure 5. Should not there be at least on entry corresponding to the planar and NVP flows for higher dimensions (>30)? Did none of these step sizes converged for flows? If so, how did we get the numbers in Figure 4? \n - At L69, what is $\\hat I$? an unbiased or biased estimator? This is more of a notational consistency issue. ", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "Yes", "review_text": "", "extra_scores": {"ethical_concerns": "None", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "3-4 hours"}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper is concerned with some broad questions in variational inference, namely what variational families should be used (e.g. light-tailed, heavy-tailed, flows) what divergence should be optimized (e.g. inclusive KL, exclusive KL) and how performance could be diagnosed after inference is complete. On a strict technical level, this paper appears to contain few contributions. Nevertheless, reviewers were overall positive about the paper's attempt to establish and support some grand (albeit somewhat informal) themes. These are 1) that mode-seeking divergences are easier to optimize than mode-spanning divergences 2) that this difficulty can be understood by considering the polynomial dependence of the divergence on importance weights and 3) that one can fit a generalized Pareto distribution and use the k statistic to diagnose inference success.

One weakness identified by several reviewers agreed on was an inadequate discussion of prior work. The paper would be stronger with a better review of prior work related to all three of the themes, here, i.e. on the difficulty of optimizing different divergences, on using different variational families for VI and on inference diagnostics. There are some citations now, but with somewhat cursory discussions. In particular, the paper would be stronger if it was self-contained so that someone not familiar with the PSIS framework can follow it. The authors have agreed to expand their discussion of prior work.

Reviewers also had some specific comments about the experimental results (both what was done and the presentation). The authors have also been receptive to this feedback.
Decision: Accept (Poster)