Heatmap and Time Series

library(data.table)
library(ggplot2)

pat_105 = fread('~/Downloads/pat_105.csv')
head(pat_105)

##    pat_id obs_week_id     tod     m_TIR same_tod_contact      m_g     sd_g
## 1:    105           2  Dinner 0.8194444                0 141.3472 29.92213
## 2:    105           2   Night 0.8256173                0 149.7546 33.88827
## 3:    105           2   Lunch 0.8566667                0 144.2279 23.13592
## 4:    105           2 Morning 0.9888889                0 146.7694 10.11993
## 5:    105           3  Dinner 0.8154762                0 149.8244 20.86109
## 6:    105           3   Lunch 0.8380952                0 154.2911 19.51535
##      p_70_g   p_30_g
## 1: 157.0833 125.5667
## 2: 171.5667 133.5333
## 3: 155.3000 139.2800
## 4: 151.2800 140.8000
## 5: 160.7143 134.9571
## 6: 158.6143 143.1714

Heatmap

pat_105[, tod := factor(tod, levels=c("Morning", "Lunch", "Dinner", "Night"))]

ggplot(pat_105,
       aes(x=tod, y=obs_week_id, fill=m_TIR)) +
    geom_tile(linewidth=0.3) + scale_fill_gradient2(low = "red", mid="white", high = "steelblue", midpoint=0.65, limits=c(0,1)) +
    scale_y_reverse() + ggtitle(paste("Pat 105 TIR by week, time of day")) + geom_point(data=pat_105[same_tod_contact==1], fill="black")

TIR Time series

ggplot(pat_105, aes(x=obs_week_id, y=m_TIR)) +
    geom_line() + ggtitle(paste("Pat 105 TIR")) + facet_grid(rows=vars(tod)) +
    theme_bw() +
    geom_point(data=pat_105[same_tod_contact==1], size=2)

Mean Glucose Time Series

ggplot(pat_105,
       aes(x=obs_week_id, y=m_g, ymin=m_g-sd_g, ymax=m_g+sd_g)) +
    geom_line() + ggtitle(paste("Patient105 Mean Glucose (interventions in green)")) + facet_grid(rows=vars(tod)) + geom_hline(yintercept=180, linetype="dashed") +
    geom_hline(yintercept=70, linetype="dashed") + 
    theme_bw() + geom_ribbon(alpha=0.3) + geom_ribbon(alpha=0.3, fill="red", aes(ymin=180, ymax=pmax(180,m_g+sd_g))) +
    geom_point(data=pat_105[same_tod_contact==1], size=4, shape=23, fill="green") + ylab("Mean Glucose") + xlab("Day in Program")

Power Simulation Results

combo_output = fread('~/Downloads/combo_output_pow_sims.csv')

head(combo_output)

##    policy n_weeks n_patients_targeted_each_week     est_ate mean_ate_se  ppos
## 1: Random      26                            15 0.009923639 0.011025308 0.148
## 2: Random      26                            20 0.009807970 0.009529174 0.190
## 3: Random      26                            25 0.009679504 0.008529523 0.208
## 4: Random      26                            30 0.009801098 0.007776239 0.242
## 5: Random      26                            35 0.009869773 0.007208651 0.278
## 6: Random      26                            40 0.009788730 0.006735489 0.290

ATE plots

ggplot(combo_output[n_weeks==52], aes(x=n_patients_targeted_each_week, y=est_ate)) + facet_wrap(~policy) + geom_line() + geom_point() + theme_bw() + coord_cartesian(ylim=c(0,0.03))

PPOS plots

ggplot(combo_output, aes(x=n_patients_targeted_each_week, y=ppos, color=factor(n_weeks))) + facet_wrap(~policy) + geom_line() + geom_point() + theme_bw() + geom_hline(yintercept=0.9, linetype='dashed')

Code used to extract power sim data for reference

Click for code

  load('~/Downloads/powerenv.RData')
  
  agg_output = function(out) {
    agg_out = out[,.(
   power = mean(pval < 0.05 & coef>0),
   m_se = mean(se),
   m_coef = mean(coef),
   m_coef_stat_sig = mean(ifelse(pval<0.05, coef, NA), na.rm=T),
   m_n_unique_pats_included_from_input_data = mean(n_unique_pats_included_from_input_data)
    ), by=c('n_weeks', 'n_patients_pop_per_week', 'n_patients_included_per_week', 'permuted')]
    agg_out[, coef_bias := m_coef_stat_sig - m_coef]
    return(agg_out)
  }
  
  combo_output = rbind(
    agg_output(many_sims_random)[,policy:='Random'],
    agg_output(many_sims_low_tir)[,policy:='Lowest TIR'],
    agg_output(many_sims_pcontactonly)[,policy:='Largest P(contact)'],
    agg_output(many_sims_p_contact_exp_dTIR)[,policy:='Largest expected effect']
  )
  
  combo_output = combo_output[permuted==FALSE, .(
    policy, n_weeks,
    n_patients_targeted_each_week = n_patients_included_per_week,
    est_ate = m_coef,
    mean_ate_se = m_se,
    ppos = power
  )]
  
  fwrite(combo_output, '~/Downloads/combo_output_pow_sims.csv')