Confidence interval using binofit

Hello all,

I was doing a test using binofit to calculate the confidence interval for a binomial distribution. I found that the confidence interval corresponds to the desired confidence level only when p is not very small. When p is small the obtained confidence level is higher than desired.

For a 95% I tried:

N=1000;
g=zeros(N,1);
Nb = 1000;
p = 1/5;
for n=1:N
    x=binornd(Nb,p);
    [phat,pci]=binofit(x,Nb);
    g(n)= (p>pci(1)) && (p<pci(2));
end
sum(g)/N

ans =

      0.9570

However, if p is decreased by a factor of 100:

N=1000;
g=zeros(N,1);
Nb = 1000;
p = 1/500;
for n=1:N
    x=binornd(Nb,p);
    [phat,pci]=binofit(x,Nb);
    g(n)= (p>pci(1)) && (p<pci(2));
end
sum(g)/N

ans =

      0.9900

Any idea of why this happens? Is there a better way to calculate the confidence interval in this case?

Thanks in advance

Answer 1

[Tom Lane]

The binofit function computes confidence intervals using the Clopper-Pearson method.

I believe what you're seeing is because the distribution is very coarse for N*p values on the order of 1. With N=1000 and p=1/500,

>> binocdf(4:5,1000,1/500)
ans =
    0.9475    0.9835

So in order for the confidence interval to have at least a 95% chance of including p=1/500, the intervals based on x=4 and x=5 would both have to include p. So the coverage probability will be over 98%.

Try increasing to N=100000 in your example and see if the results match.

Answer 2

[Javier]

Thanks!