Conditionals from a latent class mixed logit

[jazabala]

Hi everyone,

I have estimated a latent class mixed logit with normal distribution for some of the coefficients, and I would like to calculate: (1) the posterior expected values (conditionals), and their standard deviations, of random coefficients; (2) posterior class allocation probabilities for the individuals. For the first case, I have tried to apply apollo_conditionals (), but it does not work because it is a latent class model. And, for the second case, the function apollo_lcConditionals () can only be used with latent class models without continuous heterogeneity.

> conditionals=apollo_conditionals(model,apollo_probabilities,apollo_inputs)
Updating inputs...Done.
Error in apollo_conditionals(model, apollo_probabilities, apollo_inputs) :
The function 'apollo_conditionals' is not applicables for models containing latent class components!

> lc_conditionals=apollo_lcConditionals(model,apollo_probabilities,apollo_inputs)
Updating inputs...Done.
Error in apollo_lcConditionals(model, apollo_probabilities, apollo_inputs) :
apollo_lcConditionals can only be used for latent class models without continuous random heterogeneity

So, does anyone know how to do it when some coefficients follow a continuous distribution in latent class models?

Thank you.

All the best,
José Ángel

[stephanehess]

Hi Jose

when combining latent class with continuous mixtures, there are several different ways of doing so, especially a distinction between mixing inside the classes and mixing of the class allocation probabilities. A generic function like the apollo_conditionals or apollo_lcConditionals would then not work for some of these specifications.

However, it should not be difficult do this manually in Apollo. Can you describe your model specification to us, maybe by sharing the code. Then I'll make some suggestions.

Best wishes

Stephane

[jazabala]

Hi Stephane

Thanks for your quick reply and help.

Here I share the code:

# ################################################################# #
#### LOAD LIBRARY AND DEFINE CORE SETTINGS ####
# ################################################################# #

### Clear memory
rm(list = ls())

### Load Apollo library
library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
modelName = "LC-MXL-3_mm_M3",
modelDescr = "LC-MXL model_mm_3 classes",
indivID = "ID",
mixing = TRUE,
nCores = 3
)

# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS ####
# ################################################################# #

database = read.csv("Data.csv",header=TRUE, sep = ";")

### Drop protest
database = subset(database,database$protest==0)


# ################################################################# #
#### DEFINE MODEL PARAMETERS ####
# ################################################################# #

### Vector of parameters, including any that are kept fixed in estimation
apollo_beta=c(b_alto_mu_a = 0.02,
b_alto_sig_a = 0.25,
b_alto_mu_b = 0.01,
b_alto_sig_b = 0.01,
b_alto_mu_c = 0.30,
b_alto_sig_c = -0.56,
b_sal_mu_a = 0.00,
b_sal_sig_a = 0.00,
b_sal_mu_b = 0.01,
b_sal_sig_b = 0.01,
b_sal_mu_c = -0.21,
b_sal_sig_c = 0.02,
b_set_mu_a = 0.42,
b_set_sig_a = 0.98,
b_set_mu_b = 0.23,
b_set_sig_b = 0.78,
b_set_mu_c = 1,
b_set_sig_c = 0.12,
b_bpa_mu_a = 0.03,
b_bpa_sig_a = 0.03,
b_bpa_mu_b = 0.15,
b_bpa_sig_b = 0.02,
b_bpa_mu_c = 0.54,
b_bpa_sig_c = 0.08,
b_cost = -0.02,
b_sq_mu_a = -0.02,
b_sq_sig_a = 0.07,
b_sq_mu_b = -0.01,
b_sq_sig_b = -0.02,
b_sq_mu_c = -0.36,
b_sq_sig_c = 0.04,
delta_a = 0.12,
delta_b = 0.23)


### Vector with names (in quotes) of parameters to be kept fixed at their starting value in apollo_beta, use apollo_beta_fixed = c() if none
apollo_fixed = c()

# ################################################################# #
#### DEFINE RANDOM COMPONENTS ####
# ################################################################# #

### Set parameters for generating draws
apollo_draws = list(
interDrawsType = "halton",
interNDraws = 50,
interUnifDraws=c(),
interNormDraws = c("draws_alto_a","draws_sal_a", "draws_set_a", "draws_bpa_a", "draws_sq_a",
"draws_alto_b","draws_sal_b", "draws_set_b", "draws_bpa_b", "draws_sq_b",
"draws_alto_c","draws_sal_c", "draws_set_c", "draws_bpa_c", "draws_sq_c"),

intraDrawsType="mlhs",
intraNDraws=0,
intraUnifDraws=c(),
intraNormDraws=c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
randcoeff = list()

randcoeff[["b_alto_a"]] = b_alto_mu_a + b_alto_sig_a * draws_alto_a
randcoeff[["b_alto_b"]] = b_alto_mu_b + b_alto_sig_b * draws_alto_b
randcoeff[["b_alto_c"]] = b_alto_mu_c + b_alto_sig_c * draws_alto_c
randcoeff[["b_sal_a"]] = b_sal_mu_a + b_sal_sig_a * draws_sal_a
randcoeff[["b_sal_b"]] = b_sal_mu_b + b_sal_sig_b * draws_sal_b
randcoeff[["b_sal_c"]] = b_sal_mu_c + b_sal_sig_c * draws_sal_c
randcoeff[["b_set_a"]] = b_set_mu_a + b_set_sig_a * draws_set_a
randcoeff[["b_set_b"]] = b_set_mu_b + b_set_sig_b * draws_set_b
randcoeff[["b_set_c"]] = b_set_mu_c + b_set_sig_c * draws_set_c
randcoeff[["b_bpa_a"]] = b_bpa_mu_a + b_bpa_sig_a * draws_bpa_a
randcoeff[["b_bpa_b"]] = b_bpa_mu_b + b_bpa_sig_b * draws_bpa_b
randcoeff[["b_bpa_c"]] = b_bpa_mu_c + b_bpa_sig_c * draws_bpa_c
randcoeff[["b_sq_a"]] = b_sq_mu_a + b_sq_sig_a * draws_sq_a
randcoeff[["b_sq_b"]] = b_sq_mu_b + b_sq_sig_b * draws_sq_b
randcoeff[["b_sq_c"]] = b_sq_mu_c + b_sq_sig_c * draws_sq_c

return(randcoeff)
}

# ################################################################# #
#### DEFINE LATENT CLASS COMPONENTS ####
# ################################################################# #

apollo_lcPars = function(apollo_beta, apollo_inputs){
lcpars = list()
lcpars[["b_alto"]] = list(b_alto_a, b_alto_b, b_alto_c)
lcpars[["b_sal"]] = list(b_sal_a, b_sal_b, b_sal_c)
lcpars[["b_set"]] = list(b_set_a, b_set_b, b_set_c)
lcpars[["b_bpa"]] = list(b_bpa_a, b_bpa_b, b_bpa_c)
lcpars[["b_sq"]] = list(b_sq_a, b_sq_b, b_sq_c)

V=list()
V[["class_a"]] = delta_a
V[["class_b"]] = delta_b
V[["class_c"]] = 0

mnl_settings = list(
alternatives = c(class_a=1, class_b=2, class_c=3),
avail = 1,
choiceVar = NA,
V = V
)
lcpars[["pi_values"]] = apollo_mnl(mnl_settings, functionality="raw")

lcpars[["pi_values"]] = apollo_firstRow(lcpars[["pi_values"]], apollo_inputs)

return(lcpars)
}


# ################################################################# #
#### GROUP AND VALIDATE INPUTS ####
# ################################################################# #

apollo_inputs = apollo_validateInputs()

# ################################################################# #
#### DEFINE MODEL AND LIKELIHOOD FUNCTION ####
# ################################################################# #

apollo_probabilities=function(apollo_beta, apollo_inputs, functionality="estimate"){

### Function initialisation: do not change the following three commands
### Attach inputs and detach after function exit
apollo_attach(apollo_beta, apollo_inputs)
on.exit(apollo_detach(apollo_beta, apollo_inputs))

### Create list of probabilities P
P = list()

### Define settings for MNL model component
mnl_settings = list(
alternatives = c(Alt1=1, Alt2=2, Alt3=0),
avail = list(Alt1=1, Alt2=1, Alt3=1),
choiceVar = choice1
)

### Loop over classes
s=1
while(s<=3){

### Compute class-specific utilities
V=list()
V[['Alt1']] = b_alto[[s]] * (Alt1_x1 == 2) + b_sal[[s]] * (Alt1_x2 == 1) + b_set[[s]] * (Alt1_x2 == 3) + b_bpa[[s]] * (Alt1_x2 == 4) + b_cost * Alt1_x3
V[['Alt2']] = b_alto[[s]] * (Alt2_x1 == 2) + b_sal[[s]] * (Alt2_x2 == 1) + b_set[[s]] * (Alt2_x2 == 3) + b_bpa[[s]] * (Alt2_x2 == 4) + b_cost * Alt2_x3
V[['Alt3']] = b_sq[[s]]

mnl_settings$V = V

### Compute within-class choice probabilities using MNL model
P[[s]] = apollo_mnl(mnl_settings, functionality)

### Take product across observation for same individual
P[[s]] = apollo_panelProd(P[[s]], apollo_inputs ,functionality)

### Average across inter-individual draws within classes
P[[s]] = apollo_avgInterDraws(P[[s]], apollo_inputs, functionality)

s=s+1
}

### Compute latent class model probabilities
lc_settings = list(inClassProb = P, classProb=pi_values)
P[["model"]] = apollo_lc(lc_settings, apollo_inputs, functionality)

### Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}


# ################################################################# #
#### MODEL ESTIMATION ####
# ################################################################# #

model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs, estimate_settings=list(hessianRoutine="maxLik"))
_________________________________________________________________________________________________________________________

All the best,
José Ángel

[stephanehess]

Hi Jose

one initial comment. You have 15 continuous random parameters in your model, but you're using Halton draws. Please try to use some other type of draws. Beyond say 5 dimensions, Halton draws are so highly colinear that they should not be used.

In terms of what you need, I've writen the below code for you. It's not necessarily very intuitive for you to follow, but what it does is to calculate the conditionals for the latent classes first (matrix called lc_conditionals), and then the conditionals for the continuous parameters (list called conditionals). I needed to define a separate apollo_probabilities function for the within class model as the conditionals for the continuous parameters need to be calculated inside the classes where they are used given that your integration is at the class-level.

Hope this helps

Stephane

Code: Select all

model=apollo_loadModel("LC-MXL-3_mm_M3")

# ################################################################# #
#### CONDITIONALS FOR CLASS ALLOCATION                           ####
# ################################################################# #

apollo_beta  = model$estimate
class_prob     = "pi_values" 
L      = apollo_probabilities(apollo_beta, apollo_inputs, functionality="output")
apollo_attach(apollo_beta, apollo_inputs)
lcpars = apollo_lcPars(apollo_beta, apollo_inputs)
classes    = length(lcpars[[class_prob]])
post_pi = vector(mode="list", length=classes)
for(s in 1:classes) post_pi[[s]] = lcpars[[class_prob]][[s]]*L[[s]]/L[["model"]]
lc_conditionals = matrix(unlist(post_pi), ncol = length(post_pi), byrow = FALSE)
classnames   = paste("Class ",seq(1:classes),sep="")
lc_conditionals=cbind(unique(database[,apollo_control$indivID]),lc_conditionals)
colnames(lc_conditionals) = c("ID",classnames)

apollo_detach(apollo_beta, apollo_inputs)


# ################################################################# #
#### CONDITIONALS FOR CONTINUOUS PARAMS                          ####
# ################################################################# #

apollo_probabilities_class_specific=function(apollo_beta, apollo_inputs, functionality="estimate"){
  
  ### Function initialisation: do not change the following three commands
  ### Attach inputs and detach after function exit
  apollo_attach(apollo_beta, apollo_inputs)
  on.exit(apollo_detach(apollo_beta, apollo_inputs))
  
  ### Create list of probabilities P
  P = list()
  
  ### Define settings for MNL model component
  mnl_settings = list(
    alternatives  = c(Alt1=1, Alt2=2, Alt3=0),
    avail         = list(Alt1=1, Alt2=1, Alt3=1),
    choiceVar     = choice1
  )
  
  s=apollo_inputs$s
  
    ### Compute class-specific utilities
    V=list()
    V[['Alt1']]  = b_alto[[s]] * (Alt1_x1 == 2) +  b_sal[[s]] * (Alt1_x2 == 1) + b_set[[s]] * (Alt1_x2 == 3) + b_bpa[[s]] * (Alt1_x2 == 4) + b_cost * Alt1_x3 
    V[['Alt2']]  = b_alto[[s]] * (Alt2_x1 == 2) +  b_sal[[s]] * (Alt2_x2 == 1) + b_set[[s]] * (Alt2_x2 == 3) + b_bpa[[s]] * (Alt2_x2 == 4) + b_cost * Alt2_x3  
    V[['Alt3']]  = b_sq[[s]]
    
    mnl_settings$V = V
    
    ### Compute within-class choice probabilities using MNL model
    P[["model"]] = apollo_mnl(mnl_settings, functionality)
    
    ### Take product across observation for same individual
    P = apollo_panelProd(P, apollo_inputs ,functionality)
    
    ### Average across inter-individual draws within classes
    P = apollo_avgInterDraws(P, apollo_inputs, functionality)
    
  ### Prepare and return outputs of function
  P = apollo_prepareProb(P, apollo_inputs, functionality)
  return(P)
}

apollo_beta=model$estimate
P=list()
for(s in 1:3){
apollo_inputs$s=s
P[[s]]=apollo_probabilities_class_specific(apollo_beta, apollo_inputs, functionality="conditionals")
}

apollo_attach(apollo_beta, apollo_inputs)
randcoeff <- apollo_randCoeff(apollo_beta, apollo_inputs)
obsPerIndiv <- as.vector(table(database[,apollo_control$indivID]))
conditionals=list()
class1_pars=c("b_alto_a",
              "b_sal_a",
              "b_set_a",
              "b_bpa_a", 
              "b_sq_a")
class2_pars=c("b_alto_b",
              "b_sal_b",
              "b_set_b",
              "b_bpa_b", 
              "b_sq_b")
class3_pars=c("b_alto_c",
              "b_sal_c",
              "b_set_c",
              "b_bpa_c", 
              "b_sq_c")

for(j in 1:length(randcoeff)){
  b=randcoeff[[j]]
  b <- rowsum(b, group=database[,apollo_control$indivID])
  b=b/obsPerIndiv
  P_use=(names(randcoeff)[j]%in%class1_pars)*P[[1]]+(names(randcoeff)[j]%in%class2_pars)*P[[2]]+(names(randcoeff)[j]%in%class3_pars)*P[[3]]
  bn=(rowSums(b*P_use))/rowSums(P_use)
  bns=sqrt(rowSums(P_use*(b-bn)^2)/(rowSums(P_use)))
  conditionals[[names(randcoeff)[j]]]=cbind(unique(database[,apollo_control$indivID]),bn,bns)
  colnames(conditionals[[names(randcoeff)[j]]])=c("ID","post. mean","post. sd")
  rownames(conditionals[[names(randcoeff)[j]]])=c()
}

apollo_detach(apollo_beta, apollo_inputs)

[jazabala]

Hi Stephane

Ok, I'll change the type of draws, and I'll try again using MLHS.

It helps so much. I have just run the code and it works. It is what I need.

Thanks so much for your help.

All the best,
José