unbiasedKrige: Unbiased kriging
In jskoien/intamap: Procedures for Automated Interpolation

Description Usage Arguments Details Value Author(s) References Examples

unbiasedKrige is a function for modifying a kriging prediction to a prediction that can be assumed to be unbiased for a certain threshold.

1
2
3

unbiasedKrige(object, formulaString, observations, predictionLocations,
              model, outputWhat, yamamoto, iwqmaxit = 500,
              iwqCpAddLim = 0.0001, debug.level, ...)

`object`	either an object of the intamap type (see `intamap-package` for further description of the necessary content of this variable) or the output from the function `krige` in `gstat`. If the object is a result from the intamap procedure spatialPredict, the remaining arguments are not necessary.
`formulaString`	formula that defines the dependent variable as a linear model of independent variables; suppose the dependent variable has name z, for ordinary and simple kriging use the formula z~1; for universal kriging, suppose z is linearly dependent on x and y, use the formula z~x+y
`observations`	a `Spatial`*DataFrame with observations; should contain the dependent variable, independent variables, and coordinates
`predictionLocations`	the predictionLocations, only necessary if the method is "IWQSEL" and formulaString contains independent variables. Should preferentally be a grid if the method is "IWQSEL"
`model`	variogram model of dependent variable (or its residuals), defined by a call to `vgm` or `autofitVariogram`
`outputWhat`	Argument with type of unbiasedness method ("MOK" or "IWQSEL") and the thresholds.
`yamamoto`	logical describing if the yamamoto approach )is to be used in simulations. Defaults to yamamoto = FALSE when object is a `Spatial`*DataFrame.
`iwqmaxit`	maximum number of iterations in iwqsel
`iwqCpAddLim`	convergence criteria in iwqsel
`debug.level`	debug level, passed to subfunctions
`...`	other arguments that will be passed to subfunctions. These include `nsim` number of simulations necessary if the method is "IWQSEL". Defaults to nsim = 100 when object is a `Spatial`DataFrame. `maxdist` maximum number of neighbours to use in local kriging, defaults to Inf `nmax` for local kriging: the number of nearest observations that should be used in simulations for the "IWQSEL" method in terms of the space of the spatial locations. Defaults to nmax = 10 when object is a SpatialDataFrame.

It is a fact that predictions from kriging tend to be biased towards the mean of the process. The function unbiasedKrige is a function that adds one or more predictions to the original output, which are assumed to be unbiased relative to a certain threshold. The two methods supported are the IWQSEL-method (Craigmile, 2006) and MOK (Skoien et al, 2008).

an object of type intamap, as described in intamap-package, or a Spatial*DataFrame with one or more new prediction columns, representing different methods and thresholds.

Jon Olav Skoien

Craigmile, P. F., N. Cressie, T. J. Santner, and Y. Rao. 2006. A loss function approach to identifying environmental exceedances. Extremes, 8, 143-159.

Skoien, J. O., G. B. M. Heuvelink, and E. J. Pebesma. 2008. Unbiased block predictions and exceedance probabilities for environmental thresholds. In: J. Ortiz C. and X. Emery (eds). Proceedings of the eight international geostatistics congress. Gecamin, Santiago, Chile, pp. 831-840.

Pebesma, E., Cornford, D., Dubois, G., Heuvelink, G.B.M., Hristopulos, D., Pilz, J., Stohlker, U., Morin, G., Skoien, J.O. INTAMAP: The design and implementation f an interoperable automated interpolation Web Service. Computers and Geosciences 37 (3), 2011.

library(automap)
library(gstat)
data(meuse)
data(meuse.grid)
coordinates(meuse) = ~x+y
gridded(meuse.grid) = ~x+y

predictionLocations = meuse.grid[sample(1:length(meuse.grid),50),]
vmod = autofitVariogram(log(zinc)~1,meuse)$var_model
prediction = krige(log(zinc)~1,meuse,predictionLocations,vmod)
summary(prediction)

prediction <- unbiasedKrige(prediction,log(zinc)~1,
  meuse,  model = vmod, outputWhat = list(MOK = 6.0, MOK = 7.0, IWQSEL=7.0), 
  iwqmaxit = 100, iwqCpAddLim = 0.01)
summary(prediction)