Package 'marginaleffects'

Description

This function enables or disables automatic differentiation using the JAX package in Python, which can considerably speed up and increase the accuracy of standard errors when a model includes many parameters.

Usage

autodiff(autodiff = NULL, install = FALSE)

Arguments

Details

Automatic differentiation needs to be enabled once per session.

When autodiff = TRUE, this function:

• Imports the marginaleffects.autodiff Python module via reticulate::import()

• Sets the internal jacobian function to use JAX-based automatic differentiation

• Provides faster and more accurate gradient computation for supported models

• Falls back on the default finite difference method for unsupported models and calls.

Currently supports:

• Model types: lm, glm, ols

• Functions: predictions() and comparisons(), along with avg_ and plot_ variants.

• type: "response" or "link"

• by: TRUE, FALSE, or character vector.

• comparison: "difference" and "ratio"

For unsupported models or options, the function automatically falls back to the default finite difference method.

Value

When autodiff is NULL, returns TRUE if autodiff is enabled and FALSE otherwise. Otherwise called for side effects of enabling or disabling automatic differentiation or installing the Python package.

Python Configuration

By default, no manual configuration of Python should be necessary. On most machines, unless you have explicitly configured reticulate, reticulate defaults to an automatically managed ephemeral virtual environment with all Python requirements declared via reticulate::py_require().

If you prefer to use a manually managed Python installation, you can direct reticulate and specify which Python executable or environment to use. reticulate selects a Python installation using its Order of Discovery. As a convenience autodiff(install=TRUE) will install the marginaleffects Python package in a self-managed virtual environment.

To specify an alternate Python version:

library(reticulate)
use_python("/usr/local/bin/python")

To use a virtual environment:

use_virtualenv("myenv")

These configuration commands should be called before calling autodiff().

Examples

## Not run: 
# Install the Python package (only needed once)
autodiff(install = TRUE)

# Enable automatic differentiation
autodiff(TRUE)

# Fit a model and compute marginal effects
mod <- glm(am ~ hp + wt, data = mtcars, family = binomial)
avg_comparisons(mod) # Will use JAX for faster computation

# Disable automatic differentiation
autodiff(FALSE)

## End(Not run)

Description

Predict the outcome variable at different regressor values (e.g., college graduates vs. others), and compare those predictions by computing a difference, ratio, or some other function. comparisons() can return many quantities of interest, such as contrasts, differences, risk ratios, changes in log odds, lift, slopes, elasticities, etc.

• comparisons(): unit-level (conditional) estimates.

• avg_comparisons(): average (marginal) estimates.

variables identifies the focal regressors whose "effect" we are interested in. comparison determines how predictions with different regressor values are compared (difference, ratio, odds, etc.). The newdata argument and the datagrid() function control where statistics are evaluated in the predictor space: "at observed values", "at the mean", "at representative values", etc.

See the comparisons chapter on the package website for worked examples and case studies:

• https://marginaleffects.com/chapters/comparisons.html

• https://marginaleffects.com/

Usage

comparisons(
  model,
  newdata = NULL,
  variables = NULL,
  comparison = "difference",
  type = NULL,
  vcov = TRUE,
  by = FALSE,
  conf_level = 0.95,
  transform = NULL,
  cross = FALSE,
  wts = FALSE,
  hypothesis = NULL,
  equivalence = NULL,
  df = Inf,
  eps = NULL,
  numderiv = "fdforward",
  ...
)

avg_comparisons(
  model,
  newdata = NULL,
  variables = NULL,
  type = NULL,
  vcov = TRUE,
  by = TRUE,
  conf_level = 0.95,
  comparison = "difference",
  transform = NULL,
  cross = FALSE,
  wts = FALSE,
  hypothesis = NULL,
  equivalence = NULL,
  df = Inf,
  eps = NULL,
  numderiv = "fdforward",
  ...
)

Arguments

Value

A data.frame with one row per estimate. This data frame is pretty-printed by default, but users can interact with it as a regular data frame, with functions like nrow(), head(), colnames(), etc. Values can be extracted using standard ⁠[,]⁠ or $ operators, and manipulated using external packages like dplyr or data.table.

Columns may include:

• rowid: row number of the newdata data frame

• group: (optional) value of the grouped outcome (e.g., categorical outcome models)

• term: the focal variable.

• estimate: an estimate of the prediction, counterfactual comparison, or slope.

• std.error: standard errors computed via the delta method.

• p.value: p value associated to the estimate column. The null is determined by the hypothesis argument (0 by default).

• s.value: Shannon information transforms of p values. See the S values vignette at https://marginaleffects.com the marginaleffects website.

• conf.low: lower bound of the confidence (or credible) interval defined by the conf_level argument.

• conf.high: upper bound of the confidence (or credible) interval defined by the conf_level argument.

• predicted_lo: predicted outcome for the "low" value of the focal predictor in a counterfactual comparison.

• predicted_hi: predicted outcome for the "high" value of the focal predictor in a counterfactual comparison.

• p.rope.unconditional: share of posterior draws in the interval specified by the equivalence argument. This is only available for Bayesian models.

• p.rope.conditional: share of posterior draws in the interval specified by the equivalence argument, among draws in the confidence interval. This is only available for Bayesian models.

• rope: share of the posterior draws between conf.low and conf.high that are covered by the interval specified by the equivalence argument.

• statistic.noninf: test statistic for non-inferiority test (when equivalence argument is used).

• statistic.nonsup: test statistic for non-superiority test (when equivalence argument is used).

• p.value.noninf: p-value for non-inferiority test (when equivalence argument is used).

• p.value.nonsup: p-value for non-superiority test (when equivalence argument is used).

• p.value.equiv: p-value for equivalence test using Two One-Sided Tests (TOST) approach (when equivalence argument is used).

See ?print.marginaleffects for printing options.

The data.frames produced by marginaleffects stores an attribute that holds many internal objects, such as the original model, data, and much other information that can be used for post-processing. This information can be inspected using the components() function.

Warning: The internal attributes retrieved by components() are not considered part of the public API of the package. Their names and contents can change without warning or notice. Users should not rely on them.

Warning: In some cases, the internal attributes used by marginaleffects() can use up a substantial amount of memory. To clear this data, use the prune() function or set options(marginaleffects_lean=TRUE).

Functions

• avg_comparisons(): Average comparisons

Standard errors using the delta method

Standard errors for all quantities estimated by marginaleffects can be obtained via the delta method. This requires differentiating a function with respect to the coefficients in the model using a finite difference approach. In some models, the delta method standard errors can be sensitive to various aspects of the numeric differentiation strategy, including the step size. By default, the step size is set to 1e-8, or to 1e-4 times the smallest absolute model coefficient, whichever is largest.

marginaleffects can delegate numeric differentiation to the numDeriv package, which allows more flexibility. To do this, users can pass arguments to the numDeriv::jacobian function through a global option. For example:

• options(marginaleffects_numDeriv = list(method = "simple", method.args = list(eps = 1e-6)))

• options(marginaleffects_numDeriv = list(method = "Richardson", method.args = list(eps = 1e-5)))

• options(marginaleffects_numDeriv = NULL)

See the "Uncertainty" chapter on the marginaleffects website for more details on the computation of standard errors, bootstrapping, and more:

https://marginaleffects.com/chapters/uncertainty.html

Model-Specific Arguments

Some model types allow model-specific arguments to modify the nature of marginal effects, predictions, marginal means, and contrasts. Please report other package-specific predict() arguments on Github so we can add them to the table below.

https://github.com/vincentarelbundock/marginaleffects/issues

Comparison functions

Each of the quantities computed by comparisons() can be defined as a function of these quantities:

• hi: vector of predictions for the "high" side of the contrast.

• lo: vector of predictions for the "low" side of the contrast.

• y: predictions for the original data.

• x: focal predictor in the original data.

• w: weights

For example, the "lift" of a binary predictor is a popular quantity of interest, defined as the difference between predictions when the focal predictor $X = 1$, and predictions when the focal predictor is $X = 0$, normalized by the starting point. Or:

$\frac{\hat{Y}_{X=1} - \hat{Y}_{X=0}}{\hat{Y}_{X=0}}$

When, the argument is set to comparison="lift", marginaleffects will compute the quantity using this function:

function(hi, lo) { (hi - lo) / lo }

Users can supply custom functions to the comparison argument, or use one of the many shortcuts available for common quantities of interest:

Note on automatic avg switching: when the by argument is used (including by = TRUE or a character vector), marginaleffects will automatically switch comparison shortcuts to their avg variants unless they already end with avg. For example, comparison = "difference" becomes comparison = "differenceavg" when aggregating over by.

Bayesian posterior summaries

By default, credible intervals in bayesian models are built as equal-tailed intervals. This can be changed to a highest density interval by setting a global option:

options("marginaleffects_posterior_interval" = "eti")

options("marginaleffects_posterior_interval" = "hdi")

By default, the center of the posterior distribution in bayesian models is identified by the median. Users can use a different summary function by setting a global option:

options("marginaleffects_posterior_center" = "mean")

options("marginaleffects_posterior_center" = "median")

When estimates are averaged using the by argument, the tidy() function, or the summary() function, the posterior distribution is marginalized twice over. First, we take the average across units but within each iteration of the MCMC chain, according to what the user requested in by argument or tidy()/summary() functions. Then, we identify the center of the resulting posterior using the function supplied to the "marginaleffects_posterior_center" option (the median by default).

Equivalence, Inferiority, Superiority

$\theta$ is an estimate, $\sigma_\theta$ its estimated standard error, and $[a, b]$ are the bounds of the interval supplied to the equivalence argument.

Non-inferiority:

• $H_0$: $\theta \leq a$

• $H_1$: $\theta > a$

• $t=(\theta - a)/\sigma_\theta$

• p: Upper-tail probability

Non-superiority:

• $H_0$: $\theta \geq b$

• $H_1$: $\theta < b$

• $t=(\theta - b)/\sigma_\theta$

• p: Lower-tail probability

Equivalence: Two One-Sided Tests (TOST)

• p: Maximum of the non-inferiority and non-superiority p values.

Thanks to Russell V. Lenth for the excellent emmeans package and documentation which inspired this feature.

Order of operations

Behind the scenes, the arguments of marginaleffects functions are evaluated in this order:

1. newdata

2. variables

3. comparison and slope

4. by

5. vcov

6. hypothesis

7. transform

Parallel computation

The slopes() and comparisons() functions can use parallelism to speed up computation. Operations are parallelized for the computation of standard errors, at the model coefficient level. There is always considerable overhead when using parallel computation, mainly involved in passing the whole dataset to the different processes. Thus, parallel computation is most likely to be useful when the model includes many parameters and the dataset is relatively small.

Warning: In many cases, parallel processing will not be useful at all.

To activate parallel computation, users must load the future.apply package, call plan() function, and set a global option.

options(marginaleffects_parallel = TRUE): parallelize delta method computation of standard errors. options(marginaleffects_parallel_inferences = TRUE): parallelize "rsample" or "fwb" bootstrap computation in inferences(). options(marginaleffects_parallel_packages = TRUE): vector of strings with the names of modeling packages used to fit the model, ex: c("survival", "splines")

For example:

library(future.apply)
plan("multisession", workers = 4)
options(marginaleffects_parallel = FALSE)
options(marginaleffects_parallel_inferences = TRUE)
options(marginaleffects_parallel_packages = c("survival", "splines"))

slopes(model)

To disable parallelism in marginaleffects altogether, you can set a global option:

options(marginaleffects_parallel = FALSE)

Global options

The behavior of marginaleffects functions can be modified by setting global options.

Disable some safety checks and warnings:

• options(marginaleffects_startup_message = FALSE)

  • Disable the startup message printed on library(marginaleffects).

• options(marginaleffects_safe = FALSE)

  • Disable safety checks and warnings.

• options(marginaleffects_print_omit = c("p.value", "s.value"))

  • Omit some columns from the printed output.

Enforce lean return objects, sans information about the original model and data, and other ancillary attributes. Note that this will disable some advanced post-processing features and functions like hypotheses.

options(marginaleffects_lean = TRUE)

Other options:

• marginaleffects_plot_gray: logical. If TRUE, the default color of the plot is gray. Default is FALSE.

Types

The type argument determines the scale of the predictions used to compute quantities of interest with functions from the marginaleffects package. Admissible values for type depend on the model object. When users specify an incorrect value for type, marginaleffects will raise an informative error with a list of valid type values for the specific model object. The first entry in the list in that error message is the default type.

The invlink(link) is a special type defined by marginaleffects. It is available for some (but not all) models, and only for the predictions() function. With this link type, we first compute predictions on the link scale, then we use the inverse link function to backtransform the predictions to the response scale. This is useful for models with non-linear link functions as it can ensure that confidence intervals stay within desirable bounds, ex: 0 to 1 for a logit model. Note that an average of estimates with type="invlink(link)" will not always be equivalent to the average of estimates with type="response". This type is default when calling predictions(). It is available—but not default—when calling avg_predictions() or predictions() with the by argument.

Some of the most common type values are:

References

• Arel-Bundock V, Greifer N, Heiss A (2024). “How to Interpret Statistical Models Using marginaleffects for R and Python.” Journal of Statistical Software, 111(9), 1-32. doi:10.18637/jss.v111.i09 doi:10.18637/jss.v111.i09

• Arel-Bundock (2026). "Model to Meaning: How to interpret statistical models in R and Python." CRC Press. https://routledge.com/9781032908724

• Greenland S. 2019. "Valid P-Values Behave Exactly as They Should: Some Misleading Criticisms of P-Values and Their Resolution With S-Values." The American Statistician. 73(S1): 106–114.

• Cole, Stephen R, Jessie K Edwards, and Sander Greenland. 2020. "Surprise!" American Journal of Epidemiology 190 (2): 191–93. doi:10.1093/aje/kwaa136

Examples

# Linear model
tmp <- mtcars
tmp$am <- as.logical(tmp$am)
mod <- lm(mpg ~ am + factor(cyl), tmp)
avg_comparisons(mod, variables = list(cyl = "reference"))
avg_comparisons(mod, variables = list(cyl = "sequential"))
avg_comparisons(mod, variables = list(cyl = "pairwise"))

# GLM with different scale types
mod <- glm(am ~ factor(gear), data = mtcars)
avg_comparisons(mod, type = "response")
avg_comparisons(mod, type = "link")

# Contrasts at the mean
comparisons(mod, newdata = "mean")

# Contrasts between marginal means
comparisons(mod, newdata = "balanced")

# Contrasts at user-specified values
comparisons(mod, newdata = datagrid(am = 0, gear = tmp$gear))
comparisons(mod, newdata = datagrid(am = unique, gear = max))

m <- lm(mpg ~ hp + drat + factor(cyl) + factor(am), data = mtcars)
comparisons(m, variables = "hp", newdata = datagrid(FUN_factor = unique, FUN_numeric = median))

# Numeric contrasts
mod <- lm(mpg ~ hp, data = mtcars)
avg_comparisons(mod, variables = list(hp = 1))
avg_comparisons(mod, variables = list(hp = 5))
avg_comparisons(mod, variables = list(hp = c(90, 100)))
avg_comparisons(mod, variables = list(hp = "iqr"))
avg_comparisons(mod, variables = list(hp = "sd"))
avg_comparisons(mod, variables = list(hp = "minmax"))

# using a function to specify a custom difference in one regressor
dat <- mtcars
dat$new_hp <- 49 * (dat$hp - min(dat$hp)) / (max(dat$hp) - min(dat$hp)) + 1
modlog <- lm(mpg ~ log(new_hp) + factor(cyl), data = dat)
fdiff <- function(x) data.frame(x, x + 10)
avg_comparisons(modlog, variables = list(new_hp = fdiff))

# Adjusted Risk Ratio
mod <- glm(vs ~ mpg, data = mtcars, family = binomial)
avg_comparisons(mod, comparison = "lnratioavg", transform = exp)

# Adjusted Risk Ratio: Manual specification of the `comparison`
avg_comparisons(
  mod,
  comparison = function(hi, lo) log(mean(hi) / mean(lo)),
  transform = exp)
# cross contrasts
mod <- lm(mpg ~ factor(cyl) * factor(gear) + hp, data = mtcars)
avg_comparisons(mod, variables = c("cyl", "gear"), cross = TRUE)

# variable-specific contrasts
avg_comparisons(mod, variables = list(gear = "sequential", hp = 10))

# hypothesis test: is the `hp` marginal effect at the mean equal to the `drat` marginal effect
mod <- lm(mpg ~ wt + drat, data = mtcars)

comparisons(
  mod,
  newdata = "mean",
  hypothesis = "wt = drat")

# same hypothesis test using row indices
comparisons(
  mod,
  newdata = "mean",
  hypothesis = "b1 - b2 = 0")

# same hypothesis test using numeric vector of weights
comparisons(
  mod,
  newdata = "mean",
  hypothesis = c(1, -1))

# two custom contrasts using a matrix of weights
lc <- matrix(
  c(
    1, -1,
    2, 3),
  ncol = 2)
comparisons(
  mod,
  newdata = "mean",
  hypothesis = lc)

# Effect of a 1 group-wise standard deviation change
# First we calculate the SD in each group of `cyl`
# Second, we use that SD as the treatment size in the `variables` argument
library(dplyr)
mod <- lm(mpg ~ hp + factor(cyl), mtcars)
tmp <- mtcars %>%
  group_by(cyl) %>%
  mutate(hp_sd = sd(hp))
avg_comparisons(mod,
  variables = list(hp = function(x) data.frame(x, x + tmp$hp_sd)),
  by = "cyl")

# `by` argument
mod <- lm(mpg ~ hp * am * vs, data = mtcars)
comparisons(mod, by = TRUE)

mod <- lm(mpg ~ hp * am * vs, data = mtcars)
avg_comparisons(mod, variables = "hp", by = c("vs", "am"))

library(nnet)
mod <- multinom(factor(gear) ~ mpg + am * vs, data = mtcars, trace = FALSE)
by <- data.frame(
  group = c("3", "4", "5"),
  by = c("3,4", "3,4", "5"))
comparisons(mod, type = "probs", by = by)

Description

Extract components from marginaleffects objects

Usage

## S3 method for class 'marginaleffects'
components(object, component = NULL, ...)

Arguments

Details

This function provides access to the internal components stored in the mfx attribute of marginaleffects objects. The mfx attribute contains an S4 object of class "marginaleffects_internal" with various slots containing model information, data, and computational details used by the marginaleffects functions.

Warning: the internal slot names are not considered part of the public API and may change without warning in future versions of the marginaleffects package.

Value

The requested component from the mfx object

Description

Generate a data grid of user-specified values for use in the newdata argument of the predictions(), comparisons(), and slopes() functions. This is useful to define where in the predictor space we want to evaluate the quantities of interest. Ex: the predicted outcome or slope for a 37 year old college graduate.

Usage

datagrid(
  ...,
  model = NULL,
  newdata = NULL,
  by = NULL,
  grid_type = "mean_or_mode",
  response = FALSE,
  FUN = NULL,
  FUN_character = NULL,
  FUN_factor = NULL,
  FUN_logical = NULL,
  FUN_numeric = NULL,
  FUN_integer = NULL,
  FUN_binary = NULL,
  FUN_other = NULL
)

Arguments

Details

If datagrid is used in a predictions(), comparisons(), or slopes() call as the newdata argument, the model is automatically inserted in the model argument of datagrid() call, and users do not need to specify either the model or newdata arguments. The same behavior will occur when the value supplied to ⁠newdata=⁠ is a function call which starts with "datagrid". This is intended to allow users to create convenience shortcuts like:

Warning about hierarchical grouping variables: When using the default grid_type = "mean_or_mode" with hierarchical models (such as mixed models with nested grouping factors), datagrid() may create invalid combinations of grouping variables. For example, if you have students nested within schools, or countries nested within regions, the modal values of each grouping variable may not correspond to valid nested relationships in the data. This can cause prediction errors. To avoid this issue, explicitly specify valid combinations of hierarchical grouping variables in the datagrid() call, or use grid_type = "counterfactual" to preserve the original data structure.

mod <- lm(mpg ~ am + vs + factor(cyl) + hp, mtcars)
datagrid_bal <- function(...) datagrid(..., grid_type = "balanced")
predictions(model, newdata = datagrid_bal(cyl = 4))

If users supply a model, the data used to fit that model is retrieved using get_modeldata(). Consider using set_modeldata() to attach the training data explicitly when working in non-standard environments (e.g., lapply(), Shiny, nested functions).

Value

A data.frame in which each row corresponds to one combination of the named predictors supplied by the user via the ... dots. Variables which are not explicitly defined are held at their mean or mode.

Examples

# The output only has 2 rows, and all the variables except `hp` are at their
# mean or mode.
datagrid(newdata = mtcars, hp = c(100, 110))

# We get the same result by feeding a model instead of a data.frame
mod <- lm(mpg ~ hp, mtcars)
datagrid(model = mod, hp = c(100, 110))

# Use in `marginaleffects` to compute "Typical Marginal Effects". When used
# in `slopes()` or `predictions()` we do not need to specify the
# `model` or `newdata` arguments.
slopes(mod, newdata = datagrid(hp = c(100, 110)))

# datagrid accepts functions
datagrid(hp = range, cyl = unique, newdata = mtcars)
comparisons(mod, newdata = datagrid(hp = fivenum))

# The full dataset is duplicated with each observation given counterfactual
# values of 100 and 110 for the `hp` variable. The original `mtcars` includes
# 32 rows, so the resulting dataset includes 64 rows.
dg <- datagrid(newdata = mtcars, hp = c(100, 110), grid_type = "counterfactual")
nrow(dg)

# We get the same result by feeding a model instead of a data.frame
mod <- lm(mpg ~ hp, mtcars)
dg <- datagrid(model = mod, hp = c(100, 110), grid_type = "counterfactual")
nrow(dg)

# Use `by` to hold variables at group-specific values
mod2 <- lm(mpg ~ hp + cyl, mtcars)
datagrid(model = mod2, hp = mean, by = "cyl")

# Use `FUN` to apply function to all variables
datagrid(model = mod2, FUN = median)

# Use `grid_type="dataframe"` for column-wise binding instead of cross-product
datagrid(model = mod2, hp = c(100, 200), cyl = c(4, 6), grid_type = "dataframe")

Description

Downloads a dataset from the marginaleffects or the Rdatasets archives, and return it as a data frame. Opens the documentation as an HTML page. Search available datasets.

https://vincentarelbundock.github.io/Rdatasets/

Usage

get_dataset(dataset = "thornton", package = NULL, docs = FALSE, search = NULL)

Arguments

Value

A data frame containing the dataset.

Examples

dat <- get_dataset("Titanic", "Stat2Data")
head(dat)

get_dataset(search = "(?i)titanic")

# View documentation in the browser
get_dataset("Titanic", "Stat2Data", docs = TRUE)

Description

Extract Posterior Draws or Bootstrap Resamples from marginaleffects Objects

Usage

get_draws(x, shape = "long")

Arguments

Details

If DxP and PxD and the names returned by coef(x) are unique, marginaleffects sets parameter names to those names. Otherwise, it sets them to b1, b2, etc.

Value

A data.frame with drawid and draw columns.

Description

get_modeldata() retrieves the data used to fit a model. If the data was previously attached via set_modeldata(), it is returned directly. Otherwise, the data is retrieved from the model object or the environment using insight::get_data().

Usage

get_modeldata(model, ...)

set_modeldata(model, newdata)

Arguments

Details

set_modeldata() attaches a dataset to a model object as an attribute. This is useful when the training data may not be available in the environment where marginaleffects functions are called, such as inside lapply(), in Shiny apps, or in nested function calls. By attaching the data explicitly, you ensure that marginaleffects always uses the correct dataset.

Value

get_modeldata() returns a data frame of the original data used to fit the model, or NULL if not available. set_modeldata() returns the model with the data attached as an attribute.

Description

Uncertainty estimates are calculated as first-order approximate standard errors for linear or non-linear functions of a vector of random variables with known or estimated covariance matrix. In that sense, hypotheses emulates the behavior of the excellent and well-established car::deltaMethod and car::linearHypothesis functions, but it supports more models; requires fewer dependencies; expands the range of tests to equivalence and superiority/inferiority; and offers convenience features like robust standard errors.

To learn more, read the hypothesis tests vignette, visit the package website:

• https://marginaleffects.com/chapters/hypothesis.html

• https://marginaleffects.com/

Warning #1: Tests are conducted directly on the scale defined by the type argument. For some models, it can make sense to conduct hypothesis or equivalence tests on the "link" scale instead of the "response" scale which is often the default.

Warning #2: For hypothesis tests on objects produced by the marginaleffects package, it is safer to use the hypothesis argument of the original function. Using hypotheses() may not work in certain environments, in lists, or when working programmatically with *apply style functions.

Warning #3: The tests assume that the hypothesis expression is (approximately) normally distributed, which for non-linear functions of the parameters may not be realistic. More reliable confidence intervals can be obtained using the inferences() function with method = "boot".

Usage

hypotheses(
  model = NULL,
  hypothesis = NULL,
  vcov = TRUE,
  conf_level = NULL,
  df = NULL,
  equivalence = NULL,
  joint = FALSE,
  joint_test = "f",
  multcomp = FALSE,
  numderiv = "fdforward",
  ...
)

Arguments

Joint hypothesis tests

The test statistic for the joint Wald test is calculated as (R * theta_hat - r)' * inv(R * V_hat * R') * (R * theta_hat - r) / Q, where theta_hat is the vector of estimated parameters, V_hat is the estimated covariance matrix, R is a Q x P matrix for testing Q hypotheses on P parameters, r is a Q x 1 vector for the null hypothesis, and Q is the number of rows in R. If the test is a Chi-squared test, the test statistic is not normalized.

The p-value is then calculated based on either the F-distribution (for F-test) or the Chi-squared distribution (for Chi-squared test). For the F-test, the degrees of freedom are Q and (n - P), where n is the sample size and P is the number of parameters. For the Chi-squared test, the degrees of freedom are Q.

Equivalence, Inferiority, Superiority

$\theta$ is an estimate, $\sigma_\theta$ its estimated standard error, and $[a, b]$ are the bounds of the interval supplied to the equivalence argument.

Non-inferiority:

• $H_0$: $\theta \leq a$

• $H_1$: $\theta > a$

• $t=(\theta - a)/\sigma_\theta$

• p: Upper-tail probability

Non-superiority:

• $H_0$: $\theta \geq b$

• $H_1$: $\theta < b$

• $t=(\theta - b)/\sigma_\theta$

• p: Lower-tail probability

Equivalence: Two One-Sided Tests (TOST)

• p: Maximum of the non-inferiority and non-superiority p values.

Thanks to Russell V. Lenth for the excellent emmeans package and documentation which inspired this feature.

Examples

mod <- lm(mpg ~ hp + wt + factor(cyl), data = mtcars)

hypotheses(mod)

# Test of equality between coefficients
hypotheses(mod, hypothesis = "hp = wt")

# Non-linear function
hypotheses(mod, hypothesis = "exp(hp + wt) = 0.1")

# Robust standard errors
hypotheses(mod, hypothesis = "hp = wt", vcov = "HC3")

# b1, b2, ... shortcuts can be used to identify the position of the
# parameters of interest in the output of
hypotheses(mod, hypothesis = "b2 = b3")

# wildcard
hypotheses(mod, hypothesis = "b* / b2 = 1")

# term names with special characters have to be enclosed in backticks
hypotheses(mod, hypothesis = "`factor(cyl)6` = `factor(cyl)8`")

mod2 <- lm(mpg ~ hp * drat, data = mtcars)
hypotheses(mod2, hypothesis = "`hp:drat` = drat")

# predictions(), comparisons(), and slopes()
mod <- glm(am ~ hp + mpg, data = mtcars, family = binomial)
cmp <- comparisons(mod, newdata = "mean")
hypotheses(cmp, hypothesis = "b1 = b2")

mfx <- slopes(mod, newdata = "mean")
hypotheses(cmp, hypothesis = "b2 = 0.2")

pre <- predictions(mod, newdata = datagrid(hp = 110, mpg = c(30, 35)))
hypotheses(pre, hypothesis = "b1 = b2")

# The `hypothesis` argument can be used to compute standard errors for fitted values
mod <- glm(am ~ hp + mpg, data = mtcars, family = binomial)

f <- function(x) predict(x, type = "link", newdata = mtcars)
p <- hypotheses(mod, hypothesis = f)
head(p)

f <- function(x) predict(x, type = "response", newdata = mtcars)
p <- hypotheses(mod, hypothesis = f)
head(p)

# Complex aggregation
# Step 1: Collapse predicted probabilities by outcome level, for each individual
# Step 2: Take the mean of the collapsed probabilities by group and `cyl`
library(dplyr)
library(MASS)
library(dplyr)
library(magrittr)

dat <- transform(mtcars, gear = factor(gear))
mod <- polr(gear ~ factor(cyl) + hp, dat)

aggregation_fun <- function(x) {
    predictions(x, vcov = FALSE) %>%
        mutate(group = ifelse(group %in% c("3", "4"), "3 & 4", "5")) %>%
        summarize(estimate = sum(estimate), .by = c("rowid", "cyl", "group")) %>%
        summarize(estimate = mean(estimate), .by = c("cyl", "group")) %>%
        rename(term = cyl)
}

hypotheses(mod, hypothesis = aggregation_fun)

# Equivalence, non-inferiority, and non-superiority tests
mod <- lm(mpg ~ hp + factor(gear), data = mtcars)
p <- predictions(mod, newdata = "median")
hypotheses(p, equivalence = c(17, 18))

mfx <- avg_slopes(mod, variables = "hp")
hypotheses(mfx, equivalence = c(-.1, .1))

cmp <- avg_comparisons(mod, variables = "gear", hypothesis = ~pairwise)
hypotheses(cmp, equivalence = c(0, 10))

# joint hypotheses: character vector
model <- lm(mpg ~ as.factor(cyl) * hp, data = mtcars)
hypotheses(model, joint = c("as.factor(cyl)6:hp", "as.factor(cyl)8:hp"))

# joint hypotheses: regular expression
hypotheses(model, joint = "cyl")

# joint hypotheses: integer indices
hypotheses(model, joint = 2:3)

# joint hypotheses: different null hypotheses
hypotheses(model, joint = 2:3, hypothesis = 1)
hypotheses(model, joint = 2:3, hypothesis = 1:2)

# joint hypotheses: marginaleffects object
cmp <- avg_comparisons(model)
hypotheses(cmp, joint = "cyl")

# Multiple comparison adjustment
# p values and family-wise confidence intervals
cmp <- avg_comparisons(model)
hypotheses(cmp, multcomp = "hochberg")

Description

Warning: This function is experimental. It may be renamed, the user interface may change, or the functionality may migrate to arguments in other marginaleffects functions.

Apply this function to a marginaleffects object to change the inferential method used to compute uncertainty estimates.

Usage

inferences(
  x,
  method,
  R = 1000,
  conf_type = "perc",
  data_train = NULL,
  data_test = NULL,
  data_calib = NULL,
  conformal_score = "residual_abs",
  estimator = NULL,
  ...
)

Arguments

Details

When method = "simulation", we conduct simulation-based inference following the method discussed in Krinsky & Robb (1986):

1. Draw R sets of simulated coefficients from a multivariate normal distribution with mean equal to the original model's estimated coefficients and variance equal to the model's variance-covariance matrix (classical, "HC3", or other).

2. Use the R sets of coefficients to compute R sets of estimands: predictions, comparisons, slopes, or hypotheses.

3. Take quantiles of the resulting distribution of estimands to obtain a confidence interval (when conf_type = "perc") and the standard deviation of simulated estimates to estimate the standard error (which is used for a Z-test and Wald confidence intervals when conf_type = "wald").

When method = "fwb", drawn weights are supplied to the model fitting function's weights argument; if the model doesn't accept non-integer weights, this method should not be used. If weights were included in the original model fit, they are extracted by weights() and multiplied by the drawn weights. These weights are supplied to the wts argument of the estimation function (e.g., comparisons()).

Warning: custom model classes are not supported by inferences() because they are not guaranteed to come with an appropriate update() method.

Value

A marginaleffects object with simulation or bootstrap resamples and objects attached.

References

Krinsky, I., and A. L. Robb. 1986. "On Approximating the Statistical Properties of Elasticities." Review of Economics and Statistics 68 (4): 715–9.

King, Gary, Michael Tomz, and Jason Wittenberg. "Making the most of statistical analyses: Improving interpretation and presentation." American journal of political science (2000): 347-361

Dowd, Bryan E., William H. Greene, and Edward C. Norton. "Computation of standard errors." Health services research 49.2 (2014): 731-750.

Angelopoulos, Anastasios N., and Stephen Bates. 2022. "A Gentle Introduction to Conformal Prediction and Distribution-Free Uncertainty Quantification." arXiv. https://doi.org/10.48550/arXiv.2107.07511.

Barber, Rina Foygel, Emmanuel J. Candes, Aaditya Ramdas, and Ryan J. Tibshirani. 2020. "Predictive Inference with the Jackknife+." arXiv. http://arxiv.org/abs/1905.02928.

Lei, Jing, Max G'Sell, Alessandro Rinaldo, Ryan J. Tibshirani, and Larry Wasserman. 2018. "Distribution-Free Predictive Inference for Regression." Journal of the American Statistical Association 113 (523): 1094–1111.

Romano, Yaniv, Evan Patterson, and Emmanuel Candes. 2020. "Conformalized quantile regression." Advances in neural information processing systems 32.

Parallel computation

The slopes() and comparisons() functions can use parallelism to speed up computation. Operations are parallelized for the computation of standard errors, at the model coefficient level. There is always considerable overhead when using parallel computation, mainly involved in passing the whole dataset to the different processes. Thus, parallel computation is most likely to be useful when the model includes many parameters and the dataset is relatively small.

Warning: In many cases, parallel processing will not be useful at all.

To activate parallel computation, users must load the future.apply package, call plan() function, and set a global option.

options(marginaleffects_parallel = TRUE): parallelize delta method computation of standard errors. options(marginaleffects_parallel_inferences = TRUE): parallelize "rsample" or "fwb" bootstrap computation in inferences(). options(marginaleffects_parallel_packages = TRUE): vector of strings with the names of modeling packages used to fit the model, ex: c("survival", "splines")

For example:

library(future.apply)
plan("multisession", workers = 4)
options(marginaleffects_parallel = FALSE)
options(marginaleffects_parallel_inferences = TRUE)
options(marginaleffects_parallel_packages = c("survival", "splines"))

slopes(model)

To disable parallelism in marginaleffects altogether, you can set a global option:

options(marginaleffects_parallel = FALSE)

Examples

## Not run: 
library(magrittr)
set.seed(1024)
mod <- lm(Sepal.Length ~ Sepal.Width * Species, data = iris)

# bootstrap
avg_predictions(mod, by = "Species") %>%
    inferences(method = "boot")

avg_predictions(mod, by = "Species") %>%
    inferences(method = "rsample")

# Fractional (bayesian) bootstrap
avg_slopes(mod, by = "Species") %>%
    inferences(method = "fwb") %>%
    get_draws("rvar") %>%
    data.frame()

# Simulation-based inference
slopes(mod) %>%
    inferences(method = "simulation") %>%
    head()

# Two-step estimation procedure: Propensity score + G-Computation
lalonde <- get_dataset("lalonde")
estimator <- function(data) {
    # Step 1: Estimate propensity scores
    fit1 <- glm(treat ~ age + educ + race, family = binomial, data = data)
    ps <- predict(fit1, type = "response")
    # Step 2: Fit weighted outcome model
    m <- lm(re78 ~ treat * (re75 + age + educ + race),
        data = data, weight = ps
    )
    # Step 3: Compute average treatment effect by G-computation
    avg_comparisons(m, variables = "treat", wts = ps, vcov = FALSE)
}
inferences(lalonde, method = "rsample", estimator = estimator)

## End(Not run)

Description

Plot comparisons on the y-axis against values of one or more predictors (x-axis, colors/shapes, and facets).

The by argument is used to plot marginal comparisons, that is, comparisons made on the original data, but averaged by subgroups. This is analogous to using the by argument in the comparisons() function.

The condition argument is used to plot conditional comparisons, that is, comparisons made on a user-specified grid. This is analogous to using the newdata argument and datagrid() function in a comparisons() call. All variables whose values are not specified explicitly are treated as usual by datagrid(), that is, they are held at their mean or mode (or rounded mean for integers). This includes grouping variables in mixed-effects models, so analysts who fit such models may want to specify the groups of interest using the condition argument, or supply model-specific arguments to compute population-level estimates. See details below.

See the "Plots" vignette and website for tutorials and information on how to customize plots:

• https://marginaleffects.com/bonus/plot.html

• https://marginaleffects.com

Usage

plot_comparisons(
  model,
  variables = NULL,
  condition = NULL,
  by = NULL,
  newdata = NULL,
  type = NULL,
  vcov = NULL,
  conf_level = 0.95,
  wts = FALSE,
  comparison = "difference",
  transform = NULL,
  rug = FALSE,
  gray = getOption("marginaleffects_plot_gray", default = FALSE),
  draw = TRUE,
  ...
)

Arguments

Value

A ggplot2 object

Model-Specific Arguments

Some model types allow model-specific arguments to modify the nature of marginal effects, predictions, marginal means, and contrasts. Please report other package-specific predict() arguments on Github so we can add them to the table below.

https://github.com/vincentarelbundock/marginaleffects/issues

Examples

mod <- lm(mpg ~ hp * drat * factor(am), data = mtcars)

plot_comparisons(mod, variables = "hp", condition = "drat")

plot_comparisons(mod, variables = "hp", condition = c("drat", "am"))

plot_comparisons(mod, variables = "hp", condition = list("am", "drat" = 3:5))

plot_comparisons(mod, variables = "am", condition = list("hp", "drat" = range))

plot_comparisons(mod, variables = "am", condition = list("hp", "drat" = "threenum"))

# marginal comparisons
plot_comparisons(mod, variables = "hp", by = "am")

# marginal comparisons on a counterfactual grid
plot_comparisons(mod,
    variables = "hp",
    by = "am",
    newdata = datagrid(am = 0:1, grid_type = "counterfactual")
)

Description

Plot predictions on the y-axis against values of one or more predictors (x-axis, colors/shapes, and facets).

The by argument is used to plot marginal predictions, that is, predictions made on the original data, but averaged by subgroups. This is analogous to using the by argument in the predictions() function.

The condition argument is used to plot conditional predictions, that is, predictions made on a user-specified grid. This is analogous to using the newdata argument and datagrid() function in a predictions() call. All variables whose values are not specified explicitly are treated as usual by datagrid(), that is, they are held at their mean or mode (or rounded mean for integers). This includes grouping variables in mixed-effects models, so analysts who fit such models may want to specify the groups of interest using the condition argument, or supply model-specific arguments to compute population-level estimates. See details below.

See the "Plots" vignette and website for tutorials and information on how to customize plots:

• https://marginaleffects.com/bonus/plot.html

• https://marginaleffects.com

Usage

plot_predictions(
  model,
  condition = NULL,
  by = NULL,
  newdata = NULL,
  type = NULL,
  vcov = NULL,
  conf_level = 0.95,
  wts = FALSE,
  transform = NULL,
  points = 0,
  rug = FALSE,
  gray = getOption("marginaleffects_plot_gray", default = FALSE),
  draw = TRUE,
  ...
)

Arguments

Value

A ggplot2 object or data frame (if draw=FALSE)

Model-Specific Arguments

Some model types allow model-specific arguments to modify the nature of marginal effects, predictions, marginal means, and contrasts. Please report other package-specific predict() arguments on Github so we can add them to the table below.

https://github.com/vincentarelbundock/marginaleffects/issues

Types

The type argument determines the scale of the predictions used to compute quantities of interest with functions from the marginaleffects package. Admissible values for type depend on the model object. When users specify an incorrect value for type, marginaleffects will raise an informative error with a list of valid type values for the specific model object. The first entry in the list in that error message is the default type.

The invlink(link) is a special type defined by marginaleffects. It is available for some (but not all) models, and only for the predictions() function. With this link type, we first compute predictions on the link scale, then we use the inverse link function to backtransform the predictions to the response scale. This is useful for models with non-linear link functions as it can ensure that confidence intervals stay within desirable bounds, ex: 0 to 1 for a logit model. Note that an average of estimates with type="invlink(link)" will not always be equivalent to the average of estimates with type="response". This type is default when calling predictions(). It is available—but not default—when calling avg_predictions() or predictions() with the by argument.

Some of the most common type values are:

Examples

mod <- lm(mpg ~ hp + wt, data = mtcars)
plot_predictions(mod, condition = "wt")

mod <- lm(mpg ~ hp * wt * am, data = mtcars)
plot_predictions(mod, condition = c("hp", "wt"))

plot_predictions(mod, condition = list("hp", wt = "threenum"))

plot_predictions(mod, condition = list("hp", wt = range))

# marginal predictions
mod <- lm(mpg ~ hp * am, data = mtcars)
plot_predictions(mod, by = "am")

# marginal predictions on a counterfactual grid
plot_predictions(mod,
    by = "am",
    newdata = datagrid(am = 0:1, grid_type = "counterfactual")
)

Description

Plot slopes on the y-axis against values of one or more predictors (x-axis, colors/shapes, and facets).

The by argument is used to plot marginal slopes, that is, slopes made on the original data, but averaged by subgroups. This is analogous to using the by argument in the slopes() function.

The condition argument is used to plot conditional slopes, that is, slopes computed on a user-specified grid. This is analogous to using the newdata argument and datagrid() function in a slopes() call. All variables whose values are not specified explicitly are treated as usual by datagrid(), that is, they are held at their mean or mode (or rounded mean for integers). This includes grouping variables in mixed-effects models, so analysts who fit such models may want to specify the groups of interest using the condition argument, or supply model-specific arguments to compute population-level estimates. See details below. See the "Plots" vignette and website for tutorials and information on how to customize plots:

• https://marginaleffects.com/bonus/plot.html

• https://marginaleffects.com

Usage

plot_slopes(
  model,
  variables = NULL,
  condition = NULL,
  by = NULL,
  newdata = NULL,
  type = NULL,
  vcov = NULL,
  conf_level = 0.95,
  wts = FALSE,
  slope = "dydx",
  rug = FALSE,
  gray = getOption("marginaleffects_plot_gray", default = FALSE),
  draw = TRUE,
  ...
)

Arguments

Value

A ggplot2 object

Model-Specific Arguments

Some model types allow model-specific arguments to modify the nature of marginal effects, predictions, marginal means, and contrasts. Please report other package-specific predict() arguments on Github so we can add them to the table below.

https://github.com/vincentarelbundock/marginaleffects/issues

Examples

mod <- lm(mpg ~ hp * drat * factor(am), data = mtcars)

plot_slopes(mod, variables = "hp", condition = "drat")

plot_slopes(mod, variables = "hp", condition = c("drat", "am"))

plot_slopes(mod, variables = "hp", condition = list("am", "drat" = 3:5))

plot_slopes(mod, variables = "am", condition = list("hp", "drat" = range))

plot_slopes(mod, variables = "am", condition = list("hp", "drat" = "threenum"))

# marginal slopes
plot_slopes(mod, variables = "hp", by = "am")

# marginal slopes on a counterfactual grid
plot_slopes(mod,
    variables = "hp",
    by = "am",
    newdata = datagrid(am = 0:1, grid_type = "counterfactual")
)

Description

Outcome predicted by a fitted model on a specified scale for a given combination of values of the predictor variables, such as their observed values, their means, or factor levels (a.k.a. "reference grid").

• predictions(): unit-level (conditional) estimates.

• avg_predictions(): average (marginal) estimates.

The newdata argument and the datagrid() function can be used to control where statistics are evaluated in the predictor space: "at observed values", "at the mean", "at representative values", etc.

See the predictions vignette and package website for worked examples and case studies:

• https://marginaleffects.com/chapters/predictions.html

• https://marginaleffects.com/

Usage

predictions(
  model,
  newdata = NULL,
  variables = NULL,
  vcov = TRUE,
  conf_level = 0.95,
  type = NULL,
  by = FALSE,
  byfun = NULL,
  wts = FALSE,
  transform = NULL,
  hypothesis = NULL,
  equivalence = NULL,
  df = Inf,
  numderiv = "fdforward",
  ...
)

avg_predictions(
  model,
  newdata = NULL,
  variables = NULL,
  vcov = TRUE,
  conf_level = 0.95,
  type = NULL,
  by = TRUE,
  byfun = NULL,
  wts = FALSE,
  transform = NULL,
  hypothesis = NULL,
  equivalence = NULL,
  df = Inf,
  numderiv = "fdforward",
  ...
)

Arguments

Value

A data.frame with one row per estimate. This data frame is pretty-printed by default, but users can interact with it as a regular data frame, with functions like nrow(), head(), colnames(), etc. Values can be extracted using standard ⁠[,]⁠ or $ operators, and manipulated using external packages like dplyr or data.table.

Columns may include:

• rowid: row number of the newdata data frame

• group: (optional) value of the grouped outcome (e.g., categorical outcome models)

• term: the focal variable.

• estimate: an estimate of the prediction, counterfactual comparison, or slope.

• std.error: standard errors computed via the delta method.

• p.value: p value associated to the estimate column. The null is determined by the hypothesis argument (0 by default).

• s.value: Shannon information transforms of p values. See the S values vignette at https://marginaleffects.com the marginaleffects website.

• conf.low: lower bound of the confidence (or credible) interval defined by the conf_level argument.

• conf.high: upper bound of the confidence (or credible) interval defined by the conf_level argument.

• predicted_lo: predicted outcome for the "low" value of the focal predictor in a counterfactual comparison.

• predicted_hi: predicted outcome for the "high" value of the focal predictor in a counterfactual comparison.

• p.rope.unconditional: share of posterior draws in the interval specified by the equivalence argument. This is only available for Bayesian models.

• p.rope.conditional: share of posterior draws in the interval specified by the equivalence argument, among draws in the confidence interval. This is only available for Bayesian models.

• rope: share of the posterior draws between conf.low and conf.high that are covered by the interval specified by the equivalence argument.

• statistic.noninf: test statistic for non-inferiority test (when equivalence argument is used).

• statistic.nonsup: test statistic for non-superiority test (when equivalence argument is used).

• p.value.noninf: p-value for non-inferiority test (when equivalence argument is used).

• p.value.nonsup: p-value for non-superiority test (when equivalence argument is used).

• p.value.equiv: p-value for equivalence test using Two One-Sided Tests (TOST) approach (when equivalence argument is used).

See ?print.marginaleffects for printing options.

The data.frames produced by marginaleffects stores an attribute that holds many internal objects, such as the original model, data, and much other information that can be used for post-processing. This information can be inspected using the components() function.

Warning: The internal attributes retrieved by components() are not considered part of the public API of the package. Their names and contents can change without warning or notice. Users should not rely on them.

Warning: In some cases, the internal attributes used by marginaleffects() can use up a substantial amount of memory. To clear this data, use the prune() function or set options(marginaleffects_lean=TRUE).

Functions

• avg_predictions(): Average predictions

Standard errors using the delta method

Standard errors for all quantities estimated by marginaleffects can be obtained via the delta method. This requires differentiating a function with respect to the coefficients in the model using a finite difference approach. In some models, the delta method standard errors can be sensitive to various aspects of the numeric differentiation strategy, including the step size. By default, the step size is set to 1e-8, or to 1e-4 times the smallest absolute model coefficient, whichever is largest.

marginaleffects can delegate numeric differentiation to the numDeriv package, which allows more flexibility. To do this, users can pass arguments to the numDeriv::jacobian function through a global option. For example:

• options(marginaleffects_numDeriv = list(method = "simple", method.args = list(eps = 1e-6)))

• options(marginaleffects_numDeriv = list(method = "Richardson", method.args = list(eps = 1e-5)))

• options(marginaleffects_numDeriv = NULL)

See the "Uncertainty" chapter on the marginaleffects website for more details on the computation of standard errors, bootstrapping, and more:

https://marginaleffects.com/chapters/uncertainty.html

Model-Specific Arguments

Some model types allow model-specific arguments to modify the nature of marginal effects, predictions, marginal means, and contrasts. Please report other package-specific predict() arguments on Github so we can add them to the table below.

https://github.com/vincentarelbundock/marginaleffects/issues

Bayesian posterior summaries

By default, credible intervals in bayesian models are built as equal-tailed intervals. This can be changed to a highest density interval by setting a global option:

options("marginaleffects_posterior_interval" = "eti")

options("marginaleffects_posterior_interval" = "hdi")

By default, the center of the posterior distribution in bayesian models is identified by the median. Users can use a different summary function by setting a global option:

options("marginaleffects_posterior_center" = "mean")

options("marginaleffects_posterior_center" = "median")

When estimates are averaged using the by argument, the tidy() function, or the summary() function, the posterior distribution is marginalized twice over. First, we take the average across units but within each iteration of the MCMC chain, according to what the user requested in by argument or tidy()/summary() functions. Then, we identify the center of the resulting posterior using the function supplied to the "marginaleffects_posterior_center" option (the median by default).

Equivalence, Inferiority, Superiority

$\theta$ is an estimate, $\sigma_\theta$ its estimated standard error, and $[a, b]$ are the bounds of the interval supplied to the equivalence argument.

Non-inferiority:

• $H_0$: $\theta \leq a$

• $H_1$: $\theta > a$

• $t=(\theta - a)/\sigma_\theta$

• p: Upper-tail probability

Non-superiority:

• $H_0$: $\theta \geq b$

• $H_1$: $\theta < b$

• $t=(\theta - b)/\sigma_\theta$

• p: Lower-tail probability

Equivalence: Two One-Sided Tests (TOST)

• p: Maximum of the non-inferiority and non-superiority p values.

Thanks to Russell V. Lenth for the excellent emmeans package and documentation which inspired this feature.

Order of operations

Behind the scenes, the arguments of marginaleffects functions are evaluated in this order:

1. newdata

2. variables

3. comparison and slope

4. by

5. vcov

6. hypothesis

7. transform

Parallel computation

The slopes() and comparisons() functions can use parallelism to speed up computation. Operations are parallelized for the computation of standard errors, at the model coefficient level. There is always considerable overhead when using parallel computation, mainly involved in passing the whole dataset to the different processes. Thus, parallel computation is most likely to be useful when the model includes many parameters and the dataset is relatively small.

Warning: In many cases, parallel processing will not be useful at all.

To activate parallel computation, users must load the future.apply package, call plan() function, and set a global option.

options(marginaleffects_parallel = TRUE): parallelize delta method computation of standard errors. options(marginaleffects_parallel_inferences = TRUE): parallelize "rsample" or "fwb" bootstrap computation in inferences(). options(marginaleffects_parallel_packages = TRUE): vector of strings with the names of modeling packages used to fit the model, ex: c("survival", "splines")

For example:

library(future.apply)
plan("multisession", workers = 4)
options(marginaleffects_parallel = FALSE)
options(marginaleffects_parallel_inferences = TRUE)
options(marginaleffects_parallel_packages = c("survival", "splines"))

slopes(model)

To disable parallelism in marginaleffects altogether, you can set a global option:

options(marginaleffects_parallel = FALSE)

Global options

The behavior of marginaleffects functions can be modified by setting global options.

Disable some safety checks and warnings:

• options(marginaleffects_startup_message = FALSE)

  • Disable the startup message printed on library(marginaleffects).

• options(marginaleffects_safe = FALSE)

  • Disable safety checks and warnings.

• options(marginaleffects_print_omit = c("p.value", "s.value"))

  • Omit some columns from the printed output.

Enforce lean return objects, sans information about the original model and data, and other ancillary attributes. Note that this will disable some advanced post-processing features and functions like hypotheses.

options(marginaleffects_lean = TRUE)

Other options:

• marginaleffects_plot_gray: logical. If TRUE, the default color of the plot is gray. Default is FALSE.

Types

The type argument determines the scale of the predictions used to compute quantities of interest with functions from the marginaleffects package. Admissible values for type depend on the model object. When users specify an incorrect value for type, marginaleffects will raise an informative error with a list of valid type values for the specific model object. The first entry in the list in that error message is the default type.

The invlink(link) is a special type defined by marginaleffects. It is available for some (but not all) models, and only for the predictions() function. With this link type, we first compute predictions on the link scale, then we use the inverse link function to backtransform the predictions to the response scale. This is useful for models with non-linear link functions as it can ensure that confidence intervals stay within desirable bounds, ex: 0 to 1 for a logit model. Note that an average of estimates with type="invlink(link)" will not always be equivalent to the average of estimates with type="response". This type is default when calling predictions(). It is available—but not default—when calling avg_predictions() or predictions() with the by argument.

Some of the most common type values are:

References

• Arel-Bundock V, Greifer N, Heiss A (2024). “How to Interpret Statistical Models Using marginaleffects for R and Python.” Journal of Statistical Software, 111(9), 1-32. doi:10.18637/jss.v111.i09 doi:10.18637/jss.v111.i09

• Arel-Bundock (2026). "Model to Meaning: How to interpret statistical models in R and Python." CRC Press. https://routledge.com/9781032908724

• Greenland S. 2019. "Valid P-Values Behave Exactly as They Should: Some Misleading Criticisms of P-Values and Their Resolution With S-Values." The American Statistician. 73(S1): 106–114.

• Cole, Stephen R, Jessie K Edwards, and Sander Greenland. 2020. "Surprise!" American Journal of Epidemiology 190 (2): 191–93. doi:10.1093/aje/kwaa136

Examples

library("marginaleffects")
# Adjusted Prediction for every row of the original dataset
mod <- lm(mpg ~ hp + factor(cyl), data = mtcars)
pred <- predictions(mod)
head(pred)

# Adjusted Predictions at User-Specified Values of the Regressors
predictions(mod, newdata = datagrid(hp = c(100, 120), cyl = 4))

m <- lm(mpg ~ hp + drat + factor(cyl) + factor(am), data = mtcars)
predictions(m, newdata = datagrid(FUN_factor = unique, FUN_numeric = median))

# Average Adjusted Predictions (AAP)
library(dplyr)
mod <- lm(mpg ~ hp * am * vs, mtcars)

avg_predictions(mod)

predictions(mod, by = "am")

# Conditional Adjusted Predictions
plot_predictions(mod, condition = "hp")

# Counterfactual predictions with the `variables` argument
# the `mtcars` dataset has 32 rows

mod <- lm(mpg ~ hp + am, data = mtcars)
p <- predictions(mod)
head(p)
nrow(p)

# average counterfactual predictions
avg_predictions(mod, variables = "am")

# counterfactual predictions obtained by replicating the entire for different
# values of the predictors
p <- predictions(mod, variables = list(hp = c(90, 110)))
nrow(p)


# hypothesis test: is the prediction in the 1st row equal to the prediction in the 2nd row
mod <- lm(mpg ~ wt + drat, data = mtcars)

predictions(
  mod,
  newdata = datagrid(wt = 2:3),
  hypothesis = "b1 = b2")

# same hypothesis test using row indices
predictions(
  mod,
  newdata = datagrid(wt = 2:3),
  hypothesis = "b1 - b2 = 0")

# same hypothesis test using numeric vector of weights
predictions(
  mod,
  newdata = datagrid(wt = 2:3),
  hypothesis = c(1, -1))

# two custom contrasts using a matrix of weights
lc <- matrix(
  c(
    1, -1,
    2, 3),
  ncol = 2)
predictions(
  mod,
  newdata = datagrid(wt = 2:3),
  hypothesis = lc)


# `by` argument
mod <- lm(mpg ~ hp * am * vs, data = mtcars)
predictions(mod, by = c("am", "vs"))

library(nnet)
nom <- multinom(factor(gear) ~ mpg + am * vs, data = mtcars, trace = FALSE)

# first 5 raw predictions
p <- predictions(nom, type = "probs")
head(p)

# average predictions
avg_predictions(nom, type = "probs", by = "group")

by <- data.frame(
  group = c("3", "4", "5"),
  by = c("3,4", "3,4", "5"))

predictions(nom, type = "probs", by = by)

# sum of predicted probabilities for combined response levels
mod <- multinom(factor(cyl) ~ mpg + am, data = mtcars, trace = FALSE)
by <- data.frame(
  by = c("4,6", "4,6", "8"),
  group = as.character(c(4, 6, 8)))
predictions(mod, newdata = "mean", byfun = sum, by = by)

Description

This function controls the text which is printed to the console when one of the core marginalefffects functions is called and the object is returned: predictions(), comparisons(), slopes(), hypotheses(), avg_predictions(), avg_comparisons(), avg_slopes().

All of those functions return standard data frames. Columns can be extracted by name, predictions(model)$estimate, and all the usual data manipulation functions work out-of-the-box: colnames(), head(), subset(), dplyr::filter(), dplyr::arrange(), etc.

Some of the data columns are not printed by default. You can disable pretty printing and print the full results as a standard data frame using the style argument or by applying as.data.frame() on the object. See examples below.

Usage

## S3 method for class 'marginaleffects'
print(
  x,
  style = getOption("marginaleffects_print_style", default = "summary"),
  digits = getOption("marginaleffects_print_digits", default = 3),
  p_eps = getOption("marginaleffects_print_p_eps", default = 0.001),
  topn = getOption("marginaleffects_print_topn", default = 5),
  nrows = getOption("marginaleffects_print_nrows", default = 30),
  ncols = getOption("marginaleffects_print_ncols", default = 30),
  type = getOption("marginaleffects_print_type", default = TRUE),
  column_names = getOption("marginaleffects_print_column_names", default = FALSE),
  ...
)

Arguments

Examples

library(marginaleffects)
mod <- lm(mpg ~ hp + am + factor(gear), data = mtcars)
p <- predictions(mod, by = c("am", "gear"))
p

subset(p, am == 1)

print(p, style = "data.frame")

data.frame(p)

Description

Remove large attributes from marginaleffects objects to reduce memory usage. Warning: This will disable many useful post-processing features of marginaleffects

Usage

## S3 method for class 'marginaleffects'
prune(tree, component = "all", ...)

Arguments

Details

...

Value

A pruned marginaleffects object

Description

Refit a marginaleffects object with new data

Usage

refit(object, ...)

## S3 method for class 'marginaleffects'
refit(object, data = NULL, newdata = NULL, vcov = NULL, ...)

## S3 method for class 'predictions'
refit(object, data = NULL, newdata = NULL, vcov = NULL, ...)

## S3 method for class 'comparisons'
refit(object, data = NULL, newdata = NULL, vcov = NULL, ...)

## S3 method for class 'slopes'
refit(object, data = NULL, newdata = NULL, vcov = NULL, ...)

## S3 method for class 'hypotheses'
refit(object, data = NULL, newdata = NULL, vcov = NULL, ...)

Arguments

Details

If data is supplied, the underlying model is refitted using that data. If newdata is supplied, the marginaleffects call is re-evaluated with the new data. Both can be supplied together to refit the model and make predictions on new data.

Value

A marginaleffects object

Description

Partial derivative of the regression equation with respect to a regressor of interest.

• slopes(): unit-level (conditional) estimates.

• avg_slopes(): average (marginal) estimates.

The newdata argument and the datagrid() function can be used to control where statistics are evaluated in the predictor space: "at observed values", "at the mean", "at representative values", etc.

See the slopes vignette and package website for worked examples and case studies:

• https://marginaleffects.com/chapters/slopes.html

• https://marginaleffects.com/

Warning: Slopes and elasticities can only be calculated for continuous numeric variables. The slopes() functions will automatically revert to comparisons() for binary or categorical variables.

Usage

slopes(
  model,
  newdata = NULL,
  variables = NULL,
  type = NULL,
  by = FALSE,
  vcov = TRUE,
  conf_level = 0.95,
  slope = "dydx",
  wts = FALSE,
  hypothesis = NULL,
  equivalence = NULL,
  df = Inf,
  eps = NULL,
  numderiv = "fdforward",
  ...
)

avg_slopes(
  model,
  newdata = NULL,
  variables = NULL,
  type = NULL,
  by = TRUE,
  vcov = TRUE,
  conf_level = 0.95,
  slope = "dydx",
  wts = FALSE,
  hypothesis = NULL,
  equivalence = NULL,
  df = Inf,
  eps = NULL,
  numderiv = "fdforward",
  ...
)