Unbiased Approximations for Stationary Distributions of McKean-Vlasov SDEs

We consider the development of unbiased estimators, to approximate the stationary distribution of Mckean-Vlasov stochastic differential equations (MVSDEs). These are an important class of processes, which frequently appear in applications such as mathematical finance, biology and opinion dynamics. Typically the stationary distribution is unknown and indeed one cannot simulate such processes exactly. As a result one commonly requires a time-discretization scheme which results in a discretization bias and a bias from not being able to simulate the associated stationary distribution. To overcome this bias, we present a new unbiased estimator taking motivation from the literature on unbiased Monte Carlo. We prove the unbiasedness of our estimator, under assumptions. In order to prove this we require developing ergodicity results of various discrete time processes, through an appropriate discretization scheme, towards the invariant measure. Numerous numerical experiments are provided, on a range of MVSDEs, to demonstrate the effectiveness of our unbiased estimator. Such examples include the Currie-Weiss model, a 3D neuroscience model and a parameter estimation problem.

💡 Research Summary

The paper addresses the problem of estimating the invariant (stationary) distribution π of a McKean‑Vlasov stochastic differential equation (MVSDE) of the form

 dXₜ = a(Xₜ, ξ₁(Xₜ, μₜ)) dt + b(Xₜ, ξ₂(Xₜ, μₜ)) dWₜ,

where μₜ denotes the law of Xₜ and the coefficients depend on both the state and its distribution. In most applications (finance, biology, opinion dynamics) the stationary law is unknown and cannot be simulated exactly. Standard practice is to discretize the dynamics (e.g., Euler‑Maruyama) and run a long simulation, which introduces two sources of bias: (i) discretization bias, and (ii) bias due to the finite simulation horizon. The authors propose a fully unbiased Monte‑Carlo estimator that eliminates both biases.

Key ideas and methodology

1. Randomized Multilevel Monte Carlo (RMLMC) for MVSDEs – Inspired by the Rhee‑Glynn framework, the authors construct a telescoping sum over discretization levels l with mesh size Δₗ = 2⁻ˡ. At each level they run an Euler‑Maruyama particle system with Nₗ particles, producing an empirical measure μₗ,ᴺ(t). The estimator for the difference of expectations between consecutive levels is denoted ξₗ.

2. Coupling of consecutive levels – Algorithm 2.2 generates a coupled pair of particle systems at levels l and l‑1 using the same Brownian increments and a common random seed. This strong coupling is essential for variance reduction and for the telescoping identity to hold.

3. Unbiased estimator – Let L be a random level drawn from a prescribed probability mass function P_L. The final estimator is

 π̂(φ) = ξ_L / P_L(L),

where φ is any test function of interest. The authors prove that E