Fast -free Inference of Simulation Models with Bayesian Conditional Density Estimation
Many statistical models can be simulated forwards but have intractable likelihoods. Approximate Bayesian Computation (ABC) methods are used to infer properties of these models from data. Traditionally these methods approximate the posterior over parameters by conditioning on data being inside an -ball around the observed data, which is only correct in the limit → 0. Monte Carlo methods can then draw samples from the approximate posterior to approximate predictions or error bars on parameters. These algorithms critically slow down as → 0, and in practice draw samples from a broader distribution than the posterior. We propose a new approach to likelihood-free inference based on Bayesian conditional density estimation. Preliminary inferences based on limited simulation data are used to guide later simulations. In some cases, learning an accurate parametric representation of the entire true posterior distribution requires fewer model simulations than Monte Carlo ABC methods need to produce a single sample from an approximate posterior.
Conventional ABC algorithms such as the above suffer from three drawbacks. First, they only represent the parameter posterior as a set of (possibly weighted or correlated) samples. A samplebased representation easily gives estimates and error bars of individual parameters, and model predictions. However these computations are noisy, and it is not obvious how to perform some other computations using samples, such as combining posteriors from two separate analyses. Second, the parameter samples do not come from the correct Bayesian posterior, but from an approximation based on assuming a pseudo-observation that the data is within an -ball centred on the data actually observed. Third, as the -tolerance is reduced, it can become impractical to simulate the model enough times to match the observed data even once. When simulations are expensive to perform, good quality inference becomes impractical.
We propose a parametric approach to likelihood-free inference, which unlike conventional ABC does not suffer from the above three issues. Instead of returning samples from an -approximation to the posterior, our approach learns a parametric approximation to the exact posterior, which can be made as accurate as required. Preliminary fits to the posterior are used to guide future simulations, which can reduce the number of simulations required to learn an accurate approximation by orders of magnitude. Our approach uses conditional density estimation with Bayesian neural networks, and draws upon advances in parametric density estimation, stochastic variational inference, and recognition networks, as discussed in the related work section.
2 Bayesian conditional density estimation for likelihood-free inference 2.1 Simulator-based models and ABC Let θ be a vector of parameters controlling a simulator-based model, and let x be a data vector generated by the model. The model may be provided as a probabilistic program that can be easily simulated, and implicitly defines a likelihood p(x | θ), which we assume we cannot evaluate. Let p(θ) encode our prior beliefs about the parameters. Given an observation x o , we are interested in the parameter posterior p(θ | x = x o ) ∝ p(x = x o | θ) p(θ).
As the likelihood p(x = x o | θ) is unavailable, conventional Bayesian inference cannot be carried out. The principle behind ABC is to approximate p(x = x o | θ) by p( x − x o < | θ) for a sufficiently small value of , and then estimate the latter-e.g. by Monte Carlo-using simulations from the model. Hence, ABC approximates the posterior by p(θ | x − x o < ), which is typically broader and more uncertain. ABC can trade off computation for accuracy by decreasing , which improves the approximation to the posterior but requires more simulations from the model. However, the approximation becomes exact only when → 0, in which case simulations never match the observations, p( x − x o < | θ) → 0, and existing methods break down. In this paper, we refer to p(θ | x = x o ) as the exact posterior, as it corresponds to setting = 0 in p(θ | x − x o < ).
In most practical applications of ABC, x is taken to be a fixed-length vector of summary statistics that is calculated from data generated by the simulator, rather than the raw data itself. Extracting statistics is often necessary in practice, to reduce the dimensionality of the data and maintain p( x − x o < | θ) to practically acceptable levels. For the purposes of this paper, we will make no distinction between raw data and summary statistics, and we will regard the calculation of summary statistics as part of the data generating process.
The paper is on likelihood-free inference, that is on parametric inference for models where the likelihood function is too expensive to evaluate. It is proposed to obtain an approximation of the posterior distribution of the parameters by approximating the conditional distribution of data given parameters with a Gaussian mixture model (a mixture density network).  The authors see the main advantages over standard approximate Bayesian computation (ABC) in that - their approach is returning a "parametric approximation to the exact   posterior" as opposed to returning samples from an approximate posterior (line 49), - their approach is computationally more efficient. (line 55)  The paper contains a short theoretical part where the approach is shown to yield the correct posterior in the limit of infinitely many simulations if the mixture model can represent any density. The approach is verified on two toy models where the true posterior is known and two more demanding models with intractable likelihoods. 
Technical quality ---------------------  My ranking is due to the following:  1. The authors say that their approach approximates or targets the exact posterior (e.g. lines 50; 291).  I strongly disagree with the statement. The paper only concerns the approximation of the posterior after conversion of data to summary statistics. But the summary statistics are most often not sufficient so that important information is lost and targeting the "exact posterior" is not possible any more.  2. It is emphasized that the approximation can be made as accurate as required (line 51).  I suppose this becomes possible when the number of components and hidden layer in the network model is increased. Unfortunately, model choice does not seem to get discussed: How should the user choose the number of layers and components, taking in consideration that more data need to be simulated to fit more complex models?  3. While I appreciate theoretical justifications, Proposition 1 did not provide me with much information. Asymptotically, rejection ABC recovers the posterior too, under even weaker conditions, because the threshold (bandwidth) decreases to zero as the number of simulations increases (e.g. Section 3.1. in Blum2010). In that sense, also in rejection ABC, the approximation "can be made as accurate as required".  4. Regression adjustment is a classical method in ABC that deals with the increase of computational cost when the bandwidth epsilon is decreased (see e.g. Section 4.1 of the 2010 review paper by Beaumont). Some of the work is briefly mentioned toward the end of the paper, but no experimental comparison is performed. Such a comparison, however, would have been very important, at least because  - they are standard practice in ABC,  - they aim at producing (approximate) posteriors for epsilon \to zero, like the proposed method,  - there are straightforward theoretical connections to the proposed method (see below the comments on novelty).        Further comments:  1. In Fig1 (left), the more advanced algorithm 2 (MDN with proposal) seems to produce a posterior that is less accurate than the simpler algorithm 1 (MDN with prior). Would you have an explanation for that? Is it not worrisome that the more expensive algorithm 2, which builds on algorithm 1, is able to decrease the accuracy of the approximation?  2. In sections 3.3 and 3.4., the log probability (density?) that the learned posterior assigns to the true parameter value is used to assess the accuracy of the method. I think this measure could be possibly misleading because in many cases, in particular for weakly informative data, the true posterior is not centered at the data generating parameter (consider e.g. the posterior of the mean of a Gaussian). Furthermore, it does not measure the accuracy of the spread of the posterior. Getting the spread (e.g. posterior variance) right, however, is important, as pointed by the authors (line 22).  3. It is rather common to work with uniform priors on a bounded interval (e.g. [0,1]). Is it correct that the estimated posterior pHat in (3) would then be a truncated mixture of Gaussians? For truncated (mixture of) Gaussians    - exact normalization is not possible and numerically challenging in higher dimensions,    - the marginals are generally not truncated (mixtures of) Gaussians (see e.g. Horrace2005)  Can this be a problem for the proposed approach?  4. Eq (7) in the supplementary material shows that q_phi(theta|x_0)/pTilde(theta) equals the likelihood function (as the proportionality factor does not matter). Since pTilde is known and fixed in each iteration, does learning q_phi thus correspond to learning the likelihood function?   Novelty/originality ---------------------  I think the novelty claims of the paper are too strong.  Related work is mentioned in section 4 but the discussion stays on the surface and the proposed method is not placed well in the existing body of research. It is presented as a "new approach to likelihood-free inference based on Bayesian conditional density estimation". But it seems to conceptually belong to existing likelihood-free inference approaches.  1. ABC has been approached via conditional density estimation since more than 10 years (see e.g. section 4.1 of the review by Beaumont 2010). This connection resulted in the methods for regression adjustment. All these methods target the approximate posterior in the limit of epsilon \to zero.  2. When linear models are used, a pilot ABC run is needed to restrict the parameter space, which results in discarding simulations. For nonlinear models, however, all simulations can be used (see e.g. Fig 1 of Blum's 2010 paper in Statistics and Computing, reference 5 in your paper). This also corresponds to "\epsilon-free inference". Moreover, Blum used neural networks to model the relation between data (summary statistics) and parameters as well.  3. The 2012 paper by Fan, Nott, and Sisson that is cited as reference 8 is an epsilon-free inference approach based on conditional density estimation too, and it also makes use of mixture models. While Fan et al focus on approximating the likelihood, the work by Bonassi et al (2011), discussed by Fan et al, focuses on approximating the posterior.  4. Likelihood-free inference methods can be classified into parametric and nonparametric approximations (see e.g. section 3 of reference 11). Classical ABC algorithms correspond to nonparametric approximations where the "epsilon" plays the role of the bandwidth. Parametric approximations, by construction, do not require such a bandwidth (an example is reference 12). The proposed method is a parametric approach so that it is not surprising that it is "epsilon-free".  I agree that other researchers may not have used the particular neural network used in this paper. But it is unclear to me why the paper is presented as a "new approach to likelihood-free inference". It seems more like a technical difference to previous work, and its advantages relative to recent efficient likelihood-free inference methods are neither presented nor discussed.  References ------------  Blum2010: Approximate Bayesian Computation: A Nonparametric Perspective, Journal of the American Statistical Association, 2010, 105, 1178-1187   Beaumont2010: Approximate Bayesian Computation in Evolution and Ecology, Annual Review of Ecology, Evolution, and Systematics, 2010, 41, 379-406  Horrace2005: W. Some results on the multivariate truncated normal distribution, Journal of Multivariate Analysis, 2005, 94, 209-221  Bonassi2011: Bayesian Learning from Marginal Data in Bionetwork Models, Statistical Applications in Genetics and Molecular Biology, 2011, 10.   Update: ----------  Based on the explanation and promises in the author reply, I increased my scores from 2 to 3. I suggest to take the following points to heart when revising the paper.  - Please clarify that the "exact posterior" refers to the distribution of the parameters given the summary statistics and not the actually exact posterior. The expression "exact posterior" is misleading.  - Please clarify the difference to existing methods for likelihood-free inference that are based on (parametric) conditional density estimation, and update novelty claims accordingly.  - Please acknowledge that previous work on likelihood-free inference was able to greatly reduce the computational cost too (e.g. ref 11, 15).  - Please clarify why the approach is approximating the posterior rather than the likelihood function. In my opinion, Eq (3) shows that q_phi(theta|x_0)/ptilde(theta) must be proportional to the likelihood function. And thus by learning phi, you obtain an approximation of the likelihood function.