False positives in a scientist's lifetime

falsepositives.R

## store proportion of false positives
## in one lifetime of 200 experiments:
prop_fps<-rep(NA,1000)
## run k=1000 scientists, each with
## a lifetime of 200 experiments:
for(k in 1:1000){
## number of experiments for each scientist:
nexp<-200
## prob of sampling from a population
## where the null is true:
prob<- 0.20
## true mean:
## something low to reflect realistic
## studies done with low power:
effect_size<- 5
## standard deviation:
s<-50

tests_null_true<-c(0)
tests_null_false<-c(0)

for(i in 1:nexp){
   null_true<-rbinom(n=1,size=1,p=prob)
   if(null_true){
       x<-rnorm(1000,mean=0,sd=s)
       tests_null_true[i]<-t.test(x)$p.value
     } else {
         x<-rnorm(1000,mean=effect_size,sd=s)
         tests_null_false[i]<-t.test(x)$p.value
       }
   }

## Compute proportion of false positives:

## remove NAs: 
tests_null_true<-tests_null_true[!is.na(tests_null_true)]
tests_null_false<-tests_null_false[!is.na(tests_null_false)]

## Proportion of all null-true tests that ended up 
## with a Type I error:
mean(tests_null_true<0.05)

## number of true positives: 
(num_false_pos<-table(tests_null_true<0.05)[2])

## proporation of all null-false tests 
## that correctly detected effect
## ("observed" power):
mean(tests_null_false<0.05)

## number of "true positives":
(num_true_pos<-table(tests_null_false<0.05)[2])

## prop. of false positive results 
## published in lifetime:

prop_fps[k]<-num_false_pos/(num_false_pos+num_true_pos)
}

hist(prop_fps,freq=FALSE,main="Distribution of
proportion of false positives")

summary(prop_fps)