gear (version 0.1.5)

predict.geolmMan: Predict method for geostatistical models

Description

predict calculates the predicted values at specified locations. The method can additionally provide the mean square prediction error (mspe) and perform conditional simulation.

Usage

# S3 method for geolmMan
predict(object, newdata, nsim = 0, vop, vp, sp = TRUE, dmethod = "chol", ...)

Arguments

object

An object produced by the geolm function.

newdata

An optional data frame in which to look for the coordinates at which to predict. If omitted, the observed data locations are used.

nsim

A non-negative integer indicating the number of realizations to sample at the specified coordinates using conditional simulation.

vop

The cross-covariance matrix between the observed responses and the responses to predict.

vp

The covariance matrix of the responses to predict.

sp

A logical value indicating whether to object returned should be of class SpatialPointsDataFrame for easier plotting with the sp package. Default is TRUE.

dmethod

The method used to decompose the covariance matrix for conditional simulation. Valid options are "chol", "eigen", and "svd". The default is "chol".

...

Currently unimplemented.

Value

If sp = TRUE, then a SpatialPointDataFrame from the sp package is returned, with components including the prediction coordinates, the predicted responses pred, the mean square prediction error (mspe), the root mean square prediction error (rmspe), and the conditional realizations, if application (sim.1, sim.2, ...). If sp = FALSE, then a list of class gearKrige is returned, with components pred, mspe, rmspe, and sim, if relevant.

Details

The newdata data frame must include the relevant covariates for the prediction locations, where the covariates are specified on the right side of the ~ in object$formula. newdata must also include the coordinates of the prediction locations, with these columns having the names provided in object$coordnames.

Examples

# NOT RUN {
# generate response
y = rnorm(10)
# generate coordinates
x1 = runif(10); x2 = runif(10)

# data frame for observed data
data = data.frame(y, x1, x2)
coords = cbind(x1, x2)
d = as.matrix(dist(coords))
psill = 2 # partial sill
r = 4 # range parameter
evar = .1 # error variance
fvar = .1 # add finescale variance
# one can't generally distinguish between evar and fvar, but
# this is done for illustration purposes

# manually specify a an expoential covariance model 
v = psill * exp(-d/r) + (evar + fvar) * diag(10)

cmod_man = cmod.man(v = v, evar = evar)

#' # geolm for universal kriging
gearmod_uk = geolm(y ~ x1 + x2, data = data,
                 coordnames = c("x1", "x2"),
                 cmod = cmod_man)

# newdata must have columns with prediction coordinates
# add 5 unsampled sites to sampled sites
newdata = data.frame(x1 = c(x1, runif(5)), x2 = c(x2, runif(5)))
newcoords = newdata[,c("x1", "x2")]
# create vop and vp using distances
dop = sp::spDists(as.matrix(coords), as.matrix(newcoords))
dp = as.matrix(dist(newcoords))

vop = psill * exp(-dop/r) + fvar * (dop == 0)
vp = psill * exp(-dp/r) + fvar * diag(nrow(newcoords))

# prediction for universal kriging, with conditional simulation, 
# using manual covariance matrices
pred_uk_man = predict(gearmod_uk, newdata, nsim = 2, vop = vop, vp = vp,
                      dmethod = "svd")

# do the same thing, but using the std covariance function

# prediction for universal kriging, with conditional simulation
cmod_std = cmod.std("exponential", psill = psill, r = r, evar = evar, fvar = fvar)
gearmod_uk2 = geolm(y ~ x1 + x2, data = data, coordnames = c("x1", "x2"), 
                    cmod = cmod_std)
pred_uk_std = predict(gearmod_uk2, newdata, nsim = 2, dmethod = "svd")

# compare results
range(pred_uk_man$pred - pred_uk_std$pred)
range(pred_uk_man$mspe - pred_uk_std$mspe)
# }