context("mhglm")
library("lme4")
test_that("succeeds on sleepstudy", {
model <- mhglm(Reaction ~ Days + (Days | Subject), data=sleepstudy)
# fixef
fixef0 <- c("(Intercept)" = 251.4, "Days" = 10.5)
expect_that(round(fixef(model), 1), equals(fixef0))
# vcov
vcov0 <- matrix(c(44.0, -1.4, -1.4, 2.3), 2, 2)
rownames(vcov0) <- colnames(vcov0) <- c("(Intercept)", "Days")
expect_that(round(vcov(model), 1), equals(vcov0))
# VarCorr
varcor0 <- matrix(c(565.5, 11.1, 11.1, 32.7), 2, 2)
rownames(varcor0) <- colnames(varcor0) <- c("(Intercept)", "Days")
varcor <- VarCorr(model)[["Subject"]]
expect_that(attr(varcor, "stddev"), equals(sqrt(diag(varcor))))
expect_that(as.vector(attr(varcor, "correlation")),
equals(as.vector(t(varcor / attr(varcor, "stddev"))
/ attr(varcor, "stddev"))))
attr(varcor, "stddev") <- NULL
attr(varcor, "correlation") <- NULL
expect_that(round(varcor, 1), equals(varcor0))
# ranef
ranef0 <- matrix(c( 2.8, -40.0, -38.4, 22.8, 21.5, 8.8, 16.4, -7.0, -1.0,
34.7, -24.6, -12.3, 4.3, 20.6, 3.3, -24.7, 0.7, 12.1,
9.1, -8.6, -5.5, -4.7, -2.9, -0.2, -0.2, 1.0, -10.6,
8.6, 1.1, 6.5, -3.0, 3.6, 0.9, 4.7, -1.0, 1.3),
18, 2)
rownames(ranef0) <- as.character(c(308, 309, 310, 330, 331, 332, 333, 334,
335, 337, 349, 350, 351, 352, 369, 370,
371, 372))
colnames(ranef0) <- c("(Intercept)", "Days")
expect_that(round(as.matrix(ranef(model)[["Subject"]]), 1), equals(ranef0))
})
test_that("succeeds on sleepstudy in parallel", {
model <- mhglm(Reaction ~ Days + (Days | Subject), data=sleepstudy,
control=list(parallel=TRUE))
# fixef
fixef0 <- c("(Intercept)" = 251.4, "Days" = 10.5)
expect_that(round(fixef(model), 1), equals(fixef0))
# vcov
vcov0 <- matrix(c(44.0, -1.4, -1.4, 2.3), 2, 2)
rownames(vcov0) <- colnames(vcov0) <- c("(Intercept)", "Days")
expect_that(round(vcov(model), 1), equals(vcov0))
# VarCorr
varcor0 <- matrix(c(565.5, 11.1, 11.1, 32.7), 2, 2)
rownames(varcor0) <- colnames(varcor0) <- c("(Intercept)", "Days")
varcor <- VarCorr(model)[["Subject"]]
expect_that(attr(varcor, "stddev"), equals(sqrt(diag(varcor))))
expect_that(as.vector(attr(varcor, "correlation")),
equals(as.vector(t(varcor / attr(varcor, "stddev"))
/ attr(varcor, "stddev"))))
attr(varcor, "stddev") <- NULL
attr(varcor, "correlation") <- NULL
expect_that(round(varcor, 1), equals(varcor0))
# ranef
ranef0 <- matrix(c( 2.8, -40.0, -38.4, 22.8, 21.5, 8.8, 16.4, -7.0, -1.0,
34.7, -24.6, -12.3, 4.3, 20.6, 3.3, -24.7, 0.7, 12.1,
9.1, -8.6, -5.5, -4.7, -2.9, -0.2, -0.2, 1.0, -10.6,
8.6, 1.1, 6.5, -3.0, 3.6, 0.9, 4.7, -1.0, 1.3),
18, 2)
rownames(ranef0) <- as.character(c(308, 309, 310, 330, 331, 332, 333, 334,
335, 337, 349, 350, 351, 352, 369, 370,
371, 372))
colnames(ranef0) <- c("(Intercept)", "Days")
expect_that(round(as.matrix(ranef(model)[["Subject"]]), 1), equals(ranef0))
})
test_that("succeeds on cbpp", {
suppressWarnings({
model <- mhglm(cbind(incidence, size - incidence) ~ period + (period | herd),
data=cbpp, family=binomial)
})
# fixef
fixef0 <- c("(Intercept)" = -1.2, "period2" = -0.8, "period3" = -0.9,
"period4" = -1.1)
expect_that(round(fixef(model), 1), equals(fixef0))
# vcov
vcov0 <- matrix(c( 0.1, -0.1, -0.1, -0.1, -0.1, 0.2, 0.1, 0.0,
-0.1, 0.1, 0.2, 0.0, -0.1, 0.0, 0.0, 0.1), 4, 4)
rownames(vcov0) <- colnames(vcov0) <- c("(Intercept)", "period2",
"period3", "period4")
expect_that(round(vcov(model), 1), equals(vcov0))
# VarCorr
varcor0 <- matrix(c( 0.7, -0.7, -0.4, -0.5, -0.7, 1.0, 1.0, 0.3,
-0.4, 1.0, 1.4, 0.0, -0.5, 0.3, 0.0, 0.3), 4, 4)
rownames(varcor0) <- colnames(varcor0) <- c("(Intercept)", "period2",
"period3", "period4")
varcor <- VarCorr(model)[["herd"]]
expect_that(attr(varcor, "stddev"), equals(sqrt(diag(varcor))))
expect_that(as.vector(attr(varcor, "correlation")),
equals(as.vector(t(varcor / attr(varcor, "stddev"))
/ attr(varcor, "stddev"))))
attr(varcor, "stddev") <- NULL
attr(varcor, "correlation") <- NULL
expect_that(round(varcor, 1), equals(varcor0))
# ranef
ranef0 <- matrix(c(-0.1, -0.4, 0.5, 0.1, 0.3, -0.2, 0.9, 0.5, 0.0, -0.7,
-0.7, 0.0, -0.8, 1.2, -0.7,
0.8, 0.2, -0.6, 0.2, -0.8, -0.1, -0.6, -0.5, -0.2, 0.3,
1.0, 0.1, 0.7, -1.1, 0.7,
1.5, -0.1, -0.6, 0.5, -1.2, -0.5, 0.2, -0.3, -0.3, -0.2,
1.1, 0.1, 0.2, -0.7, 0.5,
-0.3, 0.3, -0.3, -0.2, 0.1, 0.3, -0.8, -0.3, 0.1, 0.6,
0.2, 0.0, 0.6, -0.7, 0.4),
15, 4)
rownames(ranef0) <- as.character(1:15)
colnames(ranef0) <- c("(Intercept)", "period2", "period3", "period4")
expect_that(round(as.matrix(ranef(model)[["herd"]]), 1), equals(ranef0))
})
test_that("succeeds on cbpp in parallel", {
suppressWarnings({
model <- mhglm(cbind(incidence, size - incidence) ~ period + (period | herd),
data=cbpp, family=binomial,
control=list(parallel=TRUE))
})
# fixef
fixef0 <- c("(Intercept)" = -1.2, "period2" = -0.8, "period3" = -0.9,
"period4" = -1.1)
expect_that(round(fixef(model), 1), equals(fixef0))
# vcov
vcov0 <- matrix(c( 0.1, -0.1, -0.1, -0.1, -0.1, 0.2, 0.1, 0.0,
-0.1, 0.1, 0.2, 0.0, -0.1, 0.0, 0.0, 0.1), 4, 4)
rownames(vcov0) <- colnames(vcov0) <- c("(Intercept)", "period2",
"period3", "period4")
expect_that(round(vcov(model), 1), equals(vcov0))
# VarCorr
varcor0 <- matrix(c( 0.7, -0.7, -0.4, -0.5, -0.7, 1.0, 1.0, 0.3,
-0.4, 1.0, 1.4, 0.0, -0.5, 0.3, 0.0, 0.3), 4, 4)
rownames(varcor0) <- colnames(varcor0) <- c("(Intercept)", "period2",
"period3", "period4")
varcor <- VarCorr(model)[["herd"]]
expect_that(attr(varcor, "stddev"), equals(sqrt(diag(varcor))))
expect_that(as.vector(attr(varcor, "correlation")),
equals(as.vector(t(varcor / attr(varcor, "stddev"))
/ attr(varcor, "stddev"))))
attr(varcor, "stddev") <- NULL
attr(varcor, "correlation") <- NULL
expect_that(round(varcor, 1), equals(varcor0))
# ranef
ranef0 <- matrix(c(-0.1, -0.4, 0.5, 0.1, 0.3, -0.2, 0.9, 0.5, 0.0, -0.7,
-0.7, 0.0, -0.8, 1.2, -0.7,
0.8, 0.2, -0.6, 0.2, -0.8, -0.1, -0.6, -0.5, -0.2, 0.3,
1.0, 0.1, 0.7, -1.1, 0.7,
1.5, -0.1, -0.6, 0.5, -1.2, -0.5, 0.2, -0.3, -0.3, -0.2,
1.1, 0.1, 0.2, -0.7, 0.5,
-0.3, 0.3, -0.3, -0.2, 0.1, 0.3, -0.8, -0.3, 0.1, 0.6,
0.2, 0.0, 0.6, -0.7, 0.4),
15, 4)
rownames(ranef0) <- as.character(1:15)
colnames(ranef0) <- c("(Intercept)", "period2", "period3", "period4")
expect_that(round(as.matrix(ranef(model)[["herd"]]), 1), equals(ranef0))
})
test_that("succeeds on sleepstudy with uncorrelated covariates assumption", {
model <- mhglm(Reaction ~ Days + (Days || Subject), data=sleepstudy)
# fixef
fixef0 <- c("(Intercept)" = 251.4, "Days" = 10.5)
expect_that(round(fixef(model), 1), equals(fixef0))
# vcov
vcov0 <- matrix(c(47.0, -2, -2.0, 2.4), 2, 2)
rownames(vcov0) <- colnames(vcov0) <- c("(Intercept)", "Days")
expect_that(round(vcov(model), 1), equals(vcov0))
# VarCorr
varcor0 <- matrix(c(583.5, 0.0, 0.0, 33.6), 2, 2)
rownames(varcor0) <- colnames(varcor0) <- c("(Intercept)", "Days")
varcor <- VarCorr(model)[["Subject"]]
expect_that(attr(varcor, "stddev"), equals(sqrt(diag(varcor))))
expect_that(as.vector(attr(varcor, "correlation")),
equals(as.vector(t(varcor / attr(varcor, "stddev"))
/ attr(varcor, "stddev"))))
attr(varcor, "stddev") <- NULL
attr(varcor, "correlation") <- NULL
expect_that(round(varcor, 1), equals(varcor0))
# ranef
ranef0 <- matrix(c( 1.9, -40.0, -38.7, 23.9, 22.4, 9.0, 16.8, -7.3, 0.1,
34.5, -25.2, -13.3, 4.7, 20.7, 3.2, -25.8, 0.8, 12.2,
9.2, -8.6, -5.4, -4.9, -3.1, -0.3, -0.2, 1.1, -10.8,
8.6, 1.2, 6.6, -3.0, 3.5, 0.9, 4.9, -1.0, 1.3),
18, 2)
rownames(ranef0) <- as.character(c(308, 309, 310, 330, 331, 332, 333, 334,
335, 337, 349, 350, 351, 352, 369, 370,
371, 372))
colnames(ranef0) <- c("(Intercept)", "Days")
expect_that(round(as.matrix(ranef(model)[["Subject"]]), 1), equals(ranef0))
})
test_that("succeeds on sleepstudy in parallel with uncorrelated covariates assumption", {
model <- mhglm(Reaction ~ Days + (Days || Subject), data=sleepstudy,
control=list(parallel=TRUE))
# fixef
fixef0 <- c("(Intercept)" = 251.4, "Days" = 10.5)
expect_that(round(fixef(model), 1), equals(fixef0))
# vcov
vcov0 <- matrix(c(47.0, -2, -2.0, 2.4), 2, 2)
rownames(vcov0) <- colnames(vcov0) <- c("(Intercept)", "Days")
expect_that(round(vcov(model), 1), equals(vcov0))
# VarCorr
varcor0 <- matrix(c(583.5, 0.0, 0.0, 33.6), 2, 2)
rownames(varcor0) <- colnames(varcor0) <- c("(Intercept)", "Days")
varcor <- VarCorr(model)[["Subject"]]
expect_that(attr(varcor, "stddev"), equals(sqrt(diag(varcor))))
expect_that(as.vector(attr(varcor, "correlation")),
equals(as.vector(t(varcor / attr(varcor, "stddev"))
/ attr(varcor, "stddev"))))
attr(varcor, "stddev") <- NULL
attr(varcor, "correlation") <- NULL
expect_that(round(varcor, 1), equals(varcor0))
# ranef
ranef0 <- matrix(c( 1.9, -40.0, -38.7, 23.9, 22.4, 9.0, 16.8, -7.3, 0.1,
34.5, -25.2, -13.3, 4.7, 20.7, 3.2, -25.8, 0.8, 12.2,
9.2, -8.6, -5.4, -4.9, -3.1, -0.3, -0.2, 1.1, -10.8,
8.6, 1.2, 6.6, -3.0, 3.5, 0.9, 4.9, -1.0, 1.3),
18, 2)
rownames(ranef0) <- as.character(c(308, 309, 310, 330, 331, 332, 333, 334,
335, 337, 349, 350, 351, 352, 369, 370,
371, 372))
colnames(ranef0) <- c("(Intercept)", "Days")
expect_that(round(as.matrix(ranef(model)[["Subject"]]), 1), equals(ranef0))
})
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.