Bayesian Confidence Interval Estimate of Theta

Question

Let $X_1,...,X_n$ be a random sample from a Bernoulli distribution with success probability $p=\theta$. Assume a prior distribution on theta; $\theta\sim UNIFORM(0, 1)$. Derive a $100(1-\alpha)\%$ Bayesian interval estimate of $\theta$.

Finding the Bayes estimator itself is not hard; after you do the calculation, you'll find that it's just simply the expected value of a $BETA(1+\sum x_i, 1 + n - \sum x_i)$ distribution, which is just $\frac{1+\sum x_i}{2+n}$.

My issue is just figuring out how to get the confidence interval from this information; I don't understand exactly how to do this, and the wording seems unclear to me. Can somebody please show me what I need to do to find the interval in question? Thank you.

Answer 1

The likelihood is $$ L(t) = \text{constant} \times t^x (1-t)^{n-x} $$ where $x$ is the observed sum of the Bernoulli random variables. The prior probability measure is $dt$ for $0<t<1.$ Thus the posterior is $$ \text{constant} \times t^x(1-t)^{n-x} \, dt \quad \text{for } 0<t<1 $$ where this "constant" may be (and in this case is) different from the other one. So this is $\operatorname{Beta}(x+1, n-x+1).$ The measure is $$ \frac{\Gamma(n+2)}{\Gamma(x+1) \Gamma(n-x+1)} t^x (1-t)^{n-x} \, dt = \frac{(n+1)!}{x!(n-x)!} t^x (1-t)^{n-x} \, dt \quad \text{for } 0<t<1. $$

Suppose you want probability $1-\alpha/2$ in each tail. You need $$ \int_0^c \frac{(n+1)!}{x!(n-x)!} t^x(1-t)^{n-x} \, dt = \frac \alpha 2. $$ Note that

• This is a polynomial function of $c.$
• Therefore you may need numerical methods to find the value of $c,$ given $\alpha.$
• But standard software packages will give you values of the inverse of the c.d.f. of the Beta distribution.

And of course you also have $$ \int_d^1 \frac{(n+1)!}{x!(n-x)!} t^x(1-t)^{n-x} \, dt = \frac \alpha 2 $$ and you similarly seek $d.$

Answer 2

Specifically, assume you have prior $\mathsf{Unif}(0,1) \equiv \mathsf{Beta}(1,1)$ as prior, and data are 20 successes in 50 trials. Then as you say the posterior dixtribution is $\mathsf{Beta}(1+20, 1+30) = \mathsf{Beta}(21, 31)$ with point estimate $21/52.$ Using R or similar software you can find a 95% Bayesian probability interval estimate $(0.276, 0.539).$ [In R, the beta quantile function (or inverse CDF) is denoted qbeta with appropriate parameters.]

qbeta(c(.025,.975), 21, 31)
[1] 0.275843 0.538859

This is a specific computed instance of @MichaelHardy's Answer (+1).

Note: This interval estimate is numerically similar to the Agresti-Coull interval using $n^+ = n + 4$ and $\tilde p = (X+2)/n^+.$ The interval is $$\tilde p \pm 1.96\sqrt{\frac{\tilde p(1-\tilde p)}{n^+}}.$$ For the example above, it computes to $(0.2763530, 0.5384619).$ [Of course, numerical similarity aside, the interpretations if the Bayesian and frequentist intervals differ.]