Problem Setup

We have two competing models:

Model M₀:
\(Y_i \mid \theta_0 \stackrel{i.i.d.}{\sim} \text{Geometric}(\theta_0)\) with PMF
\[ f_0(y \mid \theta_0) = \theta_0 (1 - \theta_0)^y, \quad y = 0,1,2,\dots \]

Model M₁:
\(Y_i \mid \theta_1 \stackrel{i.i.d.}{\sim} \text{Poisson}(\theta_1)\) with PMF
\[ f_1(y \mid \theta_1) = e^{-\theta_1} \frac{\theta_1^y}{y!}, \quad y = 0,1,2,\dots \]

Prior Distributions

  • For M₀: \(\theta_0 \sim \text{Beta}(\alpha_0, \beta_0)\)
    \[ \pi_0(\theta_0) = \frac{\theta_0^{\alpha_0-1} (1-\theta_0)^{\beta_0-1}}{B(\alpha_0, \beta_0)} \]

  • For M₁: \(\theta_1 \sim \text{Gamma}(\alpha_1, \beta_1)\) (rate parameterization)
    \[ \pi_1(\theta_1) = \frac{\beta_1^{\alpha_1}}{\Gamma(\alpha_1)} \theta_1^{\alpha_1-1} e^{-\beta_1 \theta_1} \]

Part 1: Prior Predictive Means

Under M₀ (Geometric)

By the law of total expectation: \[ E[Y \mid M_0] = E_{\theta_0 \sim \pi_0}\left[ E(Y \mid \theta_0) \right] = E_{\theta_0}\left[ \frac{1-\theta_0}{\theta_0} \right] \] \[ = \int_0^1 \frac{1-\theta_0}{\theta_0} \cdot \frac{\theta_0^{\alpha_0-1}(1-\theta_0)^{\beta_0-1}}{B(\alpha_0,\beta_0)} d\theta_0 \] \[ = \frac{1}{B(\alpha_0,\beta_0)} \int_0^1 \theta_0^{\alpha_0-2} (1-\theta_0)^{\beta_0} d\theta_0 \] \[ = \frac{B(\alpha_0-1, \beta_0+1)}{B(\alpha_0,\beta_0)} = \frac{\beta_0}{\alpha_0-1}, \quad \alpha_0 > 1 \]

Under M₁ (Poisson)

Similarly, \[ E[Y \mid M_1] = E_{\theta_1 \sim \pi_1}[\theta_1] = \frac{\alpha_1}{\beta_1} \]

Thus: \[ \boxed{E[Y \mid M_0] = \frac{\beta_0}{\alpha_0-1}, \quad E[Y \mid M_1] = \frac{\alpha_1}{\beta_1}} \]

Part 2: Bayes Factor \(BF_{01}\)

The Bayes factor comparing M₀ to M₁ is: \[ BF_{01} = \frac{m_0(\mathbf{y})}{m_1(\mathbf{y})} \] where \(m_0\) and \(m_1\) are the marginal likelihoods.

Marginal Likelihood for M₀ (Beta-Geometric)

\[ m_0(\mathbf{y}) = \int_0^1 \prod_{i=1}^n \theta_0 (1-\theta_0)^{y_i} \cdot \frac{\theta_0^{\alpha_0-1}(1-\theta_0)^{\beta_0-1}}{B(\alpha_0,\beta_0)} d\theta_0 \] \[ = \frac{1}{B(\alpha_0,\beta_0)} \int_0^1 \theta_0^{n+\alpha_0-1} (1-\theta_0)^{\sum y_i + \beta_0 - 1} d\theta_0 \] \[ = \frac{B(\alpha_0 + n,\; \beta_0 + S)}{B(\alpha_0,\beta_0)}, \quad S = \sum_{i=1}^n y_i \]

Marginal Likelihood for M₁ (Gamma-Poisson)

\[ m_1(\mathbf{y}) = \int_0^\infty \prod_{i=1}^n e^{-\theta_1} \frac{\theta_1^{y_i}}{y_i!} \cdot \frac{\beta_1^{\alpha_1}}{\Gamma(\alpha_1)} \theta_1^{\alpha_1-1} e^{-\beta_1 \theta_1} d\theta_1 \] \[ = \frac{\beta_1^{\alpha_1}}{\prod y_i! \; \Gamma(\alpha_1)} \int_0^\infty \theta_1^{\alpha_1 + S - 1} e^{-(n+\beta_1)\theta_1} d\theta_1 \] \[ = \frac{\beta_1^{\alpha_1}}{\prod y_i! \; \Gamma(\alpha_1)} \cdot \frac{\Gamma(\alpha_1 + S)}{(n+\beta_1)^{\alpha_1 + S}} \]

Bayes Factor

\[ \boxed{BF_{01} = \frac{B(\alpha_0+n,\; \beta_0+S)}{B(\alpha_0,\beta_0)} \cdot \frac{\prod y_i! \; \Gamma(\alpha_1) \; (n+\beta_1)^{\alpha_1+S}}{\beta_1^{\alpha_1} \; \Gamma(\alpha_1+S)}} \]

Part 3: Numerical Examples

We consider: - \(n = 2\) - Dataset 1: \(y_1 = y_2 = 0\) (\(S=0\)) - Dataset 2: \(y_1 = y_2 = 2\) (\(S=4\))

Prior sets:

Set \(\alpha_0\) \(\beta_0\) \(\alpha_1\) \(\beta_1\)
A 1 2 2 1
B 30 60 60 30

Numerical Values

Data Sets

  • Data 1: \(n = 2\), \(y_1 = y_2 = 0\), \(S = 0\)
  • Data 2: \(n = 2\), \(y_1 = y_2 = 2\), \(S = 4\)

Prior Sets

  • Prior set A: \(\alpha_0 = 1,\; \beta_0 = 2,\; \alpha_1 = 2,\; \beta_1 = 1\)
  • Prior set B: \(\alpha_0 = 30,\; \beta_0 = 60,\; \alpha_1 = 60,\; \beta_1 = 30\)

Detailed Calculations for Data 1 (\(S = 0\))

Prior Set A

For \(M_0\): - Beta(1,2): \(B(1,2) = \frac{\Gamma(1)\Gamma(2)}{\Gamma(3)} = \frac{1 \cdot 1}{2} = 0.5\) - \(B(\alpha_0 + n,\; \beta_0 + S) = B(3,2) = \frac{\Gamma(3)\Gamma(2)}{\Gamma(5)} = \frac{2 \cdot 1}{24} = \frac{1}{12}\) - So \(m_0 = \frac{1/12}{0.5} = \frac{1}{6}\)

For \(M_1\): - \(\alpha_1 = 2,\; \beta_1 = 1,\; S = 0\) - \(\Gamma(\alpha_1) = \Gamma(2) = 1,\; \Gamma(\alpha_1 + S) = \Gamma(2) = 1\) - \(\prod y_i! = 0! \times 0! = 1\) - \(m_1 = \frac{1^2}{1 \cdot 1} \cdot \frac{1}{(1+2)^2} = \frac{1}{9}\)

Bayes Factor: \[ BF_{01} = \frac{1/6}{1/9} = \frac{9}{6} = 1.5 \]

Prior Set B

For \(M_0\): \[ m_0 = \frac{B(32,60)}{B(30,60)} = \frac{\Gamma(32)\Gamma(60)/\Gamma(92)}{\Gamma(30)\Gamma(60)/\Gamma(90)} = \frac{\Gamma(32)}{\Gamma(30)} \cdot \frac{\Gamma(90)}{\Gamma(92)} \] \[ \frac{\Gamma(32)}{\Gamma(30)} = 31 \times 30 = 930 \] \[ \frac{\Gamma(90)}{\Gamma(92)} = \frac{1}{91 \times 92} = \frac{1}{8372} \] \[ m_0 = \frac{930}{8372} \approx 0.111085 \]

For \(M_1\): \[ m_1 = \left(\frac{30}{32}\right)^{60} = \left(\frac{15}{16}\right)^{60} \] \[ \ln(15/16) = \ln(15) - \ln(16) = 2.70805 - 2.7725887 = -0.0645385 \] \[ 60 \times (-0.0645385) = -3.87231 \] \[ m_1 = e^{-3.87231} \approx 0.020867 \]

Bayes Factor: \[ BF_{01} = \frac{0.111085}{0.020867} \approx 5.32 \]

Detailed Calculations for Data 2 (\(S = 4\))

Prior Set A

For \(M_0\): - \(B(1+2,\; 2+4) = B(3,6) = \frac{\Gamma(3)\Gamma(6)}{\Gamma(9)} = \frac{2 \cdot 120}{40320} = \frac{240}{40320} = \frac{1}{168}\) - \(B(1,2) = 0.5\) - \(m_0 = \frac{1/168}{0.5} = \frac{1}{84} \approx 0.0119048\)

For \(M_1\): - \(\prod y_i! = 2! \times 2! = 4\) - \(\Gamma(2+4) = 120,\; \Gamma(2) = 1\) - \(m_1 = \frac{1^2}{4 \cdot 1} \cdot \frac{120}{(1+2)^6} = \frac{30}{729} = \frac{10}{243} \approx 0.041152\)

Bayes Factor: \[ BF_{01} = \frac{1/84}{10/243} = \frac{243}{840} = \frac{81}{280} \approx 0.2893 \]

Prior Set B

For \(M_0\): \[ m_0 = \frac{B(32,64)}{B(30,60)} = \frac{\Gamma(32)\Gamma(64)/\Gamma(96)}{\Gamma(30)\Gamma(60)/\Gamma(90)} = \frac{\Gamma(32)}{\Gamma(30)} \cdot \frac{\Gamma(64)}{\Gamma(60)} \cdot \frac{\Gamma(90)}{\Gamma(96)} \] \[ \frac{\Gamma(32)}{\Gamma(30)} = 31 \times 30 = 930 \] \[ \frac{\Gamma(64)}{\Gamma(60)} = 63 \times 62 \times 61 \times 60 = 14,\!295,\!960 \] \[ \frac{\Gamma(90)}{\Gamma(96)} = \frac{1}{95 \times 94 \times 93 \times 92 \times 91 \times 90} = \frac{1}{625,\!757,\!605,\!200} \] \[ m_0 = 930 \times 14,\!295,\!960 \times \frac{1}{625,\!757,\!605,\!200} \approx 0.021245 \]

For \(M_1\): \[ m_1 = \frac{1}{4} \cdot \frac{\Gamma(64)}{\Gamma(60)} \cdot \left(\frac{30}{32}\right)^{60} \cdot \frac{1}{32^4} \] \[ \frac{\Gamma(64)}{\Gamma(60)} = 14,\!295,\!960,\quad \left(\frac{30}{32}\right)^{60} \approx 0.020867,\quad 32^4 = 1,\!048,\!576 \] \[ m_1 = \frac{1}{4} \times 14,\!295,\!960 \times 0.020867 \times \frac{1}{1,\!048,\!576} \approx 0.07113 \]

Bayes Factor: \[ BF_{01} = \frac{0.021245}{0.07113} \approx 0.2986 \]

Summary of Results

Dataset Prior Set A Prior Set B
(0,0) 1.5 5.32
(2,2) 0.289 0.299

R Program

# Function to compute Bayes Factor BF01
bayes_factor <- function(y, alpha0, beta0, alpha1, beta1) {
  n <- length(y)
  S <- sum(y)
  
  # Marginal likelihood M0 (Beta-Geometric)
  m0 <- exp(lbeta(alpha0 + n, beta0 + S) - lbeta(alpha0, beta0))
  
  # Marginal likelihood M1 (Gamma-Poisson)
  log_m1 <- alpha1 * log(beta1) - sum(lgamma(y + 1)) - lgamma(alpha1) + 
            lgamma(alpha1 + S) - (alpha1 + S) * log(beta1 + n)
  m1 <- exp(log_m1)
  
  return(m0 / m1)  # BF01
}

# Data sets
y1 <- c(0, 0)
y2 <- c(2, 2)

# Prior set A
bfA_1 <- bayes_factor(y1, alpha0 = 1, beta0 = 2, alpha1 = 2, beta1 = 1)
bfA_2 <- bayes_factor(y2, alpha0 = 1, beta0 = 2, alpha1 = 2, beta1 = 1)

# Prior set B
bfB_1 <- bayes_factor(y1, alpha0 = 30, beta0 = 60, alpha1 = 60, beta1 = 30)
bfB_2 <- bayes_factor(y2, alpha0 = 30, beta0 = 60, alpha1 = 60, beta1 = 30)

# Results
cat("Prior Set A (α0=1, β0=2; α1=2, β1=1):\n")
## Prior Set A (α0=1, β0=2; α1=2, β1=1):
cat("  y = (0,0): BF01 =", round(bfA_1, 4), "\n")
##   y = (0,0): BF01 = 1.5
cat("  y = (2,2): BF01 =", round(bfA_2, 4), "\n\n")
##   y = (2,2): BF01 = 0.2893
cat("Prior Set B (α0=30, β0=60; α1=60, β1=30):\n")
## Prior Set B (α0=30, β0=60; α1=60, β1=30):
cat("  y = (0,0): BF01 =", round(bfB_1, 4), "\n")
##   y = (0,0): BF01 = 5.4566
cat("  y = (2,2): BF01 =", round(bfB_2, 4), "\n")
##   y = (2,2): BF01 = 0.2995