## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

## ----setup, message = FALSE---------------------------------------------------
library(postcard)
withr::local_seed(1395878)
withr::local_options(list(postcard.verbose = 0))

## ----dat-sim------------------------------------------------------------------
n <- 1000
b0 <- 1
b1 <- 3
b2 <- 2

# Simulate data with a non-linear effect
dat_treat <- glm_data(
  Y ~ b0+b1*log(W)+b2*A,
  W = runif(n, min = 1, max = 100),
  A = rbinom(n, 1, prob = 1/2),
  family = gaussian() # Default value
)

## ----noprog-run---------------------------------------------------------------
ate <- rctglm(formula = Y ~ A * W,
              exposure_indicator = A,
              exposure_prob = 1/2,
              data = dat_treat,
              family = "gaussian") # Default value

## ----noprog-show--------------------------------------------------------------
ate

## -----------------------------------------------------------------------------
est(ate)

## ----hist-data----------------------------------------------------------------
dat_notreat <- glm_data(
  Y ~ b0+b1*log(W),
  W = runif(n, min = 1, max = 100),
  family = gaussian # Default value
)

## ----prog-run-----------------------------------------------------------------
ate_prog <- rctglm_with_prognosticscore(
  formula = Y ~ A * W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_treat,
  family = gaussian(link = "identity"), # Default value
  data_hist = dat_notreat)

## ----prog-show----------------------------------------------------------------
ate_prog

## -----------------------------------------------------------------------------
ancova <- glm(Y ~ W, data = dat_notreat)
preds_anc <- predict(ancova, dat_notreat)
lrnr <- fit_best_learner(list(mod = Y ~ W), data = dat_notreat)
preds_dsl <- dplyr::pull(predict(lrnr, new_data = dat_notreat))

power_marginaleffect(
  response = dat_notreat$Y,
  predictions = preds_anc,
  target_effect = 0.3,
  exposure_prob = 1/2
)
power_marginaleffect(
  response = dat_notreat$Y,
  predictions = preds_dsl,
  target_effect = 0.3,
  exposure_prob = 1/2
)

## -----------------------------------------------------------------------------
dat_notreat <- dplyr::mutate(dat_notreat, .pred = preds_dsl)

var_bound_ancova <- variance_ancova(Y ~ W, data = dat_notreat)
var_bound_prog <- variance_ancova(Y ~ W + .pred, data = dat_notreat)

## -----------------------------------------------------------------------------
samplesize_gs(variance = var_bound_ancova,
              power = 0.9, r = 1, ate = .8, margin = 0)
samplesize_gs(variance = var_bound_prog,
              power = 0.9, r = 1, ate = .8, margin = 0)

## -----------------------------------------------------------------------------
rpm <- repeat_power_marginaleffect(
  target_effect = 1.3, exposure_prob = 1/2,
  ns = 10:200, n_iter = 10
)

## -----------------------------------------------------------------------------
plot(rpm)