Specifying HB mixed logit with heterogeneity in mean/variance
Posted: 13 Jun 2023, 22:13
Hello:
I'm working on a mixed logit model with heterogeneity in mean and variance. The model is also estimated in WTP space. The idea is to estimate it using hierarchical Bayes. This is something I started working on last year using Stan and decided to try switching it to Apollo because the Stan model was giving me trouble. We're testing it on synthetic data at the moment.
I've been working off a few examples and gotten it to the below point. The idea is to incorporate input from [url]another forum post[viewtopic.php?p=107&hilit=heterogeneity+means#p107
] on heterogeneity in mean/variance using MSL and syntax from the HB hybrid choice example[http://www.apollochoicemodelling.com/fi ... bayesian.r]
The code runs fine with just the sigma wrt income effect (sig_tc_income), but it gives an error when I try to add a fixed mean that varies by race (mu_raceX). The error is
I'm working on a mixed logit model with heterogeneity in mean and variance. The model is also estimated in WTP space. The idea is to estimate it using hierarchical Bayes. This is something I started working on last year using Stan and decided to try switching it to Apollo because the Stan model was giving me trouble. We're testing it on synthetic data at the moment.
I've been working off a few examples and gotten it to the below point. The idea is to incorporate input from [url]another forum post[viewtopic.php?p=107&hilit=heterogeneity+means#p107
] on heterogeneity in mean/variance using MSL and syntax from the HB hybrid choice example[http://www.apollochoicemodelling.com/fi ... bayesian.r]
The code runs fine with just the sigma wrt income effect (sig_tc_income), but it gives an error when I try to add a fixed mean that varies by race (mu_raceX). The error is
Preparing user-defined functions.
Testing probability function (apollo_probabilities)
Error in apollo_probabilities(apollo_test_beta, apollo_inputs, functionality = "validate") :
object '"validate"' not found
Code: Select all
### Initialise code
apollo_initialise()
### Set core controls
apollo_control = list(
modelName = "MMNL_wtp_space_inter_HB",
modelDescr = "Mixed logit model on synthetic data using HB",
indivID = "personid",
HB = TRUE,
nCores = 4,
analyticGrad = TRUE, # Needs to be set explicitly for models with inter-intra draws (requires extra RAM).
outputDirectory = "output"
)
# ################################################################# #
#### 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 = read.csv("test_data.csv",header=TRUE)
### for data dictionary, use ?apollo_swissRouteChoiceData
# ################################################################# #
#### DEFINE MODEL PARAMETERS ####
# ################################################################# #
### Vector of parameters, including any that are kept fixed in estimation
apollo_beta = c(asc_1 = 0,
asc_2 = 0,
asc_3 = 0,
asc_4 = 0,
# mu_tc_race1 = 0,
# mu_tc_race2 = -1,
# mu_tc_race3 = -1,
# mu_tc_race4 = -1,
# mu_tc_race5 = -1,
sig_tc_income = 0,
v_tt1 = 2,
v_tt2 = 2,
v_tt3 = 2,
v_tt4 = 2)
### 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_1", "mu_tc_race1")
apollo_fixed = c("asc_1")
# ################################################################# #
#### DEFINE RANDOM COMPONENTS ####
# ################################################################# #
### Set parameters for generating draws
apollo_HB = list(
hbDist = c(asc_1 = "NR",
asc_2 = "NR",
asc_3 = "NR",
asc_4 = "NR",
# mu_tc_race1 = "NR",
# mu_tc_race2 = "NR",
# mu_tc_race3 = "NR",
# mu_tc_race4 = "NR",
# mu_tc_race5 = "NR",
sig_tc_income = "N",
v_tt1 = "LN+",
v_tt2 = "LN+",
v_tt3 = "LN+",
v_tt4 = "LN+"),
fixedA = c(sig_tc_income = 0),
gNCREP = 10000, # burn-in iterations
gNEREP = 50000, # post burn-in iterations
gINFOSKIP = 500)
# ################################################################# #
#### GROUP AND VALIDATE INPUTS ####
# ################################################################# #
apollo_inputs = apollo_validateInputs()
# ################################################################# #
#### DEFINE MODEL AND LIKELIHOOD FUNCTION ####
# ################################################################# #
apollo_probabilities=function(apollo_beta, apollo_inputs, functionality="estimate"){
### Function initialisation: do not change the following three commands
### Attach inputs and detach after function exit
apollo_attach(apollo_beta, apollo_inputs)
on.exit(apollo_detach(apollo_beta, apollo_inputs))
# Generate a random parameter for b_tc with heterogeneity in mean/variance by race/income.
# b_tc = -exp(mu_race1 * race1 + mu_race2 * race2 + mu_race3 * race3 + mu_race4 * race4 + mu_race5 * race5 + sig_tc_income * income) # Does not work
b_tc = -exp(sig_tc_income * income) # Works
### 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[["car"]] = asc_1 + b_tc*(v_tt1*travelTimeDrive + travelCostDrive)
V[["walk"]] = asc_2 + b_tc*(v_tt2*travelTimeWalk)
V[["bike"]] = asc_3 + b_tc*(v_tt3*travelTimeBike)
V[["transit"]] = asc_4 + b_tc*(v_tt4 * travelTimeTransit + travelCostTransit)
### Define settings for MNL model component
mnl_settings = list(
alternatives = c(car=1, walk=2, bike=3, transit=4),
avail = list(car=1, walk=1, bike=1, transit=1),
choiceVar = trptrans,
utilities = V
)
### Compute probabilities using MNL model
P[["model"]] = apollo_mnl(mnl_settings, functionality)
### Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}
# ################################################################# #
#### MODEL ESTIMATION ####
# ################################################################# #
### optionally load pre-estimated model from memory
# model = apollo_loadModel(apollo_control$modelName)
### Optional speedTest
#speedTest_settings=list(
# nDrawsTry = c(50, 75, 100),
# nCoresTry = 1:3,
# nRep = 10
#)
#apollo_speedTest(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs, speedTest_settings)
model = apollo_estimate(apollo_beta, apollo_fixed,apollo_probabilities, apollo_inputs)