Paper 35 (Applied Bayesian Statistics)

www.maths.cam.ac.uk/postgrad/part-iii/files/pastpapers/2014/paper_35.pdf

Table of contents
- 1 Paper 35
  - a 1
    - Solution a
  - b 1
    - Solution b
  - c 1
    - Solution c
  - d 1
    - Solution d
  - e 1
    - Solution e
  - f 1
    - Solution f
  - g 1
    - Solution g
  - h 1
    - Solution h
- 2 Paper 35
  - a 2
    - Solution a
  - b 2
    - Solution b
  - c 2
    - Solution c
  - d 2
    - Solution d
  - e 2
    - 1 e
      - Solution 1
    - 2 e
      - Solution 2
    - 3 e
      - Solution 3
    - 4 e
      - Solution 4
  - f 2
    - Solution f
  - g 2
    - Solution g
  - h 2
    - Solution h
  - i 2
    - Solution i
  - j 2
    - Solution j
- 3 Paper 35
  - a 3
    - Solution a
  - b 3
    - Solution b
  - c 3
    - Solution c
  - d 3
    - Solution d
  - e 3
    - Solution e
  - f 3
    - Solution f
  - g 3
    - Solution g
  - h 3
    - Solution h
  - i 3
    - Solution i
- 4 Paper 35
  - a 4
    - Solution a
  - b 4
    - Solution b
  - c 4
    - Solution c
  - d 4
    - Solution d
  - e 4
    - Solution e
  - f 4
    - Solution f
  - g 4
    - Solution g
  - h 4
    - Solution h
  - i 4
    - Solution i
  - j 4
    - Solution j

Solution

The likelihood function for the ordered arrival times of a Poisson process, including the absence of further arrivals before

T

, is

L (λ) = λ^{n} e^{- λ T}, 0 < t_{1} < \dots < t_{n} < T .

(1)

The arrival-time constraint does not depend on

λ

. Equivalently, the Poisson distribution of the total count gives

L (λ) = e^{- λ T} (λ T)^{n} / n!

, which differs by a parameter-independent factor. Conditional on the count, the Poisson process conditional arrival times have density

n! / T^{n}

, independent of

λ

. Thus the total count contains all the information about the rate when the exposure is known; it is a sufficient statistic.

Solution

 0  0

A conjugate prior is a family of prior distributions whose members remain in that family after updating by the likelihood function. Here multiplying a shape-rate gamma distribution density by the Poisson process likelihood changes its power of

λ

and its exponential rate, leaving a gamma distribution. The parameters change with the observations; conjugate prior does not mean that the Bayesian posterior equals the prior distribution.

Solution

 0  0

Independent increments of the Poisson process give

M ∣ λ \sim Poisson (λ)

for a one-hour interval. The Poisson distribution therefore gives

P (M = 0 ∣ λ) = e^{- λ} .

(1)

This conditional probability still depends on the unknown rate.

Solution

 0  0

The posterior predictive probability averages the conditional probability over the Bayesian posterior. Using its gamma distribution density,

p_{0} = \frac{B ^{A}}{Γ ( A )} \int_{0}^{\infty} λ^{A - 1} e^{- (B + 1) λ} d λ = (\frac{b + T}{b + T + 1})^{a + n} .

(1)

For fixed

a, b

and

n / T \to λ < \infty

lo g p_{0} = - A lo g (1 + 1/ B) = - \frac{A}{B} + O (A / B^{2}) = - \frac{n}{T} + o (1) .

(2)

Consequently

p_{0} / e^{- n / T} \to 1

, and both tend to

e^{- λ}

. This is reasonable because the posterior mean tends to the observed rate and the posterior variance tends to zero. Averaging

e^{- λ}

then approaches evaluating it at the estimated rate: with abundant observations, predictive uncertainty about the rate vanishes.

Solution

 0  0

There is a notation problem here. The posterior predictive probability

p_{0}

calculated above is already a fixed number conditional on

a, b, n, T

. Literally,

P (p_{0} < 1/2 ∣ D) = 1 {(B / (B + 1))^{A} < 1/2} .

(1)

No simulation is needed for that interpretation. The known exposure

T

is also necessary, despite its omission from this part's list of inputs.

The natural uncertain quantity is instead

q (λ) = P (M = 0 ∣ λ) = e^{- λ}

. For its posterior probability of being below one half, draw the rate from its gamma distribution Bayesian posterior and average an indicator. Rough BUGS code is

model {
  lambda ~ dgamma(a+n, b+T)
  q <- exp(-lambda)
  belowHalf <- step(lambda-log(2))
}

The monitored average of belowHalf estimates $P (λ > lo g 2 ∣ D)$ , equivalently

1 - F_{Gamma (A, B)} (lo g 2)

. The equality boundary has zero probability. This code samples the already updated Bayesian posterior; adding the count likelihood function again would count the data twice.

Solution

 0  0

For an Inhomogeneous Poisson process, the likelihood function is the product of the intensities at arrivals times the exponential of minus the integrated intensity. The integrated intensity is

θ + 2 (T - θ) = 2 T - θ

. There are

j (θ)

arrivals at rate one and

n - j (θ)

at rate two, so

L (θ) = e^{θ - 2 T} 2^{n - j (θ)} \propto e^{θ} 2^{- j (θ)}, 0 < θ < T .

(1)

Set

t_{0} = 0

and

t_{n + 1} = T

to include changes before the first or after the last arrival. An arrival exactly at the change has probability zero, so the convention for

j

at that point does not affect the Bayesian posterior.

Solution

 0  0

A uniform prior on

(0, T)

makes the posterior density proportional to the Poisson change-point posterior likelihood above. The zeros trick introduces an observed zero with Poisson distribution mean

K - lo g L (θ)

, making its likelihood function equal to

e^{- K} L (θ)

. Choose

K = T + n lo g 2 + 1

; since

lo g L (θ) \leq n lo g 2

, this mean is strictly positive. Rough BUGS code, with zero=0 supplied as data, is

model {
  theta ~ dunif(0,T)
  for (i in 1:n) {
    before[i] <- step(theta-time[i])
  }
  j <- sum(before[])
  logL <- theta-2*T+(n-j)*log(2)
  zero ~ dpois(K-logL)
}

For

n = 0

, omit the array and set j <- 0. Monitor the sampled theta to obtain its posterior mean, credible interval and interval probabilities.

There is also an exact sampling method. On

(t_{j}, t_{j + 1})

the posterior density is proportional to

2^{- j} e^{θ}

, so choose the interval with weights

w_{j} = 2^{- j} (e^{t_{j + 1}} - e^{t_{j}}), j = 0, \dots, n .

(1)

Then draw

U

from a uniform distribution on

(0, 1)

and set

θ = lo g (e^{t_{j}} + U (e^{t_{j + 1}} - e^{t_{j}}))

. These weighted interval draws sample the posterior directly, without asking a local Markov chain Monte Carlo update to cross its discontinuities.

Solution

 0  0

For a regular one-parameter sampling distribution, the Fisher information and Jeffreys prior are

I (θ) = E_{θ} [(\frac{\partial}{\partial θ} lo g p_{Y} (Y ∣ θ))^{2}], π_{J} (θ) \propto I (θ) .

(1)

Under the usual differentiation and integrability conditions,

I (θ) = - E_{θ} [\partial_{θ}^{2} lo g p_{Y} (Y ∣ θ)]

. This prior distribution transforms as a density under smooth one-to-one reparameterizations, so the rule is coordinate invariant. Its integral need not be finite; posterior propriety must still be established if it is an improper prior.

b

 0  0

Solution

 0  0

For a binomial distribution, the score function is

\frac{Y}{θ} - \frac{n - Y}{1 - θ} = \frac{Y - n θ}{θ ( 1 - θ )} .

(1)

Its squared expected value is

I (θ) = n / [θ (1 - θ)]

, using the binomial distribution variance

n θ (1 - θ)

. Hence the Jeffreys prior is

π_{J} (θ) = \frac{1}{π θ ( 1 - θ )}, 0 < θ < 1,

(2)

the Beta distribution

Beta (1/2, 1/2)

. The factor

n

is independent of

θ

and disappears on normalization.

c

 0  0

Solution

 0  0

Independent Jeffreys priors and the two independent binomial distribution likelihood factors give, by Beta-binomial conjugacy,

p_{M} ∣ D \sim Beta (7/2, 3/2), p_{T} ∣ D \sim Beta (3/2, 7/2) .

(1)

The Bayesian posteriors remain independent because each observation factor involves only its own probability. Their posterior means are respectively $0.7$ and $0.3$ . Notice that

p_{T}

is the probability of saying milk first when tea was first, rather than the probability of a correct tea identification.

d

 0  0

Solution

 0  0

In Markov chain Monte Carlo, retain the derived quantity

δ = p_{M} - p_{T}

at every iteration. Rough BUGS code using the actual observations is

model {
  pM ~ dbeta(0.5,0.5)
  pT ~ dbeta(0.5,0.5)
  milkAnswers ~ dbin(pM,4)
  teaAnswers ~ dbin(pT,4)
  delta <- pM-pT
  positive <- step(delta)
}

Supply milkAnswers=3 and teaAnswers=1. Summarize delta by its posterior mean, empirical quantiles and credible interval; the average of positive estimates

P (p_{M} > p_{T} ∣ D)

. Here $E [δ ∣ D] = 0.4$ . Check Markov chain Monte Carlo convergence diagnostics before interpreting the simulation. Since both Bayesian posteriors are independent known Beta distributions, direct independent sampling is an equally valid, simpler way to obtain the same summaries.

e

 0  0

1

 0  0

Solution

 0  0

Start with the Beta distribution draw in the Beta-binomial exchangeable coupling. Write

s = α + β

and

m = \frac{α}{s}, v = \frac{α β}{s ^{2} ( s + 1 )} .

(1)

These are the expected value and variance of

p_{M}

. The fixed

n

in this coupling is a dependence parameter chosen by the investigator, not necessarily the observed number of cups. The four numbered stages together address the single printed correlation request.

2

 0  0

Solution

 0  0

Conditional on

p_{M}

, the auxiliary binomial distribution gives

E [X / n ∣ p_{M}] = p_{M}

and

Var (X / n ∣ p_{M}) = p_{M} (1 - p_{M}) / n

. Thus large $n$ makes $X / n$ close to $p_{M}$ . The law of total variance also gives

Var (X) = nm (1 - m) + n (n - 1) v = \frac{n α β ( s + n )}{s ^{2} ( s + 1 )} .

(1)

The latent variable encodes the first draw increasingly accurately.

3

 0  0

Solution

 0  0

The next Beta distribution has conditional expected value and variance

E [p_{T} ∣ X] = \frac{α + X}{s + n}, Var (p_{T} ∣ X) = \frac{( α + X ) ( β + n - X )}{( s + n ) ^{2} ( s + n + 1 )} .

(1)

For large

n

, its expected value approaches

X / n

, while its variance is at most

1/ [4 (s + n + 1)]

. The final draw therefore stays close to the first. These are prior distributions for the cup probabilities: the auxiliary count is a device for constructing dependence, rather than additional observed tea data.

4

 0  0

Solution

 0  0

The law of iterated expectation yields

E [p_{T} ∣ p_{M}] = \frac{α + n p _{M}}{s + n}, Cov (p_{M}, p_{T}) = \frac{n}{s + n} v .

(1)

One may obtain

Var (p_{T}) = v

either from its unchanged marginal distribution proved below or directly from the law of total variance: the variance of its conditional mean is

n v / (s + n)

and its expected conditional variance is

s v / (s + n)

. Hence the correlation coefficient is exactly

Corr (p_{M}, p_{T}) = \frac{n}{n + α + β} ⟶ 1.

(2)

Also

E [(p_{T} - p_{M})^{2}] = 2 v s / (s + n) \to 0

. The two probabilities become close while keeping their original beta marginals.

f

 0  0

Solution

 0  0

Apply the law of iterated expectation to the auxiliary binomial distribution:

E [X] = E [E [X ∣ p_{M}]] = n E [p_{M}] = \frac{n α}{α + β} .

(1)

This expectation is over the full prior distribution construction, not over a particular observed auxiliary count.

g

 0  0

Solution

 0  0

Using the conditional Beta distribution expected value and the law of iterated expectation,

E [p_{T}] = E [\frac{α + X}{s + n}] = \frac{α + n α / s}{s + n} = \frac{α}{α + β} .

(1)

Equality of the expected values alone does not prove equality of the marginal distributions; the later exchangeability argument does.

h

 0  0

Solution

 0  0

Multiply the Beta distribution density of

p_{M}

, the conditional binomial distribution mass of

X

, and the conditional Beta distribution density of

p_{T}

. With

B (a, b)

the beta function, the full joint density, relative to counting measure in

x

and Lebesgue measure in the two probabilities, is

f (p_{T}, x, p_{M}) = \frac{( x n )}{B ( α , β ) B ( α + x , β + n - x )} (p_{T} p_{M})^{α + x - 1} [(1 - p_{T}) (1 - p_{M})]^{β + n - x - 1} .

(1)

Here

x = 0, \dots, n

and both probabilities lie in

(0, 1)

. The printed joint expression omits the factor $1/ B (α + x, β + n - x)$ . It is a constant when

x

is fixed and only the two probabilities vary, but is not a constant for the full joint probability distribution. Keeping it is essential when summing over the latent variable.

i

 0  0

Solution

 0  0

Two exchangeable random variables satisfy

(X, Y) = d (Y, X)

: their joint probability distribution is unchanged by swapping the coordinates. Thus for every measurable set

A

P (X \in A) = P ((X, Y) \in A \times Y) = P ((Y, X) \in A \times Y) = P (Y \in A) .

(1)

Exchangeability implies identical marginal distributions, but does not imply independence.

j

 0  0

Solution

 0  0

For each auxiliary count

x

, the correctly normalized joint expression above is symmetric in

p_{M}, p_{T}

. Summing it over

x = 0, \dots, n

preserves that symmetry, so these are exchangeable random variables. Their marginal distributions are therefore identical. Since

p_{M}

was generated from a Beta distribution,

p_{T} \sim Beta (α, β) .

(1)

The Beta-binomial exchangeable coupling changes the dependence, while leaving both prior distributions unchanged.

3

 0  0

a

 0  0

Solution

 0  0

At zero a centered normal distribution with standard deviation

σ

has density

1/ (2 π σ)

. The two prior standard deviations are

1/ n_{0}

and

c / n_{0}

, giving

p_{0} (0) / p_{1} (0) = c .

(1)

The wider prior distribution has a lower density at its center.

b

 0  0

Solution

 0  0

Let

q_{0} = n_{0}

q_{1} = n_{0} / c^{2}

and

y = \overline{y}

. These are the two prior precision parameters. Completing the square in

n (β - y)^{2} + q_{i} β^{2}

gives the Normal-normal conjugacy update

β ∣ y, H_{i} \sim N (\frac{n}{n + q _{i}} y, \frac{1}{n + q _{i}}), i = 0, 1.

(1)

Thus the posterior mean is a precision-weighted average of the observed mean and the prior center zero. The posterior variance is the reciprocal of the total precision.

c

 0  0

Solution

 0  0

Write

\overline{Y} = β + ε

, with independent

ε \sim N (0, 1/ n)

and

β ∣ H_{i} \sim N (0, 1/ q_{i})

. The convolution of independent random variables is again a normal distribution, so the prior predictive laws are

\overline{Y} ∣ H_{i} \sim N (0, V_{i}), V_{i} = 1/ n + 1/ q_{i} .

(1)

These are predictive distributions before observing

y

, hence the Bayesian model evidence for the observed mean. The residual information in the original observations is common to both models and cancels in their Bayes factor.

d

 0  0

Solution

 0  0

Divide the two normal distribution predictive densities to obtain the Bayes factor

B_{01} (y) = \frac{V _{1}}{V _{0}} exp [- \frac{y ^{2}}{2} (\frac{1}{V _{0}} - \frac{1}{V _{1}})] .

(1)

y = 0

this becomes

B_{01} (0) = \frac{n _{0} / n + c ^{2}}{n _{0} / n + 1} .

(2)

The wider prior distribution spreads its predictive mass over more possible means, giving the narrower model more Bayesian model evidence for observations very near zero. Away from zero, the exponential term opposes that factor.

e

 0  0

Solution

 0  0

Under the point hypothesis

β = 0

, the test statistic

n \overline{Y}

has a standard normal distribution. The observed statistic is three, giving a two-sided p-value

2 [1 - Φ (3)] ≃ 0.00270

. This is conventionally strong evidence against that point hypothesis. The narrow model

H_{0}

is a continuous prior distribution around zero, rather than a point hypothesis.

For

n_{0} / n = 100

and

c = 100

, Normal-normal conjugacy gives

β ∣ y, H_{0} \sim N (\frac{3}{101 n}, \frac{1}{101 n}), β ∣ y, H_{1} \sim N (\frac{300}{101 n}, \frac{100}{101 n}) .

(1)

The narrow-model posterior mean is strongly pulled toward zero and its standard deviation is about

0.0995/ n

. The wide-model posterior mean is about

2.9703/ n

, very close to the observation, with standard deviation about

0.9950/ n

. Each model produces a markedly different posterior, so choosing between them requires their predictive evidence.

f

 0  0

Solution

 0  0

The requested conclusion is false for the printed parameter values. They give

V_{0} = 1.01/ n

and

V_{1} = 101/ n

. Substituting

y = 3/ n

in the Bayes factor yields

B_{01} = 10 exp (- \frac{891}{202}) ≃ 0.12144 < 1.

(1)

Thus the Bayes factor favours $H_{1}$ by about $8.23$ to one, not

H_{0}

. Nor does the conclusion hold for arbitrary large

n_{0}, c

: both parameters enter the expression explicitly. With

n_{0} / n = 100

, taking a much wider alternative, for example

c = 1000

, would instead give

B_{01} ≃ 1.1563 > 1

. That is a different prior assumption and cannot repair the printed calculation silently.

g

 0  0

Solution

 0  0

With equal model prior probabilities, Bayes factor updating gives

w_{0} = B_{01} / (1 + B_{01})

and

w_{1} = 1/ (1 + B_{01})

. The Bayesian model averaging posterior density is therefore

p (β ∣ y) = \frac{B _{01}}{1 + B _{01}} p_{0} (β ∣ y) + \frac{1}{1 + B _{01}} p_{1} (β ∣ y) .

(1)

This is a two-component mixture model of the normal distributions already derived. For the numerical observation above, $w_{0} ≃ 0.10829$ and $w_{1} ≃ 0.89171$ . Consequently it is predominantly the wide-model posterior, not predominantly the narrow one.

h

 0  0

Solution

 0  0

Put

κ_{i} = n / (n + q_{i})

. The posterior mean of the Gaussian practical-null mixture is

E [β ∣ y] = [w_{0} (y) κ_{0} + w_{1} (y) κ_{1}] y .

(1)

Near zero,

B_{01} (0) > 1

, so the narrow component has high posterior probability and

κ_{0}

is small when

n_{0} ≫ n

. This pulls the Bayesian model averaging posterior toward zero. For the given ratios, writing

z = n y

gives

B_{01} (y) = 10 e^{- 99 z^{2} /202}, B_{01} > 1 ⟺ ∣ z ∣ < 2.16753 \dots .

(2)

Thus the order statement

y = O (n^{- 1/2})

alone does not guarantee strong pull toward zero:

z = 3

already has the opposite model preference.

∣ y ∣

grows,

V_{1} > V_{0}

makes

B_{01} (y) \to 0

, and the mixture model approaches

N (κ_{1} y, 1/ (n + q_{1}))

. Its center is exactly $κ_{1} y$ , and approximately $y$ only for a sufficiently diffuse wide prior, meaning

q_{1} ≪ n

. Here

κ_{1} = 100/101

, so the relative displacement is about one percent. Large observations alone do not remove this finite-prior shrinkage.

i

 0  0

Solution

 0  0

In a normal linear model, attach a continuous spike-and-slab prior to each regression coefficient. With suitably scaled predictors, introduce indicators

γ_{j} \sim Bernoulli (π)

and set

β_{j} ∣ γ_{j} = 0 \sim N (0, s_{0 j}^{2}), β_{j} ∣ γ_{j} = 1 \sim N (0, s_{1 j}^{2}), s_{0 j} ≪ s_{1 j} .

(1)

The narrow component describes practically negligible effects; the wide component permits substantial ones. Fit the joint Bayesian posterior of coefficients, indicators and any unknown residual variance. Bayesian model averaging gives shrinkage toward zero for poorly supported effects, while

P (γ_{j} = 1 ∣ D)

quantifies wide-component support. A shared Beta distribution prior on

π

can represent uncertainty about how many effects are substantial.

Select on scientifically meaningful effect size, for example a high

P (∣ β_{j} ∣ > δ_{j} ∣ D)

for a prechosen threshold in meaningful predictor units. Membership in the wide component alone does not imply a large realized effect: its normal distribution still permits values near zero. Correlated predictors also require interpretation of the joint Bayesian posterior, rather than treating each coefficient as an isolated test.

4

 0  0

a

 0  0

Solution

 0  0

The log odds ratio for treatment relative to control is the difference of the two log odds:

lo g \frac{θ _{T} / ( 1 - θ _{T} )}{θ _{C} / ( 1 - θ _{C} )} = (α + β /2) - (α - β /2) = β .

(1)

Thus the treatment-to-control odds ratio is $e^{β}$ . Its direction refers to whatever event the binomial count records.

b

 0  0

Solution

 0  0

The intercept is the midpoint of the two log odds:

α = \frac{1}{2} [logit θ_{C} + logit θ_{T}] .

(1)

Equivalently

e^{α}

is the geometric mean of the two odds. The logistic function evaluated at

α

represents a central event probability, and approximates the average of the two group probabilities when

β

is small. It is not generally their arithmetic average, nor is

α

specifically the control-group log odds.

c

 0  0

Solution

 0  0

Use a proper normal distribution prior centered at zero on the log odds ratio. For example,

β \sim N (0, 1)

(1)

puts approximately 95 percent of its mass between

- 2

and

2

, giving odds ratios roughly between

e^{- 2}

and

e^{2}

, close to

1/8

and

8

. Centering at zero treats reciprocal odds ratios symmetrically. If the desired central 95 percent interval is exactly

(1/8, 8)

, use standard deviation

lo g 8/ Φ^{- 1} (0.975) ≃ 1.0610

. A soft prior is appropriate for implausibility, whereas a bounded uniform prior would declare effects outside the limits impossible.

d

 0  0

Solution

 0  0

Calling the effects exchangeable random variables means that their joint prior distribution is unchanged by permuting study labels. In a hierarchical Bayesian model, conditional independent draws

β_{j} ∣ μ, τ \sim N (μ, τ^{2})

achieve this, and integrating shared hyperparameters induces dependence between studies. This permits partial pooling without asserting that all effects are identical.

The assumption is reasonable when the studies concern comparable treatments, populations, outcomes and follow-up, and no known study characteristic gives one effect a systematically different prior center. Relevant differences can instead enter a linear regression for study effects, after which the residual effects may be exchangeable. The numerical table counts deaths, so its event probabilities are mortality probabilities and

β_{j} < 0

indicates a lower mortality odds ratio. Interpreting those counts as beneficial responses would reverse the clinical meaning.

e

 0  0

Solution

 0  0

The overall mean

μ

and the between-study heterogeneity

τ

encode different information. A proper broad normal distribution prior such as

μ \sim N (0, 2^{2})

is one possible weak prior for the mean log odds ratio; information about the spread of trials alone does not determine its center.

A concrete prior calibration for normal random-effect range can make the factor-of-50 statement simultaneous across all six trials. Conditional on

τ

, each pair difference has normal distribution

β_{j} - β_{k} \sim N (0, 2 τ^{2})

. Set

R = lo g 50

m = (2 6) = 15

z = Φ^{- 1} (1 - 0.05/ (2 m))

and

τ \sim Uniform (0, A), A = \frac{R}{2 z} ≃ 0.9424.

(1)

For every

τ \leq A

, each pair exceeds

R

in absolute value with probability at most

0.05/ m

. The union bound therefore gives

P (max_{j} β_{j} - min_{j} β_{j} > R) \leq 0.05

, and integrating over the uniform prior preserves that bound. This is one explicit interpretation of “very unlikely”; a different elicited probability would change the bound. A smoother proper scale prior could be calibrated similarly.

The proposed

1/ τ

improper prior is unsuitable. The observed-data likelihood function approaches the positive common-effect likelihood as

τ ↓ 0

. After restricting

μ

and the intercepts to a compact interior region, it is bounded below there by a positive constant. Hence

\int_{0}^{ε} L (μ, τ) \frac{d τ}{τ} = \infty.

(2)

This is an improper posterior from a log-uniform random-effect scale prior. Proper conditional sampling distributions do not repair the improper joint posterior, and an arbitrary tiny cutoff would make inference depend on that cutoff.

f

 0  0

Solution

 0  0

Represent the independent locally flat intercept prior distributions by broad finite uniform priors, for example on

(- 10, 10)

; this is proper and approximately constant over plausible mortality logits. With

A

calibrated above, rough BUGS code is

model {
  mu ~ dnorm(0,0.25)
  tau ~ dunif(0,A)
  invtau2 <- pow(tau,-2)
  for (j in 1:J) {
    alpha[j] ~ dunif(-10,10)
    beta[j] ~ dnorm(mu,invtau2)
    logit(thetaC[j]) <- alpha[j]-beta[j]/2
    logit(thetaT[j]) <- alpha[j]+beta[j]/2
    rC[j] ~ dbin(thetaC[j],nC[j])
    rT[j] ~ dbin(thetaT[j],nT[j])
    oddsRatio[j] <- exp(beta[j])
  }
}

Use

J = 6

and supply treated death counts

(3, 7, 5, 102, 32, 22)

with totals

(38, 114, 69, 1533, 209, 680)

, and control death counts

(3, 14, 11, 127, 40, 39)

with totals

(39, 116, 93, 1520, 218, 674)

. In BUGS, the second dnorm argument is a precision parameter, so 0.25 corresponds to variance four. Initialize the positive scale away from zero. Monitor

μ, τ

and study odds ratios, checking Markov chain Monte Carlo convergence diagnostics and sensitivity to the finite intercept bounds and scale prior distribution. The fitted hierarchy combines binomial sampling uncertainty with between-study heterogeneity.

g

 0  0

Solution

 0  0

For the common Bayesian deviance convention

D = - 2 lo g L

used in all three models,

p_{D} = \overline{D} - D (\overline{θ}), DIC = \overline{D} + p_{D} = D (\overline{θ}) + 2 p_{D} .

(1)

Here Dhat is

D (\overline{θ})

, an at-posterior-mean fit measure, and

p_{D}

is an effective parameter count. The independent model's Dhat of 53.1 is almost identical to the exchangeable model's 53.2; both improve on the common model's 57.8. The common model's

p_{D} = 7

corresponds to six intercepts plus one shared effect. Independence uses roughly twelve effective parameters. Partial pooling reduces the exchangeable model's effective complexity to about 8.7 while retaining nearly the same fitted likelihood function as independence.

The exchangeable model has the lowest reported DIC, but the common model is competitive. Their difference is only about 1.3, whereas independence is worse by about 6.3. The deviance information criterion measures penalized fit for a predictive comparison, not model posterior probabilities, and these numbers do not establish overwhelming evidence for heterogeneity. The displayed exchangeable

\overline{D} + p_{D}

61.9 + 8.7 = 70.6

, rather than the printed 70.5; rounding of the underlying values can account for a tenth and does not change this interpretation.

h

 0  0

Solution

 0  0

Conditional on

μ, ψ > 0

, let

λ \sim χ_{4}^{2}

and independently

Z \sim N (0, 1)

. The normal distribution with the stated precision parameter can be generated as

β_{j} = μ + \frac{2 ψ Z}{λ} .

(1)

Therefore

\frac{β _{j} - μ}{ψ} = \frac{Z}{λ /4} \sim t_{4},

(2)

by the defining normal distribution and chi-squared distribution representation of Student's t-distribution. Its density is

f (β_{j} ∣ μ, ψ) = \frac{3}{8 ψ} (1 + \frac{( β _{j} - μ ) ^{2}}{4 ψ ^{2}})^{- 5/2} .

(3)

Thus $μ$ is the location and $ψ$ is the scale, not the standard deviation:

Var (β_{j} ∣ μ, ψ) = 2 ψ^{2}

. Each study gets its own independent chi-squared draw in this Student t random-effect model.

i

 0  0

Solution

 0  0

The Student t random-effect model is useful when most studies are comparable but occasional genuine departures are more frequent than a normal distribution hierarchy allows. Its heavier tails permit a study effect far from

μ

without forcing a large common between-study heterogeneity scale on every study. In the Gaussian scale mixture representation, a small study-specific

λ_{j}

lowers its precision parameter and weakens its shrinkage.

Use this as robust partial pooling when occasional atypical effects are plausible. Known systematic population or design differences should still be modeled explicitly; a heavy tail cannot identify or correct within-study bias by itself.

j

 0  0

Solution

 0  0

Fit both the normal and Student t random-effect model with comparable proper prior distributions. Compare priors on the same spread measure: a normal distribution scale

τ

is a standard deviation, whereas the

t_{4}

standard deviation is

2 ψ

Use a posterior predictive check: draw study effects and binomial counts from each fitted hierarchy and compare replicated dispersion and extreme study contrasts with the observations. For predicting a new study, generate a new effect from the hierarchy rather than reusing an existing fitted effect. A Leave-one-out cross-validation with entire studies held out can compare integrated predictive probabilities for both arms of each omitted trial, averaging over hyperparameters and its unobserved study effect. Leave-one-study-out influence analysis also reveals whether the difference is driven by a single trial.

Prefer the heavier-tailed hierarchy if it improves the relevant predictive checks and held-out study predictions robustly to reasonable prior choices. The deviance information criterion can supplement the comparison, but its effective parameter count can depend on the latent-variable representation, and six studies give limited information about tail shape. A small numerical criterion difference alone is insufficient evidence.

Past exam of the mathematics course of the University of Cambridge / 2014 / iii / Paper 35 (Applied Bayesian Statistics)

 Discussion (0)  Subscribe (1)

 Discussion (0)