Jim — Sep 29, 2013, 3:41 PM
# driving on I-75 example
# actual driving speeds
y <- c(76.2, 68.7, 77.5, 74.0, 75.0, 73.3, 76.0, 79.0, 72.5, 66.1,
69.1, 65.9, 77.1, 75.4, 75.2, 81.0, 78.0, 77.5, 58.3, 70.2)
# number of speeds smaller than 70, between 70 and 77, and over 77
f1 <- 5; f2 <- 9; f3 <- 6
# compare 2 posteriors
# NORMAL: normal sampling on y, usual noninformative prior
# GROUPED: based on observations (s, f)
library(LearnBayes)
# Use normchi2post function in LearnBayes to compute the log
# joint posterior of (mu, sigma^2) from complete data.
mycontour(normchi2post, c(65, 85, .01, 320), y,
xlab="MU", ylab="SIGMA^2")
# For grouped posterior, construct new function grouped to
# compute log posterior.
grouped <- function(theta, sf){
f1 <- sf[1]; f2 <- sf[2]; f3 <- sf[3]
mu <- theta[1]; sig <- sqrt(theta[2])
loglike <- f1 * log(pnorm(70, mu, sig)) +
f2 * log(pnorm(77, mu, sig) - pnorm(70, mu, sig)) +
f3 * log(1 - pnorm(77, mu, sig))
loglike - log(theta[2])
}
# Add the grouped posterior contour to the current graph.
mycontour(grouped, c(65, 85, .01, 320), c(f1, f2, f3),
col="red", add=TRUE)
legend(82, 300, legend=c("Normal", "Grouped"),
lty=1, col=c("black", "red"))