############################################### ## Author: Joshua M. Tebbs ## Date: 14 September 2021 ## Update: 12 October 2021 ## STAT 513 course notes: R Code ## Chapter 11 ############################################### # Figure 11.1 # Page 71 par(mar = c(4, 4.4, 4, 2)) # adjust plotting region to avoid cutting off labels y = seq(0,1,0.001) plot(y,dbeta(y,2,10),type="l",lty=1,xlab=expression(theta),ylab=expression(g(theta)),ylim=c(0,5),cex.lab=1.25) lines(y,dbeta(y,3,3),lty=4) lines(y,dbeta(y,5,1),lty=8) abline(h=0)s abline(v=1,lty=2) legend(0.4,3.5,lty=c(1,4,8),c(expression(paste(alpha,"=2, ",beta,"=10")),expression(paste(alpha,"=3, ",beta,"=3")), expression(paste(alpha,"=5, ",beta,"=1")))) # Figure 11.2 # Page 72 y = seq(0,1,0.001) plot(y,dbeta(y,2,10),type="l",lty=1,xlab=expression(theta),ylab="Prior versus posterior",ylim=c(0,9),cex.lab=1.25) lines(y,dbeta(y,42,70),lty=4) abline(h=0) legend(0.6,6,lty=c(1,4),c(expression(paste(alpha,"=2, ",beta,"=10")),expression(paste(alpha,"=42, ",beta,"=70")))) text(0.1,4.75,"Prior",cex=1.25) text(0.525,8,"Posterior",cex=1.25) # Figure 11.3 # Page 73 y = seq(0,1,0.001) plot(y,dbeta(y,42,70),type="l",lty=1,xlab=expression(theta),ylab="Posterior distribution",ylim=c(0,9),cex.lab=1.25) lines(y,dbeta(y,43,63),lty=4) lines(y,dbeta(y,45,61),lty=8) abline(h=0) legend(0.6,6,lty=c(1,4,8),c(expression(paste(alpha,"=42, ",beta,"=70")),expression(paste(alpha,"=43, ",beta,"=63")), expression(paste(alpha,"=45, ",beta,"=61")))) # Figure 11.4 # Page 76 par(mar = c(4, 4.4, 4, 2)) # adjust plotting region to avoid cutting off labels y = seq(0,5,0.01) plot(y,dgamma(y,1.25,1/0.5),type="l",lty=1,xlab=expression(theta),ylab=expression(g(theta)),xlim=c(0,4.5), ylim=c(0,1.2),cex.lab=1.25) lines(y,dgamma(y,1.25,1/0.75),lty=4) lines(y,dgamma(y,1.5,1),lty=8) abline(h=0) legend(1.75,0.75,lty=c(1,4,8),c(expression(paste(alpha,"=1.25, ",beta,"=0.5")), expression(paste(alpha,"=1.25, ",beta,"=0.75")),expression(paste(alpha,"=1.5, ",beta,"=1")))) # Figure 11.5 # Page 77 y = seq(0,5,0.01) plot(y,dgamma(y,1.5,1),type="l",lty=1,xlab=expression(theta),ylab="Prior versus posterior",ylim=c(0,3.2), xlim=c(0,4.5),cex.lab=1.25) lines(y,dgamma(y,104.5,85),lty=4) abline(h=0) legend(2.5,2,lty=c(1,4),c(expression(paste(alpha,"=1.5, ",beta,"=1")),expression(paste(alpha,"=104.5, ",beta,"=1/85")))) text(2.5,0.4,"Prior",cex=1.25) text(1.9,2.5,"Posterior",cex=1.25) # Figure 11.6 # Page 84 par(mar = c(4, 4.4, 4, 2)) # adjust plotting region to avoid cutting off labels # Left y = seq(0,1,0.001) plot(y,dbeta(y,1,1),type="l",lty=1,xlab=expression(theta),ylab=expression(g(theta)), xaxp=c(0,1,2),yaxp=c(0,1,1),ylim=c(0,2),cex.lab=1.5) abline(h=0) # Right y = seq(0,1,0.001) plot(y,dbeta(y,1/2,1/2),type="l",lty=1,xlab=expression(theta),ylab=expression(g(theta)), xaxp=c(0,1,2),yaxp=c(0,1,1),ylim=c(0,5),cex.lab=1.5) abline(h=0) # Figure 11.7 # Page 85 par(mar = c(4, 4.4, 4, 2)) # adjust plotting region to avoid cutting off labels y = seq(0,10,0.01) plot(y,dgamma(y,0.5,0.01),type="l",lty=1,xlab=expression(theta),ylab=expression(g(theta)),ylim=c(0,0.5),cex.lab=1.25) abline(h=0) # Figure 11.8 # Page 91 y = seq(0,1,0.001) plot(y,dbeta(y,9.5,28.5),type="l",lty=1,xlab=expression(p),ylab="Posterior distribution",xlim=c(0,0.6), ylim=c(0,6),cex.lab=1.25) abline(h=0) points(x=9.5/(9.5+28.5),y=0,pch=19,cex=1.25) # post.mean points(x=qbeta(0.5,9.5,28.5),y=0,pch=19,cex=1.25) # post.med points(x=8.5/36,y=0,pch=19,cex=1.25) # post.mode # Figure 11.9 # Page 94 y = seq(0,1,0.001) plot(y,dbeta(y,9.5,28.5),type="l",lty=1,xlab=expression(p),ylab="Posterior distribution",xlim=c(0,0.6), ylim=c(0,6),cex.lab=1.25) abline(h=0) points(x=qbeta(0.025,9.5,28.5),y=0,pch=19,cex=1.25) # lower quantile points(x=qbeta(0.975,9.5,28.5),y=0,pch=19,cex=1.25) # upper quantile # Figure 11.10 # Page 99 y = seq(0.5,2,0.001) plot(y,dgamma(y,104.5,85),type="l",lty=1,xlab=expression(lambda),ylab="Posterior distribution",xlim=c(0,2), ylim=c(0,3.5),cex.lab=1.25) abline(h=0) x = seq(0.5,1,0.001) y = dgamma(x,104.5,85) polygon(c(0.5,x,1),c(0,y,0),col="lightblue") points(x=1,y=0,pch=19,cex=1.25)