# Probability that a 40% shooter is better than a 30% shooter given 20 shots?

I'm working through the Regression and Other Stories (ROS) textbook and it's not clear to me how Gelman etal. would like a reader to answer question 4.3 (which I summarize below):

Two basketball shooters. Shooter A has an accuracy of 30% while shooter B has an accuracy of 40%. They each take 20 shots. What is the probability that shooter B makes more shots?

Now, I'm aware of the following question and it all makes sense. What confuses me however, is that at no point in the textbook (unless I missed it) does Gelman etal. demonstrate that approach. Therefore, I wonder if there is another way to answer the question using a "Gelman" method?

So far I have the following R code:

N = 20
p_a = 0.3  # shooter A has a 30% accuracy
se_a = sqrt(p_a * (1 - p_a) / N)
p_b = 0.4  # shooter B has a 40% accuracy
se_b = sqrt(p_b * (1 - p_b) / N)

# The probability that shooter B makes more shots than B with mean and s.e. of:
p_diff = p_b - p_a
se_diff = sqrt(se_a**2 + se_b**2)


My intuition tells me that the difference in their ability is described by a (Gaussian?) distribution with mean p_diff (0.1) and standard error of se_diff (0.15). Given that p_diff = p_b - p_a, positive values of p_diff account for instances when B "beats" A while negative values are the opposite. See below:

Therefore, the probability of B beating A is the area under the graph to the right of the red line. Which can be calculated in r with:

# p_diff is 0.1
# se_diff is 0.15
p_b_beats_a = 1 - pnorm(0, p_diff, se_diff)


Which gives the incorrect value of ~0.75 (see link above for correct value).

What step/intuition have I missed? Given that at no point in ROS has Gelman etal. mentioned continuity correction (as the above link does), I'm assuming there is some other "Gelman" method/intuition that I've missed ... Alternatively, how does one go from the method as described in the link to the "Gelman" approach?

I'd appreciate both the maths/code but also an intuitive explanation.

• This isn't a question about statistics or inference: it's about computing Binomial probabilities. The thread you reference explains this and shows how to approximate these probabilities with Normal distributions.
– whuber
Commented Nov 8, 2022 at 13:28

Since the referenced thread only includes a direct approach in a comment, here goes.

You can directly evaluate the probability of the discrete events when shooter A makes more shots than shooter B and sum them together rather than approximating a normal distribution. If $$A \sim Bin(20, 0.3)$$ and $$B \sim Bin(20, 0.4)$$ then you are looking for $$P(A < B)$$, which is the same as $$\Sigma_{n=0}^{20} P(A < n, B = n)$$

We can probably assume the shots are independent, so
$$P(A

p <- 0
for (n in 0:20) {
p <- p + dbinom(n, 20, 0.4) * pbinom(n-1, 20, 0.3)
}


If you want to approximate the Binomial distribution by a Normal one and $$A \sim Bin(n,p)$$, you need $$np$$ to be large: it might be large enough here, but it, but given that you have to do the whole continuity correction thing, it might be more hassle than its worth!

A quick search shows evaluating the Binomial distribution comes up in ROS in Ex. 3.5, but there doesn't appear to be a solid section on probability theory.

I hope this demonstrates another approach. I made the mistake of getting A and B the wrong way around initially, hopefully I haven't missed any less than/equals signs!

• +1. It's hard to see how the "continuity correction" would be a hassle: it's just a matter of characterizing the event $A\lt B$ as $A-B\le-1/2.$ The R calculation is pnorm(-1/2, 20*(0.3 - 0.4), sqrt(20*0.3*(1 - 0.3) + 20*0.4*(1 - 0.4))) and it agrees with your calculation almost to the third decimal place. You don't need $np$ to be large here, because (a) you are subtracting one skewed variable from another, yielding an almost symmetrical difference; and (b) you are not computing a probability way out in one tail. Even with $20$ replaced by $2$ the Normal approximation is excellent(!).
– whuber
Commented Nov 8, 2022 at 17:30
• @Geoffrey Liddell: Thank you for that explanation. I think I'm almost there with regards to the "intuition" however, can how you went from $\Sigma_{n=0}^{20} P(A < n, B = n)$ to $P(A<B) = \Sigma_{n=0}^{20} P(A \le n-1)P(B = n)$? Commented Nov 10, 2022 at 14:09
• @whuber: Your comment is, I think, what I had in mind. However, in "my maths" my threshold is 0 and mu/se are 0.1/0.15 whereas your threshold is -0.5 and mu/se are -2/3. Is there a way to get the correct answer using "my maths"? (Specifically, by NOT multiplying in the N). P.S. I agree your threshold of -0.5 is correct. How do I convert that to "my maths"? i.e. pnorm(???, (0.3-0.4), sqrt(0.3*(1-0.3) + 0.4*(1-0.4)) Commented Nov 10, 2022 at 14:18
• You might be overthinking this: $-1/2$ for the continuity correction will work well. Concerning your first comment: the inequality $A\lt n$ is equivalent to the inequality $A\le n-1$ when $A$ is integral. Then use the fact that $A$ and $B$ are independent.
– whuber
Commented Nov 10, 2022 at 17:28
• @user1182556 I think I might have spotted the problem: you have put down se_a = sqrt(p_a * (1 - p_a) / N), but I think you might mean sqrt(p_a * (1- p_a) * N) (check at the bottom of page 43 of ROS). The missing piece of understanding here is that standard errors and standard deviations are very different things! Here we are interested in approximating the population variance, so we want standard deviation, whereas the standard error comes in when you are looking at the sampling distribution of the mean Commented Nov 10, 2022 at 23:22

All credit must go to @Geoffrey Liddell and @whuber for their help. I post the below for completeness.

As pointed out by @Geoffrey Liddell, I was confusing standard deviations and standard errors. The solution I was looking for (as guided by @whuber) looks something like:

N = 20
p_a = 0.3  # Shooter A has a 30% shooting accuracy
mu_a = N * p_a
std_a = sqrt(N * p_a * (1 - p_a))
p_b = 0.4  # Shooter B has a 40% shooting accuracy
mu_b = N * p_b
std_b = sqrt(N * p_b * (1 - p_b))

# determine coefficients of the Gaussian that describes the difference between
# B and A:
mu_diff = mu_b - mu_a
std_diff = sqrt(std_a**2 + std_b**2)

x = seq(-20, 20, 0.001)
y = dnorm(x, mu_diff, std_diff)
threshold = 0.5  # Positive threshold because of "B - A"
{
plot(x, y, 'l', )
abline(v=threshold, col='red')
text(-10, 0.06, 'A wins')
text(10, 0.06, 'B wins')
}

p_b_beats_a = 1 - pnorm(threshold, mu_diff, std_diff)



Where p_b_beats_a is ~0.69 which is the correct answer.

Thanks again to @Geoffrey Liddell and @whuber.

ADDENDUM: In conjunction with confusing standard deviations vs standard errors, I think another issue that was confusing me was working in "proportion space" vs "actual value space". For example, originally I had p_diff = p_b - p_a. That meant I had to convert my threshold of 0.5 in "actual value space" to a value in "proportion space" (which I still didn't know how to do - assuming it's even possible?). However, by working with actual values i.e mu_diff = mu_b - mu_a, the threshold of 0.5 is now intuitive and the maths follows.