GeoModels (version 2.2.1)

GeoSimcond: Conditional (fast) simulation of Gaussian and non Gaussian Random Fields.

Description

Simulates spatial (spatio-temporal) Gaussian and non-Gaussian random fields conditioned on observed data using specified correlation models. For large conditioning datasets approximation methods like turning band can be used. For large locations to simlate local kriging can be used.

Usage

GeoSimcond(estobj = NULL, data, coordx, coordy = NULL, coordz = NULL, coordt = NULL,
           coordx_dyn = NULL, corrmodel, distance = "Eucl", grid = FALSE, loc,
           maxdist = NULL, maxtime = NULL, method = "Cholesky", model = "Gaussian",
           n = 1, nrep = 1, local = FALSE, L = 1000, neighb = NULL,
           param, anisopars = NULL, radius = 1, sparse = FALSE, time = NULL,
           copula = NULL, X = NULL, Xloc = NULL, Mloc = NULL,
           parallel=FALSE, ncores = NULL, progress=FALSE)

Value

Returns an object of class GeoSimcond containing:

• Simulated field realizations

• Model parameters and settings

• Spatial/temporal coordinates

• Covariance information

Arguments

estobj

Object of class Geofit containing model information

data

Numeric vector/matrix/array of observed data

coordx

Numeric matrix of spatial coordinates (d x 2 or d x 3)

coordy

Optional numeric vector of y-coordinates

coordz

Optional numeric vector of z-coordinates

coordt

Optional numeric vector of temporal coordinates

coordx_dyn

Optional list of dynamic spatial coordinates

corrmodel

String specifying correlation model name

distance

String specifying distance metric (default: "Eucl")

grid

Logical for regular grid (default: FALSE)

loc

Numeric matrix of prediction locations (n x 2)

maxdist

Optional maximum distance for local kriging

maxtime

Optional maximum temporal distance

method

String for decomposition method ("Cholesky", "TB", or "CE")

model

String specifying random field type (default: "Gaussian")

n

Number of trials for binomial RFs (default: 1)

nrep

Number of independent replicates (default: 1)

local

Logical for local kriging (default: FALSE)

L

Number of lines for turning bands method (default: 1000)

neighb

Optional neighborhood order for local kriging

param

List of parameter values

anisopars

List with anisotropy angle and ratio

radius

Radius for spherical coordinates (default: Earth's radius)

sparse

Logical for sparse matrix algorithms (default: FALSE)

time

Optional vector of temporal instants to predict

copula

Optional string specifying copula type

X

Optional matrix of spatio-temporal covariates

Xloc

Optional matrix of covariates for prediction locations

Mloc

Optional vector of estimated means for prediction locations

parallel

If TRUE the the computation is parallelized

ncores

Numbers of cores involved in the parallelization

progress

If TRUE then a progress bar is shown.

Author

Moreno Bevilacqua ,\ Víctor Morales Oñate ,\ Christian Caamaño-Carrillo

Details

For Gaussian RF, performs conditional simulation using three steps:

1. Unconditional simulation at observed and prediction locations

2. Simple kriging estimates at observed locations

3. Combination to produce conditional simulations

For large datasets, approximate methods ("TB" or "CE") are recommended coupled with local kriging (local=TRUE and neighb=k) and using parallelization (parallel=T).

References

Gaetan, C. and Guyon, X. (2010) Spatial Statistics and Modelling. Springer Verlag, New York.

See Also

GeoSim, GeoKrig

Examples

library(GeoModels)

##############################################
## conditional simulation of a Gaussian rf ###
##############################################
model="Gaussian"
set.seed(79)
### conditioning locations
x = runif(250, 0, 1)
y = runif(250, 0, 1)
coords=cbind(x,y)

# Set the exponential cov parameters:
corrmodel = "GenWend"
mean=0; sill=1; nugget=0
scale=0.2;smooth=0;power2=4

param=list(mean=mean,sill=sill,nugget=nugget,scale=scale,smooth=smooth,power2=power2)

# Simulation 
data = GeoSim(coordx=coords, corrmodel=corrmodel,model=model,
              param=param)$data

## estimation with pairwise likelihood
fixed=list(nugget=nugget,smooth=smooth,power2=power2)
start=list(mean=0,scale=scale,sill=1)
I=Inf
lower=list(mean=-I,scale=0,sill=0)
upper=list(mean= I,scale=I,sill=I)
# Maximum pairwise likelihood fitting :
fit = GeoFit(data, coordx=coords, corrmodel=corrmodel,model=model,
             likelihood='Marginal', type='Pairwise',neighb=3,
             optimizer="nlminb", lower=lower,upper=upper,
             start=start,fixed=fixed)

# locations to simulate 
xx=seq(0,1,0.025)
loc_to_sim=as.matrix(expand.grid(xx,xx))

# Conditional simulation
sim_result <- GeoSimcond(fit,loc = loc_to_sim,nrep=50)

cond_mean=sim_result$cond_mean  # conditional mean
cond_var =sim_result$cond_var    # conditional var

par(mfrow=c(1,3))
quilt.plot(coords, data)
quilt.plot(loc_to_sim, cond_mean)
quilt.plot(loc_to_sim, cond_var)
par(mfrow=c(1,1))

##############################################
## conditional simulation of a LogGaussian rf 
##############################################
model="LogGaussian"
set.seed(79)
### conditioning locations
x = runif(500, 0, 1)
y = runif(500, 0, 1)
coords=cbind(x,y)

# Set the exponential cov parameters:
corrmodel = "Matern"
mean=0; sill=.1; nugget=0
scale=0.2;smooth=0.5

param=list(mean=mean,sill=sill,nugget=nugget,scale=scale,smooth=smooth)

# Simulation 
data = GeoSim(coordx=coords, corrmodel=corrmodel,model=model,
              param=param)$data

## estimation with pairwise likelihood
fixed=list(nugget=nugget,smooth=smooth)
start=list(mean=0,scale=scale,sill=1)
I=Inf
lower=list(mean=-I,scale=0,sill=0)
upper=list(mean= I,scale=I,sill=I)
# Maximum pairwise likelihood fitting :
fit = GeoFit(data, coordx=coords, corrmodel=corrmodel,model=model,
             likelihood='Marginal', type='Pairwise',neighb=3,
             optimizer="nlminb", lower=lower,upper=upper,
             start=start,fixed=fixed)

# locations to simulate 
xx=seq(0,1,0.025)
loc_to_sim=as.matrix(expand.grid(xx,xx))

# Conditional simulation
sim_result <- GeoSimcond(fit,loc = loc_to_sim,nrep=50)

cond_mean=sim_result$cond_mean  # conditional mean
cond_var =sim_result$cond_var    # conditional var

par(mfrow=c(1,3))
quilt.plot(coords, data)
quilt.plot(loc_to_sim,cond_mean)
quilt.plot(loc_to_sim,cond_var)
par(mfrow=c(1,1))