Package 'medsim'

Description

Thin, documented amputer aligning mice::ampute semantics to the medsim DGM contract. Returns the input frame with NAs inserted in target column(s); column names/order are preserved (DGM data contract).

Usage

medsim_amputate(
  data,
  target,
  mechanism = c("MCAR", "MAR", "MNAR"),
  prop = 0.2,
  predictors = NULL,
  weights = NULL,
  type = c("RIGHT", "LEFT", "MID", "TAIL")
)

Arguments

Value

data with NAs inserted in target; names/order preserved.

Description

Computes comprehensive accuracy metrics for simulation results including bias, mean absolute error (MAE), root mean squared error (RMSE), and relative efficiency.

Usage

medsim_analyze(results, metrics = "all", by_scenario = TRUE)

Arguments

Details

Metrics Computed

For each estimated parameter with known ground truth:

• Bias: Mean estimation error (estimate - truth)

• MAE: Mean absolute error, robust to outliers

• RMSE: Root mean squared error, penalizes large errors

• Median AE: Median absolute error, very robust

• Max AE: Maximum absolute error, identifies worst case

• Relative Bias: Bias as percentage of true value

Interpretation

• Bias close to 0: unbiased estimator

• MAE/RMSE small: accurate estimator

• RMSE > MAE: presence of outliers/large errors

• Relative bias: standardized bias metric

Value

A list with class "medsim_analysis" containing:

• accuracy: data.frame with accuracy metrics

• by_scenario: data.frame with metrics by scenario (if by_scenario=TRUE)

• summary: Overall summary statistics

See Also

medsim_analyze_coverage(), medsim_compare_methods()

Examples

## Not run: 
# Run simulation
results <- medsim_run(method, scenarios, config)

# Analyze accuracy
analysis <- medsim_analyze(results)

# View overall metrics
print(analysis$accuracy)

# View by-scenario metrics
print(analysis$by_scenario)

# Custom metrics
analysis <- medsim_analyze(
  results,
  metrics = c("bias", "mae", "rmse"),
  by_scenario = FALSE
)

## End(Not run)

Description

Computes coverage rates for confidence intervals from simulation results. Requires that the simulation method returns confidence interval bounds (e.g., ci_lower, ci_upper).

Usage

medsim_analyze_coverage(
  results,
  ci_levels = c(0.9, 0.95, 0.99),
  ci_suffix = "_ci",
  by_scenario = TRUE
)

Arguments

Details

Coverage Interpretation

Coverage rate is the proportion of confidence intervals that contain the true parameter value. For a 95% CI, expect ~95% coverage in large samples.

• Coverage < nominal: CI too narrow (anti-conservative)

• Coverage ≈ nominal: CI has correct width

• Coverage > nominal: CI too wide (conservative)

Column Naming Conventions

The function looks for columns named:

• ⁠{parameter}_ci_lower⁠ and ⁠{parameter}_ci_upper⁠ (default)

• Or custom suffix: ⁠{parameter}{ci_suffix}_lower⁠ and ⁠{parameter}{ci_suffix}_upper⁠

Value

A list with class "medsim_coverage" containing:

• coverage: data.frame with coverage rates

• by_scenario: data.frame with coverage by scenario (if by_scenario=TRUE)

• summary: Overall summary statistics

See Also

medsim_analyze(), medsim_analyze_power()

Examples

## Not run: 
# Method that returns CIs
my_method <- function(data, params) {
  # ... estimation
  list(
    indirect = estimate,
    indirect_ci_lower = ci[1],
    indirect_ci_upper = ci[2]
  )
}

results <- medsim_run(my_method, scenarios, config)

# Analyze coverage
coverage <- medsim_analyze_coverage(results)

print(coverage$coverage)
#>   parameter  coverage_90  coverage_95  coverage_99
#>   indirect   0.902        0.951        0.991

## End(Not run)

Description

Computes empirical power (rejection rate under alternative hypothesis) from simulation results. Requires that the simulation method returns p-values.

Usage

medsim_analyze_power(
  results,
  alpha = 0.05,
  p_suffix = "_p",
  null_value = 0,
  by_scenario = TRUE
)

Arguments

Details

Power Interpretation

Power is the probability of correctly rejecting the null hypothesis when it is false (i.e., detecting a true effect).

• Power = proportion of replications where p < alpha

• Higher power = better ability to detect effects

• Power depends on: effect size, sample size, alpha level

Notes

This function computes empirical power from simulations, not theoretical power. For scenarios where the null is true (e.g., indirect effect = 0), the rejection rate represents Type I error rate, not power.

Value

A list with class "medsim_power" containing:

• power: data.frame with power rates

• by_scenario: data.frame with power by scenario (if by_scenario=TRUE)

• summary: Overall summary statistics

See Also

medsim_analyze(), medsim_analyze_coverage()

Examples

## Not run: 
# Method that returns p-values
my_method <- function(data, params) {
  # ... estimation and testing
  list(
    indirect = estimate,
    indirect_p = p_value
  )
}

results <- medsim_run(my_method, scenarios, config)

# Analyze power
power <- medsim_analyze_power(results, alpha = 0.05)

print(power$power)
#>   parameter  power    n_valid
#>   indirect   0.847    1000

## End(Not run)

Description

Deletes cache files from a directory. Can delete all files or files matching a pattern.

Usage

medsim_cache_clear(
  cache_dir,
  pattern = "*.rds",
  max_age = NULL,
  confirm = TRUE
)

Arguments

Value

Invisibly returns number of files deleted

Examples

## Not run: 
# Clear all cache
medsim_cache_clear("cache")

# Clear old cache only
medsim_cache_clear("cache", max_age = 30)

# Clear specific pattern
medsim_cache_clear("cache", pattern = "truth_*.rds")

# Without confirmation
medsim_cache_clear("cache", confirm = FALSE)

## End(Not run)

Description

Checks if a cache file exists and is readable.

Usage

medsim_cache_exists(file)

Arguments

Value

Logical: TRUE if cache exists and is readable, FALSE otherwise

Examples

## Not run: 
if (medsim_cache_exists("cache/truth.rds")) {
  truth <- medsim_cache_load("cache/truth.rds")
} else {
  truth <- list(indirect = 0.09)  # Compute your truth values
  medsim_cache_save(truth, "cache/truth.rds")
}

## End(Not run)

Description

Retrieves metadata about a cached object without loading the full object.

Usage

medsim_cache_info(file)

Arguments

Value

List with cache metadata:

• exists: Logical - file exists

• size_mb: Numeric - file size in MB

• modified: POSIXct - last modified time

• age_days: Numeric - age in days

• timestamp: POSIXct - creation time (if available)

• r_version: Character - R version used (if available)

Examples

# Save a cached object to a temp file, then inspect it
cache_file <- tempfile(fileext = ".rds")
medsim_cache_save(list(value = 42), cache_file)
info <- medsim_cache_info(cache_file)
print(info)
unlink(cache_file)

Description

Creates a cache directory with proper structure and a README file explaining the cache system.

Usage

medsim_cache_init(cache_dir, create_readme = TRUE)

Arguments

Value

Invisibly returns cache directory path

Examples

cache_dir <- file.path(tempdir(), "medsim_cache_example")
medsim_cache_init(cache_dir)
unlink(cache_dir, recursive = TRUE)

Description

Lists all cache files in a directory with their metadata.

Usage

medsim_cache_list(cache_dir, pattern = "*.rds")

Arguments

Value

data.frame with columns:

• file: File name

• path: Full path

• size_mb: Size in MB

• age_days: Age in days

• modified: Last modified time

Examples

# Create a small cache in a temp directory, then list it
cache_dir <- file.path(tempdir(), "medsim_cache_list_example")
dir.create(cache_dir, showWarnings = FALSE)
medsim_cache_save(list(a = 1), file.path(cache_dir, "result1.rds"))
medsim_cache_save(list(b = 2), file.path(cache_dir, "result2.rds"))

cache_list <- medsim_cache_list(cache_dir)
print(cache_list)

unlink(cache_dir, recursive = TRUE)

Description

Loads a previously cached R object. Returns NULL if cache file doesn't exist.

Usage

medsim_cache_load(file, check_expiry = FALSE, max_age = 30, verbose = TRUE)

Arguments

Details

Cache Validation

When loading, the cache is checked for:

• File exists

• File is readable

• Age is within max_age (if check_expiry = TRUE)

Handling Missing Cache

If cache doesn't exist:

• Returns NULL (no error)

• Calling code should compute and cache the result

Example pattern:

result <- medsim_cache_load("cache/truth.rds")
if (is.null(result)) {
  result <- expensive_computation()
  medsim_cache_save(result, "cache/truth.rds")
}

Value

Cached object if file exists and is valid, NULL otherwise

See Also

medsim_cache_save(), medsim_cache_exists()

Examples

## Not run: 
# Load from cache
truth <- medsim_cache_load("cache/truth_scenario1.rds")

if (is.null(truth)) {
  # Cache miss - compute and cache
  truth <- compute_expensive_truth()
  medsim_cache_save(truth, "cache/truth_scenario1.rds")
}

# Load with expiry check
result <- medsim_cache_load(
  "cache/old_result.rds",
  check_expiry = TRUE,
  max_age = 7  # 7 days
)

## End(Not run)

Description

Saves an R object to a cache file for later retrieval. Uses RDS format for efficient storage and compression.

Usage

medsim_cache_save(object, file, compress = TRUE, overwrite = TRUE)

Arguments

Details

Cache File Format

Cache files are stored in RDS format with metadata:

• Original object

• Timestamp of creation

• R version used

• Package version (if available)

File Organization

Recommended cache directory structure:

cache/
|-- truth_scenario_1.rds
|-- truth_scenario_2.rds
+-- ...

Value

Invisibly returns the file path

See Also

medsim_cache_load(), medsim_cache_clear()

Examples

## Not run: 
# Cache simulation truth
truth <- compute_expensive_truth()
medsim_cache_save(truth, "cache/truth_scenario1.rds")

# Cache with metadata
result <- list(value = 42, timestamp = Sys.time())
medsim_cache_save(result, "cache/my_result.rds")

## End(Not run)

Description

Checks a list of results from parallel execution for errors and provides a summary.

Usage

medsim_check_results(results, stop_on_error = FALSE)

Arguments

Value

Invisibly returns indices of failed tasks

Examples

## Not run: 
results <- medsim_run_parallel(...)
medsim_check_results(results)

## End(Not run)

Description

Creates a comparison table showing differences between test, local, and cluster configurations.

Usage

medsim_compare_configs()

Value

A data.frame with configuration comparison

Examples

comparison <- medsim_compare_configs()
print(comparison)

Description

Compares performance of multiple estimation methods across scenarios. Combines results from multiple simulation runs and computes relative performance metrics.

Usage

medsim_compare_methods(..., metrics = "all")

Arguments

Details

Comparison Metrics

• Accuracy: MAE, RMSE, bias across methods

• Timing: Mean/median execution time

• Coverage: CI coverage rates (if CIs available)

• Relative Efficiency: Ratio of MSE (lower = better)

Notes

All results must:

• Use the same scenarios (matched by name)

• Use the same configuration (n_replications, seed)

• Have ground truth available (for accuracy comparison)

Value

A list with class "medsim_comparison" containing:

• accuracy_comparison: data.frame comparing accuracy metrics

• timing_comparison: data.frame comparing timing

• coverage_comparison: data.frame comparing coverage (if applicable)

• summary: Overall comparison summary

See Also

medsim_analyze(), medsim_analyze_coverage()

Examples

## Not run: 
# Run multiple methods
results_proposed <- medsim_run(method_proposed, scenarios, config)
results_delta <- medsim_run(method_delta, scenarios, config)
results_boot <- medsim_run(method_boot, scenarios, config)

# Compare
comparison <- medsim_compare_methods(
  proposed = results_proposed,
  delta = results_delta,
  bootstrap = results_boot
)

# View accuracy comparison
print(comparison$accuracy_comparison)

# View timing comparison
print(comparison$timing_comparison)

## End(Not run)

Description

Creates a configuration object for simulation studies. Automatically detects whether running on local machine or HPC cluster and adjusts parameters accordingly.

Usage

medsim_config(
  mode = "auto",
  n_replications = NULL,
  n_cores = NULL,
  scenarios = NULL,
  output_dir = NULL,
  seed = 12345,
  ...
)

Arguments

Details

Execution Modes

Environment Detection

When mode = "auto", checks for:

• SLURM_JOB_ID (SLURM scheduler)

• PBS_JOBID (PBS/Torque scheduler)

• LSB_JOBID (LSF scheduler)

If any are found, uses "cluster" mode. Otherwise, uses "local" mode.

Custom Parameters

You can add custom parameters via ...:

config <- medsim_config(
  mode = "local",
  n_bootstrap = 5000,
  alpha = 0.05,
  custom_param = "value"
)

Value

A list with simulation configuration parameters

See Also

medsim_detect_environment(), print.medsim_config()

Examples

# Auto-detect environment
config <- medsim_config(mode = "auto")

# Explicit test mode for quick validation
config_test <- medsim_config(mode = "test")

# Local mode with custom replications
config_local <- medsim_config(
  mode = "local",
  n_replications = 500
)

# Cluster mode (auto-detects SLURM cores)
config_cluster <- medsim_config(mode = "cluster")

# Custom parameters
config_custom <- medsim_config(
  mode = "local",
  n_bootstrap = 1000,
  ci_level = 0.95
)

# Print configuration
print(config)

Description

Detects the number of CPU cores available for parallel processing. On HPC clusters, uses scheduler-allocated cores. On local machines, uses all cores minus 2 to avoid overwhelming the system.

Usage

medsim_detect_cores()

Details

On clusters, checks for:

• SLURM_CPUS_PER_TASK (SLURM)

• PBS_NUM_PPN (PBS/Torque)

• LSB_DJOB_NUMPROC (LSF)

On local machines, uses parallel::detectCores() - 2.

Value

Integer: Number of cores to use

Examples

n_cores <- medsim_detect_cores()
message(sprintf("Using %d cores", n_cores))

Description

Detects whether code is running on local machine or HPC cluster by checking for scheduler environment variables.

Usage

medsim_detect_environment()

Details

Checks for the following environment variables:

• SLURM_JOB_ID (SLURM scheduler)

• PBS_JOBID (PBS/Torque scheduler)

• LSB_JOBID (LSF scheduler)

If any are found, returns "cluster". Otherwise, returns "local".

Value

Character: "local" or "cluster"

Examples

env <- medsim_detect_environment()
if (env == "cluster") {
  message("Running on HPC cluster")
} else {
  message("Running on local machine")
}

Description

Estimates the expected speedup from parallel processing based on number of cores and task characteristics.

Usage

medsim_estimate_speedup(n_tasks, n_cores, overhead = 0.1)

Arguments

Details

Speedup is estimated using Amdahl's law with overhead:

Speedup = n_cores / (1 + overhead * n_cores)

Efficiency = Speedup / n_cores

Value

Named list with:

• speedup: Expected speedup factor

• efficiency: Parallel efficiency (0-1)

• recommendation: Character recommendation

Examples

# For 100 tasks on 10 cores
est <- medsim_estimate_speedup(100, 10)
cat(sprintf("Expected speedup: %.1fx\n", est$speedup))
cat(sprintf("Efficiency: %.0f%%\n", est$efficiency * 100))

Description

Convenience wrapper that produces the standard set of medsim plots (coverage, error boxplot, timing) and writes them to output_dir.

Usage

medsim_figures(
  results,
  output_dir = "figures",
  format = "pdf",
  width = 8,
  height = 6,
  ...
)

Arguments

Value

Invisibly, a named character vector of file paths written.

Description

Determines the optimal number of cores to use for parallel processing based on system capabilities and task requirements.

Usage

medsim_get_optimal_cores(n_tasks, reserve = 2, max_cores = NULL)

Arguments

Value

Integer: Optimal number of cores to use

Examples

# For 100 tasks
medsim_get_optimal_cores(n_tasks = 100)

# For few tasks (returns n_tasks)
medsim_get_optimal_cores(n_tasks = 3)

# With maximum limit
medsim_get_optimal_cores(n_tasks = 1000, max_cores = 8)

Description

Inverse-probability-weighting comparator. Positioned as a robustness appendix, not a co-equal main-results arm (decided 2026-06-11).

Usage

medsim_method_ipw(model, ...)

Arguments

Value

A ⁠function(data, params)⁠ returning the shared method contract (branch_switch = NA).

Description

Multiple imputation (via mice) + D4-stacked MBCO likelihood-ratio inference for the union null H0: ab = 0. Ports the validated ⁠Missing Effect/code/prototype-d4-mbco.R⁠ (phase-2, exact match vs mitml::testModels(method = "D4")): per-imputation MBCO LRTs are pooled with the Reiter/Chan-Meng D4 rule (Grund et al. 2021) into an F reference.

Usage

medsim_method_mbco_mi(model, m = 20L, ...)

Arguments

Details

The point estimate is the Rubin-pooled a*b; the CI is a Monte-Carlo product interval on the pooled paths; the p-value is the D4-pooled MBCO test. With no missingness (or without mice) it degrades to the complete-case MBCO LRT with a chi-square reference — never to a Sobel normal approximation for the test.

Value

A ⁠function(data, params)⁠ returning the shared method contract (see file header).

Description

Multiple imputation (via mice) + Monte-Carlo confidence interval for the indirect effect (via RMediation::medci() when available, else a base-R product-of-normals draw on the Rubin-pooled paths).

Usage

medsim_method_mc_ci(model, m = 20L, ...)

Arguments

Value

A ⁠function(data, params)⁠ returning the shared method contract (branch_switch = NA).

Description

Creates a combined figure with multiple panels for manuscript submission. Useful for showing multiple aspects of simulation results in one figure.

Usage

medsim_plot_combined_panel(
  results,
  panels = c("error", "timing"),
  layout = NULL,
  ...
)

Arguments

Details

Requires either patchwork or gridExtra package for combining plots.

Value

A combined plot object (patchwork or gridExtra)

Examples

## Not run: 
# Two-panel figure
p <- medsim_plot_combined_panel(
  results,
  panels = c("error", "timing"),
  layout = "1x2"
)

# Four-panel figure for manuscript
p <- medsim_plot_combined_panel(
  list(Proposed = res1, Delta = res2),
  panels = c("error", "timing", "coverage"),
  layout = "2x2"
)

## End(Not run)

Description

Creates a plot showing confidence interval coverage rates across scenarios and methods.

Usage

medsim_plot_coverage(
  coverage,
  expected_coverage = 0.95,
  tolerance = 0.02,
  title = NULL,
  ...
)

Arguments

Value

A ggplot2 object

Examples

## Not run: 
coverage <- medsim_analyze_coverage(results)
p <- medsim_plot_coverage(coverage, expected_coverage = 0.95)

## End(Not run)

Description

Creates boxplots showing the distribution of estimation errors across scenarios. Useful for comparing accuracy of different methods.

Usage

medsim_plot_error_boxplot(
  results,
  parameter = "indirect",
  log_scale = FALSE,
  color_palette = NULL,
  title = NULL,
  ...
)

Arguments

Details

Interpreting Boxplots

• Box: IQR (25th to 75th percentile)

• Line in box: Median error

• Whiskers: 1.5 × IQR

• Points: Outliers beyond whiskers

Tighter boxes = more consistent estimates

Usage Patterns

# Single method
p <- medsim_plot_error_boxplot(results)

# Compare methods
p <- medsim_plot_error_boxplot(list(
  Proposed = results1,
  Delta = results2,
  Bootstrap = results3
))

# Customize
p <- medsim_plot_error_boxplot(
  results,
  log_scale = TRUE,
  color_palette = "Set2"
)
p + ggplot2::theme_minimal()

Value

A ggplot2 object

Examples

## Not run: 
results <- medsim_run(method, scenarios, config)
p <- medsim_plot_error_boxplot(results)
print(p)

# Save to file
ggplot2::ggsave("error_boxplot.pdf", p, width = 10, height = 6)

## End(Not run)

Description

Creates bar plots comparing computational time across methods and scenarios.

Usage

medsim_plot_timing(
  results,
  metric = "mean",
  log_scale = FALSE,
  title = NULL,
  ...
)

Arguments

Value

A ggplot2 object

Examples

## Not run: 
p <- medsim_plot_timing(list(
  Proposed = results1,
  Delta = results2,
  Bootstrap = results3
))

## End(Not run)

Description

Fleishman / Vale–Maurelli power-method draws. Pure base R; no hard dependency. Warns/errors on (skew, kurtosis) outside the Fleishman-feasible region.

Usage

medsim_rnonnormal(n, mean = 0, sd = 1, skew = 0, kurtosis = 0)

Arguments

Value

Numeric vector of length n with the requested marginal moments.

Description

Executes a complete simulation study by running a user-defined method across multiple scenarios and replications. Supports parallel processing and automatic result aggregation.

Usage

medsim_run(
  method,
  scenarios,
  config,
  compute_truth = NULL,
  parallel = TRUE,
  verbose = TRUE
)

Arguments

Details

Method Function Requirements

The method function must:

• Accept data (data.frame) as first argument

• Accept params (list) as second argument

• Return a named list with at least one numeric element

Example:

my_method <- function(data, params) {
  fit_m <- lm(M ~ X, data = data)
  fit_y <- lm(Y ~ X + M, data = data)

  a <- coef(fit_m)["X"]
  b <- coef(fit_y)["M"]

  list(
    indirect = a * b,
    a_path = a,
    b_path = b,
    se_indirect = ...,  # Optional
    ci_lower = ...,     # Optional
    ci_upper = ...      # Optional
  )
}

Parallel Processing

When parallel = TRUE:

• Uses medsim_run_parallel() internally

• Number of cores from config$n_cores

• Progress bars shown (local) or suppressed (cluster)

• All necessary objects exported to workers

Ground Truth

If compute_truth is provided:

• Computed once per scenario (cached automatically)

• Used to calculate errors/bias in results

• Can use Monte Carlo or analytical methods

Value

A medsim_results object (list) containing:

• results: data.frame with all simulation results

• summary: data.frame with summary statistics

• truth: data.frame with ground truth values (if compute_truth provided)

• config: configuration used

• scenarios: scenarios used

• method_name: name of method function

• timestamp: when simulation was run

See Also

medsim_config(), medsim_scenario(), medsim_run_parallel()

Examples

## Not run: 
# Define method
my_method <- function(data, params) {
  fit_m <- lm(M ~ X, data = data)
  fit_y <- lm(Y ~ X + M, data = data)
  list(indirect = coef(fit_m)["X"] * coef(fit_y)["M"])
}

# Configure and run
config <- medsim_config("local")
scenarios <- medsim_scenarios_mediation()

results <- medsim_run(
  method = my_method,
  scenarios = scenarios,
  config = config
)

# View results
print(results)
summary(results)

## End(Not run)

Description

Executes a list of tasks in parallel using multiple CPU cores. Provides progress bars, proper error handling, and automatic cluster management.

Usage

medsim_run_parallel(
  tasks,
  fun,
  n_cores = NULL,
  progress = TRUE,
  export = NULL,
  packages = NULL,
  cluster_type = NULL
)

Arguments

Details

Cluster Types

• PSOCK (default on Windows): Creates separate R sessions. Requires explicit export of objects and packages. Works on all platforms.

• FORK (default on Unix): Copies current R session. More efficient, automatic object sharing, but Unix-only.

Progress Bars

Progress bars are shown in interactive sessions and local execution, but suppressed on HPC clusters to avoid cluttering log files.

Error Handling

If a task fails:

• Error is caught and returned as an error object

• Other tasks continue processing

• Check results with sapply(results, inherits, "error")

Memory Considerations

Each worker process requires memory. For large tasks:

• Reduce n_cores if running out of memory

• Process tasks in chunks

• Use FORK cluster (Unix) which shares memory

Value

List: Results from applying fun to each task (same length as tasks)

See Also

parallel::makeCluster(), parallel::parLapply(), pbapply::pblapply()

Examples

## Not run: 
# Simple parallel computation
results <- medsim_run_parallel(
  tasks = 1:100,
  fun = function(i) {
    Sys.sleep(0.1)  # Simulate work
    i^2
  },
  n_cores = 4
)

# With exports
my_data <- data.frame(x = 1:10, y = 11:20)

results <- medsim_run_parallel(
  tasks = 1:10,
  fun = function(i) {
    mean(my_data$x[1:i])
  },
  export = "my_data"
)

# With packages
results <- medsim_run_parallel(
  tasks = 1:10,
  fun = function(i) {
    dplyr::tibble(x = i, y = i^2)
  },
  packages = "dplyr"
)

# Check for errors
errors <- sapply(results, inherits, "error")
if (any(errors)) {
  cat("Errors occurred in:", which(errors), "\n")
}

## End(Not run)

Description

Define a custom scenario for simulation studies. A scenario consists of a name, description, data generation function, and population parameters.

Usage

medsim_scenario(name, description = "", data_generator, params = list())

Arguments

Details

Data Generator Function

The data_generator function must:

• Accept n (sample size) as first argument

• Return a data.frame with at least: X (treatment), M (mediator), Y (outcome)

• Can include additional variables or covariates

Parameters List

The params list should include known population values:

• a: X -> M path coefficient

• b: M -> Y path coefficient

• c_prime: X -> Y direct effect

• indirect: True indirect effect (a * b)

• Additional parameters as needed

Value

A scenario object (list with class "medsim_scenario")

See Also

medsim_scenarios_mediation() for standard scenarios

Examples

# Simple custom scenario
my_scenario <- medsim_scenario(
  name = "Large Effects",
  description = "Both paths have large effects",
  data_generator = function(n = 200) {
    X <- rnorm(n)
    M <- 0.7 * X + rnorm(n)
    Y <- 0.7 * M + rnorm(n)
    data.frame(X = X, M = M, Y = Y)
  },
  params = list(
    a = 0.7,
    b = 0.7,
    indirect = 0.49
  )
)

# Generate data
data <- my_scenario$data_generator(n = 100)

# Scenario with covariates
covariate_scenario <- medsim_scenario(
  name = "With Covariates",
  data_generator = function(n = 200) {
    C1 <- rnorm(n)
    C2 <- rbinom(n, 1, 0.5)
    X <- rnorm(n)
    M <- 0.3 * X + 0.2 * C1 + rnorm(n)
    Y <- 0.3 * M + 0.2 * C2 + rnorm(n)
    data.frame(X = X, M = M, Y = Y, C1 = C1, C2 = C2)
  },
  params = list(
    a = 0.3,
    b = 0.3,
    indirect = 0.09,
    gamma_m = 0.2,
    gamma_y = 0.2
  )
)

Description

Wraps medsim_scenario() with a data_generator that (1) draws complete X -> M -> Y with optional nonnormal residuals (medsim_rnonnormal()), then (2) amputes target under mechanism (medsim_amputate()).

Usage

medsim_scenario_missing(
  name,
  true_params,
  mechanism,
  prop = 0.2,
  target = "M",
  nonnormal = NULL
)

Arguments

Value

A medsim_scenario object.

Description

Convenience: expand SPEC-simulation-design cells (mechanism × prop × n × effect × nonnormality) into a list of medsim_scenario_missing() objects in one call.

Usage

medsim_scenario_missing_grid(
  true_params_list,
  mechanisms,
  props,
  nonnormal_list = list(NULL)
)

Arguments

Value

A list of medsim_scenario objects.

Description

Creates a list of standard mediation scenarios for simulation studies. These scenarios cover common patterns in mediation analysis including independent paths, various correlation structures, suppression effects, and non-standard conditions.

Usage

medsim_scenarios_mediation()

Details

Standard Scenarios

1. Independent Paths: No correlation between X, M, Y

2. Moderate Correlation: rho = 0.3 between all pairs

3. High Correlation: rho = 0.7 between all pairs

4. Suppression: Mixed positive and negative correlations

5. Non-zero Effects: Small to moderate true effects

6. Unequal Variances: Different residual variances

Each scenario generates data with:

• Sample size n (default: 200)

• Treatment X, Mediator M, Outcome Y

• Known population parameters for validation

Value

A list of scenario objects, each with name, description, and data_generator function

See Also

medsim_scenario() for creating custom scenarios

Examples

# Get all standard scenarios
scenarios <- medsim_scenarios_mediation()

# See scenario names
sapply(scenarios, function(s) s$name)

# Generate data from first scenario
data <- scenarios[[1]]$data_generator(n = 100)
head(data)

# Access scenario parameters
scenarios[[1]]$params

Description

Reads the branch_switch column emitted by medsim_method_mbco_mi() (the union-null / ab=0 diagnostic from MEMO-MBCO-MI-derivation) and returns a per-scenario summary suitable for joining as a column onto a medsim_table_coverage() table. Non-converged rows (converged == 0) are dropped before summarizing.

Usage

medsim_summarize_branch_switch(results, by = "scenario")

Arguments

Value

A data.frame: one row per group with branch_switch_rate and n_valid.

Description

Creates a publication-ready LaTeX table showing accuracy metrics (bias, MAE, RMSE) from simulation analysis results.

Usage

medsim_table_accuracy(
  analysis,
  digits = 3,
  metrics = c("bias", "mae", "rmse"),
  by_scenario = TRUE,
  caption = NULL,
  label = "tab:accuracy"
)

Arguments

Value

A medsim_table object (character vector with LaTeX code)

Examples

## Not run: 
results <- medsim_run(method, scenarios, config)
analysis <- medsim_analyze(results)

# Generate table
table <- medsim_table_accuracy(analysis)
print(table)  # Preview

# Write to file
medsim_write_table(table, "tables/accuracy.tex")

## End(Not run)

Description

Creates a publication-ready LaTeX table comparing multiple methods across accuracy, timing, and coverage metrics.

Usage

medsim_table_comparison(
  comparison,
  metrics = c("accuracy", "timing"),
  caption = NULL,
  label = "tab:comparison"
)

Arguments

Value

A medsim_table object

Description

Creates a publication-ready LaTeX table showing confidence interval coverage rates.

Usage

medsim_table_coverage(
  coverage,
  expected = 0.95,
  by_scenario = TRUE,
  caption = NULL,
  label = "tab:coverage"
)

Arguments

Value

A medsim_table object

Description

Creates a publication-ready LaTeX table showing empirical power (rejection rates) from simulation analysis.

Usage

medsim_table_power(
  power,
  by_scenario = TRUE,
  caption = NULL,
  label = "tab:power"
)

Arguments

Value

A medsim_table object

Description

Creates a publication-ready LaTeX table comparing computational time across methods.

Usage

medsim_table_timing(
  results_list,
  metric = "mean",
  include_speedup = TRUE,
  reference_method = NULL,
  by_scenario = FALSE,
  caption = NULL,
  label = "tab:timing"
)

Arguments

Value

A medsim_table object

Examples

## Not run: 
table <- medsim_table_timing(list(
  Proposed = results1,
  Delta = results2,
  Bootstrap = results3
))

## End(Not run)

Description

Convenience wrapper around medsim_tables_workflow() with sensible defaults.

Usage

medsim_tables(results, output_dir = "tables", format = "latex", ...)

Arguments

Value

Invisibly, the result of medsim_tables_workflow().

Description

Convenience function to generate all tables from simulation results and save them to an output directory.

Usage

medsim_tables_workflow(
  results,
  output_dir,
  tables = "all",
  format = "latex",
  prefix = "table_"
)

Arguments

Value

Invisibly returns a list of generated file paths

Examples

## Not run: 
results <- medsim_run(method, scenarios, config)

# Generate all tables
files <- medsim_tables_workflow(
  results,
  output_dir = "manuscript/tables",
  format = "latex"
)

## End(Not run)

Description

Checks that a scenario object is valid and can generate data correctly.

Usage

medsim_validate_scenario(scenario, n = 10)

Arguments

Value

TRUE if valid, throws error otherwise

Examples

scenarios <- medsim_scenarios_mediation()
medsim_validate_scenario(scenarios[[1]])

Description

End-to-end output pipeline: runs medsim_analyze(), medsim_figures(), and medsim_tables() and writes everything to output_dir.

Usage

medsim_workflow(
  results,
  output_dir = "results",
  figures_format = "pdf",
  tables_format = "latex",
  ...
)

Arguments

Value

Invisibly, a list with elements analysis, figures (file paths), and tables (return value of medsim_tables_workflow()).

Examples

## Not run: 
config <- medsim_config("test")
scenarios <- medsim_scenarios_mediation()
results <- medsim_run(my_method, scenarios, config)
medsim_workflow(results, output_dir = "results")

## End(Not run)

Description

Writes a medsim_table object to a file in LaTeX, CSV, or Markdown format.

Usage

medsim_write_table(table, file, format = NULL)

Arguments

Value

Invisibly returns the file path

Examples

## Not run: 
table <- medsim_table_accuracy(analysis)

# Write LaTeX
medsim_write_table(table, "tables/accuracy.tex")

# Write CSV
medsim_write_table(table, "tables/accuracy.csv", format = "csv")

## End(Not run)

Description

Print Analysis Results

Usage

## S3 method for class 'medsim_analysis'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Print Method Comparison

Usage

## S3 method for class 'medsim_comparison'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Print Configuration Summary

Usage

## S3 method for class 'medsim_config'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Print Coverage Results

Usage

## S3 method for class 'medsim_coverage'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Print Power Results

Usage

## S3 method for class 'medsim_power'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Print Simulation Results

Usage

## S3 method for class 'medsim_results'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Print Scenario Summary

Usage

## S3 method for class 'medsim_scenario'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Print medsim_table

Usage

## S3 method for class 'medsim_table'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Summarize Simulation Results

Usage

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

Arguments

Value

Summary data.frame