sirt (version 2.4-20)

mcmc_coef: Some Methods for Objects of Class mcmc.list

Description

Some methods for objects of class mcmc.list created from the coda package.

Usage

## coefficients	
mcmc_coef(mcmcobj, exclude = "deviance")
## covariance matrix 
mcmc_vcov(mcmcobj, exclude = "deviance")
## confidence interval
mcmc_confint( mcmcobj, parm, level = .95, exclude="deviance" )
## summary function
mcmc_summary( mcmcobj , quantiles=c(.025,.05,.50,.95,.975) )
## plot function
mcmc_plot(mcmcobj, ...)
## inclusion of derived parameters in mcmc object
mcmc_derivedPars( mcmcobj , derivedPars )
## Wald test for parameters
mcmc_WaldTest( mcmcobj , hypotheses )
# S3 method for mcmc_WaldTest
summary(object, digits=3, ...)

Arguments

mcmcobj

Objects of class mcmc.list as created by coda::mcmc

exclude

Vector of parameters which should be excluded in calculations

parm

Optional vector of parameters

level

Confidence level

quantiles

Vector of quantiles to be computed.

…

Parameters to be passed to mcmc_plot. See LAM::plot.amh for arguments.

derivedPars

List with derived parameters (see examples).

hypotheses

List with hypotheses of the form \(g_i( \bold{\theta})=0\).

object

Object of class mcmc_WaldTest.

digits

Number of digits used for rounding.

See Also

coda::mcmc

Examples

# NOT RUN {
#############################################################################
# EXAMPLE 1: Logistic regression in rcppbugs package 
#############################################################################


#***************************************
# (1) simulate data
set.seed(8765)
N <- 500
x1 <- stats::rnorm(N)
x2 <- stats::rnorm(N)
y <- 1*( stats::plogis( -.6 + .7*x1 + 1.1 *x2 ) > stats::runif(N) )

#***************************************
# (2) estimate logistic regression with glm
mod <- stats::glm( y ~ x1 + x2 , family="binomial" )
summary(mod)

#***************************************
# (3) estimate model with rcppbugs package
library(rcppbugs)
b <- rcppbugs::mcmc.normal( stats::rnorm(3),mu=0,tau=0.0001)
y.hat <- rcppbugs::deterministic( function(x1,x2,b){
                stats::plogis( b[1] + b[2]*x1 + b[3]*x2 ) }, 
                  x1 , x2 , b)
y.lik <- rcppbugs::mcmc.bernoulli( y , p = y.hat, observed = TRUE)
model <- rcppbugs::create.model(b, y.hat, y.lik)

#*** estimate model in rcppbugs; 5000 iterations, 1000 burnin iterations
n.burnin <- 500 ; n.iter <- 2000 ; thin <- 2
ans <- rcppbugs::run.model(model , iterations=n.iter, burn=n.burnin, adapt=200, thin=thin)
print(rcppbugs::get.ar(ans)) # get acceptance rate
print(apply(ans[["b"]],2,mean)) # get means of posterior

#*** convert rcppbugs into mcmclist object
mcmcobj <- data.frame( ans$b )
colnames(mcmcobj) <- paste0("b",1:3)
mcmcobj <- as.matrix(mcmcobj)
class(mcmcobj) <- "mcmc"
attr(mcmcobj, "mcpar") <- c( n.burnin+1 , n.iter , thin )
mcmcobj <- coda::mcmc( mcmcobj )

# coefficients, variance covariance matrix and confidence interval
mcmc_coef(mcmcobj)
mcmc_vcov(mcmcobj)
mcmc_confint( mcmcobj , level = .90 )

# summary and plot
mcmc_summary(mcmcobj)
mcmc_plot(mcmcobj, ask=TRUE)

# include derived parameters in mcmc object
derivedPars <- list( "diff12" = ~ I(b2-b1) , "diff13" = ~ I(b3-b1) )
mcmcobj2 <- sirt::mcmc_derivedPars(mcmcobj , derivedPars = derivedPars )
mcmc_summary(mcmcobj2)

#*** Wald test for parameters
 # hyp1: b2 - 0.5 = 0
 # hyp2: b2 * b3 = 0
hypotheses <- list( "hyp1" = ~ I( b2 - .5 )  , "hyp2" = ~ I( b2*b3 ) )
test1 <- sirt::mcmc_WaldTest( mcmcobj , hypotheses=hypotheses )
summary(test1)
# }