# NOT RUN {
# create artificial data
set.seed(21)
n <- 100
dat <- data.frame(
x=runif(n),
f=factor(sample(1:20, n, replace=TRUE))
)
model <- ~ reg(~ x, Q0=1, name="beta") + gen(factor=~f, name="v")
gd <- generate_data(model, data=dat)
dat$y <- gd$y
# fit a base model
model0 <- y ~ reg(~ 1, name="beta") + gen(factor=~f, name="v")
sampler <- create_sampler(model0, data=dat, block=TRUE)
sim <- MCMCsim(sampler, store.all=TRUE)
(summ0 <- summary(sim))
# fit 'true' model
model <- y ~ reg(~ x, name="beta") + gen(factor=~f, name="v")
sampler <- create_sampler(model, data=dat, block=TRUE)
sim <- MCMCsim(sampler, store.all=TRUE)
(summ <- summary(sim))
# compare random effect estimates against true parameter values
plot_coef(summ0$v, summ$v, list(est=gd$pars$v), n.se=2, offset=0.2,
maxrows=10, est.names=c("base model", "true model", "true"))
# }
# NOT RUN {
# }
Run the code above in your browser using DataLab