A Krylov Subspace Method for Covariance Approximation and Simulation of Random Processes and Fields

This paper proposes a new iterative algorithm for simultaneously computing an approximation to the covariance matrix of a random vector and drawing a sample from that approximation. The algorithm is especially suited to cases for which the elements of the random vector are samples of a stochastic process or random field. The proposed algorithm has close connections to the conjugate gradient method for solving linear systems of equations. A comparison is made between our algorithm's structure and complexity and other methods for simulation and covariance matrix approximation, including those based on FFTs and Lanczos methods. The convergence of our iterative algorithm is analyzed both analytically and empirically, and a preconditioning technique for accelerating convergence is explored. The numerical examples include a fractional Brownian motion and a random field with the spherical covariance used in geostatistics.

[1]  Thomas Kailath,et al.  Fast subspace decomposition , 1994, IEEE Trans. Signal Process..

[2]  Thomas Kailath,et al.  Application of fast subspace decomposition to signal processing and communication problems , 1994, IEEE Trans. Signal Process..

[3]  H. V. D. Vorst,et al.  An iterative solution method for solving f ( A ) x = b , using Krylov subspace information obtained for the symmetric positive definite matrix A , 1987 .

[4]  Michael K. Schneider,et al.  Krylov Subspace Estimation , 2000, SIAM J. Sci. Comput..

[5]  S. Mallat,et al.  Adaptive covariance estimation of locally stationary processes , 1998 .

[6]  B. Mandelbrot,et al.  Fractional Brownian Motions, Fractional Noises and Applications , 1968 .

[7]  C. R. Dietrich,et al.  Fast and Exact Simulation of Stationary Gaussian Processes through Circulant Embedding of the Covariance Matrix , 1997, SIAM J. Sci. Comput..

[8]  W. W. Irving,et al.  Multiscale stochastic realization and model identification with applications to large-scale estimation problems , 1995 .

[9]  Ricardo Todling,et al.  Suboptimal Schemes for Retrospective Data Assimilation Based on the Fixed-Lag Kalman Smoother , 1998 .

[10]  Austin B. Frakt.  Internal multiscale autoregressive processes, stochastic realization, and covariance extension , 1999 .

[11]  G. Evensen,et al.  Analysis Scheme in the Ensemble Kalman Filter , 1998 .

[12]  D. McLaughlin,et al.  Hydrologic Data Assimilation with the Ensemble Kalman Filter , 2002 .

[13]  Christian L. Keppenne,et al.  Data Assimilation into a Primitive-Equation Model with a Parallel Ensemble Kalman Filter , 2000 .

[14]  T. Kailath,et al.  Fast Estimation of Principal Eigenspace Using LanczosAlgorithm , 1994 .

[15]  James Demmel,et al.  Applied Numerical Linear Algebra , 1997 .

[16]  L. Knizhnerman,et al.  Extended Krylov Subspaces: Approximation of the Matrix Square Root and Related Functions , 1998, SIAM J. Matrix Anal. Appl..

[17]  C. Cabos.  Evaluation of matrix functions with the block Lanczos algorithm , 1997 .

[18]  R. Todling,et al.  Suboptimal Schemes for Atmospheric Data Assimilation Based on the Kalman Filter , 1994 .

[19]  D. Sorensen.  Numerical methods for large eigenvalue problems , 2002, Acta Numerica.

[20]  Gene H. Golub,et al.  Fast algorithms for updating signal subspaces , 1994 .

[21]  H. Vincent Poor,et al.  Signal detection in fractional Gaussian noise , 1988, IEEE Trans. Inf. Theory.

[22]  G. Evensen.  Sequential data assimilation with a nonlinear quasi‐geostrophic model using Monte Carlo methods to forecast error statistics , 1994 .

[23]  Y. Saad.  On the Rates of Convergence of the Lanczos and the Block-Lanczos Methods , 1980 .

[24]  B. Parlett,et al.  The Lanczos algorithm with selective orthogonalization , 1979 .