Learn R Programming

hbmem (version 0.2)

gammaSample: Function gammaSample

Description

Runs MCMC for the hierarchical Gamma model

Usage

gammaSample(dat, M = 10000, keep = (M/10):M, getDIC = TRUE,
freeCrit=TRUE,shape=2,jump=.005)

Arguments

dat
Data frame that must include variables cond,sub,item,lag,resp. Indexes for cond, sub, item, and respone must start at zero and have no gapes (i.e., no skipped subject numbers). Lags must be zero-centered.
M
Number of MCMC iterations.
keep
Which MCMC iterations should be included in estimates and returned. Use keep to both get ride of burn-in, and thin chains if necessary
getDIC
Logical. should the function compute DIC value? This takes a while if M is large.
freeCrit
Logical. If TRUE (default) individual criteria vary across people. If false, all participants have the same criteria (but note that overall response biases are still modeled in the means)
shape
Fixed shape across both new and studied distributuions.
jump
The criteria and decorrelating steps utilize Matropolis-Hastings sampling routines, which require tuning. All MCMC functions should self tune during the burnin perior (iterations before keep), and they will alert you to the success of tuning.

Value

  • The function returns an internally defined "uvsd" S4 class that includes the following components
  • muIndexes which element of blocks contain grand means, mu
  • alphaIndexes which element of blocks contain participant effects, alpha
  • betaIndexes which element of blocks contain item effects, beta
  • s2alphaIndexes which element of blocks contain variance of participant effects (alpha).
  • s2betaIndexes which element of blocks contain variance of item effects (beta).
  • thetaIndexes which element of blocks contain theta, the slope of the lag effect
  • estNPosterior means of block parameters for new-item means
  • estSPosterior means of block parameters for studied-item means
  • estS2Not used for gamma model.
  • estCritPosterior means of criteria
  • blockNEach iteration for each parameter in the new-item mean block. Rows index iteration, columns index parameter.
  • blockSSame as blockN, but for the studied-item means
  • blockS2Not used for gamma model.
  • s.critSamples of each criteria.
  • pDNumber of effective parameters used in DIC. Note that this should be smaller than the actual number of parameters, as constraint from the hierarchical structure decreases the number of effective parameters.
  • DICDIC value. Smaller values indicate better fits. Note that DIC is notably biased toward complexity.
  • MNumber of MCMC iterations run
  • keepMCMC iterations that were used for estimation and returned
  • b0Metropolis-Hastings acceptance rates for new-item distribution parameters. These should be between .2 and .6. If they are not, the M, keep, or jump need to be adjusted.
  • b0S2Metropolis-Hastings acceptance rates for studied-item distribution parameters.
  • b0CritMetropolis-Hastings acceptance rates for criteria.

See Also

hbmem

Examples

Run this code
#make data from gamma model
library(hbmem)
sim=gammaSim(I=30,J=200)
dat=as.data.frame(cbind(sim@subj,sim@item,sim@cond,sim@Scond,sim@lag,sim@resp))
colnames(dat)=c("sub","item","cond","Scond","lag","resp")

M=200 #set very small for demo speed
keep=50:M
gamma=gammaSample(dat,M=M,keep=keep,jump=.01)

par(mfrow=c(3,2),pch=19,pty='s')
#Look at chains of MuN and MuS
matplot(gamma@blockN[,gamma@muN],t='l',xlab="Iteration",ylab="Mu-N")
abline(h=sim@muN,col="blue")
matplot(gamma@blockS[,gamma@muS],t='l',xlab="Iteration",ylab="Mu-S")
abline(h=sim@muS,col="blue")

#Estimates of Alpha as function of true values
plot(gamma@estN[gamma@alphaN]~sim@alphaN,xlab="True
Alpha-N",ylab="Est. Alpha-N");abline(0,1,col="blue")
plot(gamma@estS[gamma@alphaS]~sim@alphaS,xlab="True
Alpha-S",ylab="Est. Alpha-S");abline(0,1,col="blue")
#Estimates of Beta as function of true values
plot(gamma@estN[gamma@betaN]~sim@betaN,xlab="True
Beta-N",ylab="Est. Beta-N");abline(0,1,col="blue")
plot(gamma@estS[gamma@betaS]~sim@betaS,xlab="True
Beta-S",ylab="Est. Beta-S");abline(0,1,col="blue")

gamma@estN[c(gamma@s2alphaN,gamma@s2betaN)]
gamma@estS[c(gamma@s2alphaS,gamma@s2betaS)]


#Look at some criteria
par(mfrow=c(2,2))
for(i in 1:4)
matplot(t(gamma@s.crit[i,,]),t='l')

Run the code above in your browser using DataLab