#R Code behind Figure 1 in
#Simonsohn, Simmons and Nelson, "P-Curve and Effect Size: Correcting for Publication Bias Using Only Significant Results"
#
#Prepared by Uri Simonsohn, uws@wharton.upenn.edu
#
##########################################################
#OVERVIEW
#This program generates Figure 1 in the paper where expected p-curves for different effect size/sample-size combinations are obtained
#It begins creating two functions
#Function 1: does all the action for the p-curve creation
#Function 2: is for % formatting in the axis of the figure
####################################################################
#Function 1 - Draw expected p-curve
drawp=function(n,d,col)
#Syntax:
#n: per-cell sample size for two-sample difference of means test
#d: standardized difference of means (Cohen's d)
#col: color of the line
{
#Degrees of freedom of the t-test
df=2*n-2
#noncentrality parameter for the student distribution
ncp=sqrt(n/2)*d
#Critical values,xc, for p=.05, .04, .03, .02 and ,01 (for 2-sided t-tests so take p/2 off the right side of the c.d.f)
x5=qt(.975,df=df);
x4=qt(.98, df=df);
x3=qt(.985,df=df);
x2=qt(.99, df=df);
x1=qt(.995,df=df)
#Probabilty of a p-value biggert han p=.04, .03, .02 and .01 given p<.05 and ncp=ncp
power=1-pt(x5,df=df,ncp=ncp) #what proportion of tests are of the predicted sign and p<.05
#The appendix in the main paper explains the calculations below
pp4=1/power*(pt(x4,df=df, ncp=ncp)-(1-power)) #p<.04
pp3=1/power*(pt(x3,df=df, ncp=ncp)-(1-power)) #p<.03
pp2=1/power*(pt(x2,df=df, ncp=ncp)-(1-power)) #p<.02
pp1=1/power*(pt(x1,df=df, ncp=ncp)-(1-power)) #p<.01
#within bins proportions (so substract the cdf up to the lower point)
prop5=pp4; #What proportion of significant p-values are .04