simulations/summarize_estimator.R

   1 library(ggdist)
   2
   3 summarize.estimator <- function(sims.df, suffix='naive', coefname='x'){
   4
   5     reported_vars <- c(
   6                        'Bxy',
   7                        paste0('B',coefname,'y.est.',suffix),
   8                        paste0('B',coefname,'y.ci.lower.',suffix),
   9                        paste0('B',coefname,'y.ci.upper.',suffix)
  10                        )
  11
  12
  13     grouping_vars <- c('N','m','B0', 'Bxy', 'Bzy', 'Bzx', 'Px', 'Py','y_explained_variance', 'prediction_accuracy','outcome_formula','proxy_formula','truth_formula','z_bias','y_bias')
  14
  15     grouping_vars <- grouping_vars[grouping_vars %in% names(df)]
  16
  17     part <- sims.df[,
  18                     unique(c(reported_vars,
  19                              grouping_vars)),
  20                     with=FALSE]
  21
  22
  23     true.in.ci <- as.integer((part$Bxy >= part[[paste0('B',coefname,'y.ci.lower.',suffix)]]) & (part$Bxy <= part[[paste0('B',coefname,'y.ci.upper.',suffix)]]))
  24     zero.in.ci <- as.integer(0 >= part[[paste0('B',coefname,'y.ci.lower.',suffix)]]) & (0 <= part[[paste0('B',coefname,'y.ci.upper.',suffix)]])
  25     bias <- part[[paste0('B',coefname,'y')]] - part[[paste0('B',coefname,'y.est.',suffix)]]
  26     sign.correct <- as.integer(sign(part$Bxy) == sign(part[[paste0('B',coefname,'y.est.',suffix)]]))
  27
  28     part <- part[,':='(true.in.ci = true.in.ci,
  29                        zero.in.ci = zero.in.ci,
  30                        bias=bias,
  31                        sign.correct =sign.correct)]
  32
  33
  34     part.plot <- part[, .(p.true.in.ci = mean(true.in.ci),
  35                           mean.bias = mean(bias),
  36                           mean.est = mean(.SD[[paste0('B',coefname,'y.est.',suffix)]],na.rm=T),
  37                           var.est = var(.SD[[paste0('B',coefname,'y.est.',suffix)]],na.rm=T),
  38                           est.upper.95 = quantile(.SD[[paste0('B',coefname,'y.est.',suffix)]],0.975,na.rm=T),
  39                           est.lower.95 = quantile(.SD[[paste0('B',coefname,'y.est.',suffix)]],0.025,na.rm=T),
  40                           mean.ci.upper = mean(.SD[[paste0('B',coefname,'y.ci.upper.',suffix)]],na.rm=T),
  41                           mean.ci.lower = mean(.SD[[paste0('B',coefname,'y.ci.lower.',suffix)]],na.rm=T),
  42                           median.ci.upper = median(.SD[[paste0('B',coefname,'y.ci.upper.',suffix)]],na.rm=T),
  43                           median.ci.lower = median(.SD[[paste0('B',coefname,'y.ci.lower.',suffix)]],na.rm=T),
  44                           ci.upper.975 = quantile(.SD[[paste0('B',coefname,'y.ci.upper.',suffix)]],0.975,na.rm=T),
  45                           ci.upper.025 = quantile(.SD[[paste0('B',coefname,'y.ci.upper.',suffix)]],0.025,na.rm=T),
  46                           ci.lower.975 = quantile(.SD[[paste0('B',coefname,'y.ci.lower.',suffix)]],0.975,na.rm=T),
  47                           ci.lower.025 = quantile(.SD[[paste0('B',coefname,'y.ci.lower.',suffix)]],0.025,na.rm=T),
  48                           N.ci.is.NA = sum(is.na(.SD[[paste0('B',coefname,'y.ci.lower.',suffix)]])),
  49                           N.sims = .N,
  50                           p.sign.correct = mean(as.integer(sign.correct & (! zero.in.ci))),
  51                           variable=coefname,
  52                           method=suffix
  53                           ),
  54                       by=grouping_vars,
  55                       ]
  56
  57     return(part.plot)
  58 }