Hearts: a mixture model for count data
The table below presents data given by Berry (1987) on the effect of a drug used to treat patients with frequent premature ventricular contractions (PVCs) of the heart.
![[hearts1]](hearts1.bmp)
Farewell and Sprott (1988) model these data as a mixture distribution of Poisson counts in which some patients are "cured" by the drug, whilst others experience varying levels of response but remain abnormal. A zero count for the post-drug PVC may indicate a "cure", or may represent a sampling zero from a patient with a mildly abnormal PVC count. The following model thus is assumed:
x
i ~ Poisson(
λi) for all patients
y
i ~ Poisson(
βλi) for all
uncured patients
P(cure) =
θ
To eliminate nuisance parameters l
i, Farewell and Sprott use the conditional distribution of y
i given t
i = x
i + y
i. This is equivalent to a binomial likelihood for y
i with denominator t
i and probability p = b /(1+b) (see Cox and Hinkley, 1974 pp. 136-137 for further details of the conditional distribution for Poisson variables). Hence the final mixture model may be expressed as follows:
P(y
i = 0 | t
i ) =
θ + (1
- θ) (1
- p) t
i P(y
i | t
i ) = (1
- θ) (t
i! / (y
i! (t
i-y
i)!)) (p
yi (1
- p)
(ti - yi) y
i = 1,2,...,t
i The BUGS code for this model is given below:
model
{
for (i in 1 : N) {
y[i] ~ dbin(P[state1[i]], t[i])
state[i] ~ dbern(theta)
state1[i] <- state[i] + 1
t[i] <- x[i] + y[i]
prop[i] <- P[state1[i]]
}
P[1] <- p
P[2] <- 0
logit(p) <- alpha
alpha ~ dnorm(0,1.0E-4)
beta <- exp(alpha)
logit(theta) <- delta
delta ~ dnorm(0, 1.0E-4)
}
Data
list(x = c(6, 9, 17, 22, 7, 5, 5, 14, 9, 7, 9, 51),
y = c(5, 2, 0, 0, 2, 1, 0, 0, 0, 0, 13, 0), N = 12)
Inits for chain 1
list(delta = 0, alpha = 0)
Inits for chain 2
list(delta = 2, alpha = 2)
Results
![[hearts2]](hearts2.bmp)