Page 1 of 2

Unconditionals

Posted: 11 Oct 2023, 21:12
by malemu588@gmail.com
Dear Stephane,

I am encountering some problems when I am trying to derive the unconditional for some parameters after successfully estimating my model (a latent class mixed logit model).

apollo_unconditionals(model, apollo_probabilities, apollo_inputs)

Error in length(apollo_draws) && is.na(apollo_draws) :
'length = 8' in coercion to 'logical(1)'

I am not quite why I am getting the above error. Could you please help. Thanks in advance.

Kind regards,
Mohammed

Re: Unconditionals

Posted: 13 Oct 2023, 10:49
by stephanehess
Hi

are you maybe using an old version of Apollo?

Stephane

Re: Unconditionals

Posted: 13 Oct 2023, 10:53
by malemu588@gmail.com
Thanks for your reply Stephane, I am using version 0.3.0.
/Mohammed

Re: Unconditionals

Posted: 13 Oct 2023, 12:06
by stephanehess
thanks for confirming

can you show the code, please?

Re: Unconditionals

Posted: 13 Oct 2023, 15:47
by malemu588@gmail.com
Yes, I am glad to share the code.

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

### Load Apollo library
library(apollo)

### Initialise code
apollo_initialise()

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

# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS ####
# ################################################################# #
database = read.csv("Data_R.csv",header=TRUE)
database <- replace(database,is.na(database),0)
head(database)
dim(database)
# ################################################################# #
#### DEFINE MODEL PARAMETERS ####
# ################################################################# #

### Vector of parameters, including any that are kept fixed in estimation
apollo_beta=c(b_ASCa = 3.5,
b_ASCb = 3.5,
b_ASCc = 3.5,
b_carba = 0.03,
b_carbb = 0.03,
b_carbc = 0.03,
b_justa = 0.3,
b_justb = 0.3,
b_justc = 0.3,
b_costa = -0.1,
b_costb = -0.1,
b_costc = -0.1,
sigma_b_asca = 2.2,
sigma_b_ascb = 2.2,
sigma_b_ascc = 2.2,
sigma_b_carba = 0.06,
sigma_b_carbb = 0.06,
sigma_b_carbc = 0.06,
sigma_b_justa = 0.4,
sigma_b_justb = 0.4,
sigma_b_justc = 0.4,
sigma_b_costa = 0.2,
sigma_b_costb = 0.2,
sigma_b_costc = 0.2,
ageea = 0,
ageeb = 0,
femalea = 0,
femaleb = 0,
#climjusta = 0,
#climjustb = 0,
objectknowa = 0,
objectknowb = 0,
#buyefa = 0,
#buyefb = 0,
seencfa = 0,
seencfb = 0,
deltaa = 0,
deltab = 0,
highedua = 0,
highedub = 0,
carbfpa = 0,
carbfpb = 0,
worka = 0,
workb = 0)
#whitepa = 0,
#whitepb = 0)

## 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 = "MLHS",
interNDraws = 4800,
interUnifDraws = c(),
interNormDraws = c("draws_carba","draws_carbb","draws_carbc","draws_justa","draws_justb","draws_justc","draws_costa","draws_costb","draws_costc","draws_asca","draws_ascb","draws_ascc"),
intraDrawsType = c(),
intraNDraws = 0,
intraUnifDraws = c(),
intraNormDraws = c()
)

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

randcoeff[["b_ascA"]] = (b_ASCa + sigma_b_asca * draws_asca)
randcoeff[["b_ascB"]] = (b_ASCb + sigma_b_ascb * draws_ascb)
randcoeff[["b_ascC"]] = (b_ASCc + sigma_b_ascc * draws_ascc)
randcoeff[["b_carbA"]] = (b_carba + sigma_b_carba * draws_carba)
randcoeff[["b_carbB"]] = (b_carbb + sigma_b_carbb * draws_carbb)
randcoeff[["b_carbC"]] = (b_carbc + sigma_b_carbc * draws_carbc)
randcoeff[["b_justA"]] = (b_justa + sigma_b_justa * draws_justa)
randcoeff[["b_justB"]] = (b_justb + sigma_b_justb * draws_justb)
randcoeff[["b_justC"]] = (b_justc + sigma_b_justc * draws_justc)
randcoeff[["b_costA"]] = -exp(b_costa + sigma_b_costa * draws_costa)
randcoeff[["b_costB"]] = -exp(b_costb + sigma_b_costb * draws_costb)
randcoeff[["b_costC"]] = -exp(b_costc + sigma_b_costc * draws_costc)

return(randcoeff)
}

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

apollo_lcPars = function(apollo_beta, apollo_inputs){
lcpars = list()
lcpars[["b_asc"]] = list(b_ascA, b_ascB, b_ascC)
lcpars[["b_carb"]] = list(b_carbA, b_carbB, b_carbC)
lcpars[["b_just"]] = list(b_justA, b_justB, b_justC)
lcpars[["b_cost"]] = list(b_costA, b_costB, b_costC)

V=list()
V[["class_a"]] = deltaa + ageea * age + femalea * gender1 + highedua * heduc + carbfpa * carbonfp + worka * workerjust + objectknowa * score + seencfa * seen_cf
V[["class_b"]] = deltab + ageeb * age + femaleb * gender1 + highedub * heduc + carbfpb * carbonfp + workb * workerjust + objectknowb * score + seencfb * seen_cf
V[["class_c"]] = 0

classAlloc_settings = list(
classes = c(class_a=1, class_b=2, class_c=3),
utilities = V
)

lcpars[["pi_values"]] = apollo_classAlloc(classAlloc_settings)

return(lcpars)
}

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

apollo_inputs = apollo_validateInputs()
# ################################################################# #
#### DEFINE MODEL AND LIKELIHOOD FUNCTION ####
# ################################################################# #

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
mnl_settings = list(
alternatives = c(alt1=1, alt2=2, alt3=3, alt4=4),
avail = list(alt1=1, alt2=1, alt3=1, alt4=1),
choiceVar = choice
)

### Loop over classes
for(s in 1:3){

### Compute class-specific utilities
V=list()
#V[['alt1']] = b_cost[[s]] * (pric1 + b_carb[[s]] * carbon1 + b_just[[s]] * justice1)
#V[['alt2']] = b_cost[[s]] * (pric2 + b_carb[[s]] * carbon2 + b_just[[s]] * justice2)
#V[['alt3']] = b_cost[[s]] * (pric3 + b_carb[[s]] * carbon3 + b_just[[s]] * justice3)
#V[['alt4']] = b_cost[[s]] * (b_asc[[s]])
V[['alt1']] = b_cost[[s]] * (pric1/10 + b_carb[[s]] * carbon1 + b_just[[s]] * justice1)
V[['alt2']] = b_cost[[s]] * (pric2/10 + b_carb[[s]] * carbon2 + b_just[[s]] * justice2)
V[['alt3']] = b_cost[[s]] * (pric3/10 + b_carb[[s]] * carbon3 + b_just[[s]] * justice3)
V[['alt4']] = b_cost[[s]] * b_asc[[s]]

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)

### Average across inter-individual draws within classes
P[[paste0("Class_",s)]] = apollo_avgInterDraws(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)

### Average across inter-individual draws in class allocation probabilities
P[["model"]] = apollo_avgInterDraws(P[["model"]], apollo_inputs, functionality)

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

### Estimate model
model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs, estimate_settings=list(maxIterations=1000))

### Show output in screen
apollo_modelOutput(model)

### Save output to file(s)
apollo_saveOutput(model)

?apollo_probabilities
apollo_probabilities(apollo_beta, apollo_inputs)

# ################################################################# #
##### POST-PROCESSING ####
# ################################################################# #

### Optional: load previously estimated model object
#model = apollo_loadModel(apollo_control$modelName)

### Print outputs of additional diagnostics to new output file (remember to close file writing when complete)
apollo_sink()

# ----------------------------------------------------------------- #
#---- CONDITIONALS AND UNCONDITIONALS ----
# ----------------------------------------------------------------- #
apollo_unconditionals(model, apollo_probabilities, apollo_inputs)

Re: Unconditionals

Posted: 17 Oct 2023, 21:30
by anders
Hi,

I think this has to do with estimating a latent class mixed logit model for attribute non-attendance. I am struggeling with the same problem of extracting unconditionals after successfully estimating a latent class mixed logit model with complete ANA structure. I suspect it is because one class is specified as full non-attendance, and hence, unconditional estimates are set to zero.

Below is my code.

Code: Select all

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

### Load Apollo library
#library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "LC-MMNL-ANA-8",
  modelDescr      = "",
  indivID         = "id",
  mixing          = TRUE, 
  nCores          =    10,
  outputDirectory = "output"
)

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

### Loading data from
database = read.csv("apollo_data_wind_main_full.csv",header=TRUE)

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

### Vector of parameters
apollo_beta = c(mu_twh2                             =              2.9828,
                mu_twh3                                     =              5.1948,
                mu_twh4                                     =              5.9565,
                mu_turb7                                    =             -4.7124,
                mu_turb14                                   =             -6.9625,
                mu_turb21                                   =             -7.8751,
                mu_cost                                     =                  -2,
                sigma_twh2                                  =                   0,
                sigma_twh3                                  =                   0,
                sigma_twh4                                  =                   0,
                sigma_turb7                                 =                   0,
                sigma_turb14                                =                   0,
                sigma_turb21                                =                   0,
                sigma_cost                                  =                   0,
                delta_a                                     =              0.6696,
                delta_b                                     =             -1.0402,
                delta_c                                     =              0.3473,
                delta_d                                     =             -2.0474,
                delta_e                                     =              0.6300,
                delta_f                                     =              1.9412,
                delta_g                                     =             -0.5162)

### 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()


### Set parameters for generating draws
apollo_draws = list(
  interDrawsType="mlhs",           
  interNDraws=500,                   
  interUnifDraws=c(),                
  interNormDraws=c("xi_twh2","xi_twh3","xi_twh4","xi_turb7","xi_turb14","xi_turb21","xi_cost"),    
  
  intraDrawsType="mlhs",
  intraNDraws=0,
  intraUnifDraws=c(),
  intraNormDraws=c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["b_twh2"]] =mu_twh2 + sigma_twh2 * xi_twh2
  randcoeff[["b_twh3"]] = mu_twh3 + sigma_twh3 * xi_twh3
  randcoeff[["b_twh4"]] = mu_twh4 + sigma_twh4 * xi_twh4
  randcoeff[["b_turb7"]] = mu_turb7 + sigma_turb7 * xi_turb7
  randcoeff[["b_turb14"]] = mu_turb14 + sigma_turb14 * xi_turb14 
  randcoeff[["b_turb21"]] = mu_turb21 + sigma_turb21 * xi_turb21
  randcoeff[["b_cost"]] =exp(mu_cost + sigma_cost * xi_cost)
  
  
  return(randcoeff)
}

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

apollo_lcPars=function(apollo_beta, apollo_inputs){
  lcpars = list()
  
  V=list()
  V[["class_a"]] = delta_a
  V[["class_b"]] = delta_b
  V[["class_c"]] = delta_c
  V[["class_d"]] = delta_d
  V[["class_e"]] = delta_e
  V[["class_f"]] = delta_f
  V[["class_g"]] = delta_g
  V[["class_h"]] = 0
  
  ### Settings for class allocation models
  classAlloc_settings = list(
    classes      = c(class_a=1, class_b=2, class_c=3, class_d=4, class_e=5, class_f=6, class_g=7, class_h=8), 
    utilities    = V  
  )
  
  lcpars[["pi_values"]] = apollo_classAlloc(classAlloc_settings)
  
  return(lcpars)
}

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

apollo_inputs = apollo_validateInputs()

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

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=1, alt2=1),
    choiceVar    = choice
  )
  
  ### Loop over classes
  for(s in 1:8){
    
    V = list()
    
    if(s==1){
      V[['alt1']]  = ((b_cost*(feeneg_1/100) + b_twh2*(twh_1==20) + b_twh3*(twh_1==30) + b_twh4*(twh_1==40)))
      V[['alt2']]  = ((b_cost*(feeneg_2/100) + b_twh2*(twh_2==20) + b_twh3*(twh_2==30) + b_twh4*(twh_2==40) + b_turb7*(turbines_2==700) + b_turb14*(turbines_2==1400) + b_turb21*(turbines_2==2100)))
    }
    
    if(s==2){
      V[['alt1']]  = ((0*(feeneg_1/100) + b_twh2*(twh_1==20) + b_twh3*(twh_1==30) + b_twh4*(twh_1==40)))
      V[['alt2']]  = ((0*(feeneg_2/100) + b_twh2*(twh_2==20) + b_twh3*(twh_2==30) + b_twh4*(twh_2==40) + b_turb7*(turbines_2==700) + b_turb14*(turbines_2==1400) + b_turb21*(turbines_2==2100)))
    }
    
    if(s==3){
      V[['alt1']]  = ((b_cost*(feeneg_1/100) + b_twh2*(twh_1==20) + b_twh3*(twh_1==30) + b_twh4*(twh_1==40)))
      V[['alt2']]  = ((b_cost*(feeneg_2/100) + b_twh2*(twh_2==20) + b_twh3*(twh_2==30) + b_twh4*(twh_2==40) + 0*(turbines_2==700) + 0*(turbines_2==1400) + 0*(turbines_2==2100)))
    }
    
    if(s==4){
      V[['alt1']]  = ((b_cost*(feeneg_1/100) + 0*(twh_1==20) + 0*(twh_1==30) + 0*(twh_1==40)))
      V[['alt2']]  = ((b_cost*(feeneg_2/100) + 0*(twh_2==20) + 0*(twh_2==30) + 0*(twh_2==40) + b_turb7*(turbines_2==700) + b_turb14*(turbines_2==1400) + b_turb21*(turbines_2==2100)))
    }
    
    if(s==5){
      V[['alt1']]  = ((b_cost*(feeneg_1/100) + b_twh2*0*(twh_1==20) + b_twh3*0*(twh_1==30) + b_twh4*0*(twh_1==40)))
      V[['alt2']]  = ((b_cost*(feeneg_2/100) + b_twh2*0*(twh_2==20) + b_twh3*0*(twh_2==30) + b_twh4*0*(twh_2==40) + b_turb7*0*(turbines_2==700) + b_turb14*0*(turbines_2==1400) + b_turb21*0*(turbines_2==2100)))
    }
    
    if(s==6){
      V[['alt1']]  = ((b_cost*0*(feeneg_1/100) + b_twh2*0*(twh_1==20) + b_twh3*0*(twh_1==30) + b_twh4*0*(twh_1==40)))
      V[['alt2']]  = ((b_cost*0*(feeneg_2/100) + b_twh2*0*(twh_2==20) + b_twh3*0*(twh_2==30) + b_twh4*0*(twh_2==40) + b_turb7*(turbines_2==700) + b_turb14*(turbines_2==1400) + b_turb21*(turbines_2==2100)))
    }
    
    if(s==7){
      V[['alt1']]  = ((b_cost*0*(feeneg_1/100) + b_twh2*(twh_1==20) + b_twh3*(twh_1==30) + b_twh4*(twh_1==40)))
      V[['alt2']]  = ((b_cost*0*(feeneg_2/100) + b_twh2*(twh_2==20) + b_twh3*(twh_2==30) + b_twh4*(twh_2==40) + b_turb7*0*(turbines_2==700) + b_turb14*0*(turbines_2==1400) + b_turb21*0*(turbines_2==2100)))
    }
    
    if(s==8){
      V[['alt1']]  = ((0*b_cost*(feeneg_1/100) + 0*b_twh2*(twh_1==20) + 0*b_twh3*(twh_1==30) + 0*b_twh4*(twh_1==40)))
      V[['alt2']]  = ((0*b_cost*(feeneg_2/100) + 0*b_twh2*(twh_2==20) + 0*b_twh3*(twh_2==30) + 0*b_twh4*(twh_2==40) + 0*b_turb7*(turbines_2==700) + 0*b_turb14*(turbines_2==1400) + 0*b_turb21*(turbines_2==2100)))
    }
    
    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)
    
    ### Average across inter-individual draws within classes
    P[[paste0("Class_",s)]] = apollo_avgInterDraws(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)
}

# ################################################################# #
#### EM ESTIMATION FOR COVARIANCE MATRIX                         ####
# ################################################################# #

model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs)

# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO SCREEN)                               ----
# ----------------------------------------------------------------- #

apollo_saveOutput(model)

unconditionals = apollo_unconditionals(model,apollo_probabilities,apollo_inputs)

Re: Unconditionals

Posted: 19 Oct 2023, 10:09
by stephanehess
Hi

in relation to the question by Mohammed, the error

Code: Select all

apollo_unconditionals(model, apollo_probabilities, apollo_inputs)

Error in length(apollo_draws) && is.na(apollo_draws) :
'length = 8' in coercion to 'logical(1)'
definitely relates to an older version of Apollo as this code no longer exists in version 0.3.0. So please double check your version, and you may also in fact now update to 0.3.1

Stephane

Re: Unconditionals

Posted: 19 Oct 2023, 10:11
by stephanehess
In relation to the post by Anders, this seems a different issue. What error are you getting. Note that you are not defining betas inside apollo_lcPars, and that will affect the functionality of apollo_conditionals and apollo_unconditionals
anders wrote: 17 Oct 2023, 21:30 Hi,

I think this has to do with estimating a latent class mixed logit model for attribute non-attendance. I am struggeling with the same problem of extracting unconditionals after successfully estimating a latent class mixed logit model with complete ANA structure. I suspect it is because one class is specified as full non-attendance, and hence, unconditional estimates are set to zero.


Re: Unconditionals

Posted: 19 Oct 2023, 13:54
by malemu588@gmail.com
Hi Stephane,

Thanks for your time and support in this forum. I see the version I am using is 0.3.0, so not quite sure when you say the error is related to an old version of Apollo. Of course, I will try to update to the newest version (if available) and see if the error can be eliminated.

/Mohammed

Re: Unconditionals

Posted: 19 Oct 2023, 20:37
by anders
Thank you very much, Stephane.

I am now defining betas inside apollo_lcPars, it provides me with the same results as my previous code.

I run the following code:

Code: Select all

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

### Load Apollo library
#library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "LC-MMNL-ANA-8",
  modelDescr      = "",
  indivID         = "id",
  mixing          = TRUE, 
  nCores          =    10,
  outputDirectory = "output"
)

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

### Loading data from
database = read.csv("apollo_data_wind_main_full.csv",header=TRUE)

#Use subset
database = subset(database,database$utvalg==2)


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

### Vector of parameters
apollo_beta = c(mu_twh2_a                                   =              2.9828,
                mu_twh3_a                                   =              5.1948,
                mu_twh4_a                                   =              5.9565,
                b_twh2_b                                    =                   0,
                b_twh3_b                                    =                   0,
                b_twh4_b                                    =                   0,
                mu_turb7_a                                  =             -4.7124,
                mu_turb14_a                                 =             -6.9625,
                mu_turb21_a                                 =             -7.8751,
                b_turb7_b                                   =                   0,
                b_turb14_b                                  =                   0,
                b_turb21_b                                  =                   0,
                b_cost_a                                    =              0.1168,
                b_cost_b                                    =                   0,
                sigma_twh2_a                                =                   0,
                sigma_twh3_a                                =                   0,
                sigma_twh4_a                                =                   0,
                sigma_turb7_a                               =                   0,
                sigma_turb14_a                              =                   0,
                sigma_turb21_a                              =                   0,
                delta_a                                     =              0.6696,
                delta_b                                     =             -1.0402,
                delta_c                                     =              0.3473,
                delta_d                                     =             -2.0474,
                delta_e                                     =              0.6300,
                delta_f                                     =              1.9412,
                delta_g                                     =             -0.5162)

### 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("b_twh2_b","b_twh3_b","b_twh4_b","b_turb7_b","b_turb14_b","b_turb21_b","b_cost_b")


### Set parameters for generating draws
apollo_draws = list(
  interDrawsType="mlhs",           
  interNDraws=500,                   
  interNormDraws=c("xi_twh2","xi_twh3","xi_twh4","xi_turb7","xi_turb14","xi_turb21")
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["b_twh2_a"]] =mu_twh2_a + sigma_twh2_a * xi_twh2
  randcoeff[["b_twh3_a"]] = mu_twh3_a + sigma_twh3_a * xi_twh3
  randcoeff[["b_twh4_a"]] = mu_twh4_a + sigma_twh4_a * xi_twh4
  randcoeff[["b_turb7_a"]] = mu_turb7_a + sigma_turb7_a * xi_turb7
  randcoeff[["b_turb14_a"]] = mu_turb14_a + sigma_turb14_a * xi_turb14 
  randcoeff[["b_turb21_a"]] = mu_turb21_a + sigma_turb21_a * xi_turb21

  
  return(randcoeff)
}

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

apollo_lcPars=function(apollo_beta, apollo_inputs){
  lcpars = list()
  lcpars[["b_twh2"  ]] = list(b_twh2_a, b_twh2_b)
  lcpars[["b_twh3"  ]] = list(b_twh3_a, b_twh3_b)
  lcpars[["b_twh4"  ]] = list(b_twh4_a, b_twh4_b)
  lcpars[["b_turb7"  ]] = list(b_turb7_a, b_turb7_b)
  lcpars[["b_turb14"  ]] = list(b_turb14_a, b_turb14_b)
  lcpars[["b_turb21"  ]] = list(b_turb21_a, b_turb21_b)
  lcpars[["b_cost"  ]] = list(b_cost_a, b_cost_b)
  
  V=list()
  V[["class_a"]] = delta_a
  V[["class_b"]] = delta_b
  V[["class_c"]] = delta_c
  V[["class_d"]] = delta_d
  V[["class_e"]] = delta_e
  V[["class_f"]] = delta_f
  V[["class_g"]] = delta_g
  V[["class_h"]] = 0
  
  ### Settings for class allocation models
  classAlloc_settings = list(
    classes      = c(class_a=1, class_b=2, class_c=3, class_d=4, class_e=5, class_f=6, class_g=7, class_h=8), 
    utilities    = V  
  )
  
  lcpars[["pi_values"]] = apollo_classAlloc(classAlloc_settings)
  
  return(lcpars)
}

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

apollo_inputs = apollo_validateInputs()

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

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=1, alt2=1),
    choiceVar    = choice
  )
  
  ### Loop over classes
  for(s in 1:8){
    
    V = list()
    
    if(s==1){
      V[['alt1']]  = ((b_cost_a*(feeneg_1/100) + b_twh2_a*(twh_1==20) + b_twh3_a*(twh_1==30) + b_twh4_a*(twh_1==40)))
      V[['alt2']]  = ((b_cost_a*(feeneg_2/100) + b_twh2_a*(twh_2==20) + b_twh3_a*(twh_2==30) + b_twh4_a*(twh_2==40) + b_turb7_a*(turbines_2==700) + b_turb14_a*(turbines_2==1400) + b_turb21_a*(turbines_2==2100)))
      
      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)
      
      ### Average across inter-individual draws within classes
      P[[paste0("Class_",s)]] = apollo_avgInterDraws(P[[paste0("Class_",s)]], apollo_inputs, functionality)
      }
    
    if(s==2){
      V[['alt1']]  = ((b_cost_b*(feeneg_1/100) + b_twh2_a*(twh_1==20) + b_twh3_a*(twh_1==30) + b_twh4_a*(twh_1==40)))
      V[['alt2']]  = ((b_cost_b*(feeneg_2/100) + b_twh2_a*(twh_2==20) + b_twh3_a*(twh_2==30) + b_twh4_a*(twh_2==40) + b_turb7_a*(turbines_2==700) + b_turb14_a*(turbines_2==1400) + b_turb21_a*(turbines_2==2100)))
      
      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)
      
      ### Average across inter-individual draws within classes
      P[[paste0("Class_",s)]] = apollo_avgInterDraws(P[[paste0("Class_",s)]], apollo_inputs, functionality)
      }
    
    if(s==3){
      V[['alt1']]  = ((b_cost_a*(feeneg_1/100) + b_twh2_a*(twh_1==20) + b_twh3_a*(twh_1==30) + b_twh4_a*(twh_1==40)))
      V[['alt2']]  = ((b_cost_a*(feeneg_2/100) + b_twh2_a*(twh_2==20) + b_twh3_a*(twh_2==30) + b_twh4_a*(twh_2==40) + b_turb7_b*(turbines_2==700) + b_turb14_b*(turbines_2==1400) + b_turb21_b*(turbines_2==2100)))
      
      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)
      
      ### Average across inter-individual draws within classes
      P[[paste0("Class_",s)]] = apollo_avgInterDraws(P[[paste0("Class_",s)]], apollo_inputs, functionality)
      }
    
    if(s==4){
      V[['alt1']]  = ((b_cost_a*(feeneg_1/100) + b_twh2_b*(twh_1==20) + b_twh3_b*(twh_1==30) + b_twh4_b*(twh_1==40)))
      V[['alt2']]  = ((b_cost_a*(feeneg_2/100) + b_twh2_b*(twh_2==20) + b_twh3_b*(twh_2==30) + b_twh4_b*(twh_2==40) + b_turb7_a*(turbines_2==700) + b_turb14_a*(turbines_2==1400) + b_turb21_a*(turbines_2==2100)))
    
      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)
      
      ### Average across inter-individual draws within classes
      P[[paste0("Class_",s)]] = apollo_avgInterDraws(P[[paste0("Class_",s)]], apollo_inputs, functionality)
      }
    
    if(s==5){
      V[['alt1']]  = ((b_cost_a*(feeneg_1/100) + b_twh2_b*(twh_1==20) + b_twh3_b*(twh_1==30) + b_twh4_b*(twh_1==40)))
      V[['alt2']]  = ((b_cost_a*(feeneg_2/100) + b_twh2_b*(twh_2==20) + b_twh3_b*(twh_2==30) + b_twh4_b*(twh_2==40) + b_turb7_b*(turbines_2==700) + b_turb14_b*(turbines_2==1400) + b_turb21_b*(turbines_2==2100)))
      
      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)
      
      ### Average across inter-individual draws within classes
      #P[[paste0("Class_",s)]] = apollo_avgInterDraws(P[[paste0("Class_",s)]], apollo_inputs, functionality)
      }
    
    if(s==6){
      V[['alt1']]  = ((b_cost_b*(feeneg_1/100) + b_twh2_b*(twh_1==20) + b_twh3_b*(twh_1==30) + b_twh4_b*(twh_1==40)))
      V[['alt2']]  = ((b_cost_b*(feeneg_2/100) + b_twh2_b*(twh_2==20) + b_twh3_b*(twh_2==30) + b_twh4_b*(twh_2==40) + b_turb7_a*(turbines_2==700) + b_turb14_a*(turbines_2==1400) + b_turb21_a*(turbines_2==2100)))
      
      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)
      
      ### Average across inter-individual draws within classes
      P[[paste0("Class_",s)]] = apollo_avgInterDraws(P[[paste0("Class_",s)]], apollo_inputs, functionality)
      }
    
    if(s==7){
      V[['alt1']]  = ((b_cost_b*(feeneg_1/100) + b_twh2_a*(twh_1==20) + b_twh3_a*(twh_1==30) + b_twh4_a*(twh_1==40)))
      V[['alt2']]  = ((b_cost_b*(feeneg_2/100) + b_twh2_a*(twh_2==20) + b_twh3_a*(twh_2==30) + b_twh4_a*(twh_2==40) + b_turb7_b*(turbines_2==700) + b_turb14_b*(turbines_2==1400) + b_turb21_b*(turbines_2==2100)))
      
      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)
      
      ### Average across inter-individual draws within classes
      P[[paste0("Class_",s)]] = apollo_avgInterDraws(P[[paste0("Class_",s)]], apollo_inputs, functionality)
      }
    
    if(s==8){
      V[['alt1']]  = ((b_cost_b*(feeneg_1/100) + b_twh2_b*(twh_1==20) + b_twh3_b*(twh_1==30) + b_twh4_b*(twh_1==40)))
      V[['alt2']]  = ((b_cost_b*(feeneg_2/100) + b_twh2_b*(twh_2==20) + b_twh3_b*(twh_2==30) + b_twh4_b*(twh_2==40) + b_turb7_b*(turbines_2==700) + b_turb14_b*(turbines_2==1400) + b_turb21_b*(turbines_2==2100)))
     
       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)
      
      ### Average across inter-individual draws within classes
      #P[[paste0("Class_",s)]] = apollo_avgInterDraws(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)
}

# ################################################################# #
#### EM ESTIMATION FOR COVARIANCE MATRIX                         ####
# ################################################################# #

model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs)

# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO SCREEN)                               ----
# ----------------------------------------------------------------- #

apollo_saveOutput(model)

unconditionals = apollo_unconditionals(model,apollo_probabilities,apollo_inputs)
It runs smoothly and I get consistent results, but when I run the unconditionals part, I get the following error message:
Error in !is.na(apollo_draws) && apollo_draws$intraNDraws == 0 :
'length = 5' in coercion to 'logical(1)'
I have also tried to run an alternative specification of the ANA model that provides the same result, but also the same error message when I run the apollo_unconditionals command. The alternative specification is:

Code: Select all

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

### Load Apollo library
#library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "LC-MMNL-ANA-8",
  modelDescr      = "",
  indivID         = "id",
  mixing          = TRUE, 
  nCores          =    10,
  outputDirectory = "output"
)

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

### Loading data from
database = read.csv("apollo_data_wind_main_full.csv",header=TRUE)

#Use subset
database = subset(database,database$utvalg==2)


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

### Vector of parameters
apollo_beta = c(mu_twh2_a                                   =              2.9828,
                mu_twh3_a                                   =              5.1948,
                mu_twh4_a                                   =              5.9565,
                mu_twh2_b                                    =                   0,
                mu_twh3_b                                    =                   0,
                mu_twh4_b                                    =                   0,
                mu_turb7_a                                  =             -4.7124,
                mu_turb14_a                                 =             -6.9625,
                mu_turb21_a                                 =             -7.8751,
                mu_turb7_b                                  =                   0,
                mu_turb14_b                                 =                   0,
                mu_turb21_b                                 =                   0,
                mu_cost_a                                   =              0.1168,
                mu_cost_b                                   =                   0,
                sigma_twh2_a                                =                   0,
                sigma_twh3_a                                =                   0,
                sigma_twh4_a                                =                   0,
                sigma_turb7_a                               =                   0,
                sigma_turb14_a                              =                   0,
                sigma_turb21_a                              =                   0,
                delta_a                                     =              0.6696,
                delta_b                                     =             -1.0402,
                delta_c                                     =              0.3473,
                delta_d                                     =             -2.0474,
                delta_e                                     =              0.6300,
                delta_f                                     =              1.9412,
                delta_g                                     =             -0.5162)

### 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("mu_twh2_b","mu_twh3_b","mu_twh4_b","mu_turb7_b","mu_turb14_b","mu_turb21_b","mu_cost_b")


### Set parameters for generating draws
apollo_draws = list(
  interDrawsType="mlhs",           
  interNDraws=500,                   
  interNormDraws=c("xi_twh2","xi_twh3","xi_twh4","xi_turb7","xi_turb14","xi_turb21","xi_cost","xi_twh2_b","xi_twh3_b","xi_twh4_b","xi_turb7_b","xi_turb14_b","xi_turb21_b","xi_cost_b")    
  
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["b_twh2_a"]] =mu_twh2_a + sigma_twh2_a * xi_twh2
  randcoeff[["b_twh3_a"]] = mu_twh3_a + sigma_twh3_a * xi_twh3
  randcoeff[["b_twh4_a"]] = mu_twh4_a + sigma_twh4_a * xi_twh4
  randcoeff[["b_turb7_a"]] = mu_turb7_a + sigma_turb7_a * xi_turb7
  randcoeff[["b_turb14_a"]] = mu_turb14_a + sigma_turb14_a * xi_turb14 
  randcoeff[["b_turb21_a"]] = mu_turb21_a + sigma_turb21_a * xi_turb21
  randcoeff[["b_cost_a"]] =mu_cost_a + 0 * xi_cost
  
  randcoeff[["b_twh2_b"]] =mu_twh2_b + 0 * xi_twh2_b
  randcoeff[["b_twh3_b"]] = mu_twh3_b + 0 * xi_twh3_b
  randcoeff[["b_twh4_b"]] = mu_twh4_b + 0 * xi_twh4_b
  randcoeff[["b_turb7_b"]] = mu_turb7_b + 0 * xi_turb7_b
  randcoeff[["b_turb14_b"]] = mu_turb14_b + 0 * xi_turb14_b 
  randcoeff[["b_turb21_b"]] = mu_turb21_b + 0 * xi_turb21_b
  randcoeff[["b_cost_b"]] =mu_cost_b + 0 * xi_cost_b
  
  
  return(randcoeff)
}

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

apollo_lcPars=function(apollo_beta, apollo_inputs){
  lcpars = list()
  lcpars[["b_twh2"  ]] = list(b_twh2_a, b_twh2_b)
  lcpars[["b_twh3"  ]] = list(b_twh3_a, b_twh3_b)
  lcpars[["b_twh4"  ]] = list(b_twh4_a, b_twh4_b)
  lcpars[["b_turb7"  ]] = list(b_turb7_a, b_turb7_b)
  lcpars[["b_turb14"  ]] = list(b_turb14_a, b_turb14_b)
  lcpars[["b_turb21"  ]] = list(b_turb21_a, b_turb21_b)
  lcpars[["b_cost"  ]] = list(b_cost_a, b_cost_b)
  
  V=list()
  V[["class_a"]] = delta_a
  V[["class_b"]] = delta_b
  V[["class_c"]] = delta_c
  V[["class_d"]] = delta_d
  V[["class_e"]] = delta_e
  V[["class_f"]] = delta_f
  V[["class_g"]] = delta_g
  V[["class_h"]] = 0
  
  ### Settings for class allocation models
  classAlloc_settings = list(
    classes      = c(class_a=1, class_b=2, class_c=3, class_d=4, class_e=5, class_f=6, class_g=7, class_h=8), 
    utilities    = V  
  )
  
  lcpars[["pi_values"]] = apollo_classAlloc(classAlloc_settings)
  
  return(lcpars)
}

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

apollo_inputs = apollo_validateInputs()

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

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=1, alt2=1),
    choiceVar    = choice
  )
  
  ### Loop over classes
  for(s in 1:8){
    
    V = list()
    
    if(s==1){
      V[['alt1']]  = ((b_cost_a*(feeneg_1/100) + b_twh2_a*(twh_1==20) + b_twh3_a*(twh_1==30) + b_twh4_a*(twh_1==40)))
      V[['alt2']]  = ((b_cost_a*(feeneg_2/100) + b_twh2_a*(twh_2==20) + b_twh3_a*(twh_2==30) + b_twh4_a*(twh_2==40) + b_turb7_a*(turbines_2==700) + b_turb14_a*(turbines_2==1400) + b_turb21_a*(turbines_2==2100)))
    
    }
    
    if(s==2){
      V[['alt1']]  = ((b_cost_b*(feeneg_1/100) + b_twh2_a*(twh_1==20) + b_twh3_a*(twh_1==30) + b_twh4_a*(twh_1==40)))
      V[['alt2']]  = ((b_cost_b*(feeneg_2/100) + b_twh2_a*(twh_2==20) + b_twh3_a*(twh_2==30) + b_twh4_a*(twh_2==40) + b_turb7_a*(turbines_2==700) + b_turb14_a*(turbines_2==1400) + b_turb21_a*(turbines_2==2100)))
      
    }
    
    if(s==3){
      V[['alt1']]  = ((b_cost_a*(feeneg_1/100) + b_twh2_a*(twh_1==20) + b_twh3_a*(twh_1==30) + b_twh4_a*(twh_1==40)))
      V[['alt2']]  = ((b_cost_a*(feeneg_2/100) + b_twh2_a*(twh_2==20) + b_twh3_a*(twh_2==30) + b_twh4_a*(twh_2==40) + b_turb7_b*(turbines_2==700) + b_turb14_b*(turbines_2==1400) + b_turb21_b*(turbines_2==2100)))
      
    }
    
    if(s==4){
      V[['alt1']]  = ((b_cost_a*(feeneg_1/100) + b_twh2_b*(twh_1==20) + b_twh3_b*(twh_1==30) + b_twh4_b*(twh_1==40)))
      V[['alt2']]  = ((b_cost_a*(feeneg_2/100) + b_twh2_b*(twh_2==20) + b_twh3_b*(twh_2==30) + b_twh4_b*(twh_2==40) + b_turb7_a*(turbines_2==700) + b_turb14_a*(turbines_2==1400) + b_turb21_a*(turbines_2==2100)))
      
    }
    
    if(s==5){
      V[['alt1']]  = ((b_cost_a*(feeneg_1/100) + b_twh2_b*(twh_1==20) + b_twh3_b*(twh_1==30) + b_twh4_b*(twh_1==40)))
      V[['alt2']]  = ((b_cost_a*(feeneg_2/100) + b_twh2_b*(twh_2==20) + b_twh3_b*(twh_2==30) + b_twh4_b*(twh_2==40) + b_turb7_b*(turbines_2==700) + b_turb14_b*(turbines_2==1400) + b_turb21_b*(turbines_2==2100)))
      
    }
    
    if(s==6){
      V[['alt1']]  = ((b_cost_b*(feeneg_1/100) + b_twh2_b*(twh_1==20) + b_twh3_b*(twh_1==30) + b_twh4_b*(twh_1==40)))
      V[['alt2']]  = ((b_cost_b*(feeneg_2/100) + b_twh2_b*(twh_2==20) + b_twh3_b*(twh_2==30) + b_twh4_b*(twh_2==40) + b_turb7_a*(turbines_2==700) + b_turb14_a*(turbines_2==1400) + b_turb21_a*(turbines_2==2100)))
      
    }
    
    if(s==7){
      V[['alt1']]  = ((b_cost_b*(feeneg_1/100) + b_twh2_a*(twh_1==20) + b_twh3_a*(twh_1==30) + b_twh4_a*(twh_1==40)))
      V[['alt2']]  = ((b_cost_b*(feeneg_2/100) + b_twh2_a*(twh_2==20) + b_twh3_a*(twh_2==30) + b_twh4_a*(twh_2==40) + b_turb7_b*(turbines_2==700) + b_turb14_b*(turbines_2==1400) + b_turb21_b*(turbines_2==2100)))
      
    }
    
    if(s==8){
      V[['alt1']]  = ((b_cost_b*(feeneg_1/100) + b_twh2_b*(twh_1==20) + b_twh3_b*(twh_1==30) + b_twh4_b*(twh_1==40)))
      V[['alt2']]  = ((b_cost_b*(feeneg_2/100) + b_twh2_b*(twh_2==20) + b_twh3_b*(twh_2==30) + b_twh4_b*(twh_2==40) + b_turb7_b*(turbines_2==700) + b_turb14_b*(turbines_2==1400) + b_turb21_b*(turbines_2==2100)))
      
    }
    
    
    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)
    
    ### Average across inter-individual draws within classes
    P[[paste0("Class_",s)]] = apollo_avgInterDraws(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)
}

# ################################################################# #
#### EM ESTIMATION FOR COVARIANCE MATRIX                         ####
# ################################################################# #

model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs)

# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO SCREEN)                               ----
# ----------------------------------------------------------------- #

apollo_saveOutput(model)

unconditionals = apollo_unconditionals(model,apollo_probabilities,apollo_inputs)