library(dplyr)
library(reshape2)
library(ggplot2)

n <- 500
null <- .75
obs_n <- seq(1,n-1)
binom_probs_low <- pbinom(obs_n, n, null)
binom_probs_upp <- pbinom(obs_n, n, null,lower.tail = FALSE)
binom_probs <- pmin(binom_probs_low, binom_probs_upp)

t_probs <- rep(0,n-2)
for (i in 1:(n-1)) {
  x = c(rep(1,i),rep(0,n-i))
  t_probs[i] <- t.test(x,mu=.75)$p.value
}

df <- data.frame(Sample = obs_n, reject = NA, Binomial = binom_probs, t = t_probs)
df <- melt(df,id.vars = c("Sample","reject"),var = "Distribution")
df$reject[df$Distribution =='Binomial'] <- ifelse(df$value[df$Distribution =='Binomial'] < .025 |df$value[df$Distribution =='Binomial'] > .975 , 1,0)
df$reject[df$Distribution =='t'] <- ifelse(df$value[df$Distribution =='t'] < .025, 1,0)
df <- df %>% group_by(Sample) %>%
  mutate(decision = if (all(reject==1)) { 1} else if (all(reject==0) ) { 0 } else { 2 })
df$decision <- factor(df$decision)
ggplot(data = df, aes(x=Sample, y = value, color = decision,shape=Distribution)) +
  geom_point() +
  coord_cartesian(xlim = c(.6*n, .9*n)) +
  ggtitle("Binomial vs. t NHST Results")

# decision = if (all(reject==1)) { 1} else if (all(reject==0) ) { 0} else { 2}