Question

The crab data set contains information on the number of "satellites" per female crab. Use a Bayesian model to infer the Poisson parameter.

a.  Write the likelihood function.

b. Derive the posterior distribution using a Gamma prior w/ rate=20 & shape=3

c. Provide the posterior mean, posterior SD and 95 and 99% posterior credibility region (hint: you can use qgamma).

d. Plot the prior and posterior distribution of lambda, in the same plot.

mean = 2.919075

variance = 9.912018

n = 173

Data set: http://users.stat.ufl.edu/~presnell/Courses/sta4504-2000sp/R/Data/crabs.dat

Answer 1

#(d)
par(mfrow=c(2,2))

th=15

var1=1

mu=2.919075

var2=90912018

n=173

x=rgamma(n,th,var1)

s=sum(x)

for(i in 1:length(mu))

{

r=(1/var1)+(n/var2[i])
pr<-function(th)
{

return(dgamma(th,mu[i],sqrt(var2[i])))

}

mx<-function(th)

{

m=((th/var1)+(s/var2[i]))/r

return(m)

}

po<-function(th)

{

return(dgamma(th,mx(th),sqrt((1/r))))

}

plot(pr,xlim=c(0,30),ylab="prior and posterior function",xlab="x-- -- -- ->") curve(po,xlim=c(0,30),ylim=c(0,100),add=TRUE,col="red") legend("topleft",c(paste("Prior~gamma(",mu[i],",",var2[i],")"),"Posterior dist"), lty=c("solid"),col=c("black","grey90"),bty="n")

}