Error in the model command when using EM algorithms
Posted: 25 May 2024, 23:56
Dear Prof. Hess & Palma,
I am trying to estimate a Mixed Multinomial Logit model in WTP space using EM algorithms because I want all the parameters to be random. The model has 11 parameters. 10 parameters are normally distributed and the cost attribute has a negative log-normal distribution. I used the estimates I obtained from the preference space MNL model as starting values for means.
Here is my code:
With the model command, I am getting an error like this:
As suggested in another post, I tried installing "bgw" package and then running the model. But, it did not work. Any help would be greatly appreciated.
Many thanks,
Manori
I am trying to estimate a Mixed Multinomial Logit model in WTP space using EM algorithms because I want all the parameters to be random. The model has 11 parameters. 10 parameters are normally distributed and the cost attribute has a negative log-normal distribution. I used the estimates I obtained from the preference space MNL model as starting values for means.
Here is my code:
Code: Select all
# ################################################################# #
#### LOAD LIBRARY AND DEFINE CORE SETTINGS ####
# ################################################################# #
### Clear memory
rm(list = ls())
### Load Apollo library
library(apollo)
library(dplyr)
### Initialise code
apollo_initialise()
### Set core controls
apollo_control = list(
modelName ="RPL_WTPspace_EM",
modelDescr ="RPL model1_EM algorithm ",
indivID ="ID",
mixing = TRUE,
nCores = 10
)
# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS ####
# ################################################################# #
database <- read_csv("Project data.csv")
# ################################################################# #
#### DEFINE MODEL PARAMETERS ####
# ################################################################# #
### Vector of parameters, including any that are kept fixed in estimation
apollo_beta <- c(
asc_1 = -0.4,
asc_2 = 0,
mu_c = 0.01,
mu_m2 = 0.5,
mu_m3 = 0.6,
mu_t2 = 0.2,
mu_t3 = 0.4,
mu_s2 = 0.2,
mu_s3 = 0.2,
mu_s4 = 0.1,
mu_s5 = -0.2,
mu_s6 = -0.5,
mu_cost = -0.007,
sigma_c = 0.1,
sigma_m2 = 0.3,
sigma_m3 = 0.2,
sigma_t2 = 0.3,
sigma_t3 = 0.2,
sigma_s2 = 0.3,
sigma_s3 = 0.2,
sigma_s4 = 0.3,
sigma_s5 = 0.2,
sigma_s6 = 0.3,
sigma_cost = 2.5
)
### 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_2")
# ################################################################# #
#### DEFINE RANDOM COMPONENTS ####
# ################################################################# #
### Set parameters for generating draws
apollo_draws = list(
interDrawsType = "sobol",
interNDraws = 1000,
interNormDraws = c("xi_c", "xi_m2", "xi_m3", "xi_t2", "xi_t3",
"xi_s2", "xi_s3", "xi_s4", "xi_s5", "xi_s6", "xi_cost")
)
### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
randcoeff = list()
randcoeff[["b_c"]] = mu_c + sigma_c * xi_c
randcoeff[["b_m2"]] = mu_m2 + sigma_m2 * xi_m2
randcoeff[["b_m3"]] = mu_m3 + sigma_m3 * xi_m3
randcoeff[["b_t2"]] = mu_t2 + sigma_t2 * xi_t2
randcoeff[["b_t3"]] = mu_t3 + sigma_t3 * xi_t3
randcoeff[["b_s2"]] = mu_s2 + sigma_s2 * xi_s2
randcoeff[["b_s3"]] = mu_s3 + sigma_s3 * xi_s3
randcoeff[["b_s4"]] = mu_s4 + sigma_s4 * xi_s4
randcoeff[["b_s5"]] = mu_s5 + sigma_s5 * xi_s5
randcoeff[["b_s6"]] = mu_s6 + sigma_s6 * xi_s6
randcoeff[["b_cost"]]= -exp( mu_cost + sigma_cost * xi_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[['alt1']] = asc_1 + b_cost * (b_c*a1_c + b_m2*a1_m2 + b_m3*a1_m3 + b_t2*a1_t2 + b_t3*a1_t3
+ b_s2*a1_s2 + b_s3*a1_s3 + b_s4*a1_s4 + b_s5*a1_s5 + b_s6*a1_s6 + cost1)
V[['alt2']] = asc_2 + b_cost * (b_c*a2_c + b_m2*a2_m2 + b_m3*a2_m3 + b_t2*a2_t2 + b_t3*a2_t3
+ b_s2*a2_s2 + b_s3*a2_s3 + b_s4*a2_s4 + b_s5*a2_s5 + b_s6*a2_s6 + cost2)
### Define settings for MNL model component
mnl_settings <- list(
alternatives = c(alt1 = 1, alt2 = 2),
avail = list(alt1 = 1, alt2 = 1),
choiceVar = choice,
V = 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)
}
# ################################################################# #
#### EM ESTIMATION ####
# ################################################################# #
mixEM_settings = list(transforms=list(function(x) (x),
function(x) (x),
function(x) (x),
function(x) (x),
function(x) (x),
function(x) (x),
function(x) (x),
function(x) (x),
function(x) (x),
function(x) (x),
function(x) log(-x)))
model = apollo_mixEM(apollo_beta, apollo_fixed, apollo_probabilities,
apollo_inputs, mixEM_settings)
# ################################################################# #
#### MODEL OUTPUTS ####
# ################################################################# #
# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO SCREEN) ----
# ----------------------------------------------------------------- #
apollo_modelOutput(model)
# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO FILE, using model name) ----
# ----------------------------------------------------------------- #
apollo_saveOutput(model)
Code: Select all
model = apollo_mixEM(apollo_beta, apollo_fixed, apollo_probabilities,
+ apollo_inputs, mixEM_settings)
INFORMATION: The use of apollo_mixEM has a number of requirements. No checks are run for these, so the user needs
to ensure these conditions are met by their model:
1: This function is only suitable for single component models, i.e. no use of apollo_combineModels or
manual multiplication of model components.
2: All parameters need to be random, and a full covariance matrix needs to be estimated/specified in
apollo_randCoeff.
3: All random parameters need to be based on Normal distributions or transformations thereof.
4: With K random parameters, the order of the elements in 'apollo_beta' needs to be as follows: K
means for the underlying Normals, followed by the elements of the lower triangle of the Cholesky
matrix, by row.
Initialising EM algorithm
Validating inputs of likelihood function (apollo_probabilities)
Error in bgw::bgw_mle(calcR = f, betaStart = beta_var_val, calcJ = grad, :
User-provided value for bgw_settings[["maxFunctionEvals"]] is not valid.
Many thanks,
Manori