data{ int K; //number of components int N; //number of samples int x[N]; //trial for each sample int y[N]; //success for each sample vector[K] phi; //hyperparameter of mixing ratio real alpha; //hyperparameter of component attached with p real beta; //hyperparameter of component attached with 1-p } parameters{ vector[K] p; simplex[K] ratio; // mixing ratio } model{ vector[K] partial_term; for(k in 1:K){ p[k] ~ beta(alpha, beta); } ratio ~ dirichlet(phi); for(i in 1:N){ for(k in 1:K){ partial_term[k] <- log(ratio[k]) + binomial_log(y[i], x[i], p[k]); } increment_log_prob(log_sum_exp(partial_term)); } }