Description Usage Arguments Details Value Author(s) References See Also Examples
This function performs posterior simulation (by MCMC) and spatial prediction in the binomiallogit spatial model.
1 2 3  binom.krige.bayes(geodata, coords = geodata$coords, data = geodata$data,
units.m = "default", locations = "no", borders,
model, prior, mcmc.input, output)

geodata 
a list containing elements 
coords 
an n x 2 matrix, each row containing Euclidean
coordinates of the n data locations. By default it takes the
element 
data 
a vector with data values. By default it takes the
element 
units.m 
ndimensional vector giving the number of trials for the data.
By default ( 
locations 
an N x 2 matrix or data frame, or a list with the 2D 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. 
model 
a list defining the components of the model. It can take an output from 
prior 
specification of priors for the model parameters.
It can take an output from 
mcmc.input 
input parameter for the MCMC algorithm. It can take an output from 
output 
parameters for controlling the output. It can take an output from 
binom.krige.bayes
is a function for Bayesian geostatistical
inference in the binomiallogit spatial model.
The Bayesian algorithm is using a discretized version of the prior distribution for the parameter phi. This means that the prior for phi is always a proper prior.
For simulating from the posterior distribution of S given y, we use a LangevinHastings type algorithm. This algorithm is a MetropolisHastings algorithm, where the proposal distribution uses gradient information from the posterior. The algorithm is described below. For shortness of presentation, we only present the MCMC algorithm for the case where beta follows a uniform prior.
When beta follows a uniform prior and the prior for sigma^2 is a scaled inversechi^2 distribution, the marginalised improper density of S is
f_I(s) propto D^T*V^{1}*D^{1/2}*V^{1/2}*{n.sigma*S2.sigma + s^T*(V^{1}V^{1}*D*(D^T*V^{1}*D)^{1}*D^T*V^{1})*s }^{(np+n.sigma)/2},
where V is the correlation matrix of S. The uniform prior for sigma^2 corresponds to S2.sigma=0 and n.sigma=2, and the reciprocal prior for sigma^2 corresponds to S2.sigma=0 and n.sigma=0.
We use the reparametrisation S = Q*Gamma, where Q is the Cholesky factorisation of V so that V=QQ^T. Posterior simulations of S are obtained by transforming MCMC simulations from the conditional distribution of Gamma given Y=y.
The log posterior density of Gamma given Y=y is
log f(gammay) = const(y)  1/2 * gamma^T*(I_n  V^{1/2}*D*(D^T*V^{1}*D)^{1}*D^T*V^{1/2})*gamma +∑_i y_i*s_i log(1+exp(s_i)),
where (s_1,...,s_n)^T = Q*gamma.
For the LangevinHastings update we use the gradient of the log target density,
nabla(gamma )^{trunc}=  (I_n  Q^{1}*D*(D^T*V^{1}*D)^{1}*D^T*(Q^{1})^T)*gamma + Q^T*{y_i\exp(s_i)/(1+exp(s_i))}_{i=1}^n .
The proposal gamma' follows a multivariate normal distribution with mean vector
xi(gamma)=gamma+(h/2)*nabla(gamma)^{trunc} and covariance matrix hI,
where h>0 is a userspecified “proposal variance” (called
S.scale
; see mcmc.control
).
When phi.prior
is not "fixed"
, we update the parameter phi by a random walk Metropolis step.
Here mcmc.input$phi.scale
(see mcmc.control
) is the proposal variance, which needs to be sufficiently large so that
the algorithm easily can move between the discrete values in prior$phi.discrete
(see prior.glm.control
).
CONTROL FUNCTIONS
The function call includes auxiliary control functions which allows
the user to specify and/or change the specification of 1) model
components
(using model.glm.control
), 2) prior
distributions (using prior.glm.control
), 3) options for the
MCMC algorithm (using mcmc.control
), and 4) options for the
output (using output.glm.control
).
Default values are available in most of the cases.
The arguments for the control functions are described in their
respective help files.
In the prediction part of the function we want to predict exp(S^*)/(1+exp(S^*)) at locations of interest. For the prediction part of the algorithm, we use the median of the predictive distribution as the predictor and 1/4 of the length of the 95 percent predictive interval as a measure of the prediction uncertainty. Below we describe the procedure for calculating these quantities.
First we perform a Bayesian Gaussian prediction with the given priors on beta
and sigma^2 on each of the simulated S“datasets” from the
posterior distribution (and in case phi is not fixed, for each sampled phi value).
This Gaussian prediction is done by an internal function which
is an extension of krige.bayes
allowing for more than one “data set”.
For calculating the probability below a threshold for the predictive distribution given the data, we first calculate this probability for each of the simulated S“datasets”. This is done using the fact that the predictive distribution for each of the simulated S“datasets” is a multivariate tdistribution. Afterwards the probability below a threshold is calculated by taking the empirical mean of these conditional probabilities.
Now the median and the 2.5 percent and 97.5 percent quantiles can be calculated by an iterative procedure, where first a guess of the value is made, and second this guess is checked by calculating the probability of being less than this value. In case the guess is too low, it is adjusted upwards, and vise verse.
A list with the following components:
posterior 
A list with results for the posterior distribution of the
model parameters and the random effects at the data locations. The components are:

predictive 
A list with results for the predictive distribution at the
prediction locations (if provided). The
components are:

model 
model components used as defined by 
prior 
priors used for the model parameters. 
mcmc.input 
input parameters used for the MCMC algorithm. 
.Random.seed 
system random seed before running the function.
Allows reproduction of results. If
the 
call 
the function call. 
Ole F. Christensen OleF.Christensen@agrsci.dk,
Paulo J. Ribeiro Jr. Paulo.Ribeiro@est.ufpr.br.
Further information about geoRglm can be found at:
http://gbi.agrsci.dk/~ofch/geoRglm.
binom.krige
for prediction with fixed parameters in the
binomiallogit normal model, pois.krige.bayes
for Bayesian prediction in the
Poisson normal model, krige.bayes
for
Bayesian prediction in the Gaussian spatial model.
1 2 3 4 5 6 7 8 9 10 11 12 13  data(b50)
if(!exists(".Random.seed", envir=.GlobalEnv, inherits = FALSE)) set.seed(1234)
## Not run:
mcmc.10 < mcmc.control(S.scale = 0.09, n.iter = 1000, phi.scale = 0.2,
phi.start = 4.5)
prior.10 < prior.glm.control(phi.discrete = seq(0.2,5,0.2))
test.10 < binom.krige.bayes(b50, locations=t(cbind(c(2.5,3.5),c(1,3.5),c(2.5,1.5),c(1,1.5))),
prior = prior.10, mcmc.input = mcmc.10)
image(test.10)
persp(test.10)
## End(Not run)

Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.