Unexpected behavior with qbinom

Question

Consider the following case:

qbinom(p=0.05,size=6040:6045,prob=0.995)
[1] 6001 6041 6042 6043 6004 6005

the quantile for n=6041,6042,6043 are 6041,6042 and 6043, respectively. This does not seems right for p=0.05.

Let's verify the result for n=6041 with pbinom:

pbinom(q=6000:6041,size=6041,prob=0.995)
 [1] 0.03490172 0.04974803 0.06943756 0.09489324 0.12695452 0.16626629
 [7] 0.21315764 0.26752726 0.32875643 0.39567139 0.46657229 0.53933683
[13] 0.61159304 0.68094122 0.74519263 0.80258630 0.85194713 0.89275986
[19] 0.92514962 0.94977961 0.96769159 0.98012377 0.98834033 0.99349835
[25] 0.99656542 0.99828757 0.99919751 0.99964817 0.99985646 0.99994584
[31] 0.99998123 0.99999408 0.99999832 0.99999957 0.99999991 0.99999998
[37] 1.00000000 1.00000000 1.00000000 1.00000000 1.00000000 1.00000000

6002 is the first quantile that is above 0.05. This contradicts with the results in qbinom.

Why is this happening?

Answer 1

The binomial distribution has two parameters, size (number of trials) and probability of success for each trial. You are holding the prob constant between the two calls, but not the size parameter.

The number of trials is different between the two functions. In the call to qbinom you asking for the quantile at given probability for trial sizes between 6040-6045. In the call to pbinom you are asking for different quantiles for a fixed trial size, 6041.

What you want to do to check is

purrr::map2_dbl(
  qbinom(p=0.05,size=6040:6045,prob=0.995),
  6040:6045, 
  \(q, size) pbinom(q, size, prob = 0.995)
)
#> [1] 0.06931045 1.00000000 1.00000000 1.00000000 0.05004380 0.05014268

Theses probabilities are all above 5%. The high probability of success of each trial makes this a VERY skewed distribution. There is a very narrow range of outcomes with any appreciable probability of occurring.

Answer 2

This is fixed in R 4.4.1.

The Bug Fixes section in Changes in 4.4.1 mentions:

qbinom() and potentially qpois(), qnbinom(), no longer sometimes fail accurate inversion (of pbinom(), etc), thanks to Christopher Chang's report and patch in PR#18711.

The detail is here and the fix is here. The issue is qbinom() jumped too far up (to 6041, 6042 and 6043 in the OP). The fix is to backtrack to the correct solution.