---
title: "Confounding and Interaction"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Confounding and Interaction}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

Confounding asks whether adjustment changes the exposure effect. Interaction
asks whether the exposure effect differs across another variable.

These checks are screening aids for viewing and organising results. Use DAGs,
subject-matter knowledge, and study design to decide which variables are
confounders or effect modifiers. Automated change-in-estimate and interaction
checks should not be used as the sole basis for model adjustment.

## Which Function Should I Use?

| Question | Use | Why |
|---|---|---|
| Could these candidate variables be confounders or effect modifiers? | `identify_confounder()` | Screens crude, adjusted, Mantel-Haenszel, and interaction signals together. |
| Does this planned interaction term improve the model? | `interaction_models()` | Compares models with and without `exposure:effect_modifier` using LRT or Wald tests. |
| Which automatically generated candidate model has the strongest fit? | `select_models()` | Fits and ranks a defined set of candidate models using model-fit statistics. |
| How do my named, fitted gtregression models differ? | `compare_models()` | Compares supplied fitted models, checks analysis-sample consistency, and reports fit statistics transparently. |
| Could part of an exposure-outcome association operate through a mediator? | `mediation_analysis()` | Estimates direct, indirect, total, and proportion mediated effects with explicit causal caveats. |
| Do I need a Mantel-Haenszel estimate? | `identify_confounder(method = "mh")` or `identify_confounder(method = "both")` | MH is a stratified pooled estimate for eligible binary/categorical settings, not a formal interaction-term test. |

A practical workflow is:

1. Use DAGs, prior literature, and study design to list important variables.
2. Use `identify_confounder()` to organise screening evidence for candidate
   confounders or effect modifiers.
3. Use `interaction_models()` when you have a planned interaction hypothesis.
4. Use `mediation_analysis()` when the mediator is part of a planned causal
   question and the temporal order is defensible.
5. If interaction is important, consider stratified reporting with
   `stratified_uni_reg()` or `stratified_multi_reg()`.

## How the Interaction Tests Differ

`identify_confounder()` is a candidate-by-candidate screening aid. For each
exposure and potential confounder, it compares `outcome ~ exposure` with
`outcome ~ exposure + potential_confounder`, examines stratum-specific
estimates, and screens the candidate interaction:

```r
outcome ~ exposure + potential_confounder + exposure:potential_confounder
```

No other covariates are included in this interaction screen. Its interaction
p-value helps organise early evidence; it is not a final covariate-adjusted
interaction analysis.

Use `interaction_models()` for a planned interaction analysis. Its `covariates`
are included in both nested models, preserving a common adjusted analysis
sample. For example, `covariates = c(age, lwt)` fits:

```r
outcome ~ exposure + effect_modifier + age + lwt
outcome ~ exposure + effect_modifier + age + lwt + exposure:effect_modifier
```

Choose covariates using the study design, a DAG, and subject-matter knowledge.
Do not automatically adjust for mediators or colliders.

The formatted interaction table labels the **effect modifier tested**, displays
the outcome, and records the two fitted model formulas in its table notes.

```{r ci-setup, message=FALSE, warning=FALSE}
library(gtregression)
library(dplyr)

data("data_birthwt", package = "gtregression")

birthwt_data <- data_birthwt |>
  mutate(
    race = factor(race, levels = c(1, 2, 3),
                  labels = c("White", "Black", "Other")),
    smoke = factor(smoke, levels = c(0, 1), labels = c("No", "Yes")),
    ht = factor(ht, levels = c(0, 1), labels = c("No", "Yes")),
    low = factor(low, levels = c(0, 1), labels = c("Normal BW", "Low BW"))
  )

attr(birthwt_data$race, "label") <- "Maternal race"
attr(birthwt_data$smoke, "label") <- "Smoking during pregnancy"
attr(birthwt_data$ht, "label") <- "Hypertension"
```

## Identify Confounders

Use `method = "change"` for the model-based change-in-estimate method. Use
`method = "mh"` or `method = "both"` when Mantel-Haenszel is appropriate.
The output is intentionally tidy and intended for viewing, not as a final
publication table. Calling the object displays the formatted table; use
`$summary` when you want the underlying tibble. Its table notes identify the
outcome and the three candidate-level models: crude, adjusted for the potential
confounder, and the interaction screen.

```{r confounder, message=FALSE, warning=FALSE}
confounder_check <- identify_confounder(
  data = birthwt_data,
  outcome = low,
  exposure = smoke,
  potential_confounder = c("race", "ht"),
  approach = logit,
  method = both,
  format = gt
)

confounder_check$table
```

The underlying summary remains available for inspection or further filtering.

```{r confounder-summary, message=FALSE, warning=FALSE}
confounder_check$summary
```

## Mantel-Haenszel Estimate

Mantel-Haenszel is useful when the question is whether a stratified pooled
estimate differs meaningfully from the crude estimate. It is available for
eligible binary/categorical settings. It is not the same as fitting an
interaction term.

```{r confounder-mh, message=FALSE, warning=FALSE}
identify_confounder(
  data = birthwt_data,
  outcome = low,
  exposure = smoke,
  potential_confounder = race,
  approach = logit,
  method = mh,
  format = flextable
)$table
```

## Test Interaction

`interaction_models()` compares models with and without the interaction term.
It is deliberately model-based and uses `LRT` or `Wald`, not
Mantel-Haenszel. Use it when the interaction term is planned or supported by
clinical, causal, or subject-matter reasoning. Covariates are included in both
models, so this is the appropriate function for an adjusted interaction test.

```{r interaction, message=FALSE, warning=FALSE}
interaction_check <- interaction_models(
  data = birthwt_data,
  outcome = low,
  exposure = smoke,
  effect_modifier = race,
  covariates = c("age", "lwt"),
  approach = logit,
  test = LRT,
  format = gt
)

interaction_check$table
```

## Survival Confounding and Interaction

The same grammar works for survival models. For `cox` and `surv`/`survreg`,
use `time` and `event` instead of `outcome`. Cox models report hazard-ratio
style estimates; parametric survival models report time-ratio style estimates.
Mantel-Haenszel screening is for binary outcome settings, so survival examples
use model-based change-in-estimate and interaction checks.

```{r survival-confounder, message=FALSE, warning=FALSE}
lung_data <- data_lungcancer |>
  dplyr::mutate(
    trt = factor(trt, levels = c(1, 2),
                 labels = c("Standard treatment", "Test treatment")),
    prior = factor(prior, levels = c(0, 10), labels = c("No", "Yes"))
  )

survival_confounder <- identify_confounder(
  data = lung_data,
  time = time,
  event = status,
  exposure = trt,
  potential_confounder = prior,
  approach = cox,
  method = change,
  format = gt
)

survival_confounder$table
```

```{r survival-interaction, message=FALSE, warning=FALSE}
survival_interaction <- interaction_models(
  data = lung_data,
  time = time,
  event = status,
  exposure = trt,
  effect_modifier = prior,
  covariates = c(age, karno),
  approach = cox,
  test = LRT,
  format = gt
)

survival_interaction$table
```

## Causal Mediation

`mediation_analysis()` asks whether part of an exposure-outcome association may
operate through a mediator. In this example, obesity is the exposure, plasma
glucose is the mediator, and diabetes is the outcome. The question is not only
"is obesity associated with diabetes?", but also "how much of that association
may operate through plasma glucose?".

The default output is a publication-ready flextable. The table reports:

- **Total effect**: the overall exposure-outcome association.
- **Direct effect**: the part not operating through the mediator.
- **Indirect effect**: the part operating through the mediator.
- **Proportion mediated**: the share of the total effect attributed to the
  indirect pathway.

For logistic outcomes, effects are reported as predicted probability
differences, not odds ratios. For example, an estimate of `0.068` means about a
6.8 percentage-point difference on the predicted probability scale.

Treat mediation output as causal only when the usual mediation assumptions are
supported by study design, DAGs, temporality, and subject-matter knowledge. In
particular, there should be no unmeasured exposure-outcome,
exposure-mediator, or mediator-outcome confounding, and the mediator should
occur before the outcome.

```{r mediation, message=FALSE, warning=FALSE}
data("data_diabetes_mediation", package = "gtregression")

diabetes_med <- mediation_analysis(
  data = data_diabetes_mediation,
  exposure = obesity,
  mediator = glucose,
  outcome = diabetes,
  covariates = c(age, blood_pressure, pregnancies, diabetes_pedigree),
  outcome_approach = logit,
  sims = 100,
  seed = 123
)

diabetes_med
```

The underlying values remain available for audit, reporting, or custom
formatting.

```{r mediation-body}
diabetes_med$table_body
```

Use `format = gt` when preparing HTML-first outputs such as websites or
teaching pages.

```{r mediation-gt, message=FALSE, warning=FALSE}
med_gt <- mediation_analysis(
  data = data_diabetes_mediation,
  exposure = obesity,
  mediator = glucose,
  outcome = diabetes,
  covariates = c(age, blood_pressure, pregnancies, diabetes_pedigree),
  outcome_approach = logit,
  format = gt,
  sims = 100,
  seed = 123
)

med_gt$table
```

`plot_mediation()` draws the same planned causal structure as a path diagram.
It is useful for teaching, presentations, or checking that the exposure,
mediator, and outcome have been specified as intended.

```{r mediation-plot, message=FALSE, warning=FALSE, fig.width=7, fig.height=4}
plot_mediation(diabetes_med)
```

If the estimates make the plot too busy, hide them and use the diagram only to
show the assumed path structure.

```{r mediation-plot-simple, message=FALSE, warning=FALSE, fig.width=7, fig.height=4}
plot_mediation(diabetes_med, show_estimates = FALSE)
```

## Overlap and Difference

| Topic | `identify_confounder()` | `interaction_models()` |
|---|---|---|
| Main purpose | Organises candidate confounder and effect-modifier screening signals. | Tests a planned exposure-by-modifier term. |
| Typical input | Exposure plus one or more candidate variables. | One exposure and one effect modifier. |
| Confounding | Crude vs adjusted change-in-estimate; optional Mantel-Haenszel comparison. | Not designed for confounder selection. |
| Effect modification | Screening signal from stratum-specific estimates and interaction p-value. | Model comparison using LRT or Wald test. |
| Best use | Early review of candidate variables, with DAGs and judgement. | Focused test of a clinically or biologically plausible interaction. |
| Output status | Viewing aid, not publication-ready evidence by itself. | Viewing aid; report with stratum-specific estimates when relevant. |

`mediation_analysis()` is different from both functions. It decomposes one
planned exposure-outcome relationship into direct and mediator-related
components; it does not select confounders or test whether effects differ across
strata.

## What To Inspect

- `identify_confounder()`: `$summary`, `$table`, `$details`, `$mh_estimate`,
  `$mh_status`, `$decision`, and `$recommendation`.
- `interaction_models()`: `$summary`, `$table`, `$p_value`, `$decision`, and
  fitted model objects.
- `mediation_analysis()`: `$table`, `$table_body`, `$models`, `$boot`,
  `$values`, `$variable_labels`, and `$complete_data`.
- Use subject-matter knowledge with these outputs. The functions support
  interpretation; they do not replace the study design.