mixed logit query

[georkapo]

Dear Stephane, Dear David,

Thank you!
I'm not sure if I need to create a new post or if I could ask here. sorry.

I also want to run a mixed logit model for the dutch sample.
The aim of my analysis is to investigate wtp for different functionalities on bicycles.
I found some deterministic heterogeneity in the data MNL model and used the same interactions, which are significant for the mixed logit model.

Here is the mixed logit model.

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

### Load Apollo library
library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
modelName = "Model_MNL_IATBR",
modelDescr = "MNL model with socios",
indivID = "ID",
mixing = TRUE,
nCores = 4
)

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

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

### Use only RP data
database = subset(database,database$country_code=='NL')

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

### Vector of parameters, including any that are kept fixed in estimation
apollo_beta=c(asc_A = 0,
asc_B = 0,
asc_C = 0,
mu_b_Assistance_emergency_unit = 0,
mu_b_Assistance_emergency_call = 0.065319,
mu_b_Assistance_smart_routes = 0.305141,
mu_b_Automatic_speed_safe_distance = 0,
mu_b_Automatic_speed_speed_limit = 0.311208,
mu_b_Automatic_speed_risky_areas = 0.175928,
mu_b_Collision_avoidance_front_side = 0,
mu_b_Collision_avoidance_rear_side = 0.144246,
mu_b_Collision_avoidance_blind_spot = 0.164900,
mu_b_cost = -2.70885 ,#-0.001955,

sigma_b_Assistance_emergency_unit = 0,
sigma_b_Assistance_emergency_call = 0,
sigma_b_Assistance_smart_routes = 0,
sigma_b_Automatic_speed_safe_distance = 0,
sigma_b_Automatic_speed_speed_limit = 0,
sigma_b_Automatic_speed_risky_areas = 0,
sigma_b_Collision_avoidance_front_side = 0,
sigma_b_Collision_avoidance_rear_side = 0,
sigma_b_Collision_avoidance_blind_spot = 0,
sigma_b_cost = 0 ,

age18_39_Automatic_speed_speed_limit = 0,
age18_39_Automatic_speed_risky_areas = 0,
age18_39_Automatic_speed_safe_distance = 0,

male_Automatic_speed_speed_limit = 0,
male_Automatic_speed_risky_areas = 0,
#male_Automatic_speed_safe_distance = 0,

male_cost = 0,
young_cost =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("asc_A",
"mu_log_b_Assistance_emergency_unit",
"mu_log_b_Automatic_speed_safe_distance",
"mu_log_b_Collision_avoidance_front_side")

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

### Set parameters for generating draws
apollo_draws = list(
interDrawsType = "mlhs", #halton
interNDraws = 500,
interUnifDraws = c(),
interNormDraws = c("draws_Ass_unit", "draws_Ass_call", "draws_Ass_routes", "draws_Aut_dist", "draws_Aut_limit", "draws_Aut_areas",
"draws_coll_front", "draws_coll_rear", "draws_coll_blind", "draws_cost"),
# intraDrawsType = "halton",
intraNDraws = 0,
intraUnifDraws = c(),
intraNormDraws = c()
)

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

randcoeff[["Assistance_emergency_unit"]] = mu_b_Assistance_emergency_unit + sigma_b_Assistance_emergency_unit * draws_Ass_unit
randcoeff[["Assistance_emergency_call"]] = mu_b_Assistance_emergency_call + sigma_b_Assistance_emergency_call * draws_Ass_call
randcoeff[["Assistance_smart_routes"]] = mu_b_Assistance_smart_routes + sigma_b_Assistance_smart_routes * draws_Ass_routes
randcoeff[["Automatic_speed_safe_distance"]] = mu_b_Automatic_speed_safe_distance + sigma_b_Automatic_speed_safe_distance * draws_Aut_dist
randcoeff[["Automatic_speed_speed_limit"]] = mu_b_Automatic_speed_speed_limit + sigma_b_Automatic_speed_speed_limit * draws_Aut_limit
randcoeff[["Automatic_speed_risky_areas"]] = mu_b_Automatic_speed_risky_areas + sigma_b_Automatic_speed_risky_areas * draws_Aut_areas
randcoeff[["Collision_avoidance_front_side"]] = mu_b_Collision_avoidance_front_side + sigma_b_Collision_avoidance_front_side * draws_coll_front
randcoeff[["Collision_avoidance_rear_side"]] = mu_b_Collision_avoidance_rear_side + sigma_b_Collision_avoidance_rear_side * draws_coll_rear
randcoeff[["Collision_avoidance_blind_spot"]] = mu_b_Collision_avoidance_blind_spot + sigma_b_Collision_avoidance_blind_spot * draws_coll_blind
randcoeff[["b_cost"]] = -exp( mu_b_cost + sigma_b_cost * draws_cost )

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: these must use the same names as in mnl_settings, order is irrelevant
V = list()
V[["A"]] = (asc_A
+ Assistance_emergency_unit * (Att1== 1)
+ Assistance_emergency_call * (Att1== 2)
+ Assistance_smart_routes * (Att1== 3)
+ (Automatic_speed_speed_limit
+ age18_39_Automatic_speed_speed_limit * (Age25_44 == 1)
+ male_Automatic_speed_speed_limit * (gender == 1) ) * (Att2== 1)
+ (Automatic_speed_risky_areas
+ age18_39_Automatic_speed_risky_areas * (Age25_44 == 1)
+ male_Automatic_speed_risky_areas * (gender == 1) ) * (Att2== 2)
+ (Automatic_speed_safe_distance
+ age18_39_Automatic_speed_safe_distance * (Age25_44 == 1)) * (Att2== 3)
+ Collision_avoidance_front_side * (Att3== 1)
+ Collision_avoidance_rear_side * (Att3== 2)
+ Collision_avoidance_blind_spot * (Att3== 3)
+ (b_cost
+ male_cost * (gender == 1)
+ young_cost * (Age25_44 == 1)) * ((Att4== 1) * 400
+(Att4== 2) * 600
+(Att4== 3) * 800
+(Att4== 4) * 1000))

V[["B"]] = (asc_B
+ Assistance_emergency_unit * (Att1.1== 1)
+ Assistance_emergency_call * (Att1.1== 2)
+ Assistance_smart_routes * (Att1.1== 3)
+ (Automatic_speed_speed_limit
+ age18_39_Automatic_speed_speed_limit * (Age25_44 == 1)
+ male_Automatic_speed_speed_limit * (gender == 1) ) * (Att2.1== 1)
+ (Automatic_speed_risky_areas
+ age18_39_Automatic_speed_risky_areas * (Age25_44 == 1)
+ male_Automatic_speed_risky_areas * (gender == 1) ) * (Att2.1== 2)
+ (Automatic_speed_safe_distance
+ age18_39_Automatic_speed_safe_distance * (Age25_44 == 1)) * (Att2.1== 3)
+ Collision_avoidance_front_side * (Att3.1== 1)
+ Collision_avoidance_rear_side * (Att3.1== 2)
+ Collision_avoidance_blind_spot * (Att3.1== 3)
+ (b_cost
+ male_cost * (gender == 1)
+ young_cost * (Age25_44 == 1)) * ((Att4.1== 1) * 400
+(Att4.1== 2) * 600
+(Att4.1== 3) * 800
+(Att4.1== 4) * 1000))

V[["C"]] = (asc_C)

### Define settings for MNL model component
mnl_settings = list(
alternatives = c(A=1, B=2, C=3),
##avail = list(alt1=Concept, alt2=Concept.1),
choiceVar = choice,
utilities = V
)

### Define settings for MNL model component
mnl_settings = list(
alternatives = c(A=1, B=2, C=3),
##avail = list(alt1=Concept, alt2=Concept.1),
choiceVar = choice,
utilities = V
)

### Compute probabilities using MNL model
P[["model"]] = apollo_mnl(mnl_settings, functionality)

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

### Average across inter-individual draws
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)


# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO FILE, using model name) ----
# ----------------------------------------------------------------- #

apollo_saveOutput(model)

I use the estimates from the MNL model as starting values for the mixed logit model. However, since the values for cost do not work I use as starting values for the mu_b_cost=log(abs(b_cost_mnl)).
Also, I use a negative lognormal distribution for the cost attribute, while for all the rest, I use a normal distribution since I assume they can take any value from positive to negative.
To calculate the wtp I take the ratio of the utilities of each attribute with the cost attribute. Is this the right way to do so? I am asking since the estimate for the cost becomes really high, affecting the willingness to pay (becomes really small), which seems meaningless.
I have also run a latent class model and there is a huge difference in wtp between these two models.

Thank you in advance for your help.

Best regards,
Georgios

[stephanehess]

Hi

I split your post and started a new topic for you as this query is about mixed logit.

Your use of starting values is correct.

in relation to WTP, how are you calculating them? Can you show us the code or an example calculation?

Stephane

[georkapo]

Hi Stephane,

Thank you.

I first calculate the unconditional and then estimate the wtp separately for each random coefficient by the cost. I also use a negative sign since I use negative lognormal distribution for the cost. is this correct?
Here is the code I use.

##uncoditionals
unconditionals <- apollo_unconditionals(model,apollo_probabilities, apollo_inputs)

wtp_routes <- -(unconditionals[["Assistance_smart_routes"]]/unconditionals[["b_cost"]])
(mean(wtp_routes)); (sd(wtp_routes))

the results look like

> (mean(wtp_routes)); (sd(wtp_routes))
[1] 140.8361
[1] 694.6622

To calculate the wtp for the rest of the attributes, since I found some deterministic heterogeneity (the interactions), which I have not defined as random coefficients, I take the ratio, e.g. male_Automatic_speed_speed_limit/male_cost. Is this the correct way to estimate them?

Thank you in advance.
Best wishes,
Georgios

[stephanehess]

Georgios

you don't need the negative sign as that is already in apollo_randCoeff and thus reflected in apollo_unconditionals.

But the bit about cost_male is not correct as you're ignoring the other random part.

So it would need to be e.g.

Code: Select all

(unconditionals[["Assistance_smart_routes"]]+model$estimate["male_Automatic_speed_speed_limit"])/(unconditionals[["b_cost"]]+model$estimate["male_cost"])

It would also be better if you included those interactions inside the random coefficient in apollo_randcoeff rather than in the utilities, as you would then guarantee the right sign for the cost coefficient

[georkapo]

Hi Stephane,

Thank you for the suggestions.

One point for the negative sign. If I do not use it, then my wtp is negative. is this correct?

Regarding the cost_male I am getting something like this
[1] 1.504081e+18
[1] 3.999951e+20

Does this make sense?

Also, when I include the interactions in the random coefficient I get NA.
Any idea why is this?

Thank you in advance.

best regards,
Georgios

[stephanehess]

Hi

for wtp and signs, you just need to interpret it on the basis of the signs of the two individual params. So the MRS between two undesirable attributes is positive, which means that the changes need to be in opposite directions (so positive WTP for decreases in time). The MRS between a desirable and an undesirable attribute is negative, which shows that increases in the desirable attribute have the same impact as decreases in cost. So the WTP for increases is positive. You can change the sign in the calculation of course if you want, but best to understand the above.

I can't answer your question on male without seeing the full estimation results

Stephane

[georkapo]

Hi Stephane,

This helped. thank you.

Here is the code when I include the interactions in the random coefficient, I get NA

Code: Select all

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

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

### Load Apollo library
library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "Model_neg_logl",
  modelDescr      = "MMNL model",
  indivID         = "ID",
  mixing          = TRUE,
  nCores          = 4

)

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

### Loading data from package
### if data is to be loaded from a file (e.g. called data.csv), 
### the code would be: database = read.csv("data.csv",header=TRUE)
#database = apollo_modeChoiceData
database = read.csv("read.csv",header=TRUE)
### for data dictionary, use ?apollo_modeChoiceData


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

### Vector of parameters, including any that are kept fixed in estimation
apollo_beta=c(asc_A = 0,
              asc_B = 0,
              asc_C = 0,
              mu_b_Assistance_emergency_unit = 0.03,
              mu_b_Assistance_emergency_call = 0.065665,
              mu_b_Assistance_smart_routes = 0.306003,
              mu_b_Automatic_speed_safe_distance = 0.1,
              mu_b_Automatic_speed_speed_limit = 0.291569,
              mu_b_Automatic_speed_risky_areas = 0.156978,
              mu_b_Collision_avoidance_front_side = 0.1,
              mu_b_Collision_avoidance_rear_side = 0.144932,
              mu_b_Collision_avoidance_blind_spot = 0.167036,
              mu_log_b_cost = -2.6704 , #-0.001955,
              
              mu_age18_39_Automatic_speed_speed_limit = 0.03,
              mu_age18_39_Automatic_speed_risky_areas = 0.03,
              mu_age18_39_Automatic_speed_safe_distance = 0.03,
              mu_male_Automatic_speed_speed_limit = 0.03,
              mu_male_Automatic_speed_risky_areas = 0.03,
              mu_male_cost    = -2.6704,

              
              sigma_b_Assistance_emergency_unit = 0,
              sigma_b_Assistance_emergency_call = 0,
              sigma_b_Assistance_smart_routes = 0,
              sigma_b_Automatic_speed_safe_distance = 0,
              sigma_b_Automatic_speed_speed_limit = 0,
              sigma_b_Automatic_speed_risky_areas = 0,
              sigma_b_Collision_avoidance_front_side = 0,
              sigma_b_Collision_avoidance_rear_side = 0,
              sigma_b_Collision_avoidance_blind_spot = 0,
              sigma_log_b_cost = 0 ,
              
              sigma_age18_39_Automatic_speed_speed_limit = 0,
              sigma_age18_39_Automatic_speed_risky_areas = 0,
              sigma_age18_39_Automatic_speed_safe_distance = 0,
              sigma_male_Automatic_speed_speed_limit = 0, 
              sigma_male_Automatic_speed_risky_areas = 0,
              sigma_male_cost = 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("asc_C",
                 "mu_b_Assistance_emergency_unit",
                 "mu_b_Automatic_speed_safe_distance",
                 "mu_b_Collision_avoidance_front_side")
# ################################################################# #
#### DEFINE RANDOM COMPONENTS                                    ####
# ################################################################# #

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "mlhs", #halton
  interNDraws    = 100,
  interUnifDraws = c(),
  interNormDraws = c("draws_Ass_unit", "draws_Ass_call", "draws_Ass_routes", "draws_Aut_dist", "draws_Aut_limit", "draws_Aut_areas",
                     "draws_coll_front", "draws_coll_rear", "draws_coll_blind", "draws_cost",
                     "draws_age18_39_Automatic_speed_speed_limit", "draws_age18_39_Automatic_speed_risky_areas",
                     "draws_male_Automatic_speed_speed_limit", "draws_ale_Automatic_speed_speed_limit", "draws_male_cost", "draws_young_cost"
                     ),
  # intraDrawsType = "halton",
  intraNDraws    = 0,
  intraUnifDraws = c(),
  intraNormDraws = c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["Assistance_emergency_unit"]] =  mu_b_Assistance_emergency_unit + sigma_b_Assistance_emergency_unit * draws_Ass_unit
  randcoeff[["Assistance_emergency_call"]] = mu_b_Assistance_emergency_call + sigma_b_Assistance_emergency_call * draws_Ass_call
  randcoeff[["Assistance_smart_routes"]] =  mu_b_Assistance_smart_routes + sigma_b_Assistance_smart_routes * draws_Ass_routes
  randcoeff[["Automatic_speed_safe_distance"]] =  mu_b_Automatic_speed_safe_distance + sigma_b_Automatic_speed_safe_distance * draws_Aut_dist
  randcoeff[["Automatic_speed_speed_limit"]] =  mu_b_Automatic_speed_speed_limit + sigma_b_Automatic_speed_speed_limit * draws_Aut_limit
  randcoeff[["Automatic_speed_risky_areas"]] =  mu_b_Automatic_speed_risky_areas + sigma_b_Automatic_speed_risky_areas * draws_Aut_areas
  randcoeff[["Collision_avoidance_front_side"]] =  mu_b_Collision_avoidance_front_side + sigma_b_Collision_avoidance_front_side * draws_coll_front
  randcoeff[["Collision_avoidance_rear_side"]] =  mu_b_Collision_avoidance_rear_side + sigma_b_Collision_avoidance_rear_side * draws_coll_rear
  randcoeff[["Collision_avoidance_blind_spot"]] =  mu_b_Collision_avoidance_blind_spot + sigma_b_Collision_avoidance_blind_spot * draws_coll_blind
  randcoeff[["b_cost"]] = -exp( mu_log_b_cost + sigma_log_b_cost * draws_cost )
  randcoeff[["age18_39_Automatic_speed_speed_limit"]] = mu_age18_39_Automatic_speed_speed_limit + sigma_age18_39_Automatic_speed_speed_limit * draws_age18_39_Automatic_speed_speed_limit
  randcoeff[["age18_39_Automatic_speed_risky_areas"]] = mu_age18_39_Automatic_speed_risky_areas + sigma_age18_39_Automatic_speed_risky_areas * draws_age18_39_Automatic_speed_risky_areas
  randcoeff[["age18_39_Automatic_speed_safe_distance"]] = mu_age18_39_Automatic_speed_safe_distance + sigma_age18_39_Automatic_speed_safe_distance * draws_male_Automatic_speed_speed_limit
  randcoeff[["male_Automatic_speed_speed_limit"]] = mu_male_Automatic_speed_speed_limit + sigma_male_Automatic_speed_speed_limit * draws_male_Automatic_speed_speed_limit
  randcoeff[["male_Automatic_speed_risky_areas"]] = mu_male_Automatic_speed_risky_areas + sigma_male_Automatic_speed_risky_areas * draws_ale_Automatic_speed_speed_limit
  randcoeff[["male_cost"]] = -exp( mu_male_cost + sigma_male_cost * draws_male_cost )
  

  
  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()
  
  
  ######create coefficients using interaction with socio

  
  ### List of utilities: these must use the same names as in mnl_settings, order is irrelevant
  V = list()
  V[["A"]]  = (asc_A
               + Assistance_emergency_unit              * (Att1== 1) 
               + Assistance_emergency_call     * (Att1== 2) 
               + Assistance_smart_routes                * (Att1== 3)   
               + (Automatic_speed_speed_limit                       
                  + age18_39_Automatic_speed_speed_limit  * (Age25_44 == 1)
                  + male_Automatic_speed_speed_limit     * (gender == 1) )       * (Att2== 1) 
               + (Automatic_speed_risky_areas 
                  + age18_39_Automatic_speed_risky_areas  * (Age25_44 == 1)
                  + male_Automatic_speed_risky_areas     *   (gender == 1) )           * (Att2== 2) 
               + (Automatic_speed_safe_distance 
                  + age18_39_Automatic_speed_safe_distance  * (Age25_44 == 1))         * (Att2== 3) 
               + Collision_avoidance_front_side         * (Att3== 1) 
               + Collision_avoidance_rear_side       *  (Att3== 2) 
               + Collision_avoidance_blind_spot  * (Att3== 3) 
               + (b_cost 
                  + male_cost      * (gender == 1)) *  ((Att4== 1) * 400
                                                                     +(Att4== 2) * 600
                                                                     +(Att4== 3) * 800
                                                                     +(Att4== 4) * 1000))

  
  V[["B"]]  = (asc_B  
               + Assistance_emergency_unit        * (Att1.1== 1) 
               + Assistance_emergency_call   * (Att1.1== 2) 
               + Assistance_smart_routes             * (Att1.1== 3) 
               + (Automatic_speed_speed_limit 
                  + age18_39_Automatic_speed_speed_limit * (Age25_44 == 1)
                  + male_Automatic_speed_speed_limit     * (gender == 1) )         * (Att2.1== 1) 
               + (Automatic_speed_risky_areas 
                  + age18_39_Automatic_speed_risky_areas * (Age25_44 == 1)
                  + male_Automatic_speed_risky_areas     * (gender == 1) )       * (Att2.1== 2) 
               + (Automatic_speed_safe_distance  
                  + age18_39_Automatic_speed_safe_distance * (Age25_44 == 1))         * (Att2.1== 3) 
               + Collision_avoidance_front_side       *   (Att3.1== 1) 
               + Collision_avoidance_rear_side      *  (Att3.1== 2) 
               + Collision_avoidance_blind_spot * (Att3.1== 3)
               + (b_cost 
                  + male_cost      * (gender == 1)) *  ((Att4.1== 1) * 400
                                                                     +(Att4.1== 2) * 600
                                                                     +(Att4.1== 3) * 800
                                                                     +(Att4.1== 4) * 1000))

  
  
  V[["C"]]  = (asc_C)
  
  ### Define settings for MNL model component
  mnl_settings = list(
    alternatives  = c(A=1, B=2, C=3), 
    ##avail         = list(alt1=Concept, alt2=Concept.1), 
    choiceVar     = choice,
    utilities     = V
  )
  
  ### Compute probabilities using MNL model
  P[["model"]] = apollo_mnl(mnl_settings, functionality)
  

  ### Take product across observation for same individual
  P = apollo_panelProd(P, apollo_inputs, functionality)
  
  ### Average across inter-individual draws
  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)

here is the estimation for the male related question (it is a different model from the one above, as I used the interaction on the utilities):

Code: Select all

Estimates:
                                          Estimate        s.e.   t.rat.(0)    Rob.s.e. Rob.t.rat.(0)
asc_A                                     2.153353     0.22820      9.4361     0.27584        7.8066
asc_B                                     2.213640     0.23576      9.3894     0.28563        7.7499
asc_C                                     0.000000          NA          NA          NA            NA
mu_b_Assistance_emergency_unit            0.030000          NA          NA          NA            NA
mu_b_Assistance_emergency_call            0.160284     0.09166      1.7487     0.08312        1.9284
mu_b_Assistance_smart_routes              0.568065     0.12157      4.6729     0.12180        4.6637
mu_b_Automatic_speed_safe_distance        0.100000          NA          NA          NA            NA
mu_b_Automatic_speed_speed_limit          0.554076     0.16194      3.4214     0.17148        3.2311
mu_b_Automatic_speed_risky_areas          0.346580     0.16550      2.0941     0.16906        2.0501
mu_b_Collision_avoidance_front_side       0.100000          NA          NA          NA            NA
mu_b_Collision_avoidance_rear_side        0.337218     0.09627      3.5027     0.09251        3.6452
mu_b_Collision_avoidance_blind_spot       0.344706     0.10559      3.2647     0.10520        3.2766
mu_log_b_cost                            -4.901457     0.08024    -61.0823     0.08429      -58.1514
sigma_b_Assistance_emergency_unit         0.626910     0.19494      3.2159     0.20613        3.0413
sigma_b_Assistance_emergency_call         0.438576     0.24915      1.7603     0.27236        1.6103
sigma_b_Assistance_smart_routes           1.380519     0.17587      7.8494     0.19638        7.0298
sigma_b_Automatic_speed_safe_distance    -0.860151     0.18737     -4.5908     0.19865       -4.3300
sigma_b_Automatic_speed_speed_limit      -0.556493     0.20036     -2.7775     0.19941       -2.7906
sigma_b_Automatic_speed_risky_areas      -0.729976     0.16935     -4.3104     0.17134       -4.2605
sigma_b_Collision_avoidance_front_side    0.764143     0.16729      4.5677     0.16398        4.6599
sigma_b_Collision_avoidance_rear_side     0.431224     0.26348      1.6366     0.26701        1.6150
sigma_b_Collision_avoidance_blind_spot    0.941558     0.16871      5.5808     0.17579        5.3561
sigma_log_b_cost                          1.158604     0.07966     14.5452     0.08257       14.0326
age18_39_Automatic_speed_speed_limit      1.409987     0.40278      3.5006     0.52463        2.6876
age18_39_Automatic_speed_risky_areas      1.720870     0.40945      4.2029     0.52805        3.2589
age18_39_Automatic_speed_safe_distance    1.817709     0.47494      3.8272     0.58679        3.0977
male_Automatic_speed_speed_limit         -0.707538     0.18182     -3.8913     0.19572       -3.6151
male_Automatic_speed_risky_areas         -0.418740     0.18599     -2.2515     0.19757       -2.1194
male_cost                                 0.001868  4.0984e-04      4.5571  4.6997e-04        3.9740

Thank you for the help.

Best regards,
Georgios

[stephanehess]

Hi

so the issue in your model is that you are explictly forcing men to be more cost sensitive than women as you are using:

Code: Select all

randcoeff[["b_cost"]] = -exp( mu_log_b_cost + sigma_log_b_cost * draws_cost )
randcoeff[["male_cost"]] = -exp( mu_male_cost + sigma_male_cost * draws_male_cost )

and then:

Code: Select all

 (b_cost 
                  + male_cost      * (gender == 1)) *  ((Att4== 1) * 400
                                                                     +(Att4== 2) * 600
                                                                     +(Att4== 3) * 800
                                                                     +(Att4== 4) * 1000))

it would be much better to use the following:

Code: Select all

randcoeff[["b_cost"]] = -exp( mu_log_b_cost + sigma_log_b_cost * draws_cost + (gender == 1) * ( mu_male_cost + sigma_male_cost * draws_male_cost))

that way, you are using a single random cost coefficient, but allow for differences between men and women in the mean and in the heterogeneity

Stephane

[georkapo]

Dear Stephane,

Thank you for your recommendation.
I created a single random coefficient you suggested.

randcoeff[["b_cost"]] = -exp( mu_log_b_cost + sigma_log_b_cost * draws_cost + (gender == 1) * ( mu_male_cost + sigma_male_cost * draws_male_cost))

Now, I am trying to get the full distribution of the random coefficient in order to estimate the WTP. I use apollo_unconditionals for this, however, I am not sure if the code below is the correct one.

wtp_male_Automatic_limit = (unconditionals[["Automatic_speed_speed_limit"]]/(unconditionals[["b_cost"]]+(model$estimate["mu_male_cost"]/model$estimate["sigma_male_cost"])))

Could you please help me with this?

Thank you in advance.
Best regards,
Georgios

[georkapo]

Dear Stephane,

One additional question related to the previous one.
When I use a single random coefficient allowing for differences between men and women, there is a higher sd error for the attributes e.g.

> (mean(wtp_unit)); (sd(wtp_unit))
[1] 12.42137
[1] 587.157

while when I do not allow for differences between men and women, the sd is lower

> (mean(wtp_unit)); (sd(wtp_unit))
[1] 8.141867
[1] 360.8153

Thank you very much for your support.
Best regards,
Georgios