Lakens

library("meta")
library("metafor")
#EXPLANATION:
#This script is adapted from Simonsohn, Nelson, & Simmons, 2014 - the original (2- Fig 2 - File-drawering and effect size estimation.R) can be downloaded here: http://www.p-curve.com/Supplement/Rcode/   
#The original created simulations of many different meta-analyses - I stripped everything to look at a single meta-analysis

#SIMPLIFIED LOSS FUNCTION, USED THROUGHOUT BELOW
#see appendix in paper for explanations and more robust version (robust to user input with t<0 and/or p>.05)
loss=function(t_obs,df_obs,d_est) {  
  ncp_est=sqrt((df_obs+2)/4)*d_est

Lakens

library(MBESS)
library(pwr)

nSims <- 100000 #number of simulated experiments
p <-numeric(nSims) #set up empty container for all simulated p-values
obs_pwr <-numeric(nSims) #set up empty container
t <-numeric(nSims) #set up empty container
d_all<-numeric(nSims) 
N<-33 #number of participants

Lakens

require(car) #Car package required for Levene's test

n1<-38 #size condition x
n2<-22 #size condition y
sd1<-1.11 #sd condition x
sd2<-1.84 #sd condition y 
m1<-0
m2<-0
trueD<-(m2-m1)/(sqrt((((n1 - 1)*((sd1^2))) + (n2 - 1)*((sd2^2)))/((n1+n2)-2)))
trueD

Lakens

nSims <- 100000 #number of simulated experiments (the more, the more accurate the numbers you get, but the longer it takes. I used 1000000 simulations for my blog)
N<-32 #number of participants
lowp<-0.04
highp<-0.05

#set up some variables
p<-numeric(nSims) 
p2<-numeric(nSims) 
obs_pwr<-numeric(nSims) 
t<-numeric(nSims) 

Lakens

#The necessary packages (please ensure that these have been #installed before running the following code). 
#Many of the analyses are performed with both packages to double-check the calculations
library(meta)
library(metafor)

#The data for the full sample. All outputs that end in ".all" #are based on the full sample. 
r=0.5 #set correlation between dv's to 0.5, using formula's for variance from dependent test for one-sample t-test
#effect sizes
d.all<-c(-0.054458115, 0.093112835, 0.206774549, 0.478004184, 0.141213676, 0.272767619, 0.251, 0.195142683, 0.258913272, 0.204022676, 0.223, 0.145336391, 0.092630988, 0.192, 0.418607214, 0.290090298, 0.108461538, 0.316418861, 0.417283465, -0.043799784, 0.147113469, 0.139283883, 0.060796002, -0.113389342, 0, 0.114947393, 0.169535252, -0.048878416, -0.029138576, -0.004628788, 0.086060345, 0.231428571, 0.161220346, -0.050911688, 0.044376016, 0.185846777, 0.251, 0.52, 0.355321158, 0.02665009, 0.052177581, -0.011716899, -0.023772371, 0.252476027, 0.210898339, 0.067797354,

Lakens

set.seed(2)
options(scipen=20) #disable scientific notation for numbers

nSim<-10 #numbber of simulated studies

library(pwr)
library(MBESS)
library(gsDesign) # The group sequential design package
library(BayesFactor)

Lakens

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

Lakens

#Simulate bias in eta-squared, omega-squared, and epsilon-squared.
#R script by Kensuke Okada from: Okada, K. (2013). Is omega squared less biased? A comparison of three major effect size indices in one-way ANOVA. Behaviormetrika, 40(2), 129-147.

muvec <- c(0.00,0.00,0.8,0.8)
meanmu <- mean(muvec)
sigb <- sum((muvec-meanmu)^2)/4
eta2p <- sigb/(sigb+1)
k <- length(muvec)
nsim <- 1000000 #Bias should decrease as sample size increases. Set simulations to 1000000 for best results
njs <- c(10,20,30,40,50,60,70,80,90,100)

Lakens

# # # # # # # # # # #
#Initial settings----
# # # # # # # # # # #
if(!require(ggplot2)){install.packages('ggplot2')}
library(ggplot2)
if(!require(MBESS)){install.packages('MBESS')}
library(MBESS)
if(!require(pwr)){install.packages('pwr')}
library(pwr)
if(!require(meta)){install.packages('meta')}

Lakens

#Additional Analyses of Nuijten et al: https://mbnuijten.files.wordpress.com/2013/01/nuijtenetal_2015_reportingerrorspsychology1.pdf
#First run the original script to read in the data: https://osf.io/e9qbp/


#Select only errors.
subdata<-subset(data, data$Error == TRUE) 
subdata$pdif<-subdata$Reported.P.Value-subdata$Computed #Compute difference in p-values.
#Plot differences in reported and computed p-values for all errors
ggplot(as.data.frame(subdata$pdif), aes(subdata$pdif))  + 
  geom_histogram(colour="black", fill="grey", binwidth = 0.01) + ggtitle("All Errors") + xlab("Reported P-value minus Computed P-value") + ylab("Frequency")  + theme_bw(base_size=20)