Can someone help with JAGS code for a Bayesian multinomial logistic model with one categorical X variable (Dirichlet prior)? My representative sample is the matrix "z" in the code below that represents the 3 outcomes and "site", in the bottom line of code, is the categorical x variable.
I can get the code that estimates each of the 3 outcomes but I'm stumped on how to add a categorical X (hospital site).
I would like to use the first outcome, z[, 1], as the reference and 'a' as the reference for 'site'.
Here is example code that estimates the outcomes (NO categorical X). This is what I have so far. Any help on extending this model with X would be appreciated.
library('rjags')
z <- matrix(c(rep(1,70), rep(0,30),
rep(0,70), rep(1,22), rep(0,8),
rep(0,92), rep(1,8)),
nrow=100, ncol=3)
## The model ##
modelString = "
model
{
for (j in 1:K)
{
alpha[j] <- 1
}
#Prior
pi ~ ddirch(alpha[1:K])
for (i in 1:N)
{
z[i, 1:K] ~ dmulti(pi, 1)
}
}
"
writeLines( modelString , con="TEMPmodel.txt" )
#Run jags
jags <- jags.model('TEMPmodel.txt',
data = list('z' = z,
'N' = nrow(z),
'K' = ncol(z)),
n.chains = 4,
n.adapt = 1000)
mcmc.samples <- coda.samples(jags,
c('pi'),
2500)
#Estimates are similar to observed data
summary(mcmc.samples)
#5 category predictor hospital site to add to model
set.seed(1)
site <- sample(rep(letters[1:5], 20), 100)
x
in your sample data. Can you please edit your question to include representative sample data? – Maurits Everssite
is supposed to do/denote. Perhaps I'm misunderstanding. Generally when including covariates in a multinomial logistic model (it doesn't matter whether they are numeric or categorical, in case of the latter we can simply introduce index covariates) you'd need a measurement for every outcome of your multinomial response. So how do I know with which outcomesite
is associated? – Maurits Evers