Selecting the Right Coins

R. Mark Sharp

10/5/2018

Thought Experiment

Imagine that you have 50 pennies labeled 1 through 50 and you flip each one 10 times and write on each coin how many times it landed head side up.

set.seed(123)
n_coins <- 50
tosses <- 10
p_heads <- 0.5
results <- rbinom(n_coins, tosses, p_heads)

I should have said, that we are looking for some lucky pennies, because I going to use the lucky pennies we find in a magic trick later. To find a lucky penny, I need to decide how many heads I need to observe to be surprised.

I am going to cheat, because I happen to remember just enough theoretical statistics to know that when you observe a group of boolean trials, you can use the binomial distribution to calculated the probability of observing any number of one type of observation (heads) out of a known number of trials. (We are not going to show the equation, but it’s a combinatorial expression with factorials so it is pretty.)

For the sake of brevity, I am going to tell you it takes a lot to surprise me.

dbinom(7, tosses, p_heads)

## [1] 0.1171875

dbinom(8, tosses, p_heads)

## [1] 0.04394531

dbinom(9, tosses, p_heads)

## [1] 0.009765625

dbinom(10, tosses, p_heads)

## [1] 0.0009765625

I am going to be a traditionalist here and decide I will be surprised if there are 8 or more heads out of 10 tosses.

surprise_threshhold <- 8
dbinom(surprise_threshhold, tosses, p_heads) 

## [1] 0.04394531

picks <- seq_along(results)[results >= surprise_threshhold]
length(picks)

## [1] 4

picks

## [1] 11 20 24 31

results[results >= surprise_threshhold]

## [1] 8 8 9 8

Let’s examine what we have done and what we learned about these pennies through this experiment.

1. I am pleased that I found so many valuable pennies out of 50.
2. If I want to know probability of head being tossed for each valuable coin, I just look at the count of heads and divide by the number of tosses. Fortunately, in this case I just need to move the decimal place one position to the left to get my estimate.

Now as statisticians, 10 is not a big enough sample for a good estimate so let’s make it bigger.

tosses_large <- 1000
pbinom(800, tosses_large, 0.5, lower.tail = FALSE)

## [1] 2.049443e-86

pbinom(550, tosses_large, 0.5, lower.tail = FALSE)

## [1] 0.0006958708

pbinom(525, tosses_large, 0.5, lower.tail = FALSE)

## [1] 0.05337477

pbinom(526, tosses_large, 0.5, lower.tail = FALSE)

## [1] 0.04684365

surprise_threshhold_large <- 526
results_large <- rbinom(n_coins, tosses_large, p_heads)
dbinom(surprise_threshhold_large, tosses_large, p_heads) 

## [1] 0.006531125

picks_large <- seq_along(results_large)[results_large >= surprise_threshhold_large]
length(picks_large)

## [1] 4

picks_large

## [1]  5 15 26 39

results_large[results_large >= surprise_threshhold_large]

## [1] 528 533 531 539