pois.krige: Conditional Simulation and Prediction for the Poisson Spatial...
In geoRglm: A Package for Generalised Linear Spatial Models

Description Usage Arguments Details Value Author(s) References See Also Examples

This function performs conditional simulation (by MCMC) and spatial prediction in the Poisson normal model with link function from the Box-Cox class for fixed covariance parameters. Available types of prediction are: sk (simple kriging; fixed beta), ok (ordinary kriging; uniform prior on beta).

1
2
3

pois.krige(geodata, coords = geodata$coords, data = geodata$data,
               units.m = "default", locations = NULL, borders,
               mcmc.input, krige, output)

`geodata`	a list containing elements `coords` and `data` as described next. Typically an object of the class `"geodata"` - a geoR data set. If not provided the arguments `coords` and `data` must be provided instead. The list may also contain an argument `units.m` as described below.
`coords`	an n x 2 matrix, each row containing Euclidean coordinates of the n data locations. By default it takes the element `coords` of the argument `geodata`.
`data`	a vector with data values. By default it takes the element `data` of the argument `geodata`.
`units.m`	n-dimensional vector of observation times for the data. By default (`units.m = "default"`), it takes `geodata$units.m` in case this exist and else a vector of 1's.
`locations`	an N x 2 matrix or data frame, or a list with the 2-D coordinates of the N prediction locations.
`borders`	optional. If a two column matrix defining a polygon is provided the prediction is performed only at locations inside this polygon.
`mcmc.input`	input parameter for the MCMC algorithm. It can take an output from `mcmc.control` or a list with elements as for the arguments in `mcmc.control`. See documentation for `mcmc.control`. ATTENTION: the argument `S.scale` is necessary while all the others have default values.
`krige`	defines the model components and the type of kriging. It can take an output from `krige.glm.control` or a list with elements as for the arguments in `krige.glm.control`. See documentation for `krige.glm.control`.
`output`	parameters for controlling the output. It can take an output from `output.glm.control` or a list with elements as for the arguments in `output.glm.control`. See documentation for `output.glm.control`.

For simulating the conditional distribution of S given y, the Langevin-Hastings algorithm with the parametrisation in Papaspilliopoulus, Roberts and Skold (2003) is used. This algorithm is a Metropolis-Hastings algorithm, where the proposal distribution uses gradient information from the log-posterior distribution.

The proposal variance (called S.scale; see mcmc.control) for the algorithm needs to be scaled such that approximately 60 percent of the proposals are accepted. We also recommend that the user to studies plots of the autocorrelations.

The prediction part of the program consist of performing trans-Gaussian kriging on each of the simulated g^{-1}(S)-“datasets” from the conditional distribution. Afterwards the predictor is obtained by taking the mean of prediction means, and the prediction variance is obtained by taking the mean of the prediction variances plus the variance of the prediction means. The trans-Gaussian kriging is done by calling an internal function which is an extension of krige.conv allowing for more than one “data set”, and using a second order Taylor approximation of the inverse Box-Cox transformation function g^{-1} when the transformation parameter lambda>0; for the exponential function, i.e. logarithmic link and lambda=0, an exact formula is used instead of the Taylor approximation.

A list with the following components:

`predict`	a vector with predicted values.
`krige.var`	a vector with predicted variances.
`mcmc.error`	estimated Monte Carlo errors on the predicted values.
`beta.est`	estimate of the beta parameter. Not included in the output if `type.krige = "sk"`.
`intensity`	an n x n.sim matrix with n.sim being the number of MCMC simulations. containing g^{-1}(S_i). Each column corresponds to a conditional simulation of the conditional distribution of g^{-1}(S_i) at the data locations. Only returned when no prediction locations are given.
`acc.rate`	matrix with acceptance rates from MCMC. Only returned when no prediction locations are given.
`simulations`	an ni x n.sim matrix where ni is the number of prediction locations and n.sim is the number of MCMC simulations. Each column corresponds to a conditional simulation of the predictive distribution g^{-1}(S^{})*. Only returned if `sim.predict = TRUE`.
`call`	the function call.

Ole F. Christensen OleF.Christensen@agrsci.dk,
Paulo J. Ribeiro Jr. Paulo.Ribeiro@est.ufpr.br.

O. Papaspiliopoulus and G. O. Roberts and M. Skold (2003). Non-centered parameterizations for hierarchical models and data augmentation. Bayesian statistics 7 (eds. J. M. Bernardo, S. Bayarri, J. O. Berger, A. P. Dawid, D. Heckerman, A. F. M. Smith and M. West), Oxford University Press, 307-326.

Further information about geoRglm can be found at:
http://gbi.agrsci.dk/~ofch/geoRglm.

pois.krige.bayes for Bayesian prediction in the Poisson-normal model, binom.krige for prediction with fixed parameters in the binomial-logit normal model, krige.conv for prediction in the linear Gaussian model.

if(!exists(".Random.seed", envir=.GlobalEnv, inherits = FALSE)) set.seed(1234)
data(p50)
# First we scale the algorithm, and study how well the chain is mixing.
test <- pois.krige(p50, krige = list(cov.pars = c(1,1), beta = 1),
      mcmc.input = mcmc.control(S.scale = 0.2, thin = 1))
plot(log(test$intensity[45,]), type = "l")
acf(log(test$intensity[45,]), type = "correlation", plot = TRUE)
## Not run: # Now we make prediction (we decide to thin to every 10, which is the default),
# where we now use S.scale = 0.55.
test2 <- pois.krige(p50, locations = cbind(c(0.5,0.5), c(1,0.4)), 
      krige = krige.glm.control(cov.pars = c(1,1), beta = 1), 
      mcmc.input = mcmc.control(S.scale = 0.55)) 
## End(Not run)