# run the CRP Gibbs function in Appendix B.
alpha <- 0.01
mu0 <- matrix(rep(0, 2), ncol = 2, byrow = TRUE)
sigma0 <- diag(2) * 3^2
sigma_y <- diag(2) * 1
c_init <- rep(1, nrow(datc))
results <- crp_gibbs(data = datc, alpha = alpha, mu0 = mu0,
sigma0 = sigma0, sigma_y = sigma_y, c_init = rep(1, nrow(datc)))
#####
# Gibbs sampling iterations are saved in
# `results'. Each row contains the latent class
# membership distribution of one person over
# MCMC iterations. A person is deemed a member
# of c=2 if, e.g., s/he gets assigned most
# often to c2. So we tabulate the frequency
# counts and whichever cluster with highest
# frequency is that person's latent class.
#####
tab <- apply( results, 1, FUN = function(x) {
tab <- table(x)
ans <- names(tab[which.max(tab)])
return(ans) } )
table(tab)