# Find more examples and instructions in the package vignette:
# vignette('hglm')
require(hglm)
# --------------------- #
# semiconductor example #
# --------------------- #
data(semiconductor)
m11 <- hglm(fixed = y ~ x1 + x3 + x5 + x6,
random = ~ 1|Device,
family = Gamma(link = log),
disp = ~ x2 + x3, data = semiconductor)
summary(m11)
plot(m11, cex = .6, pch = 1,
cex.axis = 1/.6, cex.lab = 1/.6,
cex.main = 1/.6, mar = c(3, 4.5, 0, 1.5))
# ------------------- #
# redo it using hglm2 #
# ------------------- #
m12 <- hglm2(y ~ x1 + x3 + x5 + x6 + (1|Device),
family = Gamma(link = log),
disp = ~ x2 + x3, data = semiconductor)
summary(m12)
# -------------------------- #
# redo it using matrix input #
# -------------------------- #
attach(semiconductor)
m13 <- hglm(y = y, X = model.matrix(~ x1 + x3 + x5 + x6),
Z = kronecker(diag(16), rep(1, 4)),
X.disp = model.matrix(~ x2 + x3),
family = Gamma(link = log))
summary(m13)
# --------------------- #
# verbose & likelihoods #
# --------------------- #
m14 <- hglm(fixed = y ~ x1 + x3 + x5 + x6,
random = ~ 1|Device,
family = Gamma(link = log),
disp = ~ x2 + x3, data = semiconductor,
verbose = TRUE, calc.like = TRUE)
summary(m14)
# --------------------------------------------- #
# simulated example with 2 random effects terms #
# --------------------------------------------- #
## Not run:
# set.seed(911)
# x1 <- rnorm(100)
# x2 <- rnorm(100)
# x3 <- rnorm(100)
# z1 <- factor(rep(LETTERS[1:10], rep(10, 10)))
# z2 <- factor(rep(letters[1:5], rep(20, 5)))
# Z1 <- model.matrix(~ 0 + z1)
# Z2 <- model.matrix(~ 0 + z2)
# u1 <- rnorm(10, 0, sqrt(2))
# u2 <- rnorm(5, 0, sqrt(3))
# y <- 1 + 2*x1 + 3*x2 + Z1%*%u1 + Z2%*%u2 + rnorm(100, 0, sqrt(exp(x3)))
# dd <- data.frame(x1 = x1, x2 = x2, x3 = x3, z1 = z1, z2 = z2, y = y)
#
# (m21 <- hglm(X = cbind(rep(1, 100), x1, x2), y = y, Z = cbind(Z1, Z2),
# RandC = c(10, 5)))
# summary(m21)
# plot(m21)
#
# (m22 <- hglm2(y ~ x1 + x2 + (1|z1) + (1|z2), data = dd, vcovmat = TRUE))
# image(m22$vcov, main = 'Variance-covariance Matrix')
# summary(m22)
# plot(m22)
#
# m31 <- hglm2(y ~ x1 + x2 + (1|z1) + (1|z2), disp = ~ x3, data = dd)
# print (m31)
# summary(m31)
# plot(m31)
# ## End(Not run)
Run the code above in your browser using DataLab