I expect identical results calculating a CI for a proportion using qbinom
and binom.test
, but they are actually slightly different:
success <- 360
n <- 1226
lci <- qbinom(0.025, n, success/n)/n
uci <- qbinom(0.975, n, success/n)/n
c(lci, uci)
[1] 0.2683524 0.3189233
binom.test(success, n, success/n)$conf.int
[1] 0.2682571 0.3200123
What am I missing here?
The source code for the binom.test
function uses qbeta
rather than qbinom
, as this is the accepted formula for exact binomial confidence intervals. Where x
is the number of successes, the lci and uci as given by binom.test
are:
p.L <- function(x, alpha) {
if (x == 0)
0
else qbeta(alpha, x, n - x + 1)
}
p.U <- function(x, alpha) {
if (x == n)
1
else qbeta(1 - alpha, x + 1, n - x)
}
alpha <- (1 - 0.95)/2
lci <- p.L(x, alpha)
uci <- p.U(x, alpha)