Source-linked AI summary
Adaptive approximate Bayesian computation
Mark A. Beaumont, Jean-Marie Cornuet, Jean-Michel Marin, Christian P. Robert
TL;DR
Approximate Bayesian computation addresses inference when likelihoods are unavailable, but sequential partial rejection control can bias its posterior approximation. This paper introduces a genuine-importance-sampling population Monte Carlo correction and reports improved performance in a population genetics example.
Problem
When likelihoods lack closed forms, approximate Bayesian computation replaces likelihood evaluation with simulation and tolerance-based acceptance; sequential partial rejection control can nevertheless produce biased posterior weights.
Method
The paper replaces partial rejection-control weighting with a population Monte Carlo correction based on genuine importance sampling and adaptive forward kernels.
Results
The proposed correction improves on the partial rejection-control version and is efficient relative to two standard alternatives in a population genetics example.
Takeaways & Limitations
Sequential approximate Bayesian computation can use progressively improved proposals, but its weighting scheme must be justified by genuine importance sampling arguments.
Takeaways & Limitations
The partial rejection-control weight is incorrect except under a special backward-kernel condition, and random-walk proposals do not satisfy that exception.
Abstract
from arXiv · showhide
Sequential techniques can enhance the efficiency of the approximate Bayesian computation algorithm, as in Sisson et al.'s (2007) partial rejection control version. While this method is based upon the theoretical works of Del Moral et al. (2006), the application to approximate Bayesian computation results in a bias in the approximation to the posterior. An alternative version based on genuine importance sampling arguments bypasses this difficulty, in connection with the population Monte Carlo method of Cappe et al. (2004), and it includes an automatic scaling of the forward kernel. When applied to a population genetics example, it compares favourably with two other versions of the approximate algorithm.
1. INTRODUCTION
Approximate Bayesian computation avoids evaluating unavailable likelihoods by accepting simulations near the observed data, while sequential variants improve proposals across iterations. The paper identifies bias in the partial rejection-control weighting and motivates a corrected population Monte Carlo approach.
- Approximate Bayesian computation: Approximate Bayesian computation bypasses unavailable likelihood calculations by simulating parameters and data, accepting draws whose discrepancy from observations is below a tolerance.For continuous observations, exact equality is replaced by ϱ(x, y) ≤ϵ.
- Approximate Bayesian computation: The algorithm’s output targets a tolerance-conditioned distribution, π{θ | ϱ(x, y) < ϵ}, rather than generally the exact posterior.
- Sequential methods: Marjoram et al.’s Markov chain version has the posterior as its stationary distribution, with exact acceptance replaced by a discrepancy tolerance for continuous data.
- Illustrative example: In the toy mixture example, ϵ = 0 · 025 produced a density indistinguishable from the exact posterior, separating tolerance approximation from proposal-mixing problems.The random-walk scale τ controls convergence, whereas ϵ controls posterior approximation.
- Sequential methods: Sequential approximate Bayesian computation generates populations across iterations, reuses earlier samples for later proposals, and decreases ϵ over time.The tolerance sequence may be deterministic or based on earlier-iteration quantiles.
- Paper contribution: The proposed correction uses genuine importance sampling and population Monte Carlo ideas, rather than the backward-kernel construction used in the partial rejection-control method.The paper reports theoretical and experimental evidence that the original weight is biased and presents the correction as an improvement.
Distribution of the partial rejection control sample
The partial rejection-control sample has a joint density involving the previous and current parameters, but its proposed weight is generally biased. The bias disappears only under a special backward-kernel condition.
- Under ϵ = 0 and a correctly sampled previous parameter, an accepted pair has joint density proportional to π(θ(t−1) | y)Kt(θ(t) | θ(t−1))f(y | θ(t)).
- The weight is biased unless its inner function integrates over θ(t−1) to the same constant for every θ(t).
- The weighting is incorrect for random-walk backward kernels because the weighted output is not distributed from π(θ | y).An exception is possible when Lt−1(θ(t−1) | θ(t)) = g(θ(t−1)).
- The partial rejection-control weight omits f(y | θ(t−1)) from its denominator, so replacing the posterior with the prior there is not justified by approximate Bayesian computation’s accept-reject principle.
A mixture illustration
The partial rejection control algorithm remains biased relative to the tolerance target across tolerance values, whereas the population Monte Carlo version recovers that target in the mixture example.
- A mixture illustration: The bias in the standard algorithm's weights has no reason to vanish under the tolerance approximation, although it generally decreases as ϵ increases.The decrease is consistent with the limiting tolerance target becoming the prior.
- A mixture illustration: At τ = 0 · 15, five consecutive iterations of the partial rejection control algorithm poorly cover the tails for every tested ϵt.The tested sequence decreases from ϵ1 = 2 to ϵ5 = 0 · 01.
- A mixture illustration: The figure compares weighted samples from five consecutive iterations of the standard and population Monte Carlo algorithms across ϵt = 2, 1.5, 1, 0 · 5, 0 · 01.Both rows use τ = 0 · 15; dotted curves show the exact posterior and full lines show the tolerance target.
Population Monte Carlo
The population Monte Carlo correction replaces the biased weighting with genuine importance sampling and adaptively scales the forward kernel. In the mixture example, it recovers the target under fixed, large, or adaptive scales, while its dominant computational cost remains simulation.
- Population Monte Carlo: The correction uses a standard importance sampling perspective rather than estimating the unknown likelihood from earlier samples.The missing factor in the biased weight is the unknown likelihood f(x | θ(t−1)).
- Population Monte Carlo: The population Monte Carlo scheme assigns importance weights to accepted simulations generated from the iteration-specific proposal distribution.Its unbiasedness follows from an importance-sampling identity.
- Population Monte Carlo: The method remains valid when Kt depends on simulations from earlier iterations, but asymptotic Kullback–Leibler improvement requires modifying Kt across iterations.If the random-walk variance is unchanged, the target approximation also does not change, making one iteration with twice as many points more efficient.
- Population Monte Carlo: The scheme adaptively selects component-wise random-walk scales by optimizing a Kullback–Leibler criterion.For each component, the optimal scale is linked to posterior variability and can be updated using empirical variance.
- Population Monte Carlo: The algorithm uses a decreasing sequence of tolerance thresholds ϵ1 ≥ · · · ≥ ϵT and iterates over N particles at each stage.The implementation begins at t = 1 and repeats the particle procedure for 2 ≤ t ≤ T.
- Population Monte Carlo: O(T N^2) is the computational cost of the weight-update step, while repeated sampling dominates realistic applications.In the population genetics example, at least 95% of computing time is spent in the repeat loop and less than 5% elsewhere.
- Population Monte Carlo: The population Monte Carlo version recovers the mixture target with τ = 0 · 15, τ = 1/0 · 15, or adaptive τt’s.The estimated variance stabilizes very quickly, and graphical differences among these implementations are difficult to detect.
A population genetics example
The population genetics experiment compares standard approximate Bayesian computation, tempered Markov chain Monte Carlo, and population Monte Carlo under matched computing times and shared priors. Figure 2 evaluates posterior-density variability across five independent runs against a reference posterior.
- Experimental setting: The experiment models two populations that diverged from a common ancestor using genotypes from 50 diploid individuals at five microsatellite loci.The loci follow a strict stepwise mutation model, and the populations do not exchange genes after divergence.
- Experimental setting: Simulated datasets use N1 = Nanc = 10,000, N2 = 2,000, tdiv = 1,000, and µ = 0.0005, with identical conclusions across datasets.
- Experimental design: The tolerance region uses twelve summary statistics, with Euclidean distance between observed and simulated summary-statistic vectors.Statistics include diversity, variance, Garza-Williamson, FST, (δµ)2, and assignment measures for the two populations.
- Experimental design: Three analyses use the same prior distributions: standard approximate Bayesian computation, tempered Markov chain Monte Carlo, and population Monte Carlo.The priors cover effective sizes, divergence time, and mutation rate.
- Results: Figure 2 compares five independent runs of the three methods across posterior densities for four identifiable parameters against a reference posterior.Rows represent population Monte Carlo, tempered Markov chain Monte Carlo, and Beaumont et al.’s version; vertical lines mark true parameter values.
- Experimental design: 4.5 minutes is the average computing time required by all three versions in the matched comparison.Reference posteriors use 3 × 10^4 simulated datasets with ϵ set to the 0.01 distance quantile.
4. CONCLUSION
The conclusion identifies biased weights as a problem for Sisson et al.’s scheme and presents an importance-sampling population Monte Carlo alternative. This alternative uses the same transition kernels and computing power while adapting proposals toward the tolerance target.
- Conclusion: Sisson et al.’s algorithm induces biased weights that visibly reduce approximation quality.
- Conclusion: The same Markov transition kernels and computing power can instead produce an unbiased scheme.
- Conclusion: Population Monte Carlo uses an importance argument that avoids the backward kernel required in the earlier formulation.
- Conclusion: The adaptive scheme adjusts the proposal kernel toward the target π{θ | ϱ(x, y) < ϵ}.
- Conclusion: Iterations can stop when proposal parameters settle, while tolerance-sequence selection can rely on stabilised posterior-quantity estimators.