Alt

Step 1

Let \(X_i\) represent the number of customers who visit the webpage in a given hour.

Assume \(X_1, \ X_2, \ ..., \ X_n \overset{iid}\sim Poisson(\lambda)\). Where the prior distribution of \(\lambda \overset{iid}\sim Exponential(\beta)\).

Now we have:

\(p(x_i|\lambda)=\frac{\lambda^{\sum x_i}e^{-n\lambda}}{\prod {x!}}\) and

\(p(\lambda|x_i)=\frac{\left(\frac{\lambda^{\sum x_i}e^{-n\lambda}}{\prod {x!}}\right)\left(\frac{1}{\beta}e^{-\frac{\lambda}{\beta}}\right)}{\int_0^\infty \left(\left(\frac{\lambda^{\sum x_i}e^{-n\lambda}}{\prod {x!}}\right)\left(\frac{1}{\beta}e^{-\frac{\lambda}{\beta}}\right)\right)d\lambda}\)

The denominator simplifies to

\(\int_0^\infty \left(\left(\frac{\lambda^{\sum x_i}e^{-n\lambda}}{\prod {x!}}\right)\left(\frac{1}{\beta}e^{-\frac{\lambda}{\beta}}\right)\right)d\lambda=\frac{1}{\beta\prod x!}\int_0^\infty \lambda^{\sum x_i }e^{-\left(\frac{n\beta+1}{\beta}\right)\lambda}d\lambda\)

\(\frac{1}{\beta \prod x!}\int_0^\infty \lambda^{\sum x_i}e^{-\left(\frac{n\beta+1}{\beta}\right)\lambda}d\lambda=\frac{\Gamma\left(\sum x_i+1\right)}{\beta \prod x!\left(\frac{\beta+1}{\beta}\right)^{\sum x_i+1}}\int_0^\infty \frac{\left(\frac{\beta+1}{\beta}\right)^{\sum x_i+1}}{\Gamma\left(\sum x_i+1\right)}\lambda^{\sum x_i}e^{-\left(\frac{n\beta+1}{\beta}\right)\lambda}d\lambda\)

\(\frac{1}{\beta \prod x!}\int_0^\infty \lambda^{\sum x_i}e^{-\left(\frac{n\beta+1}{\beta}\right)\lambda}d\lambda=\frac{\Gamma\left(\sum x_i+1\right)}{\beta \prod x!\left(\frac{\beta+1}{\beta}\right)^{\sum x_i+1}}\int_0^\infty Gamma\left(\left(\sum x_i+1\right),\frac{n\beta+1}{\beta}\right)d\lambda\)

which simplifies to

\(\int_0^\infty \left(\left(\frac{\lambda^{\sum x_i}e^{-n\lambda}}{\prod {x!}}\right)\left(\frac{1}{\beta}e^{-\frac{\lambda}{\beta}}\right)\right)d\lambda=\frac{\Gamma\left(\sum x_i+1\right)}{\beta \prod x!\left(\frac{\beta+1}{\beta}\right)^{\sum x_i+1}}\)

Now,

\(p(\lambda|x)=\frac{\left(\frac{\lambda^{\sum x_i}e^{-n\lambda}}{\prod {x!}}\right)\left(\frac{1}{\beta}e^{-\frac{\lambda}{\beta}}\right)}{\frac{\Gamma\left(\sum x_i+1\right)}{\beta \prod x!\left(\frac{\beta+1}{\beta}\right)^{\sum x_i+1}}}\)

\(p(\lambda|x)=\frac{\left(\frac{n\beta+1}{\beta}\right)^{\sum x_i+1}}{\Gamma\left(\sum x_i+1\right)}\lambda^{\sum x_i}e^{-\left(\frac{n\beta+1}{\beta}\right)\lambda}\)

But \(\bar{x}=\frac{\sum x_i}{n}\) implying \(\sum x_i=n\bar{x}\)

Therefore, the posterior distribution is

\(p(\lambda|x)\sim gamma\left(n\bar{x}+1,\ \frac{n\beta+1}{\beta}\right)\)

Step 2

Figure 1 shows plot of prior and three posteriors under \(\bar{x}=0.30, 0.45, 0.49, \beta = 0.05\) and \(n=100\).

n <- 100
Beta <- 0.05
xbar1 <- 0.3
xbar2 <- 0.45
xbar3 <- 0.49
prior <- function(x){
  exp(-x/Beta)/Beta
}
posterior1 <- function(x){
  ((((n*Beta+1)/Beta)^(n*xbar1+1))/gamma(n*xbar1+1))*x^(n*xbar1)*exp(-((n*Beta+1)/Beta)*x)
}
posterior2 <- function(x){
  ((((n*Beta+1)/Beta)^(n*xbar2+1))/gamma(n*xbar2+1))*x^(n*xbar2)*exp(-((n*Beta+1)/Beta)*x)
}
posterior3 <- function(x){
  ((((n*Beta+1)/Beta)^(n*xbar3+1))/gamma(n*xbar3+1))*x^(n*xbar3)*exp(-((n*Beta+1)/Beta)*x)
}
curve(prior, from = 0, to = 1, col = 1, ylim = c(0,9),
      main = 
        "Figure 1: Plot of prior and posterior (n = 100, beta = 0.05)",
      xlab = "lambda",
      ylab = "Prior/Posterior")
curve(posterior1, from = 0, to = 1, col = 2, 
      add = TRUE, ylim = c(0,9))
curve(posterior2, from = 0, to = 1, col = 3, 
      add = TRUE, ylim = c(0,9))
curve(posterior3, from = 0, to = 1, col = 4, 
      add = TRUE, ylim = c(0,9))
legend("topright", c("Prior", 
                     "Posterior mean = 0.30", 
                     "Posterior mean = 0.45",
                     "Posterior mean = 0.49"),
       col = c(1,2,3,4),
       lty = c(1,1), 
       lwd = c(1,1))

Figure 2 shows plot of prior and three posterior plots under \(\bar{x}=0.030, 0.045, 0.049, \beta = 0.05\) and \(n=1000\).

n <- 1000
Beta <- 0.05
xbar1 <- 0.030
xbar2 <- 0.045
xbar3 <- 0.049
prior <- function(x){
  exp(-x/Beta)/Beta
}
posterior1 <- function(x){
  ((((n*Beta+1)/Beta)^(n*xbar1+1))/gamma(n*xbar1+1))*x^(n*xbar1)*exp(-((n*Beta+1)/Beta)*x)
}
posterior2 <- function(x){
  ((((n*Beta+1)/Beta)^(n*xbar2+1))/gamma(n*xbar2+1))*x^(n*xbar2)*exp(-((n*Beta+1)/Beta)*x)
}
posterior3 <- function(x){
  ((((n*Beta+1)/Beta)^(n*xbar3+1))/gamma(n*xbar3+1))*x^(n*xbar3)*exp(-((n*Beta+1)/Beta)*x)
}
curve(prior, from = 0, to = 0.08, col = 1, ylim = c(0,90),
      main = 
        "Figure 2: Plot of prior and posterior (n = 1000, beta = 0.05)",
      xlab = "lambda",
      ylab = "Prior/Posterior")
curve(posterior1, from = 0, to = 0.08, col = 2, 
      add = TRUE, ylim = c(0,100))
curve(posterior2, from = 0, to = 0.08, col = 3, 
      add = TRUE, ylim = c(0,100))
curve(posterior3, from = 0, to = 0.08, col = 4, 
      add = TRUE, ylim = c(0,100))
legend("topright", c("Prior", 
                     "Posterior mean = 0.30", 
                     "Posterior mean = 0.45",
                     "Posterior mean = 0.49"),
       col = c(1,2,3,4),
       lty = c(1,1), 
       lwd = c(1,1))

As the sample size increases the prior becomes more flat. However, as the beta approaches the sample mean, the distribution becomes more symetric. mean.

References

  1. Bibliography Ashutosh Nandeshwar. 2019. “How to Plot Multiple Lines on the Same Graph Using R.” YouTube Video. YouTube. https://www.youtube.com/watch?v=_gl-CNA7wck.

  2. Data Mentor. 2017. “R Multiple Plot Using Par() Function.” DataMentor. November 23, 2017. https://www.datamentor.io/r-programming/subplot/.Gelman,

  3. Andrew, Abdel El-Shaarawi, and Walter Piegorsch. 2002. “Posterior Distribution.” http://www.stat.columbia.edu/~gelman/research/published/p032-_o.pdf.