Nothing
library("JointAI")
library("splines")
# Sys.setenv(IS_CHECK = "true")
skip_on_cran()
if (identical(Sys.getenv("NOT_CRAN"), "true")) {
set_seed(1234)
longDF <- JointAI::longDF
# gamma variables
longDF$L1 <- rgamma(nrow(longDF), 2, 4)
longDF$L1mis <- rgamma(nrow(longDF), 2, 4)
longDF$L1mis[sample.int(nrow(longDF), 20)] <- NA
# beta variables
longDF$Be1 <- plogis(longDF$time - longDF$C1)
longDF$Be2 <- plogis(longDF$y + longDF$c1)
longDF$Be2[c(1:20) * 5] <- NA
run_glmm_models <- function() {
sink(tempfile())
on.exit(sink())
invisible(force(suppressWarnings({
models = list(
# no covariates
m0a1 = lme_imp(y ~ 1 + (1 | id), data = longDF, n.adapt = 5, n.iter = 10,
seed = 2020, warn = FALSE, mess = FALSE),
m0a2 = glme_imp(y ~ 1 + (1 | id), family = gaussian(link = "identity"),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0a3 = glme_imp(y ~ 1 + (1 | id), family = gaussian(link = "log"),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0a4 = glme_imp(y ~ 1 + (1 | id), family = gaussian(link = "inverse"),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0b1 = glme_imp(b1 ~ 1 + (1 | id), family = binomial(link = "logit"),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0b2 = glme_imp(b1 ~ 1 + (1 | id), family = binomial(link = "probit"),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0b3 = glme_imp(b1 ~ 1 + (1 | id), family = binomial(link = "log"),
data = longDF, n.adapt = 50, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0b4 = glme_imp(b1 ~ 1 + (1 | id), family = binomial(link = "cloglog"),
data = longDF, n.adapt = 50, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0c1 = glme_imp(L1 ~ 1 + (1 | id), family = Gamma(link = "inverse"),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0c2 = glme_imp(L1 ~ 1 + (1 | id), family = Gamma(link = "log"),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0d1 = glme_imp(p1 ~ 1 + (1 | id), family = poisson(link = "log"),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0d2 = glme_imp(p1 ~ 1 + (1 | id), family = poisson(link = "identity"),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0e1 = lognormmm_imp(L1 ~ 1 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m0f1 = betamm_imp(Be1 ~ 1 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
# only complete
m1a = lme_imp(y ~ C1 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m1b = glme_imp(b1 ~ C1 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020, family = binomial(),
warn = FALSE, mess = FALSE),
m1c = glme_imp(L1 ~ C1 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020, family = Gamma(),
warn = FALSE, mess = FALSE),
m1d = glme_imp(p1 ~ C1 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020, family = poisson(),
warn = FALSE, mess = FALSE),
m1e = lognormmm_imp(L1 ~ C1 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m1f = betamm_imp(Be1 ~ C1 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
# only incomplete
m2a = lme_imp(y ~ c2 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020, warn = FALSE, mess = FALSE),
m2b = glme_imp(b2 ~ c2 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020, family = binomial(),
warn = FALSE, mess = FALSE),
m2c = glme_imp(L1mis ~ c2 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020, family = Gamma(),
warn = FALSE, mess = FALSE),
m2d = glme_imp(p2 ~ c2 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020, family = poisson(),
warn = FALSE, mess = FALSE),
m2e = lognormmm_imp(L1mis ~ c2 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m2f = betamm_imp(Be2 ~ c2 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
# no intercept
m3a = lme_imp(y ~ 0 + C2 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m3b = glme_imp(b2 ~ 0 + C2 + (1 | id), data = longDF,
family = binomial(), n.adapt = 5, n.iter = 10,
seed = 2020, warn = FALSE, mess = FALSE),
m3c = glme_imp(L1mis ~ 0 + C2 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020, family = Gamma(),
warn = FALSE, mess = FALSE),
m3d = glme_imp(p2 ~ 0 + C2 + (1 | id), data = longDF, n.adapt = 5,
n.iter = 10, seed = 2020, family = poisson(),
warn = FALSE, mess = FALSE),
m3e = lognormmm_imp(L1mis ~ 0 + C2 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m3f = betamm_imp(Be2 ~ 0 + C2 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
# as covariate
m4a = lme_imp(c1 ~ c2 + B2 + p2 + L1mis + Be2 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
models = c(p2 = "glmm_poisson_log",
L1mis = "glmm_gamma_inverse",
Be2 = "glmm_beta"),
warn = FALSE, mess = FALSE),
m4b = lme_imp(c1 ~ c2 + b2 + p2 + L1mis + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
models = c(c2 = "glmm_gaussian_inverse",
p2 = "glmm_poisson_identity",
b2 = "glmm_binomial_probit",
L1mis = "glmm_lognorm"),
warn = FALSE, mess = FALSE),
m4c = lme_imp(c1 ~ c2 + b2 + p2 + L1mis + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020, no_model = "time",
models = c(c2 = "glmm_gaussian_log",
p2 = "glmm_poisson_identity",
L1mis = "glmm_gamma_log",
b2 = "glmm_binomial_log"),
warn = FALSE, mess = FALSE),
m4d = lme_imp(c1 ~ c2 + b2 + p2 + L1mis + Be2 + (1 | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
trunc = list(Be2 = c(0, 1)),
shrinkage = "ridge",
models = c(c2 = "glmm_gaussian_log",
p2 = "glmm_poisson_identity",
L1mis = "glmm_gamma_log",
b2 = "glmm_binomial_log"),
warn = FALSE, mess = FALSE),
# complex structures
m5a = lme_imp(y ~ M2 + o2 * abs(C1 - c2) + log(C1) + time + I(time^2) +
(time | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020,
warn = FALSE, mess = FALSE),
m5b = glme_imp(b1 ~ L1mis + abs(c1 - C2) + log(Be2) + time +
(time + I(time^2) | id),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
models = c(C2 = "glm_gaussian_log",
L1mis = "glmm_gamma_inverse",
Be2 = "glmm_beta"), shrinkage = "ridge",
family = binomial(),
warn = FALSE, mess = FALSE),
# no random intercept
m6a = lme_imp(y ~ b2 + C1 + C2 + time + (0 + time | id), data = longDF,
n.adapt = 5, n.iter = 10, seed = 2020, no_model = "time",
warn = FALSE, mess = FALSE),
m6b = glme_imp(b1 ~ c1 + C2 + B1 + time + (0 + time + I(time^2) | id),
data = longDF, n.adapt = 5, n.iter = 10, seed = 2020,
shrinkage = "ridge",
family = binomial(),
warn = FALSE, mess = FALSE),
# spline random effects
m7a = lme_imp(y ~ ns(time, df = 2), random = ~ ns(time, df = 2)|id,
data = longDF, n.iter = 10, adapt = 5, seed = 2020),
m7b = lme_imp(y ~ bs(time, df = 3), random = ~ bs(time, df = 3)|id,
data = longDF, n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m7c = lme_imp(y ~ C1 + c1 + ns(time, df = 3),
random = ~ ns(time, df = 3)|id, data = longDF,
n.iter = 10, nadapt = 5, seed = 2020),
m7d = lme_imp(y ~ C1 + C2 + c1 + ns(time, df = 3), random = ~ time|id,
data = longDF, n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m7e = lme_imp(y ~ C1 + C2 + c1 + ns(time, df = 3),
random = ~ ns(time, df = 3)|id, data = longDF,
no_model = "time", n.iter = 10, n.adapt = 5, seed = 2020),
m7f = lme_imp(y ~ C1 + C2 + c1 + ns(time, df = 3), random = ~ time|id,
data = longDF,n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
# complex random effects
m8a = lme_imp(y ~ c1 + c2 + time, random = ~ time + c2|id,
no_model = "time", data = longDF,
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8b = lme_imp(y ~ c1 + c2 + time, random = ~ time + c2|id,
data = longDF,
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8c = lme_imp(y ~ B2 * c1 + c2 + time, random = ~ time + c1|id,
data = longDF, no_model = "time",
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8d = lme_imp(y ~ B2 * c1 + c2 + time, random = ~ time + c1|id,
data = longDF, n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8e = lme_imp(y ~ C1 + B2 * c1 + c2 + time, random = ~ time + c2|id,
data = longDF, n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8f = lme_imp(y ~ C1 + B2 * c1 + c2 + time, random = ~ time + c2|id,
data = longDF, no_model = "time",
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8g = lme_imp(y ~ C1 + B2 * c1 + c2 + time, random = ~ time + c2|id,
data = longDF, no_model = c("time", "c1"),
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8h = lme_imp(y ~ C1 + B2 * c2 + c1 + time, random = ~ time + c1|id,
data = longDF, n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8i = lme_imp(y ~ C1 + B2 * c2 + c1 + time, no_model = "time",
random = ~ time + c1|id, data = longDF,
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8j = lme_imp(y ~ C1 + B2 * c2 + c1 + time, random = ~ time + c2 | id,
data = longDF, n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8k = lme_imp(y ~ C1 + B2 * c2 + c1 + time, random = ~ time + c2 | id,
data = longDF, n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
# 3-fold interaction on different levels with random slope variable
m8l = lme_imp(y ~ C1 + B2 * c1 * time, random = ~ time + I(time^2) | id,
data = longDF, n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8m = lme_imp(y ~ c1 * b1 + o1, random = ~b1|id, data = longDF,
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
m8n = lme_imp(y ~ c1 + C1 * time + b1 + B2,
random = ~C1 * time|id, data = longDF,
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
# only lvl-1 variables, crossed random effects
m9a = lme_imp(y ~ c1 + b1 + time + (1|id) + (1 | o1), data = longDF,
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE),
# bugfix: random slope and monitor analysis_random (issue with
# re-scaling)
m9b = lme_imp(y ~ C1 + C2 + B1 + time + (time | id), data = longDF,
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE,
monitor_params = c(analysis_random = TRUE)),
# bugfix: random intercept only (gave warning about failing to trace
# monitor for RinvD)
m9c = lme_imp(y ~ C1 + C2 + B1 + (1 | id), data = longDF,
n.iter = 10, n.adapt = 5, seed = 2020,
warn = FALSE, mess = FALSE,
monitor_params = c(analysis_random = TRUE))
)
}
)))
}
models <- run_glmm_models()
models0 <- set0_list(models)
test_that("models run", {
for (k in seq_along(models)) {
expect_s3_class(models[[k]], "JointAI")
}
})
test_that("models have the correct model and analysis model type", {
for (i in seq_along(models)) {
expect_false(any(sapply(compare_modeltype(models), isFALSE)))
}
})
test_that("there are no duplicate betas/alphas in the jagsmodel", {
expect_null(unlist(lapply(models, find_dupl_parms)))
})
test_that("MCMC is mcmc.list", {
for (i in seq_along(models)) {
expect_s3_class(models[[i]]$MCMC, "mcmc.list")
}
})
test_that("MCMC samples can be plottet", {
for (k in seq_along(models)) {
expect_silent(traceplot(models[[k]]))
expect_silent(densplot(models[[k]]))
expect_silent(plot(MC_error(models[[k]])))
}
})
test_that("data_list remains the same", {
skip_on_os("linux")
expect_snapshot(lapply(models, "[[", "data_list"))
})
test_that("jagsmodel remains the same", {
expect_snapshot(lapply(models, "[[", "jagsmodel"))
})
test_that("GRcrit and MCerror give same result", {
expect_snapshot(lapply(models0, GR_crit, multivariate = FALSE))
expect_snapshot(lapply(models0, MC_error))
})
test_that("summary output remained the same", {
expect_snapshot(lapply(models0, print))
expect_snapshot(lapply(models0, coef))
expect_snapshot(lapply(models0, confint))
expect_snapshot(lapply(models0, summary, missinfo = TRUE))
expect_snapshot(lapply(models0, function(x) coef(summary(x))))
})
test_that("prediction works", {
for (k in seq_along(models)) {
expect_s3_class(predict(models[[k]], type = "link", warn = FALSE)$fitted,
"data.frame")
expect_s3_class(predict(models[[k]], type = "response",
warn = FALSE)$fitted,
"data.frame")
}
# prediction without specifying the type
expect_s3_class(predict(models$m5a, warn = FALSE)$fitted, "data.frame")
# prediction with newdata
ndf <- predDF(models$m5a, vars = ~ c2)
expect_s3_class(ndf, "data.frame")
expect_s3_class(predict(models$m5a, newdata = ndf, warn = FALSE)$fitted,
"data.frame")
})
test_that("residuals", {
for (k in seq_along(models)[8]) {
expect_type(residuals(models[[k]], type = "response", warn = FALSE),
"double")
if (models[[k]]$analysis_type %in% c("beta", "glmm_beta")) {
expect_error(residuals(models[[k]], type = "working", warn = FALSE))
expect_error(residuals(models[[k]], type = "pearson", warn = FALSE))
} else {
expect_type(residuals(models[[k]], type = "working", warn = FALSE),
"double")
expect_type(residuals(models[[k]], type = "pearson", warn = FALSE),
"double")
}
}
expect_type(residuals(models$m5a, warn = FALSE), "double")
})
test_that("model can be plottet", {
for (i in seq_along(models)) {
if (models[[i]]$analysis_type %in% c("lognorm", "glmm_lognorm",
"beta", "glmm_beta")) {
expect_error(plot(models[[i]]))
} else {
expect_silent(plot(models[[i]]))
}
}
})
}
# test_that("wrong models give errors", {
#
# # cauchit link is not implemented
# expect_error(glm_imp(B1 ~ 1, family = binomial(link = "cauchit"),
# data = wideDF, seed = 2020))
#
# # sqrt is not an allowed link (not implemented)
# expect_error(glm_imp(P1 ~ 1, family = poisson(link = "sqrt"), data = wideDF,
# seed = 2020))
#
# # gives JAGS model error (no dedicated error message)
# # glm_imp(time ~ 1, family = Gamma(link = "identity"), data = wideDF,
# # n.adapt = 5, n.iter = 10, seed = 2020)
#
# # no family specified
# expect_error(glm_imp(B1 ~ C1, data = wideDF, seed = 2020))
#
# # unknown covariate model type
# expect_error(lm_imp(C1 ~ C2 + Be2, data = wideDF,
# n.adapt = 5, n.iter = 10, seed = 2020,
# models = c(Be2 = "betareg")))
# })
# Sys.setenv(IS_CHECK = "")
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.