Question 1

Suppose the globe tossing data from Chapter 2 had turned out to be 9 water in 16 tosses. Using grid approximation, construct the posterior distribution. Use the same flat prior as in the book. Plot your posterior.

# define grid
p_grid <- seq(from = 0, to = 1, length.out = 16)

# define flat prior from book
prior <- rep(1, 16)

# compute likelihood at each value in grid
likelihood <- dbinom(9, size = 16, prob = p_grid)

# compute product of likelihood and prior
unstd.posterior <- likelihood * prior

# standardize the posterior, so it sums to 1
posterior <- unstd.posterior / sum(unstd.posterior)

#plot posterior
plot(p_grid, posterior, type = "b",
      xlab = "probability of water", 
      ylab = "posterior probability") + 
  mtext("16 points")

## integer(0)

Question 2

Redo number 1, but change your prior to be zero below p = .5 and a constant above .5 (there is code shown in the chapter for how to set this prior). First, describe in a sentence what this prior is stating conceptually and why it is better than our flat prior. Next, construct your posterior and graph it. What difference has the prior made?

# define grid
p_grid <- seq(from = 0, to = 1, length.out = 16)

# prior < .5 = 0; prior > .5 = 1
prior <- ifelse(p_grid < 0.5 , 0 , 1 )

# compute likelihood at each value in grid
likelihood <- dbinom(9, size = 16, prob = p_grid)

# compute product of likelihood and prior
unstd.posterior <- likelihood * prior

# standardize the posterior, so it sums to 1
posterior <- unstd.posterior / sum(unstd.posterior)

#plot posterior
plot(p_grid, posterior, type = "b",
      xlab = "probability of water", 
      ylab = "posterior probability") + 
  mtext("16 points")

## integer(0)

ANSWER

This prior is stating conceptually…

It is better than a flat prior because it actually influences the posterior. Other reasons too??

Question 3

Now redo number 1, including the flat prior, using quadratic approximation. What are the mean and standard deviation of your posterior?

globe.qa <- quap(
  alist(
    W ~ dbinom(W + L, p), # binomial likelihood
    p ~ dunif(0, 1) #uniform prior
    ),
  data = list(W = 6, L = 3))

# display summary of quadratic approximation
precis(globe.qa)
##        mean        sd      5.5%     94.5%
## p 0.6666663 0.1571339 0.4155361 0.9177966

ANSWER mean: .67 sd: .16