inst/doc/model-fit.R

In postcard: Estimating Marginal Effects with Prognostic Covariate Adjustment

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

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

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

dat_gaus <- glm_data(
  Y ~ b0+b1*sin(W)^2+b2*A,
  W = runif(n, min = -2, max = 2),
  A = rbinom(n, 1, prob = 1/2)
)

dat_gaus_hist <- glm_data(
  Y ~ b0+b1*sin(W)^2,
  W = runif(n, min = -2, max = 2)
)

dat_pois <- glm_data(
  Y ~ b0+b1*sin(W)^2+b2*A,
  W = runif(n, min = -2, max = 2),
  A = rbinom(n, 1, 1/2),
  family = poisson(link = "log")
)

## -----------------------------------------------------------------------------
 # Default amount of printing
ate <- rctglm(
  formula = Y ~ A + W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_gaus,
  verbose = 2)
ate_prog <- rctglm_with_prognosticscore(
  formula = Y ~ A + W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_gaus,
  data_hist = dat_gaus_hist,
  verbose = 2)
# At little less printing
ate <- rctglm(
  formula = Y ~ A + W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_gaus,
  verbose = 1)
ate_prog <- rctglm_with_prognosticscore(
  formula = Y ~ A + W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_gaus,
  data_hist = dat_gaus_hist,
  verbose = 1)
# No printing
ate <- rctglm(
  formula = Y ~ A + W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_gaus,
  verbose = 0)
ate_prog <- rctglm_with_prognosticscore(
  formula = Y ~ A + W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_gaus,
  data_hist = dat_gaus_hist,
  verbose = 0)

## ----rate-ratio-run-show------------------------------------------------------
rate_ratio <- rctglm(
  formula = Y ~ A + W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_pois,
  family = "poisson",
  estimand_fun = "rate_ratio",
  verbose = 1)
rate_ratio$estimand_funs

## ----nonsense-run-show--------------------------------------------------------
nonsense_estimand_fun <- function(psi1, psi0) {
  psi1 / sqrt(psi0) * 2 - 1
}

nonsense_estimand <- rctglm(
  formula = Y ~ A * W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_pois,
  family = poisson(),
  estimand_fun = nonsense_estimand_fun,
  verbose = 1)
nonsense_estimand$estimand_funs

## ----default-learners-show----------------------------------------------------
default_learners

## ----fit-own-learners---------------------------------------------------------
learners <- list(
  rf = list(
    model = parsnip::rand_forest(
      mode = "regression",
      trees = 500,
      min_n = parsnip::tune("min_n")
    ) %>% 
      parsnip::set_engine("ranger"),
    grid = data.frame(
      min_n = 1:10
    )
  ),
  svm.linear = list(
    model = parsnip::svm_linear(
      mode = "regression",
      cost = parsnip::tune("cost"),
      margin = parsnip::tune("margin")) %>% 
      parsnip::set_engine("LiblineaR"),
    grid = data.frame(
      cost = 1:5,
      margin = seq(0.1, 0.5, 0.1)
    )
  )
)

model_own_learners <- rctglm_with_prognosticscore(
  formula = Y ~ A * W,
  exposure_indicator = A,
  exposure_prob = 1/2,
  data = dat_gaus,
  data_hist = dat_gaus_hist,
  learners = learners)

## ----prognostic-info-show-----------------------------------------------------
prog_info <- prog(model_own_learners)
prog_info$data <- head(prog_info$data)
prog_info