R/plot_posteriors2_v7.R

In jean997/sherlockAsh: CAUSE: Causal Analysis Using Summary Effect Estimates

plot_posteriors <- function(fits, elpd){
  plts <- lapply(fits, function(fit){
        post_marge <- lapply(seq_along(fit$post_marge), function(i){
                     p <- fit$post_marge[[i]]
                     p$param <- fit$params[i]
                     return(p)})
        post_marge <- do.call(rbind, post_marge)
        post_marge <- select(post_marge, mid, width, param, post, prior) %>%
              gather("dist", "pdf", -mid, -param, -width)
        post_marge$param <- factor(post_marge$param, levels =c("gamma", "eta", "q"))
        plt <- ggplot(post_marge) + geom_line(aes(x=mid, y=pdf/width, linetype=dist)) +
          xlab("parameter value") +
          theme_bw() + theme(legend.position="none",
                             axis.title.y=element_blank(),
                             axis.title.x=element_blank()) +
          facet_wrap(~param, scale="free")
        return(plt)})

  tab <- sapply(fits, function(fit){
            gamma <- eta <- q <- NA
            full_params <- c("gamma", "eta", "q")
            tt <- c()
            for(j in seq_along(full_params)){
                if(!full_params[j] %in% fit$params){
                    tt[j] <- NA
                }else{
                    ix <- which(fit$params == full_params[j])
                    med <- with(fit$post_marge[[ix]], step_quantile(0.5, begin, end, post))
                    x025 <- with(fit$post_marge[[ix]], step_quantile(0.025, begin, end, post))
                    x975 <- with(fit$post_marge[[ix]], step_quantile(0.975, begin, end, post))
                    tt[j]  <- paste0(round(med, digits=2),
                                " (", round(x025, digits=2), ", ",
                                round(x975, digits=2), ")")
                }
            }
            names(tt) <- full_params
            return(tt)
       })
  tab <- data.frame(t(tab))
  plts[[4]] <- tableGrob(tab)
  elpd <- elpd %>% mutate(z_score = delta_elpd/se_delta_elpd,
                               delta_elpd = round(delta_elpd, digits=2),
                               se_delta_elpd = round(se_delta_elpd, digits=2),
                               z_score = round(z_score, digits=2))
  plts[[5]] <- tableGrob(elpd)
  h <- arrangeGrob(grobs = plts, widths=c(1, 1, 1),
                layout_matrix = rbind(c(1, 4, 4), c(5, 2, 2), c(3, 3, 3)))
  return(h)
}