Learn R Programming

DiffusionRgqd (version 0.1.2)

GQD.dic: Summarize MCMC Selection Output for a List of GQD.mcmc or BiGQD.mcmc objects.

Description

GQD.dic() summarizes the MCMC output from a list of GQD.mcmc() objects. This may be used to neatly summarize the MCMC output of various models fitted to a given dataset.

Usage

GQD.dic(model.list, type = "col")

Arguments

model.list
A list of GQD.mcmc() objects.
type
Shoould output be of row ('row') or column ('col') format.

Value

A data frame with summary of model output. See Details.

Details

GQD.dic() summarizes the output from various models fitted via GQD.mcmc(). By ranking them according to DIC. [=] indicates which model has the minimal DIC.
Elapsed_Time Time_Homogeneous p DIC pD
N Model 1 00:00:47 No 5.00 389.30
3.92 201 Model 2 00:01:00 No 5.00
[=]386.45 4.13 201 Model 3 00:02:50 No
5.00 391.37 3.94 201 Elapsed_Time

References

Updates available on GitHub at https://github.com/eta21.

See Also

GQD.mcmc

Examples

Run this code

#===============================================================================
# Simulate a time inhomogeneous diffusion.
#-------------------------------------------------------------------------------

 data(SDEsim1)
 attach(SDEsim1)
 par(mfrow=c(1,1))
 expr1=expression(dX[t]==2*(5+3*sin(0.5*pi*t)-X[t])*dt+0.5*sqrt(X[t])*dW[t])
 plot(Xt~time,type='l',col='blue',xlab='Time (t)',ylab=expression(X[t]),main=expr1)


 #------------------------------------------------------------------------------
 # Define coefficients of model 1
 #------------------------------------------------------------------------------

  # Remove any existing coeffients. If none are pressent NAs will be returned, but
  # this is a safeguard against overlapping.
    GQD.remove()

  # Define time dependant coefficients. Note that all functions have a single argument.
  # This argument has to be `t' in order for the dependancy to be recognized.
  # theta does not have to be defined as an argument.

  G0 <- function(t){theta[1]*(theta[2]+theta[3]*sin(0.25*pi*t))}
  G1 <- function(t){-theta[1]}
  Q0 <- function(t){theta[4]*theta[4]}

  theta.start   <-  c(1,1,1,1)                     # Starting values for the chain
  proposal.sds  <-  c(0.8,0.1,0.1,0.1)             # Std devs for proposal distributions
  mesh.points   <-  10                             # Number of mesh points
  updates       <-  50000                          # Perform 50000 updates

  priors=function(theta){dnorm(theta[1],5,5)}
  m1 <- GQD.mcmc(Xt,time,mesh=mesh.points,theta=theta.start,sds=proposal.sds,
                 updates=updates)

 
 #------------------------------------------------------------------------------
 # Define coefficients of model 2
 #------------------------------------------------------------------------------

  GQD.remove()
  G0 <- function(t){theta[1]*(theta[2]+theta[3]*sin(0.25*pi*t))}
  G1 <- function(t){-theta[1]}
  Q1 <- function(t){theta[4]*theta[4]}

  proposal.sds  <-  c(0.8,0.1,0.1,0.1)
  m2 <- GQD.mcmc(Xt,time,mesh=mesh.points,theta=theta.start,sds=proposal.sds,
                 updates=updates)

  # Checkthe parameter estimates:
  GQD.estimates(m2,thin = 200)
 #------------------------------------------------------------------------------
 # Summarize the output from the models.
 #------------------------------------------------------------------------------

 GQD.dic(list(m1,m2))
 
#===============================================================================

Run the code above in your browser using DataLab