Gaussian Process simulation and circulant embeddings

I might shoehorn Whittle likelihoods in here too

2017-11-10 – 2017-11-10
Maybe related 🤷:

An converse problem to covariance estimation. Related: phase retrieval. Gaussian process regression.

TODO

Discuss in the context of

  • deterministic covariance
  • expected autocorrelation in covariance

Simulating Gaussian fields with the desired covariance structure

Following the introduction in (Dietrich and Newsam 1993):

Let’s say we wish to generate a stationary Gaussian process \(Y(x)\) on a points \(\Omega\). \(\Omega=(x_0, x_1,\dots, x_m)\).

Stationary in this context means that the covariance function \(r\) is translation-invariance and depend only on distance, so that it may be given \(r(|x|)\). Without loss of generality, we assume that \(\bb E[Y(x)]=0\) and \(\var[Y(x)]=1\).

The problem then reduces to generating a vector \(\vv y=(Y(x_0), Y(x_1), \dots, Y(x_m) )\sim \mathcal{N}(0, R)\) where \(R\) has entries \(R[p,q]=r(|x_p-x_q|).\)

Note that if \(\bb \varepsilon\sim\mathcal{N}(0, I)\) is an \(m+1\)-dimensional normal random variable, and \(AA^T=R\), then \(\vv y=\mm A\bb \varepsilon\) has the required distribution.

The circulant embedding trick

But what if that’s too slow? If we have additional structure, we can do better.

Suppose further that our points form a grid, \(\Omega=(x_0, x_0+h,\dots, x_0+mh)\); specifically, equally-spaced-points on a line.

We know that \(R\) has a Toeplitz structure. Moreover it is non-negative definite, with \(\vv x^t\mm R \vv x \geq 0\forall \vv x.\) (Why?) 🏗

Simulating point processes with the desired covariance structure

For now, see spatial point processes. 🏗

References

Alexanderian, Alen. 2015. “A Brief Note on the Karhunen-Loève Expansion.” October 26, 2015. http://arxiv.org/abs/1509.07526.
Chan, G., and A. T. A. Wood. 1999. “Simulation of Stationary Gaussian Vector Fields.” Statistics and Computing 9 (4): 265–68. https://doi.org/10.1023/A:1008903804954.
Choromanski, Krzysztof, and Vikas Sindhwani. 2016. “Recycling Randomness with Structure for Sublinear Time Kernel Expansions.” May 29, 2016. http://arxiv.org/abs/1605.09049.
Cotter, S. L., G. O. Roberts, A. M. Stuart, and D. White. 2013. “MCMC Methods for Functions: Modifying Old Algorithms to Make Them Faster.” Statistical Science 28 (3): 424–46. https://doi.org/10.1214/13-STS421.
Davies, Tilman M., and David Bryant. 2013. “On Circulant Embedding for Gaussian Random Fields in R.” Journal of Statistical Software 55 (9). https://doi.org/10.18637/jss.v055.i09.
Dietrich, C. R., and G. N. Newsam. 1993. “A Fast and Exact Method for Multidimensional Gaussian Stochastic Simulations.” Water Resources Research 29 (8): 2861–69. https://doi.org/10.1029/93WR01070.
Ellis, Robert L., and David C. Lay. 1992. “Factorization of Finite Rank Hankel and Toeplitz Matrices.” Linear Algebra and Its Applications 173 (August): 19–38. https://doi.org/10.1016/0024-3795(92)90420-F.
Graham, Ivan G., Frances Y. Kuo, Dirk Nuyens, Rob Scheichl, and Ian H. Sloan. 2017a. “Analysis of Circulant Embedding Methods for Sampling Stationary Random Fields.” October 2, 2017. http://arxiv.org/abs/1710.00751.
———. 2017b. “Circulant Embedding with QMC – Analysis for Elliptic PDE with Lognormal Coefficients.” October 25, 2017. http://arxiv.org/abs/1710.09254.
Gray, Robert M. 2006. Toeplitz and Circulant Matrices: A Review. Vol. 2. https://doi.org/10.1561/0100000006.
Guinness, Joseph, and Montserrat Fuentes. 2016. “Circulant Embedding of Approximate Covariances for Inference From Gaussian Data on Large Lattices.” Journal of Computational and Graphical Statistics 26 (1): 88–97. https://doi.org/10.1080/10618600.2016.1164534.
Heinig, Georg, and Karla Rost. 2011. “Fast Algorithms for Toeplitz and Hankel Matrices.” Linear Algebra and Its Applications 435 (1): 1–59. https://doi.org/10.1016/j.laa.2010.12.001.
Powell, Catherine E. 2014. “Generating Realisations of Stationary Gaussian Random Fields by Circulant Embedding.” Matrix 2 (2): 1.
Stroud, Jonathan R., Michael L. Stein, and Shaun Lysen. 2017. “Bayesian and Maximum Likelihood Estimation for Gaussian Processes on an Incomplete Lattice.” Journal of Computational and Graphical Statistics 26 (1): 108–20. https://doi.org/10.1080/10618600.2016.1152970.
Whittle, P. 1953a. “The Analysis of Multiple Stationary Time Series.” Journal of the Royal Statistical Society: Series B (Methodological) 15 (1): 125–39. https://doi.org/10.1111/j.2517-6161.1953.tb00131.x.
———. 1953b. “Estimation and Information in Stationary Time Series.” Arkiv för Matematik 2 (5): 423–34. https://doi.org/10.1007/BF02590998.
Whittle, Peter. 1952. “Some Results in Time Series Analysis.” Scandinavian Actuarial Journal 1952 (1-2): 48–60. https://doi.org/10.1080/03461238.1952.10414182.
Ye, Ke, and Lek-Heng Lim. 2016. “Every Matrix Is a Product of Toeplitz Matrices.” Foundations of Computational Mathematics 16 (3): 577–98. https://doi.org/10.1007/s10208-015-9254-z.