Confidence intervals of the binomial distribution

conf_ints.py

import numpy as np
import scipy.stats
from scipy.interpolate import interp1d
#import matplotlib as mpl
from matplotlib import pyplot as plt


def cints(n, levels):    
    
    ci_out = np.zeros((len(levels), 1 + n))

    k = np.arange(n + 2) - 1
    
    for m in range(1, n):
        p = m / n   # in python2 use float(m) / n
        cdf = scipy.stats.binom.cdf(k,n,p)
        f = interp1d(cdf, k)
        ci_out[:, m] = f(levels)

    ci_out[:, n] = n    
    return ci_out


n = 1320
levels = [0.05, 0.95]
ci = cints(n, levels)

plt.plot(ci[0,:])
plt.plot(ci[1,:])
plt.show()

• probability = k / N where k is the horizontal axis value plotted (ranging from 0 to N)
• number of samples = N

Just linear interpolation used here when interpolating CDF onto required confidence intervals - maybe something fancier could be used