Q1: Single Parameter Poisson

Consider the count of airline crashes per year over a 10-year period from 1976 to 1985: \(y = (24, 25, 31, 31, 22, 21, 26, 20, 16, 22)\).

(a) Show that the Gamma distribution is the conjugate prior

Solution

To prove Gamma distribution is conjugate to poisson ,we will let \[\boxed{\mathbf{y} \sim \text{Pois}(\lambda)} \], it follows through Bayes’ theorem that \[ f(\theta|\mathbf{y}) \propto \text{likelihood of } \mathbf{\theta} \times \text{prior } \theta \].

Therefore, the likelihood of \(\theta\) assuming \(y \sim \text{Pois}(\lambda)\) is as follows:

\[ L(\lambda|y) \sim \prod_{i=1}^{n} e^{-\lambda} \lambda^{y} \]

NB: We will remove constants and introduce the proportionality symbol \(\sim\).

\[ L(\lambda|y) \sim \prod_{i=1}^n e^{-\lambda} \lambda^{y} \]

Hence:

\[ L(\lambda | y) \sim e^{-n\lambda} \lambda^{n\bar{y}} \]

Assuming \(\lambda \sim \gamma(\alpha, \beta)\):

Then from Bayes:

\[ \text{Posterior}(\lambda|y_{observed}) \sim f(y_{observed}|\lambda)\times \text{prior}(\lambda)=L(\lambda|y_{observed}) \times \text{prior}(\lambda) \]

From above:

\[ \text{Posterior}(\lambda|y) \sim e^{-n\lambda} \lambda^{n\bar{y}} \times \gamma(\alpha, \beta) = (e^{-\beta\lambda}) \lambda^{\alpha-1} \]

\[ \text{Posterior}(\lambda|y) \sim e^{-n\lambda} \lambda^{n\bar{y}} \times e^{-\beta\lambda} \lambda^{\alpha-1} \]

The above can be expressed as:

\[ \text{Posterior}(\lambda | y) \sim e^{-\lambda(n+\beta)} \cdot \lambda^{n\bar{y} + \alpha - 1} \]

The above posterior is therefore a \(\gamma(n+\beta, n\bar{y} + \alpha)\). Hence, Gamma is conjugate for Poisson.

(b) Show graphics of the resulting posterior distribution for choices of various gamma priors.

Solution

From the data, the average rate for airline crashes is

y <- c(24, 25, 31, 31, 22, 21, 26, 20, 16, 22)
average_crash <- mean(y)
average_crash
## [1] 23.8
Aircrash distribution with Kernel Density Estimate

Aircrash distribution with Kernel Density Estimate 

(c) Show a table of the prior, MLE and posterior estimates of the poison mean (θ) under different choices of the gamma priors in (a) above.

The mean for \(\lambda \sim \Gamma(\alpha,\beta)\) is given by: \(\frac{\alpha}{\beta}\)

Prior MLE Posterior Mean
Non-informative Gamma(5,0.5) 23.8 \(5/0.5=10\)
Non-informative Gamma(11,0.5) 23.8 \(11/0.5=22\)
Informative Gamma(249,10.5) 23.8 \(249/10.5=23.7\)

MLE = \(\bar{y} = 23.8\)

Q2: CR and HPD (a) Differentiate between Credible Intervals and the Highest Posterior Density (HPD) in Bayesian analysis Credible Intervals (CI) and Highest Posterior Density (HPD) are both statistical measures used in Bayesian analysis to quantify uncertainty in parameter estimates, but they have distinct interpretations and characteristics.

  1. Credible Interval (CI):
    • A Credible Interval is a range of values for a parameter of interest.
    • It represents a probability mass, typically centered around the posterior median or mean.
    • The interval contains a specified probability (credibility level) of the posterior distribution.
    • Common choices for the credibility level include 95%, 90%, etc.
    • The width of the interval may vary depending on the shape of the posterior distribution.
  2. Highest Posterior Density (HPD):
    • The Highest Posterior Density interval is the narrowest interval that contains a specified probability mass.
    • It is defined by the region where the posterior density is maximized.
    • HPD intervals are useful for capturing the most plausible values for the parameter.
    • The width of the HPD interval is not fixed and depends on the shape of the posterior distribution.
    • Unlike CIs, HPD intervals do not necessarily correspond to a specific credibility level.

In summary, while both Credible Intervals and Highest Posterior Density intervals provide information about parameter uncertainty, CIs are based on probability mass, often centered around the median or mean, and have a fixed credibility level. On the other hand, HPD intervals focus on the region of maximum posterior density and do not adhere to a specific credibility level, providing a more localized measure of uncertainty.

  1. Graphically demonstrate this for a disease prevalence given that out of 150 individuals 18 had disease. Hint: assume Binomial distribution

From the above, we can see that \(\hat{p} = \frac{18}{150} = 0.12\).

Now assuming a \(\beta\) prior for p it follows that the posterior for p will be beta with parameters().

Using Bayes’ theorem:

\[ \text{Posterior}(p|y) \sim L(p|y) \times \text{prior}(p) \] \[L(p|y_{observed}) \sim \prod_{i=1}^{n} \binom{n}{y} p^{y} (1-p)^{n-y}\times y \] Taking out constants that dont involve p leaves us with the following expression: \[L(p|y_{observed}) \sim \prod_{i=1}^{n} \ p^{y} (1-p)^{n-y}\times y \]

\[L(p|y_{observed}) \sim \prod_{i=1}^{n} \ p^{y} (1-p)^{n-y}\times y \] The above then becomes

\[L(p|y_{\text{observed}}) \sim p^{\sum_{i=1}^n y} (1-p)^{n-\sum_{i=1}^n y} \times y\]

The conjugate prior for Binomial is the beta: So this implies that if \[p\sim \beta\space distribution \space with \space parameters (\alpha,\beta)\] Then \[f(p|y)\sim p^{\alpha -1}\times q^{\beta-1}\]

(p|y) = \[f(p|y) \sim p^{\alpha - 1} \times q^{\beta - 1} \cdot p^{\sum_{i=1}^n y} (1-p)^{n - \sum_{i=1}^n y} \]

The above posterior will then reduce to:

(p|y) = \[f(p|y) \sim p^{\alpha + \sum_{i=1}^n y} - 1 \times q^{\beta+n - \sum_{i=1}^y} - 1 \]

It follows that the posterior for P is \[\text{Beta}(\alpha + \sum_{i=1}^n y, \beta + n - \sum_{i=1}^ny)\]

Where P \[\sum_{i=1}^n y=n\times p=150*0.12=18 \]

Lets generate non -informative priors for p
>shape 1 <- 0.12
>shape 2 <- 3
>mean_beta_non_informative <- shape1 / (shape1 + 
shape2)

Now through inspection, the Credible intervals and HPD will be as follows:

Through visual inspection, the HPD is within \(x=15\) and \(x=21\). Hence, this implies \(p\) will fall between \(\frac{15}{180}=0.08\) or \(\frac{21}{180}=0.12\).

For Credible intervals: At a 95% confidence level, \(p\) will fall between \(0.12 \pm 1.96 \sqrt{0.12 \times 0.88/180}\). This will be \(0.12 \pm 0.024\)=[0.1,0.14]