Source-linked AI summary
Control functionals for Monte Carlo integration
Chris J. Oates, Mark Girolami, Nicolas Chopin
TL;DR
Monte Carlo integration can be too slow for complex models because root-n convergence is problematic when simulations or evaluations are expensive. The paper introduces control functionals, a gradient-based non-parametric extension of control variates that is unbiased, compatible with un-normalised densities, and post-hoc while achieving super-root-n convergence. Theoretical and empirical results support its use for integration in hierarchical models and nonlinear differential-equation models.
Problem
Complex models make root-n Monte Carlo estimation costly, motivating methods that combine unbiasedness, un-normalised-density compatibility, super-root-n convergence, and post-hoc use.
Method
Control functionals use gradient information from the sampling density within a reproducing-kernel, gradient-based function space to construct non-parametric control-variate estimators.
Results
Control functionals provide guaranteed super-root-n convergence and empirical state-of-the-art estimation in hierarchical models and models based on nonlinear differential equations.
Takeaways & Limitations
Variance reduction can be applied after sampling without modifying sampling or model code, with kernel choices deferred until after expensive simulations.
Takeaways & Limitations
The presented control-functional version is limited in the dimensionality of problems for which it is effective and is not parameterisation-invariant.
Abstract
from arXiv · showhide
A non-parametric extension of control variates is presented. These leverage gradient information on the sampling density to achieve substantial variance reduction. It is not required that the sampling density be normalised. The novel contribution of this work is based on two important insights; (i) a trade-off between random sampling and deterministic approximation and (ii) a new gradient-based function space derived from Stein's identity. Unlike classical control variates, our estimators achieve super-root-$n$ convergence, often requiring orders of magnitude fewer simulations to achieve a fixed level of precision. Theoretical and empirical results are presented, the latter focusing on integration problems arising in hierarchical models and models based on non-linear ordinary differential equations.
1 Introduction
Complex models make root-n Monte Carlo estimation impractical when simulations or function evaluations are expensive. The paper introduces control functionals to provide unbiased, super-root-n, post-hoc variance reduction compatible with un-normalised densities.
- Motivation: Complex models can make root-n Monte Carlo convergence too slow when simulations or function evaluations are expensive.A recent climate model required approximately 10^6 core-hours per simulation.
- Motivation: The paper evaluates variance-reduction methods against unbiasedness, un-normalised-density compatibility, super-root-n convergence, and post-hoc applicability.Post-hoc methods separate sampling from variance reduction and require no sampling-code modification.
- Contribution: No existing technique satisfies all four desiderata, whereas control functionals are presented as satisfying them simultaneously.The method applies retrospectively to samples and to densities known only up to proportionality.
- Contribution: Control functionals extend classical control variates non-parametrically by constructing surrogate functions from basis functions with zero expectation.Classical control variates solve a generally misspecified linear regression and typically provide at most a constant-factor variance reduction.
- Contribution: The proposed gradient-based function-space approach yields closed-form estimators with guaranteed faster-than-root-n convergence and empirical support in hierarchical and nonlinear differential-equation models.The paper reports state-of-the-art estimation in these applications.
2.1 Set-up and notation
The setup considers a positive density on an open subset of R^d and a square-integrable test function whose expectation is the integration target. The method uses cached function values and density-gradient information at sampled states for post-hoc estimation.
- 2.1 Set-up and notation: X is a random vector in an open set Ω⊆R^d with positive density π, and f is assumed square-integrable under π.The target is the expectation μ(f) under the sampling density.
- 2.1 Set-up and notation: At sampled states x_i, the function values f(x_i) and gradients ∇_x log π(x_i) are pre-computed and cached.These quantities support later estimation without further model or function evaluations.
- 2.1 Set-up and notation: After caching, the method requires no further evaluations of π or f, making it a post-hoc scheme.The sampling states and previously computed quantities are reused for variance reduction.
2.2 From control variates to control functionals
Control functionals extend control variates by combining deterministic approximation with residual random sampling. A Stein-identity-based function space provides analytically tractable corrections whose increasing approximation accuracy can yield super-root-n convergence.
- Deterministic approximation: The method splits the data into an approximation subset D0 and an independent sampling subset D1, trading deterministic approximation against random sampling.For m < n, D1 supplies the unbiased Monte Carlo correction after fitting the surrogate on D0.
- Deterministic approximation: The surrogate s_f,D0 approximates f using D0, has analytically tractable expectation, and produces a residual estimator with the same target mean.Its residual variance satisfies σ2(f_D0) = σ2(f − s_f,D0).
- From control variates to control functionals: An increasingly accurate infinite-dimensional approximation space, rather than a finite-dimensional control-variate family, enables super-root-n convergence.The paper identifies this expansion of the approximation space as the key insight beyond classical control variates.
- Theoretical rate: If m = O(n^γ) and the expected functional approximation error decays with exponent δ, the mean-squared error is O(n^-1-γδ).The result gives faster-than-root-n convergence when the approximation error decreases sufficiently with m.
- Control variates based on Stein’s identity: The construction assumes a differentiable sampling density and boundary or tail conditions that make the relevant surface integrals vanish.For unbounded domains, the boundary requirement is interpreted as a tail condition.
- Control variates based on Stein’s identity: The Stein-based correction uses ψ(x) = ∇x · φ(x) + φ(x) · u(x), with u(x) = ∇x log π(x), and has zero expectation under boundary conditions.Consequently, for s_f,D0(x) = c + ψ(x), the analytically tractable expectation is c.
2.3 Theory
Control functionals construct a gradient-based Hilbert space of zero-mean functions from an RKHS, enabling regularized approximation of integrands. Under stated assumptions, the resulting estimators are unbiased and can achieve super-root-n convergence, with explicit error bounds and practical computational advantages.
- A Hilbert space of control functionals: The control-functional space H0 is generated from RKHS functions through ψ = ∇φ + φ∇log π and contains only zero-mean functions.The construction remains valid for possibly un-normalized π under the paper’s assumptions.
- A Hilbert space of control functionals: The gradient-based kernel k0 defines H0, while H+ = C + H0 augments it with constant functions to represent integrands as f = c + ψ.The analysis assumes f ∈ H+, equivalently that the fundamental equation has a solution in the chosen RKHS.
- Consistent approximation and asymptotics: O(n^-7/6) mean-squared error is obtained for IID D0 samples with λ = O(m^-1/2), while preserving unbiasedness under the theorem’s assumptions.The result improves on the usual Monte Carlo rate, whose baseline is root-n convergence.
- Non-asymptotic bounds: The method’s additional RLS cost is O(m^3), but this is treated as negligible when simulation from π or evaluation of f dominates computation.The paper also notes that multi-splitting remains unbiased, whereas the simplified m = n estimator reduces variance when consistency is sufficient.
- Non-asymptotic bounds: Theorem 3 supplies an explicit error bound and a discrepancy-based design principle, recovering the root-n rate when kernel correlations vanish between distinct points.The discrepancy framework parallels deterministic approximation ideas, while the extreme zero-correlation case reduces to ordinary Monte Carlo behavior.
2.4 Illustration
In a one-dimensional synthetic integration problem, control functionals showed lower variance and super-root-n convergence, while performance remained favorable as dimension and problem structure became more challenging.
- Synthetic illustration: The simplified CF estimator had negligible bias (<∼10^-3) compared with Monte Carlo error (∼10^-2).This supports using the simplified estimator when reduced variance is preferred and exact unbiasedness is not essential.
- Synthetic illustration: CF estimators had the lowest variance across all sample sizes considered, although their precision gains carried comparable computational costs in this unconstrained synthetic example.In realistic applications, CF computation was expected to be negligible relative to simulation from the model.
- Asymptotic assessment: CFs achieved super-root-n convergence and significantly outperformed piecewise-linear Riemann sums in the synthetic experiment.The arithmetic mean and control variates showed constant-factor variance reduction, whereas Riemann sums and CFs achieved super-root-n convergence.
- Generality: Increasing dimensionality to d = 3 and d = 5 challenged all variance-reduction methods, but CFs continued to outperform alternatives.Across additional test functions and densities, CFs could dramatically reduce variance when low-degree polynomial control variates were inadequate or the density was complex.
- Generality: CF performance degraded when theoretical assumptions (A2) and (A5) were violated.The reported examples marked “CF ×” showed poor estimator performance under those assumption violations.
3 Applications
The applications examine control functionals in hierarchical Gaussian-process prediction and nonlinear ODE normalising-constant estimation. Across both settings, control functionals reduce estimator variance, including relative to established control-variate approaches.
- 3.1 Hierarchical models: Control functionals address hyper-parameter marginalisation in hierarchical models, where posterior prediction can require prohibitive computation.The SARCOS application uses hierarchical GP regression for robot-arm inverse dynamics.
- 3.1 Hierarchical models: Each integrand evaluation in the GP prediction problem requires O(N 3) operations because of matrix inversion.A subset-of-regressors approximation reduces this cost, but realistic evaluations remain expensive.
- 3.1 Hierarchical models: For the vast majority of 4,449 SARCOS integration problems, CF has lower estimator variance than both the arithmetic mean and ZV control variates.Estimator standard deviations were computed over 10 independent Monte Carlo realisations for each test sample.
- 3.1 Hierarchical models: CFs combined with RQMC outperform RQMC alone in the evaluated extension.The extension uses RQMC to draw the design points D0.
- 3.2 Non-linear ODE models: In the van der Pol oscillator, controlled TI improves on standard TI, while CF produces a more substantial reduction in estimator variance.The comparison concerns distributions of 100 independent realisations of estimators for log p(y).
- 3.2 Non-linear ODE models: CF adds essentially zero cost relative to population MCMC, whose repeated ODE-system solutions dominate computation.The CF experiment uses a simplified but biased estimator because numerical quadrature already makes TI biased.
4 Discussion
The discussion presents control functionals as a general, post-hoc integration method with super-root-n convergence and identifies methodological, theoretical, and application boundaries. The paper highlights usefulness when gradient information is available and dimensionality is not too large.
- 4 Discussion: Control functionals provide a post-hoc variance-reduction approach requiring no changes to sampling or model code.Kernel and other implementation choices can also be made after expensive simulations.
- 4 Discussion: The method achieves super-root-n convergence and has demonstrated practical utility when gradient information is available.The discussion also limits effectiveness to problems whose dimensionality is not too large, such as dimensions ≤10.
- 4 Discussion: The presented estimates are not parameterisation-invariant, and effective parametrisations remain an open post-hoc design question.The paper notes that both f and π can be transformed through importance sampling.
- 4 Discussion: The presented CF version is limited in the dimensionality of problems for which it is effective.High-dimensional functional-approximation techniques are identified as a direction for future research.
- 4 Discussion: For unbounded domains, existing theoretical analysis relies on a uniform ergodicity assumption that is too strong.Future work aims to relax this constraint for MCMC-based sampling with un-normalised densities.
- 4 Discussion: The proposed applications target un-normalised Bayesian densities, while extensions to unknown parameter-dependent normalising constants are suggested.Examples include Markov random fields and random network models.
A Proofs
The proofs establish unbiasedness, RKHS structure, variance decomposition, and super-root-n mean-squared-error rates for control-functional estimators. They also derive the estimator through kernel interpolation and quantify the role of design-point approximation.
- Variance analysis: The estimator’s mean-squared error equals (n−m)^−1 times the expected residual variance after fitting the control functional.This separates random-sampling error from approximation error at the design points.
- Convergence rates: With m=O(n^γ) and residual variance O(m^−δ), the overall rate is O(n^−1−γδ).The proof substitutes the residual approximation rate into the variance decomposition.
- Unbiasedness: Stein’s identity gives µ(ψ)=0, so the control-functional estimator preserves the target expectation through its constant component.The divergence theorem and boundary condition establish the zero-mean property.
- Theorem 1: The control-functional construction uses gradient-based functions in an RKHS and defines a vector space H0 through the canonical feature map.The proof verifies differentiability under kernel smoothness assumptions.
- Convergence rates: Under the stated well-posedness conditions, residual variance O(m^−1/6) yields MSE O(n^−7/6).This uses λ=O(m^−1/2), γ=1, and δ=1/6.
- Estimator construction: The RKHS interpolation problem is reduced to solving for a constant and kernel coefficients, with fitted values evaluated at the held-out points.The representer theorem and Woodbury identity produce the coefficient formulas and estimator weights.
- Estimator construction: The resulting CF weights sum to one, and the discrepancy term is D(D0,D1)=w^T K w.This expression appears after substituting the estimator representation into the error calculation.