Trouble turning Bayesian finite mixture model into hierarchical BFMM

specification

#1

I want to fit a mixture model hierarchically but have hit a wall. Currently, I have the parameters of a Dirichlet distribution describing the probability for a model. The 4 clusters that are thought to compose my data are 3 von Mises distributions (one centered at mu, one centered at swapLoc, and one centered at swapCLoc), and 1 uniform guessing distribution. The problem I have is that I want to get both individual-level and group-level estimates for these mixing proportions. The below code is a working non-hierarchical mixture model (the transformed parameters section I took from here). Note that all the data is is circular (in radians), which is why I use von Mises distributions.

I assume to implement hierarchy, I need to have different p[c]'s for each individual and v needs to be v ~ beta(a,b) instead of beta(1,1), where a and b are drawn from some larger distribution. But turning p[c] into a matrix like p[i,c] has presented me with a flurry of errors, and I’m not sure how to specify this larger distribution for v ~ beta(a,b) (if that’s even what I need to do to get this working). I am at a loss for how to proceed. All help appreciated.


data {
  int NumSubj; //44
  int NumData; //9635
  int C; //num of clusters
  vector[NumData] ID; // subject ID (1-44)
  int NumResp[NumSubj]; // subject's num of responses
  int CumResp[NumSubj]; // cumulative num of responses
  vector[NumData] errors;
  vector[NumData] swapLoc;
  vector[NumData] swapCLoc;
}
parameters {
  real<lower=-pi()/6,upper=pi()/6> mu;
  real<lower=0,upper=100> kappa;
  real<lower=0,upper=100> kappa2;
  real<lower=0,upper=1> v[C];
}
transformed parameters {
  simplex[C] p; // mixing proportions for each of C clusters (sums to 1)
  p[1] = v[1];
  // stick-breaking Dirichlet Process based on The BUGS book Chapter 11 (p.294)
  for (l in 2:(C-1)) {
    p[l]= v[l]*(1-v[l-1])*p[l-1]/v[l-1]; 
  }
  p[C]=1-sum(p[1:(C-1)]); // to make a simplex
}
model {
  real ps[C]; 
  v ~ beta(1,1);
  
  kappa ~ normal(10,100);
  kappa2 ~ normal(10,100);
  mu ~ normal(0,5);

  for (i in 1:NumSubj) {
    for (j in 1:NumResp[i]) {
      for (c in 1:C) {
        if (c == 1)       
          ps[c] = log(p[c])+von_mises_lpdf(errors[CumResp[i]-(j-1)] | swapLoc[CumResp[i]-(j-1)], kappa2); // Swap
        else if (c == 2)  
          ps[c] = log(p[c])+von_mises_lpdf(errors[CumResp[i]-(j-1)] | swapCLoc[CumResp[i]-(j-1)], kappa2); // Swap-Comparison
        else if (c == 3)  
          ps[c] = log(p[c])+uniform_lpdf(errors[CumResp[i]-(j-1)] | -pi(), pi()); // Guess
        else              
          ps[c] = log(p[c])+von_mises_lpdf(errors[CumResp[i]-(j-1)] | mu, kappa); // Target
        }
      target += log_sum_exp(ps);
    }
  }
}