Uncertainty in a fractional count

Question

What is the uncertainty (68% confidence level) of $N/M$, where $N$ is the number of entries that pass a cut and $M$ is the total number of entries? ($N$ and $M$ are both integers, and I'm interested in the extreme where $N$ or $M - N$ is a small integer, maybe zero.)

In the past, I've always assumed a Binomial model, where $N$ is the number of coin tosses that come up heads and $M$ is the total number of coin tosses. Following this logic, I've used the variance $Mp(1-p)$ to conclude that the uncertainty is $\sqrt{\frac{p(1-p)}{M}}$ (with $p=N/M$). However, I'm beginning to think this is flawed: this "uncertainty" is exactly zero if $p=0$ or $p=1$. Getting five heads in a row shouldn't lead one to conclude that the coin will always yield heads with perfect certainty.

In general, I think the upper uncertainty will be different from the lower uncertainty; that it should come from some integration that has a hard cut-off at $p=0$ and $p=1$, introducing an asymmetry close to the border. Should this come from a Bayesian formalism because I'm making inferences about an unknown distribution from measurement?

(I'm also surprised that I haven't found an answer online: I would have thought it to be a very common problem. Putting uncertainties on trigger efficiencies in physics, for instance.)

Answer 1

Apparently, there are many answers to this question: it has its own Wikipedia page and R package. The uncertainty range I described above is the "normal approximation interval":

$\displaystyle p + z \sqrt{\frac{p(1 - p)}{M}}$

where $z$ is signed, $z=0$ is the central value ($p = N/M$), $z=1$ is "one sigma" (68% confidence level) above the central value and $z=-1$ is "one sigma" below the central value. It has the failures I discussed, and they stem from the fact that the distribution of true values around an observation are not themselves binomial.

There are many ways to correct the problem; this blog shows an overlay plot of the different estimators (drawn from the R package). Each has different properties and a different justification.

The simplest correction (if you're not running R) is the Wilson score interval:

$\displaystyle \frac{1}{1+z^2/M} \left(p + \frac{z^2}{2M} + z \sqrt{\frac{p(1-p)}{M} + \frac{z^2}{4M^2}} \right)$

Note that $z \to 0$ in the above yields $p$, so the central value is unchanged, but the positive and negative errors ($z = 1$ and $z = -1$) are now asymmetric and not equal to each other as $p \to 0$ or $1$.

I would take this as a simple formula to use for putting error bars on plots.

Answer 2

I'd like to point out that the "exact" solution to this problem in the frequentist sense is the Clopper-Pearson interval. As the author noted, there are a number of nice approximate intervals with different characteristics, and the Wilson score interval is often a good choice. However I thought I would also point out that Bayesian inference with a uniform prior on $p$ yields a very simple formula for the posterior distribution for $p$

$$P(p|N,M) = \frac{(M+1)!}{N!(M-N)!} p^N (1-p)^{M-N}$$

with the following features:

$${\rm mode}(p|N,M) = \frac{N}{M}$$ $${\rm mean}(p|N,M) = \frac{N+1}{M+2}$$ $${\rm variance}(p|N,M)= \frac{(N+1) (N+2) }{(M+2)(M+3) } -\frac{(N+1)^2}{(M+2)^2} $$

Using [mean $\pm \sqrt{\rm variance}$] gives easy-to-compute 1$\sigma$ confidence intervals with very nice properties even when $N=0$ or $N=M$. One can even show that these intervals have descent frequentist coverage properties, especially at larger $M$, and only undercover when the "true" value of $p$ is unmeasurably close to 0 or 1.