p-values within boundaries

p-values_within_boundaries.R

nSims <- 1000000 #number of simulated experiments
p <-numeric(nSims) #set up empty container for all simulated p-values
d <-numeric(nSims) #set up empty container for all simulated d's

n=20

for(i in 1:nSims){ #for each simulated experiment
  x<-rnorm(n = n, mean = 0.68, sd = 1) #produce n simulated participants
  #with mean=100 and SD=20
  y<-rnorm(n = n, mean = 0, sd = 1) #produce n simulated participants
  #with mean=100 and SD=20
  z<-t.test(x,y) #perform the t-test
  p[i]<-z$p.value #get the p-value and store it
  d[i] <-  (((mean(x)-mean(y))) / sqrt((((n - 1)*((sd(x)^2))) +  ((n - 1)*((sd(y)^2))))/((n * 2) -2)))  
}

#now plot the histogram
hist(p, main="Histogram of p-values", xlab=("Observed p-value"), breaks=20)
hist(d, main="Histogram of Cohen's d", xlab=("Observed Cohen's d"))

#effect size
mean(d)

#Observed Power
power<-p[p<0.05]
length(power)/nSims

#Count how often 2 p-values occur in a row
runs <- rle(p>0.025&p<0.05) #set boundaries
result<-with(runs, table(values, lengths)) 
#Return total number of 2 p-values in a row.
(result[2,2]+((result[2,3])*2)+((result[2,4]*3))+((result[2,5]*4))+((result[2,6]*5)))/nSims