VISA: Variational Inference with Sequential Sample-Average Approximations
Authors: Heiko Zimmermann
Christian A. Naesseth
Jan-Willem van de Meent
Venue: NeurIPS 2024 poster (2024)
Keywords: Variational Inference
Sample Average Approximations
Importance Sampling
Primary area: probabilistic_methods
TL;DR: Forward-KL Variational Inference with Sequential Sample-Average Approximations

Abstract: We present variational inference with sequential sample-average approximations (VISA), a method for approximate inference in computationally intensive models, such as those based on numerical simulations. VISA extends importance-weighted forward-KL variational inference by employing a sequence of sample-average approximations, which are considered valid inside a trust region. This makes it possible to reuse model evaluations across multiple gradient steps, thereby reducing computational cost. We perform experiments on high-dimensional Gaussians, Lotka-Volterra dynamics, and a Pickover attractor, which demonstrate that VISA can achieve comparable approximation accuracy to standard importance-weighted forward-KL variational inference with computational savings of a factor two or more for conservatively chosen learning rates.

VISA: Variational Inference with Sequential Sample-Average Approximations
Heiko Zimmermann h. wuw.merraannuua....mmmmmmm nl
Christian A. Naesseth C. a. nawsseth@uva.nn ml
Jan- Willem van de Meent w. vandemeent@uva nl
Amsterdam Machine Learning Lab University of Amsterdam Amsterdam, The Netherlands
Abstract
We present variational inference with sequential sample-average approximations (VISA), a method for approximate inference in computationally intensive models, such as those based on numerical simulations. VISA extends importance-weighted forward-KL. variational inference (IWFVI) by employing a sequence of sample- average approximations, which are considered valid inside trust region. This makes it possible to reuse model evaluations across multiple gradient steps, thereby reducing computational cost. We perform experiments on high-dimensional Gaus- sians, Lotka- Volterra dynamics, and a Pickover attractor. We demonstrate that VISA achieve comparable approximation accuracy standard importance- weighted forward-KI. variational inference while requirering significantly fewer samples for conservatively chosen learning rates.
Introduction
Gradient-based methods have become the workhorse for inference in simulation-based models. When a model defines fully differentiable density, methods based on Hamiltonian Monte Carlo (Duane et al., 1987; Neal al., 2011) and reparameterized variational inference (Kingma & Welling, 2013; Rezende et al., 2014) are often considered the gold standard for generating high quality samples from the posterior distribution. However, it is not always practical possible to use differentiable model. The implementation of the simulator may not support differentiation, or the model itself may not be differentiable, for example because it employs discrete random variables stochastic control flow. In such cases, inference falls back on methods based on score- core-funniion estimators (Glynn, 1990; Wingate & Weber, 2013; Ranganath et al., 2014) IWFVI, which derives from reweighted wake-sleep methods (Bornschein & Bengio, 2015; Le et al., 2018). These methods are often less efficient, as larger number of samples is required compensate for higher variance gradient estimates. This can be problematic when model is expensive to evaluate, for example because evaluation involves numerical simulations. Nonetheless, these methods remain the most viable option in a substantial number of use cases.
In this paper, present VISA, method that substantially improve the computational efficiency of variational inference for models that are non- non-differentiable and computationally intensive. VISA is designed for applications where evaluation of the variational approximation is cheap relative that of the model. Here we can save computation by reusing model evaluations across multiple updates of the variational posterior. To this end, we adapt IWFVI to employ series of sample-average approximations (SAAs) (Nemirovski et al., 2009; Kim et al., 2015), which use a fixed set of samples
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

that defines deterministic surrogate to the objective, rather than generating fresh set of samples at each gradient step.
SAA-based methods were recently studied in the context of reparameterized black-box variational inference (Giordano et al., 2024; Burroni et al., 2024) which optimizes the reverse Kl.-divergence. These methods fix samples from parameter-free distribution, which are transformed to samples from the approximate posterior using a differentiable map, whose parameters are optimized maximize the variational bound. VISA differs from these methods in that it optimizes a forward-KL divergence and does not require differentiable model. Concretely, VISA fixes samples from parameterized variational distribution, rather than samples from parameter-fie free distribution. Since the variational distribution will change during optimization, we construct new SAA whenever the optimization procedure leaves trust region, which we define in terms of the effective sample size. The result is a drop-in replacement for IWFVI that re-uses samples as much as possible, thereby saving computation.
We evaluate VISA in the context of three experiments. We first consider high-dimensional Gaussians, where the approximation error can be computed exactly. We then consider inference in Lotka- Volterra system and a Pickover Attractor, where numerical integration is performed as part of the forward simulation. Our results show that VISA with a conservative (i.e. smaller than needed) step size can converge in smaller number of model evaluations than IWFVI with more carefully tuned step size. These results come with the caveat that VISA is more susceptible bias than IWFVI, especially when used with a low effective sample size threshold. Our experiments show, that VISA converges substantially faster than IWFVI for conservartively chosen learning rates, while VISA performs par with IWFVI for more carefully tuned step sizes.
2
Background
We first briefly review variational inference (VI) with stochastic gradient descent (SGD) and SAAs, before we introduce VISA in Section 3. Readers familiar with these topics can safely skip ahead.
2.1
Variational Inference
VI approximates intractable target density with a tractable variational distribution by solving an optimization problem. The objective is typically to minimize a divergence measure D between the variational approximation 4s with parameters & E & and the target density
min {L(Ø) = D(gg. T)}- Φ
(1)
We assume that the target density is the posterior of probabilistic model (2) = p(z y) for which we are able to point-wise evaluate the joint density p(y. The two most common approaches to VI are to minimize the reverse forward KL divergence, for which objectives can be defined in terms of a lower bound -- Lg and upper bound LE on the log marginal log p(y)(
P(U, 2) au (2(2)
p(y, x)
L = - E log
2 - p(v),
Lx(Ø) = E log
> log p(v).
p(-1x)
(2)
We will briefly discuss standard reparameterized VI, which maximizes the lower bound Lng, and IWFVI, which minimizes the upper bound Lg.
Reparameterized VI. When maximizing a lower bound with stochastic gradient descent, we can either employ score-function estimators (Ranganath et al., 2014; Paisley et al., 2012), which tend to exhibit a high degree of variance, or make use of the reparameterization trick to compute pathwise gradient estimates (Kingma & Welling, 2013; Titsias & Lázaro-Gredilla 2014; Rezende et al., 2014).
Reparameterization generates samples from a pushforward density = = Ta(E) ~ 2ø by fist sampling E ~ a from a distribution independent of and then transforming the these samples using a differentiable function Tg. Under conditions that allow to exchange the order of differentiation and integration we can compute an unbiased estimate of the gradient
d P(y,To(E))
N
ply. TT((())) log
(o( =
log
g(i gg.
do
95 d (T(((E))
N
dø
4s (T (T)()))))
i=

---

Importantly, reparameterized VI requires model p(y, 2) that is differentiable with respect to & in order compute the pathwise derivative.
Forward-KL VI. Approximating the gradient of the forward KL- divergence does not require differentiability of the model p(y. 2) with respect to 2, but does require approximate inference to generate samples from the posterior p(z y). IWFVI uses self-normalized importance sampling to propose samples from the variational distribution and reweights them accord- ing to the posterior by introducing an importance weights we ((v,z)/ug(z)
Importance-Weigee
d dΦ
Lo)
BEE pl-bv) de
log c0(2)
Σ
( d dø
log (=)0
Σ/=
N
0)
where 20) ~ 4s- The resulting estimate is biased but consistent, meaning that it converges almost surely the true gradient as N DC.
2.2 VI with Sample-Average Approximations
Sample-average approximations (see Kim et al. (2015) for a review) approximate an expected loss with surrogate loss in the form of Monte Carlo estimate. In contrast standard VI, the samples that the SAA is based on remain fixed throughout the optimization process. This means that the surrogate objective can be treated like any other deterministic functiin, which can be optimized using standard optimization tools. Concretely, a sample-average approximation applies to an optimization problem of the form
[(0): : E(((),))
min $E
in which the density p(z) does not depend on the parameters Φ. This means that we can compute surrogate loss (()) that is unbiased estimate of the original loss L(Ø) by averaging samples from p.
N Σ 1(200,))
1 L(0) = N
20 ~ p.
Under mild conditions, the number of samples N - the minimizer I arg ming (()) and the minimal value L()) converge almost surely the minima 6" arg min (Ø) and minimal value ((o") of the original problem.
In the context of reparameterized VI, sample-average approximation be constructed by fixing set of samples {g()) ) from distribution independent of o,
p(u. (g())))
(0)
log g (T
In an SAA-based approach to reparameterized VI (Giordano et al., 2024; Burroni et al., 2024), optimization of the parameters o will move the transformed samples 20) = To(g()) to match the posterior density. whilst keeping the noise realizations ((() fixed. Empirical evaluations show that combining the SAA approximation with an off-the-sheif second-order optimizer can result in substantial computational gains as well as more reliable convergence to the optimum.
3 SAA for Forward-KL Variational Inference
The primary motivation behind existing SAA-based methods for reparameterized VI (Giordano et al., 2024; Burroni et al., 2024) is that fixing the noise realizations defines a completely deterministic surrogate objective (o)),. which is compatible with standard second-order optimization methods. The main requirement from an implementation point of view is that the model density p(y. 2) is differentiable with respect to In this setting, it is also necessary to evaluate the model for every update, since any change to o also changes the values of the transformed samples 24)) (9))

---

Algorithm 1 VISA
Input: Initial param. Po. trust region threshold data y
Zz cr (Pt)
- 0u + (2 ~ =9g1==1
b Initialize proposal parameter
40
(1) 900 100
1
Initialize samples
[ 3) - v log (1)
(o(=(())
Initialize SAA
8 2g (8a)
for I,...., do optimizer-step( Ly
43
23 (4g (3)
if di & Sz.a (6) then á - di
Not inside trust region
Update proposal
8 $zg.n
(1( 973
Refresh samples
46
Zo 166 m 11m1
[ (b) Σi=
th log pla
Refresh SAA
end end for
Figure I: Visualization of parameter traces and trust regions corresponding to different SAAs. after an update @ $ S 2,0(,) we set á + & to construct new SAA and corresponding trust region.
In developing VISA, both the motivation and implementation requirements are somewhat different. Our primary interest in minimizing the total number of model evaluations at convergence. We also wish to develop a method that is applicable when the model density p(y, 2) is not differentiable, either because the implementation simply does not support (automatic) derivatives, or because the model incorporates discrete variables stochastic control flow, which introduce discontinuities in the density p(y. 2). To this end, we propose a method that optimizes a forward KL with an importance weighted objective that incorporates ideas from SAA-based approaches.
In a setting where we already have access to samples from the posterior, we can trivially define a SAA for the upper bound Lr(g) (Ø))
1 N log (0(2())
N
P(U, z(0)
(o
p(- 11)-
In practice, this naive approach is unlikely to be usefuu in setting where evaluation of p(y, 2) is computationally expensive, since we would still need to carry out approximate inference to generate set of samples from the posterior. We therefore adopt the approach used in IWFVI, which uses the variational distribution as a proposal in a self ceff-normalized importance-sample the objective in this setting, we express the objective parameters & in terms of an expectation with approximate
To define an SAA for
respect distribution from the same family with fixed parameters, Φ. This allows the objective by means of a sample-average appoximation (d;;; Φ) with samples 2(6) ~ 93-
log p(U,2) (2)
p(U.2) 2) (2)
P(U, 2(0) (o) (x(21))
[+(#)=
E 25
1 wz(2)) log
log
p(-lv)
The quality of the approximation depends on how closely the proposal with parameters e matches the posterior. Since our approximation of the posterior will improve during optimization, we will update é to the current parameter values & at some interval, resulting in an approach that we will refer to as a sequential SAA.
To determine when we need to generate fresh SAA, we will introduce a notion of a trust region. This serves define an optimization process in which the SAA is refreshed whenever the optimization trajectory leaves the current trust region. This optimization process is illustrated in Figure I and described schematically in Algorithm 1. We begin by setting the proposal parameters Ø = du to the initial variational parameters, generating a set of samples 2 = {2(0) N i= and defining an SAA of the objective ( Φ) and a trust region Sz,((o) based on these samples. We then repeatedly update o, using an optimizer until the value Φ; no longer lies in the trust region. At this point, we update the proposal parameters & < Φe and generate a fresh sample set, which we then use to update the SAA and the trust region.

---

Defining Trust regions. The importance weight in the objective can decomposed into two parts, (1) the ratio of the variational density and the trust region density, and (2) the ratio between posterior and variational density,
(0( 93
w wz(2) - (.))() (s(2),
2 d (z)
The variance of wa is independent of the fixed proposal parameters and decreases as approaches the posterior p(= y) during optimization. Similarly, the variance of Vo.ū measures the mismatch between 93 is to qø- but can be set to zero by updating the proposal parameters to the parameters of the current variational approximation. The variance of the importance weights contributes to the variance of the gradient updates, which we generally expect to increase with increased variance of reuse samples long
the weights. Accordingly, we want trust-region criteria that allows as the proposal distribution does not differ too much from the variational distribution, to limit the itroduction of additional variance and bias compared IWFVI.
We formalize the notion of a trust region Sz,a by defining a set-valued function which, for a given threshold a and samples Z = [, maps each parameter to a corresponding trust region based on scoring function sz
(2 = to E &b ss((),)) > a}.
In other words, for given threshold a we can verify & E (z) by checking 82 $2(Ø, Φ) > a. Intuitively, the scoring function function should measure how well the proposal- and variational approximations, 93 and gø. match on particular set of samples Z. We visualize how new trust regions, corresponding to different SAAs, are constructed sequentially during optimization in Figure 1.
In this work, we propose to use the effective sample size (ESS) (Kong, 1992) as a trust-region criteria which can motivated a proxy measure of the variance the importance weight or from the perspective of minimizing a x² --divergence. More specifically. let qa be the current variational approximation and 93 the current proposal, then the following relationships hold,
(EN " (z1)
(2) (3(2)
(2(2) (x(=)
(go (g)) Vargg
1, Reff E (2)22
U "o)(())
Teff
where 2 ~ 9% (see Appendix for a derivation). We see that a higher ESS indicates a lower x² -ivveeegee and lower variance of the importance weights. Consequently, chose our scoring function
EEE2PO0(1) E.Ez (2)
8 $2 Φ) Pafff [2, 61 x))
Deff (2)
where we express the ESS as a function of the sample set and distribution parameters, and normalize such that we can set the threshold parameter a E [0. 1].
Note that, assuming sampling from the proposal distribution is cheap, we could evaluate the ESS on fresh set of samples from gz. instead of using the sample set Z. This would allow us use potentially much larger set of samples and hence obtain a better approximation of the variance of the importance weight. However, because we do not refresh samples at every iteration, we need to take into account that the actual sample set 2 is not a perfect representation of the current proposal parameter (especially in the small-sample domain). As a consequence, we want our scoring function to not only represent how well d3 matches qa. but how well the distributions match for the particular set of samples Z.
Effect of the Trust-Region on Convergence. The gradient estimator of VISA is expected to exhibit higher variance and bias than IWFVI as go is generally closer to the posterior than qz. As result, decreasing a impedes the improvement of the gradient estimates of VISA compared to IWFVI as it delays updating the proposal distribution. As a - the frequency of sample acquisition increases VISA reduces to standard IWFVI.
and, for
5

---

For a low a, the algorithm might converge to locally optimal parameter & within the trust region of the current SAA that does not yet satisfy our global convergence criteria, i.c. minimizing the forward -i--eiveeeece. In these cases, if the variational approximation is not degenerate, we can recover by increasing such that à & (ó) and continue optimization. In this work we propose to choose high enough that the algorithm does not converge prematurely to an optima of an intermediate SAA. In these cases we find that convergence of the training loss can be used indicator for convergence of the forward KL-divergence or corresponding upper bound, which we verify in our experiments. We also experimented with caching past sample sets to compute secondary loss based on the last M SAAs. While we did not find it to add additional value in our experiments, it be a useful tool to assess convergence in other settings.
Efficient Implementation. To avoid recomputing density values for old sample locations, we cache both sample location =00) and the corresponding log-joint density of the model pgg((()) y). If sampling from the proposal is cheap and memory is of concern, e.g. for large samples set or if past sample sets are stored to compute a validation loss, we can store the random seed instead of the sample and rematerialize the sample when needed.
4 Related Work
VI with SAAS. While there exist several works that incorporate SAAS in VI (Gianniotis et al., 2016; Giordano et al., 8; Wycoff et al., 2022), only recent work (Giordano et al., 2024; Burroni et al., 2024) has been focused on theoretically and empirically evaluating SAAs in the context of reparameterized black-box VI. These methods optimize reverse KL-divergence and rely on reparameterization to move samples to arcas of high posterior density while keeping a fixed set of noise realizations from the base distribution, which does not depend of the variational parameters. Optimizing the resulting deterministic objective allows the use of second-order optimization and linear response methods (Giordano al., 2015) fit covariances. This allows for substantial gains in terms of inference quality and efficiency but requires differentiability of the model.
Stochastic second-order optimization. There is also work outside of the context of SAAs that aims to incorporate second order information to improve stochastic optimization and variational inference. Byrd et al. (2016) propose batched-L-BFGG by sub-sampled Hessian-vector products instead of computing gradient differences at every iteration. This work has been further adopted to the variational inference setting by Liu & Owen (2021). Pathfinder (Zhang et al., 2022) uses a quasi-Newton method to find the mode of the a target density and construct normal approximations normal approximations are consequently used to define variational approximation that minimizes an evidence lower bound. Similar to SAA-based methods, pathfinder can significantly reduce the number of model evaluations, but requires a differentiable model.
which computes stable curvature estimates
the density along the optimization path. The intermediate
VI and trust-region optimization. There are various works that combine VI with trust region optimization, e.g. with natural gradient descent (Theis & Hoffman, 2015), policy search (Arenz al., 2018), automatic differentiation VI to leverage second-order information within trust regions (Regier et al., 2017). However, in contrast to VISA, these methods do not optimize a forward KL-divergence and impose additional constraints on the variational distribution.
VI with forward KI -ivdiieeeeeeee VISA is also related to other methods that aim to optimize a forward KL-divergence or its stochastic upper bound. This includes reweighted-wake sleep (and wake-wake) methods (Bornschein & Bengio, 2015; Le et al., 2018) to which we compare VISA in the experiment section, well as their doublyyrepprrmmetrrieed variants (Tucker et al., 2018; Finke & Thiery, 2019; Bauer & Mnih, 2021), which are not directly comparable as they require a differentiable model. While the methods above use single importance sampling step using the variational approximation as a proposal, other methods use more complex proposal including MCMC proposals (Naesseth et al., 2020: Zhang et al., 2023), approximate Gibbs kernels (Wu et al., 2020), or proposal defined by probabilitic programs (Stites et al., 2021; Zimmermann et al., 2021). While these methods do not necessarily require differentiable proposal they are not designed to be sample efficient but to approximate complex target densities.

---

Symmetric KI-IHvergenoe
VI9A a-0.99
VISA a-0.9
VISA a-0.6
VISA a-0.3
IWPVI
DDVI-SF
DEVLRP
Demmon Clos. Ding, Coo..
200
2500
5000
7500
10000
0 2500
5000
7500
10000
2500
5000
1500
10000
2500
5000
7500
10000
00
3000
1000
1000
8000
1000
cvaltuation
evaluations
coaluntions
evaluations
Ir = 0.001
Ir = 0.005
I: = 0.01
Ir = 0.05
Figure 2: Symmetric KL-divergence as a function of the number of model evaluations for a Gaussian target with diagonal covariance matrix (top row) and dense covariance matrix (bottom row). For small learning rates (0.001, 0.005, 0.01) IWFVI and BBVI-SE, need larger number of model evaluations converge. VISA converges much faster it compensates for the small step size by reusing samples. For learning rate of 0.05 VISA becomes unstable and fails to reliably converge, while IWFVI still converges. For even higher learning rate all methods but BBVI-RP fail converge.
5 Experiments
We numerically evaluate VISA on three tasks, targeting medium- to high- dimensional Gaussians, Lotka-Voltent dynamics, and Pickover Attractor model. For all experiments we use Adam (Kingma Ba, 2015) as an optimizer with the learning rates as indicated in the experiments. For the Gaussian experiment report additional results with alternative optimizers in Appendix B.
5.1 Gaussians
To study the effect of different learning rates and ESS threshold parameters, we first evaluate VISA on approximating medium- to high- dimensional Gaussians and compare the inference performance over the number of model evaluations to IWFVI, standard reparameterized variational inference (BBVI-RP) and variational inference using ecore-function gradient estimator (BBVI-SF). Notably, we include BBVI-RP as a reference only, showcasing that faster convergence can achieved by leveraging the differentiability of the model. To allow for fair comparison between methods that optimize forward KL-divergence (VISA, IWFVI) and methods that optimize a reverse KL-divergence (BBVI), evaluate the inference qualify in terms of the symmetric KL-divergence.
We study two target densities, a (1) 128 dimensional Gaussian with a diagonal covariance matrix Cdiag and (2) D = 32 dimensional Gaussian with a dense covariance matrix Cdense-
Cding = diag (Comin +
D T max Tmin
AAT (IIAAT + +0.1T)
Cdeese
D -
i=1
with Omin 0.1, o max = 1, and Aij ~ (1..111 1). Figure 2 shows the results for different learning rates Ir E (0.001, 0.005, 0.01, 0.05) and ESS thresholds E {0.3, 0.6, 0.9, 0.99). We compute gradient estimates for VISA, IWFVI and BBVI-SF with N = 10 samples, and gradient estimates for BBVI-RP using single sample. Each setting is evaluated based on 10 independent runs.
We observe that VISA converges substantially faster than IWFVI and BBVI-SF at lower learning rates (0.001, 0.005, 0.01), both for the diagonal- and dense covariance matrix. The difference in the convergence rate becomes less pronounced as the learning rate increases. For large learning rates (0.05) VISA fail to converge reliably, while IWFVI still converges. For even higher learning rates all methods but BBVI-RP fail to converge. Moreover, we observe that VISA with small converges faster in the early stages of training but can fail converge fully in the later stages of training the thresholds is too low (see variance as samples are not refreshed frequently enough to prevent overfitting to high weight samples. Comparing convergence across learning rates and ESS thresholds, we observe that VISA with a 0.3 converges faster or at same rate as IWFVI with the same or higher learning rate.
0.3 in Figure 2). This is a result of underestimating posterior
Python code for the experiments is available ht https:: : pttp:://giii //gtth///./ttttt// com/zmheiko/visa
7

---

VISA 0-0.99
VISA 0-0.95
VISA a-0.9
IWEVI
VISA 0-0.999
bas IWFVI
That Lowe Objjectiion
100
-00
-0
-20
-20
20
1000
ovslusttioss 0.001)
exsilaation
0.01)
Figure 3: Results for Voltera model with different learning rates. (Top row) Training objective over number of model evaluations. (Middle row) Approximate forward KL--dvergence computed on reference samples obtained by MCMC. For smaller step sizes (0.001, 0.005) VISA achieves comparable forward KL-divergence to IWFVI while requiring significantly less model evaluations converge (see vertical lines). For larger step sizes (0.01) VISA only converges with a high ess threshold (0.99) and requires more evaluations than IWFVI and VISA with a smaller step size (0.005).
5.2 Lotka-Volterry
The Lotka-Von Vorterra predator- prey population dynamics (Lotka, 1925; Volterra, 1927) are modeled by a pair of first-order ordinary differential equations (ODEs),
du dt (a Bu)u,
dv (-7 ōu)c,
dt
where denotes the predator- and u denotes the prey population. We will in the following denote the pair predator-prey populations at time with 2 = (ue, Uz). The dynamics of the ODE are governed by its population growth and shrinkage parameters 0 = (ax, 7. 6), which we want to infer together with the initial conditions of the system given noisy observations yır = (y1,..., ..r) Following Carpenter (2018), we place priors over the initial population sizes 20 and system parameters,
=prey
= -rprrd [ogNormal(log(()), 1), o, Normal 0.5), 3,5 N0.0al 0.05),
and assume a fractional observation error,
prey pred
LogNormaa(log Zt, Oz),
0, LogNormal(1111)
Given an initial population 20, system parameters and observations y1:r we can solve the ODE numerically obtain approximate population sizes 21.7 for time steps 1, ...... T which we to compute the likelihood of the observed predator-prey populations, p(hi: 20, 0) II= =))
Our goal is to learn an approximation to the posterior p(0, =0 y) by minimizing the upper bound
(Ø( = (2.,0.))pp( log p(z0: y, 0)
go (2. 0)
We model the variational approximation qa for the interaction parameters 0 and initial population sizes zn as jointly log-normal and initialize ø such that the the marginal over 2g matches the prior and the marginal over has similar coverage to the prior. We approximate the objective and its gradient using N = 100 samples. To specify common convergence criteria, we compute the highest common test loss value for VISA with a 0.99 and IWFVI that is not exceeded by more than I nat by all consecutive test loss values. The convergence threshold computed this way is 7712.7 nats.
To evaluate the inference performance, we first generate N 100 000 approximate posterior samples using a No-U-Turn Sampler (NUTS) (Hoffman & Gelman, 2014) with 100 000 burn-in steps and window-adaption, which generally provides good performance of the box. We use the approximate posterior samples to approximate an "oracle" for the upper bound in Equation (3),
N
p(=6 y. g(0)
L NUTS (6) =
log
N
qs (=(=))) 0(0)
NUTS is an adaptive Hamiltonian Monte Carlo sampler and uses the gradient information of the model to guide the generation of proposals. As such it requires the log-joint density model to be differentiable which is required by VISA IWFVI
8

---

which we evaluate along with the training objective during optimization to assess convergence.
We find that VISA is able obtain variational distribution of similar quality IWFVI while requiring significantly fewer model evaluations for smaller learning rates (see Figure 3). Interestingly, IWFVI requires significantly fewer model evaluations per gradient step during the early stages of training, while requiring slightly more evaluations per gradient step thereafter. We hypothesise that this is again result of under approximating posterior variance in the later stages of training. As result, even small changes in the variational distribution can lead to big changes in the ESS, which triggers the drawing of fresh samples. For IWFVI, also find more pronounced difference between different ESS thresholds and their influence on convergence. Runs with higher ESS thresholds converge more stably and able achieve lower test loss in the final stages of training. The reported means and standard are computed based on 10 independent runs.
5.3 Pickover Attractor
Following Rainforth et al. (2016), we model 3D Pickover attractor (Pickover, 1995) with system parameters (3, n) (where 22 IR3 is vector and the superscript denotes the component),
1 sin(3x?) (2s(2.5x)
# sin(1.5x) )a} cos(nzt), ) r sin(x)).
Due to its chaotic nature the system is sensitive to small perturbations in its initial state, i.c. small variations in the initial state lead to exponentially fast diverging trajectories. Therefore, to track the evolution of the system, noisy observations == yı.r and introduces auxiliary variables = 21.r to model the latent state of the system. We define the prior system parameters
employ bootstrap particle filter (Gordon et al., 1993) which assumes
{ 1/18 0
-3 0 S 3
p(0)
otherwise
and model the transition and observation as,
20 ~ W- 03, lla),
2 ~~~ (),
Yt ~ N 22, oy)
where a₂ 0.01, o, = 0.2, and h evolves the system by one time step using the equations of the Pickover attractor described above. The particle filter is used to simulate T = 100 time steps with M 500 particles, which renders evaluating the model expensive. To restrict the proposal to the same domain as the prior by first sampling from a Normal, which we denote Ts and parameterize by mean and lower Cholesky factor, and then transform the sample by
f(0) (tanh(3 01), tanh 02 + 1.5)).
The density of the transformed samples is (0) (f)]dct of 1-1 which appropriately restricted to the domain of the prior.
initialioos.om VAr, ipprixx
[log Pla.
VISA
IWFVI
-2000
100
200
300
2
2
ground Truth attractor (-231 1,25)
average AEtractor (- 2.33, 1.29)
evaluations
Figure 4: Results for Pickover attractor. (a) Approximate log-joint density over number batch- evaluations of model. (b) Log-joint approximation plotted over domain of prior. The variational approximation capture the high density area containing the data. (c) Visualization of pickover attractor with ground truth parameters system parameters computed over 10.000 samples from the learned variational approximation. Each evaluation in the plot corresponds to evaluating a batch of N 10 samples.
|-2.3, 1.25). (d) Visualization of attractor with average

---

We are interested in approximating the marginal posterior p(0 y) over system parameters by optimizing the evidence upper bound
Po (y, (v,0) (o(0)
By pt
PPy [0]] (o(
log ily 00
0
E (0(0))
log
E (o(0lv)
log
E
Ep
(e(0lv)
(5(
To obtain tractable objective we replace the intractable marginal likelihood ply 0) ~ P 8) with the marginal likelihood estimate obtained by running the particle filter (Naesseth et al., 2019), similar to pseudo-marginal methods (Andrieu et al., 2010) and approximate the gradient with N = 10 samples. As the likelihood estimate is non-differentiable due to the discrete ancestor choices made inside the particle filter, we cannot run NUTS to obtain approximate posterior sample as before, but instead report the log-joint density of the variational distribution.
We observe that VISA converges more stably with fewer samples compared to IWFVI and find that attractors corresponding to samples from the variational approximation look qualitatively similar to those based on the true parameters. We summarize the result in Figure 4.
6 Discussion and Limitations
In this paper we developed VISA, method for approximate inference for expensive to evaluate models that optimizes the forward KL-divergence through a sequence of sample-average approxima- tions. Each SAA is optimized deterministically and fixes a single set of samples, hereby requiring new model evaluations only when the SAA is refreshed. To track the approximation quality of the current SAA, VISA computes the ESS of the ratio between the current variational distribution and the proposal distribution that was used to construct the SAA. If the ESS falls below a predefined thresh- old, a new SAA approximation is constructed based on fresh samples from the current variational distribution. We observe a significant reduction of required model evaluations for conservatively chosen step sizes, while achieving similar posterior approximation accuracy as IWFVI, the equivalent method that does not employ the sequential sample-average approximation. In the following we are discussing the limitations of our method.
Underapproximation of posterior variance. Both reparameterized VI, which optimizes the reverse KL-divergence, and IWFVI, which optimizes the forward KL-divergence via importance sampling, are prone to under approximating posterior variance. In the case of reparameterized VI, this can often be attributed to the mode seeking behaviour of the reverse KL-divergence, while in IWFVI the low effective samples sizes can lead to over fitting to a small number of hight samples. We found that keeping the samples fixed for too long, i.c. using ESS threshold that is too low, can exacerbate this problem, as the optimizer can take multiple steps towards the same neighh samples.
Giordano et al. (2024) and Burroni et al. (2024) showed that when applying SAA to reparameterized VI, is possible to make use of second-order methods. We experimented with optimizing SAAS with L-BGFS, which is a quasi-Newton method with line search. However, found that in the setting of optimizing forward KL with relatively few samples, L-BGFS can amplify the problem of overfitting, often leading to instabilities and collapsed variational distributions.
Number of latent variables and parameters. Because VISA employs a relatively small number of samples and does not refresh samples at every iteration, we found that it is not well-suited to models with a large number of latent variables or large number of parameters. This agrees with theoretical findings by Giordano al. (2024), who show that SAAS for full covariance Gaussian fail if the number of samples is not at least in the same regime as the number of latent dimensions. Burroni et al. (2024) manage to train full covariance normal approximation by using sequence of SAAS using increasingly large sample sizes, however, this is directly opposed to our goal of reducing the number of model evaluations in expensive to evaluate models.
Impact Statement
This paper presents work whose goal is to advance the field of Machine Learning. While we recognize that the advancement of machine learning methods can have broad societal impact, we do not foresee specific consequences of this particular research.
10

---

Acknowledgement
The authors would like to thank Tamara Broderick for helpful discussion about sample-average approximations in the context of reparameterized variational inference.
References
Andrieu, C., Doucet, A., and Holenstein, R. Particle Markov Chain Monte Carlo Methods. Journal of the Royal Statistical Society Series B: Statistical Methodology. 72(3):269 342, June 2010.
Arenz, O., Neumann, G., and Zhong, M. Efficient gradient-free variational inference using policy search. In International conference machine learning, pp. 204 243, PMLR, 2018.
Bauer, M. and Mnih, A. Generalized Doubly Reparameterized Gradient Estimators. International Conference on Machine Learning, pp. 738-747. PMLR, 2021.
Bornschein, J. and Bengio, Y. Reweighted wake-sleep. In International Conference on Learning Representations, 2015.
Burroni, J., Domke, J., and Sheldon, D. Sample Average Approximation for Black-Box VI. In Uncertainty in Artificial Intelligence, 2024.
Byrd, R. H., Hansen, S. L.,, Nocedal, J., and Singer, Stochastic Quasi-Newton Method for Large-Scale Optimization. SIAM Journal Optimization, 26(2):1008-1033 January 2016. ISSN 1052-6234. doi: 10. 137/14495944...
Carpenter, B. Predator-Prey Population Dynamics: The Lotka- Volterra model in Stan. https://mc- setaaarrrratttrrreeeeee January 2018.
Duane, S., Kennedy, A. D., Pendleton, B. J., and Roweth, D. Hybrid monte carlo. Physics letters B, 195(2): 222, 1987.
Finke, A. and Thiery, A. H. On importannce-wiighted autoencoders. arXiv preprint arX1v:1907.10477. 2019.
Gianniotis, N., Schnörr, C., Molkenthin, C., and Bora, S. S. Approximate variational inference based on finite sample of gaussian latent variables. Pattern Analysis and Applications, 19:475-485. 2016.
Giordano, R., Broderick, T., and Jordan, M. Covariances, robustness, and variational bayes. Journal of machine learning research, 19(51):1-49 2018.
Giordano, R., Ingram, M., and Broderick, T. Black box variational inference with deterministic objective: Faster, more accurate, and even more black box. Journal of Machine Learning Research, 25(18):1-39, 2024.
Giordano, R. J., Broderick, T., and Jordan, M. I. Linear Response Methods for Accurate Covariance Estimates from Mean Field Variational Bayes, In Advances in Neural Information Processing Systems, volume 28. Curran Associates, Inc., 2015.
Glynn, W. Likelihood ratio gradient estimation for stochastic systems. Communications of the ACM, 33(10):75-84, October 1990. ISSN 0001-0782. doi: 10.1145/84537.844 84455..
Gordon, N. J., Salmond, D. J., and Smith, A. F. M. Novel approach to nonlinearnoon-Gaussian Bayesian state estimation. IEE Proceedings F (Radar and Signal Processing). 140(2):107-113, April 1993. ISSN 2053-9045. doi: 10. (0299/
Hoffman, M. D. and Gelman, A. The No-U-Turn Sampler: Adaptively Setting Path Lengths in Hamiltonian Monte Carlo. Journal of Machine Learning Research, 15(47):1593-1623, 2014. ISSN 1533-7928.
Kim, S., Pasupathy, R., and Henderson, S. Guide to Sample Average Approximation. Handbook of simulation optimization, pp. 207-243, 2015.

---

Kingma, D. P. and Ba, J. Adam: A method for stochastic optimization. In Bengio, Y. and LeCun, Y. (eds.), 3rd International Conference on Learning Representations, ICLR 2015, San Diego, CA, USA, May 7-9, 2015, Conference Track Proceedings, 2015.
Kingma, D. P. and Welling, M. Auto-Encoding Variational Bayes. ICLR 2014 conference submission, December 2013.
Kong, A. note importance sampling using standardized weights. University of Chicago, Dept. of Statistics, Tech. Rep, 348:14, 14, 1992.
Le, T. A., Kosiorek, A. R., Siddharth, N., Teh, Y. W., and Wood, F. Revisiting Reweighted Wake- Wake-Sleep. International Conference on Learning Representations, September 2018.
Liu, S. and Owen, A. B. Quasi-Monte Carlo Quasi-Newtrr in Variational Bayes. Journal of Machine Learning Research, 22(243):1 23,, 2021. ISSN 1533-7928.
Lotka, A. J. Principles of physical biology. 1925.
Naesseth, C., Lindsten, F., and Blei, Markovian Score Climbing: Variational Inference with KL(p I1 q). In Advances in Neural Information Processing Systems, volume 33, pp. Associates, Inc., 2020.
Curran
Naesseth, C. A., Lindsten, F., and Thomas B., S. Elements of Sequential Monte Carlo. Foundations and Trends in Machine Learning. 2019.
Neal, R. M. et al. Meme using hamiltonian dynamics. Handbook of markov chain monte carlo, 2 (11):2, 2011.
Nemirovski, A., Juditsk,, A., Lan, G., and Shapiro, A. Robust stochastic approxiittinn approach to stochastic programming. SIAM Journal on optimization, 19(4):1574-100, 2009.
Paisley, J., Blei, D., and Jordan, M. Variational bayesian inference with stochastic search. In Proceedings of the 37th International Conference Machine Learning, 2012.
Pickover, C. A. The Pattern Book: Fractals, Art And Nature. World Scientific, June 1995. ISBN 978-981-4504-03-4
Rainforth, T., Le, T. A., van de Meent, J., Osborne, M. A., and Wood, F. Bayesian optimization for probabilistic programs. Advances in Neural Information Processing Systems, 29, 2016.
Ranganath, R., Gerrish, S., and Blei, D. Black Box Variational Inference. In Proceedings of the Seventeenth International Conference on Artificial Intelligence and Statistics, pp. 814-822. PMIR, April 2014.
Regier, J., Jordan, M. I., and McAuliffe, J. Fast black-box variational inference through stochastic trust-region optimization. Advances in Neural Information Processing Systems, 30, 2017.
Rezende, D. J., Mohamed, S., and Wierstra, D. Stochastic Backpropagation and Approximate Inference in Deep Generative Models. In Proceedings of the 31st International Conference on Machine Learning, pp. 1278-1286. PMLR, June 2014.
Stites, S., Zimmermann, H., Wu, H., Sennesh, E.., and de Meent, J.-W. Learning proposals for probabilistic programs with inference combinators. Proceedings of the Thirty-Seventh Conference on Uncertainty Artificial Intelligence, pp. 1056-1066. PMLR, December 2021.
Theis, L. and Hoffman, M. A trust-region method for stochastic variational inference with applications to streaming data. In International conference machine learning, pp. 2503-2511, PMLR, 2015.
Titsias, M. and Lázaro-Gredilla, M. Doubly stochastic variational bayes for non-conjugate inference. In International conference machine learning, pp. 1971-1979, PMLR, 2014.
Tucker, G., Lawson, D., Gu, S., and Maddison, C. J. Doubly Reparameterized Gradient Estimators for Monte Carlo Objectives. In International Conference Learning Representations, 2018.
12

---

Volterra, V. Fluctuations in the Abundance of a Species considered Mathematically. Nature, 119 (2983):12-13, January 1927. ISSN 1476-4687. doi: 10.1038/1 19012b0.
Wingate, D. and Weber, T. Automated variational inference probabilistic programming. arXiv preprint arXiv:1301. 1299, 2013.
Wu, H., Zimmermann, H., Sennesh, E., Le, T. A., and Meent, J.-W. V. D. Amortized Population Gibbs Samplers with Neural Sufficient Statistics. In Proceedings of the 37th International Conference on Machine Learning, pp. (0421-10431 PMLR, November 2020.
Wycoff, N., Arab, A., Donato, M., and Singh, L. O. Sparse bayesian lasso via variable-cooeffcieent & penalty, arXiv preprint arXiv:2211.05082 2022.
Zhang, L., Carpenter, B., Gelman, A., and Vehtari, A. Pathfinder: Parallel quasi-Newiin variational inference. Journal of Machine Learning Research, 23(306): 1-49, 2022. ISSN 1533-7928.
Zhang, L., Blei, D., and Naesseth, C. A. Transport Score Climbing: Variational Inference Using Forward KL and Adaptive Neural Transport. Transactions on Machine Learning Research, May 2023. ISSN 2835-8856.
Zimmermann, H., Wu, H., Esmaeili, B., and van de Meent, J.-W. Nested Variational Inference. In Advances in Neural Information Processing Systems, volume 34, pp. 20423-20435. Curran Associates, Inc., 2021.
13

---

A Motivation of trust region.
Using the ESS as trust-region criteria can be motivated from a perspective of minimizing x² divergence or the variance of the corresponding importance weight. Let be the current variational approximation and q the current proposal, then the following relationship holds,
[E) E== W)
((2)
xx² q)
Vara
Reff
((2)
Refff
2i=1 w
Note that the ESS criteria defined in section 3 uses the normalized ESS and hence differs from Neff by a factor n,
1
1.
(4)
Teff
82(0.9)
Proof.
Dx=(0 () = dz (q(z) - {(z))²
d2 {(2)2
q(2)
((2)
dz
dz
((z)
g(z)
= BB [(()
(())
27
- 1 = Eg
-1
[(2)
((2)
((()
(())
7(2))
((=)
27
Var;
((=)
= E,
((z))
7(z)
7(2)
-
9(2))
27
E5
((x)
((=)
Σ 12
w
Eg
1
-1
((z)
Eg
[(x))
E w/
2
71eff
((())
Kish design effect
B Additional results for Gaussians
We evaluate the performance of VISA with two additional optimizers, R.SPrrp and SGD. We use the implementation and standard parameters provided by Optax. We find that VISA converges faster with Adam and RMSProp compared to SGD. We also observe that, for a small learning rate, SGD can get stuck in local optima, while RMSProp gets unstable more quickly with higher learning rates.
Symmetric KL-Divergence
VIBA-SGD a-0.99
VISA-Adum a-0.99
VISA-RMSPPos a=0.99
VISA-SGD a-0.8
VISA-Adam
VISA-RM5Prop
VIBA-SCD a-0.6
VISA-Adam a-0.6
VISA-RMSPPro a-0..6
Donse Cov. Ding, Coo..
200
2000
4000
6000
8000
10000
2000
600D
8000
10000
0000
8000
10000
50
1000
2003
3000
1000
2000
5000
1000
2000
3000
ovsllsstions
evalustions
evaluations
Ir
0.001
Ir 0.005
Ir 0.01
14

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: Claims are supported by experiments in Section 5.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how be expected to generalize other settings.
much the results
is fine to include aspirational goals as motivation as long is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: Limitations are discussed in Section 6.
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those not discussed in the paper.
The authors are encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification. asymptotic approximations only holding locally). The authors should reflect how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low images are taken in low lighting. Or speech-to-text system might be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach address problems of privacy and faimess.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory Assumptions and Proofs
Question:
For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [NA]
15

---

Justification:
Paper does not include theoretical results.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but they appear in the supplemental material, the authors are encouraged to provide short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental Result Reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless whether the code and data are provided or not)?
Answer: [Yes]
Justification:
Information reproduce results is provided in Section 5.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make clear how to reproduce that algorithm. (b) the contribution is primarily a new model architecture, the paper should describe the architccture clearly and fully.
(c) If the contribution is new model (e.g., a large language model), then there should either way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reprodueibility. the case of closed-source models, may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers have some path to reproducing or verifying the results.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
16

---

Answer:
[Yes]
Justification:
link to the code is provided in the paper.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (https://nips.. cc/ Sguuides/oodeSummt
plyccy for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https: / //iips.cooc wttpcccccppgpppaaaeeece /ttpa//ggiieeeeeeeubbiiiiiiii
/Codeeubeeeeeeeeeee cettuubmissiiii sesiisssisnnoooii seiiinyylcc)) for more details.
The authors should provide instructions on data access and preparation, including how access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: Training and test details are provided in Section 5. Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, or as supplemental material.
7. Experiment Statistical Significance
Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer: [Yes]
Justification: We report error bars.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, statistical significance tests, at least for the experiments that support main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
should be clear whether the error bar is the standard deviation the standard error of the mean.
17

---

It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of enrors is not verified.
For asymmetric distributions, the authors should be careful to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative eтror rates).
If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of exccution) needed to reproduce the experiments?
Question:
Answer: [No]
Justification: Experiments do note require specialized hardware or extensive computational resources.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported the paper (e.g., preliminary or failed experiments that didn't make into the paper).
9. Code Of Ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https http///eerrppppeeeee wwwpuurpss.cppuuuue http://blic//hhicccuueeeeeeeeeee
Answer: [Yes]
Justification:
Our work conforms wit the NeurIPS Code of Ethics.
Guidelines:
The answer NA means that the authors have not reviewed the NeurlPS Code of Ethics.
If the authors answer No, they should explain the special cireumstances that require deviation from the Code of Ethics.
The authors should make to preserve anonymity (e.g.. if there is special consid- cration due to laws or regulations their jurisdiction..
10. Broader Impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [NA]
Justification:
No societal impact.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
18

---

The conference expects that many papers will be foundational research and not tied particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology being used intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology.
If there negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms monitor how system learns from feedback over time, improving the efficiency and accessibility of MI..
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: The paper does not pose such risks.
Guidelines:
The answer NA means that the paper poses such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Lieenses for existing assets
Are the creators or original owners of assets (e.g., code, data, models), used in the paper. properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [NA]
Justification: No existing assets used.
Guidelines:
The answer NA means that the paper does existing assets.
The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset is used and, if possible, include URL.
The name of the license (c.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (c.g., website), the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswi thcode. com com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided.
19

---

If this information is not available online, the authors are encouraged to reach out to the asset's 's creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification: We will release documented code upon publication.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/moodl as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsoureing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According the NeurIPS Code of Ethies, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurlPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
20


Reviews:
[{"review_id": "ZHofV9KDR5", "reviewer": "Reviewer_snbt", "summary": "The authors propose to reduce the number of model evaluations during the optimization of the variational lower bound. To this end, they integrate sample avarage approximations (SAAs) into IWFVI framework, which allows updates of variational parameters while keeping approximate samples generated from a variational distribution obtained a few steps before during optimization. The method was evaluated on a few synthetic experiments.", "questions": "See weakness.", "limitations": "The authors thinks there is no social impact. However, no justifications are given.", "rating": 4, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "1. The method is effective in the experiments considered (it requires less model evaluations to converge).\n2. The clarity of the paper is good and is easy to follow.", "weaknesses": "1. The method seems to be a combination of a few engineering tricks. It is unclear what the effect of inconsistency between approximate samples ($\\tilde{\\phi}$) and variational distribution ($\\phi$) during optimization is.\n2. Only synthetic experiments are considered.\n3. In the paper, the model is assumed to have no tubable parameters, while many applications requires tunable model as well (e.g. VAE).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "**Clarification regarding our methodology.**\nVISA does not optimize the reverse KL-divergence or corresponding variational lower bound but a forward KL-divergence or corresponding variational upper bound. This is a crucial difference to previous related work, which studies SAAs in the context of the reparameterized reverse KL-divergence.\n\n**Presentation.**\nGiven the clarity of our paper was highlighted as a strength, we are slightly confused by the low score for the presentation. Please let us know in which ways we can further improve our presentation.\n\n**Inconsistency between proposal and variational distribution.**\nThe discrepancy between the proposal and variational distribution is captured by the importance weight $q(z \\mid \\phi) / q(z \\mid \\tilde \\phi)$, which contributes to the overall importance weight used in the sample average approximation and corrects for the error introduced by not sampling directly from $q(z \\mid \\phi)$. In general, the bias and variance of the gradient estimates will be higher if the proposal has not been updated recently and the discrepancy to the variational distribution is large. To quantify the effect that this discrepancy has on the overall performance of the algorithm numerically we conducted experiments with various ESS thresholds $\\alpha$, which are a way to formalize an upper bound on the maximum allowed discrepancy between the proposal and variational distribution.\n\n**Only synthetic experiments.**\nSee general response.\n\n**No tunable model parameters.**\nVAEs are trained by jointly optimizing a generative model and a variational posterior approximation by maximizing an ELBO in the context of amortized variational inference. In contrast, our work considers posterior inference for a given model, a setting that is frequently encountered when working with simulation-based models, and treats parameter estimation of that model as an orthogonal problem.", "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": "", "extra_scores": {}}, {"review_id": "RVOiXfKNzM", "reviewer": "Reviewer_kZyU", "summary": "The paper proposes to use a sample-average approximation for variational inference. In contrast to previous works, the method uses the forward-KL which does not require differentiability of the joint likelihood. In order to sample from an approximate posterior instead of the exact posterior, a sequential trust-region approach is introduced. The method is evaluated on small-scale problems and some benefits of the forward-KL SAA approach compared to BBVI or reverse-KL SAA are demonstrated.", "questions": "- When compared to other SAA methods (Giordano et al 2024), (Burroni et al 2023), the novelty of the proposed method seems to be in the use of forward KL rather than reverse KL.  I can see the computational advantages of using forward-KL, but does the mode-covering property of forward-KL pose problems sometimes? For example, when the approximate posterior is simple but the true posterior is multimodal with high loss inbetween? \n- Natural-gradient VI methods can also be used in black-box settings -- see for instance \"Variational Adaptive-Newton Method for Explorative Learning\", Khan et al., 2017. Related natural evolution strategies (such as CMA-ES) have shown state-of-the-art performance in optimizing loss functions with only gradient evaluations. Have you considered comparing to these methods, or would you expect them to perform similary as BBVI-SF?", "limitations": "Limitations are adequately discussed in the last section of the paper.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "- The derivatinos and mathematical formulation are clear, the paper is easy to read. \n- The proposed trust-region approach is rigorously evaluated on small-scale examples with various ablation studies.", "weaknesses": "- The proposed VISA method (and other baselines as well) seem quite sensitive to the learning rate and other tuning parameters. Therefore it is hard to judge whether the proposed approach has significant benefits over other methods, or whether the benefits are due to a tuning of the learning rate. \n- The relevance of the considered applications to the wider machine-learning audience seems unclear. The work is motivated by non-differentiability but then only toy examples are shown. Some real-world applications (for example, reinforcement learning) would make the paper stronger.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "**Sensitivity to the learning rate.**\nSee general response.\n\n**No real-world experiments.**\nSee general response.\n\n**Failure cases of the Forward KL-divergence.**\nWhile we consider the mode-covering behavior of the forward KL-divergence as a positive, there can be instances where a mode of the variational approximation covers multiple modes of the target density and as a result assigns not insignificant mass to a low density area. To prevent these scenarios it is important to choose an appropriate variational approximation that allows to model the expected degree of multimodality. Moreover, in practice, optimizing an importance-weighted approximation to the forward KL-divergence often results in an under approximation of posterior variance similar to optimizing a reverse KL-divergence, especially when targeting well separated modes in the small sample regime.\n\n**Adding more specialized baselines.**\nWe chose to compare VISA to IWFV and BBVI-SF because these methods are equally back-box, in the sense that they do not assume differentiability of the model and are not restrictive in the variational families that can be used. However, in settings where more specialized methods can be applied, we expect these methods to outperform our baselines and potentially VISA. One instance of this can be seen in the first experiment where we included a BBVI-RP baseline, which is generally considered a black-box method but leverages the differentiability of the model. \n\nThe suggested baseline, VAN and CMA-ES, specifically use a multivariate Gaussian as a variational approximation and VAN-D makes use of a reparameterization trick to approximate the diagonal of the Hessian, which requires a differentiable model. That being said, we agree that it would be insightful to compare to and potentially combine VISA with a black-box (in the sense described above)  stochastic second-order method. We will work to include an appropriate method in the final manuscript.", "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": "", "extra_scores": {}}, {"review_id": "fHoUWxAzMU", "reviewer": "Reviewer_Bkjs", "summary": "This paper introduces a new method for variational inference called VISA, which stands for Variational Inference Using Sequential Sample-Average Approximations. VISA is based on importance-weighted forward-KL variational inference and allows for reusing model evaluations across multiple gradient steps. This makes it suitable for non-differentiable and computationally intensive models. The authors conducted experiments using simulated and real data to validate VISA, showing that it can achieve similar approximation accuracy to standard IWFVI while requiring significantly fewer samples.", "questions": "1.\tCould the algorithm's performance be enhanced by incorporating adaptive adjustments for the ESS threshold $\\alpha$ and the learning rate? \n2.\tIt appears there may be a typo in the caption for Figure 3 when referring to the methods VISA and IWFVI. Could you please verify and correct this if necessary to ensure clarity and accuracy in the presentation of your results?", "limitations": "The authors have adequately discussed the limitations of their approach.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1. The paper is well-written and organized, making complex concepts easy to understand through a logical flow of information.\n2. VISA effectively addresses the challenge of minimizing the number of model evaluations, thereby achieving significant computational efficiency. \n3. The experiments conducted provide robust evidence that VISA converges faster, requiring fewer model evaluations compared to existing methods.", "weaknesses": "1. While the concept of trust regions based on Effective Sample Size (ESS) is intuitive, the paper would benefit from a more formal theoretical justification for the specific choice of trust region. \n2. The paper currently lacks a discussion on the robustness of the proposed method to various optimization methods. While it touches on the use of L-BGFS, including an analysis of how VISA performs with other optimization techniques such as classical stochastic gradient descent (SGD) and adaptive methods like RMSProp and Adam would enhance the comprehensiveness of the study.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "**Motivation for ESS as trust-region criteria.**\nSee general response.\n\n**Robustness to optimization methods.**\nWe indeed use Adam (optax implementation) for all experiments and will clarify this in the final version of the manuscript. Based on your suggestion, we conducted additional experiments for the Gaussian target densities to demonstrate that VISA can be used with SGD, Adam, and RMSProp for appropriately selected stepsizes. Our evaluations shows that plain SGD (without momentum or Nesterov acceleration) might run into local optima more easily as Adam and RMSProp for small step sizes, while RMSProp might not be as robust to larger step sizes as Adam. Overall, the experiment suggests that Adam is a good standard choice for VISA. We have included the corresponding convergence plots in the PDF file attached to the general rebuttal and will add a polished version to the final manuscript.\n\n**Adaptive ESS threshold.**\nWhile we did not test adaptive strategie, we experimented with using a schedule for the ess threshold. However, we did not find consistent improvement over using a fixed threshold and hence opted for the simplest possible presentation of the algorithm. That being said, we expect that problem specific scheduling and adaptation strategies are able to outperform fixed ess threshold. However, determining a performant schedule or adaption hyperparameters might require running inference multiple times, which contributes to the overall number of evaluations and complicates a fair comparison.\n\n**Typo in Figure 3.**\nThere was indeed a typo, thanks for pointing it out! The caption should read: \"For smaller step sizes (0.001, 0.005) VISA achieves\ncomparable forward KL-divergence to IWFVI while requiring significantly less model evaluations\nto converge (see vertical lines). For larger step sizes (0.01) VISA only converges with a high ess\nthreshold (0.99) and requires more evaluations than IWFVI and VISA with a smaller step size (0.005).\"", "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": "", "extra_scores": {}}, {"review_id": "XFu8GoA4Pl", "reviewer": "Reviewer_sEay", "summary": "This paper presents a method for approximate inference in expensive to evaluate models where the gradients may not be available. The presented approach proceeds in a trust-region-optimization fashion where a series of deterministic surrogates objectives are optimized. Each objective is optimized until a trust region condition is met, indicating that the objective needs to refreshed. The proposed approach shows promising results when compared to other competing approaches.", "questions": "Please, see the Weaknesses section.", "limitations": "The scalability of the proposed approach is rather limited as pointed out by the authors. The proposed method seems extremely sensitive to the choice of learning rate and the alpha. It is unclear how to chose them properly.", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "The paper presents a simple procedure that is clearly laid out. It is well written for most parts and seems to be promising direction for future research.", "weaknesses": "- I think the use of ESS for determining the trust region is rather arbitrary. I did not find a sound explanation for why this should be the goto way for determining the trust region. \n- Also, the analysis around the choice of alpha seems rather experimental. The performance seems extremely critical to the choice of this parameter and there doesn't seem to be any good way to select this parameter. \n- Moreover, rerunning the optimization to select this parameter can result in a lot of evaluations. Since we wanted to avoid this in the first place, would a user not be better of using the IWFVI approach?\n- Also, what exactly is the difference between the IWFVI approach and the standard Re-weighted wake sleep? Why rename an already named algorithm? \n- I did not find enough details about the exact optimization procedure used to optimize the VISA approach. Which optimizer was used? \n- I am also not sure how are the number of evaluations calculated for the Figure 2. Are we only counting the evaluations during the refresh step? If so, should we not see step decays in the metrics for VISA approaches?", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "**Motivation for ESS as trust-region criteria.** \nSee general response.\n\n**Choice of ESS threshold.** \nSee general response.\n\n**IWFVI vs Reweighted Wake-Sleep (RWS).** \nRWS is used in the context of amortized variational inference and optimizes both the parameters of the model and the parameters of the (amortized) variational approximation. Moreover standard RWS optimizes different objectives to update the variational approximation for the sleep- and wake-phases of the algorithm, which use simulated and real data for training respectively. IWFVI only optimize the variational distribution and does not use simulated data from the model (sleep-phase in RWS). Hence, IWFVI updates correspond to the wake-phase updates to the variational parameters in RWS only; both are based on optimizing an importance-sampling based estimate of the forward KL-divergence.\n\n**Details on exact optimization procedure.**\nWe use Adam with standard parameters (as implemented in optax) for all experiments. We will clarify this in the final version of the manuscript.\n\n**Number of evaluations.**\n It is correct that only freshly drawn samples are counted towards the number of evaluations, hence the number of evaluations only increases after a “refresh step”. The reason why there is no visible staircase pattern in the corresponding plots is that VISA refreshes after a relatively small amount of evaluations (for conservatively chosen ess thresholds) and that the plots are averaged over multiple runs. For a single VISA run with a low enough ess threshold we would indeed see a visible staircase pattern.", "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": "", "extra_scores": {}}]

Author rebuttal:
We would like to thank the reviewers for their detailed reviews! We are delighted to see the clarity of our manuscript and rigor of our experiments listed among the strengths of our work. We address the point raised by multiple reviewers below and respond to the individual concerns and questions in dedicated rebuttals under the corresponding reviews.

**Sensitivity to learning rate and choice of ESS threshold (kZyU, sEay).**
We believe that one of the strengths of VISA is its robustness to the choices of learning rate for conservatively chosen thresholds $\alpha$. As seen in the first two experiments, the rate of convergence of VISA is relatively stable across all but the largest presented learning rates, which are chosen too large on purpose. In contrast, the convergence rates of IWFVI and BBVI are more sensitive to the learning rate and only competitive when chosen close to the upper limit in terms of training stability. In general, while finding the optimal training hyperparameters for either method might require careful tuning, we found that a simple recipe for selecting the learning rate and ess-threshold for VISA is often enough to outperform IWFVI:
1) Choose the same learning rate as for IWFVI and a conservative ess-threshold, e.g. $\alpha=0.99$. 
2) If convergence is stable, one might choose to reduce $\alpha$ to achieve even faster convergence, otherwise decrease the learning rate slightly.


**Motivation for ESS as a trust-region criteria (Bkjs, sEay).**
There is in fact a good reason to use the ESS to define the trust region, and we will clarify this in the manuscript. The Kish effective sample size is directly related to the variance of the importance ratio between the current variational distribution and the trust region distribution, which can in turn be expressed as a Chi-Square divergence between these distributions.

$
\chi^2(q(\cdot; \phi) \mid q(\cdot; \tilde{\phi})) =
 \mathrm{Var}_{q(\cdot; \tilde{\phi})} \left[\frac{q(z; \phi)}{q(z; \tilde{\phi})}\right] 
\approx \frac{n}{\mathrm{ess}} - 1
$

As a result, the ESS-based trust region can be motivated from (1) a divergence perspective. i.e. we want the sampling distribution close to the variational distribution as measured by the chi-square divergence, or from (2) the perspective of the variance of the corresponding importance weight (or ESS as a proxy measure), which we want to be low (or high ESS respectively). The second perspective also motivates the use of the ESS to assess the quality of importance samplers or adaptive resampling strategies based on an ESS criteria. 

We will expand upon the motivation of the effective sample size as a trust-region criteria in the manuscript and will add deviations for the above equalities to the appendix.

**Only synthetic experiments (snbt, kZyU).**
We see the contribution of this work in proposing VISA and studying its behavior by exploring a wide spectrum of hyperparameters. While this setting makes it somewhat challenging to evaluate VISA on truly expensive real-world models, we want to highlight that the Lotka-Volterra model and Pickover attractor capture the properties of many real-world scientific simulation models. The Pickover attractor model in particular can exhibit extremely chaotic behavior, is not fully differentiable, and is fairly expensive to evaluate. To be able to track the attractor state we run SMC with 500 particles for 100 time steps for each of the 100 samples that make up a single batch evaluation of the model.

**VISA runs with different optimizers (Bkjs).**
Please find convergence plots in the attached PDF file.

Meta-review: This is a borderline paper.

Given an expensive-to-evaluate $p(y, z)$, not-differentiable wrt $z$ and known $y$, the goal of this paper is to speed up importance-weighted forward-KL variational inference (IWFVI), which optimizes a variational posterior without requiring differentiability wrt $z$.

This work stands on top two previously existing works:
   1. VI with Sample-Average Approximations: These use Reverse KL + SAA. But they require differentiability of $p(y, z)$ wrt $z$.
   2. Importance-weighted forward-KL VI (IWFVI): Does not require differentiability wrt $z$, but requires many samples of $p(y, z)$ to reduce the variance of the gradients.

The novelty of this work is the introduction of a mechanism to speed up IWFVI by reusing the samples of $p(y, z)$ as the optimization of the variational distribution $q(z)$ proceeds. $p(y, z)$ is evaluated at locations sampled from an initial $q(z)$. As we optimize $q$, if it doesn't change much, optimizations steps can proceed without re-evaluating the expensive $p(y, z)$. Whenever $q$ departs too much (above a threshold) from the last sampled one, a resampling is triggered and reevaluation of  $p(y, z)$ is triggered.

This idea accelerates IWFVI, by reducing the number of $p(y, z)$ evaluations required for a given performance. 

On the negative side, this paper only includes synthetic experiments, and I agree with the critical reviewers that this paper would be strengthened by the inclusion of real-world use cases where $p(y, z)$ is very expensive to evaluate. Fortunately, some of the selected synthetic experiments are not trivial and do stress-test the proposed model, so it is reasonable to expect the conclusions of the current work to carry over to at least some real-world use cases.

Given this and the potential value of the novel algorithm, I rule in favor of acceptance.

Authors: T
- Take into account the reviewer's comments to strengthen the paper and improve its clarity. 
- Clarify the use of ESS to define the trust region and include other points of contention raised in the discussion.
- Include a real-world experiment if possible. This can be taken from prior work on Forward-KL Variational Inference, such as Le et al. "Revisiting Reweighted Wake-Sleep". It should be easy to show a reduced number of evaluations to achieve the same performance.

Minor comments:

Under line 81, there's a sign error. Additionally, a $\phi$ is missing in $\hat{w}^{(i)}$.

Line 161 "we can recover by decreasing α" -> Increasing α.
Decision: Accept (poster)