Fast simulation of Gaussian random fields with flexible correlation models in Euclidean spaces
Abstract
The efficient simulation of Gaussian random fields with flexible correlation structures is fundamental in spatial statistics, machine learning, and uncertainty quantification. In this work, we revisit the \emph{spectral turning-bands} (STB) method as a versatile and scalable framework for simulating isotropic Gaussian random fields with a broad range of covariance models. Beyond the classical Matérn family, we show that the STB approach can be extended to two recent and flexible correlation classes that generalize the Matérn model: the Bummer-Tricomi model, which allows for polynomially decaying correlations and long-range dependence, and the Gauss-Hypergeometric model, which admits compactly supported correlations, including the Generalized Wendland family as a special case. We derive exact stochastic representations for both families: a Beta-prime mixture formulation for the Kummer-Tricomi model and complementary Beta- and Gasper-mixture representations for the Gauss-Hypergeometric model. These formulations enable exact, numerically stable, and computationally efficient simulation with linear complexity in the number of spectral components. Numerical experiments confirm the accuracy and computational stability of the proposed algorithms across a wide range of parameter configurations, demonstrating their practical viability for large-scale spatial modeling. As an application, we use the proposed STB simulators to perform parametric bootstrap for standard error estimation and model selection under weighted pairwise composite likelihood in the analysis of a large climate dataset.
Summary
This paper addresses the problem of efficiently simulating Gaussian random fields (GRFs) with flexible correlation structures, a crucial task in spatial statistics, machine learning, and uncertainty quantification. The authors revisit and extend the spectral turning-bands (STB) method, a scalable framework for simulating isotropic GRFs. Their key contribution is expanding the STB approach beyond the classical Matérn family to encompass two more general correlation classes: the Kummer-Tricomi (K) model, which allows for polynomially decaying correlations and long-range dependence, and the Gauss-Hypergeometric (GH) model, which admits compactly supported correlations. The paper derives exact stochastic representations for both the K and GH families, specifically a Beta-prime mixture for the K model and Beta/Gasper mixtures for the GH model. These representations enable numerically stable and computationally efficient simulation with linear complexity in the number of spectral components. The authors demonstrate the accuracy and computational stability of their proposed STB algorithms through numerical experiments across a wide range of parameter configurations, confirming their practical viability for large-scale spatial modeling. They further showcase the utility of these simulators by performing parametric bootstrap for standard error estimation and model selection using weighted pairwise composite likelihood in the analysis of a large climate dataset. This application highlights the importance of fast and accurate GRF simulation for statistical inference in real-world scenarios. The algorithms have been implemented in the R package GeoModels.
Key Insights
- •The paper provides exact stochastic representations for the Kummer-Tricomi (K) and Gauss-Hypergeometric (GH) covariance models, enabling efficient simulation using the spectral turning-bands (STB) method.
- •A Beta-prime mixture formulation is derived for the K model, while Beta and Gasper-mixture representations are developed for the GH model. These mixture representations allow for direct sampling from the spectral density, avoiding numerical integration.
- •The proposed STB algorithms achieve linear complexity O(NL) in the number of spatial locations (N) and spectral components (L), enabling large-scale simulations.
- •For the GH model, a universal auxiliary distribution is found that encapsulates the spectral behavior, independently of some model parameters. This allows decoupling sampling into a universal draw and a local rescaling.
- •Adaptive rejection sampling is used to efficiently sample from the universal auxiliary variable, achieving acceptance rates of 50-60% across typical parameter ranges.
- •The paper identifies simulatable and non-simulatable regions in the parameter space of the GH model, and provides alternative sampling strategies (Beta vs. Gasper mixtures) to address these limitations.
- •Numerical experiments validate the accuracy of the proposed STB algorithms by comparing empirical and theoretical semivariograms.
Practical Implications
- •The developed algorithms are directly applicable to spatial statistics, machine learning, and uncertainty quantification problems where efficient simulation of Gaussian random fields with flexible correlation structures is required.
- •Practitioners and engineers can use the GeoSimapprox function in the R package GeoModels to simulate Gaussian random fields with Matérn, Kummer-Tricomi, and Gauss-Hypergeometric covariance models.
- •The parametric bootstrap application demonstrates how these simulators can be used for standard error estimation and model selection in the analysis of large spatial datasets, particularly in climate science and environmental modeling.
- •The proposed algorithms can be extended to simulate non-Gaussian random fields derived from latent Gaussian random fields, opening up possibilities for modeling spatial extremes, environmental risk assessment, and subsurface modeling.
- •Future research could focus on extending these methods to non-isotropic covariance models and exploring applications in other domains, such as image processing and finance.