Estimation of MNL and MMNL in WTP space issues
Posted: 31 Jan 2024, 16:06
Hello, I am new to DCE estimation using Apollo.
I am trying to estimate data from a labelled DCE using both a MNL and a Mixed MNL in WTP space, but I am ecountering some issues. The DCE cards included 4 alternatives (differing for the labelling attribute – the country of origin (COO), taking 4 levels – and the price – including three levels), and 1 opt-out alternative.
I will start with the MNL in WTP space issue. I want to estimate the WTP for the different COOs, hence the WTP of the ASCs. The issue is that when I specify utilities scaling the ASCs by the price to get WTP values, I get the following WARNING: Iteration limit exceeded. No covariance matrix will be computed. You may wish to use the current estimates as starting values for a new estimation with a higher limit for iterations.
Here is the code for the MNL:
Regarding the MMNL I would like to allow random preferences for all the alternatives (hence allowing the ASCs to randomly vary across the population) and to allow random preferences for the price (however constraining it to a negative distribution). I think I might have misspecified the utilities, since the results are not as expected (when I tried with a different estimation program, asc_italy_mu and asc_eu_mu were positive, asc_extraeu_mu and asc_nobuy_mu were negative).
Thank you for any clarification you will be able to give me.
Giuditta
I am trying to estimate data from a labelled DCE using both a MNL and a Mixed MNL in WTP space, but I am ecountering some issues. The DCE cards included 4 alternatives (differing for the labelling attribute – the country of origin (COO), taking 4 levels – and the price – including three levels), and 1 opt-out alternative.
I will start with the MNL in WTP space issue. I want to estimate the WTP for the different COOs, hence the WTP of the ASCs. The issue is that when I specify utilities scaling the ASCs by the price to get WTP values, I get the following WARNING: Iteration limit exceeded. No covariance matrix will be computed. You may wish to use the current estimates as starting values for a new estimation with a higher limit for iterations.
Here is the code for the MNL:
Code: Select all
rm(list = ls())
apollo_initialise()
database <- read_excel("Dati_final.xlsx")
database <- subset(database, database$pasta==1)
apollo_initialise()
### Set core controls
apollo_control = list(
modelName ="MNL",
modelDescr ="MNL model COO in WTP space",
indivID ="Individual"
)
apollo_beta=c(asc_italy = 0,
asc_eu = 0,
asc_extraeu = 0,
asc_nolabel = 0,
asc_nobuy = 0,
b_price = 0,
cet_italy = 0,
cet_eu = 0,
cet_extraeu = 0,
cet_nolabel = 0,
cet_nobuy = 0,
look_italy = 0,
look_eu = 0,
look_extraeu = 0,
look_nolabel = 0,
look_nobuy = 0,
knowledge_italy = 0,
knowledge_eu = 0,
knowledge_extraeu = 0,
knowledge_nolabel = 0,
knowledge_nobuy = 0)
apollo_fixed = c("asc_nobuy", "cet_nobuy", "look_nobuy", "knowledge_nobuy")
apollo_inputs = apollo_validateInputs()
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[['italy']] = b_price*(asc_italy + cet_italy * CET_high + look_italy * origin_look_always +
knowledge_italy * knowledge + price_italy)
V[["eu"]] = b_price*(asc_eu + cet_eu * CET_high + look_eu * origin_look_always +
knowledge_eu * knowledge + price_eu)
V[["extraeu"]] = b_price*(asc_extraeu + cet_extraeu * CET_high + look_extraeu * origin_look_always +
knowledge_extraeu * knowledge + price_extraeu)
V[["nolabel"]] = b_price*(asc_nolabel + cet_nolabel * CET_high + look_nolabel * origin_look_always +
knowledge_nolabel * knowledge + price_nolabel)
V[["nobuy"]] = b_price*(asc_nobuy + cet_nobuy * CET_high + look_nobuy * origin_look_always +
knowledge_nobuy * knowledge + price_nobuy)
### Define settings for MNL model component
mnl_settings = list(
alternatives = c(nolabel=1, italy=2, eu=3, extraeu=4, nobuy=5),
choiceVar = answer,
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)
### Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}
model_mnl_wtp_pasta = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs)
apollo_modelOutput(model_mnl_wtp_pasta)
rm(list = ls())
apollo_initialise()
database <- read_excel("Dati_final.xlsx")
database <- subset(database, database$pasta==1)
apollo_initialise()
### Set core controls
apollo_control = list(
modelName ="MNL",
modelDescr ="MNL model COO in WTP space",
indivID ="Individual"
)
apollo_beta=c(asc_italy = 0,
asc_eu = 0,
asc_extraeu = 0,
asc_nolabel = 0,
asc_nobuy = 0,
b_price = 0,
cet_italy = 0,
cet_eu = 0,
cet_extraeu = 0,
cet_nolabel = 0,
cet_nobuy = 0,
look_italy = 0,
look_eu = 0,
look_extraeu = 0,
look_nolabel = 0,
look_nobuy = 0,
knowledge_italy = 0,
knowledge_eu = 0,
knowledge_extraeu = 0,
knowledge_nolabel = 0,
knowledge_nobuy = 0)
apollo_fixed = c("asc_nobuy", "cet_nobuy", "look_nobuy", "knowledge_nobuy")
apollo_inputs = apollo_validateInputs()
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[['italy']] = asc_italy + b_price*(cet_italy * CET_high + look_italy * origin_look_always +
knowledge_italy * knowledge + price_italy)
V[["eu"]] = asc_eu + b_price*(cet_eu * CET_high + look_eu * origin_look_always +
knowledge_eu * knowledge + price_eu)
V[["extraeu"]] = asc_extraeu + b_price*(cet_extraeu * CET_high + look_extraeu * origin_look_always +
knowledge_extraeu * knowledge + price_extraeu)
V[["nolabel"]] = asc_nolabel + b_price*(cet_nolabel * CET_high + look_nolabel * origin_look_always +
knowledge_nolabel * knowledge + price_nolabel)
V[["nobuy"]] = asc_nobuy + b_price*(cet_nobuy * CET_high + look_nobuy * origin_look_always +
knowledge_nobuy * knowledge + price_nobuy)
### Define settings for MNL model component
mnl_settings = list(
alternatives = c(nolabel=1, italy=2, eu=3, extraeu=4, nobuy=5),
choiceVar = answer,
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)
### Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}
model_mnl_wtp_pasta = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs)
apollo_modelOutput(model_mnl_wtp_pasta)
# The following are the estimated parameters
Estimated parameters:
Estimate
asc_italy 3.696e+04
asc_eu 3462.774
asc_extraeu -3.319e+04
asc_nolabel -2.640e+04
asc_nobuy 0.000
b_price 4.281e-05
cet_italy 2.184e+04
cet_eu 1.081e+04
cet_extraeu -8550.641
cet_nolabel 5200.902
cet_nobuy 0.000
look_italy -2.212e+04
look_eu -7736.336
look_extraeu 1.240e+04
look_nolabel -1.501e+04
look_nobuy 0.000
knowledge_italy 3.790e+04
knowledge_eu 4.253e+04
knowledge_extraeu 5.930e+04
knowledge_nolabel 6.054e+04
knowledge_nobuy 0.000
Regarding the MMNL I would like to allow random preferences for all the alternatives (hence allowing the ASCs to randomly vary across the population) and to allow random preferences for the price (however constraining it to a negative distribution). I think I might have misspecified the utilities, since the results are not as expected (when I tried with a different estimation program, asc_italy_mu and asc_eu_mu were positive, asc_extraeu_mu and asc_nobuy_mu were negative).
Code: Select all
################################################################################
> ################################## RPL Pooled ##################################
> ################################################################################
> ### Clear memory
> rm(list = ls())
> database <- read_excel("Dati_final.xlsx")
> apollo_initialise()
Apollo ignition sequence completed
> ### Set core controls
> apollo_control = list(
+ modelName = "RPL_pooled_WTP_space",
+ modelDescr = "RPL model pooled in WTP space",
+ indivID = "Individual",
+ mixing = TRUE, #true for models including random parameters
+ nCores = 3,
+ outputDirectory = "output"
+ )
> # ################################################################# #
> #### LOAD DATA AND APPLY ANY TRANSFORMATIONS ####
> # ################################################################# #
> database <- read_excel("Dati_final.xlsx")
> ### Vector of parameters, including any that are kept fixed in estimation
> apollo_beta = c(asc_italy_mu = 0,
+ asc_italy_sig = 0,
+ asc_eu_mu = 0,
+ asc_eu_sig = 0,
+ asc_extraeu_mu = 0,
+ asc_extraeu_sig = 0,
+ asc_nolabel_mu = 0,
+ asc_nolabel_sig = 0,
+ asc_nobuy_mu = 0,
+ asc_nobuy_sig = 0,
+ b_price_mu = 0,
+ b_price_sig = 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_nolabel_mu", "asc_nolabel_sig")
> ### Set parameters for generating draws
> apollo_draws = list(
+ interDrawsType = "sobol",
+ interNDraws = 1000,
+ interNormDraws = c("draws_asc_italy","draws_asc_eu", "draws_asc_extraeu",
+ "draws_asc_nolabel", "draws_asc_nobuy", "draws_b_price")
+ )
> ### Create random parameters
> apollo_randCoeff = function(apollo_beta, apollo_inputs){
+ randcoeff = list()
+
+ randcoeff[["asc_italy"]] = asc_italy_mu + asc_italy_sig * draws_asc_italy
+ randcoeff[["asc_eu"]] = asc_eu_mu + asc_eu_sig * draws_asc_eu
+ randcoeff[["asc_extraeu"]] = asc_extraeu_mu + asc_extraeu_sig * draws_asc_extraeu
+ randcoeff[["asc_nolabel"]] = asc_nolabel_mu + asc_nolabel_sig * draws_asc_nolabel
+ randcoeff[["asc_nobuy"]] = asc_nobuy_mu + asc_nobuy_sig * draws_asc_nobuy
+ randcoeff[["b_price"]] = b_price_mu + b_price_sig * draws_b_price
+
+ return(randcoeff)
+ }
> apollo_inputs = apollo_validateInputs()
apollo_draws and apollo_randCoeff were found, so apollo_control$mixing was set
to TRUE
Several observations per individual detected based on the value of Individual.
Setting panelData in apollo_control set to TRUE.
All checks on apollo_control completed.
All checks on database completed.
Generating inter-individual draws ...... Done
> 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()
+
+ ### List of utilities: these must use the same names as in mnl_settings, order is irrelevant
+ V = list()
+ V[["italy"]] = b_price * ( asc_italy + price_italy )
+ V[["eu"]] = b_price * ( asc_eu + price_eu )
+ V[["extraeu"]] = b_price * ( asc_extraeu + price_extraeu )
+ V[["nolabel"]] = b_price * ( asc_nolabel + price_nolabel )
+ V[["nobuy"]] = b_price * ( asc_nobuy + price_nobuy )
+
+ ### Define settings for MNL model component
+ mnl_settings = list(
+ alternatives = c(nolabel=1, italy=2, eu=3, extraeu=4, nobuy=5),
+ choiceVar = answer,
+ 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 = apollo_estimate(apollo_beta, apollo_fixed,
+ apollo_probabilities, apollo_inputs)
Preparing user-defined functions.
Testing likelihood function...
Setting "avail" is missing, so full availability is assumed.
Overview of choices for MNL model component :
nolabel italy eu extraeu nobuy
Times available 864.00 864.00 864.00 864.00 864.00
Times chosen 42.00 512.00 145.00 27.00 138.00
Percentage chosen overall 4.86 59.26 16.78 3.12 15.97
Percentage chosen when available 4.86 59.26 16.78 3.12 15.97
Pre-processing likelihood function...
Creating cluster...
Preparing workers for multithreading...
Testing influence of parameters
Starting main estimation
BGW using analytic model derivatives supplied by caller...
Iterates will be written to:
output/RPL_pooled_WTP_space_iterations.csv it nf F RELDF PRELDF RELDX MODEL stppar
0 1 1.390554356e+03
1 2 1.342786415e+03 3.435e-02 1.198e-03 1.00e+00 G -5.60e-04
2 3 1.312630929e+03 2.246e-02 2.976e-02 8.37e-01 G 3.15e-02
3 4 1.006640954e+03 2.331e-01 3.221e-02 1.00e+00 G 1.06e-04
4 5 9.368733968e+02 6.931e-02 4.176e-02 3.55e-01 G 0.00e+00
5 6 8.726781237e+02 6.852e-02 4.183e-02 2.53e-01 G 0.00e+00
6 7 8.432197868e+02 3.376e-02 9.594e-02 5.53e-01 S 0.00e+00
7 8 7.732461899e+02 8.298e-02 6.083e-02 2.84e-01 S 0.00e+00
8 9 7.355182646e+02 4.879e-02 3.943e-02 2.14e-01 S 0.00e+00
9 11 7.173275081e+02 2.473e-02 1.743e-02 5.33e-02 S 1.09e-01
10 13 6.998903603e+02 2.431e-02 1.328e-02 3.26e-01 S-G 0.00e+00
11 14 6.908877066e+02 1.286e-02 9.946e-03 6.60e-02 G 0.00e+00
12 15 6.827644060e+02 1.176e-02 6.844e-03 3.72e-02 G 0.00e+00
13 16 6.769681613e+02 8.489e-03 5.156e-03 3.11e-02 G 0.00e+00
14 17 6.695499781e+02 1.096e-02 1.044e-02 1.01e-01 S 0.00e+00
15 18 6.689020071e+02 9.678e-04 2.972e-03 2.45e-02 S 0.00e+00
16 19 6.672236127e+02 2.509e-03 1.702e-03 1.52e-02 S 0.00e+00
17 20 6.661568301e+02 1.599e-03 2.266e-03 6.31e-02 G 0.00e+00
18 21 6.654997320e+02 9.864e-04 9.121e-04 1.77e-02 S 0.00e+00
19 22 6.652625255e+02 3.564e-04 3.280e-04 1.10e-02 S 0.00e+00
20 23 6.650346093e+02 3.426e-04 2.063e-04 6.80e-03 S 0.00e+00
21 25 6.649653261e+02 1.042e-04 4.520e-04 1.34e-02 G-S 0.00e+00
22 26 6.648268274e+02 2.083e-04 3.053e-04 1.12e-02 S 0.00e+00
23 27 6.648218495e+02 7.487e-06 9.895e-05 4.11e-03 S 0.00e+00
24 28 6.647741132e+02 7.180e-05 7.302e-05 2.23e-03 S 0.00e+00
25 30 6.647695561e+02 6.855e-06 1.192e-05 1.14e-03 S 2.05e-01
26 31 6.647657737e+02 5.690e-06 4.384e-06 6.60e-04 S 0.00e+00
27 32 6.647645355e+02 1.863e-06 2.304e-06 8.67e-04 S 0.00e+00
28 33 6.647643813e+02 2.319e-07 8.480e-07 2.87e-04 S 0.00e+00
29 34 6.647641596e+02 3.336e-07 4.840e-07 1.57e-04 S 0.00e+00
30 35 6.647641130e+02 7.009e-08 1.026e-07 7.61e-05 S 0.00e+00
31 36 6.647641047e+02 1.236e-08 1.268e-08 3.35e-05 S 0.00e+00
32 37 6.647641046e+02 1.673e-10 1.715e-10 4.73e-06 S 0.00e+00
33 38 6.647641046e+02 1.631e-11 1.668e-11 1.33e-06 S 0.00e+00
***** Relative function convergence *****
Additional convergence test using scaled estimation. Parameters will be scaled
by their current estimates and additional iterations will be performed.
BGW using analytic model derivatives supplied by caller...
Iterates will be appended to:
output/RPL_pooled_WTP_space_iterations.csv it nf F RELDF PRELDF RELDX MODEL stppar
0 1 6.647641046e+02
***** Relative function convergence *****
Estimated parameters with approximate standard errors from BHHH matrix:
Estimate BHHH se BHH t-ratio (0)
asc_italy_mu -0.58641 0.03916 -14.9743
asc_italy_sig -0.24828 0.04138 -5.9998
asc_eu_mu -0.22368 0.03231 -6.9240
asc_eu_sig -0.14508 0.04683 -3.0981
asc_extraeu_mu 0.07850 0.07254 1.0821
asc_extraeu_sig 0.11142 0.11185 0.9961
asc_nolabel_mu 0.00000 NA NA
asc_nolabel_sig 0.00000 NA NA
asc_nobuy_mu 2.48685 0.19652 12.6544
asc_nobuy_sig 0.97919 0.17659 5.5449
b_price_mu -10.03383 1.25862 -7.9721
b_price_sig 4.86259 1.23623 3.9334
Final LL: -664.7641
Calculating log-likelihood at equal shares (LL(0)) for applicable models...
Calculating log-likelihood at observed shares from estimation data (LL(c)) for
applicable models...
Calculating LL of each model component...
Calculating other model fit measures
Computing covariance matrix using analytical gradient.
0%....25%....50%....75%.100%
Negative definite Hessian with maximum eigenvalue: -0.533522
Computing score matrix...
Your model was estimated using the BGW algorithm. Please acknowledge this by
citing Bunch et al. (1993) - DOI 10.1145/151271.151279
> # ----------------------------------------------------------------- #
> #---- FORMATTED OUTPUT (TO SCREEN) ----
> # ----------------------------------------------------------------- #
> apollo_modelOutput(model)
Model run by giudi using Apollo 0.3.1 on R 4.2.3 for Darwin.
Please acknowledge the use of Apollo by citing Hess & Palma (2019)
DOI 10.1016/j.jocm.2019.100170
www.ApolloChoiceModelling.com
Model name : RPL_pooled_WTP_space
Model description : RPL model pooled in WTP space
Model run at : 2024-01-31 17:03:01
Estimation method : bgw
Model diagnosis : Relative function convergence
Optimisation diagnosis : Maximum found
hessian properties : Negative definite
maximum eigenvalue : -0.533522
reciprocal of condition number : 0.000150828
Number of individuals : 96
Number of rows in database : 864
Number of modelled outcomes : 864
Number of cores used : 3
Number of inter-individual draws : 1000 (sobol)
LL(start) : -1390.55
LL at equal shares, LL(0) : -1390.55
LL at observed shares, LL(C) : -1000.42
LL(final) : -664.76
Rho-squared vs equal shares : 0.5219
Adj.Rho-squared vs equal shares : 0.5148
Rho-squared vs observed shares : 0.3355
Adj.Rho-squared vs observed shares : 0.3295
AIC : 1349.53
BIC : 1397.14
Estimated parameters : 10
Time taken (hh:mm:ss) : 00:00:47.05
pre-estimation : 00:00:9.82
estimation : 00:00:15.1
initial estimation : 00:00:14.06
estimation after rescaling : 00:00:1.04
post-estimation : 00:00:22.13
Iterations : 34
initial estimation : 33
estimation after rescaling : 1
Unconstrained optimisation.
Estimates:
Estimate s.e. t.rat.(0) Rob.s.e. Rob.t.rat.(0)
asc_italy_mu -0.58641 0.04942 -11.866 0.06929 -8.463
asc_italy_sig -0.24828 0.03087 -8.043 0.03397 -7.309
asc_eu_mu -0.22368 0.03295 -6.789 0.04105 -5.449
asc_eu_sig -0.14508 0.02845 -5.099 0.02185 -6.641
asc_extraeu_mu 0.07850 0.06137 1.279 0.06565 1.196
asc_extraeu_sig 0.11142 0.06575 1.695 0.05520 2.018
asc_nolabel_mu 0.00000 NA NA NA NA
asc_nolabel_sig 0.00000 NA NA NA NA
asc_nobuy_mu 2.48685 0.09584 25.948 0.05875 42.332
asc_nobuy_sig 0.97919 0.08122 12.056 0.04464 21.937
b_price_mu -10.03383 1.18306 -8.481 1.69716 -5.912
b_price_sig 4.86259 0.82722 5.878 0.93897 5.179
Giuditta