Important: Read this before posting to this forum

  1. This forum is for questions related to the use of Apollo. We will answer some general choice modelling questions too, where appropriate, and time permitting. We cannot answer questions about how to estimate choice models with other software packages.
  2. There is a very detailed manual for Apollo available at http://www.ApolloChoiceModelling.com/manual.html. This contains detailed descriptions of the various Apollo functions, and numerous examples are available at http://www.ApolloChoiceModelling.com/examples.html. In addition, help files are available for all functions, using e.g. ?apollo_mnl
  3. Before asking a question on the forum, users are kindly requested to follow these steps:
    1. Check that the same issue has not already been addressed in the forum - there is a search tool.
    2. Ensure that the correct syntax has been used. For any function, detailed instructions are available directly in Apollo, e.g. by using ?apollo_mnl for apollo_mnl
    3. Check the frequently asked questions section on the Apollo website, which discusses some common issues/failures. Please see http://www.apollochoicemodelling.com/faq.html
    4. Make sure that R is using the latest official release of Apollo.
  4. If the above steps do not resolve the issue, then users should follow these steps when posting a question:
    1. provide full details on the issue, including the entire code and output, including any error messages
    2. posts will not immediately appear on the forum, but will be checked by a moderator first. We check the forum at least twice a week. It may thus take a couple of days for your post to appear and before we reply. There is no need to submit the post multiple times.

Latent Class model with EM algorithm

Ask questions about how to estimate models and how to change your settings for estimation.
Post Reply
Andrew
Posts: 3
Joined: 10 Mar 2025, 14:56

Latent Class model with EM algorithm

Post by Andrew »

Dear team and Apollo users,

I have only recently started working with Apollo. I also use Stata. I was able to estimate the first models and also reproduce an LC model that I had previously estimated with Stata. Below is the code of my LC model. My dataset is based on a traditional DCE with 2 unlabeled choice alternatives, 4 attributes x1 to x4 (linear-coded in the dataset), no ASC, no opt-out.

Code: Select all

# ################################################################# #
#### INITIALISATION                                              ####
# ################################################################# #

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

### Load Apollo library
library(apollo)

### Initialise Apollo
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "LC_DCE_Model_linear",
  modelDescr      = "LC model for DCE data",
  indivID         = "respid",
  outputDirectory = "output_LC_linear"
)

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

### Load data from CSV/Excel
library(readxl)
database = read_excel("data.xlsx")

### Reshaping the data set from long (one alt per row) to wide (one obs per row)
settings<-list(
  alternative_column="alt",
  alternative_specific_attributes= c("x1", "x2", "x3", "x4"),
  choice_column="choice",
  ID_column="respid",
  observation_column="setid"
)

database = apollo_longToWide(database, longToWide_settings = settings)

# Recoding of column choice_new as vector
database$choice_new <- as.vector(database$choice_new)

# ################################################################# #
#### MODEL PARAMETERS                                            ####
# ################################################################# #

### Define model parameters
apollo_beta = c(
  beta_x1_a = 0,
  beta_x1_b = 0,
  beta_x2_a = 0,
  beta_x2_b = 0,
  beta_x3_a = 0,
  beta_x3_b = 0,
  beta_x4_a = 0,
  beta_x4_b = 0,
  delta_a = 0,
  delta_b = 0
)

### No fixed parameters
apollo_fixed = c("delta_b")

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

apollo_lcPars = function(apollo_beta, apollo_inputs) {
  
  lcpars = list()
  lcpars[["beta_x1"]] = list(beta_x1_a, beta_x1_b)
  lcpars[["beta_x2"]] = list(beta_x2_a, beta_x2_b)
  lcpars[["beta_x3"]] = list(beta_x3_a, beta_x3_b)
  lcpars[["beta_x4"]] = list(beta_x4_a, beta_x4_b)

  ### Utilities of class allocation model
  V = list()
  V[["class_a"]] = delta_a
  V[["class_b"]] = delta_b
  
  ### Settings for class allocation models
  classAlloc_settings = list(
    classes = c(class_a=1, class_b=2) ,
    utilities = V
  )
  
  lcpars[["pi_values"]] = apollo_classAlloc(classAlloc_settings)
  
  return(lcpars)
}

# ################################################################# #
#### RUN apollo_validateInputs FUNCTION          VALIDATE INPUTS ####
# ################################################################# #

apollo_inputs = apollo_validateInputs()

# ################################################################# #
#### MODEL DEFINITION                                            ####
# ################################################################# #

apollo_probabilities = function(apollo_beta, apollo_inputs, functionality = "estimate") {
  
  ### 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 that are generic across classes
  mnl_settings = list(
    alternatives = c(alt1 = 1, alt2 = 2),
    avail = list(alt1 = avail_1, alt2 = avail_2),
    choiceVar = choice_new
  )
    
  ### Loop over classes
  for (s in 1:2) {    
    
    ### Define utilities for each alternative (remove "database$" prefix)
    V = list()
    V[["alt1"]] = beta_x1[[s]]*x1_1 + beta_x2[[s]]*x2_1 + beta_x3[[s]]*x3_1 + beta_x4[[s]]*x4_1
    V[["alt2"]] = beta_x1[[s]]*x1_2 + beta_x2[[s]]*x2_2 + beta_x3[[s]]*x3_2 + beta_x4[[s]]*x4_2
    
    mnl_settings$utilities     = V
    mnl_settings$componentName = paste0("Class_",s)
    
    ### Compute within-class choice probabilities using MNL model
    P[[paste0("Class_", s)]] = apollo_mnl(mnl_settings, functionality)
    
    ### Take product across observation for same individual
    P[[paste0("Class_", s)]] = apollo_panelProd(P[[paste0("Class_", s)]], apollo_inputs, functionality)
  }
  
  ### Compute latent class model probabilities
  lc_settings = list(inClassProb = P, classProb = pi_values)
  
  ### Compute latent class model output
  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)

# ################################################################# #
#### MODEL OUTPUTS                                               ####
# ################################################################# #

### Print results
apollo_modelOutput(model)

### Save results
apollo_saveOutput(model)

Here the output. The result matches the Stata LC model which was estimated via EM algorithm.

Code: Select all

Model name                                  : LC_DCE_Model_linear
Model description                           : LC model for DCE data
Model run at                                : 2025-03-11 11:38:41.380007
Estimation method                           : bgw
Model diagnosis                             : Relative function convergence
Optimisation diagnosis                      : Maximum found
     hessian properties                     : Negative definite
     maximum eigenvalue                     : -56.994819
     reciprocal of condition number         : 0.0325824
Number of individuals                       : 300
Number of rows in database                  : 4500
Number of modelled outcomes                 : 4500

Number of cores used                        :  1 
Model without mixing

LL(start)                                   : -3119.16
LL (whole model) at equal shares, LL(0)     : -3119.16
LL (whole model) at observed shares, LL(C)  : -3115.64
LL(final, whole model)                      : -2548.56
Rho-squared vs equal shares                  :  0.1829 
Adj.Rho-squared vs equal shares              :  0.18 
Rho-squared vs observed shares               :  0.182 
Adj.Rho-squared vs observed shares           :  0.1798 
AIC                                         :  5115.12 
BIC                                         :  5172.83 

LL(0,Class_1)                    : -3119.16
LL(final,Class_1)                : -3557.47
LL(0,Class_2)                    : -3119.16
LL(final,Class_2)                : -3215.11

Estimated parameters                        : 9
Time taken (hh:mm:ss)                       :  00:00:3.14 
     pre-estimation                         :  00:00:0.44 
     estimation                             :  00:00:1.03 
     post-estimation                        :  00:00:1.67 
Iterations                                  :  19  

Unconstrained optimisation.

Estimates:
             Estimate        s.e.   t.rat.(0)    Rob.s.e. Rob.t.rat.(0)
beta_x1_a     0.40473     0.04095       9.883     0.04292         9.430
beta_x1_b    -0.16036     0.03718      -4.313     0.05585        -2.871
beta_x2_a     0.44770     0.04723       9.479     0.05483         8.165
beta_x2_b    -0.90816     0.04320     -21.020     0.04438       -20.462
beta_x3_a    -0.19622     0.04148      -4.730     0.04706        -4.169
beta_x3_b    -0.59027     0.03628     -16.270     0.02992       -19.731
beta_x4_a    -0.09546     0.02897      -3.296     0.02833        -3.370
beta_x4_b    -0.38496     0.02688     -14.323     0.03007       -12.804
delta_a      -0.57324     0.13137      -4.364     0.13511        -4.243
delta_b       0.00000          NA          NA          NA            NA


Summary of class allocation for model component :
         Mean prob.
Class_1      0.3605
Class_2      0.6395


Overview of choices for MNL model component Class_1:
                                    alt1    alt2
Times available                  4500.00 4500.00
Times chosen                     2339.00 2161.00
Percentage chosen overall          51.98   48.02
Percentage chosen when available   51.98   48.02



Overview of choices for MNL model component Class_2:
                                    alt1    alt2
Times available                  4500.00 4500.00
Times chosen                     2339.00 2161.00
Percentage chosen overall          51.98   48.02
Percentage chosen when available   51.98   48.02



Classical covariance matrix:
            beta_x1_a   beta_x1_b   beta_x2_a   beta_x2_b   beta_x3_a   beta_x3_b   beta_x4_a   beta_x4_b     delta_a
beta_x1_a    0.001677   6.892e-05  4.0294e-04  1.8008e-04   3.485e-07   4.745e-05  -1.847e-05   6.571e-05 -4.4254e-04
beta_x1_b   6.892e-05    0.001383  3.6344e-04  4.5968e-04  1.6278e-04  1.3440e-04   5.327e-05  1.1449e-04 -7.5228e-04
beta_x2_a  4.0294e-04  3.6344e-04    0.002231  4.4771e-04  1.5542e-04  1.1316e-04  -1.900e-05  1.9090e-04   -0.001290
beta_x2_b  1.8008e-04  4.5968e-04  4.4771e-04    0.001867  2.1360e-04  5.9414e-04   6.376e-05  4.1918e-04   -0.001054
beta_x3_a   3.485e-07  1.6278e-04  1.5542e-04  2.1360e-04    0.001721   5.139e-06   5.233e-05   9.582e-05 -5.9470e-04
beta_x3_b   4.745e-05  1.3440e-04  1.1316e-04  5.9414e-04   5.139e-06    0.001316   2.433e-05  2.3959e-04 -2.5699e-04
beta_x4_a  -1.847e-05   5.327e-05  -1.900e-05   6.376e-05   5.233e-05   2.433e-05  8.3912e-04  -1.740e-05 -1.5977e-04
beta_x4_b   6.571e-05  1.1449e-04  1.9090e-04  4.1918e-04   9.582e-05  2.3959e-04  -1.740e-05  7.2238e-04 -3.9763e-04
delta_a   -4.4254e-04 -7.5228e-04   -0.001290   -0.001054 -5.9470e-04 -2.5699e-04 -1.5977e-04 -3.9763e-04    0.017258

Robust covariance matrix:
            beta_x1_a   beta_x1_b   beta_x2_a   beta_x2_b   beta_x3_a   beta_x3_b   beta_x4_a   beta_x4_b     delta_a
beta_x1_a    0.001842  1.9789e-04  5.9266e-04  2.5317e-04   3.986e-05   4.099e-05   6.373e-05  1.1930e-04 -6.3217e-04
beta_x1_b  1.9789e-04    0.003119  9.8382e-04  9.2295e-04  5.0050e-04   3.842e-05   4.843e-05  7.3538e-04   -0.001857
beta_x2_a  5.9266e-04  9.8382e-04    0.003006  7.3263e-04  6.1693e-04   8.468e-05 -2.1467e-04  4.5132e-04   -0.002027
beta_x2_b  2.5317e-04  9.2295e-04  7.3263e-04    0.001970  4.0418e-04  3.8064e-04   2.399e-05  6.5701e-04   -0.001438
beta_x3_a   3.986e-05  5.0050e-04  6.1693e-04  4.0418e-04    0.002215   1.130e-05  -3.621e-06  2.4249e-04   -0.001079
beta_x3_b   4.099e-05   3.842e-05   8.468e-05  3.8064e-04   1.130e-05  8.9493e-04   4.530e-06  1.5919e-04 -1.5585e-04
beta_x4_a   6.373e-05   4.843e-05 -2.1467e-04   2.399e-05  -3.621e-06   4.530e-06  8.0259e-04  -1.776e-05  -5.543e-05
beta_x4_b  1.1930e-04  7.3538e-04  4.5132e-04  6.5701e-04  2.4249e-04  1.5919e-04  -1.776e-05  9.0397e-04 -8.3494e-04
delta_a   -6.3217e-04   -0.001857   -0.002027   -0.001438   -0.001079 -1.5585e-04  -5.543e-05 -8.3494e-04    0.018256

Classical correlation matrix:
            beta_x1_a   beta_x1_b   beta_x2_a   beta_x2_b   beta_x3_a   beta_x3_b   beta_x4_a   beta_x4_b     delta_a
beta_x1_a     1.00000     0.04526     0.20833     0.10178  2.0517e-04    0.031936    -0.01557     0.05970    -0.08226
beta_x1_b     0.04526     1.00000     0.20695     0.28613    0.105535    0.099628     0.04946     0.11456    -0.15400
beta_x2_a     0.20833     0.20695     1.00000     0.21941    0.079331    0.066046    -0.01388     0.15038    -0.20793
beta_x2_b     0.10178     0.28613     0.21941     1.00000    0.119188    0.379061     0.05095     0.36098    -0.18568
beta_x3_a  2.0517e-04     0.10554     0.07933     0.11919    1.000000    0.003415     0.04355     0.08595    -0.10913
beta_x3_b     0.03194     0.09963     0.06605     0.37906    0.003415    1.000000     0.02315     0.24572    -0.05392
beta_x4_a    -0.01557     0.04946    -0.01388     0.05095    0.043550    0.023147     1.00000    -0.02235    -0.04198
beta_x4_b     0.05970     0.11456     0.15038     0.36098    0.085948    0.245719    -0.02235     1.00000    -0.11262
delta_a      -0.08226    -0.15400    -0.20793    -0.18568   -0.109132   -0.053922    -0.04198    -0.11262     1.00000

Robust correlation matrix:
            beta_x1_a   beta_x1_b   beta_x2_a   beta_x2_b   beta_x3_a   beta_x3_b   beta_x4_a   beta_x4_b     delta_a
beta_x1_a     1.00000     0.08256     0.25184     0.13291    0.019732    0.031926    0.052411     0.09245    -0.10901
beta_x1_b     0.08256     1.00000     0.32128     0.37236    0.190418    0.022999    0.030611     0.43795    -0.24612
beta_x2_a     0.25184     0.32128     1.00000     0.30106    0.239072    0.051625   -0.138197     0.27377    -0.27365
beta_x2_b     0.13291     0.37236     0.30106     1.00000    0.193501    0.286689    0.019076     0.49237    -0.23980
beta_x3_a     0.01973     0.19042     0.23907     0.19350    1.000000    0.008029   -0.002716     0.17137    -0.16971
beta_x3_b     0.03193     0.02300     0.05163     0.28669    0.008029    1.000000    0.005345     0.17699    -0.03856
beta_x4_a     0.05241     0.03061    -0.13820     0.01908   -0.002716    0.005345    1.000000    -0.02085    -0.01448
beta_x4_b     0.09245     0.43795     0.27377     0.49237    0.171366    0.176993   -0.020847     1.00000    -0.20553
delta_a      -0.10901    -0.24612    -0.27365    -0.23980   -0.169707   -0.038558   -0.014481    -0.20553     1.00000

 20 most extreme outliers in terms of lowest average per choice prediction:
  ID Avg prob per choice
 253           0.3760466
 279           0.4003265
 220           0.4301473
  89           0.4343926
 291           0.4351980
  85           0.4402339
 232           0.4412028
 226           0.4437410
 237           0.4448000
 248           0.4454061
 245           0.4516578
 225           0.4519976
 240           0.4567660
 288           0.4580246
 219           0.4581241
 145           0.4586679
 130           0.4586679
 280           0.4589208
 289           0.4589612
 276           0.4597902

Changes in parameter estimates from starting values:
              Initial    Estimate  Difference
beta_x1_a       0.000     0.40473     0.40473
beta_x1_b       0.000    -0.16036    -0.16036
beta_x2_a       0.000     0.44770     0.44770
beta_x2_b       0.000    -0.90816    -0.90816
beta_x3_a       0.000    -0.19622    -0.19622
beta_x3_b       0.000    -0.59027    -0.59027
beta_x4_a       0.000    -0.09546    -0.09546
beta_x4_b       0.000    -0.38496    -0.38496
delta_a         0.000    -0.57324    -0.57324
delta_b         0.000     0.00000     0.00000

Settings and functions used in model definition:

apollo_control
--------------
                       Value                  
modelName              "LC_DCE_Model_linear"  
modelDescr             "LC model for DCE data"
indivID                "respid"               
outputDirectory        "output_LC_linear/"    
debug                  "FALSE"                
nCores                 "1"                    
workInLogs             "FALSE"                
seed                   "13"                   
mixing                 "FALSE"                
HB                     "FALSE"                
noValidation           "FALSE"                
noDiagnostics          "FALSE"                
calculateLLC           "TRUE"                 
analyticHessian        "FALSE"                
memorySaver            "FALSE"                
panelData              "TRUE"                 
analyticGrad           "TRUE"                 
analyticGrad_manualSet "FALSE"                
overridePanel          "FALSE"                
preventOverridePanel   "FALSE"                
noModification         "FALSE"                

Hessian routines attempted
--------------------------
numerical jacobian of LL analytical gradient

Scaling used in computing Hessian
---------------------------------
               Value
beta_x1_a 0.40473131
beta_x1_b 0.16036290
beta_x2_a 0.44770498
beta_x2_b 0.90815967
beta_x3_a 0.19622323
beta_x3_b 0.59026911
beta_x4_a 0.09546364
beta_x4_b 0.38496044
delta_a   0.57323632


apollo_lcPars
---------------
function(apollo_beta, apollo_inputs) {
  
  lcpars = list()
  lcpars[["beta_x1"]] = list(beta_x1_a, beta_x1_b)
  lcpars[["beta_x2"]] = list(beta_x2_a, beta_x2_b)
  lcpars[["beta_x3"]] = list(beta_x3_a, beta_x3_b)
  lcpars[["beta_x4"]] = list(beta_x4_a, beta_x4_b)

  ### Utilities of class allocation model
  V = list()
  V[["class_a"]] = delta_a
  V[["class_b"]] = delta_b
  
  ### Settings for class allocation models
  classAlloc_settings = list(
    classes = c(class_a=1, class_b=2) ,
    utilities = V
  )
  
  lcpars[["pi_values"]] = apollo_classAlloc(classAlloc_settings)
  
  return(lcpars)
}
<bytecode: 0x0000028f8afe1b78>


apollo_probabilities
----------------------
function(apollo_beta, apollo_inputs, functionality = "estimate") {
  
  ### 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 that are generic across classes
  mnl_settings = list(
    alternatives = c(alt1 = 1, alt2 = 2),
    avail = list(alt1 = avail_1, alt2 = avail_2),
    choiceVar = choice_new
  )
    
  ### Loop over classes
  for (s in 1:2) {    
    
    ### Define utilities for each alternative (remove "database$" prefix)
    V = list()
    V[["alt1"]] = beta_x1[[s]]*x1_1 + beta_x2[[s]]*x2_1 + beta_x3[[s]]*x3_1 + beta_x4[[s]]*x4_1
    V[["alt2"]] = beta_x1[[s]]*x1_2 + beta_x2[[s]]*x2_2 + beta_x3[[s]]*x3_2 + beta_x4[[s]]*x4_2
    
    mnl_settings$utilities     = V
    mnl_settings$componentName = paste0("Class_",s)
    
    ### Compute within-class choice probabilities using MNL model
    P[[paste0("Class_", s)]] = apollo_mnl(mnl_settings, functionality)
    
    ### Take product across observation for same individual
    P[[paste0("Class_", s)]] = apollo_panelProd(P[[paste0("Class_", s)]], apollo_inputs, functionality)
  }
  
  ### Compute latent class model probabilities
  lc_settings = list(inClassProb = P, classProb = pi_values)
  
  ### Compute latent class model output
  P[["model"]] = apollo_lc(lc_settings, apollo_inputs, functionality)
  
  ### Prepare and return outputs of function
  P = apollo_prepareProb(P, apollo_inputs, functionality)
  return(P)
}

About my request. When I try to implement the EM algorithm using apollo_lcEM, I get error messages and a completely wrong result. I used the example of the EM_LC_no_covariates model on the Apollo website as a reference. Below is the code and corresponding output.

Code: Select all

# ################################################################# #
#### INITIALISATION                                              ####
# ################################################################# #

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

### Load Apollo library
library(apollo)

### Initialise Apollo
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "LC_DCE_Model_linear",
  modelDescr      = "LC model for DCE data",
  indivID         = "respid",
  outputDirectory = "output_LC_linear"
)

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

### Load data from CSV/Excel
library(readxl)
database = read_excel("data.xlsx")

### Reshaping the data set from long (one alt per row) to wide (one obs per row)
settings<-list(
  alternative_column="alt",
  alternative_specific_attributes= c("x1", "x2", "x3", "x4"),
  choice_column="choice",
  ID_column="respid",
  observation_column="setid"
)

database = apollo_longToWide(database, longToWide_settings = settings)

# Recoding of column choice_new as vector
database$choice_new <- as.vector(database$choice_new)

### create weights column for use in EM algorithm
database$weights = 1

# ################################################################# #
#### MODEL PARAMETERS                                            ####
# ################################################################# #

### Define model parameters
apollo_beta = c(
  beta_x1_a = 0,
  beta_x1_b = 0,
  beta_x2_a = 0,
  beta_x2_b = 0,
  beta_x3_a = 0,
  beta_x3_b = 0,
  beta_x4_a = 0,
  beta_x4_b = 0,
  delta_a = 0,
  delta_b = 0
  )

### No fixed parameters
apollo_fixed = c("delta_b")

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

apollo_lcPars = function(apollo_beta, apollo_inputs) {
  
  lcpars = list()
  lcpars[["beta_x1"]] = list(beta_x1_a, beta_x1_b)
  lcpars[["beta_x2"]] = list(beta_x2_a, beta_x2_b)
  lcpars[["beta_x3"]] = list(beta_x3_a, beta_x3_b)
  lcpars[["beta_x4"]] = list(beta_x4_a, beta_x4_b)

  lcpars[["pi_values"]] = list(class_a = 1/(1 + exp(delta_b - delta_a)),
                               class_b = 1/(1 + exp(delta_a - delta_b)))

    return(lcpars)
}

# ################################################################# #
#### RUN apollo_validateInputs FUNCTION          VALIDATE INPUTS ####
# ################################################################# #

apollo_inputs = apollo_validateInputs()

# ################################################################# #
#### MODEL DEFINITION                                            ####
# ################################################################# #

apollo_probabilities = function(apollo_beta, apollo_inputs, functionality = "estimate") {
  
  ### 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 that are generic across classes
  mnl_settings = list(
    alternatives = c(alt1 = 1, alt2 = 2),
    avail = list(alt1 = avail_1, alt2 = avail_2),
    choiceVar = choice_new
  )
  
  ### Loop over classes
  for (s in 1:length(pi_values)) {
    
    ### Define utilities for each alternative (remove "database$" prefix)
    V = list()
    V[["alt1"]] = beta_x1[[s]]*x1_1 + beta_x2[[s]]*x2_1 + beta_x3[[s]]*x3_1 + beta_x4[[s]]*x4_1
    V[["alt2"]] = beta_x1[[s]]*x1_2 + beta_x2[[s]]*x2_2 + beta_x3[[s]]*x3_2 + beta_x4[[s]]*x4_2
    
    mnl_settings$utilities     = V
    mnl_settings$componentName = paste0("Class_",s)
    
    ### Compute within-class choice probabilities using MNL model
    P[[paste0("Class_", s)]] = apollo_mnl(mnl_settings, functionality)
    
    ### Take product across observation for same individual
    P[[paste0("Class_", s)]] = apollo_panelProd(P[[paste0("Class_", s)]], apollo_inputs, functionality)
  }
  
  ### 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_lcEM(apollo_beta, apollo_fixed, apollo_probabilities,
                    apollo_inputs,lcEM_settings = list(EMmaxIterations=100))

# ################################################################# #
#### MODEL OUTPUTS                                               ####
# ################################################################# #

### Print results
apollo_modelOutput(model)

### Save results
apollo_saveOutput(model)

And the output. I get a warning that says: Some eigenvalues of the Hessian are positive, which indicates convergence to a saddle point! Warning: In sqrt(diag(varcov)): NaNs were generated.

Code: Select all

Model name                                  : LC_DCE_Model_linear
Model description                           : LC model for DCE data
Model run at                                : 2025-03-11 11:42:46.09308
Estimation method                           : EM algorithm (bfgs) -> Maximum likelihood (bfgs)
Model diagnosis                             : successful convergence
Optimisation diagnosis                      : Saddle point found
     hessian properties                     : Some eigenvalues are positive and others negative
     maximum eigenvalue                     : 2727.476417
     reciprocal of condition number         : not calculated (Hessian is not negative definite)
Number of individuals                       : 300
Number of rows in database                  : 4500
Number of modelled outcomes                 : 4500

Number of cores used                        :  1 
Model without mixing

LL(start)                                   : -3119.16
LL (whole model) at equal shares, LL(0)     : -3119.16
LL (whole model) at observed shares, LL(C)  : -3115.64
LL(final, whole model)                      : -2825.25
Rho-squared vs equal shares                  :  0.0942 
Adj.Rho-squared vs equal shares              :  0.0913 
Rho-squared vs observed shares               :  0.0932 
Adj.Rho-squared vs observed shares           :  0.091 
AIC                                         :  5668.5 
BIC                                         :  5726.21 

LL(0,Class_1)                    : -3119.16
LL(final,Class_1)                : -2825.25
LL(0,Class_2)                    : -3119.16
LL(final,Class_2)                : -2825.25

Estimated parameters                        : 9
Time taken (hh:mm:ss)                       :  00:00:5.36 
     pre-estimation                         :  00:00:0.52 
     estimation                             :  00:00:1.83 
     post-estimation                        :  00:00:3.01 
Iterations                                  :  3 (EM) & 2 (bfgs)  

Unconstrained optimisation.

Estimates:
             Estimate        s.e.   t.rat.(0)    Rob.s.e. Rob.t.rat.(0)
beta_x1_a     0.07806     0.01948       4.008     0.03081         2.534
beta_x1_b     0.07806     0.01948       4.008     0.03083         2.532
beta_x2_a    -0.30321     0.01879     -16.138     0.04055        -7.478
beta_x2_b    -0.30321     0.01879     -16.140     0.04054        -7.479
beta_x3_a    -0.34679     0.04925      -7.042     0.02304       -15.049
beta_x3_b    -0.34679     0.04925      -7.041     0.02307       -15.035
beta_x4_a    -0.21865     0.01648     -13.266     0.01859       -11.762
beta_x4_b    -0.21865     0.01643     -13.305     0.01852       -11.804
delta_a       0.00000         NaN         NaN     3.05278         0.000
delta_b       0.00000          NA          NA          NA            NA


Summary of class allocation for model component :
         Mean prob.
Class_1      0.5000
Class_2      0.5000


Overview of choices for MNL model component Class_1:
                                    alt1    alt2
Times available                  4500.00 4500.00
Times chosen                     2339.00 2161.00
Percentage chosen overall          51.98   48.02
Percentage chosen when available   51.98   48.02



Overview of choices for MNL model component Class_2:
                                    alt1    alt2
Times available                  4500.00 4500.00
Times chosen                     2339.00 2161.00
Percentage chosen overall          51.98   48.02
Percentage chosen when available   51.98   48.02



Classical covariance matrix:
            beta_x1_a   beta_x1_b   beta_x2_a   beta_x2_b   beta_x3_a   beta_x3_b   beta_x4_a   beta_x4_b     delta_a
beta_x1_a  3.7936e-04  6.3457e-04  1.7396e-04 -1.8959e-04 -3.1268e-04  2.6417e-04 -3.4755e-04  3.0262e-04     -67.647
beta_x1_b  6.3457e-04  3.7938e-04 -1.8921e-04  1.7434e-04  2.6448e-04 -3.1236e-04  3.0318e-04 -3.4701e-04      67.646
beta_x2_a  1.7396e-04 -1.8921e-04  3.5300e-04  6.8427e-04 -5.3040e-04  6.2414e-04  -4.640e-06   6.201e-05       8.280
beta_x2_b -1.8959e-04  1.7434e-04  6.8427e-04  3.5291e-04  6.2379e-04 -5.3013e-04   6.249e-05  -5.254e-06      -8.323
beta_x3_a -3.1268e-04  2.6448e-04 -5.3040e-04  6.2379e-04    0.002425   -0.001374  7.6342e-04 -6.9069e-04     -55.122
beta_x3_b  2.6417e-04 -3.1236e-04  6.2414e-04 -5.3013e-04   -0.001374    0.002426 -6.8978e-04  7.6340e-04      55.086
beta_x4_a -3.4755e-04  3.0318e-04  -4.640e-06   6.249e-05  7.6342e-04 -6.8978e-04  2.7167e-04  2.8596e-04      97.663
beta_x4_b  3.0262e-04 -3.4701e-04   6.201e-05  -5.254e-06 -6.9069e-04  7.6340e-04  2.8596e-04  2.7007e-04     -97.725
delta_a        -67.65       67.65       8.280      -8.323  -55.122250   55.086089       97.66      -97.73  -7.621e+06

Robust covariance matrix:
            beta_x1_a   beta_x1_b   beta_x2_a   beta_x2_b   beta_x3_a   beta_x3_b   beta_x4_a   beta_x4_b     delta_a
beta_x1_a  9.4914e-04  9.4989e-04  3.3156e-04  3.3146e-04   6.187e-05   6.248e-05  2.1870e-04  2.1762e-04    -0.04228
beta_x1_b  9.4989e-04  9.5065e-04  3.3321e-04  3.3311e-04   6.249e-05   6.311e-05  2.1956e-04  2.1847e-04    -0.04244
beta_x2_a  3.3156e-04  3.3321e-04    0.001644    0.001644  4.0614e-04  4.0748e-04  2.4873e-04  2.4635e-04    -0.09291
beta_x2_b  3.3146e-04  3.3311e-04    0.001644    0.001644  4.0606e-04  4.0740e-04  2.4862e-04  2.4624e-04    -0.09289
beta_x3_a   6.187e-05   6.249e-05  4.0614e-04  4.0606e-04  5.3102e-04  5.3153e-04   8.034e-05   7.943e-05    -0.03530
beta_x3_b   6.248e-05   6.311e-05  4.0748e-04  4.0740e-04  5.3153e-04  5.3204e-04   8.104e-05   8.013e-05    -0.03543
beta_x4_a  2.1870e-04  2.1956e-04  2.4873e-04  2.4862e-04   8.034e-05   8.104e-05  3.4556e-04  3.4432e-04    -0.04819
beta_x4_b  2.1762e-04  2.1847e-04  2.4635e-04  2.4624e-04   7.943e-05   8.013e-05  3.4432e-04  3.4309e-04    -0.04796
delta_a      -0.04228    -0.04244   -0.092910   -0.092890    -0.03530    -0.03543    -0.04819    -0.04796     9.31944

Classical correlation matrix:
            beta_x1_a   beta_x1_b   beta_x2_a   beta_x2_b   beta_x3_a   beta_x3_b   beta_x4_a   beta_x4_b     delta_a
beta_x1_a      1.0000      1.6727     0.47537    -0.51814     -0.3260      0.2754    -1.08262     0.94544         NaN
beta_x1_b      1.6727      1.0000    -0.51703     0.47646      0.2757     -0.3256     0.94439    -1.08410         NaN
beta_x2_a      0.4754     -0.5170     1.00000     1.93870     -0.5732      0.6745    -0.01498     0.20084         NaN
beta_x2_b     -0.5181      0.4765     1.93870     1.00000      0.6742     -0.5729     0.20182    -0.01702         NaN
beta_x3_a     -0.3260      0.2757    -0.57321     0.67424      1.0000     -0.5663     0.94047    -0.85339         NaN
beta_x3_b      0.2754     -0.3256     0.67445    -0.57294     -0.5663      1.0000    -0.84967     0.94313         NaN
beta_x4_a     -1.0826      0.9444    -0.01498     0.20182      0.9405     -0.8497     1.00000     1.05571         NaN
beta_x4_b      0.9454     -1.0841     0.20084    -0.01702     -0.8534      0.9431     1.05571     1.00000         NaN
delta_a           NaN         NaN         NaN         NaN         NaN         NaN         NaN         NaN         NaN

Robust correlation matrix:
            beta_x1_a   beta_x1_b   beta_x2_a   beta_x2_b   beta_x3_a   beta_x3_b   beta_x4_a   beta_x4_b     delta_a
beta_x1_a     1.00000     1.00000      0.2654      0.2654     0.08714     0.08792      0.3819      0.3813     -0.4495
beta_x1_b     1.00000     1.00000      0.2665      0.2665     0.08796     0.08874      0.3831      0.3825     -0.4509
beta_x2_a     0.26541     0.26652      1.0000      1.0000     0.43465     0.43567      0.3300      0.3280     -0.7506
beta_x2_b     0.26537     0.26648      1.0000      1.0000     0.43462     0.43564      0.3299      0.3279     -0.7505
beta_x3_a     0.08714     0.08796      0.4346      0.4346     1.00000     1.00000      0.1875      0.1861     -0.5018
beta_x3_b     0.08792     0.08874      0.4357      0.4356     1.00000     1.00000      0.1890      0.1875     -0.5032
beta_x4_a     0.38188     0.38307      0.3300      0.3299     0.18755     0.18899      1.0000      1.0000     -0.8493
beta_x4_b     0.38135     0.38254      0.3280      0.3279     0.18610     0.18755      1.0000      1.0000     -0.8481
delta_a      -0.44950    -0.45091     -0.7506     -0.7505    -0.50177    -0.50320     -0.8493     -0.8481      1.0000

 20 most extreme outliers in terms of lowest average per choice prediction:
  ID Avg prob per choice
 226           0.3592952
 279           0.3598091
 265           0.3731624
 225           0.3922842
 237           0.3944891
 229           0.3951572
 298           0.3959765
 278           0.4027529
 270           0.4027822
 286           0.4030448
 223           0.4155056
 268           0.4155056
 264           0.4157283
 228           0.4161612
 248           0.4166112
 214           0.4207322
 232           0.4207365
 202           0.4239741
 239           0.4241082
 206           0.4242185

Changes in parameter estimates from starting values:
              Initial    Estimate  Difference
beta_x1_a     0.07806     0.07806   2.491e-06
beta_x1_b     0.07806     0.07806   2.491e-06
beta_x2_a    -0.30321    -0.30321  -4.400e-07
beta_x2_b    -0.30321    -0.30321  -4.400e-07
beta_x3_a    -0.34679    -0.34679   1.367e-06
beta_x3_b    -0.34679    -0.34679   1.367e-06
beta_x4_a    -0.21864    -0.21865  -3.771e-06
beta_x4_b    -0.21864    -0.21865  -3.771e-06
delta_a       0.00000     0.00000       0.000
delta_b       0.00000     0.00000       0.000

Settings and functions used in model definition:

apollo_control
--------------
                       Value                  
modelName              "LC_DCE_Model_linear"  
modelDescr             "LC model for DCE data"
indivID                "respid"               
outputDirectory        "output_LC_linear/"    
debug                  "FALSE"                
nCores                 "1"                    
workInLogs             "FALSE"                
seed                   "13"                   
mixing                 "FALSE"                
HB                     "FALSE"                
noValidation           "TRUE"                 
noDiagnostics          "TRUE"                 
calculateLLC           "TRUE"                 
analyticHessian        "FALSE"                
memorySaver            "FALSE"                
panelData              "TRUE"                 
analyticGrad           "FALSE"                
analyticGrad_manualSet "FALSE"                
overridePanel          "FALSE"                
preventOverridePanel   "FALSE"                
noModification         "FALSE"                

Hessian routines attempted
--------------------------
numerical second derivative of LL (using numDeriv)


apollo_lcPars
---------------
function(apollo_beta, apollo_inputs) {
  
  lcpars = list()
  lcpars[["beta_x1"]] = list(beta_x1_a, beta_x1_b)
  lcpars[["beta_x2"]] = list(beta_x2_a, beta_x2_b)
  lcpars[["beta_x3"]] = list(beta_x3_a, beta_x3_b)
  lcpars[["beta_x4"]] = list(beta_x4_a, beta_x4_b)

  lcpars[["pi_values"]] = list(class_a = 1/(1 + exp(delta_b - delta_a)),
                               class_b = 1/(1 + exp(delta_a - delta_b)))

    return(lcpars)
}
<bytecode: 0x0000028f893dcc78>


apollo_probabilities
----------------------
function(apollo_beta, apollo_inputs, functionality = "estimate") {
  
  ### 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 that are generic across classes
  mnl_settings = list(
    alternatives = c(alt1 = 1, alt2 = 2),
    avail = list(alt1 = avail_1, alt2 = avail_2),
    choiceVar = choice_new
  )
  
  ### Loop over classes
  for (s in 1:length(pi_values)) {
    
    ### Define utilities for each alternative (remove "database$" prefix)
    V = list()
    V[["alt1"]] = beta_x1[[s]]*x1_1 + beta_x2[[s]]*x2_1 + beta_x3[[s]]*x3_1 + beta_x4[[s]]*x4_1
    V[["alt2"]] = beta_x1[[s]]*x1_2 + beta_x2[[s]]*x2_2 + beta_x3[[s]]*x3_2 + beta_x4[[s]]*x4_2
    
    mnl_settings$utilities     = V
    mnl_settings$componentName = paste0("Class_",s)
    
    ### Compute within-class choice probabilities using MNL model
    P[[paste0("Class_", s)]] = apollo_mnl(mnl_settings, functionality)
    
    ### Take product across observation for same individual
    P[[paste0("Class_", s)]] = apollo_panelProd(P[[paste0("Class_", s)]], apollo_inputs, functionality)
  }
  
  ### 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)
}

The code seems to calculate an MNL, but not for the individual classes. I have also searched the manual and forum for possible solutions. I may be missing something. I would appreciate any advice.

Thank you very much.
Andrew
stephanehess
Site Admin
Posts: 1264
Joined: 24 Apr 2020, 16:29

Re: Latent Class model with EM algorithm

Post by stephanehess »

Hi there

it looks like your EM estimation converged to a pretty bad solution. This is somewhat surprising, but EM does not in fact guarantee convergence to a better solution than classical estimation. Did you try some other starting values?

Stephane
--------------------------------
Stephane Hess
www.stephanehess.me.uk
Andrew
Posts: 3
Joined: 10 Mar 2025, 14:56

Re: Latent Class model with EM algorithm

Post by Andrew »

Dear Stephane,

thank you very much!

If I change the delta_a in the model parameters of the EM model from 0 to 0.5 or another number, I get a result that is identical to the other model. Since I'm just starting with Apollo, I'm a little unsure which model I should use to estimate. Perhaps the classical estimation is already what I need.

Andrew
dpalma
Posts: 215
Joined: 24 Apr 2020, 17:54

Re: Latent Class model with EM algorithm

Post by dpalma »

Hi Andrew,

Glad to hear that using different starting values led you to the expected solution. In general, when using latent classes we advice against using the same starting values for all classes. In fact, the latest version of Apollo (released recently) will complain if you use the same starting values for two or more classes.

Best wishes,
David
Andrew
Posts: 3
Joined: 10 Mar 2025, 14:56

Re: Latent Class model with EM algorithm

Post by Andrew »

Hi David,

Thank you very much! So far I have used Stata (lclogit, lclogit2, gllamm) for latent class models. I used these to generate hypotheses and hardly ever changed the starting values. Could you share more about starting values or recommend a paper on the topic? Or is there a new Apollo manual available?

Also, is it possible to estimate Scale-adjusted Latent Class Models with Apollo? I know these are often done with Latent Gold, but I believe I read that Apollo might support them too. I couldn't find anything about it in the forum.

Kind regards,
Andrew
dpalma
Posts: 215
Joined: 24 Apr 2020, 17:54

Re: Latent Class model with EM algorithm

Post by dpalma »

Hi Andrew,

Out of the top of my mind, I don't have a reference for the issue with using the same starting values for latent classes (LC). The literature on identification of LC models is spread across multiple research areas and it is tricky business. However, imagine a model with two classes, only constants in the class allocation function, and analogous utility functions for the choice but with different parameters. The class allocation probability at the starting values is 50% for each class. Then the initial gradients for both classes will be the same, and there is no reason to expect that their path towards convergence may differ at any point. If you include sociodemographics in the class allocation function, however, gradients will be different.

We are preparing a new Apollo manual for version 0.3.5. We will release it soon.

About the scale-adjusted LC model, you can estimate one in Apollo. If you want something like what's proposed in Burke et al. (2010), it would look somewhat similar to the DM example in the Apollo website, but combining only two class allocation probabilities (one for the scale class, and another for the preferences class). However, Hess & Train (2017) state that its is not really possible to separate scale heterogeneity from heterogeneity in preferences, so at least from a theoretical perspective it would not be possible to state that your model can differentiate between the two.

Best wishes,
David
Andrew
Posts: 3
Joined: 10 Mar 2025, 14:56

Re: Latent Class model with EM algorithm

Post by Andrew »

David, thank you very much for the detailed reply. And the paper by Burke et al. was indeed new to me - very interesting! When reviewers asked whether the heterogeneity in my LC model was due to preferences or survey behavior, I often responded with a reference to Hess & Train. I know from colleagues that they had difficulties interpreting scale and preference classes resulting from SALC models. Which I can relate, since the results in Burke are not easy for me to fully comprehend. However, I'll gladly take a closer look at the DM model, and I’m looking forward to the new Apollo manual! Great forum, by the way!
Andrew
Last edited by Andrew on 09 Apr 2025, 11:27, edited 2 times in total.
Post Reply