grf {geoR}R Documentation

Simulation of Gaussian Random Fields

Description

Generates simulations of Gaussian random fields for given covariance parameters.

Usage

grf(n, grid = "irreg", nx, ny, xlims = c(0, 1), ylims = c(0, 1),
    nsim = 1, cov.model = "matern",
    cov.pars = stop("cov. parameters (sigmasq and phi) needed"), 
    kappa = 0.5, nugget = 0, lambda = 1, aniso.pars,
    method = c("cholesky", "svd", "eigen", "circular.embedding"),
    messages.screen = TRUE)

Arguments

n number of points (spatial locations) in each simulations.
grid optional. An n x 2 matrix with coordinates of the simulated data.
nx optional. Number of points in the X direction.
ny optional. Number of points in the Y direction.
xlims optional. Limits of the area in the X direction. Defaults to [0,1].
ylims optional. Limits of the area in the Y direction. Defaults to [0,1].
nsim Number of simulations. Defaults to 1.
cov.model correlation function. See cov.spatial for further details. Defaults to the exponential model.
cov.pars a vector with 2 elements or an n x 2 matrix with values of the covariance parameters sigma^2 (partial sill) and phi (range parameter). If a vector, the elements are the values of sigma^2 and phi, respectively. If a matrix, corresponding to a model with several structures, the values of sigma^2 are in the first column and the values of phi are in the second.
kappa additional smoothness parameter required only for the following correlation functions: "matern", "powered.exponential", "cauchy" and "gneiting.matern". More details on the documentation for the function cov.spatial.
nugget the value of the nugget effect parameter tau^2.
lambda value of the Box-Cox transformation parameter. The value lambda = 1 corresponds to no transformation, the default. For any other value of lambda Gaussian data is simulated and then transformed.
aniso.pars geometric anisotropy parameters. By default an isotropic field is assumed and this argument is ignored. If a vector with 2 values is provided, with values for the anisotropy angle psi_A (in radians) and anisotropy ratio psi_A, the coordinates are transformed, the simulation is performed on the isotropic (transformed) space and then the coordinates are back-transformed such that the resulting field is anisotropic. Coordinates transformation is performed by the function coords.aniso.
method simulation method. Defaults to the Cholesky decomposition. See section DETAILS below.
messages.screen logical, indicating whether or not status messages are printed on the screen (or output device) while the function is running. Defaults to TRUE.

Details

For the methods "cholesky", "svd" and "eigen" the simulation consists of multiplying a vector of standardized normal deviates by a square root of the covariance matrix. The square root of a matrix is not uniquely defined. The three available methods differs in the way they compute the square root of the (positive definite) covariance matrix.

For method = "circular.embedding" the algorithm implements the method described by Wood & Chan (1994) which is based on Fourier transforms. Only regular and equally spaced grids can be generated using this method.
The code for the "circular.embedding" method was provided by Martin Schlather, University of Bayreuth
(http://btgyn8.geo.uni-bayreuth.de/~martin/).

WARNING: The code for the "circular.embedding" method is no longer being maintained. Martin has released a package called RandomFields (available on CRAN) for simulation of random fields. We strongly recommend the use of this package for simulations on fine grids with large number of locations.

Value

A list with the components:

coords an n x 2 matrix with the coordinates of the simulated data.
data a vector (if nsim = 1) or a matrix with the simulated values. For the latter each column corresponds to one simulation.
cov.model a string with the name of the correlation function.
nugget the value of the nugget parameter.
cov.pars a vector with the values of sigma^2 and phi, respectively.
kappa value of the parameter kappa.
lambda value of the Box-Cox transformation parameter lambda.
aniso.pars a vector with values of the anisotropy parameters, if provided in the function call.
method a string with the name of the simulation method used.
sim.dim a string "1d" or "2d" indicating the spatial dimension of the simulation.
.Random.seed the random seed by the time the function was called.
messages messages produced by the function describing the simulation.
call the function call.

Author(s)

Paulo Justiniano Ribeiro Jr. Paulo.Ribeiro@est.ufpr.br,
Peter J. Diggle p.diggle@lancaster.ac.uk.

References

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

Schlather, M. (1999) Introduction to positive definite functions and to unconditional simulation of random fields. Tech. Report ST–99–10, Dept Maths and Stats, Lancaster University.

Further information on the package geoR can be found at:
http://www.est.ufpr.br/geoR.

See Also

plot.grf and image.grf for graphical output, coords.aniso for anisotropy coordinates transformation and, chol, svd and eigen for methods of matrix decomposition.

Examples

#
sim1 <- grf(100, cov.pars = c(1, .25))
# a display of simulated locations and values
points(sim1)   
# empirical and theoretical variograms
plot(sim1)
#
# a "smallish" simulation
sim2 <- grf(441, grid = "reg", cov.pars = c(1, .25)) 
image(sim2, xlab="X Coord", ylab="Y Coord")
#
# a "bigger" one
sim3 <- grf(40401, grid = "reg", cov.pars = c(10, .2), met = "circ") 
image(sim3, xlab="X Coord", ylab="Y Coord")
##
## 1-D simulations using the same seed and different noise/signal ratios
##
sim11 <- grf(100, ny=1, cov.pars=c(1, 0.25), nug=0)
assign(".Random.seed",  sim11$.Random.seed, envir=.GlobalEnv)
sim12 <- grf(100, ny=1, cov.pars=c(0.75, 0.25), nug=0.25)
assign(".Random.seed",  sim11$.Random.seed, envir=.GlobalEnv)
sim13 <- grf(100, ny=1, cov.pars=c(0.5, 0.25), nug=0.5)
##
par.ori <- par(no.readonly = TRUE)
par(mfrow=c(3,1), mar=c(3,3,.5,.5))
yl <- range(c(sim11$data, sim12$data, sim13$data))
image(sim11, type="l", ylim=yl)
image(sim12, type="l", ylim=yl)
image(sim13, type="l", ylim=yl)
par(par.ori)

[Package Contents]