---
title: "Cookbook: common tasks"
format: html
vignette: >
  %\VignetteIndexEntry{Cookbook}
  %\VignetteEngine{quarto::html}
  %\VignetteEncoding{UTF-8}
---

```{r}
#| label: setup
#| include: false
knitr::opts_chunk$set(collapse = TRUE, comment = "#>", eval = FALSE)
```

Short recipes for common tasks. Every recipe uses the simulated data below (continuous
treatment `X`, mediator `M`, outcome `Y`, covariate `C`). Chunks are not evaluated when the
site is built; the code was run against version 0.4.1.

```{r}
library(probmed)
set.seed(123)
n <- 300
d <- data.frame(X = rnorm(n), C = rnorm(n))
d$M <- 0.5 * d$X + 0.3 * d$C + rnorm(n)
d$Y <- 0.4 * d$M + 0.2 * d$X + 0.2 * d$C + rnorm(n)
```

## Point estimate only

```{r}
fit <- pmed(Y ~ X + M + C, formula_m = M ~ X + C, data = d,
            treatment = "X", mediator = "M", method = "plugin", seed = 1)
fit@estimate
```

`"plugin"` returns the point estimate with no interval. Use it for quick checks.

## Interval by method

All methods share one call; only `method` changes.

```{r}
args <- list(Y ~ X + M + C, formula_m = M ~ X + C, data = d,
             treatment = "X", mediator = "M", seed = 1)
pb   <- do.call(pmed, c(args, method = "parametric_bootstrap", n_boot = 1000))
npb  <- do.call(pmed, c(args, method = "nonparametric_bootstrap", n_boot = 1000))
mbco <- do.call(pmed, c(args, method = "mbco"))
```

See [Choosing a confidence interval](confidence-intervals.html) for when to use which.

## Read the result

```{r}
mbco@estimate; mbco@ci_lower; mbco@ci_upper   # P_med and its interval
mbco@ie_estimate                               # indirect effect a*b
mbco@diagnostics$pmed                          # interval status, near_double_null flag
```

For `method = "mbco"`, check that the printed interval status reads `ok`. A missing bound means
a verified likelihood-ratio endpoint was not found; the method reports that instead of a guess.

## Binary mediator or outcome with covariates

Use the model-based bootstrap, which refits both equations on simulated data.

```{r}
pmed(Y_bin ~ X + M + C, formula_m = M ~ X + C, data = d_bin,
     treatment = "X", mediator = "M",
     method = "model_bootstrap", n_boot = 500, seed = 1)
```

A rare outcome (about 5%) at small N often returns `bootstrap_incomplete` instead of an
interval. That is the method declining to report, not an error. Poisson, ordinal and other
families, interactions and transformed treatment or mediator terms are rejected.

## From a fitted lavaan or mediation object

`extract_mediation()` returns a `MediationData` object, and `pmed()` has a method for it. See
[lavaan integration](lavaan-integration.html) and [mediation integration](mediation-integration.html).

## Two mediators

For parallel mediators use the joint P_med; see [Parallel mediation](parallel-mediation.html).

## Tipping point for unmeasured M-Y confounding (mean-based shares)

```{r}
pmed_sensitivity(indirect = 0.30, total = 0.50, threshold = 0, evalue = TRUE)
```

More in [Sensitivity analysis](sensitivity-analysis.html).

## Weak identification check

```{r}
ward_residual(data.frame(A = as.integer(d$X > 0), M = d$M, Y = d$Y, C = d$C))
```

`ward_residual()` expects a binary 0/1 treatment `A` and covariates named in `covars` (default `"C"`). It reports gates for weak identification;
see [Gauge residual](gauge-residual.html).

## Missing values

`pmed()` stops if any variable named in either formula has a missing value, before fitting
anything. Drop incomplete rows first so both models use the same sample.

```{r}
d_complete <- d[stats::complete.cases(d[, c("X", "M", "Y", "C")]), ]
```

Columns that neither formula uses may contain missing values. Since 0.4.1 this check runs on
the data you pass, so it also holds with medfit 0.10 and later.

## Reproducibility

Pass `seed` to every call. Plug-in calls honor it, MBCO is deterministic, and bootstrap
methods draw from `seed`.