Source-linked AI summary

Sequential Entropic Optimal Transport for Bayesian Inverse Problems via Multilevel Monte Carlo

Jason Beh, Pierre Del Moral, Ajay Jasra

arXiv:2610.29546v1math.NAcs.NAstat.CO

TL;DR

PDE-based Bayesian inverse problems require numerical posterior approximations, and multilevel estimation depends on coupling posteriors across solver resolutions. This paper learns entropic couplings with light Schrödinger bridges and corrects transport with Metropolis–Hastings steps; its numerical results show faster performance than MLSMC at tight accuracies in the reported tests.

  • Problem

    PDE-based Bayesian inverse problems require numerical solvers, while multilevel Monte Carlo needs coupled posterior samples across increasingly precise solves to reduce computational effort at a prescribed mean square error.

    Source: §1 · p. 1, §1.3 · p. 3
  • Method

    The method learns entropic couplings between consecutive posterior targets with light Schrödinger bridges and follows transport with Metropolis–Hastings steps that preserve pairing.

    Source: §2.3 · p. 5
  • Results

    TransportMH reaches a relative RMSE of 3×10^-4 in 3.5 s against 5.9 s in test case 1, and 10^-2 in 13 s against 40 s in test case 2.

    Source: §4.5 · p. 13
  • Takeaways & Limitations

    The reported gains are greatest at coarse transitions where incremental weights degenerate; the authors suggest using transport selectively at low effective sample size and reducing fine-level sweeps.

    Source: §5 · p. 16

    On problems where all transitions have high ESS, the MLSMC variants are faster at every tested accuracy; the transport advantage is specific to transitions with degenerate incremental weights.

    Source: §4.5 · p. 16

Abstract

from arXiv · show

We consider the problem of Bayesian estimation associated to inverse problems governed by partial differential equations (PDEs). In most cases of practical interest, one must use numerical solvers to implement computational methods based on Monte Carlo methods. This inevitably leads to an adaptation of multilevel Monte Carlo (MLMC) methods to reduce computational efforts to achieve a pre-specified mean square error (MSE). MLMC is based on sampling coupled versions of Bayesian posteriors associated to increasingly precise numerical solves of the PDE. In this article we present a novel method to implement MLMC based upon sequential entropic optimal transport couplings learned with light Schr\"odinger bridges, and correct the transport with Metropolis--Hastings steps that preserve the pairing. We bound the mean square error associated to a version of our new estimation method and show, for an elliptic inverse problem, that to achieve a MSE of $\mathcal{O}(\delta^2)$ ($\delta >0$) one can have an optimal cost of $\mathcal{O}(\delta^{-2})$. The results feature explicit bounds on the effect of the fitted transport and of the Metropolis--Hastings correction. From a numerical perspective, on elliptic, Darcy-flow and Navier-Stokes problems, the method reduces the variance of the most precise solver by a factor of five to ten and is faster than multilevel sequential Monte Carlo [arXiv:1503.07259] at tight accuracies.

1 Introduction

PDE-based Bayesian inverse problems require numerical posterior approximations, motivating computational methods that couple posteriors across solver resolutions. The paper introduces MLMC methods using learned entropic couplings, with an analysis and numerical evaluation of the approach.

1 Introduction

  • Bayesian inverse problems with PDEs require numerical solves, so posterior densities are available only point-wise and with solver error.This motivates computational methods such as Markov chain simulation.Source: §1 · p. 1

1.3 Related work, contributions and paper structure

  • MLMC estimates a finest-level posterior expectation as a telescoping sum of corrections between consecutive levels, using coupled samples to reduce variance and cost.The paper places entropic optimal transport couplings within this framework.Source: §1.3 · pp. 2–3

1.2 Background on entropic optimal transport

  • Entropic optimal transport provides stochastic couplings of two distributions by minimizing expected squared distance with an entropy penalty; light Schrödinger bridges make these couplings cheap to fit and exact to simulate.The bridge parametrizes the coupling potential with a Gaussian mixture.Source: §1.2 · p. 2, §1.3 · p. 3

1.3 Related work, contributions and paper structure

  • The paper presents two MLMC methods using entropic couplings between consecutive posteriors: one with light Schrödinger bridges and another augmented by Metropolis–Hastings steps.It bounds the MSE for a version of the estimator and reports an optimal O(δ^-2) cost for MSE O(δ^2) in an elliptic inverse problem.Source: §1.3 · p. 3

2 Multilevel Estimation with Learned Entropic Couplings

The method learns entropic couplings between consecutive posterior targets and uses the resulting paired samples for multilevel estimation. TransportMH follows transport with invariant Metropolis–Hastings steps, retaining parent-child pairing while reducing transport error.

2.2 Multilevel transport with light Schrödinger bridges

  • At each level, the method learns an entropic coupling of consecutive targets with a light Schrödinger bridge and uses it to move particles for multilevel estimation.The estimator is based on paired differences along the transport.Source: §2.2 · p. 4

2.3 Transport followed by Metropolis–Hastings

  • TransportMH applies a target-invariant Markov kernel for a finite number of steps, starting from each transported particle and keeping it paired with its parent.The implementation uses a preconditioned Crank–Nicolson Metropolis–Hastings kernel.Source: §2.3 · pp. 5–6
  • Finite Metropolis–Hastings sweeps reduce but do not eliminate the distance between the particle law and the target, so the estimator remains biased.The step count and proposal scale adapt to transition difficulty measured by the effective sample size fraction.Source: §2.3 · pp. 5–6
  • Larger Metropolis–Hastings moves trade weaker parent-child pairing and higher variance for smaller bias.The transition cost includes target-density evaluations for training weights and transported particles and sweeps, plus the stochastic-gradient fit.Source: §2.3 · p. 6

3 Mathematical Analysis

The analysis bounds the paired estimator’s MSE through transport costs between levels, then extends the bounds to bounded-misfit posteriors and fitted kernels with or without Metropolis–Hastings correction. It identifies conditions for multilevel complexity rates and explains how fitting error and correction-induced displacement affect those rates.

3.1 MSE under exact couplings

  • The oracle estimator’s MSE is controlled by transport costs between levels, with paired-level covariances bounded by the square of the sum of level standard deviations.The analysis isolates the effect of transport cost by replacing fitted kernels with exact couplings.Source: §3 · p. 6, §3.1 · p. 7
  • With the stated bias, variance, and sample-cost rates, classical sample allocation yields multilevel complexity bounds even though paired level corrections are dependent.The rates depend on the problem and solver, while the coupling choice affects the variance constant.Source: §3.1 · p. 7

3.2 Posteriors with bounded misfit

  • The bounded-misfit analysis covers the elliptic test case, while the Darcy-flow and Navier–Stokes priors fall outside its compact-support assumption.Those latter problems are tested numerically despite falling outside (A2).Source: §3.2 · p. 8
  • Under bounded-misfit assumptions, consecutive posteriors are close in W2 at the forward-map-error rate, supporting α = 1 and β = 2 for Lipschitz observables.The corresponding rate conclusions require couplings whose costs are close to optimal.Source: §3.2 · p. 8
  • For the one-dimensional elliptic problem, α = 1, β = 2, and ε_l ≲ h_l^2/(d log(1/h_l)) give MSE at most δ^2 at cost O(δ^-2).For solver cost growth ζ ≥ β, the stated alternatives are O(δ^-2 log^2 δ) or O(δ^-2-(ζ-β)/α).Source: §3.2 · p. 9
  • Fitting adds bias of order (dΔ_L)^1/2 and a term of order dΔ_l^1/2 to each level variance; preserving β = 2 requires every fit to reach the finest-level accuracy scale.Without correction, accumulated fitting error can leave a bias that does not decrease with L, and LightSB stalls.Source: §3.2 · p. 9
  • If the Metropolis–Hastings kernels contract uniformly and fitting errors are bounded, the accumulated bias is bounded uniformly in the number of levels.The correction kernels leave each π_l invariant; the stated bound assumes τ_l ≤ τ̄ < 1 and E_l ≤ Ē.Source: §3.2 · pp. 9–10
  • More Metropolis–Hastings sweeps reduce inherited bias but increase displacement, loosening the pairing; at fine levels, three sweeps multiply paired-difference variance by about four.The reported fine-level trade-off occurs when the weights are nearly uniform.Source: §3.2 · p. 10

4 Numerical Experiments

The numerical experiments compare TransportMH and related samplers on elliptic, Darcy-flow, and Navier–Stokes inverse problems. Their setup uses common initial samples and sample counts across methods, with TransportMH incurring greater cost at equal level count.

4 Numerical Experiments

  • The experiments cover elliptic, Darcy-flow, and Navier–Stokes inverse problems; the first two use finite-element forward maps, while the third uses a nonlinear, time-dependent model.The Navier–Stokes task recovers initial vorticity from noisy later-time point observations.Source: §4 · p. 10
  • The first transition is difficult across all three tests, with ESS fractions of 3%, 1.2%, and 1.8%; later transitions are easier in tests 2 and 3.At the second transition, the ESS fractions are about 60% and 73% in test cases 2 and 3.Source: §4 · p. 10

4.1 Common setup

  • All four methods share initial samples and level sample counts; TransportMH differs from the ESS-sweep MLSMC variant only by transport and costs more at equal L.TransportMH applies transport evaluations and Metropolis–Hastings sweeps at every level, including L, where MLSMC does not mutate.Source: §4.1 · pp. 10–11

LightSB fits.

The experiments fit light Schrödinger bridge transports across three PDE inverse problems and compare their accuracy, runtime, and level-correction variance with MLSMC. TransportMH performs best at tight accuracy, with its variance gains concentrated at coarse, difficult transitions.

LightSB fits.

  • Sample sizes follow Corollary 3.2, with N_l ≥ 100; the geometric entropic schedule ε_l = 0.05·4^−(l−1), floored at 10−5, follows shrinking distances between levels.Given the stated rates, β > ζ for all three problems, and the allocation targets MSE δ^2 at cost O(δ^−2).Source: LightSB fits · p. 11

4.2 Test case 1: one-dimensional elliptic problem

  • The one-dimensional elliptic case uses piecewise-linear finite elements on 2^l+3 uniform elements and observes the solution at a fixed point using the continuum solution.A finite-element solve is a tridiagonal system; the observable is held constant across levels.Source: §4.2 · p. 12

4.3 Test case 2: two-dimensional Darcy flow

  • The Darcy case uses P1 elements on successively refined square meshes, 25 pressure sensors, and mean permeability as its observable, which is not globally Lipschitz.Consequently, Section 3 does not cover this observable.Source: §4.3 · p. 12

4.4 Test case 3: two-dimensional Navier–Stokes

  • The Navier–Stokes case uses finite differences and Runge–Kutta time stepping; finite differences preserve nonzero level corrections for this smooth flow, and the observable is linear and Lipschitz.Vorticity is observed at 32 points and two times; the observable is mean initial vorticity over a quarter-domain.Source: §4.4 · p. 12

4.5 Results

  • TransportMH reaches relative RMSE 3×10−4 in 3.5 s versus 5.9 s (case 1), 10−2 in 13 s versus 40 s (case 2), and 2×10−2 in 5.0 s versus 59 s (case 3).In cases 1 and 2, MLSMC variants are cheaper at loose accuracy; TransportMH alone reaches the smallest targets. In case 3, its L = 4 MSE is 4.5 times smaller than MLSMC-pCN with ESS sweeps for 1.8 times the wall time.Source: §4.5 · p. 13, Table 2 · p. 13, Table 3 · p. 13
  • Uncorrected LightSB stalls at relative bias 10−3–10−2; MH removes accumulated transport error, while level 1 dominates variance and TransportMH cuts its standard deviation about 6×, 5×, and 10×.The reductions are for test cases 1–3, respectively, against MLSMC-pCN; the level-1 correction dominates the variance for all methods.Source: §4.5 · pp. 13, 15, 16
  • TransportMH and MLSMC-pCN with ESS sweeps use nearly the same number of PDE solves, but TransportMH adds solves where they cost most; three fine-level MH sweeps roughly quadruple paired-difference variance.Skipping those sweeps could reduce the case-2 L = 6 runtime from 119 s to roughly 60–65 s, but its effect on bias was not tested.Source: §4.5 · p. 16

5 Discussion

The discussion proposes applying learned transport selectively at difficult transitions and reducing fine-level sweeps. It also identifies limits in the fitted-kernel analysis and the potential need for structured potentials in higher dimensions.

5 Discussion

  • A hybrid scheme could use transport only where effective sample size is low and fewer or no fine-level sweeps, with adaptive ε_l and mixture-component counts.The numerical results place TransportMH's gain at coarse, hard transitions.Source: §5 · p. 16
  • The analysis bounds fitted-kernel effects only crudely; sharper bounds require controlling stochastic-gradient fit error, Gaussian-mixture approximation error, and pCN-kernel contraction.The theoretical analysis otherwise covers exact couplings and the fitted kernels, with and without Metropolis–Hastings correction.Source: §5 · pp. 16–17
  • Higher-dimensional parameters such as discretized random fields may require more structured potentials, for example low-rank covariances.Source: §5 · p. 17

A Proofs for Section 3

The proofs combine posterior-distance bounds, bias and variance estimates, and the standard multilevel complexity argument. Under the stated level-variance and cost assumptions, the resulting sample allocations are non-increasing across levels.

Proofs for Section 3

  • The multilevel complexity argument uses the bias, cost, and variance bounds to choose the terminal level and level sample counts; decreasing mesh sizes make those counts non-increasing.The proof invokes the multilevel complexity theorem of [9] and notes as the original result.Source: A Proofs for Section 3 · p. 17
  • The mean square error is decomposed into squared discretization bias and variance, with the variance bounded using independent, identically distributed lineages.The bias is bounded by the Lipschitz constant of the test function and the posterior approximation error.Source: A Proofs for Section 3 · pp. 17–19
  • A transport inequality and a bound on the posterior density ratio yield a bound on the distance between consecutive level posteriors.The density-ratio argument bounds the change in potential; the transport inequality converts the resulting divergence bound into a Wasserstein bound.Source: A Proofs for Section 3 · pp. 17–18
  • A discretized coupling approximates the optimal transport cost within 8dε + 4dε^2, while its mutual-information penalty is bounded by d log(1 + 2/ε).The construction partitions the domain into cubes of side at most ε and resamples within paired cubes.Source: A Proofs for Section 3 · p. 18
  • The fitted transport and Metropolis–Hastings correction are controlled through total-variation bounds, a contraction recursion, and a bound on particle displacement.The recursion uses the MH kernel's invariance and Dobrushin coefficient; the corrected lineages remain independent and identically distributed.Source: A Proofs for Section 3 · pp. 18–19

B.1 LightSB fits

LightSB fits a Gaussian-mixture transport at each transition using the preceding level's weighted particle cloud, with transition-specific training effort.

B.1 LightSB fits

  • The transport is a diagonal-covariance Gaussian mixture fitted by minimizing the self-normalized objective on the level-(l−1) cloud.The cloud is whitened, and the fit starts from an identity-coupling initialization using anchors drawn from the weighted cloud.Source: §B.1 · p. 19
  • Training uses Adam with minibatches capped at 1024, gradient-norm clipping at 10, and a cosine learning-rate schedule.The number of iterations is transition-specific because most fitting error occurs in the first few transitions.Source: §B.1 · p. 19

B.2 PDE solves per method

The methods differ in their level-l PDE solve counts, with LightSB and TransportMH also requiring transport fitting and evaluation at transported particles.

B.2 PDE solves per method

  • TransportMH uses Nl−1 + Nl + nlNl level-l solves, adding MH sweeps to LightSB's Nl−1 + Nl solve count.MLSMC-pCN uses Nl−1 + nlNl solves; its mutation is absent at l = L.Source: §B.2 · p. 19
  • Transported particles require γl evaluations for the MH ratio and next-transition training weights, while LightSB and TransportMH incur stochastic-gradient fitting time without PDE solves.The PDE discretization and solver differ by test case; a level-4 solve in test case 3 takes about 0.4 s.Source: §B.2 · p. 19

B.3 Sample sizes and reference values

The appendix specifies level-dependent sample allocations and reference values for three test cases, including cases where higher-level corrections are below reference uncertainty.

B.3 Sample sizes and reference values

  • For test case 1 at L = 9, the allocation is (100000, 35356, 12500, 4420, 1563, 553, 196, 100, 100, 100), and the continuum reference is 50.33335 ± 0.00010.The reference uses self-normalized importance sampling with 2 × 10^8 continuum-model samples.Source: §B.3 · p. 20
  • Test case 2 scores L = 6 runs against 0.99780 ± 0.00112 because its level-6 correction is (1 ± 1) × 10^-5, far below the reference standard error.The reference comes from 40 independent runs of method (ii) at ten times the L = 5 allocation.Source: §B.3 · p. 20
  • For test case 3, the target is π3(φ) with reference 0.7981 ± 0.0029; its level-4 correction is at most 6 × 10^-4, below that standard error.The reference uncertainty contributes about a sixth of TransportMH's mean squared error at L = 4.Source: §B.3 · p. 20

B.4 Breakdown of the computing time

The runtime breakdown separates shared costs at the first transition from TransportMH’s extra fine-level solves. Fine-level MH sweeps have high acceptance but add substantial variance and expensive PDE work.

B.4 Breakdown of the computing time

  • 1.06 × 10^6 of 1.37 × 10^6 solves in test case 1 and 2.5 × 10^5 of 4.25 × 10^5 in test case 2 come from the 30 first-transition sweeps.These sweeps dominate the solve count of both ESS-adaptive methods; plain MLSMC-pCN avoids them but has a larger level-2 correction standard deviation in test case 1: 0.0061 versus 0.0011.Source: §B.4 · p. 20
  • TransportMH and MLSMC-pCN with ESS sweeps make similar solve counts, but TransportMH adds solves at the finest, most expensive levels.The counts are 1.43 versus 1.37 × 10^6 in test case 1 and 4.35 versus 4.25 × 10^5 in test case 2.Source: §B.4 · p. 20
  • TransportMH’s level-5 and level-6 transported-particle and MH solves take 74 s of its 119 s at L = 6 in test case 2.LightSB fitting costs 3.0 s (2.5%) in test case 2 and 3.0 s (2%) in test case 3; in test case 1 it costs 2.7 s, about a fifth of total time.Source: §B.4 · p. 21
  • At fine levels, MH acceptance is 94–99%, yet three sweeps multiply paired-difference variance by about four: 2.3 × 10^-6 to 9.4 × 10^-6 in test case 2 and 2.5 × 10^-6 to 9.7 × 10^-6 in test case 1.At coarse transitions, the sweeps improve accuracy; without them, LightSB has relative bias of −2 × 10^-3 and −3 × 10^-2.Source: §B.4 · p. 21
Loading 2610.29546v1…