This is the code for Alexander Etz's blog at http://wp.me/p4sgtg-SQ

MaxBF.R

#This is the code for Alexander Etz's blog at http://wp.me/p4sgtg-SQ

#Code to make the figure and find max oracle BF for normal data
maxL <- function(mean=1.96,se=1,h0=0){
        L1 <-dnorm(mean,mean,se)
        L2 <-dnorm(mean,h0,se)
        Ratio <- L1/L2
        curve(dnorm(x,mean,se), xlim = c(-2*mean,2.5*mean), ylab = "Likelihood", xlab = "Population mean", las=1,
              main = "Likelihood function for the mean", lwd = 3)
        points(mean, L1, cex = 2, pch = 21, bg = "cyan")
        points(h0, L2, cex = 2, pch = 21, bg = "cyan")
        lines(c(mean, h0), c(L1, L1), lwd = 3, lty = 2, col = "cyan")
        lines(c(h0, h0), c(L1, L2), lwd = 3, lty = 2, col = "cyan")
        return(Ratio) ## Returns the Bayes factor for oracle hypothesis vs null
}

#Code to find the max subtle oracle prior BF for normal data
t = 1.96                                #Critical value corresponding to the p value
maxBF = 1/(sqrt(exp(1))*t*exp(-t^2/2))  #Formula from Berger and Sellke 1987
#multiply by 2*(1-p/2) if using a one-sided prior