# Mixed Models (lme4, glmmTMB, MixMod, ...) -----------------------------
# =======================================================================
#' @rdname get_predicted
#' @export
get_predicted.lmerMod <- function(x,
data = NULL,
predict = "expectation",
ci = NULL,
ci_method = NULL,
include_random = "default",
iterations = NULL,
verbose = TRUE,
...) {
# Sanitize input
my_args <- .get_predicted_args(
x,
data = data,
predict = predict,
ci = ci,
ci_method = ci_method,
include_random = include_random,
verbose = verbose,
...
)
# Make prediction only using random if only random
if (all(names(my_args$data) %in% find_random(x, flatten = TRUE))) {
random.only <- TRUE
} else {
random.only <- FALSE
}
# Prediction function
predict_function <- function(x, ...) {
stats::predict(
x,
newdata = my_args$data,
type = my_args$type,
re.form = my_args$re.form,
random.only = random.only,
allow.new.levels = my_args$allow_new_levels,
...
)
}
# 1. step: predictions
if (is.null(iterations)) {
predictions <- predict_function(x)
} else {
predictions <- .get_predicted_boot(
x,
data = my_args$data,
predict_function = predict_function,
iterations = iterations,
verbose = verbose,
...
)
}
# 2. step: confidence intervals
ci_data <- get_predicted_ci(x, predictions,
data = my_args$data, ci = ci,
ci_method = ci_method, ci_type = my_args$ci_type,
...
)
# 3. step: back-transform
out <- .get_predicted_transform(x, predictions, my_args, ci_data, verbose = verbose)
# 4. step: final preparation
.get_predicted_out(out$predictions, my_args = my_args, ci_data = out$ci_data)
}
#' @export
get_predicted.merMod <- get_predicted.lmerMod
# glmmTMB ---------------------------------------------------------------
# =======================================================================
#' @export
get_predicted.glmmTMB <- function(x,
data = NULL,
predict = "expectation",
ci = NULL,
include_random = "default",
iterations = NULL,
verbose = TRUE,
...) {
# validation checks
if (!is.null(predict) && predict %in% c("prediction", "predicted", "classification")) {
predict <- "expectation"
if (verbose) {
format_warning(
"\"prediction\" and \"classification\" are currently not supported by the `predict` argument for `glmmTMB` models.",
"Changing to `predict=\"expectation\"`."
)
}
}
# TODO: prediction intervals
# https://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#predictions-andor-confidence-or-prediction-intervals-on-predictions
# Sanitize input
my_args <- .get_predicted_args(
x,
data = data,
predict = predict,
ci = ci,
include_random = include_random,
verbose = verbose,
...
)
# Prediction function
predict_function <- function(x, data, ...) {
stats::predict(
x,
newdata = data,
type = my_args$type,
re.form = my_args$re.form,
allow.new.levels = my_args$allow_new_levels,
...
)
}
# 1. step: predictions
rez <- predict_function(x, data = my_args$data, se.fit = TRUE)
if (is.null(iterations)) {
predictions <- as.numeric(rez$fit)
} else {
predictions <- .get_predicted_boot(
x,
data = my_args$data,
predict_function = predict_function,
iterations = iterations,
verbose = verbose,
...
)
}
# "expectation" for zero-inflated? we need a special handling
# for predictions and CIs here.
if (my_args$scale == "response" && my_args$info$is_zero_inflated) {
# intermediate step: prediction from ZI model, for non-truncated families!
# for truncated family, behaviour in glmmTMB changed in 1.1.5 to correct
# conditional and response predictions
if (!my_args$info$is_hurdle) {
zi_predictions <- stats::predict(
x,
newdata = data,
type = "zprob",
re.form = my_args$re.form,
...
)
predictions <- link_inverse(x)(predictions) * (1 - as.vector(zi_predictions))
}
# 2. and 3. step: confidence intervals and back-transform
ci_data <- .simulate_zi_predictions(
model = x,
newdata = data,
predictions = predictions,
nsim = iterations,
ci = ci
)
out <- list(predictions = predictions, ci_data = ci_data)
} else {
# 2. step: confidence intervals
ci_data <- .get_predicted_se_to_ci(
x,
predictions = predictions,
se = rez$se.fit,
ci = ci,
verbose = verbose,
...
)
# 3. step: back-transform
out <- .get_predicted_transform(x, predictions, my_args, ci_data, verbose = verbose)
}
# 4. step: final preparation
.get_predicted_out(out$predictions, my_args = my_args, ci_data = out$ci_data)
}
# GLMMadaptive: mixed_model (class MixMod) ------------------------------
# =======================================================================
#' @export
get_predicted.MixMod <- function(x,
data = NULL,
predict = "expectation",
ci = NULL,
include_random = "default",
iterations = NULL,
verbose = TRUE,
...) {
# validation checks
if (!is.null(predict) && predict %in% c("prediction", "predicted", "classification")) {
predict <- "expectation"
if (verbose) {
format_warning(
"\"prediction\" and \"classification\" are currently not supported by the `predict` argument for `GLMMadaptive` models.",
"Changing to `predict=\"expectation\"`."
)
}
}
# Sanitize input
my_args <- .get_predicted_args(
x,
data = data,
predict = predict,
ci = ci,
include_random = include_random,
verbose = verbose,
...
)
# Prediction function
predict_function <- function(x, data, ...) {
stats::predict(
x,
newdata = data,
type_pred = my_args$type,
type = ifelse(isTRUE(my_args$include_random), "subject_specific", "mean_subject"),
...
)
}
# 1. step: predictions
rez <- predict_function(x, data = my_args$data)
if (is.null(iterations)) {
predictions <- as.numeric(rez)
} else {
predictions <- .get_predicted_boot(
x,
data = my_args$data,
predict_function = predict_function,
iterations = iterations,
verbose = verbose,
...
)
}
# "expectation" for zero-inflated? we need a special handling
# for predictions and CIs here.
if (my_args$scale == "response" && my_args$info$is_zero_inflated) {
# 2. and 3. step: confidence intervals and back-transform
ci_data <- .simulate_zi_predictions(model = x, newdata = data, predictions = predictions, nsim = iterations, ci = ci)
out <- list(predictions = predictions, ci_data = ci_data)
} else {
# 2. step: confidence intervals
ci_data <- get_predicted_ci(
x,
predictions,
data = my_args$data[colnames(my_args$data) != find_response(x)],
ci = ci,
ci_type = my_args$ci_type,
...
)
# 3. step: back-transform
out <- .get_predicted_transform(x, predictions, my_args, ci_data)
}
# 4. step: final preparation
.get_predicted_out(out$predictions, my_args = my_args, ci_data = out$ci_data)
}
# HGLM: mixed_model (class hglm) ------------------------------
# =============================================================
#' @export
get_predicted.hglm <- function(x, verbose = TRUE, ...) {
# hglm only provide fitted values
x$fv
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.