## Not run: ------------------------------------
# ## Simulated fusion Lasso example.
# bmu <- c(0,0,0,2,2,2,4,4,4)
# bsigma <- c(0,0,0,-2,-2,-2,-1,-1,-1)
# id <- factor(sort(rep(1:length(bmu), length.out = 300)))
#
# ## Response.
# set.seed(123)
# y <- bmu[id] + rnorm(length(id), sd = exp(bsigma[id]))
#
# ## Estimate model:
# ## fuse=1 -> nominal fusion,
# ## fuse=2 -> ordinal fusion,
# ## first, unpenalized model to be used for adaptive fusion weights.
# f <- list(y ~ la(id,fuse=2,fx=TRUE), sigma ~ la(id,fuse=1,fx=TRUE))
# b0 <- bamlss(f, sampler = FALSE)
#
# ## Model with single lambda parameter.
# f <- list(y ~ la(id,fuse=2), sigma ~ la(id,fuse=1))
# b1 <- bamlss(f, sampler = FALSE, optimizer = lasso,
# criterion = "BIC", zeromodel = b0)
#
# ## Plot information criterion and coefficient paths.
# lasso.plot(b1, which = 1)
# lasso.plot(b1, which = 2)
# lasso.plot(b1, which = 2, name = "mu.s.la(id).id")
# lasso.plot(b1, which = 2, name = "sigma.s.la(id).id")
#
# ## Extract coefficients for optimum Lasso parameter.
# coef(b1, mstop = lasso.stop(b1))
#
# ## Predict with optimum Lasso parameter.
# p1 <- predict(b1, mstop = lasso.stop(b1))
#
# ## Full MCMC, needs lasso.transform() to assign the
# ## adaptive weights from unpenalized model b0.
# b2 <- bamlss(f, optimizer = FALSE, transform = lasso.transform,
# zeromodel = b0, nobs = length(y), start = coef(b1, mstop = lasso.stop(b1)))
# summary(b2)
# plot(b2)
#
# ci <- confint(b2, model = "mu", pterms = FALSE, sterms = TRUE)
# lasso.plot(b1, which = 2, name = "mu.s.la(id).id", spar = FALSE)
# for(i in 1:8) {
# abline(h = ci[i, 1], lty = 2, col = "red")
# abline(h = ci[i, 2], lty = 2, col = "red")
# }
## ---------------------------------------------
Run the code above in your browser using DataLab