[Paper Review] Theoretical analysis of a Stochastic Approximation approach for computing Quasi-Stationary distributions

José Blanchet, Peter W. Glynn+1 more|arXiv (Cornell University)|Jan 2, 2014

Scientific Research and Discoveries5 references4 citations

TL;DR

This paper proposes a stochastic approximation framework for computing quasi-stationary distributions in Markov chains without explicit matrix manipulation, proving convergence and deriving a central limit theorem (CLT) for the estimator. It identifies scenarios where convergence is slow, introduces a projected and averaged variant using Polyak-Ruppert averaging, and demonstrates via numerical experiments that the improved algorithm achieves significantly faster convergence rates, especially in high-dimensional or poorly conditioned systems.

ABSTRACT

This paper studies a method, which has been proposed in the Physics literature by [8, 7, 10], for estimating the quasi-stationary distribution. In contrast to existing methods in eigenvector estimation, the method eliminates the need for explicit transition matrix manipulation to extract the principal eigenvector. Our paper analyzes the algorithm by casting it as a stochastic approximation algorithm (Robbins-Monro) [23, 16]. In doing so, we prove its convergence and obtain its rate of convergence. Based on this insight, we also give an example where the rate of convergence is very slow. This problem can be alleviated by using an improved version of the algorithm that is given in this paper. Numerical experiments are described that demonstrate the effectiveness of this improved method.

Motivation & Objective

• To provide a theoretical foundation for a heuristic sampling-based method proposed in physics literature to estimate quasi-stationary distributions in interacting particle systems.
• To establish convergence and a general central limit theorem (CLT) for the estimator, extending beyond specific functionals considered in prior work.
• To identify conditions under which the original algorithm suffers from very slow convergence, particularly when eigenvalues are poorly separated.
• To develop and analyze an improved algorithm using projection and iterate averaging (Polyak-Ruppert) that ensures a valid CLT under all scenarios.
• To validate the improved method through numerical experiments on loopy chains, M/M/1 queues, and the contact process.

Proposed method

• The algorithm is recast as a Robbins-Monro stochastic approximation process, where the update rule is derived from the Kolmogorov forward equation and empirical averages of pathwise observations.
• Convergence is proven using martingale limit theory and the ODE method, showing almost sure convergence of the empirical measure to the true quasi-stationary distribution.
• A central limit theorem (CLT) is established for the estimator under general conditions, not restricted to functionals related to non-principal eigenvectors.
• The improved algorithm applies projection to ensure boundedness and uses iterate averaging (Polyak-Ruppert) to reduce asymptotic variance and accelerate convergence.
• The continuous-time version is analyzed by deriving a corresponding ODE system with rate matrix $ Q $, and proving that $ \exp(-Q)^{-1} $ has non-negative entries using series expansion.
• Theoretical results rely on the Perron-Frobenius theorem for substochastic or rate matrices, ensuring existence and uniqueness of the principal left eigenvector.

Experimental results

Research questions

• RQ1 Under what conditions does the original stochastic approximation algorithm for quasi-stationary distributions converge, and what is its rate of convergence?
• RQ2 Why does the original algorithm exhibit very slow convergence in certain cases, and what structural properties of the transition matrix cause this?
• RQ3 Can a modified version of the algorithm be constructed to ensure a valid central limit theorem and faster convergence across all scenarios?
• RQ4 How does the performance of the improved algorithm compare to the original in high-dimensional or poorly conditioned systems?
• RQ5 Does the theoretical CLT hold for the improved algorithm in cases where it fails for the original, particularly in systems with eigenvalues near the spectral gap?

Key findings

• The original algorithm is formally recognized as a stochastic approximation process, enabling rigorous convergence and CLT analysis beyond prior heuristic or urn-process-based results.
• A general central limit theorem is established for the estimator, valid for all functionals, not just those tied to non-principal eigenvectors, significantly broadening applicability.
• The paper identifies specific scenarios—particularly when the spectral gap is small or eigenvalues are clustered—where the original algorithm's convergence rate becomes extremely slow.
• The improved algorithm, using projection and Polyak-Ruppert iterate averaging, achieves a faster rate of convergence and ensures a valid CLT even in problematic cases.
• Numerical experiments on loopy chains, M/M/1 queues, and the contact process confirm that the improved algorithm converges orders of magnitude faster than the original, especially in heavy-traffic or high-dimensional settings.
• The continuous-time version of the algorithm is shown to converge under similar conditions, with $ \exp(-Q)^{-1} $ proven to be non-negative, ensuring the validity of the ODE approximation.