Help in Bayesian Bernoulli-Beta Model (solution verification)

Question

I have some doubts that -for many of you- might seem basic and I am basically looking for some guidance verifying my attempt of solution.

The Problem

Let $p(x\mid\theta)=Bernoulli(\theta)$, and let's suppose the Prior Distribution is given by $p(\theta)=Beta(\theta \mid \alpha,\beta)$. ($\theta$ is unknown, and $\alpha, \beta$ are known).

Assume independet observations and obtain: 1) The Posterior Distribution of $\theta$, 2)The Prior Predictive Distribution, and 3)The Posterior Predictive Distribution.

Attempt of Solution

Let $I=\{1,2,...,n\}$, $\mathrm{X}=\{X_i\}_{i \in I}$ the independent observatrions above mentioned, and $\mathbb{I}$ be the Indicator function.

---

1)

Given the information of the excercise,

$p(\theta \mid \mathrm{x})=\frac{[\prod_{i \in I}\theta^{x_i}(1-\theta^{1-{x_i}})\mathbb{I}_{\{0,1\}}(x_i)][\theta^{\alpha-1}(1-\theta)^{1-\beta} \mathbb{I}_{[0,1]}(\theta)]}{\int_\Theta [\prod_{i \in I}\theta^{x_i}(1-\theta^{1-{x_i}})\mathbb{I}_{\{0,1\}}(x_i)][\theta^{\alpha-1}(1-\theta)^{1-\beta} \mathbb{I}_{[0,1]}(\theta)]d\theta}$

$p(\theta \mid \mathrm{x})=\frac{\theta^{\alpha + r-1}(1-\theta)^{\beta+n-r-1}\mathbb{I}_{\{0,1\}}(x_i)\mathbb{I}_{[0,1]}(\theta)}{\int_0^1\theta^{\alpha + r-1}(1-\theta)^{\beta+n-r-1}\mathbb{I}_{\{0,1\}}(x_i)d\theta}$

$p(\theta \mid \mathrm{x})=Beta(\theta \mid \alpha+r,\beta+n-r),$ $r:=\sum_{i \in I}x_i$.

---

2) The Prior Predictive Distribution $p(x)$ is given then by

$p(x)=\frac{1}{Beta(\theta \mid \alpha,\beta)}\int_0^1\theta^{\alpha+x-1}(1-\theta)^{\beta-x}\mathbb{I}_{\{0,1\}}(x)d\theta$

$p(x)=\frac{1}{Beta(\theta \mid \alpha,\beta)}\int_0^1(\theta^{\alpha+x-1}-\theta^{\alpha+\beta-1})\mathbb{I}_{\{0,1\}}(x)d\theta$

$p(x)=\frac{1}{Beta(\theta \mid \alpha,\beta)}(\frac{1}{\alpha+x}-\frac{1}{\alpha+\beta})\mathbb{I}_{\{0,1\}}(x)$

---

3) The Posterior Predictive Distribution $p(x \mid \mathrm{x})$ is given then by

$p(x \mid \mathrm{x})=\frac{1}{Beta(\theta \mid \alpha+r,\beta+n-r)}\int_0^1\theta^{\alpha+x+r-1}(1-\theta)^{\beta+n-x-r}\mathbb{I}_{\{0,1\}}(x)d\theta$

$p(x \mid \mathrm{x})=\frac{1}{Beta(\theta \mid \alpha+r,\beta+n-r)}\int_0^1(\theta^{\alpha+x+r-1}-\theta^{\alpha+\beta+n-1})\mathbb{I}_{\{0,1\}}(x)d\theta$

$p(x \mid \mathrm{x})=\frac{1}{Beta(\theta \mid \alpha+r,\beta+n-r)}(\frac{1}{\alpha+x+r}-\frac{1}{\alpha+\beta+n})\mathbb{I}_{\{0,1\}}(x)$

---

Answer 1

I have not verified every line of your Question, but it seems your initial formulation for each part is on the right track. In this particular situation, you may be able to save some work by noticing that the beta prior and the binomial likelihood are 'conjugate', as discussed below. [I will not comment on the posterior predictive distribution because you have not defined it or explained the meaning of $r.]$

Prior. As a specific example, suppose you are trying to interpret a public opinion poll to judge which of two candidates for mayor of a large city is currently ahead. Many American elections result in something near a 50:50 vote count. So it might be reasonable to use the prior distribution $\mathsf{Beta}(5, 5),$ which puts nearly half of its probability between 40% and 60% of the vote. [Computation in R.]

diff(pbeta(c(.4,.6), 5, 5))
[1] 0.4668646

If $\theta$ is the percentage in favor of Candidate A, then we might write this as $p(\theta) = K\theta^{5-1}(1-\theta)^{5-1},$ where $K$ is the norming constant (constant of integration). Sometimes we write this as $p(\theta) \propto \theta^{5-1}(1-\theta)^{5-1},$ where the symbol $\propto$ (read 'proportional to') acknowledges than the factor $K$ is absent.

Likelihood. Then suppose results from a reliable pollster based on interviewing $n = 700$, reports $x = 372$ in favor of Candidate A. Then the binomial likelihood function for the poll would be $p(x|\theta)\propto \theta^{372}(1-\theta)^{700-372}.$

Posterior. Then by Bayes' Theorem we have that the posterior distribution $p(\theta|x) = p(\theta) \times p(x|\theta).$ In our case this gives $$p(\theta|x) = \theta^{5-1}(1-\theta)^{5-1} \times \theta^{372}(1-\theta)^{700-372} \propto \theta^{377-1}(1-\theta)^{333-1},$$ where we recognize the right-hand side as proportional to the density function of $\mathsf{Beta}(377, 333).$

Conjugacy. Because the beta and binomial distributions are 'conjugate' (mathematically compatible), we do not need to do further computation to find the norming constant for the posterior distribution.

Posterior probability interval. A 95% posterior probability interval for the population proportion in favor of Candidate A would be $(0.494, 0.568).$ Disappointingly, this interval includes 0.5, so the candidate should maintain his/her campaign, but chances of success seem good.

qbeta(c(.025, .975), 377, 333) 
[1] 0.4942327 0.5675738

Note: For the poll discussed above, a frequentist Agresti-Coull 95% confidence interval is $(0.4944, 0.5681),$ numerically about the same as the Bayesian posterior (or 'credible') interval, but subject to somewhat different interpretations.

Addendum per comments on shortest posterior probability interval.

The length of the standard (probability-symmetric) interval above, which puts probability $0.025$ in each tail, is:

std.ci = qbeta(c(.025,.975), a1, a2); diff(std.ci)
[1] 0.07334116

Your posterior distribution with shape parameters 377 and 333, which have a ratio of about 1, is very nearly symmetrical, so the shortest interval will not be much different from the standard interval:

a1 = 377; a2 = 333; cred = .95 
ltp = seq(.0001, 1-cred, by=.00001) # vector of lower tail probs
ll = qbeta(ltp, a1,a2); ul=qbeta(cred+ltp, a1,a2)
len = ul-ll
ltp.short = ltp[len==min(len)];  ltp.short
[1] 0.02518    # not much different from 0.025
short.ci = qbeta(c(ltp.short, cred+ltp.short),a1,a2)
diff(short.ci)
[1] 0.07334081 # not much shorter than standard int

As @Xi'an comments, this shortest interval is called the HPD interval, for 'highest probability density', for reasons covered in various comments.

Answer 2

The answer to (2) is incorrect since $$\theta^{\alpha+x-1}(1-\theta)^{\beta-x}\ne\theta^{\alpha+x-1}-\theta^{\alpha+\beta-1}$$ Similarly, the answer to (3) is incorrect due to the same mistake.