13
$\begingroup$

I'm trying to fit a hierarchical model using jags, and the rjags package. My outcome variable is y, which is a sequence of bernoulli trials. I have 38 human subjects which are performing under two categories: P and M. Based on my analysis, every speaker has a probability of success in category P of $\theta_p$ and a probability of success in category M of $\theta_p\times\theta_m$. I'm also assuming that there is some community level hyperparameter for P and M: $\mu_p$ and $\mu_m$.

So, for every speaker: $\theta_p \sim beta(\mu_p\times\kappa_p, (1-\mu_p)\times\kappa_p)$ and $\theta_m \sim beta(\mu_m\times\kappa_m, (1-\mu_m)\times\kappa_m)$ where $\kappa_p$ and $\kappa_m$ control how peaked the distribution is around $\mu_p$ and $\mu_m$.

Also $\mu_p \sim beta(A_p, B_p)$, $\mu_m \sim beta(A_m, B_m)$.

Here's my jags model:

model{
## y = N bernoulli trials
## Each speaker has a theta value for each category
for(i in 1:length(y)){
    y[i] ~ dbern( theta[ speaker[i],category[i]])
}

## Category P has theta Ptheta
## Category M has theta Ptheta * Mtheta
## No observed data for pure Mtheta
##
## Kp and Km represent how similar speakers are to each other 
## for Ptheta and Mtheta
for(j in 1:max(speaker)){
    theta[j,1] ~ dbeta(Pmu*Kp, (1-Pmu)*Kp)
    catM[j] ~ dbeta(Mmu*Km, (1-Mmu)*Km)
    theta[j,2] <- theta[j,1] * catM[j]
}

## Priors for Pmu and Mmu
Pmu ~ dbeta(Ap,Bp)
Mmu ~ dbeta(Am,Bm)

## Priors for Kp and Km
Kp ~ dgamma(1,1/50)
Km ~ dgamma(1,1/50)

## Hyperpriors for Pmu and Mmu
Ap ~ dgamma(1,1/50)
Bp ~ dgamma(1,1/50)
Am ~ dgamma(1,1/50)
Bm ~ dgamma(1,1/50)
}

The issue I have is that when I run this model with 5000 iterations for adapting, then take 1000 samples, Mmu and Km have converged to single values. I've been running it with 4 chains, and each chain doesn't have the same value, but within each chain there is just a single value.

I'm pretty new to fitting hierarchical models using MCMC methods, so I'm wondering how bad this is. Should I take this as a sign that this model is hopeless to fit, that something is wrong with my priors, or is this par for the course?

Edit: In case it matters, the value for $\mu_m$ it converged to (averaged across chains) was 0.91 and $\kappa_m$ was 1.78

Improve this question
$\endgroup$
9
  • $\begingroup$ If I'm understanding you correctly, these parameters "converge" on one fixed value in each chain (after some iterations it doesn't change at all), but that value is different for every chain you run? That sounds bad, like maybe a really crappy Metropolis Hastings step. It could be your model, it could be JAGS, it could be a combination of both. Presumably this model doesn't take too long to fit, so I'd try running (much) longer chains first, especially for the adapting period. $\endgroup$
    – JMS
    Commented Jun 13, 2011 at 18:11
  • $\begingroup$ So, I updated the model with 5000 more iterations, and the parameters in question didn't budge. I didn't realize they could fall into local minima like this. $\endgroup$
    – JoFrhwld
    Commented Jun 13, 2011 at 19:36
  • 1
    $\begingroup$ quick sugesstions:1. Try to use dbin, with n=1. And use bounds to limit the value of p. Something like this: p.bound[i] <- max(0, min(1, p[i])) $\endgroup$
    – Manoel Galdino
    Commented Jun 13, 2011 at 20:25
  • 1
    $\begingroup$ A couple of clarifying questions: 1. you have 38 subjects under category P and 38 subjects under category M, such as length(y) = 76? 2. Could you give more background information about the rationale for the hyperparamters and the experiment? It's a bit confusing to me. $\endgroup$
    – Manoel Galdino
    Commented Jun 13, 2011 at 20:35
  • 1
    $\begingroup$ I'd probably just fix all the hyperparameters, as in theta[j,1] ~ dbeta(1.1, 1.1) or something, and see what results you get before trying to move on to a hyperprior. Also, theta[j,2]'s prior is the product of two betas, which is not, in general, a beta itself, and of course theta[j,2] < theta[j,1] as a result. It looks like you intended this; just out of curiosity, why? $\endgroup$
    – jbowman
    Commented Dec 4, 2011 at 3:07

2 Answers 2

Reset to default
2
$\begingroup$

This is more a comment, but as I do not have enough reputation I might as well answer.

From my limited experience with MCMC samplers, what I have observed is that the parameters tend to stay fixed when the hyperparameters are too narrow. As they control the spreading of the parameters, they prevent the solution space to be efficiently sampled.

Try to take larger values for the hyperparameters and see what happens.

This technical paper helped me a lot to understand the MCMC samplers. It is composed of two samplers, Gibbs (the one you are using), and Hybrid Monte Carlo, and rapidly explaining how to choose priors, hyperpriors and values for the parameters and hyperparameters.

Improve this answer
$\endgroup$
0
$\begingroup$

This could be a problem of the structure of the chain. Where you end up depends on where you start. To use MCMC you want the chain to be recurrent which means that no matter where you start you can get to every other state in the state space. If the chain is not recurrent you can be trapped in a subset of the state space. The idea of MCMC is to have an existing stationary distribution that the chain will eventually wind up in. This stationary distribution usually has a positive probability for being in any of the states in the chain and not wind up trapped at a single point as you described. I can't check your algorithm but maybe you have a mistake in it. It is also possible that you have defined a problem where your Markov chain is not recurrent.

If you want to get knowledgeable in MCMC I recommend that you take a look at the Handbook of Markov Chain Monte Carlo which has articles that describe every aspect of MCMC. It was published by CRC Press in 2011.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.