################################################### ### code chunk number 1: binom-s.rnw:15-23 ################################################### theta <- c(2,4,6,8)/10 prior <- c(1,1,1,1)/4 x <- 1 n <- 1 like <- dbinom( x, n, theta ) like.pr <- prior * like post <- like.pr / sum( like.pr ) round( cbind( theta, prior, like, like.pr, post ), 3 ) ################################################### ### code chunk number 2: binom-s.rnw:32-40 ################################################### theta <- c(2,4,6,8)/10 prior <- c(1,1,1,1)/4 x <- 15 n <- 20 like <- dbinom( x, n, theta ) like.pr <- prior * like post <- like.pr / sum( like.pr ) round( cbind( theta, prior, like, like.pr, post ), 3 ) ################################################### ### code chunk number 3: binom-s.rnw:44-45 ################################################### round( cbind( theta, prior, like, like.pr, post ), 17 ) ################################################### ### code chunk number 4: binom-s.rnw:60-68 ################################################### theta <- c(0,2,4,6,8,10)/10 prior <- c(1,1,1,1,1,1)/6 x <- 15 n <- 20 like <- dbinom( x, n, theta ) like.pr <- prior * like post <- like.pr / sum( like.pr ) round( cbind( theta, prior, like, like.pr, post ), 7 ) ################################################### ### code chunk number 5: binom-s.rnw:73-81 ################################################### theta <- c(0,2,4,6,8,10)/10 prior <- c(1,1,1,1,1,1)/6 x <- 1 n <- 1 like <- dbinom( x, n, theta ) like.pr <- prior * like post <- like.pr / sum( like.pr ) round( cbind( theta, prior, like, like.pr, post ), 3 ) ################################################### ### code chunk number 6: binom-s.rnw:104-110 ################################################### m <- 0.4 s <- 0.1 a.plus.b <- m*(1-m)/s^2 - 1 a <- m * a.plus.b b <- a.plus.b - a c(m,s,a,b) ################################################### ### code chunk number 7: prior1 ################################################### # Points where we plot: p <- seq(from=0,to=1,length=100) # Graph of the prior plot( p, dbeta( p, a, b ), lwd=4, bty="n", type="l" ) ################################################### ### code chunk number 8: lik1 ################################################### x <- 15 n <- 20 plot( p, dbinom( x, n, p ), lwd=4, bty="n", type="l" ) ################################################### ### code chunk number 9: post1 ################################################### # Graph of the likelihood plot( p, dbeta( p, a+x, b+n-x ), lwd=4, bty="n", type="l" ) ################################################### ### code chunk number 10: three ################################################### par( mfcol=c(3,1) ) plot( p, dbeta( p, a, b ), lwd=4, bty="n", type="l" ) plot( p, dbinom( x, n, p ), lwd=4, bty="n", type="l" ) plot( p, dbeta( p, a+x, b+n-x ), lwd=4, bty="n", type="l" ) ################################################### ### code chunk number 11: three-x ################################################### par( mfcol=c(3,1), mar=c(3,3,0,0) ) plot( p, dbeta( p, a, b ), lwd=4, bty="n", type="l" ) text( par("usr")[1], par("usr")[4], "\n Prior", adj=c(0,1) ) plot( p, dbinom( x, n, p ), lwd=4, bty="n", type="l" ) text( par("usr")[1], par("usr")[4], "\n Likelihood", adj=c(0,1) ) plot( p, dbeta( p, a+x, b+n-x ), lwd=4, bty="n", type="l" ) text( par("usr")[1], par("usr")[4], "\n Posterior", adj=c(0,1) ) ################################################### ### code chunk number 12: binom-s.rnw:174-192 ################################################### Bayes.ill <- function( m, s, x, n, ... ) { p <- seq(0,1,,1000) a.plus.b <- m*(1-m)/s^2 - 1 a <- m * a.plus.b b <- a.plus.b - a plot( p, dbeta( p, a, b ), lwd=4, bty="n", type="l", ... ) text( par("usr")[1], par("usr")[4], paste("\n Prior\n m=", m, ",s=", s, "\n a=", a,", b=", b), adj=c(0,1) ) plot( p, dbinom( x, n, p ), lwd=4, bty="n", type="l", ... ) text( par("usr")[1], par("usr")[4], paste("\n Likelihood\n n=", n,", x=",x), adj=c(0,1) ) plot( p, dbeta( p, a+x, b+n-x ), lwd=4, bty="n", type="l", ... ) text( par("usr")[1], par("usr")[4], paste("\n Posterior\n Beta(", a+x, ",", b+n-x, ")"), adj=c(0,1) ) } ################################################### ### code chunk number 13: twelve-a ################################################### par( mfcol=c(3,2), mar=c(2,4,0,0) ) Bayes.ill( 0.4, 0.2, 15, 20, col=gray(0.5) ) Bayes.ill( 0.4, 0.1, 15, 20, col=gray(0.5) ) ################################################### ### code chunk number 14: twelve-b ################################################### par( mfcol=c(3,2), mar=c(2,4,0,0) ) Bayes.ill( 0.4, 0.2, 55, 100, col=gray(0.5) ) Bayes.ill( 0.4, 0.1, 75, 100, col=gray(0.5) ) ################################################### ### code chunk number 15: binom-s.rnw:235-236 ################################################### pbeta( 0.6, 100, 100 ) - pbeta( 0.4, 100, 100 ) ################################################### ### code chunk number 16: binom-s.rnw:239-241 ################################################### zz <- rbeta( 10000, 100, 100 ) mean( zz<0.6 & zz>0.4 ) ################################################### ### code chunk number 17: births ################################################### a <- b <- 100 m <- a/(a+b) s <- sqrt(m*(1-m)/(a+b+1)) par( mfcol=c(3,1), mar=c(4,2,0,0) ) Bayes.ill( m, s, 511, 1000, xlim=c(0.4,0.6), xlab="% male births" ) abline(v=0.5) ################################################### ### code chunk number 18: binom-s.rnw:272-273 ################################################### pbeta(0.5,611,589) ################################################### ### code chunk number 19: binom-s.rnw:279-280 ################################################### pbeta(0.5,512,490)