First, my apologies for asking yet another question. Thanks for all your help thus far with my learning of DCE estimation.
I'm having difficulties running a code for mixed logit. I would like to make the preferences for only one attribute random. When I start with the basic form with only attributes, I get NA for the error terms. I fixed this by setting workInLogs=TRUE in apollo_control. That worked and I was able to get estimated results. Now, I would like to add interaction terms. Specifically, I am thinking of adding interactions of four respondent characteristics with all attributes (and levels for categorical attributes). However, it seems that once I add above a certain number of interactions, it no longer computes.
My code is:
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 = "MMNL_preference_space",
modelDescr = "Mixed logit model, uncorrelated Lognormal of system cost in preference space",
indivID = "id",
mixing = TRUE,
nCores = 4,
outputDirectory = "output",
workInLogs=TRUE
)
# ################################################################# #
#### 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:
library(readxl)
database = read_excel("")
### for data dictionary, use ?apollo_modeChoiceData
# ################################################################# #
#### DEFINE MODEL PARAMETERS ####
# ################################################################# #
### Vector of parameters, including any that are kept fixed in estimation
apollo_beta=c(asc_alt1 = 1.352605913,
asc_alt2 = 1.393373359,
asc_alt3 = 1.183729019,
asc_out = 0,
b_tenant = 0,
b_propertyinvestor = -0.19286773,
b_thirdparty = -0.700911191,
b_upfrontyes = 0,
b_upfrontno = 0.049493571,
b_upfrontlinkedprop = 0.100254433,
mu_log_b_systemcost = -3,
sigma_log_b_systemcost = -0.01,
b_tenant_highhhincome = 0,
b_propertyinvestor_highhhincome = 0,
b_thirdparty_highhhincome = 0,
b_upfrontyes_highhhincome = 0,
b_upfrontno_highhhincome = 0,
b_upfrontlinkedprop_highhhincome = 0,
b_systemcost_highhhincome = 0,
b_tenant_rentpriceagree = 0,
b_propertyinvestor_rentpriceagree = 0,
b_thirdparty_rentpriceagree = 0,
b_upfrontyes_rentpriceagree = 0,
b_upfrontno_rentpriceagree = 0,
b_upfrontlinkedprop_rentpriceagree = 0,
b_systemcost_rentpriceagree = 0,
b_tenant_rentimp = 0,
b_propertyinvestor_rentimp = 0,
b_thirdparty_rentimp = 0,
b_upfrontyes_rentimp = 0,
b_upfrontno_rentimp = 0,
b_upfrontlinkedprop_rentimp = 0,
b_systemcost_rentimp = 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_out", "b_tenant", "b_upfrontyes", "b_tenant_highhhincome", "b_upfrontyes_highhhincome", "b_tenant_rentpriceagree", "b_upfrontyes_rentpriceagree", "b_tenant_rentimp", "b_upfrontyes_rentimp")
# ################################################################# #
#### DEFINE RANDOM COMPONENTS ####
# ################################################################# #
### Set parameters for generating draws
apollo_draws = list(
interDrawsType = "halton",
interNDraws = 500,
interUnifDraws = c(),
interNormDraws = c("draws_systemcost"),
intraDrawsType = "halton",
intraNDraws = 0,
intraUnifDraws = c(),
intraNormDraws = c()
)
### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
randcoeff = list()
randcoeff[["b_systemcost"]] = -exp( mu_log_b_systemcost + sigma_log_b_systemcost * draws_systemcost )
return(randcoeff)
}
# ################################################################# #
#### 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))
### Create list of probabilities P
P = list()
### Create coefficients using interactions with socio-demographics
b_tenants_highhhincome = b_tenant + b_tenant_highhhincome * highhhincome
b_propertyinvestors_highhhincome = b_propertyinvestor + b_propertyinvestor_highhhincome * highhhincome
b_thirdpartys_highhhincome = b_thirdparty + b_thirdparty_highhhincome * highhhincome
b_upfrontyess_highhhincome = b_upfrontyes + b_upfrontyes_highhhincome * highhhincome
b_upfrontnos_highhhincome = b_upfrontno + b_upfrontno_highhhincome * highhhincome
b_upfrontlinkedprops_highhhincome = b_upfrontlinkedprop + b_upfrontlinkedprop_highhhincome * highhhincome
b_systemcosts_highhhincome = b_systemcost + b_systemcost_highhhincome * highhhincome
b_tenants_rentpriceagree = b_tenant + b_tenant_rentpriceagree * rentpriceagree
b_propertyinvestors_rentpriceagree = b_propertyinvestor + b_propertyinvestor_rentpriceagree * rentpriceagree
b_thirdpartys_rentpriceagree = b_thirdparty + b_thirdparty_rentpriceagree * rentpriceagree
b_upfrontyess_rentpriceagree = b_upfrontyes + b_upfrontyes_rentpriceagree * rentpriceagree
b_upfrontnos_rentpriceagree = b_upfrontno + b_upfrontno_rentpriceagree * rentpriceagree
b_upfrontlinkedprops_rentpriceagree = b_upfrontlinkedprop + b_upfrontlinkedprop_rentpriceagree * rentpriceagree
b_systemcosts_rentpriceagree = b_systemcost + b_systemcost_rentpriceagree * rentpriceagree
b_tenants_rentimp = b_tenant + b_tenant_rentimp * rentimp
b_propertyinvestors_rentimp = b_propertyinvestor + b_propertyinvestor_rentimp * rentimp
b_thirdpartys_rentimp = b_thirdparty + b_thirdparty_rentimp * rentimp
b_upfrontyess_rentimp = b_upfrontyes + b_upfrontyes_rentimp * rentimp
b_upfrontnos_rentimp = b_upfrontno + b_upfrontno_rentimp * rentimp
b_upfrontlinkedprops_rentimp = b_upfrontlinkedprop + b_upfrontlinkedprop_rentimp * rentimp
b_systemcosts_rentimp = b_systemcost + b_systemcost_rentimp * rentimp
### List of utilities: these must use the same names as in mnl_settings, order is irrelevant
V = list()
V[["alt1"]] = asc_alt1 + b_tenants_highhhincome*(meb1 == 1) + b_propertyinvestors_highhhincome*(meb1 == 2) + b_thirdpartys_highhhincome*(meb1 == 3) + b_upfrontyess_highhhincome*(pu1 == 1) + b_upfrontnos_highhhincome*(pu1 == 2) + b_upfrontlinkedprops_highhhincome*(pu1 == 3) + b_systemcosts_highhhincome*sc1 + b_tenants_rentpriceagree*(meb1 == 1) + b_propertyinvestors_rentpriceagree*(meb1 == 2) + b_thirdpartys_rentpriceagree*(meb1 == 3) + b_upfrontyess_rentpriceagree*(pu1 == 1) + b_upfrontnos_rentpriceagree*(pu1 == 2) + b_upfrontlinkedprops_rentpriceagree*(pu1 == 3) + b_systemcosts_rentpriceagree*sc1 + b_tenants_rentimp*(meb1 == 1) + b_propertyinvestors_rentimp*(meb1 == 2) + b_thirdpartys_rentimp*(meb1 == 3) + b_upfrontyess_rentimp*(pu1 == 1) + b_upfrontnos_rentimp*(pu1 == 2) + b_upfrontlinkedprops_rentimp*(pu1 == 3) + b_systemcosts_rentimp*sc1
V[["alt2"]] = asc_alt2 + b_tenants_highhhincome*(meb2 == 1) + b_propertyinvestors_highhhincome*(meb2 == 2) + b_thirdpartys_highhhincome*(meb2 == 3) + b_upfrontyess_highhhincome*(pu2 == 1) + b_upfrontnos_highhhincome*(pu2 == 2) + b_upfrontlinkedprops_highhhincome*(pu2 == 3) + b_systemcosts_highhhincome*sc2 + b_tenants_rentpriceagree*(meb2 == 1) + b_propertyinvestors_rentpriceagree*(meb2 == 2) + b_thirdpartys_rentpriceagree*(meb2 == 3) + b_upfrontyess_rentpriceagree*(pu2 == 1) + b_upfrontnos_rentpriceagree*(pu2 == 2) + b_upfrontlinkedprops_rentpriceagree*(pu2 == 3) + b_systemcosts_rentpriceagree*sc2 + b_tenants_rentimp*(meb2 == 1) + b_propertyinvestors_rentimp*(meb2 == 2) + b_thirdpartys_rentimp*(meb2 == 3) + b_upfrontyess_rentimp*(pu2 == 1) + b_upfrontnos_rentimp*(pu2 == 2) + b_upfrontlinkedprops_rentimp*(pu2 == 3) + b_systemcosts_rentimp*sc2
V[["alt3"]] = asc_alt3 + b_tenants_highhhincome*(meb3 == 1) + b_propertyinvestors_highhhincome*(meb3 == 2) + b_thirdpartys_highhhincome*(meb3 == 3) + b_upfrontyess_highhhincome*(pu3 == 1) + b_upfrontnos_highhhincome*(pu3 == 2) + b_upfrontlinkedprops_highhhincome*(pu3 == 3) + b_systemcosts_highhhincome*sc3 + b_tenants_rentpriceagree*(meb3 == 1) + b_propertyinvestors_rentpriceagree*(meb3 == 2) + b_thirdpartys_rentpriceagree*(meb3 == 3) + b_upfrontyess_rentpriceagree*(pu3 == 1) + b_upfrontnos_rentpriceagree*(pu3 == 2) + b_upfrontlinkedprops_rentpriceagree*(pu3 == 3) + b_systemcosts_rentpriceagree*sc3 + b_tenants_rentimp*(meb3 == 1) + b_propertyinvestors_rentimp*(meb3 == 2) + b_thirdpartys_rentimp*(meb3 == 3) + b_upfrontyess_rentimp*(pu3 == 1) + b_upfrontnos_rentimp*(pu3 == 2) + b_upfrontlinkedprops_rentimp*(pu3 == 3) + b_systemcosts_rentimp*sc3
V[["out"]] = asc_out
### Define settings for MNL model component
mnl_settings = list(
alternatives = c(alt1=1, alt2=2, alt3=3, out=4),
avail = list(alt1=1, alt2=1, alt3=1, out=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) ----
# ----------------------------------------------------------------- #
modelOutput_settings = list(printPVal=2)
apollo_modelOutput(model,modelOutput_settings)
# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO FILE, using model name) ----
# ----------------------------------------------------------------- #
saveOutput_settings = list(printPVal=2)
apollo_saveOutput(model, saveOutput_settings)
Thanks and regards,
Mara