#DataColada.org/74
#Written by Uri Simonsohn
#Last update: 2018 11 27
#
#Please contact Uri privately and directly if you see any problems: urisohn@gmail.com
#
###########################################

#Description:
#The expression of concern proposes that using one scale in the morning and
#another in the afernoon, with different precision (e.g., to 2 or to 5 grams) 
#may explain the prevalence of 3s and 7s. The code below shows that a uniform is expected
#anyway,
###############################################################################


#Function to extract last digit
  last.digit =function(x1) as.numeric(substr(x1,nchar(x1),nchar(x1)))

#Generate data
  n=500000
  x1=round(rnorm(n,mean=50,sd=5))       #Variable precise to 1 gram 
  x2=round(rnorm(n,mean=50,sd=5)/2)*2   #Variable precise to 5 gram 
  x5=round(rnorm(n,mean=100,sd=5)/5)*5  #Variable precise to 5 gram 
  x2.100=round(rnorm(n,mean=100,sd=5)/2)*2  #Variable precise to 5 gram 
#Extract last digits of differences
  ld15=last.digit(x5-x1)
  ld25=last.digit(x5-x2)
  ld12=last.digit(x2.100-x1)
  
#Printout to show mechanics
  cbind(x5,x2,x5-x2,ld25)[1:30,]

#Plot of the distributions
  par(mfrow=c(1,3))
    barplot(table(x1),ylab='',xlab='grams', col='dodgerblue',main='1 gm scale')
    barplot(table(x2),ylab='',xlab='grams', col='dodgerblue',main='2 gm scale')
    barplot(table(x5),ylab='',xlab='grams', col='orange',main='5 gm scale')
  
#Plot of the last digits
  par(mfrow=c(1,2))
    barplot(table(ld15),ylab='',xlab='last digit', col='purple1',main='5gm scale (M=100,SD=5) -\n 1 gm scale (M=50,SD=5)')
    barplot(table(ld25),ylab='',xlab='last digit', col='purple3',main='5gm scale(M=100,SD=5) -\n 2 gm scale (M=50,SD=5)')
    barplot(table(ld12),ylab='',xlab='last digit', col='purple4',main='2gm scale (M=100,SD=5) -\n 1 gm scale (M=50,SD=5)')