The Binomial Distribution

The binomial formula gives the probability of k successes in n independent trials each with k probability of success:

\[ P(k | n,p) = {{n}\choose{k}} \cdot p^k (1-p)^{n-k} \]

where \({{n}\choose{k}} = \frac{n!}{k!(n-k)!} =\) “n choose k” = “the number of ways to choose k things from n”

dbinom

The dbinom function simply plugs in to the binomial formula. If you want to know the probability of exactly 5 successes in 10 trials if each trial has a 70% chance of success you can type:

dbinom(x=5, size=10, prob=0.7)

or simply:

dbinom(5, 10, 0.7)

Better yet, dbinom can do operations on vectors of numbers. If I wanted to know the probability of any number of successes, 0 through 10, in 10 trials each with a 70% chance of success we could do:

dbinom(0:10, 10, 0.7)

and we could even plot these probabilities versus the number of successes:

plot(x = 0:10, y=dbinom(0:10, 10, 0.7), type="h", xlab="Successes", ylab="Probability")

To check your understanding:

  1. Calculate the probability of exactly 10 successes in 30 independent trials, if there is a 30% chance of success in each trial.
  2. Calculate the probability of every number of successes, 0 through 30.
  3. Plot these probabilities.

A Simplified Basketball Problem

Imagine that there are only three types of basketball players when it comes to choosing corner threes, 20% of players have a “true” 30% success rate (meaning that if they kept shooting corner’s threes forever, they’d hit 30% of them), 60% of players have a “true” 40% success rate and 20% of players have a “true” 50% success rate.

We can write these probabilities and player qualities as vectors:

probs <- c(0.2, 0.6, 0.2)
qualities <- c(0.30, 0.40, 0.50)

and we could calculate the average player quality in this league by weighing the qualities by the probabilities:

Average Shooting Ability \(= \frac{\Sigma \ x \cdot P(x)}{\Sigma P(x)}\) [Note that because the probabilities add to 1, dividing by the sum of the probabilities isn’t necessary here BUT later when we’re dealing with relative probabilities which don’t add to 1, this will be important.]

sum(probs*qualities)/sum(probs)

Now, let’s imagine that we see a player hit 12 out of 20 corner threes. How likely is she to below to each of the three player quality groups?

Bayes’ Theorem says that we can get relative “posterior probabilities” by multiplying our “prior probabilities” (the chance that she belonged to each group before we saw any shooting) by the relative chance that members of each group would make 12 out of 20 corner threes.

We can use our friend dbinom to calculate the probability that 30%, 40% and 50% shooters would make exactly 12 out of 20 corner threes:

dbinom(12, 20, qualities)

Now let’s multiply these by the prior probabilities (which are just the percentages of players of each quality):

dbinom(12, 20, qualities)*probs

The line of code above gave us a vector of relative probabilities which clearly don’t add to one. To translate these into probabilities, we can just divide each relative probability by the sum of the relative probabilities:

rel.probs <- dbinom(12, 20, qualities)*probs

rel.probs/sum(rel.probs)

Now, we know how likely she is to belong to each “quality” group. If we’ve done everything correctly we should have found that she has roughly a 52% chance of being a 50% shooter and only a 17% chance of being a 30% shooter – two possibilities which were equally likely before we saw her shoot.

We can find her “average” shooting ability (really our best estimate or her true shooting ability) in the same way we found the average for all players, by weighing qualities by probabilities:

posterior.probs <- rel.probs/sum(rel.probs)

sum(posterior.probs*qualities)/sum(posterior.probs)

The Normal Distribution

The question posed in class, was more challenging and at least somewhat more realistic. Instead of only three possible shooting qualities, we said that players were, on average 40% shooters with a standard deviation of 5% (and shooting ability is normally distributed). We can use this to determine relative (“prior”) likelihoods of any shooting ability using the dnorm function.

The following code find the relative frequencies of 30%, 40% and 50% shooters given a normal distribution centered on 40% with a standard deviation of 5%:

dnorm(c(30, 40, 50), 40, 5)

We could also find the probability of any shooting ability 0% to 100% given our distribution that’s centered on 40% with a standard deviation of 5.

dnorm(0:100, 40, 5)

and we can plot these:

plot(0:100, dnorm(0:100, 40, 5), type="h", xlab="shooting ability", ylab="relative frequency", main="prior probabilities")

Now, let’s calculate the probability that shooters of each of these possible qualities would make exactly 12 of 20 shots:

dbinom(12, 20, (0:100)/100)

We can plot this too:

plot(0:100, dbinom(12, 20, (0:100)/100), type="h", xlab="shooting ability", ylab="likelihood of making 12 of 20", main="likelihoods")

As you may have guessed, we can now calculate our relative posterior probabilities by multiplying our priors by our likelihoods:

posteriors <- dnorm(0:100, 40, 5)*dbinom(12, 20, (0:100)/100)

Let’s plot these too!

plot(0:100, posteriors, type="h", xlab="shooting ability", ylab="relative posterior prob", main="posteriors")

I can get the average of our posterior, our best guess as the the chance that she makes her next shot just as we did before:

sum(posteriors*(0:100))/sum(posteriors)

We should find that she has a 43.4% chance of making her next shot.

Challenges:

  1. Some Saint Ann’s Students are more likely to arrive on time than others, and in secret teacher meetings (where grading is rampant) we typical group them as follows:

If a student arrives on time for 40 of the first 60 days of school, what is the probability that they will arrive on time for the 61st day of school?

  1. At a rival school, a median student is twice as likely to arrive on time as they are to be late but, there too, not all students are the same. While the typical student has \(log_2{odds} = 1\) of timeliness the standard deviation in \(log_2{odds}\) is also 1 (and \(log_2{odds}\) are normally distributed).

We can calculate, and plot, the relative probabilities that the \(log_2{odds}\) of a student’s timeliness has different possible values (between -4 and 4) by doing:

possible.logodds.values <- seq(-4,4, 0.1)

dnorm(possible.logodds.values, 1, 1)
plot(possible.logodds.values, dnorm(possible.logodds.values, 1, 1), type="h", xlab="log2 odds of timeliness", ylab="relative frequency")

We could also plot the relative chances of one of these students having different odds of timeliness as:

possible.odds.values <- 2^possible.logodds.values

plot(possible.odds.values, dnorm(possible.logodds.values, 1, 1), type="h", xlab="odds of timeliness", ylab="relative frequency")
possible.prob.values <- possible.odds.values/(possible.odds.values+1)

plot(possible.prob.values, dnorm(possible.logodds.values, 1, 1), type="h", xlab="probability of timeliness", ylab="relative frequency")

If a student at this rival school arrives on time for 40 of the first 60 days of school, what is the probability that they will arrive on time for the 61st day of school?