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

## ----setup--------------------------------------------------------------------
library(delaydiscount)
library(dplyr)

## -----------------------------------------------------------------------------
set.seed(8678)
sim_data <- simulate_dataset(groups = c("EFT", "HIT", "NCC"),
                             num_subj = c(120, 142, 172),
                             time_points = c(30, 90, 180, 365, 1095, 1825, 3650),
                             mean_ln_k = c(-6.877, -5.861, -6.078),
                             sigma_sq = 1.978, g = 10.407)

## -----------------------------------------------------------------------------
# These methods can be run on the output from simulate_dataset
prep_data <- prepare_data_frame(sim_data)
rule_check_results <- jb_rule_check(sim_data)

# These methods need to be run on the output from prepare_data_frame
subj_est_ln_k <- get_subj_est_ln_k(prep_data)
sim_model <- dd_hyperbolic_model(prep_data)

## -----------------------------------------------------------------------------
# Get a data frame consisting of one observation per subject, containing each
#  subject's true ln(k)
subj_true_ln_k <- sim_data %>%
  select(subj, true_ln_k, group) %>%
  unique()

# Merge the data frames containing the true subject ln_k values, estimated
#  subject ln_k values, and rule check results.
subj_ln_k_comp <- merge(subj_est_ln_k, subj_true_ln_k)
subj_ln_k_comp_wrc <- merge(subj_ln_k_comp, rule_check_results)

# Get the overall rule check result for each subject
subj_ln_k_comp_wrc$rc_pass <- subj_ln_k_comp_wrc$C1 & subj_ln_k_comp_wrc$C2

# Make a plot showing true vs estimated ln(k) for each subject,
#  coloring the points based on whether or not the corresponding subject
#  passed the rule check.
plot(x = subj_ln_k_comp_wrc$true_ln_k, y = subj_ln_k_comp_wrc$ln_k,
     col = ifelse(subj_ln_k_comp_wrc$rc_pass, "green", "blue"),
     main = "Subject True vs Maximum Likelihood Estimated ln(k)",
     xlab = "True subject ln(k)", ylab = "Estimated subject ln(k)")
abline(a=0, b=1)
legend("topleft", 
       legend = c(
         "Passes rule check",
         "Fails rule check"),
       col = c("green", "blue"), pch = 1)