Generalized Bayesian Inference for Scientific Simulators via Amortized Cost Estimation
Authors: Richard Gao
Michael Deistler
Jakob H. Macke
Venue: NeurIPS 2023 poster (2023)
Keywords: simulation-based inference
generalized bayesian inference
neural network
machine learning for science
TL;DR: We train a neural network to estimate a cost function for Generalized Bayesian Inference, improving parameter fitting for (misspecified) scientific simulator models.

Abstract: Simulation-based inference (SBI) enables amortized Bayesian inference for simulators with implicit likelihoods. But when we are primarily interested in the quality of predictive simulations, or when the model cannot exactly reproduce the observed data (i.e., is misspecified), targeting the Bayesian posterior may be overly restrictive. Generalized Bayesian Inference (GBI) aims to robustify inference for (misspecified) simulator models, replacing the likelihood-function with a cost function that evaluates the goodness of parameters relative to data. However, GBI methods generally require running multiple simulations to estimate the cost function at each parameter value during inference, making the approach computationally infeasible for even moderately complex simulators. Here, we propose amortized cost estimation (ACE) for GBI to address this challenge: We train a neural network to approximate the cost function, which we define as the expected distance between simulations produced by a parameter and observed data. The trained network can then be used with MCMC to infer GBI posteriors for any observation without running additional simulations. We show that, on several benchmark tasks, ACE accurately predicts cost and provides predictive simulations that are closer to synthetic observations than other SBI methods, especially for misspecified simulators. Finally, we apply ACE to infer parameters of the Hodgkin-Huxley model given real intracellular recordings from the Allen Cell Types Database. ACE identifies better data-matching parameters while being an order of magnitude more simulation-efficient than a standard SBI method. In summary, ACE combines the strengths of SBI methods and GBI to perform robust and simulation-amortized inference for scientific simulators.

Generalized Bayesian Inference for Scientific Simulators via Amortized Cost Estimation
Richard Gao 21
Michael Deistler
.g..gg.a aig.gga@ggaill com
mi michael.co michael.deistlorre malsstllrrumiituebngeeeee
Jakob H. Macke
jakob. baackkeeuni-uubbbngenn.
Machine Learning in Science, Excellence Cluster Machine Learning, University of Tübingen
Tübingen Al Center
DDepartment Empirical Inference, Max Planck Institute for Intelligent Systems Tübingen, Germany 'qqual contributions.
Abstract
Simulation-base inference (SBI) enables amortized Bayesian inference for sim- ulators with implicit likelihoods. But when are primarily interested in the quality of predictive simulations, or when the model cannot exactly reproduce the observed data (i.e., misspecified), targeting the Bayesian posterior may be overly restrictive. Generalized Bayesian Inference (GBI) aims to robustify infer- ence for (misspecified)) simulator models, replacing the fikelihodd-function with cost function that evaluates the goodness of parameters relative to data. However, GBI methods generally require running multiple simulations to estimate the cost function at each parameter value during inference, making the approach compu- tationally infeasible for even moderately complex simulators. Here, we propose amortized cost estimation (ACH) for GBI to address this challenge: We train a neural network approximate the cost function, which we define the expected distance between simulations produced by parameter and observed data. The trained network can then be used with MCMC to infer GBI posteriors for any obser- vation without running additional simulations. We show that, on several benchmark tasks, ACE accurately predicts cost and provides predictive simulations that are closer synthetic observations than other SBI methods, especially for misspecified simulators. Finally, we apply ACE to infer parameters of the Bodggkin-Hulle model given real intracellular recordings from the Allen Cell Types Database. ACE identifies better data-matching parameters while being an order of magnitude more simulation-efficieer than standard SBI method. In summary, ACE combines the strengths of SBI methods and GBI to perform robust and simulation cotilatioonmmotizeed inference for scientific simulators.
Introduction
Mechanistic models expressed as computer simulators are used in wide range of scientific domains, from astronomy, geophysics, to ncurobiology. The parameters of the simulator, 0, encode mechanisms of interest, and simulating different parameter values produces different outputs, i.c., sim(0,) X,, where each model-simulation X; can be compared to experimentally observed data, x,- Using such simulators, quantitatively reason about the contribution of mechanisms behind experimental
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

measurements. But to do so, a key objective is often to find all those parameter values that can produce simulations consistent with observed data.
One fruitful approach towards this goal is simulation- inference (SBI) 1. which makes it possible to perform Bayesian inference on such models by interpreting simulator outputs as samples
from an implicit likelihood 12]. x ~ p(x 0). Standard Bayesian inference targets the parameter p(x. 0)p(0), where p(0) captures prior p(xo)
posterior distribution given observed data, i..., p(0)x,)
knowledge and constraints over model parameters, and the likelihood function p(x, ()) is evaluated as a function of 0 for fixed x SBI methods can differ in whether they aim to approximate the
likelihood B 4 5 6 or the posterior directly 7 g 91, and can be amortized, i..., do not require new simulations and retraining for new data [00 m.. In the end, each method provides samples from the posterior, which are all, in theory, capable of producing simulations that are identical to the observation we condition on. Furthermore, by definition, the posterior probability of drawing sample scales as the product of its prior probability and, critically, the likelihood that this sample can produce simulation that is exactly equal to the observation.
However, targeting the exact posterior may be overly restrictive. In many inference scenarios, modelers are primarily interested in obtaining a diverse collection of parameter values that can explain the observed data. This desire is also reflected the common usage of posterior predictive checks, where seeing predictive simulations that resemble the data closely (in some specific aspects) is used to gauge the success of the inference process. In particular, it is often clear that the scientific model is only coarse approximation to the data- generating process, and in some cases even cannot generate data-matching simulations, i.e., is misspecified [12]. For example, in the life-sciences, it is uncommon use idealized, theoretically motivated models with few parameters, and it would be unrealistic to expect that they precisely capture observations of highly complex biological systems. In such cases, or in cases where the model is fully deterministic, it is nonsensical use the probability of exactly reproducing the data. In contrast, it would still be useful to find parameter values that produce simulations which are *close enough', or as close possible to the data. Therefore, instead of sampling parameters according to how often they produce simulations that match the data exactly, many use cases call for sampling parameters according to how closely their corresponding simulations reproduce the observed data.
Generalized Bayesian Inference (GBI) 13 of- fers a principled way to do so by replacing the (log) likelihood function with a cost function that
parameter space (9)
data space (x)
simulate: P019)
simply scores a parameter given an observation, such the expected distance between x,, and all possible simulations x produced by 0: (Fig. D. Several recent works have leveraged this frame- work to perform inference for models with implicit or intractable likelihoods, especially to tackle model misspecification: Matsubara et al. [14] use a Stein Discrepancy as the cost function (which requires the evaluation of an unnormalized likeli- hood and multiple i.i.d. data samples), and Cherief- Abdellatif et al. and Dellaporta et al. [15 [6] use simulator samples to estimate maximum mean discrepancy and directly optimize over this cost function via stochastic gradient descent (which requires differentiable simulator). More broadly, have been used to generalize approximate Bayesian computation (ABC) [17] and synthetic likelihood approaches [8 where the Monte Carlo estimate requires (multiple) simulations from p(x(0).
x/
dis
cost function: N[XxK] = E. [d00x
Karget
Figure I:
Estimating cost from simulations. Using the expected distance between simulated and target data as the cost function, GBI assigns high probability to parameter values that, on av- erage. produce simulations that are close-but not necessarily equal to the observation.
cost functions such as scoring rule estimators
Thus, existing GBI approaches for SBI either require many simulations to be run during MCMC sampling of the posterior (similar to classical ABC methods), or are limited to differentiable simula- tors. Moreover, performing inference for new observations requires re-rumning simulations, rendering such methods cimulattoo-inefffcien and expensive at inference-time, and ultimately impractical for scientific simulators with even moderate computational burden.
We here propose to perform GBI for scientific simulators with amortized cost estimation (ACE). which inherits the flexibility of GBI but amortizes the overhead of simulations by training a neural
2

---

network to predict the cost function for any parameter-observation pair. We first outline the GBI formalism in Sec. 2 then introduce ACE in Sec. 3 In Sec. we show that ACE provides GBI posterior predictive simulations that are close to synthetic observations for a variety of benchmark tasks, especially when the simulator misspecified. We showcase its real-world applicability in Sec. [ using experimental data from the Allen Cell Types Database, ACE successfully infers parameters of the Hodgkin-Huxley, single-neuron simulator with superior predictive performance and an order of magnitude higher simulation efficiency compared to neural posterior estimation (10). Finally, we discuss benefits and limitations of GBI and ACE, and related work (Sec. (
Background
To construct the GBI posterior, the likelihood, p(), is replaced by a "generalized likelihood function', L(0; xo), which does not need to be a probabilistic model of the data-generating process, as long as it can be evaluated for any pair of 0 and X,. Following the convention in Bissiri et al. [3. we define T.(0;x, xa) = e -e((0x)), where E(O; x,) is cost function that encodes the quality of 0 relative to an observation x.. and 3 a scalar inverse temperature hyperparameter that controls how much the posterior weighs the cost relative to the prior. Thus, the GBI posterior can be written as
p p())),) X -E((:xx)) p(0).
As noted previously [13] if we define x(0;x.) = - log () (i.e., self-information) and / = 1, then we recover the standard Bayesian posterior, and "tempered" or "power" posterior for 3 f 1 [19 20) The advantage of GBI is that, instead of adhering strictly to the (implicit) likelihood, the user allowed choose arbitrary cost functions to rate the goodness of relative to an observation X.. which is particularly useful when the simulator is misspecified. Previous works have referred to ((0:x,) as risk function (1]) loss function [3] or (proper) scoring rules when they satisfy certain properties 22 [8] (further discussed in Section 6. Here we adopt 'cost' to avoid overloading the terms 'loss' and 'score' in the deep learning context.
3 Amortized Cost Estimation for GBI
3.1 Estimating cost function with neural networks
In this work, we consider cost functions that can be written as an expectation over the likelihood:
x()xx)) = Ex(x1o) Id(x, x.)] = / d(x, xx))(x
Many popular cost functions and scoring rules can be written in this form, including the average mean- geeansquared error (MSE) [13]. the maximum mean discrepancy (MMD²) [23), and the energy score (ES) [22 (details in Appendix A2). While ((0;x, can be estimated via Monte Carlo sampling, doing so in SBI setting is simulation--mefficient and time-intensive, since the inference procedure must be repeated for every observation x., and simulations must be run in real-time during MCMC sampling of the posterior. Furthermore, this does not take advantage of the structure in parameter-space or data-space around neighboring points that have been simulated.
We propose to overcome these limitations by training regression neural network (NN) to learn E(0; x.). Our first insight is that cost functions of the form of Eq. 2 can be estimated from a dataset consisting of pairs of parameters and outputs in particular, from single simulation run per 0 for MSE, and finitely many simulation runs for MMD² and ES. We leverage the well-known property that NN regression converges to the conditional expectation of the data labels given the data: If we compute the distances d(x, x,) between a single observation x, and every x in our dataset, then neural network fo(+) trained to predict the distances given parameters will denoise the noisy distance labels d(x, x,) and converge onto the desired cost fo + E(0: x,,) = [ppxjey]dd [d(x, x.)], approximating the cost of any 0 relative to x. (see Appendix A3.1 for formal statement and proof).
3.2 Amortizing over observations
As outlined above, a regression NN will converge onto the cost function ((0: x.) for a particular observation x,. However, naively applying this procedure would require retraining of the network
3

---

for any new observation Xo. which prevents application of this method in timeccriical or high- throughput scenarios. We therefore propose amortize cost estimation over target distribution p(x+). Specifically, NN which receives as input parameter A and an independently sampled target datapoint x will converge to ((0: x.) for all x, on the support of the target distribution (Fig. 21.b), enabling estimation of the cost function for any pair of (0, x,) (see Appendix A3.2 for formal statement and proof).
Naturally, we use the already simulated x ~ pxx) as target data during training, and therefore do not require further simulations. In order have good accuracy on potentially misspecified observations, however, such datapoints should be within the support of the target distribution. Thus, in practice, we augment this target dataset with noisy simulations to broaden the support of p(xx). Furthermore, if the set of observations (i.e., real data) is known upfront, they can also be appended to the target dataset during training. Lastly, to keep training efficient and avoid quadratic scaling the number of simulations, we randomly subsample a small number of per 0 in each training epoch (2 in our experiments), thus ensuring linear scaling as function of simulation budget. Fig. 2b.b summarizes dataset construction and network training for ACE (details in Appendix A4.1
3.3 Sampling from the generalized posterior
Given a trained cost estimation network fe('s -), an observed datapoint x., and a user-selected inverse temperature 3, the generalized posterior probability (Eq. [T) can be computed up to proportionality for any 0: p(0) ()) & exp( (3 fo xx) (p()) and, thus, this term can be sampled with MCMC
(Fig. 2 ). The entire algorithm summarized in Algorithm [
Algorithm 1: Generalized Bayesian Inference with Amortized Cost Estimation (ACE)
Inputs: prior p(0), simulator with implicit likelihood p(x(0) (0). number of simulations N, feedforward NN fe with parameters Φ, NN learning rate v, distance function d(., -), noise level 0, number of noise-augmmeted samples S, inverse temperature 3, number of target datapoints per 0 Nurges K observations x (I.....)_
Outputs: M samples from generalized posteriors given K observations.
Generate dataset:
sample prior and simulate: 0,x + {0: -=p()). x; ~ p(xx(0;))
add noise and concatenate:
Xiarget (x, X1...s + e, xol. (00,011)
Training:
while not converged do
for (0.x) in batch do
xiedd +- sample Nuargtt datapoints from Xtarget for in .ised do
< L + (0, d(xx x:)) 7 Adam Adam(V.CC) Vv.C)
// and reset zero
Sampling:
for k E [L,..., K] do
Draw M samples, with MCMC, from: exp(-B (0(0 p(0)
3.4
Considerations for choosing the value of
We note that the choice of value for / is an important decision, though this is an issue not only for ACE but GBI methods in general [33. A good "baseline" value for is such that the average distance across a subset of the training data is scaled to be in the same range as the (log) prior probability, both of which can be computed on prior simulations. From there, increasing 3 sacrifices sample diversity for predictive distance, and as 3 approaches infinity, posterior samples converge onto the minimizer of the cost function. In practice, we recommend experimenting with a range of (roughly) log-spaced values since, as we show below, predictive sample quality tend to improve with increasing 3.

---

training data construction
c
inferance
specilled
misspecifed
example 1D simulalor
ak
AGE networ
d(-.x)
0
(f'aa
b
network training
MCMC
d(xx)
wall-specffed X, misspaciled X,
0
Figure 2: Schematic of dataset construction, network training, and inference. (a-b) The neural network is trained to predict the distance between pairs of x (red) and x, (green), as a noisy sample of the cost function (i.e., expected distance) evaluated on 0 (grey) and x, (c) At inference time, the trained ACE network predicts the cost for any parameter 0 given observation x. (top row), which is used to evaluate the GBI posterior under different / (bottom row, darker for larger () for MCMC sampling without running additional simulations. The distance is well-defined and can be approximated even when the simulator is misspecified (dashed lines).
Benchmark experiments
4.1 Experiment setup
Tasks We first evaluated ACE on four benchmark tasks (modified from Lueckmann et al. (24) with a variety of parameter- and nuta-dimeesionnlltyy as well as choice of distance measure: (1) Uniform ID: ID 0 and the simulator implements even polynomial with uniform noise likelihood, uniform prior (Fig. ) (2) 2 Moons: 2D and x, simulator produces a half-circle with constant mean radius and radially uniform noise of constant width, translated as a function of 9, uniform prior; (3) Linear Gaussian: 10D O and x, Gaussian model with mean and fixed covariance, Gaussian prior; (4) Caussian Mixture: 2D 0 and x, simulator returns five i.i.d. samples from mixture of two Gaussians, both with mean 0, and fixed covariances, one with broader covariance than the other, uniform prior.
For the first three tasks, we the mean-squared error between simulation and observation the distance function. For the Gaussian Mixture task, we use maximum mean discrepancy (MMD?) to measure the statistical distance between the four tasks, we can compute the integral in Eq. either analytically or accurately capture it with quadrature over Hence, obtain the true cost [(0; and subsequently the 'ground- GBI posterior (with Eq. [ and use that to draw, for each value of B and x., 5000 samples (GT-GBI, black in Fig. 3. See Appendix A4.2 for more detailed descriptions of tasks and distance functions.
sets of five i.idd. samples. Importantly, for each of
Training data For each task, we simulate 10,000 pairs of (0, x) and construct the target dataset as in Fig. 2a with 100 additional noise-augmented targets and 20 synthetic observations well- specified and 10 misspecified-foo additional prior predictive samples, while misspecified observations are created by moving prior predictive samples outside the boundaries defined by the minimum and maximum of 100,000 prior predictive simulations (e.g.. by successively adding Gaussian noise). See Appendix A4.3 for details.
total of 10120 x, data points. Well-specified observations
Test data To evaluate inference performance, we use ACE to sample approximate GBI posteriors conditioned on 40 different synthetic observations. 20 of which were included in the target dataset x and 10 additional well-specified and misspecified observations which were not included the target dataset. We emphasize that including observations in the target data is not case of test data leakage, but represents real use case where some experimental data which one wants to perform inference on already available, while the network should also be amortized for unseen observations measured after training. Nevertheless, we report in Fig. results for 'unscen' observations, i.c., not in the target dataset. Results are almost identical for those that were in the target dataset (Appendix AT). We drew 5000 posterior samples per observation, for 3 different / values for each task.
We are primarily interested in two aspects of performance: approximate posterior pre- dictive distance and cost estimation accuracy, First, motivated above, want to find parameter
Metrics
5

---

posierior distance well-specified xo
posterior preciclive dislance misspecified
estimation accuracy
GBI accuracy ACE
Uniform
MSE
0.1
proc clod
0.5
C2ST
0.5
0.0
0
N/A
4
20
100
N/A
4
20
100
0.0
0.5
20
100
0.10
Moons
1.0
prodicted
MSE
0.2
C2ST
0.05
0.5
0.00
0.0
0.0
NA
10
100
1000
N/A
10
100
1000
0.0
0.2
100
1000
10
Gauzsion
MSE
20
predicted
C2ST
0.5
0
NA
10
100
NA
10
100
0
5
10
10
100
Gaussian
Mixture
0.10
predicted
0.50
MMD
0.2
C29T
0.05
0.25
0.5
0.00
0.0
0.00
N/A
10
25
75
N/
10
25
75
0.0
0.2
0.4
10
25
75
beta
beta
true cost
bata
GT-GBI
ABC
NLE
-
wall-specifid
AGE-GBI
NPE
NLE (tempered)
misspocffied
Performance benchmark tasks. ACE obtains posterior samples with low average distance to observations, and accurately estimates cost function. Rows: results for each task. Columns: average predictive distance compared to SBI methods and GT and 2nd), cost estimation accuracy evaluated on ACE posterior samples for different 3 (lighter blue shades are lower values of (3) (3rd), and C2ST accuracy relative to GT GBI posterior (4th, lower is better).
Figure 3:
configurations which produce simulations that are as close as possible to the observation, as mea- sured by the task-specific distance function. Therefore, we simulate using each of the 5000 ACE GBI posterior samples, and compute the average distance between predictive simulations and the observation. Mean and standard deviation are shown for well-specified and misspecified observations separately below (Fig. 3 stt and 2nd columns). Second, we want to confirm that ACE accurately approximates (0; ).). which is prerequisite for correctly inferring the GBI posterior. Therefore, we compare the ACE- ACE-predicted and true cost across 5000 samples from each GBI posterior, as well as the classifier 2-sample test (C2ST, 25 24) score between the ACE approximate and ground-truth GBI posterior (Fig. 3 3rd and 4th columns). Note that cost estimation accuracy can be evaluated for parameter values sampled in any way (e.g., from the prior), but here we evaluate accuracy as samples become more concentrated around good parameter values, i.e., from GBI posteriors with increasing 3. We expect that these tasks become increasingly challenging with higher values 3, since these settings require the cost estimation network to be highly accurate in tiny regions of parameter-space.
Other algorithms As a comparison against SBI methods that target the standard Bayesian posterior (but which nevertheless might produce good predictive samples), we also tested approximate Bayesian computation (ABC), neural posterior estimation (NPE, [7]), and neural likelihood estimation (NLE, 4 (6) on the same tasks. NPE and NLE were trained on the same 10,000 simulations, and 5000 approximate posterior samples were obtained for each x,- We used the amortized (single-round) variants of both as fair comparison against ACE, Additionally, to rule out the possibility that ACE is simply benefiting from the additional inverse temperature hyperparameter, we implemented "tempered" version NLE by scaling the NLE-approximated log-likelihood term (by the same (s) during MCMC sampling. For ABC, we used the 10,000 training samples as reference set, from which 50 were drawn as posterior samples with probability scaling inversely with the distance between their corresponding simulation and the observation, i.c., ABC with acceptance kernel [27].

---

4.2 Benchmark results
Overall, we see that for well-specifien x, (i.e., observations for which the simulator specifieedddd ACE obtains GBI posterior samples that achieve low average posterior predictive simulation distance across all four tasks, especially at high values of 3 (Fig. E 1st column). In comparison, ABC is worse for the Linear Gaussian task (which has a higher parameter dimensionality than all other tasks), whereas NPE, NLE, and tempered NLE achieve similarly low posterior predictive distances.
On misspecified observations, across all tasks and simulation-budgets (with the exception of Gaussian mixture on 10k simulations) we see that ACE achieves lower or equally low average posterior predictive simulation distance as neural SBI methods, even at moderate values of 8 (Fig. 3 2nd column, Figs. A4 A5). This is in line with our intuition that ACE returns a valid and accurate cost even if the simulator is incapable of producing data anywhere near the observation, while Bayesian likelihood and posterior probabilities estimated by NLE and NPE are in these cases nonsensical 28 22 30, BT. Furthermore, simply concentrating the approximate Bayesian posterior via tempering does not lead to more competitive performance than ACE on such misspecified observations, and is sometimes even detrimental (e.g., tempered NLF at high 8 on Uniform 1D and Linear Gaussian tasks). Therefore, we see that ACE can perform valid inference for a broad range of simulators, obtaining a distribution of posterior samples with predictive simulations close to observations, and is automatically robust against model-misspeciiiaation as a result of directly targeting the cost function.
For both well-specified and misspecified observations, ACE-GBI samples achieve posterior predictive distance very close to ground-truth (GT)-GBI samples, all values of (Fig. 3] 1st and 2nd column), suggesting that ACE is able to accurately predict the expected distance. Indeed, especially for low to moderate values of 3. the ACE-predicted cost closely matches the true cost (Fig. 3 3rd column, light blue for well-specified x.. Fig. (A2 for misspecified). For higher values of B, ACE- E-prediiceed cost is still similar to true cost, although the error is, as expected, larger for very large 3 (Fig. 3] 3rd column, dark blue).
As result, highly concentrated generalized posteriors are estimated with larger (relative) discrepan- cies, which is reflected in the classifier 2-sample score (C2ST) between ACE and GT GBI posteriors (Fig. 3 4th column): ACE posterior samples are indistinguishable from GT samples low B, even for the 10D Linear Gaussian task, but becomes less accurate with increasing B. Nevertheless, predictive simulation distance dramatically increases with 3 even when ACE is less accurate, suggesting that sampling to explicitly minimize cost function which targets parameters with data-similar simula- tions is productive goal. Relative performance results across algorithms are qualitatively similar when using a training simulation budget of 200 (Fig. A4) and 1000 (Fig. A5), but ABC required a sufficiently high simulation budget and performed poorly for 1000 training simulations or less.
5
Houely inference from Allen Cell Types Database recordings
Finally, we applied ACE to a commonly used scientific simulator and real data: we used single- compartment Hodgkin-Huxley (HH) simulator from neuroscience and aimed to infer eight parameters of the simulator given electrophysiologiaal recordings from the Allen Cell Types Database 33 33 34].
While this simulator can generate a broad range of voltage traces, it is still crude approximation to real neurons: it models only subset of ion channels, it ignores the spatial structure of neurons, and it ignores many intracellular mechanisms (35). It has been demonstrated that parameters of the HH-model given synthetic recordings can be efficiently estimated with standard NPE [00, but estimating parameters given experimental recordings has been challenging (36) and has required ad-hoc changes to the inference procedure (e.g., Bernaerts et al. (37] added noise to the summary statistics, and Gonçalves et al. [10] used a custom multi-round scheme with a particular choice of density estimator). We will demonstrate that ACE can successfully perform simulation-amortized inference given experimental recordings from the Allen Cell Types Database (Fig.
We trained NPE and ACE given 100K prior-sampled simulations (details in Appendix A4.4). After training. ACE accurately predicts the true cost of parameters given experimental observations (Fig. We then used slice sampling to draw samples from the GBI posterior for three different values of (25, 50, 100) and for ten observations from the Allen Cell Types Database. Interestingly, the marginal distributions between NPE and ACE posteriors are very similar, especially for rather low values of 3 (Fig. c. cornerplot in Appendix Fig. A7 p The quality of posterior predictive samples, however, strongly differs between NPE and ACE: across the ten observations from the Allen
7

---

observatins
b
2.5
predicee
ACE predictives 100)
0.5
9Na
80.0
box
3000
0.0
0.0
2.5
true oost
94
15
VV
90
0.4
NPE predictives
0.3
0.0 Qwxx
0.6 0.00 nakes 0.15
0.1
0.0
NPE ACE ACE
ACE
25
50
100
0.0 0.6
35
-F 100
Figure 4: Applicaton of ACE to Allen data. (a) Three observations from the Allen Cell Types Database. (b) True cost (evaluated as Monte-Carlo average over 10 simulations) per vs ACE-
predicted cost. Colors are different observations. (c) Marginals of posterior distributions for NPE (orange) and ACE (shades of blue. Light blue: 5 = 25, medium blue: 3 = 50, dark blue: 3 = 100). (d) Top: Two GBI predictive samples for each observation. Bottom: Two NPE predictive samples. Additional samples in Appendix (88 All (e) Average predictive distance to observation for NPE and
ACE with B {25, 50, 100}.
Cell Types database, only 35.6% of NPE posterior predictives produced more than five spikes (all observations have at least 12 spikes), whereas the ACE posterior predictives closely match the data, even for low values of 3 (Fig. 4d, samples for all observations and all /3 values in Figs. A AAAATO and for NPE in Fig. AIT 66% [B 25). 87% (5 50), 96% (B = 100) of samples have more than five spikes). Indeed, across all ten observations, the average posterior predictive distance of ACE was significantly smaller than that of NPE, and for large values of ß the distance is even less than half (Fig- 44). Finally, for rejection-ABi ABC, only the top samples (out of the full training budget of 100K simulati lauiiostis) had a distance that is less than the average posterior predictive distance achieved by ACE.
To investigate these differences between NPE and ACE, we also evaluated NPE posterior predictive performance on synthetic data (prior predictives) and found that it had an average predictive distance 0.19, which roughly matches the performance of ACE on the experimental observations (0. 174 for 8=50). This suggest that. in line with previous results (37,, NPE indeed struggles with experimental observations. for which the simulator is inevitably imperfect. We then trained NPE with 10 times more simulations (IM in total). With this increased simulation budget, NPE performed significantly better than with 100K simulations, but still produced poorer predictive samples than ACE trained with 100K simulations (for B = (50, 100)), although the marginals were similar between NPE (1M) and ACE (100K) (Fig. A6 samples for all observations in Appendix Fig. A12).
Overall, these results demonstrate that ACE can successfully be applied real-world simulators on which vanilla NPE fails. On the Hodgkin-Huxley simulator, ACE generates samples with improved predictive accuracy despite order of magnitude fewer simulations and despite the marginal distributions being similar to those of NPE.
Discussion
We presented ACE, method to perform distance-aware inference for scientific simulators within the Generalized Bayesian Inference (GBI) framework. Contrary to 'standard' simulation-based inference (SBI), our method does not target the Bayesian posterior, but replaces the likelihood function with cost function. For real-world simulators, doing so can provide practical advantages over standard Bayesian inference:
8

---

First, the likelihood function quantifies the probability that parameter generates data which exactly matches the data. However, cases where the model is rather crude approximation the real system being studied, scientists might well want to include parameters that can generate data that sufficiently close (but not necessarily identical) subsequent analyses. Our method makes this possible, and is advantageous over other GBI-based methods since it is amortized over observations and the inverse temperature 3. Second, many simulators are formulated as noise-free models, and it can be hard to define appropriate stochastic extensions (e.g.. (38). In these cases, the likelihood function ill-defined and, in practice, this setting would require 'standard' SBI methods, whose generally built to model continuous distributions, to model discrete jumps
density estimators in the posterior density. In contrast, our method can systematically and easily deal with noise-free simulators, and in such situations more closely resembles varameter-fitiingg algorithms. Lastly. standard Bayesian inference is challenging when the model is misspecified, and the performance of neural nework-based SBI methods can suffer drastically in this scenario 33)
6.1
Related work
GBI for Approximate Bayesian Computation Several studies have proposed methods that perform GBI on simulators with either an implicit (i.e., comulation-basede likelihood or an unnormalized likelihood. Wilkinson et al. [39] argued that rejection-ABC performs exact inference for modified model (namely, one that appends an additive uniform error) instead of approximate inference for the original model. Furthermore, ABC with arbitrary probabilistic acceptance kernels can also be interpreted having different error models, and Schmon et al. [17] integrate this view to introduce generalized posteriors for ABC, allowing the user to replace the hard-ttreeshold kernel (i.c., e-ball of acceptance) with arbitrary loss function that measures the discrepancy between x and x, for MCMC- -CMC-sampling of the approximate generalized posterior.
Other recent GBI methods require differentiable simulator 16 15 or build tractable cost functions that can be sampled with MCMC [14 18, but this still requires running simulations at inference time (i.e., during MCMC) and does not amortize the cost of simulations and does not reuse already simulated datapoints.
Finally, Bayesian Optimization for Likelihood-free Inference (BOLFI, 40)) and error-guided LFI- MCMC 4T] are not cast as generalized Bayesian inference approaches, but are related ACE, Similarly as in ACE, they train models (for BOLFI, a Gaussian process and, for error-guided LFI- estimate the discrepancy between observation and simulation. BOLFI, the estimator is then used to iteratively select new locations at which to simulate. However, contrary to ACE, neither of these two methods amortizes the cost of simulations over observations.
MCMC, a classifier)
Misppecification-avaai SBI Several other methods have been proposed to overcome the problem of misspecification SBI: For example, Bernaerts et al. (37 add noise to the summary statistics in the training data, Ward et al.
42 use MCMC to make the misspecified data well-specified, and Kelly et al. (43) introduce auxiliary variables to shift the (misspecified) observation towards being well-specified. All of these methods, however, maintain that the inference result should be 'as close as possible' version of the posterior distribution. Contrary to that, our method does not aim to obtain the Bayesian posterior distribution (which, for misspecified models, often be nonsensical even undefined if the evidence is zero), but is specifically targeted towards parameter regions that are a specified distance from the observation. More broadly, recent advances in uncertainty quantification in deep neural networks, where standard preeiccre regression networks are supplemented with uncertainty- pariance-predictin network [44 45 may serve further connect loss-minimizing deep learning with (misspecification-awaal SBI.
6.2 Limitations
While our method amortizes the cost of simulations and of training, it still requires another method to sample from the posterior distribution. We used multi-chain slice-sampling [46) for efficiency, but any other MCMC algorithm, as well variational inference, could also be employed 47 48). While sampling incurs an additional cost, this cost generally small in comparison to potentially expensive simulations.

---

In addition, our method can perform inference for distance functions which can be written as expectations over the likelihood. As we demonstrated, this applies to many popular and widely used distances. Our method can, however, not be applied to arbitrary distance functions (e.g., the minimum distance between simulator samples and the observation). While the distances we investigated here are certainly useful to practioners, they do not necessarily fulfill the criterion of being 'proper' scoring rules [22] [8]. Furthermore, we note that the cost functions considered here by default give rise unnormalized generalized likelihoods. Therefore, depending on whether the user aims to approximate the generalized posterior given the normalized or unnormalized likelihood, different MCMC schemes should be used in conjunction with ACE (standard MCMC vs. doubly-intractabe intractable MCMC, e.g., the Exchange algorithm [49
Compared to 'standard' SBI, GBI introduces an additional hyperparameter to the inference procedure, the inverse temperature 3. This hyperparameter has to be set by the user and its choice strongly affects inference behaviour: low values of 3 will include regions of parameter-space whose data do not necessarily match the observation closely, whereas high values of 3 constrain the parameters to only the best-fitting parameter values. While provide heuristics for selecting (3, we acknowledge the inconvenience of an additional hyperparameter. However, our method amortized over 3, which makes exploration of different 33 values possible, and which could simplify automated methods for setting B, similar to works where 3 is taken as the exponent the likelihood function [50].
Finally, as with any method leveraging deep neural networks, including neural density estimator- based SBI methods (such as NPE), sensitivity the number of training samples and the dimensionality of the task should always be considered. As we demonstrate above, increasing simulation budget improves the performance of any algorithm, and reasonable number of training simulations yielded improved performance on real-world neuroscience application, while the amortization property shifts the cost of simulations up front. addition, we consider tasks up to 10 dimensions here, as most existing SBI methods have been benchmarked and shown to perform adequately on such tasks [24], though it remains to be seen how ACE can extend to higher dimensional parameter-space and data-space and whether embedding networks will be similarly helpful.
7 Conclusion
We presented a method that performs generalized Bayesian inference with amortized cost estimation. Our method produces good predictive samples on several benchmark tasks, especially in the case of misspecified observations, and we showed that it allows amortized parameter estimation of Hodgkin-Huxley models given experimental recordings from the Allen Cell Types Database.
8
Acknowledgements
RG is supported by the European Union's Horizon 2020 research and innovation program under the Marie Sklodowska-Curdt grant agreement No. 101030918 (AutoMIND). MD supported by the International Max Planck Research School for Intelligent Systems (IMPRS-IS), The authors are funded by the Machine Learning Cluster of Excellence, EXC number (366//1390777655. This work was supported by the Tübingen AI Center (Agile Research Funds), the German Federal Ministry of Education and Research (BMBF): Tübingen AI Center, FKZ: OIIS 18099A, and the German Research Foundation (DFG): SFB 1233, Robust Vision: Inference Principles and Neural Mechanisms, project number: 276693517
We would like to thank Jan Boelts, Janne Lappalainen, and Auguste Schulz for feedback on the manuscript, Julius Vetter for feedback and discussion on proper scoring rules, Poornima Ramesh and Mackelab members for extensive discussions throughout the project, as well as Francois-Xavier Briol for suggestions on sampling doubly-intractable posteriors, and the reviewers for their constructive comments on the readability of the manuscript and suggestions for additional analyses.
10

---

References
[1] Kyle Cranmer, Johann Brehmer, and Gilles Louppe. The frontier of simulation-based inference. Proceedings of the National Academy of Sciences, 117(48):30555-3006) 2020.
[2] Peter Diggle and Richard J Gratton. Monte carlo methods of inference for implicit statistical models. Journal of the Royal Statistical Society. Series B, Statistical methodology, 46(2): 227, 1984.
[3] Simon N Wood. Statistical inference for noisy nonlinear ecological dynamic systems. Nature, (1667310):1102-1110 2010.
[4] George Papamakarios, David Sterratt, and lain Murray. Sequential neural likelihood: Fast likelihood- free inference with autoregressive flows. In The 22nd International Conference on Artificial Intelligence and Statistics, pages 837-848. PMLR, 2019.
[5] Joeri Hermans, Volodimir Begy, and Gilles Louppe. Likelihood-fire free meme with amortized approximate ratio estimators. In International Conference on Machine Learning, pages 4239- 4248. PMLR, 2020.
[6] Owen Thomas, Ritabrata Dutta, Jukka Corander, Samuel Kaski, and Michael U Gutmann. Likelihood- frre inference by ratio estimation. Bayesian Analysis, 17(1):1-31, 2022.
[7] George Papamakarios and lain Murray, Fast e-free inference of simulation models with bayesian conditional density estimation. Advances in neural information processing systems, 29, 2016.
[8] David Greenherg, Marcel Nonnenmacher, and Jakob Macke. Automatic posterior transformation for likelihood-fred inference. In International Conference on Machine Learning, pages 2404- 2414. PMLR, 2019.
[9] Katalin Csilléry, Michael GB Blum, Oscar E Gaggiotti, and Olivier François. Approximate bayesian computation (abc) in practice. Trends in ecology & evolution, 25(7):410-418, 2010.
[10] Pedro Gonçalves, Jan-Matthis Lueckmann, Michael Deistler, Marcel Nonnenmacher, Kaan Öcal, Giacomo Bassetto, Chaitanya Chintaluri, William F Podlaski, Sara A Haddad, Tim P Vogels, et al. Training deep neural density estimators to identify mechanistic models of neural dynamics. Elife, 9:e56261, 2020.
[11] Stefan T Radev, Ulf K Mertens, Andreas Voss, Lynton Ardizzone, and Ullrich Köthe. Bayesflow: Learning complex stochastic models with invertible neural networks. IEEE transactions on neural networks and learning systems, 2020.
[12] Stephen G Walker. Bayesian inference with misspecified models. Journal of statistical planning and inference, 143(10):1621-1633 October 2013.
[13] Pier Giovanni Bissiri, Chris Holmes, and Stephen Walker. A general framework for updating belief distributions. Journal of the royal statistical society, series b, statistical methodology. 18(5):1103, 2016.
[14] Takuo Matsubara, Jeremias Knoblauch, François- Xavier Briol, Chris Oates, et al. Robust generalised bayesian inference for intractable likelihoods. arXiv preprint arXiv:2104.07359, 2021.
[15] Badr-Eddine Chericf-Abdellatif and Pierre Alquier. MMD-Bayes: Robust bayesian estimation via maximum mean discrepancy. In Cheng Zhang, Francisco Ruiz, Thang Bui, Adji Bousso Dieng, and Dawen Liang, editors, Proceedings of The 2nd Symposium on Advances in Approxi- mate Bayesian Inference, volume 1 118 of Proceedings of Machine Learning Research, pages 1-21. PMLR, December 2020.
[16] Charita Dellaporta, Jeremias Knoblauch, Theodoros Damoulas, and François- -aanee Robust bayesian inference for simulator-basse models via the mmd posterior bootstrap. In International Conference on Artificial Intelligence and Statistics, pages 943-970. PMLR, 2022.
Xavier Briol.
[17] Sebastian M Schmon, Patrick W Cannon, and Jeremias Knoblauch. Generalized posteriors in approximate bayesian computation. arXiv preprint arX1w:2011.08644, 2020.
11

---

[18] Lorenzo Pacchiardi and Ritabrata Dutta. Generalized bayesian likelihood free inference using scoring rules estimators. arXiv preprint arXiv:2104.01889 2021.
[19] Enzo Marinari and Giorgio Parisi. Simulated tempering: A new monte carlo scheme. May 1992.
[20] N Friel and N Pettitt. Marginal likelihood estimation via power posteriors. Journal of the Royal Statistical Society. Series B, Statistical methodology, 7063):389-607, 2008.
[21] Wenxin Jiang and Martin A Tanner. Gibbs posterior for variable selection in High-Dimensional classification and data mining. Annals of statistics, 36(5):2207-22)) 2008.
[22] Tilmann Gneiting and Adrian E Raftery. Strictly proper scoring rules, prediction, and estimation. Journal of the American statistical Association, 102(477):359-373 2007.
[23] Arthur Gretton, Karsten M Borgwardt, Malte J Rasch, Bernhard Schölkopf, and Alexander Smola. A kernel two-sample test. The Journal of Machine Learning Research, 13(1):723-773, 2012.
[24] Jan-Matthis Lueckmann, Jan Boelts, David Greenberg, Pedro Goncalves, and Jakob Macke. Benchmarking simulation- inference. In International Conference Artificial Intelligence and Statistics, pages 343-351. PMLR, 2021.
[25] David Lopez-Paz and Maxime Oquab. Revisiting classifier two-sample tests. arXiv preprint 1rXiv:1610.06545, 2016.
[26] Mat-Misthis Lueckmann, Giacomo Bassetto, Theofanis Karaletsos, and Jakob H Macke. Likelihood-free inference with emulator networks. In Symposium Advances in Approx- imate Bayesian Inference, pages 32-53. PMLR, 2019.
[27] Scott A Sisson and Yanan Fan. Likelihood-fre markov chain monte carlo. arXiv preprint arXiv:1001.2058. 2010.
[28] Peter Grünwald and John Langford. Suboptimal behavior of bayes and mdl in classification under misspecification. Machine Learning, 66:1 19-119 2007.
[29] Peter Grünwald and Thijs Van Ommen. Inconsistency of bayesian inference for misspecified linear models, and a proposal for repairing it. Bayesian Analysis, 2017.
[30] Patrick Cannon, Daniel Ward, and Sebastian M Schmon. Investigating the impact of model misspecification in neural simulation-basedt inference. arXiv preprint arXiv:2209.00 209001149 2022.
[31] Daniel Ward, Patrick Cannon, Mark Beaumont, Matteo Fasiolo, and Sebastian Schmon. Robust neural posterior estimation and statistical model criticism. Advances in Neural Information Processing Systems, 35:33845-33859, 2022.
[32] Allen Institute for Brain Science. Allen cell types database. http://cellxxpesbbbaaimmm.ppcp.p 2016.
[33] Corinne Teeter, Ramakrishnan Iyer, Vilas Menon. Nathan Gouwens, David Feng, Jim Berg, Aaron Szafer, Nicholas Cain, Hongkui Zeng, Michael Hawrylycz, Christof Koch, and Stefan Mihalas. Generalized leaky integrate-ane and-ate--nn-finn and-fire models classify multiple neuron types. Nature communications, (1): February 2018
[34] Martin Pospischil. Maria foledo-Rodriguce Cyril Monier, Zuzanna Piwkowska, Thierry Bal, Yves Frégnac, Henry Markram, and Alain Destexhe. Minimal hodgkin-huxley type models for different classes of cortical and thalamic neurons. Biological cybernetics, 99:427 441,, 2008
[35] Romain Brette. What is the most realistic single-compartment model of spike initiation? PLoS computational biology, 11(4)cc10041 14, 2015.
[36] Nicholas Tolley, Pedro LC Rodrigues, Alexandre Gramfort, and Stephanie R Jones. Methods and considerations for estimating parameters in biophysically detailed neural models with simulation based inference. bioRxiy, pages 2023-04, 2023.
12

---

[37] Yves Bernaerts, Michael Deistler, Pedro J Goncalves, Jonas Beck, Marcel Stimberg. Federico Scala, Andreas Tolias, Jakob Macke, Dmitry Kobak, and Philipp Berens. Combined statistical-mechansttii modeling links ion channel genes to physiology of cortical neuron types. bioRviv, pages 2023-03, 2023.
[38] Joshua H Goldwyn and Eric Shea-Brown The what and where of adding channel noise to the hodgkin-huxley equations. PLoS computational biology, 7(11):e1002247 2011.
[39] Richard David Wilkinson. Approximate bayesian computation (ABC) gives exact results under the assumption of model error. Statistical applications in genetics and molecular biology, 12(2):129-141) May 2013.
[40] Michael U Gutmann and Jukka Corander. Bayesian optimization for likelihood inference of simulator-based statistical models. Journal of Machine Learning Research, 2016.
[41] Volodimir Begy and Erich Schikuta. Error-guided likelihood-free MCMC. arXiv,, October 2020.
[42] Daniel Ward, Patrick Cannon, Mark Beaumont, Matteo Fasiolo, and Sebastian M Schmon. Robust neural posterior estimation and statistical model criticism. In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho, editors, Advances Neural Information Processing Systems, 2022.
[43] Ryan P Kelly, David J Nott, David T Frazier, David J Warne, and Chris Drovandi. Misspecificcation-rouus sequential neural likelihood. arXiv preprint arXiv:2301. 13368, 2023.
[44] Ian Osband, Zheng Wen, Seyed Mohammad Asghari, Vikranth Dwaracherla, Morteza Ibrahimi, Xiuyuan Lu, and Benjamin Van Roy. Epistemic neural networks. July 2021.
[45] Salem Lahlou, Moksh Jain, Hadi Nekoei, Victor I Butoi, Paul Bertin, Jarrid Rector- Rector-Brooks, Maksym Korablyov, and Yoshua Bengio. DEUP: Direct epistemic uncertainty prediction. October 2022.
[46] Radford M Neal. Slice sampling. The annals of statistics, 31(3):705-767, 2003.
[47] Samuel Wiqvist, Jes Frellsen, and Umberto Picchini. Sequential neural posterior and likelihood approximation. arXiv preprint arXiv:2102.06652, 2021.
[48] Manuel Glöckler, Michael Deistler, and Jakob H. Macke. Variational methods for simulation- based inference. In International Conference on Learning Representattions, 2022.
[49] lain Murray, Zoubin Ghahramani, and David MacKay, MCMC for doubly-intractable distribu- tions. June 2012.
[50] Pei-Shien Wu and Ryan Martin. A comparison of learning rate selection methods in generalized bayesian inference. Bayesian Analysis, 18(1): 2023.
[51] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al. Pytorch: An imperative style, high-performance deep learning library. Advances in neural information processing systems, 32, 2019.
[52] Omry Yadan. Hydra framework for elegantly configuring complex applications. Github, 2019.
[53] Scott A Sisson, Yanan Fan, and Mark M Tanaka. Sequential monte carlo without likelihoods. 2007.
Proceedings of the National Academy of Sciences,
[54] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Deep residual learning for image recognition. In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 770-778, 2016.
[55] Manzil Zaheer, Satwik Kottur, Siamak Ravanbakhsh, Barnabas Poczos, Russ R Salakhutdinov, and Alexander Smola. Deep sets. Advances in neural information processing systems, 30, 2017.
13

---

[56] Alvaro Tejero- Caneter, Jan Boelts, Michael Deistler, Jan- Matthis Lueckmann, Conor Durkan, Pedro J. Gonçalves, David S. Greenberg, and Jakob H. Macke. sbi: A toolkit for simulation- based inference. Journal of Open Source Software, 5(52):2505, 2020.
[57] Conor Durkan, Artur Bekasov, lain Murray, and George Papamakarios, Neural spline flows. Advances in neural information processing systems, 32, 2019.
[58] Michael Deistler, Pedro J. Goncalves, and Jakob H. Macke. Truncated proposals for scalable and hassle-free simulation- based inference. In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho, editors, Advances in Neural Information Processing Systems, 2022.
14

---

Appendix
A1 Software and Data
We used PyTorch for all neural networks [5T] and hydra to track all configurations reproduce results is available at /tttts:://gthhhe www.coo/cccaa/// comm/mccceelabbbnuurlllg..
52
Code to
A2 Common cost functions as expectations of the likelihood
As described Sec. 3 our framework includes all cost functions which can written as expectations over the likelihood: 2 x.) [p(x [d(xx, x.)]. Below, we demonstrate that the average mean- squared error (MSE), the maximum mean discrepancy (MMD²) 23 and the energy score [22 fall into this category.
Average Mean-squared error The average mean trann-quared error between samples from the likeli- hood p(x () and an observation x, is defined as:
MSE(p(x(O), x.) = Ex Ex-p(xxo)] - xo.1²]
which is trivially of the form required for our framework.
Maximum mean discrepancy Throughout this paper, the maximum mean discrepancy squared (MMD²) [23). The MMD between the likelihood p(x(0) and the data distribution p(x.) is defined as:
MMD"(p(x(00,,oo p(xx))
E Ex,x'<p(x10) ( x')] xx wp(x.) (k(x:: x()] (xxxp(xkl))xx.xp(xx. (K(x, xo)]
Assume we are given K iid samples as observations x!" ().... Thus, the MMDD becomes:
MMD p(x.)))))()))))) p(xx))
K
K
1 K
Ex bx +
k(x), x²) + Ex-
Σ K
k(x,
K²
==]
which can be written single expectation over p MMD3(p(x10),
1 K
E Exxp(x)0) F [xx,,
k(x
x +
Σ x(xxxxx xx)
k
1=1
and follows the form required for our framework.
Energy Following Gneiting al. (22), the (negative, in order to fit notation) energy score (ES) is defined as
(S(p(x()\xx.) -x x'Il Ex 0) |||x xx11%1.
This can straight-forwardly be written single expectation over the likelihood
ES(p(x(0))x.) Exxx(xx0) Fx x'11a] ||x x=
and, thus, falls within our framework.
A3 Convergence proofs
The proofs closely follow standard proofs that regression converges the conditional expectation.
15

---

A3.1 Proposition and proof
Let p(O.x x) be the joint distribution over parameters and data. Let ((0:x.) be cost function that can be written as ((0;xo) = Ep(x(o) [d(x, x.)) [x)))x)))) f be function parameterized by 6. Then, the loss function C = (p(e (f - d(x, x,))a] is minimized if and only if. for all 0 E supp(p(O)), fa(0) = [pixio) x.)].
Proposition 1.
= Ix d(x, xx)p(x 0) and let
Proof. We aim to prove that
xx xe) [d((x, Xo) - ((0)) 2 xxxx=x(xx))(((( xo) - Ex p(x)))x))) xo)1)²)
for every function g0)). We begin by rearranging terms:
10,x-p(0 [d( xo) - 9(0))²]=
_0,x~p(0xx) [[d x,) (x [d(x', xx)x + Ex oxx xo)] -
500xxp(00xx [d((x, xx) =px100) [d(x', xx) 2 + (Exx xx)] ((0))²
with
X = B _xxx~p(0xxi (d (x, x.) Ex ~p(x0) [d(x', X., x)) x xp(x [d(x' x,)] ,()))
By the law of iterated expectations, one can show that X 0: X Eg Ee-p(o) p(x [(d(x, x.) p(xx)( [d(x', xx (Bxx [x)x -(())
The first term in the products reads a difference of the same term, and, thus, is zero
Thus, since X = 0 and since Ex' ~p( [d(x', x-)] 9(0) O, we have
Eg xx ((dx) [([(x, x.) -(0))²² 2 E. (xxx [d((x, xx) Ex xx [d(x', xx)))
A3.2 Proposition and proof
Proposition 2. Let p(0, x) be the joint distribution over parameters and data and let p(x;) be a distribution of larget samples. Let e(0;x, x.) be cost function and (-)-) parameterized function as
in proposition [ Then, the loss function L = xp((x))(xx)) [. xt) dd(x, x1)) is minimized if and only if. for all O E supp(p(O)) and all xx E supp(p(xx) (xx)) we have (0, xx) = [Ep(x00] [d(x, x:)]
Proof. We aim to prove that
Eg Ex p(xx)) [(((x 9(0,x,)) >
Es Exxxpie x p(x,) [[((x, x,)
Ex (pxx)))))) [d(x', x. )x)
for every function 0 x.). We begin as in proposition]
Eg xxxx=p(0xx) x) xxxx=pppxxxxxxxxxxx)())) x(xx) [dd(x, x.)
g(0, x:))
E Ep(x.) [[0x ((d(xx - x.))
Below, we prove that, for any Xt, the optimal g(0, x+) is the conditional expectation Exixo(xio) [d(x', x:)]:
Ex (d(xxxxx)) 9 x1
Ex [d((x, x.) Ex 0) id(x', xx)] + Exx xx)] g(0, x1))
Eg Baxxpp(0x0 [d((x x:) - Excxpix 0) i(x'' x:)))² + (Ex ld(x', x:)] g (,, xx + X
with
X= Eg [[d(x, xx) Ex ~n(())) [d(x', (Bx (Bx) [d(x xx)] - g(0, x, x1)
which, as above, is X = 0 (proof is identical to proposition D. Thus, since Ex'<p(x(0) [d(x', xx)]
g(0, x2 > 0, we have:
Eg,x-p(0xx) [[d(x, x.) g(O, (10x 2 Eg xxxxpp((xx) [d((x x.) [x)) [d(x', [d(x'xx))))' {1
Because this inequality holds for any xt, the average over p(x,) will also be minimized if and only
g(0.x:)) matches the conditional expectation Exr=pix [d(x', x.)] for any x, within the support of p(xx).
16

---

A4 Further experimental details
A4.1 Details on training procedure
Dataset construction We generate samples from the prior 0; ~ p(0) and corresponding simulations
X; ~ (((00)). ) to obtain paramoter-simulatoon pairs, θ, X = (0,,,x):=1...h)) Next, a dataset of larget data points, Xarget, is constructed by concatenating in the batch-dimension: 1) X, 2) random subset of X augmented with Gaussian noise with a specified variance, i..., X; + C,, C; NN0, I2), and 3) any number of real experimental observations the user wishes to include, Xeat:: total, Nurget N + Nuoised + Nreat- Note that Xaaget does not technically need to include any simulated or real data, only "realistic" targets in the neighhorhood of the observed data that one eventually wants to perform inference for. Practically, 1) has already been simulated, and 2) and 3) does not require additional simulations while improving performance through noise augmentation and having access to real data targets. Finally, a pairwise distance matrix D is computed with elements d,, d(x,,x/) X2 E X, X₁ E Xumga- D can either be pre-computed in full, or partially within the training loop.
Network optimization and convergence to GBI loss Given dataset θ X, Xarget D, fully connected feed-forward deep neural network with weights & is trained minimize the mean squared other words, for a parameter configuration 0; and a target data point X,, the network is trained to predict the distance di,, between the corresponding single (stochastic) simulation x; and the target xN. In every training epoch, Daagee target simulations (usually 1-10, compared to Naaget > 10000) are randomly sub-sampled per O,, drastically reducing training time and allowing the relevant di,, to be computed on the fly.
error loss: Nampes Nan Eeer (NN, (0., x.) d(x,, x,))².
Note that we do not wish to accurately predict the distance d(x;;x, x,) in data-space for individual simulations X,, but rather the loss X.) Conveniently, with mean squared error the objective function for network training, for fixed pair of i,, x,, the trained network predicts the mean distance, Eptxle)) [d(x, x,)], precisely the loss function target in Eq 2 (proof in Appendix).
A4.2
Description of benchmark tasks, distance and cost functions
A4.2.1
Uniform ID
A one-dimensional task with uniform noise to illustrate how expected distance from observation can be used as cost function for inference, especially in the case of model misspecification:
Prior
u(-1.5,,.15
Simulator
x = g(z) + e, where U(-0.25, 0.25), 0.8 (0 + 0.25), and
Ppoiseless g(2) 0.1627 0.90732 - 1.21972 - 1.4639a3 + 1.4381
Dimensionality
0 E 11 x E IR
Cost function
d(x, x.): MSE; p() computed exactly given uniform noise. For obtaining the true cost, Eq. is analytically integrated over Tnoiseless + 0.25
Posterior
Ground GBI posterior samples are obtained via rejection sampling. We used the prior as proposal distribution.
A4.2.2 Two Moons
two-dimensional task with a posterior that exhibits both global (bimodality) and local (crescent shape) structure to illustrate how algorithms deal with multimodality:
17

---

Prior
U(-1,1)
r cos((a) + 0.25
-|0: + 021l/22
Simulator
æ/0
+
where a ~ -((-1/2,7/2)
sin((a)
[--1 + 02)/v2
and (0001,00001))
Dimensionality
0 E E R
Cost function
d(x x,): MSE; and p(x computed exactly. To obtain the true cost, Eq. is numerically integrated over a 2D grid of dx with 500 equal bins in both dimensions, where r(1) E [-11.2, 0.4] (2) E 1-1.6, 1.6].
Posterior
Ground-truth GBI posterior samples are obtained via rejection sampling with the prior as proposal distribution.
References
24 g
A4.2.3 Linear Gaussian
Inference of the mean of a 10-d Gaussian model, in which the covariance is fixed. The (conjugate) prior is Gaussian:
Prior
00,00.10)
Simulator
x ~ ((xmmm =0
Dimensionality
0 E 10 1
Cost function
d(x x.): MSE; and p(x computed exactly. The true cost (Eq. be computed analytically.
Posterior
Ground-truth GBI posterior samples are obtained following the procedure in Lueckmann al. [24]: We first ran MCMC to generate 10k samples from the GBI posterior. We then trained neural spline flow on these 10k samples and, finally, used the trained flow as proposal distribution for rejection sampling.
References
24
A4.2.4
Gaussian Mixture
The single-trial version of this task is common in the ABC literature. It consists of inferring the common mean of mixture of two wwo-dimensional Gaussian distributions, one with much broader covariance than the other. this study, we used this task infer the common mean of the two distributions from five i.i.d. simulations:
Prior
u1-10,
Simulator
x 0.5 (x/mma 0,S 1)+ 0.5 (x 0,S 0.01 I)
Dimensionality
O R. E R
Cost function
d(x, x.): MMD The true cost (Eq. 2 is obtained by integrating the distance function grid for every trial independently and multiplying over the trials.
Posterior
Ground-truth GBI posterior samples obtained via rejection sampling. As proposal, we used Normal distribution centered around the ground truth parameter and with variance 50
3
References
24 53
18

---

A4.3 Training, inference, and evaluation for benchmark tasks
For the benchmark tasks, the cost estimation network is a residual network with 3 hidden layers of 64 units [54]. The training dataset was split 90:10 into training and validation sets, with Pugget = 2 and 5 x, randomly sampled each epoch for evaluating training and validation loss, respectively. We used batchsize of 500, i.c., 500 0. 2 xg. for 1000 cost targets per training batch. Networks usually converge within 500 epochs, with 100 epochs of non-decreasing validation loss as the convergence criterion. For inference, we ran multi-chain slice sampling with 100 chains sample the potential function. For the Gaussian Mixture task with 5 i.d.d. samples per simulation/observan we appended to the cost estimation network fully connected, permutation-invarian embedding network (55) (2 layers. 100 units each), which preprocessed the i.d.d. datapoints into a single 20-dimensional vector (but was trained end-to-end).
ACE
Kernel ABC We compared our algorithm with version of kernel ABC on the benchmark tasks. For kernel ABC, we accepted samples (from the fixed set of 10000 prior simulations) with probability exp(-B (d)--,d((x,, xo)). where d(.,-) is the distance function. In many cases, and especially for large B, this yielded very few or no accepted samples. We, therefore, reduced / until at least 50 samples were accepted.
We used the implementation in the shi toolbox [56], with neural spline flow [57 as density estimator (five transformation layers) to approximate the posterior. For tasks with bounded priors (all except Linear Gaussian), we appended sigmoidal transformation to the flow such that its support matches the support of the prior [58] For the Gaussian Mixture task, the same permutation-invariat embedding net was used as for ACE. Posterior samples were directly obtained from the trained flow. All other hyperparameters were the same as in the shi toolbox, version 09. [56].
NPE
NLE and tempered NLE We used the implementation in the shi toolbox [56]. with neural spline flow 57 as density estimator (five transformation layers) to approximate the likelihood. For the Gaussian Mixture task with 5 i.d.d. simulations per samplc/observation, they were split up 5 independent simulations with the same repeated 5 times, effectively having times the training simulation budget. We ran multi-chain slice sampling with 100 chains to sample the potential function as the sum of log-prior probability and approxppo log-likelihood. Tempered NLE uses the same learned likelihood estimator, but the likeliidddo term in the potential function is scaled by B during MCMC sampling (same values as used for ACE). All other byperparameters were the same as in the sbi toolbox, version 0.19.2 [56].
Noise augmentation For each task, randomly subsampled 100 of all simulated datapoints (x), and added Gaussian noise with zero-mean and standard deviation equaling two times the standard deviation of all prior predictives, i.e., /O, (2). We additionally varied noise amplitude (0, 20,, 50,) and found little to no effect on benchmark performance (see Appendix A3).
Synthetic misspecified observations For the Uniform ID, 2 Moons, and Linear Gaussian task, we generated an additional 100,000 prior predictive simulations and defined the prior predictive bound as the minimum and maximum (of each dimension) of those simulations and a₂ as their standard deviation. Then, to create 20 misspecified observations, generated 20 model-simulations and added Gaussian noise with zero-mean and standard deviation equaling 0.5 a. (i.c., ~ (O (0.50,,) (1)) to the 20 observations iteratively until they are outside of the prior predictive bounds in every dimension. For the Gaussian Mixture task, we replaced the second Gaussian above by (11.5. x sign(0), 0.5²I), ie. displacing the mean the corner of the quadrant whose signs match 0, and slightly inereasing the variance.
Metrics For the metrics in Fig. 3. in columns 1, 2, and 4 and all corresponding figures in the Appendix, results were aggregated across the 10 well rell-specified and misspecified samples separately, where marker and error bars represent mean and standard deviation over 10 observations. Columns I and 2 show average distance between observation and 5000 posterior predictive simulations from each algorithm (except ABC, for which there were 50 posterior samples). Column 3 shows, for each of the 10 well-specified observations and 3 B values, 3 random posterior samples for which we compare the ACE-estimated and ground-truth cost, i.c., 10 x 3 x 3 = 90 points are shown for each
19

---

task. Column 4 shows C2ST score between 5000 ground-truth GBI posterior samples and ACE GBI posterior samples.
A4.4 Training procedure for Hodgkin-Huxxey model
For NPE, we used the implementation in the sbi toolbox [56]. As density estimator, we used a neural spline flow [57] with five transformation layers. We used batchsize of 5,000, We appended sigmoidal transformation to the flow such that its support matches the support of the prior [58) All other hyperparameters were the same in the shi toolbox, version 0.19.2 [56].
For ACE, we used a residual neural network [54] with 7 layers and 100 hidden units each. We used 10% of all simulations as held-out validation set and stopped training when the validation set did not decrease for 100 epochs. We used batchsize of 5000. For sampling from the generalized posterior, multi-chain slice sampling with 100 chains. As distance function, error, where the distance in each summary statistic was z-scored with the standard deviation of the prior predictives. To generate Xs, we used all simulations which generated between 5 and 40 spikes (since any reasonable experimental recording would fall within this range) and we appended 1000 augmented simulations which we added Gaussian noise with times the standard deviation of prior predictives that produce between and 40 spikes. We did not use the experimental recordings from the Allen Cell Types database in X,.
used mean-absolute-
A5 Hodgkin-Huxley model
We used the same model as Gonçalves et al. [1] which follows the model proposed in Pospischil et al. (34). Briefly, the model contains four types conductances (sodium, delayed-rectiffer potassium, slow voltage-dependent potassium, and leak) and has total of eight parameters that generate time series which we reduce to seven summary statisties (spike count, mean resting potential, standard deviation of the resting potential, and the first four voltage moments mean, standard deviation, skew, curtosis, same as in Gonçalves et al. 10
20

---

A6 Supplementary figures
possarior predictive distance well-specified xo
posterior predictive distance miespecified xo
cost matioon accuracy
postarior accuracy ACE vs. GT
Uniform 1D
1.0
0.10
predicted
0.05
0.5
100
Moons
1.0
predicted
MSE
C2ST
0.5
0.5
0.00
0.0
0.0
N/A
10
100
1000
N/A
10
100
1000
0.0
0.2
10
100
1000
Linear
MSE
10
predicted
C2ST
2
0.5
100
100
Gaussian
MWD
0.10
0.4
predicted
0.2
C2ST
0.2
0.00
0.0
0.0
0.5
N/A
10
25
75
N/A
10
25
75
0.0
0.2
10
25
75
bata
GT-GBI
ACE-CBI
NLE
well-specifed
misspectfied
Figure A1: Benchmark performance for 'seen' observations Panels are the same as in Fig. but
for the 20 To that were included in the target dataset during training.
21

---

well-specifid,
well-specifidd,
misspacified, unseen
misspecifidd,
Uniform 1D
prodicted
0.2
0.2
0.0
0.0
0.0
0.0
0.0
0.5
0.0
0.2
0.0
0.5
0.0
Moons
predicted
0.2
0.2
0.2
0.2
0.0
0.0
0.0
0.0
0.0
0.2
0.0
0.2
0.0
0.2
0.0
0.2
10
10
predicted
Linear
5
0
Gaussian
Mixture
predicted
0.50
0.50
0.2
0.25
0.25
0.00
0.0
0.00
0.0
0.0
0.2
0.0
0.2
0.0
0.2
0.4
0.0
0.2
Figure A2: Cost estimation accuracy for all observations Columns are same as 3rd column in Fig. 3 but for all combinations of (unseen, seen) and (well-specified, misspecified) observations (10 each, 40 total) for each task. 1st column here is identical to Fig. 33 3rd column.
uniform two moons
linear gaussian
gaussian mixture
vell-spacifie xo
0.10
0.10
a=5
0.1
0.05
0.05
0=0
no augm.
0.0
0.00
0.00
4
20
100
10
100
1000
10
100
10
25
75
misspecified Xo
0.5
0.2
10
0.2
0.0
0.0
0.0
20
100
100
1000
10
100
10
25
75
beta
beta
beta
beta
Figure A3: Benchmark performance for ACE trained with various augmentation noise levels. Results based on 10k training budget plus 100 noise augmented samples with varying noise amplitudes (blue, orange, green), well as removing data augmentation altogether (i.c., noise augmented simulated nor observed data, red). Overall, changing o does not significantly impact performance compared to the original results (2az, blue).
22

---

posserior prediclive distance well-specified xo
posteriar predictive distance misspecified xo
cosl mstiion accuracy
GBI postarior accuracy ACE vs. GT
Uniform 1D
MSE
0.2
predicted
0.4
C2ST
0.2
0.5
100
20
0.0
0.5
20
100
Moons
MSE
0.2
predicted
0.2
C2ST
0.0
0.0
0.5
NAA
100
1000
N/A
10
100
1000
0.0
0.1
0.2
100
1000
5.0
5.0
Linear
MSE
predicted
C2ST
2.5
2.5
0.0
0.0
NA
10
100
N/A
1
10
100
0.0
2.5
5.0
10
100
1.0
1.0
Caussian Mixture
predicted
0.4
MM00
0.5
0.5
0.2
C2ST
0.5
0.0
beta
beta
true cost
beta
GT-GBI
ACE GBI
ABC
NPE
NLE
well specified
-D-
misspecified
Figure A4:
Benchmark performance with simulation budget of 200 Panels are the same as in but all algorithms trained with simulation budget of 200.
Fig. 3
postarior predictive distance well-specified xo
posterlor predictive distance misspecified xo
ccst est mation accurasy
GBI postarior accuracy ACE. vs. GT
Uniform 1D
gredicted
MSE
0.1
0.5
C2ST
0.5
100
Moons
predicted
MSE
C2ST
0.05
0.1
0.00
0.0
0.5
NAA
10
100
1000
N/A
10
100
1000
0.0
0.2
10
100
1000
Linear
MSE
10
predicted
C2ST
N/
10
100
N/A
10
100
5
100
Caussian Mixture
0.50
Redicied
0.50
MMD
C2ST
0.1
0.25
0.5
0.0
0.00
0.00
beta
bets
true cost
beta
GT-GBI
ACE-
GBI
ABC
NPE
NLE
well- specfepedifeed
-B- misspecified
Figure A5: Benchmark performance with simulation budget of 1000 Panels are the same as in Fig. 3 but all algorithms trained with simulation budget of 1000.
23

---

observations
0.4
pradictive distance
0.3
0.2
0.5 She
80.0
50
3000
NPE (1M sims) prediciives
0.1
0.3
100
9k
40
-V4
0.0 9w<
0.00
0.15
9m
35 -Enu
100
Figure A6: NPE performance with 1 million simulations. Panels are the same as in Fig. but NPE was run with million simulations. ACE is still run with 100K simulations and, thus, the ACE data the same as Fig.
0.5
9
0.0001
9k
0.0001
0.6
3000
Arow
40
90
-V
0.0001
0.16
-E
Figure A7: Cornerplot of posterior distributions. Diagonals ID-marginals, upper diagonals are 2D-marginals (e8th-percentile contour from 5000 posterior samples). Orange: NPE with 100K simulations. Shades of blue: ACE with 100K simulations with 25 (light blue), 50 (medium blue), 8 = 100 (dark blue).
24

---

ACE (B = 25)
Figure A8: Posterior predictive samples of ACE (with 100K simulations) with B = 25. Top row (black): 10 experimental recordings from the Allen Cell Types database. Below: nine predictive samples given each of the ten observations.
25

---

ACE (B 50)
Figure A9: Posterior predictive samples of ACE (with 100K simulations) with B = 50. Top row (black): 10 experimental recordings from the Allen Cell Types database. Below: Nine predictive samples given each of the ten observations.
26

---

ACE (B = 100)
Figure A10: Posterior predictive samples of ACE (with 100K simulations) with 3 100. Top row (black): 10 experimental recordings from the Allen Cell Types database. Below: Nine predictive samples given each of the ten observations.
27

---

NPE (100K simulations)
Figure All: Posterior predictive samples of NPE with 100K simulations. Top row (black); 10 experimental recordings from the Allen Cell Types database. Below: Nine predictive samples given each of the ten observations.
28

---

NPE (1M simulations)
Figure A12: Posterior predictive samples of NPE with 1M simulations. Top row (black): 10 experimental recordings from the Allen Cell Types database. Below: Nine predictive samples given each of the ten observations.
29


Reviews:
[{"review_id": "vKq8gRzuse", "reviewer": "Reviewer_bf2N", "summary": "The paper presents a new technique - amortized cost estimation (ACE) - that, as stated on the can, amortizes a broad class of loss functions used in generalized Bayesian inference (GBI) in place of the (log) likelihood. \nAfter training on a moderate-to-large number of model simulations (10K-100K in the examples of the paper), the amortized loss can then be used in place of the real loss to perform GBI (e.g., via MCMC). Crucially, this step does not require further simulations from the model.\nThe authors show that the amortized loss generally matches the ground-truth loss, and show that their amortized method inherits the advantages of GBI, in that it is typically more robust to misspecification than standard and amortized Bayesian inference, as in neural likelihood estimation and neural posterior estimation", "questions": "- Can you expand on the role of $\\beta$, e.g. how you chose (and a practitioner should choose) reasonable ranges of the hyperparameter?\n- What about having an inverse temperature to the log likelihood of standard Bayesian inference (as per standard tempering, also seen helping in Bayesian deep learning applications)? Can you show that in your examples (for NLE, and possibly NPE)?\n- Consider adding a lesion study or analysis about the role of $\\sigma$ and $S$.\n- Not strictly necessary, but it would be nice to see the method at work on another (simple) example with real data.", "limitations": "The authors address the limitations of their method (and I already asked above to expand on the role of $\\beta$).\nThe work has no particular potential negative social impact (at least not above any other general method for statistical inference).", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- **Quality:** The general polish and quality of the paper is high.\n- **Clarity:** The paper is well written and generally clear. The aim and methodology are very well explained and motivated, with some points that would benefit from further expansion.\n- **Significance:** The idea of the paper (amortizing GBI for simulator-based inference, with the additional goal of dealing with misspecification) is very timely and there are surely many applications.\n- **Originality:** The idea at the core is not particularly original - amortizing the loss in GBI seems a natural direction for the field at the moment - but clearly worth pursuing. Similarly, the execution seems fairly straightforward (which is not a bad thing).\n\n\n### Post-rebuttal:\n\nThanks for having addressed most of my comments. I appreciated the additional experiments with tempered NLE, with varying noise augmentation, and the recommendations for the choice of $\\beta$. The proposed heuristic seems both reasonable and practical. I understand that given the limited time, it was not feasible to conduct an additional set of experiments.\n\nOverall, I am satisfied with the changes and I increased my score accordingly from 6 to 7.", "weaknesses": "Most of my concerns have to do with the treatment of the hyperparameter $\\beta$ (inverse temperature), a key element of GBI.\n\nThe role of $\\beta$ could be explained and discussed more (it is mentioned in the **Limitations**).\n\nMany experiments show the results for a range of $\\beta$ (which also changes from experiment to experiment), and it is unclear how the practitioner should choose its value. For example, $\\beta \\rightarrow \\infty$, (generalized) Bayesian inference tends to maximum-a-posteriori (or minimum-generalized-loss) point estimation. The paper is ambiguous on what the user should be trying to achieve. I appreciate that some of these issues are not of ACE (the authors' method) but of GBI, but the paper would benefit from directly addressing these points.\n\nThe comparison between GBI and Bayesian inference is somewhat unfair in that Bayesian inference is not allowed an inverse temperature hyperparameter (i.e., likelihood tempering), but in principle it could be easily applied (at least for NLE). Is ACE+GBI truly better here, or is the advantage just given by the fact that ACE/GBI has an extra free hyperparameter? In some cases, standard Bayesian inference can do better with large $\\beta$. We know this can be the case for example with neural network posteriors (possibly a case of prior misspecification).\n- With NLE, just run MCMC with a scaling factor $\\beta$ on the amortized log likelihood. \n- Amortizing different temperatures with NPE would be more complex. Naively, one would need to retrain the network for different values of $\\beta$. One could also subtract the log prior from the NPE posterior to get the log likelihood, and then rescale it, but that might lead to instabilities.\n\nSimilarly to $\\beta$ (but less importantly), the paper could explore more the role of the noise $\\sigma$ added to $S$ observations. This parameter seems to be less relevant though - it is very nice that the authors use a fixed value throughout the paper, and this value can be determined from prior-predictive checks. Nonetheless, a bit of exploration / insight could help.\n\nFor the rest, the empirical evaluation is acceptable, but it would have been nice to see at least a couple of applications to real-world data (even another simple one), since this is arguably (also according to the authors) where GBI shines.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate the reviewer’s encouraging remarks regarding the timeliness and broad applicability of our work in leveraging GBI for (misspecified) SBI problems, as well as noting the high quality of execution and presentation in our paper. Furthermore, they raised several interesting questions and comments regarding the effect of the beta and sigma hyperparameters in our original results, for which we conducted additional experiments and discuss below.\n\n**Experiments with tempered NLE**:\nWe very much agree with the reviewer that the standard Bayesian posterior from neural SBI methods (e.g., NLE) could in theory benefit from a temperature hyperparameter, and therefore could constitute a more fair experiment vs. ACE. Therefore, we implemented “tempered NLE” by adding a similar beta parameter (inverse temperature) to the log-likelihood term for MCMC sampling, i.e., beta*NLE_loglikelihood + logprior. \n\nIn short, we find that increasing beta for tempered NLE does not improve posterior predictive distance as it does for GT/ACE-GBI, and is more often detrimental. Results on all tasks are presented in Fig. 1 of the attached pdf (10k training budget). The one exception is the 2 Moons misspecified observation task, where higher beta improves NLE, but ACE still systematically outperforms tempered NLE at all betas.\n\nIt’s an open question whether tempering NLE is detrimental due to targeting of entirely different objectives (likelihood vs. predictive distance), or if tempering exacerbates errors in the normalizing flow’s approximation of the likelihood. We will include the new results and discussions for potential future investigations in the paper, and we thank the reviewer for this interesting suggestion.\n\n**Experiments with varying noise augmentation (sigma)**:\nSimilarly, we agree that the role of sigma could be of interest to explore further. We therefore conducted additional experiments by varying our data augmentation procedure under a simulation budget of 10k: we trained ACE with sigma of 0, 2, and 5 (original results with sigma=2), as well as entirely removing data augmentation, i.e., no noised simulations nor real observations seen during training. \n\nOverall, varying noise augmentation during training has barely any effect on predictive distance across all tasks (Fig. 2 in attached pdf). This may be due to the fact that the number of augmented samples (100) is small compared to the simulation budget (10k), providing reassurance regarding the robustness of our main results. We also studied the effect of removing training augmentation with a training budget of 1000 (not shown), which slightly decreased performance in two tasks where the observation is misspecified. This supports our original motivation for noise and real data augmentation as expanding prior simulation range to combat misspecification, though its effect is more pronounced with smaller training simulation budgets.\n\nAs the reviewer pointed out, we kept sigma constant throughout our original experiments, and it is straightforward to determine a good value from prior predictive checks. These new experiments nevertheless provide interesting insights and sanity checks, which will be included in the results and discussion sections of the paper, and we thank the reviewer for their suggestions.\n\n**Discussion of beta and guideline for practitioners**:\nWe agree that the non-trivial role of beta can be mentioned earlier in the work, and will discuss it in more detail in the methods section in the context of GBI. As the reviewer pointed out, while the issue is not unique to ACE but GBI in general, using our proposed method nevertheless requires the practitioner to choose a beta value suitable for their goals:\n\nTaking the log of Eq. 1, we see that beta weighs the cost function against the impact of the log-prior regularization. In this particular case, a Monte Carlo estimate of the cost function (i.e., distances) is a quantity that can be computed on simulated data alone (unlike log-likelihood in e.g., tempered NLE). As such, a practitioner should choose a value of beta that scales the cost function relative to the log prior probability, both of which can be computed on prior simulations before training ACE. In addition, the choice of beta should consider how broadly the cost function is distributed, which can be straightforwardly estimated with a 1D histogram. Since the cost is always greater or equal to 0, larger betas penalize larger costs more heavily due to the exponential. \n\nConsidering the above, one reasonable heuristic is to choose beta such that it scales the mean or median of the empirical distribution of distances (computed on random pairs of prior simulations) to be in the same range as the mean or median of log-prior probabilities of those same prior samples. As a concrete example, for the Uniform 1D problem, the log-prior is -1.1 (uniform), while the median of prior simulation distances from 10k random pairs is 0.077. Therefore, a beta of 10-20 is a good start, and increasing by a factor of 5-10 further trades off posterior sample diversity for predictive simulation distance / performance.\n\nLastly, many GBI methods for setting beta have been developed (see, e.g., Wu, Martin, Bayesian Analysis 2023). Similar methods could be applied to ACE and, since ACE amortizes inference over beta, it is even amenable to methods which require repeated posterior sampling for different values of beta.\n\nThanks to the reviewer for raising this important practical issue, and we will include a summary of the above in the method and discussion sections.\n\n**Another application with real data**:\nGiven the time constraint of the rebuttal phase, we are not able to conduct an additional set of experiments on a real problem, which requires a suitable dataset and simulator model, in-house expertise, and comparison with existing methods. However, we agree that further applications on real problems are important, and are currently being pursued."}, {"review_id": "KJIzRF3mhX", "reviewer": "Reviewer_CFba", "summary": "This paper proposes amortized cost estimation (ACE) for generalized Bayesian inference (GBI) for SBI. The paper trains a neural network to approximate the cost function. The paper demonstrates results on baseline synthetic SBI examples followed by a real-world application using experimental data from the Allen Cell Types Database. ", "questions": "* Compared to BGI, equation (2) seems to not include the KL divergence term. Can this approach really be called GBI if it is not being regularized by some prior? Why was the KL term not included?\n* Following on from the question regarding the prior, is the addition of the noise ($\\epsilon$) implicitly adding a prior? It might be the case that adding this noise is also enforcing some implicit distance metric between nearby parameters. I.e. does it assume that there is a local Euclidean distance metric?\n* Why does this approach require multi-chain slice-sampling? Is random-walk Metropolis-Hastings not sufficient?\n", "limitations": "* The idea behind the paper is interesting, but the experimental results are difficult to interpret. Further clarifications on the above would help increase the score.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "* The paper is written clearly and is structured well.\n* The idea appears novel. While, GBI was previously defined in [13], its application within the SBI literature is interesting and is therefore of value to the community.\n* Algorithm box makes the implementation of the algorithm clear to the reader.\n", "weaknesses": "* The main weakness of the paper seems to be the experimentation. While the figures are neatly presented, it is not that clear what the reader is supposed to conclude from the results:\n    * For example, it was not clearly defined what the difference between specified and misspecified samples means within Figure 3. It would be helpful to mathematically define these differences rather than describe them within the text. \n    * What is the definition of the GT-GBI? It first appears on line 173 and it is not clear what the true distance is and how it is obtained. This misunderstanding makes it further difficult to understand Figure 3.\n    * Another question is why are the first three rows only showing MSE and the last row showing MMD? This is confusing and makes one question why the Figure is not consistent with the metric. \n    * Additionally, for the MSE/MMD it seems that as $\\beta$ is increased, the MSE/MMD goes down, but the C2ST seems to go up. Why is this the case? What happens in the limit that $\\beta$ goes to infinity?\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for noting the novelty of our work and the clear structure of the paper, and their positive score of 3s (good) in soundness, presentation, and contribution. We appreciate their requests for further clarifications and the opportunity to increase the score. We apologize for the confusion caused by the lack of explicitly presented information in these cases. We address each concern below and will modify the paper accordingly. We hope this clarifies the reviewer’s concerns and that it allows them to recommend our paper for acceptance.\n\n**Difference between well-specified and misspecified observations**:\nThis was a clear oversight on our part, as we failed to refer to Appendix A4.3 on line 615 in the main text, which contains the mathematical definitions for the synthetic well-specified and misspecified observations. We do so now in the main text and expand on the original description, and summarize A4.3 here for convenience:\n\nWell-specified observations for all tasks were prior simulations. For the Uniform 1D, 2 Moons, and Linear Gaussian task, misspecified observations were prior simulations that were successively perturbed with Gaussian noise with fixed variance (i.e., Gaussian random walk) until the sample was outside the range (i.e., [min, max]) defined by 100k prior simulations in all dimensions. For the Gaussian Mixture task, misspecified samples were generated by replacing the second of the two Gaussians in the simulator with N(12.5 * sign(theta), 0.5^2*eye).\n\n**Definition of GT-GBI and true distance**:\nWe apologize for, and will correct the mistake on line 172: It should have been “true cost”, not “true distance”. Throughout the paper, we use “distance” to refer to the output of the distance function (e.g., MSE) between two points in data space, while “cost” refers to the _expectation of the distance_ between all simulations generated by a parameter (theta) and an observation in data space (x). We compute the ground-truth cost function (Eq. 2) for all benchmark tasks, either computing the integral analytically, or numerically via quadrature (i.e., summing over a fine grid). This is then used to define the unnormalized GBI posterior in Eq. 1 and sampled using MCMC or rejection sampling, referred to as GT-GBI samples. This will be clarified in the main text.\n\n**MSE vs. MMD in Fig. 3**:\n3 of the benchmark tasks use MSE as the distance function in data space, since Euclidean distance between two points are easily measured. In the fourth task (Gaussian Mixture), each observation is a set of 5 independently sampled data points, therefore the distance function must measure the statistical distance between two distributions (which MSE cannot, but MMD can). In addition, users may want to use different distance functions, and it is a feature of ACE that it works for a diverse set of distance functions. The first and second columns of Fig. 3 quantify the average distance the posterior predictive simulations achieve for each of the algorithms, hence the first 3 rows have the y-axis labeled as MSE, and the last as MMD. More details in Appendix A4.2.\n\n**MSE/MMD goes down while C2ST goes up with increasing beta**:\nDue to the exponential in Eq. 1, as beta increases, the GBI posterior becomes more concentrated near parameter regions with low cost. Therefore, high-beta posterior predictive simulations have lower distance to the observation (thus lower MSE & MMD). At the same time, since the cost estimation network only learns to approximate cost using finite simulations, if the ground truth posterior is very narrow, then small errors in cost estimation can lead to large changes in C2ST, resulting in generally higher C2ST with increasing beta.\n\nIn the limit of beta approaching infinity, the true GBI posterior would collapse onto the minimizer of the cost function (i.e., the parameter that produces simulations with the lowest average distance to the observation). Samples from this ground-truth posterior would achieve the lowest cost by definition, while samples from the ACE posterior would be close, but slightly higher in cost, and the two sets of samples coming from two different delta functions would be perfectly separable, resulting in a C2ST score of 1 but comparably small predictive distances. We will expand on the current discussion of this in the results section (line 235).\n\n**Prior regularization and addition of noise**:\nOur Eq. 2 does not contain the KL term because it is only a loss to train the neural network to learn an approximation of the cost function. Prior regularization in our GBI posterior happens in Eq. 1, via standard Bayesian updating that balances the (generalized) likelihood and the prior probabilities. In ref 13 (Bissiri et al 2016), Eq. 7 is a loss over the entire posterior approximation, such that it can be optimized directly, thus including a prior KL term. In section 2.3, those authors note that when the generalized likelihood is defined by a cost function, the solution to the posterior approximation problem exactly follows our Eq. 1, thus the two formulations are equivalent.\n\nThe addition of noise to the training samples in data space during learning of the cost function does not have any effect on regularizing the posterior density in parameter space. Instead, it expands the range of prior simulations such that the simulation support may cover observed data for which the simulator is otherwise misspecified, and does not assume local Euclidean distance.\n\n**Alternate MCMC methods**:\nMulti-chain slice-sampling was used for computational efficiency in our experiments. Any MCMC algorithm, such as Metropolis-Hastings, is applicable, since the cost estimation network simply provides the (generalized) log-likelihood. We will note this in the discussion of limitations (line 335)."}, {"review_id": "KM6d8zW4dj", "reviewer": "Reviewer_HKbp", "summary": "This paper studies the problem of simulation-based inference - which can is encountered in a wide range of scientific problems - where one is interested in performing Bayesian inference using simulators with implicit likelihoods. Scientific problems can have two unique properties - a) the predictive quality of the simulation is more important than the posterior b) the model can be misspecified. These can be handled though Generalized Bayesian Inference - which generalizes traditional Bayesian Inference by admitting generalized likelihood functions which are defined by some general loss / cost function defined on the parameters. A challenge with employing Generalized Bayesian Inference in the simulation-based case is that estimating the cost function can be computationally expensive. The authors propose learning an amortized neural network estimator for the cost function to reduce the overhead for enabling GBI in simulation-based inference. This amortized cost estimator can be combined with standard inference algorithms (MCMC) to sample from the Generalized Bayesian posterior. The overall proposed approach consists of 3 stages: 1) Collecting data 2) Training ACE 3) Sampling using ACE. The authors present several experiments on a variety of benchmark simulation-based inference tasks and present a case study with the Hodgkin-Huxley simulator. ", "questions": "* There are already some results with different simulation budgets - but a more careful study of the effect of the number of samples available to train the ACE might be quite useful. \n* Could you comment on the connections to the work on uncertainty estimation mentioned above? \n* Do the authors have any thoughts on extending the approach to higher dimensional settings? ", "limitations": "The authors already discuss several limitations in the paper including a) need for another method to sample from the posterior induced by ACE, b) accommodation of arbitrary distances, c) additional hyper parameters. \n\nIn addition to these I would also add access to sufficient samples for training ACE and being limited to low dimensional settings as limitations. ", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "* The problem of simulation-based inference is an important one for various scientific problems and misspecification is an often under-studied but important aspect in practice. The paper attempts to address this important issue relevant to various communities. \n* The proposed approach is novel within the context of the simulation based inference. (Although as I mention below, similar ideas have been employed is some related areas)\n* Even despite some caveats that I discuss below, the idea is quite neat, and relatively simple (which is a good thing!) Replacing the estimation of the cost with an amortized predictor seems like a natural extension to enable generalized posteriors in the simulation-based inference setting. The simpler algorithm also makes it easier to adopt by practitioners with little overhead. \n* The experiments cover a fairly wide variety of tasks used in prior work in the SBI setting. The case study on the HH model was a nice example of leveraging the method on realistic tasks. \n* The presentation in the paper was generally quite clear. The paper is well written and easy to follow.\n* I also appreciate the authors open-sourcing the code to aid reproducibility. ", "weaknesses": "* A fundamental weakness I see with the method is that it relies on the neural estimator to generalize well from finitely many simulations. This might not be true in practice (certainly isn’t guaranteed). Put differently, an implicit assumption for the approach to work is that the ACE is trained with _enough_ data to be able to provide informative guidance to the sampler learned in the next stages. As with other work relying on NN estimators this is a fundamental failure mode. \n* On a similar note, the method appears to be practical only for relatively low-dimensional problems. The NN learning issues would crop up if the dimensionality of the problem is increased. \n* Similar ideas of learning amortized neural estimators have been recently explored in the uncertainty estimation literature [1,2]. In particular, the secondary predictor learned in [1] amortizes a similar quantity. \n\n[1] DEUP: Direct Epistemic Uncertainty Prediction. Lahlou et al., TMLR 2023.\n\n[2] Epistemic Neural Networks. Osband et al., 2023. ", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their positive assessment of our work, who was very kind to note the importance, novelty, and simplicity of the contribution, as well as our effort on presenting the idea clearly and with included code for reproducibility. We enjoyed the concise and accurate summary of the problem setting surrounding scientific simulators, and we agree that simpler algorithms (with good performance) are more likely to be adopted by practitioners. \n\nWe address the reviewer’s questions and suggestions point by point here:\n\n**Problem of neural network estimation with finite training data**:\nWe agree wholeheartedly with the reviewer’s concern that the reliance of neural networks on having enough training data is a fundamental limitation of our (and all other) neural network-based SBI algorithms. We will explicitly acknowledge this as a limitation in the discussion. We also agree that a careful characterization of performance vs. simulation budget is a critical sensitivity to measure for any algorithm. We demonstrate precisely the reviewer’s point in our experiments: with a simulation budget of 200 (Suppl. Fig. A3), ACE’s cost estimation accuracy is poor, even though posterior predictive distance is still on par or better compared to other neural methods in this low-training data regime. With a moderate increase in simulation budget (to 1000, Suppl. Fig. A4), we already see a marked improvement in cost estimation. This provides guidance on the simulation budget for problems with similar parameter and data dimensionality. We additionally characterize this dependence in our real-world problem, noting that ACE achieves very good performance compared to NPE even with 10 times less training data (100k vs. 1 million), and is further improved with 1 million training samples. \n\n**Extending to higher-dimensional problems**:\nWe acknowledge that here we only provide evidence for good performance in low to moderate dimensionality problems (up to 10 parameters), and we agree that the finite data issue will be further exacerbated in high-dimensional problems. We will include these points in the discussion of limitations as well. However, we note that thus far most successful applications of simulation-based inference are in regimes of parameter dimensionality around 10. For example, the benchmark tasks presented in Lueckmann et al. 2019 all have parameter dimensionality smaller or equal to 10.\n\nIn addition, compared to methods such as NLE and NPE, ACE casts the density estimation problem into a regression problem, where the output of the neural network is always kept as 1-dimensional. Therefore, on a task with high dimensionality, e.g., 30 parameters and 30 data dimensions, NPE needs to learn a rather complex 30-D to 30-D transformation, which can be challenging for a normalizing flow with invertibility constraints, whereas ACE only performs a 60-D to 1-D regression, drastically simplifying the problem. Other approaches employed by existing methods, such as adding a preprocessing/embedding network to first reduce parameter and/or data dimensionality, or extending the algorithm to perform sequential inference that concentrates on a specific observation, can be applied here as well. We will include the above in the discussion section, and we thank the reviewer for raising the possibility of future works along these lines.\n\n\n**Relations to work on uncertainty estimation**:\nThanks for bringing to our attention a related area of literature. From our understanding (and correct us if we’re wrong), both referenced papers train neural networks to predict the _predictive_ uncertainty of a neural network. The “standard” deep NN solves a regression task by predicting the expected (mean) output given an input, while the uncertainty networks in the cited works additionally predict deviations from that, noting various sources of uncertainty (e.g., approximation uncertainty vs. misspecification/epistemic uncertainty). As noted in Lahlou et al 2023, this encapsulates Bayesian NNs, which produce samples from the posterior predictive distribution directly. Furthermore, since the predictive uncertainty can be computed per input without retraining, it is amortized. \n\nConditional density estimation indeed tackles a similar problem, i.e., predicting not just the mean output, but (samples from) a conditional distribution. Furthermore, standard methods in SBI (such as NLE) currently do not naturally account for misspecification. Therefore, both uncertainty-aware / epistemic neural networks and ACE-GBI are relevant for the high-level problem of producing a diversity of good predictive samples under possible misspecification.\n\nOn the other hand, they have very different motivations and regimes of operation, i.e., inverting (misspecified) stochastic mechanistic simulator models vs. giving black-box neural networks awareness of sources of uncertainty. As such, our overall goal, the implementation of amortized estimation of the predictive distance, and its usage in SBI are all quite different from the referenced works on amortized prediction of uncertainty in deep NNs (predictive uncertainty vs. parameter uncertainty, targeting cost vs. posterior variance, using MCMC for sampling vs. directly predicting the full distribution). \n\nNevertheless, there are interesting potential connections between predictive uncertainty, Bayesian NNs, and GBI, which can be explored in future works for robustifying both deep NN and SBI approaches to misspecification and other sources of uncertainty. We thank the reviewer for raising these potential connections, and will note the above discussion in the paper as well."}, {"review_id": "zqKxrjPfqX", "reviewer": "Reviewer_9qMX", "summary": "The authors present a novel approach for SBI by utilizing neural networks (NN) to estimate a generalized cost in GBI. Their proposed method, ACE, demonstrates superior computational efficiency compared to previous approaches, without compromising competitive performance across various evaluation metrics. The authors conducted extensive benchmarking of their method under diverse experimental settings, ranging from 1D to 10D. Furthermore, they conducted experiments using real intracellular recordings, adding an additional dimension to their research.", "questions": "nothing, great work!", "limitations": "The authors do address the limitations of their work.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "The paper is well written and maintains a smooth flow throughout. Furthermore, the contribution made by the authors appears to be a novel and valuable improvement upon existing ideas in the field. Moreover, the authors have thoughtfully acknowledged and addressed the limitations of their work, placing their findings within the broader context of the field. It is worth noting that their code repository is well-organized, adding to the overall credibility and reproducibility of their research.", "weaknesses": "1. To enhance clarity, it would be beneficial to improve the visual presentation of Figure 2. When referring to subfigure C, the arrangement could be adjusted to prevent the automatic focus on the third pane.\n2. Propositions 1 and 2 appear somewhat forced, as their results do not seem entirely compelling or noteworthy. It may be worth revisiting these propositions. In the appendix the authors state, \"The proofs closely follow standard proofs that regression converges to the conditional expectation.\" Thus highlighting the lack of importance of these propositions.\n3. A minor typo on line 116, \"parameter- or\".", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We really appreciate the reviewer’s encouraging and positive comments about our work, in particular highlighting the novelty and value of our proposed contribution to SBI, our extensive benchmark experiments, and the quality of the writing and code for better communication and reproducibility. \n\nIn response to their suggestions regarding clarity, we will:\n\n- Rearrange Figure 3C such that the first panel is bigger and more prominent to prevent automatic focus on the third panel,\n- Move the propositions to the Appendix, such that in the main text we keep the presentation more concise and just state the well-known convergence proof, and refer to the Appendix for more details.\n- Fix the noted typo (thanks for the detailed read!)\n"}]

Author rebuttal:
We would like to thank all the reviewers for their constructive and detailed engagement with our work, resulting in many helpful comments, questions, and opportunities for clarification, as well as ideas for future work. We are especially grateful for several reviewers’ acknowledgement that the problem we tackle is “well-motivated” (bf2N), “important and relevant” (HKbp), that our contribution to the SBI & GBI literature is novel (CFba, HKbp, 9aMX), “timely” (bf2N), and will be of value to the community (CFba, HKbp, 9aMX), while also pointing out that the method itself being simple is a good thing (bf2N, HKbp). We are happy to hear that all 4 reviewers thought the paper was “well-written” and clearly structured, and appreciated the open-source code (HKbP, 9aMX).

At the same time, we agree with many of the concerns and questions raised by the reviewers, and have responded individually to each review in detail. We summarize here the main points and results from new experiments, referring to individual reviewers’ comments whenever applicable. We hope that these follow-up discussions and results clarify their concerns, allowing all reviewers to recommend our work for acceptance.

**Tempered NLE (rev:bf2N)**:
We implemented the tempered Bayesian posterior with NLE (neural likelihood estimation), and tested it on all tasks with the same betas as used for GBI. Overall, we find that tempering NLE does not result in lower posterior predictive distances that we aim to minimize, and is in fact often detrimental compared to standard NLE. “Colder” posteriors (higher beta) slightly improved NLE in just one of the eight scenarios (2 Moons, misspecified observations), but is still outperformed by ACE in this (and all other) tasks. Results in attached pdf, Fig. 1. We will update Fig. 3, method and discussion sections in the paper.

**Noise augmentation lesion study (rev:bf2N)**:
We varied the augmentation noise variance (sigma=0,2,5) during training for ACE in the 10k simulation budget experiment and observed no impact in subsequent performance. We further removed all training data augmentation, including observed data points, and saw barely any effects (attached pdf, Fig. 2). We therefore conclude that at sufficient simulation budgets, removing augmentation would not degrade performance. We did see a slight decrease in performance for two tasks with misspecified observations when all augmentation was removed for ACE trained with 1k simulation budget, supporting the idea that data augmentation aids learning outside of regions produced by the simulator. These results will be discussed in the paper.

**On choosing beta (rev:bf2N) and when beta approaches infinity (rev:bf2N & CFba)**:
We agree that earlier acknowledgement of the additional hyperparameter increases transparency, and will do so in the method section as well. We also agree that choosing beta is an important practical concern, and provide here a heuristic for choosing starting values: a good “baseline value” is such that it scales the average distance across a subset of the training data (precomputed on prior simulations) to be in the same range as the (log) prior probability, both of which can be computed on prior simulations. From there, increasing beta sacrifices sample diversity for predictive distance, as was pointed out that when beta approaches infinity, posterior samples converge onto the minimizer of the cost function. We also note that, since our method is amortized over beta, beta-selection methods which require posterior sampling for several different values can be performed at low computational cost. See individual responses for more details.

**Prior regularization via dKL (rev:CFba)**:
In ref 13, the GBI loss is for optimizing the entire posterior approximation directly, hence prior regularization is included explicitly as the KL term. In contrast, our loss function (Eq. 2) is only for training the cost estimation network, while prior regularization occurs via the standard Bayesian updating in Eq. 1. These two views are consistent, as is also stated in ref 13, Sec. 2.3. 

**Clarification on experimental details and results (rev:CFba)**:
We apologize for erroneous or missing references to methodological details and definitions. We now refer to the definition of well-specified and misspecified observations in the Appendix (A4.3), fixed inconsistencies when referring to distance vs. cost functions and how ground-truth was obtained, and make further clarifications regarding the MSE vs. MMD, and their differing directions of change compared to C2ST with increasing beta. Full details in the individual response.

**Relationship to uncertainty estimation and Bayesian NN (rev:HKbp)**:
We discuss how posterior predictive samples and their distance to the observation (our optimized quantity) conceptually relate to learning predictive uncertainty in works on uncertainty estimation (e.g., epistemic neural networks) and Bayesian NNs, while noting the very different goals and implementations between these two areas of literature.

**Acknowledgement of limitations regarding finite data, high-dimensional problems, hyperparameter sensitivities**:
We thank the reviewers for noting further limitations in our work that we had failed to discuss, and will explicitly acknowledge in the discussion section, including: errors in neural network approximation when learning from finite data (rev:HKbp), potential issues and strategies when applying to higher dimensional problems (rev:HKbp), and sensitivity to various hyperparameters,e.g., beta, sigma, and training budget. Detailed discussions in individual responses.

**Clarification and improvements on presentation and visualization**:
We make a number of clarifications or changes following questions and suggestions raised by all reviewers, including adjustment to Fig. 2C (rev:9qMX), moving the propositions to Appendix for clarity (rev:9qMX), and further explanations on results in Fig. 3 (rev:CFba).

Meta-review: The reviewers unanimously recommend acceptance (7-6-6-7). Several minor issues have been raised and discussed during the author-reviewer discussion period. The authors are encouraged to take this feedback into account in the final version of the paper.
Decision: Accept (poster)