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

# Use a wider console width for readable tabular output, while preserving and
# restoring the user's original option at the end of the vignette.
old_options <- options(width = 140)

library(twoPhaseSens)

## ----example-data-------------------------------------------------------------------------------------------------------------------------
dat <- twophase_example_data(
  n = 300,
  seed = 2026
)

head(dat)

## ----phase2-structure---------------------------------------------------------------------------------------------------------------------
table(dat$phase2)

with(
  dat,
  table(
    phase2 = phase2,
    complete_marker_block = complete.cases(marker1, marker2, marker3)
  )
)

## ----one-outcome--------------------------------------------------------------------------------------------------------------------------
fit_one <- twophase_sensitivity(
  data = dat,
  outcome = "outcome1",
  exposure = "exposure",
  covariates = c("age", "sex"),
  phase2_covariates = c("marker1", "marker2", "marker3"),
  phase2 = "phase2",
  methods = c("naive", "cca", "ipw"),
  seed = 1001
)

fit_one

## ----one-outcome-summary------------------------------------------------------------------------------------------------------------------
summary(fit_one)

## ----one-outcome-table--------------------------------------------------------------------------------------------------------------------
tbl_sensitivity(
  fit_one,
  outcome_labels = c(outcome1 = "Outcome 1"),
  caption = "**Sensitivity analysis for Outcome 1**"
)

## ----one-outcome-data-frame---------------------------------------------------------------------------------------------------------------
tbl_sensitivity(
  fit_one,
  outcome_labels = c(outcome1 = "Outcome 1"),
  output = "data.frame"
)

## ----decimals-----------------------------------------------------------------------------------------------------------------------------
tbl_sensitivity(
  fit_one,
  outcome_labels = c(outcome1 = "Outcome 1"),
  estimate_digits = 2,
  p_digits = 4,
  output = "data.frame"
)

## ----se-column----------------------------------------------------------------------------------------------------------------------------
tbl_sensitivity(
  fit_one,
  outcome_labels = c(outcome1 = "Outcome 1"),
  include_se = TRUE,
  output = "data.frame"
)

## ----multiple-outcomes--------------------------------------------------------------------------------------------------------------------
fit_two <- twophase_sensitivity(
  data = dat,
  outcome = c("outcome1", "outcome2"),
  exposure = "exposure",
  covariates = c("age", "sex"),
  phase2_covariates = c("marker1", "marker2", "marker3"),
  phase2 = "phase2",
  methods = c("naive", "cca", "ipw"),
  seed = 2001
)

## ----inspect-seeds------------------------------------------------------------------------------------------------------------------------
unique(fit_two$results[c("outcome", "analysis_seed")])

## ----multiple-outcome-table---------------------------------------------------------------------------------------------------------------
tbl_sensitivity(
  fit_two,
  outcome_labels = c(
    outcome1 = "Outcome 1",
    outcome2 = "Outcome 2"
  )
)

## ----explicit-seeds-----------------------------------------------------------------------------------------------------------------------
fit_seeded <- twophase_sensitivity(
  data = dat,
  outcome = c("outcome1", "outcome2"),
  exposure = "exposure",
  covariates = c("age", "sex"),
  phase2_covariates = c("marker1", "marker2", "marker3"),
  phase2 = "phase2",
  methods = c("naive", "cca"),
  seed = c(
    outcome2 = 9002,
    outcome1 = 9001
  )
)

unique(fit_seeded$results[c("outcome", "analysis_seed")])

## ----method-subset------------------------------------------------------------------------------------------------------------------------
fit_subset <- twophase_sensitivity(
  data = dat,
  outcome = "outcome1",
  exposure = "exposure",
  covariates = c("age", "sex"),
  phase2_covariates = c("marker1", "marker2", "marker3"),
  phase2 = "phase2",
  methods = c("cca", "fcs_mi", "ipw"),
  n_imputations = 5,
  mice_maxit = 5,
  seed = 3001
)

tbl_sensitivity(
  fit_subset,
  outcome_labels = c(outcome1 = "Outcome 1")
)

## ----mi-settings, eval=FALSE--------------------------------------------------------------------------------------------------------------
# fit_full <- twophase_sensitivity(
#   data = dat,
#   outcome = "outcome1",
#   exposure = "exposure",
#   covariates = c("age", "sex"),
#   phase2_covariates = c("marker1", "marker2", "marker3"),
#   phase2 = "phase2",
#   n_imputations = 50,
#   mice_maxit = 20,
#   jomo_nburn = 2000,
#   jomo_nbetween = 1000,
#   seed = 4001
# )

## ----confidence-level---------------------------------------------------------------------------------------------------------------------
fit_90 <- twophase_sensitivity(
  data = dat,
  outcome = "outcome1",
  exposure = "exposure",
  covariates = c("age", "sex"),
  phase2_covariates = c("marker1", "marker2", "marker3"),
  phase2 = "phase2",
  methods = c("naive", "cca"),
  conf_level = 0.90
)

tbl_sensitivity(
  fit_90,
  outcome_labels = c(outcome1 = "Outcome 1"),
  output = "data.frame"
)

## ----diagnostics--------------------------------------------------------------------------------------------------------------------------
diagnostics(fit_one)

## ----diagnostics-status-------------------------------------------------------------------------------------------------------------------
diagnostics(fit_one, type = "status")

## ----diagnostics-weights------------------------------------------------------------------------------------------------------------------
diagnostics(fit_one, type = "weights")

## ----diagnostics-mi, eval=FALSE-----------------------------------------------------------------------------------------------------------
# diagnostics(fit_full, type = "mi")

## ----diagnostics-filter-------------------------------------------------------------------------------------------------------------------
diagnostics(
  fit_two,
  outcome = "outcome2",
  method = "ipw"
)

## ----one-outcome-forest-------------------------------------------------------------------------------------------------------------------
plot_sensitivity(
  fit_one,
  outcome_labels = c(outcome1 = "Outcome 1")
)

## ----one-outcome-forest-generic, eval=FALSE-----------------------------------------------------------------------------------------------
# plot(fit_one)

## ----multiple-outcome-forest, fig.width=12.5, fig.height=5.4------------------------------------------------------------------------------
plot_sensitivity(
  fit_two,
  outcome_labels = c(
    outcome1 = "Outcome 1",
    outcome2 = "Outcome 2"
  ),
  annotate = TRUE
)

## ----multiple-outcome-forest-fixed, fig.width=12.5, fig.height=5.4------------------------------------------------------------------------
plot_sensitivity(
  fit_two,
  outcome_labels = c(
    outcome1 = "Outcome 1",
    outcome2 = "Outcome 2"
  ),
  facet_scales = "fixed",
  annotate = TRUE
)

## ----save-forest-examples, eval=FALSE-----------------------------------------------------------------------------------------------------
# p_one <- plot_sensitivity(
#   fit_one,
#   outcome_labels = c(outcome1 = "Primary outcome"),
#   estimate_label = "Adjusted beta"
# )
# 
# ggplot2::ggsave(
#   "sensitivity_single_outcome.png",
#   p_one,
#   width = 9.5,
#   height = 6,
#   dpi = 320,
#   bg = "white"
# )
# 
# p_two <- plot_sensitivity(
#   fit_two,
#   outcome_labels = c(
#     outcome1 = "Outcome 1",
#     outcome2 = "Outcome 2"
#   ),
#   estimate_label = "Adjusted beta",
#   facet_scales = "fixed"
# )
# 
# ggplot2::ggsave(
#   "sensitivity_two_outcomes.png",
#   p_two,
#   width = 14,
#   height = 6,
#   dpi = 320,
#   bg = "white"
# )

## ----custom-labels------------------------------------------------------------------------------------------------------------------------
tbl_sensitivity(
  fit_one,
  method_labels = c(
    Naive = "Phase-1 naive",
    CCA = "Complete case"
  ),
  outcome_labels = c(
    outcome1 = "Primary outcome"
  ),
  output = "data.frame"
)

## ----restore-options, include=FALSE-------------------------------------------
options(old_options)