Package {causalsim}

causalsim: Simulation-Ready Causal Data Generating Processes

Description

Construct, simulate, and evaluate causal data generating processes (DGPs) with known ground truth. Covariate roles (confounder, effect modifier, noise) and heterogeneous treatment effects are first-class concepts, and estimator performance is summarised with bias, RMSE, confidence-interval coverage, and power. Results are returned as plain data frames for further analysis or plotting.

Author(s)

Maintainer: Chayce Reed Chayce.Reed.HSE@dartmouth.edu

See Also

Useful links:

• https://chaycereed.github.io/causalsim/

• https://github.com/chaycereed/causalsim

• Report bugs at https://github.com/chaycereed/causalsim/issues

Simulate a causal dataset with known ground truth

Description

A convenience wrapper around causalsim_dgp() and causalsim_draw(). Defines a causal data generating process and returns one simulated dataset in a single call.

Usage

causalsim(
  n,
  effect = 1,
  propensity = "moderate",
  baseline = "moderate",
  sigma = 1,
  covariates = list(),
  n_confounders = 0L,
  n_effect_modifiers = 0L,
  n_noise = 0L,
  mc_draws = 10000L,
  seed = NULL
)

Arguments

Details

The returned data frame includes all covariate columns, the treatment indicator A, the outcome Y, and two ground-truth columns:

.tau

Individual-level causal effect (CATE).

.p

True propensity score.

For full control over the DGP (multiple draws, grid evaluation, estimator benchmarking), use causalsim_dgp() and causalsim_draw() directly.

Value

A data frame with covariate columns, A, Y, .tau, and .p.

Examples

# One confounder, constant effect of 2, moderate confounding
data <- causalsim(n = 500, n_confounders = 1, effect = 2, seed = 1L)
head(data)

# Heterogeneous effect with an explicit covariate spec
data2 <- causalsim(
  n = 500,
  covariates = list(
    W = causalsim_covar("normal", role = "confounder"),
    V = causalsim_covar("binary", role = "effect_modifier", prob = 0.4)
  ),
  effect = function(V) 2 + 1.5 * V,
  propensity = function(W) plogis(0.5 * W),
  seed = 42L
)
head(data2)

Define a covariate for a causal DGP

Description

Constructs a fully-specified covariate descriptor for use in causalsim_dgp(). Every covariate has a distribution family, optional distribution parameters, and one or more causal roles that control how it enters the automatically-generated treatment and outcome models.

Usage

causalsim_covar(dist = "normal", role = "confounder", ...)

Arguments

Value

An S3 object of class causalsim_covar.

Examples

# Standard normal confounder
causalsim_covar("normal", role = "confounder", mean = 0, sd = 1)

# Binary effect modifier
causalsim_covar("binary", role = "effect_modifier", prob = 0.4)

# Variable that both confounds and moderates the effect
causalsim_covar("normal", role = c("confounder", "effect_modifier"))

Create a causal data generating process

Description

Defines a causal DGP with known ground truth. Covariates can be specified via shorthand count arguments (Option B), an explicit named list of causalsim_covar() objects (Option A), or both combined.

Usage

causalsim_dgp(
  n,
  effect = 1,
  propensity = "moderate",
  baseline = "moderate",
  sigma = 1,
  covariates = list(),
  n_confounders = 0L,
  n_effect_modifiers = 0L,
  n_noise = 0L,
  mc_draws = 10000L
)

Arguments

Details

Structural model

W ~ covariate_spec
A ~ Bernoulli(propensity(W))
Y  = baseline(W) + effect(W) * A + N(0, sigma^2)

A is used for treatment throughout to avoid collision with R's built-in T alias.

Function calling convention

The effect, propensity, and baseline functions are called with named arguments matching covariate names, not a data frame. Write:

effect = function(W) 2 + 1.5 * W
propensity = function(W1, W2) plogis(0.3 * W1 + 0.5 * W2)

Every argument name is validated against the DGP's covariate spec at construction time, so mismatches surface immediately rather than at draw time. The propensity function is additionally evaluated on a small test draw to confirm it returns values in [0, 1].

True ATE

For scalar effect the true ATE is exact. For function effect it is approximated via Monte Carlo over mc_draws draws from the covariate distribution.

Value

An S3 object of class causalsim_dgp with components:

n

Sample size (integer)

covar_spec

Named list of causalsim_covar() objects

effect_fn

Normalized effect function

propensity_fn

Normalized propensity function

baseline_fn

Normalized baseline function

sigma

Outcome noise standard deviation

true_ate

True ATE: exact for scalar, MC approximation otherwise

heterogeneous

Logical; TRUE if effect was a function

mc_draws

Monte Carlo draws used (integer)

Examples

# Minimal: one confounder, constant effect, moderate confounding
dgp <- causalsim_dgp(n = 500, n_confounders = 1, effect = 2)
dgp

# Heterogeneous effect, explicit covariate spec
dgp2 <- causalsim_dgp(
  n = 500,
  covariates = list(
    W = causalsim_covar("normal", role = "confounder"),
    V = causalsim_covar("binary", role = "effect_modifier", prob = 0.4)
  ),
  effect = function(V) 2 + 1.5 * V,
  propensity = function(W) plogis(0.5 * W),
  baseline = function(W) 1.5 * W
)

# Mixed: shorthand confounders + explicit noise covariate + RCT propensity
dgp3 <- causalsim_dgp(
  n = 1000,
  n_confounders = 2,
  covariates = list(X = causalsim_covar("normal", role = "noise")),
  effect = 1,
  propensity = 0.5
)

Draw a dataset from a causal DGP

Description

Simulates one dataset from a causalsim_dgp() object using the structural model:

Usage

causalsim_draw(dgp, seed = NULL)

Arguments

Details

W ~ covariate_spec
A ~ Bernoulli(propensity(W))
Y  = baseline(W) + effect(W) * A + N(0, sigma^2)

Value

A data frame with dgp$n rows and the following columns:

Covariate columns

One column per covariate in dgp$covar_spec, named to match the spec (e.g. W, W1, Z).

A

Binary treatment indicator (0/1).

Y

Observed outcome.

.tau

Individual treatment effect (CATE). Ground truth.

.p

Individual propensity score. Ground truth.

The .tau and .p columns carry individual-level ground truth and are prefixed with . to distinguish them from observed variables. causalsim_eval() uses them directly without re-estimation.

Examples

dgp <- causalsim_dgp(n = 500, n_confounders = 1, effect = 2)
d   <- causalsim_draw(dgp, seed = 1)
head(d)

# Reproducible draws
d1 <- causalsim_draw(dgp, seed = 42)
d2 <- causalsim_draw(dgp, seed = 42)
identical(d1, d2)  # TRUE

Evaluate a causal estimator against a known DGP

Description

Repeatedly draws datasets from a causalsim_dgp() object, applies a user-supplied estimator to each draw, and returns a tidy summary of estimator performance against the known ground truth.

Usage

causalsim_eval(
  dgp,
  estimator,
  reps = 200L,
  metrics = c("bias", "rmse", "coverage", "power"),
  seed = NULL
)

Arguments

Details

Estimator convention

The estimator receives the full data frame returned by causalsim_draw(), including ground-truth columns .tau and .p. It should return a named numeric vector or one-row data frame with at minimum:

my_estimator <- function(data) {
  fit <- lm(Y ~ A + W, data = data)
  est <- coef(fit)["A"]
  se <- sqrt(vcov(fit)["A", "A"])
  c(estimate = est, ci_lower = est - 1.96 * se, ci_upper = est + 1.96 * se)
}

Metrics

"bias"

Monte Carlo estimate of bias: mean(estimate) - true_ate.

"rmse"

Root mean squared error: sqrt(mean((estimate - true_ate)^2)).

"coverage"

Proportion of replications where the CI brackets true_ate. Requires ci_lower and ci_upper.

"power"

Proportion of replications where the CI excludes zero, i.e., the rejection rate of H0: ATE = 0. When true_ate != 0, this is power; when true_ate == 0, it is the Type I error rate. Requires ci_lower and ci_upper.

Monte Carlo standard errors (MCSE) are reported alongside each metric. For bias: sd(estimates) / sqrt(reps). For RMSE: delta-method approximation. For coverage and power: Bernoulli SE.

Raw draws

Per-replication estimates are stored in result$draws and can be used for custom analysis or plotting.

Value

An S3 object of class causalsim_eval with components:

metrics

Tidy data frame: metric, value, se.

draws

Data frame of per-replication estimates, one row per rep.

true_ate

The DGP's true ATE.

reps

Number of replications run.

Examples

dgp <- causalsim_dgp(n = 300, n_confounders = 1, effect = 2,
                     propensity = 0.5)

# OLS estimator (unbiased under RCT)
ols_estimator <- function(data) {
  fit <- lm(Y ~ A + W, data = data)
  est <- coef(fit)[["A"]]
  se <- sqrt(vcov(fit)["A", "A"])
  c(estimate = est, ci_lower = est - 1.96 * se, ci_upper = est + 1.96 * se)
}


result <- causalsim_eval(dgp, ols_estimator, reps = 100L, seed = 1L)
result

Evaluate an estimator across a parameter grid

Description

Runs causalsim_eval() over the Cartesian product of one or more DGP parameter values, returning a tidy data frame of performance metrics for every cell. Designed for studying how estimator behavior changes with sample size, confounding strength, effect size, or noise level.

Usage

causalsim_eval_grid(
  dgp,
  estimator,
  vary,
  reps = 200L,
  metrics = c("bias", "rmse", "coverage", "power"),
  seed = NULL,
  verbose = FALSE
)

Arguments

Details

How it works

For each cell in expand.grid(vary), causalsim_eval_grid() takes dgp's stored original parameters, overrides the cell's values, reconstructs a new causalsim_dgp(), runs causalsim_eval(), and tags the resulting metrics with the cell's parameter values.

Constraints on vary

Elements of vary must be atomic vectors (character, numeric, integer). Functions cannot be varied via this interface. covariates (a list of causalsim_covar() objects) is excluded. For complex covariate variation, construct DGPs manually and use causalsim_eval() directly.

Note that varying n_confounders from 1 to 2 changes auto-generated covariate names from W to ⁠W1, W2⁠. If the base DGP's effect function references W, the reconstructed DGP will error at construction time. Design the base DGP accordingly.

Value

An S3 object of class causalsim_eval_grid with components:

results

Tidy data frame: one row per (cell, metric). Columns: grid parameter names, metric, value, se.

grid

Data frame of grid cells, one row per cell.

vary

Character vector of varied parameter names.

reps

Replications per cell.

metrics

Metrics evaluated.

Examples

dgp <- causalsim_dgp(n = 500, n_confounders = 1, effect = 2,
                     propensity = 0.5)

ols_estimator <- function(data) {
  fit <- lm(Y ~ A + W, data = data)
  est <- coef(fit)[["A"]]
  se  <- sqrt(vcov(fit)["A", "A"])
  c(estimate = est, ci_lower = est - 1.96 * se, ci_upper = est + 1.96 * se)
}


grid_result <- causalsim_eval_grid(
  dgp = dgp,
  estimator = ols_estimator,
  vary = list(n = c(250L, 500L, 1000L)),
  reps = 50L,
  seed = 1L
)
grid_result

Plot the distribution of estimates from a causalsim_eval result

Description

Draws a histogram of per-replication estimates with vertical lines marking the true ATE (solid) and the mean estimate (dashed).

Usage

## S3 method for class 'causalsim_eval'
plot(x, ...)

Arguments

Value

x, invisibly.

Summarise a causalsim_eval result

Description

Returns a structured summary of estimator performance including the metrics table and the distribution of per-replication estimates (mean, SD, median, and 10th / 90th percentiles).

Usage

## S3 method for class 'causalsim_eval'
summary(object, ...)

Arguments

Value

An S3 object of class causalsim_eval_summary with components:

metrics

The metrics data frame from the original result.

mean_estimate

Mean of per-replication estimates.

sd_estimate

Standard deviation of per-replication estimates.

median_estimate

Median of per-replication estimates.

p10_estimate

10th percentile of per-replication estimates.

p90_estimate

90th percentile of per-replication estimates.

true_ate

True ATE from the DGP.

reps

Number of replications.