## ----setup, include=FALSE-----------------------------------------------------
knitr::opts_chunk$set(
  collapse  = TRUE,
  comment   = "#>",
  fig.width = 7,
  fig.height = 4
)

## ----load, message=FALSE, warning=FALSE---------------------------------------
library(causalMed)
library(data.table)

## ----data---------------------------------------------------------------------
data("nonsurvivaldata")
head(nonsurvivaldata, 10)

## ----spec-model-args, eval=FALSE----------------------------------------------
# spec_model(
#   formula,           # response ~ predictors
#   var_type,          # how values are simulated
#   mod_type,          # role in the causal structure
#   subset = NULL,     # optional: restrict the model to some rows
#   recode = NULL      # optional: recodes() applied before this model
# )

## ----recodes-example----------------------------------------------------------
init_rc <- recodes(lag1_A = 0, lag1_L1 = 0, lag1_L2 = 0)   # first time step
in_rc   <- recodes(lag1_A = A, lag1_L1 = L1, lag1_L2 = L2) # later steps

## ----quickstart-gformula------------------------------------------------------
models_te <- list(
  spec_model(A     ~ V + lag1_A + lag1_L1 + lag1_L2 + time,
             var_type = "binary", mod_type = "exposure"),
  spec_model(L1    ~ V + A + lag1_L1 + time,
             var_type = "normal", mod_type = "covariate"),
  spec_model(L2    ~ V + A + lag1_L2 + time,
             var_type = "binary", mod_type = "covariate"),
  spec_model(Y_bin ~ V + A + L1 + L2,
             var_type = "binary", mod_type = "outcome")
)

fit_te <- gformula(
  data         = nonsurvivaldata,
  id_var       = "id",
  time_var     = "time",
  base_vars    = "V",
  exposure     = "A",
  models       = models_te,
  intervention = list(natural = NULL, always_treat = 1, never_treat = 0),
  ref_int      = "natural",
  init_recode  = init_rc,
  in_recode    = in_rc,
  mc_sample    = 10000,
  R            = 1,        # set R > 1 for bootstrap confidence intervals
  quiet        = TRUE,
  seed         = 2025
)

fit_te$effect_size   # mean outcome under each strategy
fit_te$estimate      # risk difference and ratio against the natural course

## ----quickstart-mediation-----------------------------------------------------
init_med <- recodes(lag1_A = 0, lag1_L1 = 0, lag1_L2 = 0, lag1_M = 0)
in_med   <- recodes(lag1_A = A, lag1_L1 = L1, lag1_L2 = L2, lag1_M = M)

models_med <- list(
  spec_model(A     ~ V + lag1_A + lag1_L1 + lag1_L2 + time,
             var_type = "binary", mod_type = "exposure"),
  spec_model(L1    ~ V + A + lag1_L1 + time,
             var_type = "normal", mod_type = "covariate"),
  spec_model(L2    ~ V + A + lag1_L2 + time,
             var_type = "binary", mod_type = "covariate"),
  spec_model(M     ~ V + A + L1 + L2 + lag1_M + time,
             var_type = "normal", mod_type = "mediator"),
  spec_model(Y_bin ~ V + A + M + L1 + L2,
             var_type = "binary", mod_type = "outcome")
)

fit_med <- mediation(
  data        = nonsurvivaldata,
  id_var      = "id",
  time_var    = "time",
  base_vars   = "V",
  exposure    = "A",
  outcome     = "Y_bin",
  models      = models_med,
  init_recode = init_med,
  in_recode   = in_med,
  mc_sample   = 10000,
  R           = 1,
  quiet       = TRUE,
  seed        = 2025
)

fit_med$estimate

## ----parallel, eval=FALSE-----------------------------------------------------
# future::plan(future::multisession)
# fit <- gformula(..., R = 500)
# future::plan(future::sequential)