Skip to content

Instantly share code, notes, and snippets.

@ericnovik
Last active July 31, 2019 20:58
Show Gist options
  • Save ericnovik/20e8f292c5ddb125c79bff4bc8e8f77a to your computer and use it in GitHub Desktop.
Save ericnovik/20e8f292c5ddb125c79bff4bc8e8f77a to your computer and use it in GitHub Desktop.
library(dplyr)
library(tidyr)
library(magrittr)
library(bayesplot)
plot_areas <- function(d, n, mu, s, cutoff, left = TRUE, h, v) {
if (left) {
Pr <- apply(d, 2, function (x) mean(x < cutoff))
} else {
Pr <- apply(d, 2, function (x) mean(x >= cutoff))
}
key <- colnames(d)
d <- gather(d)
d %<>%
group_by(key) %>%
mutate(x = density(value, bw = "SJ", n = n)$x,
y = density(value, bw = "SJ", n = n)$y,
x_pos = mean(x),
y_pos = mean(y))
if (left) {
d <- mutate(d, area = x < cutoff)
} else {
d <- mutate(d, area = x >= cutoff)
}
label <- paste("Pr =", as.character(round(Pr, 2)))
dat_text <- data.frame(label = label, key = key)
p <- ggplot(d, aes(x, y)) +
geom_ribbon(aes(ymin = 0, ymax = y, fill = area), alpha = 1/2) +
facet_wrap(~key) +
geom_vline(xintercept = cutoff, size = 0.2, linetype = 2) +
ylab("") +
geom_text(
data = dat_text,
mapping = aes(x = -Inf, y = -Inf, label = label),
hjust = h,
vjust = v
) +
bayesplot::theme_default() +
theme(legend.position="none",
axis.title.y = element_blank(),
axis.text.y = element_blank(),
axis.ticks.y = element_blank(),
axis.line.y = element_blank())
return(p)
}
n <- 1e3
mu <- c(0.69, 0.76, 0.81, 0.78)
s <- c(0.1, 0.12, 0.15, 0.17)
cutoff <- 0.6
d <- tibble(A = rnorm(n, mu[1], s[1]),
B = rnorm(n, mu[2], s[2]),
C = rnorm(n, mu[3], s[3]),
D = rnorm(n, mu[4], s[4]))
p <- plot_areas(d, n, mu, s, cutoff, TRUE, -0.5, -3) + xlab("Hazard Ratio")
p
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment