CS计算机代考程序代写 rm(list=ls()); set.seed(20210127);

rm(list=ls()); set.seed(20210127);

# Data

n1 <- 10; y1 <- 5 # Data from experiment 1 n2 <- 10; y2 <- 8 # Data from experiment 2 S <- 1000 # Number of MC samples theta1 <- rbeta(S, y1+1, n1-y1+1) theta2 <- rbeta(S, y2+1, n2-y2+1) # Compare histograms for posteriors of theta1 and theta2 op <- par(mfrow=c(1,2)) hist(theta1, freq=F, right=F, xlab=expression(theta[1]), xlim=0:1, breaks=25, main="MC approximate posterior") hist(theta2,, freq=F, right=F, xlab=expression(theta[2]), xlim=0:1, breaks=25, main="MC approximate posterior") par(op) # Scatterplot of Monte Carlo sample # Can use to approximate posterior prob that theta2>theta1

(y1+1)/(n1+2) # True post mean

mean(theta1) # MC estimate

plot(theta1, theta2, xlim=c(0,1), ylim=c(0,1),
xlab=expression(theta[1]), ylab=expression(theta[2]))

abline(0,1)

mean(theta2>theta1)