Reproducing a didactic example of Lindley (1993)

Question

Lindley (1993) discusses the following mixed discrete and continuous prior for the tea tasting lady experiment, where $\pi$ is probability of a correct classification: $p(\pi=0.5) = 0.8$ (discrete part) and $p(\pi | \pi > 0.5) = 1.6 (1-\pi)$ (continuous part). The prior reflects a skeptic stance towards Muriel Bristol's tea tasting ability, but not entirely discards the possibility pf $\pi > 0.5$. Lindley states the posterior probability for $P(\pi = 0.5)$ for several combinations of correct and wrong classifications, e.g., 5 right, one wrong: $P(\pi = 0.5 | 5R, 1W) = 0.59$.

The math for this uncommon mixed prior is straightforward: By Bayes' theorem $$ P(\pi = 0.5 | 5R, 1W) = \frac{P(5R,1W|\pi=0.5)p(\pi = 0.5)}{P(5R,1W)} $$ and the denominator leads to an integral of the continuous part of the posterior, which I evaluate numerically. However, I failed to reproduce the numbers, coming up with $P(\pi = 0.5 | 5R, 1W) = 0.605$. The discrepancy for $P(\pi = 0.5 | 6R)$ is even higher. Also for the other example in Lindley (1993) [wine tasting] I find the maximum-a-posteriori estimate closer to 0.8 than to 0.75 (what Lindley states).

prior_tea_discrete <- function(p){(p==0.5) * 0.8}
prior_tea_cont <- function(p){1.6*(1-p)}
y1 <- c(1,1,1,1,1,0)
lik1 <- function(p){choose(length(y1),sum(y1)) * p^sum(y1) * (1-p)^{length(y1)-sum(y1)}}
lik1 <- Vectorize(lik1)
y2 <- c(1,1,1,1,1,1)
lik2 <- function(p){choose(length(y2),sum(y2)) * p^sum(y2) * (1-p)^{length(y2)-sum(y2)}}
lik2 <- Vectorize(lik2)

post_tea_cont_1 <- function(p)lik1(p)*prior_tea_cont(p)
v1 <- integrate(post_tea_cont_1, 0.5, 1)$value
post_tea_discrete_1 <- lik1(0.5)* prior_tea_discrete(0.5)
post_prob_discr1 <- post_tea_discrete_1/(post_tea_discrete_1 + v1) 
# 0.605, slightly higher than Lindley's (1993) 0.59

post_tea_cont_2 <- function(p)lik2(p)*prior_tea_cont(p)
v2 <- integrate(post_tea_cont_2, 0.5, 1)$value
post_tea_discrete_2 <- lik2(0.5)* prior_tea_discrete(0.5)
post_prob_discr2 <- post_tea_discrete_2/(post_tea_discrete_2 + v2) 
# 0.312, substantially higher than Lindley's 0.23

Reference

Lindley, D.V. (1993). The analysis of experimental data: The appreciation of tea and wine. Teaching Statistics, 15(1), 22-25.

Answer 1

Yes, and like you, for five right and one wrong, I get (using Lindley's $P$ rather than your $\pi$) $$\mathbb P(P=\tfrac12 \mid RRRRRW) = \frac{ 0.8 P^5(1-P)\Big|_{P=1/2}}{0.8 P^5(1-P)\Big|_{P=1/2} + \int\limits_{1/2}^1 1.6(1-P)P^5(1-P)\, dP} \\= \frac{112}{185} \approx 0.6054$$

and for six right and none wrong, I get $$\mathbb P(P=\tfrac12 \mid RRRRRR) = \frac{ 0.8 P^6\Big|_{P=1/2}}{ 0.8 P^6\Big|_{P=1/2} + \int\limits_{1/2}^1 1.6(1-P)P^6\, dP} \\= \frac{112}{359} \approx 0.3120.$$

---

For the wine tasting, Lindley used a prior of $48(1-P)(P-\frac12)$ on $[\frac12,1]$. By assigning a zero prior probability to $P<\frac12$ he seems to have forgotten his separate advocacy of "Cromwell's rule".

In the 5R,1W case, I get a posterior density of $\frac{258048}{275}(1-P)^2(P-\frac12)P^5$ on $[\frac12,1]$ with a maximum, i.e. a mode or MAP, at $P=\frac{\sqrt{41}+19}{32} \approx 0.794$, presumably the same as your calculation, and a mean of $\frac{1058}{1375} \approx 0.769$.

In the 6R,0W case, I get a posterior density of $\frac{258048}{1291}(1-P)(P-\frac12)P^6$ on $[\frac12,1]$ with a maximum at $P=\frac{\sqrt{57}+21}{32} \approx 0.892$, higher than the roughly $0.87$ suggested in Lindley's Figure 1, and a mean of $\frac{5397}{6455} \approx 0.836$.