Information From Randomness?

info_puzzle.r

# https://www.quantamagazine.org/20150722-solution-information-from-randomness/

n = 1000
vals = runif(n) # values shown by opened hand
g = runif(n) # g values
highest = sample(c(TRUE,FALSE),n,replace=TRUE) # was the value in the opened hand the higher value of the two hands?
# one of two things happens
# 1) g is less than the value you chose, you stay, and you are right
# 2) g is greater than the value you chose, you switched, and you were right
# other cases are failures (add 0 to s) so don't care about them since they are recoverable via n-s
s = sum(ifelse((g < vals & highest) | (g > vals & !highest),1,0))

# plot curve of beta distr
f = n-s
mean = s/(s + f)
var = (s*f)/((s+f)^2*(s+f+1))
sd = sqrt(var)
curve(dbeta(x,s,f),xlim=(mean+c(-1,1)*4*sd))


###############
## generate a and b instead of randomly choosing value shown in hand
## gives 66% answer (if you set g = rep(0.5,n))
n = 1000
a = runif(n)
b = runif(n,a,1)
g = runif(n)
highest = sample(c(TRUE,FALSE),n,replace=TRUE)

s = 0
for(i in 1:n){
  high = highest[i]
  if(!high) shown = a[i]
  if(high) shown = b[i]
  if(g[i] < shown & high) s = s + 1 # you stayed and you were right
  if(g[i] > shown & !high) s = s + 1 # you switched and you were right 
}

f = n-s
print(qbeta(c(0.025,0.975),s,f))

## A more optimized approach for R
n = 1000
g = runif(n)
highest = sample(c(TRUE,FALSE),n,replace=TRUE)
a = runif(n)
b = runif(n,a,1)
v = ifelse(highest,b,a)

s = sum(ifelse((g < v & highest) | (g > v & !highest),1,0))
f = n-s
print(qbeta(c(0.025,0.975),s,f))