Page 1 of 1

Correlations across alternatives

Posted: 21 Jan 2025, 08:57
by ifurqanbhat
Hello Prof. Stephane,

I am trying to estimate a joint model with the data from two stated preference choice surveys. In my data collection, the respondents were first asked choose their preferred vehicle (between electric and conventional). Further in the same questionnaire, the same respondents were asked to choose their preferred location of charging, if they were to purchase an electric vehicle in future. I am trying to estimate a joint model (I have some variables that are same across both the games) to understand if the sensitivities to the common variables remains the same for the two types of decisions. I also wanted to introduce correlation terms to see if there is a correlation between choosing a type of vehicle and a type of location. For example, I wanted to see if there is correlation between choosing an electric vehicle and a particular type of charging location and similarly with the choice of conventional vehicle and choice of location.

What is the correct way to introduce these correlation terms?

Code: Select all

### Clear memory
rm(list = ls())
# install.packages("apollo")
library(apollo)
database = read.csv(file.choose(),header=TRUE)
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName  ="V1_Joint",
  modelDescr ="Joint SP1/SP2 Model Estimation",
  indivID    ="Person_ID",
  nCores = 32
)

apollo_beta=c(asc_CV                  =  0,
              asc_EV                  =  0,
              asc_L                   =  0,
              asc_W                   =  0,
              asc_H                   =  0,
              
              b_pp                    =  0,
              b_oc                    =  0,
              b_ct                    =  0,
              b_rn                    =  0,
              b_cf                    =  0,
              b_em                    =  0,
              
              b_cc                    =  0,
              b_wt                    =  0,
              
              mu_SP1                       =  1,
              mu_SP2                       =  1,
              sigma                        =  0
)

apollo_fixed = c("asc_CV", "asc_H",
                 "mu_SP1")

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

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "mlhs",
  interNDraws    = 1000,
  interUnifDraws = c(),
  interNormDraws = c("draws_alt1", "draws_alt2", "draws_alt3", "draws_alt4", "draws_alt5"),
  
  intraDrawsType = "halton",
  intraNDraws    = 0,
  intraUnifDraws = c(),
  intraNormDraws = c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["ec1"]]  = sigma*draws_alt1
  randcoeff[["ec2"]]  = sigma*draws_alt2
  randcoeff[["ec3"]]  = sigma*draws_alt3
  randcoeff[["ec4"]]  = sigma*draws_alt4
  randcoeff[["ec5"]]  = sigma*draws_alt5
  
  return(randcoeff)
}

# ################################################################# #
#### 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()
  
  ### List of utilities (before applying scales): these must use the same names as in mnl_settings, order is irrelevant
  V = list()
  
  #######################################################################################
  
  V[['EV']]         = asc_EV  +  b_pp*ppEV* +  b_oc*ocEV  +  b_ct*ctEV + b_rn*rEV  + b_cf*(cfEV) + b_em*(eEV) + ec1
  V[['CV']]         = asc_CV  +  b_pp*ppCV* +  b_oc*ocCV               + b_rn*rCV                             + ec2

  V[['leisure']]    = asc_L  +  b_ct*ctLeisure     +  b_cc*ccLeisure    +  b_wt*wtLeisure   + b_em*eLeisure   + ec3
  V[['workplace']]  = asc_W  +  b_ct*ctWorkplace   +  b_cc*ccWorkplace  +  b_wt*wtWorkplace + b_em*eWorkplace + ec4
  V[['highway']]    = asc_H  +  b_ct*ctHighway     +  b_cc*ccHighway    +  b_wt*wtHighway   + b_em*eHighway   + ec5
  
  ### Compute probabilities for the SP1 part of the data using MNL model
  mnl_settings_SP1 = list(
    alternatives  = c(EV=1, CV=2), 
    avail         = list(EV=avEV, CV=avCV), 
    choiceVar     = choice, 
    utilities     = list(EV  = mu_SP1*V[["EV"]],
                         CV  = mu_SP1*V[["CV"]]
    ),
    rows          = (SP1==1)
  )
  
  P[["SP1"]] = apollo_mnl(mnl_settings_SP1, functionality)
  
  ### Compute probabilities for the SP2 part of the data using MNL model
  mnl_settings_SP2 = list(
    alternatives  = c(leisure=3, workplace=4, highway=5), 
    avail         = list(leisure=avLeisure, workplace=avWorkplace, highway=avHighway), 
    choiceVar     = choice, 
    utilities     = list(leisure  = mu_SP2*V[["leisure"]],
                         workplace  = mu_SP2*V[["workplace"]],
                         highway  = mu_SP2*V[["highway"]]),
    rows          = (SP2==1)
  )
  
  P[["SP2"]] = apollo_mnl(mnl_settings_SP2, functionality)
  
  ### Combined model
  P = apollo_combineModels(P, apollo_inputs, functionality)
  
  ### Take product across observation for same individual
  P = apollo_panelProd(P, apollo_inputs, functionality)
  
  ### Average across inter-individual draws
  P = apollo_avgInterDraws(P, 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                                               ####
# ################################################################# #

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

apollo_modelOutput(model)
apollo_saveOutput(model)






Re: Correlations across alternatives

Posted: 23 Jan 2025, 17:26
by stephanehess
Hi

my suggested approach here would be to have a joint model (you can look at the RP-SP example on the website) and to add error components to create the correlations you want.

Let me know if this is clear or if you need more help once you've set up an initial attempt

Stephane

Re: Correlations across alternatives

Posted: 26 Jan 2025, 08:35
by ifurqanbhat
Thank you for your reply. I have set up the problem as a joint model for the two SP surveys. I have added error components to take care of panel effects first.

Code: Select all

### Clear memory
rm(list = ls())
# install.packages("apollo")
library(apollo)
database = read.csv(file.choose(),header=TRUE)
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName  ="V1_Joint",
  modelDescr ="Joint SP1/SP2 Model Estimation",
  indivID    ="Person_ID",
  nCores = 32
)

apollo_beta=c(asc_CV                  =  0,
              asc_EV                  =  0,
              asc_L                   =  0,
              asc_W                   =  0,
              asc_H                   =  0,
              
              b_pp                    =  0,
              b_oc                    =  0,
              b_ct                    =  0,
              b_rn                    =  0,
              b_cf                    =  0,
              b_em                    =  0,
              
              b_cc                    =  0,
              b_wt                    =  0,
              
              mu_SP1                       =  1,
              mu_SP2                       =  1,
              sigma_panel                =  0
)

apollo_fixed = c("asc_CV", "asc_H",
                 "mu_SP1")

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

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "mlhs",
  interNDraws    = 1000,
  interUnifDraws = c(),
  interNormDraws = c("draws_alt1", "draws_alt2", "draws_alt3", "draws_alt4", "draws_alt5"),
  
  intraDrawsType = "halton",
  intraNDraws    = 0,
  intraUnifDraws = c(),
  intraNormDraws = c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["ec1"]]  = sigma_panel*draws_alt1
  randcoeff[["ec2"]]  = sigma_panel*draws_alt2
  randcoeff[["ec3"]]  = sigma_panel*draws_alt3
  randcoeff[["ec4"]]  = sigma_panel*draws_alt4
  randcoeff[["ec5"]]  = sigma_panel*draws_alt5
  
  return(randcoeff)
}

# ################################################################# #
#### 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()
  
  ### List of utilities (before applying scales): these must use the same names as in mnl_settings, order is irrelevant
  V = list()
  
  #######################################################################################
  
  V[['EV']]         = asc_EV  +  b_pp*ppEV* +  b_oc*ocEV  +  b_ct*ctEV + b_rn*rEV  + b_cf*(cfEV) + b_em*(eEV) + ec1
  V[['CV']]         = asc_CV  +  b_pp*ppCV* +  b_oc*ocCV               + b_rn*rCV                             + ec2

  V[['leisure']]    = asc_L  +  b_ct*ctLeisure     +  b_cc*ccLeisure    +  b_wt*wtLeisure   + b_em*eLeisure   + ec3
  V[['workplace']]  = asc_W  +  b_ct*ctWorkplace   +  b_cc*ccWorkplace  +  b_wt*wtWorkplace + b_em*eWorkplace + ec4
  V[['highway']]    = asc_H  +  b_ct*ctHighway     +  b_cc*ccHighway    +  b_wt*wtHighway   + b_em*eHighway   + ec5
  
  ### Compute probabilities for the SP1 part of the data using MNL model
  mnl_settings_SP1 = list(
    alternatives  = c(EV=1, CV=2), 
    avail         = list(EV=avEV, CV=avCV), 
    choiceVar     = choice, 
    utilities     = list(EV  = mu_SP1*V[["EV"]],
                             CV  = mu_SP1*V[["CV"]]
    ),
    rows          = (SP1==1)
  )
  
  P[["SP1"]] = apollo_mnl(mnl_settings_SP1, functionality)
  
  ### Compute probabilities for the SP2 part of the data using MNL model
  mnl_settings_SP2 = list(
    alternatives  = c(leisure=3, workplace=4, highway=5), 
    avail         = list(leisure=avLeisure, workplace=avWorkplace, highway=avHighway), 
    choiceVar     = choice, 
    utilities     = list(leisure        = mu_SP2*V[["leisure"]],
                              workplace  = mu_SP2*V[["workplace"]],
                                 highway  = mu_SP2*V[["highway"]]),
    rows          = (SP2==1)
  )
  
  P[["SP2"]] = apollo_mnl(mnl_settings_SP2, functionality)
  
  ### Combined model
  P = apollo_combineModels(P, apollo_inputs, functionality)
  
  ### Take product across observation for same individual
  P = apollo_panelProd(P, apollo_inputs, functionality)
  
  ### Average across inter-individual draws
  P = apollo_avgInterDraws(P, 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                                               ####
# ################################################################# #

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

apollo_modelOutput(model)
apollo_saveOutput(model)
Further to measure correlations between alternatives at two stages of the survey (SP1 and SP2), should I introduce terms between each pair of alternatives? For example, should I introduce one error component to measure correlation between EV (choice of vehicle in the first SP) and Workplace (choice of location in the second SP) and another for EV and leisure place and so on for other pairs as well?

Re: Correlations across alternatives

Posted: 27 Jan 2025, 10:02
by stephanehess
Yes, that's what I would suggest

Re: Correlations across alternatives

Posted: 29 Jan 2025, 04:12
by ifurqanbhat
Dear Prof. Hess,

Thank you for your help. Below is how I have setup the problem. Please check if it is correct.

Code: Select all

### Clear memory
rm(list = ls())
# install.packages("apollo")
library(apollo)
database = read.csv(file.choose(),header=TRUE)
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName  ="V1_Joint",
  modelDescr ="Joint SP1/SP2 Model Estimation",
  indivID    ="Person_ID",
  nCores = 32
)

apollo_beta=c(asc_CV                  =  0,
              asc_EV                  =  0,
              asc_L                   =  0,
              asc_W                   =  0,
              asc_H                   =  0,
              
              b_pp                    =  0,
              b_oc                    =  0,
              b_ct                    =  0,
              b_rn                    =  0,
              b_cf                    =  0,
              b_em                    =  0,
              
              b_cc                    =  0,
              b_wt                    =  0,
              
              mu_SP1                       =  1,
              mu_SP2                       =  1,
              sigma_panel                =  0,
              
              sigma_cor_EV_Work            =  0,
              sigma_cor_EV_Leis            =  0,
              sigma_cor_EV_High            =  0,
              sigma_cor_CV_Work            =  0,
              sigma_cor_CV_Leis            =  0,
              sigma_cor_CV_High            =  0
              
)

apollo_fixed = c("asc_CV", "asc_H",
                 "mu_SP1")

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

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "mlhs",
  interNDraws    = 1000,
  interUnifDraws = c(),
  interNormDraws = c("draws_alt1", "draws_alt2", "draws_alt3", "draws_alt4", "draws_alt5", "draws_cor_EV", "draws_cor_CV"),
  
  intraDrawsType = "halton",
  intraNDraws    = 0,
  intraUnifDraws = c(),
  intraNormDraws = c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["ec1"]]  = sigma_panel*draws_alt1
  randcoeff[["ec2"]]  = sigma_panel*draws_alt2
  randcoeff[["ec3"]]  = sigma_panel*draws_alt3
  randcoeff[["ec4"]]  = sigma_panel*draws_alt4
  randcoeff[["ec5"]]  = sigma_panel*draws_alt5
  
  randcoeff[["cor_EV_Work"]]  = sigma_cor_EV_Work*draws_cor_EV
  randcoeff[["cor_EV_Leis"]]  = sigma_cor_EV_Leis*draws_cor_EV
  randcoeff[["cor_EV_High"]]  = sigma_cor_EV_High*draws_cor_EV
  randcoeff[["cor_CV_Work"]]  = sigma_cor_CV_Work*draws_cor_CV
  randcoeff[["cor_CV_Leis"]]  = sigma_cor_CV_Leis*draws_cor_CV
  randcoeff[["cor_CV_High"]]  = sigma_cor_CV_High*draws_cor_CV
  
  return(randcoeff)
}

# ################################################################# #
#### 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()
  
  ### List of utilities (before applying scales): these must use the same names as in mnl_settings, order is irrelevant
  V = list()
  
  #######################################################################################
  
  V[['EV']]         = asc_EV  +  b_pp*ppEV* +  b_oc*ocEV  +  b_ct*ctEV + b_rn*rEV  + b_cf*(cfEV) + b_em*(eEV) + ec1 + cor_EV_Leis + cor_EV_Work + cor_EV_High
  V[['CV']]         = asc_CV  +  b_pp*ppCV* +  b_oc*ocCV               + b_rn*rCV                             + ec2 + cor_CV_Leis + cor_CV_Work + cor_CV_High

  V[['leisure']]    = asc_L  +  b_ct*ctLeisure     +  b_cc*ccLeisure    +  b_wt*wtLeisure   + b_em*eLeisure   + ec3 + cor_EV_Leis + cor_CV_Leis
  V[['workplace']]  = asc_W  +  b_ct*ctWorkplace   +  b_cc*ccWorkplace  +  b_wt*wtWorkplace + b_em*eWorkplace + ec4 + cor_EV_Work + cor_CV_Work 
  V[['highway']]    = asc_H  +  b_ct*ctHighway     +  b_cc*ccHighway    +  b_wt*wtHighway   + b_em*eHighway   + ec5 + cor_EV_High + cor_CV_High
  
  ### Compute probabilities for the SP1 part of the data using MNL model
  mnl_settings_SP1 = list(
    alternatives  = c(EV=1, CV=2), 
    avail         = list(EV=avEV, CV=avCV), 
    choiceVar     = choice, 
    utilities     = list(EV  = mu_SP1*V[["EV"]],
                             CV  = mu_SP1*V[["CV"]]
    ),
    rows          = (SP1==1)
  )
  
  P[["SP1"]] = apollo_mnl(mnl_settings_SP1, functionality)
  
  ### Compute probabilities for the SP2 part of the data using MNL model
  mnl_settings_SP2 = list(
    alternatives  = c(leisure=3, workplace=4, highway=5), 
    avail         = list(leisure=avLeisure, workplace=avWorkplace, highway=avHighway), 
    choiceVar     = choice, 
    utilities     = list(leisure        = mu_SP2*V[["leisure"]],
                              workplace  = mu_SP2*V[["workplace"]],
                                 highway  = mu_SP2*V[["highway"]]),
    rows          = (SP2==1)
  )
  
  P[["SP2"]] = apollo_mnl(mnl_settings_SP2, functionality)
  
  ### Combined model
  P = apollo_combineModels(P, apollo_inputs, functionality)
  
  ### Take product across observation for same individual
  P = apollo_panelProd(P, apollo_inputs, functionality)
  
  ### Average across inter-individual draws
  P = apollo_avgInterDraws(P, 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                                               ####
# ################################################################# #

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

apollo_modelOutput(model)
apollo_saveOutput(model)

Re: Correlations across alternatives

Posted: 17 Feb 2025, 04:57
by stephanehess
Hi

sorry for the slow reply.

That doesn't seem quite correct. The correlation is introduced by having shared draws. In your case, you are sharing draws between the EV options (and between the CV), not between EV and work, for example

Stephane

Re: Correlations across alternatives

Posted: 03 Mar 2025, 06:37
by ifurqanbhat
Thank you Prof. Hess for your help.
I have incorporated the shared draws to introduce correlated across alternatives. Can you please confirm if it is correct now?

Code: Select all

### Clear memory
rm(list = ls())
# install.packages("apollo")
library(apollo)
database = read.csv(file.choose(),header=TRUE)
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName  ="V1_Joint",
  modelDescr ="Joint SP1/SP2 Model Estimation",
  indivID    ="Person_ID",
  nCores = 32
)

apollo_beta=c(asc_CV                  =  0,
              asc_EV                  =  0,
              asc_L                   =  0,
              asc_W                   =  0,
              asc_H                   =  0,
              
              b_pp                    =  0,
              b_oc                    =  0,
              b_ct                    =  0,
              b_rn                    =  0,
              b_cf                    =  0,
              b_em                    =  0,
              
              b_cc                    =  0,
              b_wt                    =  0,
              
              mu_SP1                       =  1,
              mu_SP2                       =  1,
              sigma_panel                =  0,
              
              sigma_cor_EV_Work            =  0,
              sigma_cor_EV_Leis            =  0,
              sigma_cor_EV_High            =  0,
              sigma_cor_CV_Work            =  0,
              sigma_cor_CV_Leis            =  0,
              sigma_cor_CV_High            =  0
              
)

apollo_fixed = c("asc_CV", "asc_H",
                 "mu_SP1")

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

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "mlhs",
  interNDraws    = 1000,
  interUnifDraws = c(),
  interNormDraws = c("draws_alt1", "draws_alt2", "draws_alt3", "draws_alt4", "draws_alt5", "draws_cor_EV_leis", "draws_cor_EV_work", "draws_cor_EV_high", "draws_cor_CV_leis", "draws_cor_CV_work", "draws_cor_CV_high"),
  
  intraDrawsType = "halton",
  intraNDraws    = 0,
  intraUnifDraws = c(),
  intraNormDraws = c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["ec1"]]  = sigma_panel*draws_alt1
  randcoeff[["ec2"]]  = sigma_panel*draws_alt2
  randcoeff[["ec3"]]  = sigma_panel*draws_alt3
  randcoeff[["ec4"]]  = sigma_panel*draws_alt4
  randcoeff[["ec5"]]  = sigma_panel*draws_alt5
  
  randcoeff[["cor_EV_Leis"]]  = sigma_cor_EV_Leis*draws_cor_EV_leis
  randcoeff[["cor_EV_Work"]]  = sigma_cor_EV_Work*draws_cor_EV_work
  randcoeff[["cor_EV_High"]]  = sigma_cor_EV_High*draws_cor_EV_high
  randcoeff[["cor_CV_Leis"]]  = sigma_cor_CV_Leis*draws_cor_CV_leis
  randcoeff[["cor_CV_Work"]]  = sigma_cor_CV_Work*draws_cor_CV_work
  randcoeff[["cor_CV_High"]]  = sigma_cor_CV_High*draws_cor_CV_high
  
  return(randcoeff)
}

# ################################################################# #
#### 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()
  
  ### List of utilities (before applying scales): these must use the same names as in mnl_settings, order is irrelevant
  V = list()
  
  #######################################################################################
  
  V[['EV']]         = asc_EV  +  b_pp*ppEV* +  b_oc*ocEV  +  b_ct*ctEV + b_rn*rEV  + b_cf*(cfEV) + b_em*(eEV) + ec1 + cor_EV_Leis + cor_EV_Work + cor_EV_High
  V[['CV']]         = asc_CV  +  b_pp*ppCV* +  b_oc*ocCV               + b_rn*rCV                             + ec2 + cor_CV_Leis + cor_CV_Work + cor_CV_High
  
  V[['leisure']]    = asc_L  +  b_ct*ctLeisure     +  b_cc*ccLeisure    +  b_wt*wtLeisure   + b_em*eLeisure   + ec3 + cor_EV_Leis + cor_CV_Leis
  V[['workplace']]  = asc_W  +  b_ct*ctWorkplace   +  b_cc*ccWorkplace  +  b_wt*wtWorkplace + b_em*eWorkplace + ec4 + cor_EV_Work + cor_CV_Work 
  V[['highway']]    = asc_H  +  b_ct*ctHighway     +  b_cc*ccHighway    +  b_wt*wtHighway   + b_em*eHighway   + ec5 + cor_EV_High + cor_CV_High
  
  ### Compute probabilities for the SP1 part of the data using MNL model
  mnl_settings_SP1 = list(
    alternatives  = c(EV=1, CV=2), 
    avail         = list(EV=avEV, CV=avCV), 
    choiceVar     = choice, 
    utilities     = list(EV  = mu_SP1*V[["EV"]],
                         CV  = mu_SP1*V[["CV"]]
    ),
    rows          = (SP1==1)
  )
  
  P[["SP1"]] = apollo_mnl(mnl_settings_SP1, functionality)
  
  ### Compute probabilities for the SP2 part of the data using MNL model
  mnl_settings_SP2 = list(
    alternatives  = c(leisure=3, workplace=4, highway=5), 
    avail         = list(leisure=avLeisure, workplace=avWorkplace, highway=avHighway), 
    choiceVar     = choice, 
    utilities     = list(leisure        = mu_SP2*V[["leisure"]],
                         workplace  = mu_SP2*V[["workplace"]],
                         highway  = mu_SP2*V[["highway"]]),
    rows          = (SP2==1)
  )
  
  P[["SP2"]] = apollo_mnl(mnl_settings_SP2, functionality)
  
  ### Combined model
  P = apollo_combineModels(P, apollo_inputs, functionality)
  
  ### Take product across observation for same individual
  P = apollo_panelProd(P, apollo_inputs, functionality)
  
  ### Average across inter-individual draws
  P = apollo_avgInterDraws(P, 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                                               ####
# ################################################################# #

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

apollo_modelOutput(model)
apollo_saveOutput(model)
The other question I have is about the interpretation of results of correlation coefficients. The results from the model are:

Code: Select all

Estimates:
                          Estimate        s.e.   t.rat.(0)    Rob.s.e. Rob.t.rat.(0)
sigma_panel              -0.250732    0.256126     -0.9789    0.135326       -1.8528
sigma_cor_EV_Work        -0.686001    0.133476     -5.1395    0.108626       -6.3153
sigma_cor_EV_Leis         0.265219    0.274686      0.9655    0.155272        1.7081
sigma_cor_EV_High        -0.431194    0.244756     -1.7617    0.199271       -2.1639
sigma_cor_CV_Work        -0.764715    0.127185     -6.0126    0.111881       -6.8351
sigma_cor_CV_Leis         0.657683    0.160999      4.0850    0.149798        4.3905
sigma_cor_CV_High        -0.132411    0.310501     -0.4264    0.105109       -1.2598

Re: Correlations across alternatives

Posted: 10 Mar 2025, 14:39
by stephanehess
Hi

what you really need to do here is a Cholesky transformation, as that would allow you to capture positive as well as negative correlation. Your specification is simply telling you there is correlation as the same sigma parameters are used in both utilities, so if multipled together they will be positive.

The details for this are in Kenneth Train's book, but you basically need something like this:

U1 = .... + sig1 * draws1
U2 = .... + sig21 * draws1 + sig22 * draws2
...

then the covariance between U1 and U2 is sig1*sig22

Re: Correlations across alternatives

Posted: 21 Mar 2025, 11:26
by ifurqanbhat
Thank you Prof. Hess for your response.

Is this approach correct?

Code: Select all

V[['EV']]      = .......... + sigma_EV_Work_1*draws_EV_Work + sigma_EV_Leis_1*draws_EV_Leis + sigma_EV_High_1*draws_EV_High
V[['CV']]      = .......... + sigma_CV_Work_1*draws_CV_Work + sigma_CV_Leis_1*draws_CV_Leis + sigma_CV_High_1*draws_CV_High
V[['Work']]  = .......... + sigma_EV_Work_2*draws_EV_Work + sigma_CV_Work_2*draws_CV_Work + sigma_Work*draws_Work
V[['Leis']]    = .......... + sigma_EV_Leis_2*draws_EV_Leis + sigma_CV_Leis_2*draws_CV_Leis + sigma_Leis*draws_Leis
V[['High']]   = .......... + sigma_EV_High_2*draws_EV_High + sigma_CV_High_2*draws_CV_High + sigma_High*draws_High
And will the covariance between EV and work then be: sigma_EV_Work_1*sigma_Work
And the covariance between CV and work be: sigma_CV_Work_1*sigma_Work ?

Re: Correlations across alternatives

Posted: 03 Apr 2025, 23:57
by dpalma
Hi,

I would recommend the approach below. Once estimated, you can calculate the covariances between error components based on the elements of the Cholesky decomposition (you can use apollo_deltaMethod to get the s.e. for them).

Best wishes,
David

Code: Select all

### Clear memory
rm(list = ls())
# install.packages("apollo")
library(apollo)
database = read.csv(file.choose(),header=TRUE)
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName  ="V1_Joint",
  modelDescr ="Joint SP1/SP2 Model Estimation",
  indivID    ="Person_ID",
  nCores = 32
)

apollo_beta=c(asc_CV                  =  0,
              asc_EV                  =  0,
              asc_L                   =  0,
              asc_W                   =  0,
              asc_H                   =  0,
              
              b_pp                    =  0,
              b_oc                    =  0,
              b_ct                    =  0,
              b_rn                    =  0,
              b_cf                    =  0,
              b_em                    =  0,
              
              b_cc                    =  0,
              b_wt                    =  0,
              
              mu_SP1                  =  1,
              mu_SP2                  =  1,
              
              s11 = 0, 
              s21 = 0, s22 =2, 
              s31 = 0, s32 = 0, s33 = 0, 
              s41 = 0, s42 = 0, s43 = 0, s44 = 0, 
              s51 = 0, s52 = 0, s53 = 0, s54 = 0, s55 =0
)

apollo_fixed = c("asc_CV", "asc_H",
                 "mu_SP1")

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

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "mlhs",
  interNDraws    = 1000,
  interUnifDraws = c(),
  interNormDraws = c("draws1", "draws2", "draws3", "draws4", "draws5"),
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["ec1"]]  = s11*draws1
  randcoeff[["ec2"]]  = s21*draws1 + s22*draws2
  randcoeff[["ec3"]]  = s31*draws1 + s32*draws2 + s33*draws3
  randcoeff[["ec4"]]  = s41*draws1 + s42*draws2 + s43*draws3 + s44*draws4
  randcoeff[["ec5"]]  = s51*draws1 + s52*draws2 + s53*draws3 + s44*draws4 + s55*draws5
  
  return(randcoeff)
}

# ################################################################# #
#### 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()
  
  ### List of utilities (before applying scales): these must use the same names as in mnl_settings, order is irrelevant
  V = list()
  
  #######################################################################################
  
  V[['EV']]         = asc_EV  +  b_pp*ppEV* +  b_oc*ocEV  +  b_ct*ctEV + b_rn*rEV  + b_cf*(cfEV) + b_em*(eEV) + ec1
  V[['CV']]         = asc_CV  +  b_pp*ppCV* +  b_oc*ocCV               + b_rn*rCV                             + ec2
  
  V[['leisure']]    = asc_L  +  b_ct*ctLeisure     +  b_cc*ccLeisure    +  b_wt*wtLeisure   + b_em*eLeisure   + ec3
  V[['workplace']]  = asc_W  +  b_ct*ctWorkplace   +  b_cc*ccWorkplace  +  b_wt*wtWorkplace + b_em*eWorkplace + ec4
  V[['highway']]    = asc_H  +  b_ct*ctHighway     +  b_cc*ccHighway    +  b_wt*wtHighway   + b_em*eHighway   + ec5
  
  ### Compute probabilities for the SP1 part of the data using MNL model
  mnl_settings_SP1 = list(
    alternatives  = c(EV=1, CV=2), 
    avail         = list(EV=avEV, CV=avCV), 
    choiceVar     = choice, 
    utilities     = list(EV  = mu_SP1*V[["EV"]],
                         CV  = mu_SP1*V[["CV"]]
    ),
    rows          = (SP1==1)
  )
  
  P[["SP1"]] = apollo_mnl(mnl_settings_SP1, functionality)
  
  ### Compute probabilities for the SP2 part of the data using MNL model
  mnl_settings_SP2 = list(
    alternatives  = c(leisure=3, workplace=4, highway=5), 
    avail         = list(leisure=avLeisure, workplace=avWorkplace, highway=avHighway), 
    choiceVar     = choice, 
    utilities     = list(leisure        = mu_SP2*V[["leisure"]],
                         workplace  = mu_SP2*V[["workplace"]],
                         highway  = mu_SP2*V[["highway"]]),
    rows          = (SP2==1)
  )
  
  P[["SP2"]] = apollo_mnl(mnl_settings_SP2, functionality)
  
  ### Combined model
  P = apollo_combineModels(P, apollo_inputs, functionality)
  
  ### Take product across observation for same individual
  P = apollo_panelProd(P, apollo_inputs, functionality)
  
  ### Average across inter-individual draws
  P = apollo_avgInterDraws(P, 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                                               ####
# ################################################################# #

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

apollo_modelOutput(model)
apollo_saveOutput(model)