Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
Type: Package
Package: modelbased
Title: Estimation of Model-Based Predictions, Contrasts and Means
Version: 0.17.0.10
Version: 0.17.0.11
Authors@R:
c(person(given = "Dominique",
family = "Makowski",
Expand Down
9 changes: 9 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,15 @@

## Bug fixes

* With `backend = "marginaleffects"`, `estimate_means()`, `estimate_slopes()`
and `estimate_contrasts()` no longer drop or rename a `by` variable with a
name such as `type`, `n` or `value`. They also no longer mislabel other
columns because of such a name. Before, for example, `estimate_means()`
dropped a `by` variable named `type` and renamed `n` to `n_Obs`. It also
labeled a variable named `Median` as `Mean`, and the means as `Coefficient`.
This fix does not yet apply to Bayesian models or to inequality comparisons
(for example `comparison = "inequality"`).

* `estimate_means()` with `backend = "emmeans"` now shows the response level
names for models from `ordinal::clmm()`. Before, the `Response` column showed
`1`, `2`, `3`, and the `Threshold` column showed `1|2`, `2|3`.
Expand Down
69 changes: 66 additions & 3 deletions R/format.R
Original file line number Diff line number Diff line change
Expand Up @@ -722,6 +722,16 @@ equivalence_columns <- c(
is_contrast_analysis = FALSE,
...
) {
# focal terms can have names that `model_parameters()` renames or that we
# remove below (e.g. "type", "n" or "p"). We give these columns temporary
# names and restore the original names after all renaming is done.
focal_terms <- attr(x, "focal_terms", exact = TRUE)
protected_focal <- .protected_focal_columns(x, focal_terms)
if (length(protected_focal)) {
colnames(x)[match(protected_focal, colnames(x))] <- names(protected_focal)
focal_terms[match(protected_focal, focal_terms)] <- names(protected_focal)
}

# tidy output - we want to tidy the output, using `model_parameters()` or
# `describe_posterior()` for Bayesian models. We also need to know how the
# coefficient column is named, because we replace that column name with an
Expand Down Expand Up @@ -783,9 +793,8 @@ equivalence_columns <- c(
params <- params[c(setdiff(colnames(params), relocate_columns), relocate_columns)]

# relocate focal terms to the beginning
by <- attr(x, "focal_terms", exact = TRUE)
if (!is.null(by) && all(by %in% colnames(params))) {
params <- datawizard::data_reorder(params, by, verbose = FALSE)
if (!is.null(focal_terms) && all(focal_terms %in% colnames(params))) {
params <- datawizard::data_reorder(params, focal_terms, verbose = FALSE)
}

# rename coefficient name and statistics columns
Expand Down Expand Up @@ -829,6 +838,12 @@ equivalence_columns <- c(
params <- .safe(datawizard::data_rename(params, "group", "Class"), params)
}

# restore original names of protected focal terms
if (length(protected_focal)) {
restore <- match(names(protected_focal), colnames(params))
colnames(params)[restore[!is.na(restore)]] <- protected_focal[!is.na(restore)]
Comment thread
jmgirard marked this conversation as resolved.
}

# extract information for potential inequality analysis
by_terms <- attributes(x)$hypothesis_by
comparison_hypothesis <- list(...)$hypothesis
Expand Down Expand Up @@ -919,6 +934,54 @@ equivalence_columns <- c(
}


# Returns the focal terms that are columns of the marginaleffects output `x`
# and hold the values of the data grid, named with temporary column names.
# Columns that marginaleffects writes itself (e.g. "hypothesis" or "s.value")
# do not hold the data grid values, so these are not protected. Bayesian
# models are tidied by `describe_posterior()`, which may drop unknown columns,
# so we do not protect focal terms there.

#' @keywords internal
.protected_focal_columns <- function(x, focal_terms) {
if (is.null(focal_terms) || .is_bayesian_marginaleffects(x)) {
return(character(0))
}
# the grid used for predictions is the "newdata" component, which differs
# from the "datagrid" attribute when users provide their own `newdata`
grids <- list(
.safe(marginaleffects::components(x, "newdata")),
attr(x, "datagrid", exact = TRUE)
)
grids <- grids[!vapply(grids, is.null, logical(1))]
if (!length(grids)) {
return(character(0))
}
protected <- vapply(
focal_terms,
function(i) {
i %in% colnames(x) &&
any(vapply(
grids,
function(g) {
i %in% colnames(g) && all(as.character(x[[i]]) %in% as.character(g[[i]]))
},
logical(1)
))
},
logical(1)
)
protected <- unname(focal_terms[protected])
if (!length(protected)) {
return(character(0))
}
# temporary names must not match any existing column
temporary <- make.unique(
c(colnames(x), paste0(".modelbased_focal_", seq_along(protected)))
)[-seq_along(colnames(x))]
stats::setNames(protected, temporary)
}


# This function ensures that the formatted object still has all relevant
# information saved as attributes

Expand Down
11 changes: 10 additions & 1 deletion tests/testthat/test-brms.R
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,16 @@ test_that("estimate_means - brms", {
estim2 <- estimate_means(model, iterations = 200, backend = "marginaleffects")
expect_identical(attributes(estim2)$iterations, 200)

expect_true(all(abs(estim1$Median - estim2$Median) > 0.001))
# a random subsample of 200 draws gives other medians than all draws. The
# differences can be smaller than 0.001, so we only test that they differ.
expect_false(isTRUE(all.equal(estim1$Median, estim2$Median, tolerance = 0)))
estim3 <- estimate_means(
model,
iterations = 200,
keep_iterations = TRUE,
backend = "marginaleffects"
)
expect_length(grep("^iter_", colnames(estim3)), 200)
})

test_that("estimate_relation - brms", {
Expand Down
228 changes: 228 additions & 0 deletions tests/testthat/test-focal-name-clashes.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,228 @@
skip_on_cran()
skip_if_not_installed("marginaleffects", minimum_version = "0.29.0")

# Focal variables can have names that insight, parameters or modelbased use
# for their own columns. Each test renames the focal variable and compares the
# output with the output for a neutral name (`safe_var`).

# names that `insight::standardize_names()` changes, or that modelbased drops
# by name
focal_reported_names <- c(
"type", "p", "CI", "S", "s.value", "predicted", "conf.level", "df_residual",
"df_error", "response", "n", "value", "component", "effect", "weight",
"weights", "sumsq", "meansq", "std.dev", "y.level", "hypothesis", "Median",
"Mean", "MAP", "Dxy", "rho", "r", "tau", "Difference", "Estimate", "F",
"Chi2", "chisq", "Chisq", "W"
)

focal_data <- function(name, numeric = FALSE) {
d <- datasets::mtcars
d$am_f <- factor(d$am)
d[[name]] <- if (numeric) d$wt else factor(d$cyl)
d
}

# the formula environment is the calling function, so that the data `d` is
# found there and not in the global environment
focal_formula <- function(rhs, env = parent.frame()) {
stats::as.formula(paste("mpg ~", rhs), env = env)
}

focal_calls <- list(
means = function(name) {
d <- focal_data(name)
m <- stats::lm(focal_formula(paste0("`", name, "` + hp")), data = d)
estimate_means(m, name)
},
means_two = function(name) {
d <- focal_data(name)
m <- stats::lm(focal_formula(paste0("`", name, "` * am_f")), data = d)
estimate_means(m, c("am_f", name))
},
slopes_by = function(name) {
d <- focal_data(name)
m <- stats::lm(focal_formula(paste0("`", name, "` * wt")), data = d)
estimate_slopes(m, "wt", by = name)
},
contrasts_by = function(name) {
d <- focal_data(name)
m <- stats::lm(focal_formula(paste0("`", name, "` * am_f")), data = d)
estimate_contrasts(m, "am_f", by = name)
},
contrasts = function(name) {
d <- focal_data(name)
m <- stats::lm(focal_formula(paste0("`", name, "` + hp")), data = d)
estimate_contrasts(m, name)
}
)

# Returns NULL when `out` equals `ref` with "safe_var" replaced by `name`,
# otherwise a short description of the first difference.
focal_difference <- function(ref, out, name) {
if (inherits(out, "error")) {
return(paste("error:", conditionMessage(out)))
}
swap <- function(x) {
if (is.character(x)) gsub("safe_var", name, x, fixed = TRUE) else x
}
ref_names <- swap(colnames(ref))
if (!identical(ref_names, colnames(out))) {
return(paste0(
"columns [", toString(colnames(out)), "], expected [",
toString(ref_names), "]"
))
}
for (i in seq_along(ref)) {
r <- ref[[i]]
o <- out[[i]]
if (!identical(class(r), class(o))) {
return(paste0("column ", i, " has class ", class(o)[1]))
}
if (is.factor(r)) {
levels(r) <- swap(levels(r))
}
r <- swap(r)
if (!isTRUE(all.equal(r, o, check.attributes = FALSE))) {
return(paste0("column ", i, " has other values"))
}
}
for (a in c("focal_terms", "by", "trend", "contrast")) {
if (!identical(swap(attr(ref, a, exact = TRUE)), attr(out, a, exact = TRUE))) {
return(paste0("attribute ", a, " differs"))
}
}
NULL
}

expect_focal_names_match <- function(call_name, names) {
call <- focal_calls[[call_name]]
ref <- call("safe_var")
for (name in names) {
out <- tryCatch(suppressWarnings(call(name)), error = function(e) e)
diff <- focal_difference(ref, out, name)
expect_null(diff, label = paste0(call_name, "(\"", name, "\"): ", diff))
}
}


test_that("estimate_means() keeps focal variables named like renamed columns", {
names <- setdiff(focal_reported_names, c("Mean", "s.value"))
expect_focal_names_match("means", names)
expect_focal_names_match("means_two", names)
})


test_that("estimate_slopes() keeps by-variables named like renamed columns", {
names <- setdiff(focal_reported_names, c("p", "s.value", "predicted"))
expect_focal_names_match("slopes_by", names)
})


test_that("estimate_contrasts() keeps by-variables named like renamed columns", {
names <- setdiff(
focal_reported_names,
c("p", "s.value", "hypothesis", "Difference")
)
expect_focal_names_match("contrasts_by", names)
})


test_that("focal variable names that worked before still work", {
# names that matched on upstream main b964ea75, computed once from a list
# of 99 candidate names (the reported names and other column, argument and
# output names)
matched_before <- list(
means = c(
"z", "rowidcf", "wts", "predicted_lo", "predicted_hi", "emmean",
"prob", "rate", "lower.CL", "upper.CL", "asymp.LCL", "asymp.UCL",
"Level1", "Level2", "Parameter", "Comparison", "Response", "Threshold",
"Row", "Residuals", "Group", "Level", "Slope", "Coefficient",
"Probability", "Contrast", "Component", "Effects", "Type", "Ratio",
"iter_1", "x", "model", "data", "at", "ci", "trend", "predict",
"comparison", "length", "range", "fun", "i", "out", "means", "grid"
),
means_two = c(
"z", "rowidcf", "wts", "predicted_lo", "predicted_hi", "emmean",
"prob", "rate", "lower.CL", "upper.CL", "asymp.LCL", "asymp.UCL",
"Level1", "Level2", "Parameter", "Comparison", "Response", "Threshold",
"Row", "Residuals", "Group", "Level", "Slope", "Coefficient",
"Probability", "Contrast", "Component", "Effects", "Type", "Ratio",
"iter_1", "x", "model", "data", "at", "ci", "trend", "predict",
"comparison", "length", "range", "fun", "i", "out", "means", "grid"
),
slopes_by = c(
"rowidcf", "wts", "emmean", "prob", "rate", "lower.CL", "upper.CL",
"asymp.LCL", "asymp.UCL", "Level1", "Level2", "Response", "Threshold",
"Row", "Residuals", "Group", "Level", "Probability", "Contrast",
"Component", "Effects", "Type", "Ratio", "iter_1", "x", "model",
"data", "at", "ci", "trend", "predict", "comparison", "length",
"range", "fun", "i", "out", "means", "grid"
),
contrasts_by = c(
"z", "rowidcf", "wts", "predicted_lo", "predicted_hi", "emmean",
"prob", "rate", "lower.CL", "upper.CL", "asymp.LCL", "asymp.UCL",
"Level1", "Level2", "Comparison", "Response", "Threshold", "Row",
"Residuals", "Group", "Level", "Probability", "Contrast", "Component",
"Effects", "Type", "Ratio", "iter_1", "x", "model", "data", "at", "ci",
"trend", "predict", "comparison", "length", "range", "fun", "i", "out",
"means", "grid"
),
contrasts = c(
"z", "F", "Chi2", "chisq", "Chisq", "W", "Median", "Mean", "MAP",
"Dxy", "rho", "r", "tau", "df_residual", "df_error", "std.dev",
"s.value", "component", "effect", "predicted", "response", "y.level",
"conf.level", "n", "type", "weight", "weights", "hypothesis", "value",
"sumsq", "meansq", "Estimate", "rowidcf", "wts", "predicted_lo",
"predicted_hi", "emmean", "prob", "rate", "lower.CL", "upper.CL",
"asymp.LCL", "asymp.UCL", "Predicted", "CI", "Level1", "Level2",
"Comparison", "Response", "Threshold", "Row", "Residuals", "Group",
"Level", "Slope", "Coefficient", "Probability", "Contrast",
"Component", "Effects", "Type", "Ratio", "S", "iter_1", "x", "model",
"data", "at", "ci", "trend", "predict", "comparison", "length",
"range", "fun", "i", "out", "means", "grid"
)
)
for (call_name in names(matched_before)) {
expect_focal_names_match(call_name, matched_before[[call_name]])
}
})


test_that("estimate_means() keeps numeric focal variables for glm", {
names <- setdiff(focal_reported_names, c("Mean", "s.value"))
call <- function(name) {
d <- focal_data(name, numeric = TRUE)
d$vs <- datasets::mtcars$vs
m <- stats::glm(
stats::as.formula(paste0("vs ~ am_f * `", name, "`")),
family = stats::binomial(),
data = d
)
estimate_means(m, c("am_f", name))
}
ref <- suppressWarnings(call("safe_var"))
for (name in names) {
out <- tryCatch(suppressWarnings(call(name)), error = function(e) e)
diff <- focal_difference(ref, out, name)
expect_null(diff, label = paste0("glm(\"", name, "\"): ", diff))
}
})


test_that("estimate_means() keeps focal variables with custom newdata", {
# the values in `newdata` are not part of the generated data grid
names <- setdiff(focal_reported_names, c("Mean", "s.value"))
call <- function(name) {
d <- focal_data(name, numeric = TRUE)
m <- stats::lm(focal_formula(paste0("`", name, "` + hp")), data = d)
newdata <- data.frame(x = c(2.123, 4.567), hp = 150)
colnames(newdata)[1] <- name
estimate_means(m, by = name, newdata = newdata)
}
ref <- call("safe_var")
for (name in names) {
out <- tryCatch(suppressWarnings(call(name)), error = function(e) e)
diff <- focal_difference(ref, out, name)
expect_null(diff, label = paste0("newdata(\"", name, "\"): ", diff))
}
})
3 changes: 3 additions & 0 deletions tests/testthat/test-visualisation_recipe.R
Original file line number Diff line number Diff line change
Expand Up @@ -434,6 +434,9 @@ test_that("visualization_recipe, emmeans, recipes that do not change", {
expect_null(out$color)
.expect_same_print(x)
}
})


test_that("visualization_recipe, informative error when nothing goes on x-axis", {
# intercept-only model
model <- lm(qsec ~ 1, data = mtcars)
Expand Down
Loading