Latent class (nested) model - elasticity

[Ying Lu]

Hi prof Hess,

I estimated a latent class nested logit model by using apollo_modeChoiceData.csv. I attempted to estimate the elasticity with regards to travel time for each of the four travel options across the two classes. The following is my syntax of elasticity calculations:
############ Class 1 elasticity (car)
predictions_base_car_c1 = apollo_prediction(model,apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=1))
database$time_car_scaled = 1.01*database$time_car_scaled
predictions_new_car_c1 = apollo_prediction(model, apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=1))
database$time_car_scaled = 1/1.01*database$time_car_scaled
log(sum(predictions_new_car_c1[,"car"])/sum(predictions_base_car_c1[,"car"]))/log(1.1)
############# Class 2 Elasticity (car)
predictions_base_car_c2 = apollo_prediction(model,apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=2))
database$time_car_scaled = 1.01*database$time_car_scaled
predictions_new_car_c2 = apollo_prediction(model, apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=2))
database$time_car_scaled = 1/1.01*database$time_car_scaled
log(sum(predictions_new_car_c2[,"car"])/sum(predictions_base_car_c2[,"car"]))/log(1.1)

(1) Could I kindly know is this the correct way to calculate the elasticity of car with regards to travel time across two classes (PS. I run this syntax and it works. The own elasticity for car in class one is -0.1652, and -0.2201 in class two)? (2) When I add runs=30 in prediction_settings, an error occurred showing that: Error in predicitons_base_car_c1[,"car"]: incorrect number of dimensions.

predictions_base_car_c1 = apollo_prediction(model,apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=1,runs=30))

Any suggestions are appreciated! Thank you.

Kind regards,
Ying
----------------------------------------------------------------------------------------------------
The latent class nested logit model syntax is also attached below:
### Clear memory
rm(list = ls())

### Load Apollo library
library(apollo)

### Initialise code
apollo_initialise()

apollo_control = list(
modelName ="Apollo_example_20",
modelDescr ="LC model with class allocation model",
indivID ="ID",
nCores = 1
)
database = read.csv("apollo_modeChoiceData.csv",header = TRUE)
database = subset(database,database$SP==1)
## scaled variables
database$income_scaled = database$income/1000
database$time_car_scaled = database$time_car/100
database$time_bus_scaled = database$time_bus/100
database$time_air_scaled = database$time_air/100
database$time_rail_scaled = database$time_rail/100
database$cost_car_scaled = database$cost_car/100
database$cost_bus_scaled = database$cost_bus/100
database$cost_air_scaled = database$cost_air/100
database$cost_rail_scaled = database$cost_rail/100
database$access_car_scaled = database$access_car/100
database$access_bus_scaled = database$access_bus/100
database$access_air_scaled = database$access_air/100
database$access_rail_scaled = database$access_rail/100
apollo_beta = c(asc_car_c1 = 0,
asc_bus_c1 = 0.0453,
asc_air_c1 = 1.0239,
asc_rail_c1 = -0.8988,
b_tt_car_c1 = -0.9114,
b_tt_bus_c1 = -1.3664,
b_tt_air_c1 = -1.7559,
b_tt_rail_c1 = -0.3124,
b_cost_c1 = -4.8636,
b_access_c1 = -2.1111,
lambda_PT_c1 = 0.8054,
lambda_MOT_c1 = 0.7314,
gamma_gender_c1 = 0,
gamma_income_c1 = 0,
delta_c1 = 0,
asc_car_c2 = 0,
asc_bus_c2 = 0.032,
asc_air_c2 = 0.72,
asc_rail_c2 = -0.63,
b_tt_car_c2 = -0.64,
b_tt_bus_c2 = -0.96,
b_tt_air_c2 = -1.23,
b_tt_rail_c2 = -0.21,
b_cost_c2 = -3.405,
b_access_c2 = 1.478,
lambda_PT_c2 = 0.95,
lambda_MOT_c2 = 0.95,
gamma_gender_c2 = 0,
gamma_income_c2 = 0,
delta_c2 = 0
)
apollo_fixed = c("asc_car_c1","asc_car_c2","gamma_gender_c1","gamma_income_c1","delta_c1")

## Define latent class components
apollo_lcPars = function(apollo_beta, apollo_inputs){
lcpars = list()
lcpars[["b_tt_car"]] = list(b_tt_car_c1,b_tt_car_c2)
lcpars[["b_tt_bus"]] = list(b_tt_bus_c1,b_tt_bus_c2)
lcpars[["b_tt_air"]] = list(b_tt_air_c1,b_tt_air_c2)
lcpars[["b_tt_rail"]] = list(b_tt_rail_c1,b_tt_rail_c2)
lcpars[["b_cost"]] = list(b_cost_c1,b_cost_c2)
lcpars[["b_access"]] = list(b_access_c1, b_access_c2)
lcpars[["asc_car"]] = list(asc_car_c1,asc_car_c2)
lcpars[["asc_bus"]] = list(asc_bus_c1,asc_bus_c2)
lcpars[["asc_air"]] = list(asc_air_c1,asc_air_c2)
lcpars[["asc_rail"]] = list(asc_rail_c1,asc_rail_c2)
V = list()
V[["class_1"]] = delta_c1 + gamma_gender_c1*female + gamma_income_c1*income_scaled
V[["class_2"]] = delta_c2 + gamma_gender_c2*female + gamma_income_c2*income_scaled
mnl_settings = list(
alternatives = c(class_1=1, class_2=2),
avail = 1,
choiceVar = NA,
V = V
)
lcpars[["pi_values"]] = apollo_mnl(mnl_settings, functionality = "raw")
lcpars[["pi_values"]] = apollo_firstRow(lcpars[["pi_values"]], apollo_inputs)
return(lcpars)
}

### Group and validate inputs
apollo_inputs = apollo_validateInputs()
### Define model and likelihood function

apollo_probabilities = function(apollo_beta, apollo_inputs, functionality = "estimate"){
## 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()
### Specify nests for NL model
nlNests = list()
nlNests[[1]] = list(root=1, PT_c1=lambda_PT_c1, MOT_c1 = lambda_MOT_c1)
nlNests[[2]] = list(root=1,PT_c2=lambda_PT_c2, MOT_c2 = lambda_MOT_c2)
nlStructure = list()
nlStructure[[1]] = list()
nlStructure[[2]] = list()
nlStructure[[1]][["root"]] = c("PT_c1","MOT_c1")
nlStructure[[2]][["root"]] = c("PT_c2","MOT_c2")

nlStructure[[1]][["PT_c1"]] = c("air","rail")
nlStructure[[1]][["MOT_c1"]] = c("car","bus")
nlStructure[[2]][["PT_c2"]] = c("air","rail")
nlStructure[[2]][["MOT_c2"]] = c("car","bus")

nl_settings <- list(
alternatives = c(car=1, bus=2, air=3, rail=4),
avail = list(car=av_car, bus=av_bus, air=av_air, rail=av_rail),
choiceVar = choice,
nlNests = nlNests,
nlStructure = nlStructure
)
### Loop over classes
s=1
while(s<=2){
### Compute class-specific utilities
V=list()
V[['car']] = asc_car[[s]] + b_tt_car[[s]] * time_car_scaled + b_cost[[s]] * cost_car_scaled
V[['bus']] = asc_bus[[s]] + b_tt_bus[[s]]*time_bus_scaled + b_cost[[s]]*cost_bus_scaled + b_access[[s]]*access_bus_scaled
V[['air']] = asc_air[[s]] + b_tt_air[[s]]*time_air_scaled + b_cost[[s]]* cost_air_scaled + b_access[[s]]*access_air_scaled
V[['rail']] = asc_rail[[s]] + b_tt_rail[[s]]*time_rail_scaled + b_cost[[s]]*cost_rail_scaled + b_access[[s]]*access_rail_scaled
nl_settings$V = V
nl_settings$nlNests = nlNests[[s]]
nl_settings$nlStructure = nlStructure[[s]]
P[[s]] = apollo_nl(nl_settings, functionality)
P[[s]] = apollo_panelProd(P[[s]],apollo_inputs,functionality)

s=s+1
}
lc_settings = list(inClassProb = P, classProb=pi_values)
P[["model"]] = apollo_lc(lc_settings,apollo_inputs,functionality)

P = apollo_prepareProb(P,apollo_inputs,functionality)
}

model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities,apollo_inputs)
apollo_modelOutput(model,modelOutput_settings = list(printPVal=TRUE))
############ Class 1 elasticity (car)
predictions_base_car_c1 = apollo_prediction(model,apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=1))
database$time_car_scaled = 1.01*database$time_car_scaled
predictions_new_car_c1 = apollo_prediction(model, apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=1))
database$time_car_scaled = 1/1.01*database$time_car_scaled
log(sum(predictions_new_car_c1[,"car"])/sum(predictions_base_car_c1[,"car"]))/log(1.1)
############# Class 2 Elasticity (car)
predictions_base_car_c2 = apollo_prediction(model,apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=2))
database$time_car_scaled = 1.01*database$time_car_scaled
predictions_new_car_c2 = apollo_prediction(model, apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=2))
database$time_car_scaled = 1/1.01*database$time_car_scaled
log(sum(predictions_new_car_c2[,"car"])/sum(predictions_base_car_c2[,"car"]))/log(1.1)

[stephanehess]

Dear Ying

please have a look at the manual for further details on this.

Your code is correct but when you include runs, predictions becomes a list rather than just a matrix.

So now, predictions_base_car_c1 has two components, namely:

predictions_base_car_c1[["at_estimates"]]

and

predictions_base_car_c1[["draws"]]

What you're looking for is

predictions_base_car_c1[["at_estimates"]][,"car"]

Best wishes

Stephane

[Ying Lu]

Dear Prof Hess,

Thank you very much for your reply. I still have one more question about this issue.

When I put modelComponent=3 within prediction_settings, I think Apollo will show an error because I do not have three classes according to the code. However, it works and I get a value which is -0.175933. An error happens until I set modelComponent=4. Could I kindly know why this happens?

predictions_base_car_c1 = apollo_prediction(model,apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=3))
database$time_car_scaled = 1.01*database$time_car_scaled
predictions_new_car_c1 = apollo_prediction(model, apollo_probabilities, apollo_inputs,prediction_settings = list(modelComponent=3))
database$time_car_scaled = 1/1.01*database$time_car_scaled
log(sum(predictions_new_car_c1[,"car"])/sum(predictions_base_car_c1[,"car"]))/log(1.1)

Thank you for your helps.

Kind regards,
Ying

[stephanehess]

Ying

you have 3 components. The third component is your actual latent class model, so the full model.

Stephane

[Ying Lu]

Dear Prof Hess,

Now I understand. Thank you for your answers!

Ying