| pois.krige {geoRglm} | R Documentation |
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).
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 |
an |
data |
a vector with data values. By default it takes the
element |
units.m |
|
locations |
an |
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 |
krige |
defines the model components and the type of
kriging. It can take an output from |
output |
parameters for controlling the output. It can take an output from |
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 |
intensity |
an |
acc.rate |
matrix with acceptance rates from MCMC. Only returned when no prediction locations are given. |
simulations |
an |
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)