Error in estimte_slopes() for emmeans-backend when trend is a factor (#619)

* Error in `estimte_slopes()` for emmeans-backend when `trend` is a factor
Fixes #305

* address comments

* minor formatting

* add tests

* Fixes #618

* conditional on emmeans

* docs

Author: strengejacke

DESCRIPTION

@@ -1,7 +1,7 @@
 Type: Package
 Package: modelbased
 Title: Estimation of Model-Based Predictions, Contrasts and Means
-Version: 0.15.0
+Version: 0.15.0.1
 Authors@R:
     c(person(given = "Dominique",
              family = "Makowski",

NEWS.md

@@ -1,3 +1,10 @@
+# modelbased (devel)
+
+## Changes
+
+* Informative error message when the `trend` variable in `estimate_slopes()` is
+  not numeric and `backend = "emmeans"`.
+
 # modelbased 0.15.0
 
 ## Breaking Changes

R/estimate_contrasts.R

@@ -218,7 +218,7 @@
 #'   confidence bands: Theory, implementation, and an application to SVARs.
 #'   Journal of Applied Econometrics, 34(1), 1–17. \doi{10.1002/jae.2656}
 #'
-#' @examplesIf all(insight::check_if_installed(c("lme4", "marginaleffects", "parameters", "datawizard", "rstanarm"), quietly = TRUE))
+#' @examplesIf all(insight::check_if_installed(c("lme4", "emmeans", "marginaleffects", "parameters", "datawizard", "rstanarm"), quietly = TRUE))
 #' \dontrun{
 #' # Basic usage --------------------------------
 #' # --------------------------------------------
@@ -259,8 +259,15 @@
 #' # "time" only has integer values and few values, so it's treated like a factor
 #' estimate_contrasts(model, "time", by = "education")
 #'
-#' # we set `integer_as_continuous = TRUE` to treat integer as continuous
-#' estimate_contrasts(model, "time", by = "education", integer_as_continuous = 1)
+#' # Setting `integer_as_continuous = TRUE` treats "time" as a continuous
+#' # variable. This allows us to compare its average slope rather than making
+#' # discrete pairwise comparisons at each time point.
+#' estimate_contrasts(
+#'   model,
+#'   contrast = "time",
+#'   by = "education",
+#'   integer_as_continuous = TRUE
+#' )
 #'
 #' # pairwise comparisons for multiple groups
 #' estimate_contrasts(
@@ -290,6 +297,17 @@
 #' model <- lm(happiness ~ puppy_love * dose, data = puppy_love)
 #' estimate_slopes(model, "puppy_love", by = "dose", comparison = cond_tx)
 #'
+#' # Note: for the emmeans-backend, we need to use `estimate_contrasts()` for
+#' # the above example:
+#' cond_tx <- list(`no treatment` = c(1, 0, 0), treatment = c(0, 0.5, 0.5))
+#' estimate_contrasts(
+#'   model,
+#'   contrast = "puppy_love",
+#'   by = "dose",
+#'   comparison = cond_tx,
+#'   backend = "emmeans"
+#' )
+#'
 #' # Other models (mixed, Bayesian, ...) --------
 #' # --------------------------------------------
 #' data <- iris

R/estimate_means.R

@@ -117,7 +117,7 @@
 #' @param ... Other arguments passed, for instance, to [insight::get_datagrid()],
 #' to functions from the **emmeans** or **marginaleffects** package, or to process
 #' Bayesian models via [bayestestR::describe_posterior()]. Examples:
-#' - `insight::get_datagrid()`: Argument such as `length`, `digits` or `range`
+#' - `insight::get_datagrid()`: Arguments such as `length`, `digits` or `range`
 #'   can be used to control the (number of) representative values. For integer
 #'   variables, `protect_integers` modulates whether these should also be
 #'   treated as numerics, i.e. values can have fractions or not.
@@ -134,8 +134,8 @@
 #' - **emmeans**: Internally used functions are `emmeans()` and `emtrends()`.
 #'   Additional arguments can be passed to these functions.
 #' - Bayesian models: For Bayesian models, parameters are cleaned using
-#'   `describe_posterior()`, thus, arguments like, for example, `centrality`,
-#'   `rope_range`, or `test` are passed to that function.
+#'   `bayestestR::describe_posterior()`, thus, arguments like, for example,
+#'   `centrality`, `rope_range`, or `test` are passed to that function.
 #' - Especially for `estimate_contrasts()` with integer focal predictors, for
 #'   which contrasts should be calculated, use argument `integer_as_continuous`
 #'   to set the maximum number of unique values in an integer predictor to treat

R/get_emtrends.R

@@ -13,13 +13,15 @@
 #' get_emtrends(model)
 #' get_emtrends(model, by = "Sepal.Width")
 #' @export
-get_emtrends <- function(model,
-                         trend = NULL,
-                         by = NULL,
-                         predict = NULL,
-                         keep_iterations = FALSE,
-                         verbose = TRUE,
-                         ...) {
+get_emtrends <- function(
+  model,
+  trend = NULL,
+  by = NULL,
+  predict = NULL,
+  keep_iterations = FALSE,
+  verbose = TRUE,
+  ...
+) {
   # check if available
   insight::check_if_installed("emmeans")
 
@@ -38,7 +40,11 @@ get_emtrends <- function(model,
   ))
 
   # handle distributional parameters
-  if (!is.null(predict) && inherits(model, "brmsfit") && predict %in% .brms_aux_elements(model)) {
+  if (
+    !is.null(predict) &&
+      inherits(model, "brmsfit") &&
+      predict %in% .brms_aux_elements(model)
+  ) {
     fun_args$dpar <- predict
   } else {
     fun_args$type <- predict
@@ -70,11 +76,13 @@ get_emtrends <- function(model,
 # =========================================================================
 
 #' @keywords internal
-.guess_emtrends_arguments <- function(model,
-                                      trend = NULL,
-                                      by = NULL,
-                                      verbose = TRUE,
-                                      ...) {
+.guess_emtrends_arguments <- function(
+  model,
+  trend = NULL,
+  by = NULL,
+  verbose = TRUE,
+  ...
+) {
   # Gather info
   model_data <- insight::get_data(model, verbose = FALSE)
   predictors <- intersect(
@@ -86,19 +94,37 @@ get_emtrends <- function(model,
   if (is.null(trend)) {
     trend <- predictors[sapply(model_data[predictors], is.numeric)][1]
     if (!length(trend) || is.na(trend)) {
-      insight::format_error("Model contains no numeric predictor. Please specify `trend`.")
+      insight::format_error(
+        "Model contains no numeric predictor. Please specify `trend`."
+      )
     }
     if (verbose) {
-      insight::format_alert(paste0("No numeric variable was specified for slope estimation. Selecting `trend = \"", trend, "\"`.")) # nolint
+      insight::format_alert(paste0(
+        "No numeric variable was specified for slope estimation. Selecting `trend = \"",
+        trend,
+        "\"`."
+      ))
     }
   }
   if (length(trend) > 1) {
     trend <- trend[1]
     if (verbose) {
-      insight::format_alert(paste0("More than one numeric variable was selected for slope estimation. Keeping only ", trend[1], ".")) # nolint
+      insight::format_alert(paste0(
+        "More than one numeric variable was selected for slope estimation. Keeping only `",
+        trend[1],
+        "`."
+      ))
     }
   }
 
+  if (trend %in% colnames(model_data) && !is.numeric(model_data[[trend]])) {
+    insight::format_error(paste0(
+      "Variable `",
+      trend,
+      "` is not numeric. Slopes can only be estimated for numeric variables when `backend = \"emmeans\"`. Please use `estimate_contrasts()` instead."
+    ))
+  }
+
   my_args <- list(trend = trend, by = by)
   .process_emmeans_arguments(model, args = my_args, data = model_data, ...)
 }

man/estimate_contrasts.Rd

@@ -131,6 +131,30 @@ test_that("estimate_slopes", {
 })
 
 
+test_that("estimate_slopes, errors for emmeans when trend is non-numeric", {
+  skip_if_not_installed("emmeans")
+  data(iris)
+  model <- lm(Sepal.Length ~ Species * Sepal.Width, data = iris)
+  expect_error(
+    estimate_slopes(model, trend = "Species", by = "Sepal.Width", backend = "emmeans"),
+    regex = "Variable `Species` is not numeric",
+    fixed = TRUE
+  )
+  model <- lm(Sepal.Length ~ Species, data = iris)
+  expect_error(
+    estimate_slopes(model, backend = "emmeans"),
+    regex = "Model contains no numeric predictor",
+    fixed = TRUE
+  )
+  model <- lm(Sepal.Length ~ Species * Sepal.Width, data = iris)
+  expect_message(
+    estimate_slopes(model, backend = "emmeans"),
+    regex = "No numeric variable was specified",
+    fixed = TRUE
+  )
+})
+
+
 test_that("estimate_slopes, johnson-neyman p-adjust", {
   data(iris)
   model <- lm(Sepal.Width ~ Petal.Width * Petal.Length, data = iris)