Identification issues with Latent Class Analysis/Mixture model

Modeling
, ,
17

@Kevin_Linares i did a quick test, and if you use the modified code shared by @spinkney with polytomous items its works as long as you separate them into all possible binary/dummy variables

Something like this example

data {
  int<lower=1> I;                                     // number of items
  int<lower=1> J;                                     // number of respondents
  matrix[J, I] y;                                     // score for obs n
  int<lower=1> C;           // # of attribute profiles (latent classes) 
}
parameters {
  simplex[C] nu;
  // average intercept and main effect
  real mean_intercept;
  // item level deviations from average effects
  vector[I] dev_intercept;
  vector<lower=-7, upper=7>[I] class_interaction[C];
  ordered[C] class_effect_non;
}

transformed parameters {
  matrix[I, C] log_prob;    // Probability of correct response for each class on each item
  // Probability of correct response for each class on each item
  for (c in 1:C) 
     log_prob[, c] = log_inv_logit(mean_intercept + dev_intercept + class_interaction[c] + class_effect_non[c]);
}
model{
  // Priors
  mean_intercept ~ normal(0, 5);
  
  for (c in 1:C)
    class_interaction[c] ~ std_normal();
  
  class_effect_non ~ normal(0, 3);
  dev_intercept ~ normal(0, 3);
 
  {
    for (j in 1:J) {
      row_vector[C] ps = log(nu') + y[j] * log_prob + (1 - y[j]) *  log1m_exp(log_prob);
      target += log_sum_exp(ps);
    }
  }
}
generated quantities {
  matrix[I, C] probs;
  
  for (c in 1:C)
    probs[, c] = inv_logit(mean_intercept + dev_intercept + class_interaction[c] + class_effect_non[c]);
}

library(psych)
library(poLCA)
data(election)


f2a <- cbind(MORALG,CARESG,KNOWG)~1
nes2a <- poLCA(f2a,election,nclass=2,nrep=10,  maxiter=20000)
nes2a


d2 <- election[,c("MORALG","CARESG","KNOWG")]
head(d2)
apply(d2, 2, table)
summary(d2)

table(d2[,1],dummy.code(d2[,1])[,1] )

d3 <- cbind(dummy.code(d2[,1]),
            dummy.code(d2[,2]),
            dummy.code(d2[,3]) )

head(d3)
d4 <- na.omit(d3)
dim(d4)
head(d4)

wide_data <- list(
  I = ncol(d4),
  J = nrow(d4),
  y = d4,
  C = 2
)

pars <- c("nu","probs")
iter <- 1000
warmup <- iter - 500
wide_model <- stan(model_code=lca_wide, 
                   data = wide_data,
                   pars=pars,
                   chains = 3, 
                   iter = iter,
                   warmup = warmup)
wide_model
2 Likes