Learn R Programming

bayesics (version 2.0.2)

mediate_b: Mediation using Bayesian methods

Description

Mediation analysis done in the framework of Imai et al. (2010). Currently only applicable to linear models.

Usage

mediate_b(
  model_m,
  model_y,
  treat,
  control_value,
  treat_value,
  n_draws = 500,
  ask_before_full_sampling = TRUE,
  CI_level = 0.95,
  seed = 1,
  mc_error = ifelse("glm_b" %in% model_y, 0.01, 0.002),
  batch_size = 500
)

Value

A list with the following elements:

  • summary - tibble giving results for causal mediation quantities

  • posterior_draws (of counterfactual expectations)

  • mc_error absolute error used, including any rescaling to match the scale of the outcome

  • other inputs to mediate_b

Arguments

model_m

a fitted model object of class lm_b for mediator.

model_y

a fitted model object of class lm_b for outcome.

treat

a character string indicating the name of the treatment variable used in the models. NOTE: Treatment variable must be numeric (even if it's 1's and 0's).

control_value

value of the treatment variable used as the control condition. Default is the 1st quintile of the treat variable.

treat_value

value of the treatment variable used as the treatment condition. Default is the 4th quintile of the treat variable.

n_draws

Number of preliminary posterior draws to assess final number of posterior draws required for accurate interval estimation

ask_before_full_sampling

logical. If FALSE, the user will not be asked if they want to complete the full sampling. Defaults to TRUE, as this can be a computationally intensive procedure.

CI_level

numeric. Credible interval level.

seed

integer. Always set your seed!!!

mc_error

positive scalar. The number of posterior samples will, with high probability, estimate the CI bounds up to \(\pm\)mc_error\(\times\)sd(y).

batch_size

positive integer. Number of posterior draws to be taken at once. Higher values are more computationally intensive, but values which are too high might take up significant memory (allocates on the order of batch_size\(\times\)nrow(model_y$data)).

Details

The model is the same as that of Imai et al. (2010): $$ M_i(X) = w_i'\alpha_m + X\beta_m + \epsilon_{m,i}, \\ y_i(X, M(\tilde X)) = w_i'\alpha_y + X\beta_y + M(\tilde X)\gamma + \epsilon_{y,i}, \\ \epsilon_{m,i} \overset{iid}{\sim} N(0,\sigma^2_m), \\ \epsilon_{y,i} \overset{iid}{\sim} N(0,\sigma^2_y), \\ $$ where \(M_i(X)\) is the mediator as a function of the treatment variable \(X\), and \(w_i\) are confounder covariates.

Unlike the mediation R package, the estimation in mediate_b is fully Bayesian (as opposed to "quasi-Bayesian").

References

Imai, Kosuke, et al. “A General Approach to Causal Mediation Analysis.” Psychological Methods, vol. 15, no. 4, 2010, pp. 309–34, https://doi.org/10.1037/a0020761.

Examples

Run this code
# \donttest{
# Simplest case
## Generate some data
set.seed(2025)
N = 500
test_data = 
  data.frame(tr = rnorm(N),
             x1 = rnorm(N))
test_data$m = 
  rnorm(N, 0.4 * test_data$tr - 0.25 * test_data$x1)
test_data$outcome = 
  rnorm(N,-1 + 0.6 * test_data$tr + 1.5 * test_data$m + 0.25 * test_data$x1)

## Fit the mediator and outcome models
m1 = 
  lm_b(m ~ tr + x1,
       data = test_data)
m2 = 
  lm_b(outcome ~ m + tr + x1,
       data = test_data)
## Estimate the causal mediation quantities
m3 <-
  mediate_b(m1,m2,
            treat = "tr",
            control_value = -2,
            treat_value = 2,
            n_draws = 500,
            mc_error = 0.05,
            ask_before_full_sampling = FALSE)
m3
summary(m3,
        CI_level = 0.9)

# More complicated scenario
## Generate some data
set.seed(2025)
N = 500
test_data = 
  data.frame(tr = rep(0:1,N/2),
             x1 = rnorm(N))
test_data$m = 
  rnorm(N, 0.4 * test_data$tr - 0.25 * test_data$x1)
test_data$outcome = 
  rpois(N,exp(-1 + 0.6 * test_data$tr + 1.5 * test_data$m + 0.25 * test_data$x1))

## Fit the mediator and outcome models
m1 = 
  lm_b(m ~ tr + x1,
       data = test_data)
m2 = 
  glm_b(outcome ~ m + tr + x1,
        data = test_data,
        family = poisson())

##  Estimate the causal mediation quantities
m3 <-
  mediate_b(m1,m2,
            treat = "tr",
            control_value = 0,
            treat_value = 1,
            n_draws = 500,
            mc_error = 0.05,
            ask_before_full_sampling = FALSE)
summary(m3)
# }



Run the code above in your browser using DataLab