Using linear solvers to sample large Gaussians

Document Type

Presentation Abstract

Presentation Date

10-19-2009

Abstract

Generating samples from large multivariate normal (Gaussian) distributions is useful for simulating Gaussian Processes and Gauss-Markov Random Fields, which are commonly used in conjunction with Markov Chain Monte Carlo methods. Techniques from numerical linear algebra to solve linear systems are the same methods which are used to produce samples from Gaussian distributions. For example, the Cholesky factorization, the preferred method of solving linear systems for moderately sized problems, is the conventional way to produce samples from a Gaussian. For linear systems and Gaussian models with dimension 106 or more (eg models of global total column ozone, tropical ocean surface winds, or the structure of the Earth's mantel and outer core), iterative linear solvers and iterative samplers are the only feasible option due to their inexpensive cost per iteration, and small computer memory requirements. Motivated by ample examples, I will provide an explicit recipe which shows how to convert any stationary linear solver into an iterative sampler of a multivariate normal distribution, which has the same convergence properties as the corresponding linear solver.

The last fifty years has seen an explosion of theoretical results and algorithmic development making linear solvers faster and more efficient, so that for large problems, stationary processes are used as pre-conditioners at best, while the method of conjugate gradients and polynomial accelerators are the current state-of-the-art for solving linear systems in a finite number of steps. Perhaps less well known is that iterative samplers are merely stationary processes, which were used as very slow (ie geometrically converging) linear solvers in the 1950's.

I will show how iterative samplers can be sped up appreciably by using common acceleration techniques from numerical linear linear algebra such as successive over-relaxation (SOR) and polynomial preconditioning. Some samplers which will be presented are derived from Chebyshev accelerated symmetric-SOR, the method of conjugate gradients, and the Lanczos method for estimating the extreme eigenvalues of a matrix.

Additional Details

Monday, 19 October 2009
3:10 p.m. in Math 103
4:00 p.m. Refreshments in Math Lounge 109

This document is currently not available here.

Share

COinS