Learn R Programming

SpatioTemporal (version 0.9.2)

run.MCMC: MCMC Inference of Parameters in the Spatio-Temporal Model.

Description

Estimates parameters and parameter uncertainties for the spatio-temporal model using a Metropolis-Hastings based Markov Chain Monte Carlo (MCMC) algorithm.

Usage

run.MCMC(x, mesa.data.model, type = "f", N = 1000, 
         Hessian.prop = NA, Sigma.prop = NA, silent = TRUE)

Arguments

x
Point at which to start the MCMC. Could be either only log-covariance parameters or regression and log-covariance parameters. If regression parameters are given but not needed they are dropped, if they are needed but not give
mesa.data.model
Data structure holding observations, and information regarding which geographic and spatio-temporal covariates to use when fitting the model. See create.data.model and
type
A single character indicating the type of log-likelihood to compute. Valid options are "f" or "r", for full, or restricted maximum likelihood (REML).

Since profile is not a proper likelihood type="p"

N
Number of MCMC iterations to run.
Hessian.prop
Hessian (information) matrix for the log-likelihood, can be used to create a proposal matrix for the MCMC.
Sigma.prop
Proposal matrix for the MCMC.
silent
If FALSE outputs status information and brief progress information every 100:th iteration. If TRUE no output.

Value

  • Returns a list containing:
  • parA N - by - (number of parameters) matrix with the MCMC paths of the estimated parameters.
  • log.likeA vector of length N with the log-likelihood values at each iteration.
  • acceptanceA vector of length N with the acceptance probability of the MCMC proposal at each iteration. Can be used to study the overall acceptance rate.
  • chol.propCholeskey factor of the proposal matrix.

encoding

latin1

Details

The function runs uses a Metropolis-Hastings algorithm (Hastings, 1970) to sample from the parameters of the spatio-temporal model, assuming flat priors for all the parameters (flat on the log-scale for the covariance parameters).

At each iteration of the MCMC new parameters are proposed using a random-walk with a proposal covariance matrix. The proposal matrix is determined as: [object Object],[object Object],[object Object] The resulting proposal matrix is checked to ensure that it is positive definite before proceeding, all(eigen(Sigma.prop)$value > 1e-10).

References

W.K. Hastings. (1970) Monte Carlo sampling methods using Markov chains and their applications, Biometrika, 57, 97--109.

G.O. Roberts, A. Gelman, W.R. Gilks. (1997) Weak convergence and optimal scaling of random walk Metropolis algorithms, Annals of Probability, 7, 110--120.

See Also

See fit.mesa.model for ML-parameter estimation, both functions use log-likelihood given by loglike.

Expected names for x are given by loglike.var.names. For further optimization functions see loglike, loglike.var.names, create.data.model, and cond.expectation.

Examples

Run this code
##load data
data(mesa.data.model)
data(mesa.data.res)

##Extract parameters,
par <- mesa.data.res$par.est$res.best$par.all
##... and the Hessian.
H <- mesa.data.res$par.est$res.best$hessian.all
##run the MCMC, this may take a while...
MCMC.res <- run.MCMC(par, mesa.data.model, N = 5000, 
                     Hessian.prop = H, silent = FALSE)
##Get the precomputed results instead.
MCMC.res <- mesa.data.res$MCMC.res

##components of the MCMC results
names(MCMC.res)

##The acceptance probability (alpha) for each step 
##in the Metropolis-Hastings algorithm.
summary(MCMC.res$acceptance)

##The MCMC-estimated parameters
summary(MCMC.res$par)

##MCMC tracks for four of the parameters
par(mfrow=c(4,1),mar=c(2,2,2.5,.5))
for(i in c(4,9,13,15)){
  plot(MCMC.res$par[,i], ylab="", xlab="", type="l",
       main=colnames(MCMC.res$par)[i])
}

Run the code above in your browser using DataLab