Given N samples drawn from a Bernoulli distribution with parameter $p$, so that $P(X=1)=p$ and $P(X=0)=1-p$, what method is most suitable to provide an estimate of the underlying parameter $p$ together with a measure of the uncertainty in this estimate?
I am aware that the proportion $\frac{x}{N}$ of outcomes for which $X=1$ is the maximum likelihood estimate for $p$, but this does not provide any information about how certain this estimate $\hat{p}$ is to coincide with the true $p$. If I take 5 samples and all 5 of them happen to have the outcome $X=1$, for example, this does not mean that $p$ has to be 1. Sure, $p=1$ is the one for which the probability of obtaining these observed sample values is highest, but other values of $p$ can also have the same result. Can someone provide a detailed discussion or point to references with a rigorous treatment of what methods can be used to quantify the uncertainty in this estimate in one sense or another? In particular, I am interested in the general case where $p$ can be very close to 0 or 1.