Become an expert in R — Interactive courses, Cheat Sheets, certificates and more!
Get Started for Free

RPcirculant

Circulant Embedding methods


Description

Circulant embedding is a fast simulation method for stationary (possibly anisotropic) Gaussian random fields on regular grids based on Fourier transformations. It is guaranteed to be an exact method for covariance functions with finite support, e.g. the spherical model. The method is admissible for any dimension apart from memory restrictions.
The simulation is performed on a torus which represents the bended grid. To remove wrong dependencies occuring at different borders of the grid which would be close on the torus, the simulation area is multiplied by a natural number. There is also a multivariate version of circulant embedding.

Cut-off embedding is a fast simulation method for stationary, isotropic Gaussian random fields on square lattices based on the standard RPcirculant method, so that exact simulation is garantueed for further covariance models, e.g. the RMwhittle model.

In fact, the circulant embedding is called with the cutoff hypermodel, see RMcutoff. Cutoff halfens the maximum number of elements models used to define the covariance function of interest (from 10 to 5).

Here, multiplicative models are not allowed (yet).
For details see RMcutoff.

Intrinsic embedding is a fast simulation method for intrinsically stationary, isotropic Gaussian random fields on square lattices based on the standard RPcirculant method, for further variogram models, e.g. RMfbm.

Note that the simulated random field is always non-stationary. In fact, the circulant embedding is called with the Intrinsic hypermodel, see RMintrinsic.

Here, multiplicative models are not allowed (yet).
For details see RMintrinsic.

Usage

RPcirculant(phi, boxcox, force, mmin, strategy,
 maxGB, maxmem,  tolIm, tolRe, trials, useprimes, dependent,
 approx_step, approx_maxgrid)

RPcutoff(phi, boxcox, force, mmin, strategy,
 maxGB, maxmem, tolIm, tolRe, trials, useprimes, dependent,
 approx_step, approx_maxgrid, diameter, a) 
 

RPintrinsic(phi, boxcox, force, mmin, strategy,
 maxGB, maxmem, tolIm, tolRe, trials, useprimes, dependent,
 approx_step, approx_maxgrid, diameter, rawR)

Arguments

phi

See RPgauss.

boxcox

the one or two parameters of the box cox transformation. If not given, the globally defined parameters are used. See RFboxcox for details.

force

Logical. Circulant embedding does not work if the constructed circulant matrix has negative eigenvalues. Sometimes it is convenient to replace all the negative eigenvalues by zero (force=TRUE) after trials number of trials. Default: FALSE.

mmin

Scalar or vector, integer if positive. CE.mmin determines the initial size of the circulant matrix. If CE.mmin=0 the minimal starting size is determined automatically according to the dimensions of the grid. If CE.mmin>0 then the absolute starting size is given. If CE.mmin<0 then the automatically determined matrix size is multiplied by |\code{CE.mmin}|; here, CE.mmin must be smaller than -1; the value -1 takes over the minimal starting size.
Note: in any cases, the initial size might be increased according to CE.useprimes.
Default: 0.

strategy

Logical. 0: If the circulant matrix has negative eigenvalues then the size in each direction is doubled;
1: The size is enhanced only in one direction, namely that one where the covariance function has the largest value at the end point of the grid — note that the default value of trials is probably too small in that case.

In some cases strategy=0 works better, in other cases strategy=1. Just try. Clearly, if the field is isotropic and a square grid should be simulated, then strategy=0 is the better choice.

Default: 0.

maxGB

Maximal memory used for the circulant matrix in units of GB. If this argument is set then maxmem is set to MAXINT.

Default: 1.

maxmem

Integer. maximal number of entries in a row of the circulant matrix. The total amount of memory needed for the internal calculations is 32 (= 4 * sizeof(double)) time as large (factor 2 is needed as complex numbers must be considered for calculating the fft of the covariance matrix; another factor 2 is needed for storing the simulated result).

The value of maxmem must be at least 2^d times as large as the number of points to be simulated. Here, d is the space dimension. In some cases even much larger.

Note that maxmem can be used to control the automatic choice of the simulation algorithm. Namely, in case of huge circulant matrices, other simulation methods (TBM) might be faster and might be preferred by the user.

If this argument is set then maxGB is set to Inf.

Default: MAXINT.

tolIm

If the modulus of the imaginary part is less than tolIm then the eigenvalue is always considered as real (independently of the value of force).

Default: 1E-3.

tolRe

Eigenvalues between tolRe and 0 are always considered as 0 and set 0 (independently of the value of force).

Default: -1E-7.

trials

Integer. A larger circulant matrix is likely to make more eigenvalues non-negative. If at least one of the thresholds tolRe and tolIm are missed then the matrix size is doubled according to strategy, and the matrix is checked again. This procedure is repeated up to trials - 1 times. If there are still negative eigenvalues, the simulation method fails if force=FALSE.

Default: 3.

useprimes

Logical. If FALSE the columns of the circulant matrix have length 2^k for some k. Otherwise the algorithm tries to find a nicely factorizable number close to the size of the given matrix.

Default: TRUE.

dependent

Logical. If FALSE then independent random fields are created. If TRUE then at least 4 non-overlapping rectangles are taken out of the the expanded grid defined by the circulant matrix. These simulations are dependent. See RFoptionsAdvanced for an example. See trials for some more information on the circulant matrix.

Default: FALSE.

approx_step

Real value. It gives the grid size of the approximating grid in case circulant embedding is used although the points do not lie on a grid.

If NA then approx_step is chosen such that approx_maxgrid is nearly reached.

Default: NA.

approx_maxgrid

It defaults to maxmem.

diameter

See RMcutoff or RMintrinsic.

a

See RMcutoff.

rawR

See RMintrinsic.

Details

Here, the algorithms by Dietrich and Newsam are implemented.

Value

An object of class RMmodel.

Author(s)

References

Circulant Embedding

  • Chan, G. and Wood, A.T.A. (1997) An Algorithm for Simulating Stationary Gaussian Random Fields. Journal of the Royal Statistical Society: Series C (Applied Statistics), 46, 171–181.

  • Dietrich, C. R. and G. N. Newsam (1993) A fast and exact method for multidimensional gaussian stochastic simulations. Water Resour. Res. 29(8), 2861–2869.

  • Dietrich, C. R. and G. N. Newsam (1996) A fast and exact method for multidimensional Gaussian stochastic simulations: Extension to realizations conditioned on direct and indirect measurements Water Resour. Res. 32(6), 1643–1652.

  • Dietrich, C. R. and Newsam, G. N. (1997) Fast and Exact Simulation of Stationary Gaussian Processes through Circulant Embedding of the Covariance Matrix. SIAM J. Sci. Comput. 18, 1088–1107.

  • Wood, A. T. A. and Chan, G. (1994) Simulation of Stationary Gaussian Processes in [0, 1]^d. Journal of Computational and Graphical Statistics 3, 409–432.

Cutoff and Intrinsic

  • Gneiting, T., Sevecikova, H, Percival, D.B., Schlather M., Jiang Y. (2006) Fast and Exact Simulation of Large Gaussian Lattice Systems in $R^2$: Exploring the Limits. J. Comput. Graph. Stat. 15, 483–501.

  • Stein, M.L. (2002) Fast and exact simulation of fractional Brownian surfaces. J. Comput. Graph. Statist. 11, 587–599

See Also

Examples

RFoptions(seed=0) ## *ANY* simulation will have the random seed 0; set
##                   RFoptions(seed=NA) to make them all random again

model <- RMstable(s=1, alpha=1.8)
x <- seq(-3,3,0.1)

z <- RFsimulate(model=RPcirculant(model), x=x, y=x, n=1)
plot(z)

model <- RMexp(var=10, s=2)
z <- RFsimulate(model=RPcirculant(model), 1:10)
plot(z)

model <- RMfbm(Aniso=diag(c(1,2)), alpha=1.5)
z <- RFsimulate(model, x=1:10, y=1:10)
plot(z)

RandomFields

Simulation and Analysis of Random Fields

v3.3.10
GPL (>= 3)
Authors
Martin Schlather [aut, cre], Alexander Malinowski [aut], Marco Oesting [aut], Daphne Boecker [aut], Kirstin Strokorb [aut], Sebastian Engelke [aut], Johannes Martini [aut], Felix Ballani [aut], Olga Moreva [aut], Jonas Auel[ctr], Peter Menck [ctr], Sebastian Gross [ctr], Ulrike Ober [ctb], Paulo Ribeiro [ctb], Brian D. Ripley [ctb], Richard Singleton [ctb], Ben Pfaff [ctb], R Core Team [ctb]
Initial release

We don't support your browser anymore

Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.