MCMC.STmodel: 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.

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).

Usage

# S3 method for STmodel
MCMC(object, x, x.fixed = NULL, type = "f", N = 1000,
  Hessian.prop = NULL, Sigma.prop = NULL, info = min(ceiling(N/50), 100),
  ...)
MCMC(object, ...)

Arguments

object

STmodel for which to run MCMC.

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 given they are inferred by calling predict.STmodel with only.pars=TRUE.

x.fixed

Vector with parameter to be held fixed; parameters marked as NA will still be estimated.

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" will revert (with a warning) to using the full log-likelihood.

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.

info

Outputs status information every info:th iteration. If info=0 no output.

...

ignored additional arguments.

Value

mcmcSTmodel object with elements:

par

A N - by - (number of parameters) matrix with trajectories of the parameters.

log.like

A vector of length N with the log-likelihood values at each iteration.

acceptance

A vector of length N with the acceptance probability for each iteration.

Sigma.prop, chol.prop

Proposal matrix and it's Choleskey factor.

x.fixed

Any fixed parameters.

Details

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:

1: If Sigma.prop is given then this is used.
2: If Sigma.prop=NULL then we follow Roberts et.al. (1997) and compute c <- 2.38*2.38/dim(Hessian.prop)[1] Sigma.prop <- -c*solve(Hessian.prop).
3: If both Sigma.prop=NULL and Hessian.prop=NULL then the Hessian is computed using loglikeSTHessian and Sigma.prop is computed according to point 2.

The resulting proposal matrix is checked to ensure that it is positive definite before proceeding, all(eigen(Sigma.prop)$value > 1e-10).

Examples

Run this code

# NOT RUN {
##load data
data(mesa.model)
##and results of estimation
data(est.mesa.model)

##strating point
x <- coef(est.mesa.model)
##Hessian, for use as proposal matrix
H <- est.mesa.model$res.best$hessian.all
# }
# NOT RUN {
  ##run MCMC
  MCMC.mesa.model <- MCMC(mesa.model, x$par, N = 2500, Hessian.prop = H)
# }
# NOT RUN {
##lets load precomputed results instead
data(MCMC.mesa.model)

##Examine the results
print(MCMC.mesa.model)

##and contens of result vector
names(MCMC.mesa.model)
 
##Summary
summary(MCMC.mesa.model)

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

Run the code above in your browser using DataLab

Get 50% off unlimited learning