SpatioTemporal (version 1.1.2)

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 class 'STmodel':
MCMC(object, x, x.fixed = NULL,
    type = "f", N = 1000, Hessian.prop = NULL,
    Sigma.prop = NULL, info = min(ceiling(N/50), 100), ...)

Arguments

object

STmodel for which to run MCMC.

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 given th

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 re

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.

info

Outputs status information every silent: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: [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).

See Also

Other mcmcSTmodel methods: print.mcmcSTmodel, print.summary.mcmcSTmodel, summary.mcmcSTmodel

Other STmodel methods: createSTmodel, c.STmodel, estimateCV.STmodel, estimate.STmodel, plot.STdata, plot.STmodel, predictCV.STmodel, predict.STmodel, print.STmodel, print.summary.STmodel, simulate.STmodel, summary.STmodel

Examples

##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
##run MCMC
  MCMC.mesa.model <- MCMC(mesa.model, x$par, N = 5000, Hessian.prop = H)
##lets load precomputed results instead
data(MCMC.mesa.model)

####TODO: NEEDS S3-update
##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$l, ylab="", xlab="", type="l",
     main="log-likelihood")
for(i in c(4,9,13,15)){
  plot(MCMC.mesa.model$par[,i], ylab="", xlab="", type="l",
       main=colnames(MCMC.mesa.model$par)[i])
}

##And estimated densities for the log-covariance parameters.
##The red line is the approximate normal distribution given by
##the maximum-likelihood estimates, e.g. ML-estimate and standard 
##deviation from the observed information matrix.
par(mfrow=c(3,3),mar=c(4,4,2.5,.5))
for(i in 9:17){
  xd <- sort(unique(MCMC.mesa.model$par[,i]))
  yd <- dnorm(xd, mean=x$par[i],sd=x$sd[i])
  dens <- density(MCMC.mesa.model$par[,i])
  plot(dens, ylim=c(0,max(c(dens$y,yd))), main =
       colnames(MCMC.mesa.model$par)[i])
  lines(xd,yd,col=2)
}