Example: Bayesian analysis of a bioassay experiment (logistic, nonconjugate)

1. Data and scientific setup

We have a classic toxicity/bioassay experiment: give animals different doses of a compound and record how many die.

Data (Table 3.1):

dose $x_i$​ (log g/ml)	animals $n_i$	deaths $y_i$
$-0.86$	5	0
$-0.30$	5	1
$-0.05$	5	3
$0.73$	5	5

So $k = 4$ dose levels, $n_i = 5$ each, total 20 animals.

Within each dose group it’s reasonable to treat animals as exchangeable → binomial model.

$y_i \mid θ_i \sim \text{Binomial}(n_i, θ_i)$

where $θ_i$​ = probability of death at dose $x_i$​.

But we don’t want four unrelated $θ_i$​’s — death probability should increase with dose. So we model $θ_i$​ as a function of $x_i$​.

---

2. Dose–response model (logistic regression)

A linear model $θ_i = α + β x_i$ would go outside $[0,1]$, so we use the logit link:

$\text{logit}(θ_i) = α + β x_i, \quad \text{that is} \quad θ_i = \text{logit}^{-1}(α + β x_i) = \frac{1}{1 + e^{-(α + β x_i)}}$

So the likelihood for group $i$ is $p(y_i \mid α, β) \propto \bigl(\text{logit}^{-1}(α + β x_i)\bigr)^{y_i} \bigl(1 – \text{logit}^{-1}(α + β x_i)\bigr)^{n_i – y_i}.$

Assuming independence across i:

$p(\mathbf{y} \mid α, β) = \prod_{i=1}^k p(y_i \mid α, β)$

---

3. Prior

They choose the simplest thing:

$p(α, β) \propto 1$

(a flat prior in $(α, β)$). This gives a posterior $p(α, β \mid \mathbf{y}) \propto \prod_{i=1}^k p(y_i \mid α, β).$

This model is not conjugate (logistic + flat prior ≠ nice closed-form posterior), so we need a numerical method.

---

4. Rough location of the posterior

They first run an ordinary logistic regression (i.e. MLE) to see roughly where the posterior mode is:

$(\hat α, \hat β) \approx (0.8, 7.7)$

with big standard errors (1.0 and 4.9), so the posterior will be fairly spread out. This helps choose a grid region.

---

5. Grid approximation to the posterior

Since it’s just 2 parameters, we can evaluate the unnormalized posterior on a grid:

• pick a rectangle in $(α, β)$-space, e.g. $α \in [-5, 10], \quad β \in [-10, 40],$ big enough to contain essentially all posterior mass
• at each grid point $(α, β)$, compute $\tilde p(α, β) = \prod_{i=1}^k p(y_i \mid α, β)$ (log-likelihood first, subtract max, then exponentiate → avoids under/overflow)
• normalize over the grid so that total mass = 1

From this discrete posterior they draw contour lines → Below image shows a slanted blob of high density.

---

6. Sampling from the discrete posterior

They want actual draws $(α^{(s)}, β^{(s)})$, so they do:

1. From the grid, marginalize over β to get a discrete $p(α \mid y)$.
2. Sample $α^{(s)}$ from that 1D discrete distribution.
3. Then, conditional on that sampled $α^{(s)}$, sample $β^{(s)}$ from the discrete conditional $p(β \mid α^{(s)}, y)$.
4. Add a tiny jitter (uniform over the grid cell) so points don’t all sit on grid centers.

Repeat 1000 times → Below image (scatter of 1000 posterior draws).

That procedure is just doing $p(α, β \mid y) = p(α \mid y) \, p(β \mid α, y)$

in discrete form — exactly the marginal/conditional idea from earlier.

---

7. LD50: posterior for a derived quantity

In toxicity, LD50 = dose at which death probability is 0.5.

In this model, $\text{logit}(θ) = α + β x.$

Set $θ = 0.5 \Rightarrow \text{logit}(0.5) = 0$, so $0 = α + β x \quad \Rightarrow \quad x = -\frac{α}{β}.$

So LD50 on the log-dose scale is $\text{LD}_{50} = -\frac{α}{β}.$

Now the nice Bayesian part: since we already have 1000 draws $(α^{(s)}, β^{(s)})$, we just compute $\text{LD}_{50}^{(s)} = -\frac{α^{(s)}}{β^{(s)}}$

for each draw → histogram

---

8. What if $β \le 0$?

If $β \le 0$, the dose–response curve is flat or decreasing → the notion “dose that kills 50%” stops making sense.

So they:

1. look at $P(β > 0 \mid y)$,
2. then condition on $β > 0$ when summarizing LD50.

In this dataset, all 1000 draws had $β > 0$ → posterior probability the drug is harmful is > 0.999.

Then the LD50 histogram is just from those (all) draws.

This also shows why posterior mean of LD50 can be a bad summary: if some draws had $β \le 0$, the ratio $-α/β$ could blow up or be meaningless, so you’d report a conditional posterior instead.

---

9. Practical notes from the example

• Grid methods work nicely in 2D, but not in high dimensions.
• You must choose grid bounds and spacing carefully — too narrow and you miss posterior mass; too wide and you waste resolution.
• Always work on the log scale for the posterior to avoid numerical underflow.
• This example is an instance of a generalized linear model (binomial + logit), but handled here with brute-force Bayesian computation.

---

Bottom line

• Model: $y_i \sim \text{Bin}(n_i, \text{logit}^{-1}(α + β x_i))$
• Prior: $p(α, β) \propto 1$
• Posterior: not conjugate → approximate on a grid
• Inference for LD50: compute $-α/β$ on posterior draws, condition on $β > 0$
• Posterior says: very high probability the compound is harmful, and we can quantify the dose at which it kills half the animals.