Important: Read this before posting to this forum

  1. This forum is for questions related to the use of Apollo. We will answer some general choice modelling questions too, where appropriate, and time permitting. We cannot answer questions about how to estimate choice models with other software packages.
  2. There is a very detailed manual for Apollo available at http://www.ApolloChoiceModelling.com/manual.html. This contains detailed descriptions of the various Apollo functions, and numerous examples are available at http://www.ApolloChoiceModelling.com/examples.html. In addition, help files are available for all functions, using e.g. ?apollo_mnl
  3. Before asking a question on the forum, users are kindly requested to follow these steps:
    1. Check that the same issue has not already been addressed in the forum - there is a search tool.
    2. Ensure that the correct syntax has been used. For any function, detailed instructions are available directly in Apollo, e.g. by using ?apollo_mnl for apollo_mnl
    3. Check the frequently asked questions section on the Apollo website, which discusses some common issues/failures. Please see http://www.apollochoicemodelling.com/faq.html
    4. Make sure that R is using the latest official release of Apollo.
  4. If the above steps do not resolve the issue, then users should follow these steps when posting a question:
    1. provide full details on the issue, including the entire code and output, including any error messages
    2. posts will not immediately appear on the forum, but will be checked by a moderator first. This may take a day or two at busy times. There is no need to submit the post multiple times.

Unreasonable estimate for cost

Ask questions about the results reported after estimation. If the output includes errors, please include your model code if possible.
Post Reply
Nishihara
Posts: 9
Joined: 19 May 2023, 02:14

Unreasonable estimate for cost

Post by Nishihara »

Dear Sir,

I am running a mixed logit model and getting a very unreasonable estimate for the log-normalized cost: -1066.11898

Please see the code and results. I cannot figure out what is the issue here because I already run the similar codes with different variable names in different datasets (same project).

It would be highly appreciated if you could let me know the reason behind it and how to fix it.

Thank you very much in advance.

Hara

Code:

Code: Select all

library(readr)
library(support.CEs)
library(survival)
library(openxlsx)
library(tidyr)
library(tidymodels)
library(apollo)
library(chatgpt)
library(dplyr)


setwd("~")

### Clear memory
rm(list = ls())

### Load Apollo library
library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "Bangladesh CE",
  modelDescr      = "Mixed logit model on Bangladesh CE data",
  indivID         = "ID",  
  mixing          = TRUE,
  # workInLogs = TRUE,
  nCores          = 15,
  outputDirectory = "~"
)


des1 <- rotation.design(attribute.names = list(Labor_intensity = c ("L0_labor", "L1_labor", "L2_labor", "L3_labor"), Ecosystem_integration = c ("L0_ecosystem", "L1_ecosystem", "L2_ecosystem", "L3_ecosystem"), Government_Support = c ("L0_support", "L1_support", "L2_support", "L3_support"), Management_system = c ("L0_system", "L1_system", "L2_system", "L3_system"), Fingerling_cost = c ("500", "1000", "1500", "2000")), nalternatives = 2, nblocks = 1, row.renames = FALSE, randomize = TRUE, seed = 4592)



des2 <- rotation.design(attribute.names = list(Labor_intensity = c ("0", "1", "2", "3"), Ecosystem_integration = c ("0", "1", "2", "3"), Government_Support = c ("0", "1", "2", "3"), Management_system = c ("0", "1", "2", "3"), Fingerling_cost = c ("500", "1000", "1500", "2000")), nalternatives = 2, nblocks = 1, row.renames = FALSE, randomize = TRUE, seed =4592)

alt <- des2$alternatives  
alt1 <- alt$alt.1 
# Add "1" after variables' name
names(alt1)[names(alt1) == "Labor_intensity"] <- paste0("Labor_intensity1")
names(alt1)[names(alt1) == "Ecosystem_integration"] <- paste0("Ecosystem_integration1")
names(alt1)[names(alt1) == "Government_Support"] <- paste0("Government_Support1")
names(alt1)[names(alt1) == "Management_system"] <- paste0("Management_system1")
names(alt1)[names(alt1) == "Fingerling_cost"] <- paste0("Fingerling_cost1")

alt2 <- alt$alt.2
names(alt2)[names(alt2) == "Labor_intensity"] <- paste0("Labor_intensity2")
names(alt2)[names(alt2) == "Ecosystem_integration"] <- paste0("Ecosystem_integration2")
names(alt2)[names(alt2) == "Government_Support"] <- paste0("Government_Support2")
names(alt2)[names(alt2) == "Management_system"] <- paste0("Management_system2")
names(alt2)[names(alt2) == "Fingerling_cost"] <- paste0("Fingerling_cost2")

alt.all <- cbind(alt1,alt2)

alt.all <- alt.all[,-c(1,3,9:11)]


input1=read.csv("559.csv")

num <- length(unique(input1$ID))
alt.all.expand <- alt.all[rep(seq_len(nrow(alt.all)), times = num), ]

desmat1 <- make.design.matrix(choice.experiment.design = des1, optout = TRUE, categorical.attributes = c("Labor_intensity", "Ecosystem_integration", "Government_Support", "Management_system"), continuous.attributes = c("Fingerling_cost"), unlabeled = TRUE)

dataset <- make.dataset(respondent.dataset = input1, choice.indicators = c("q1", "q2", "q3", "q4", "q5", "q6", "q7", "q8", "q9", "q10", "q11", "q12", "q13", "q14", "q15", "q16"), design.matrix = desmat1)    
choice <- dataset[,c(1,3:5)]
choice <- subset(choice,RES=="TRUE")

database <- cbind(choice,alt.all.expand)
database$Fingerling_cost1 <- as.numeric(as.character(database$Fingerling_cost1))
database$Fingerling_cost2 <- as.numeric(as.character(database$Fingerling_cost2))


### Vector of parameters, including any that are kept fixed in estimation
apollo_beta = c(asc = 0,
                sigma_asc = 0,
                #mu_Labor_intensity_L0    = 0,
                #sigma_Labor_intensity_L0 = 0,
                mu_Labor_intensity_L1    = 0,
                sigma_Labor_intensity_L1 = 0,
                mu_Labor_intensity_L2    = 0,
                sigma_Labor_intensity_L2 = 0,
                mu_Labor_intensity_L3    = 0,
                sigma_Labor_intensity_L3 = 0,
                #mu_Ecosystem_integration_L0    = 0,
                #sigma_Ecosystem_integration_L0 = 0,
                mu_Ecosystem_integration_L1    = 0,
                sigma_Ecosystem_integration_L1 = 0, 
                mu_Ecosystem_integration_L2    = 0,
                sigma_Ecosystem_integration_L2 = 0, 
                mu_Ecosystem_integration_L3    = 0,
                sigma_Ecosystem_integration_L3 = 0, 
                #mu_Government_Support_L0    = 0,
                #sigma_Government_Support_L0 = 0,
                mu_Government_Support_L1    = 0,
                sigma_Government_Support_L1 = 0,
                mu_Government_Support_L2    = 0,
                sigma_Government_Support_L2 = 0,
                mu_Government_Support_L3    = 0,
                sigma_Government_Support_L3 = 0,
                #mu_Management_system_L0    = 0,
                #sigma_Management_system_L0 = 0,
                mu_Management_system_L1    = 0,
                sigma_Management_system_L1 = 0,
                mu_Management_system_L2    = 0,
                sigma_Management_system_L2 = 0,
                mu_Management_system_L3    = 0,
                sigma_Management_system_L3 = 0,
                mu_log_Fingerling_cost    = -10,
                sigma_log_Fingerling_cost = 1
)

### 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()

# ################################################################# #
#### DEFINE RANDOM COMPONENTS                                    ####
# ################################################################# #

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "mlhs",
  interNDraws    = 500,
  interUnifDraws = c(),
  interNormDraws = c("draws_Labor_intensity_L0","draws_Labor_intensity_L1","draws_Labor_intensity_L2","draws_Labor_intensity_L3","draws_Ecosystem_integration_L0","draws_Ecosystem_integration_L1","draws_Ecosystem_integration_L2","draws_Ecosystem_integration_L3","draws_Government_Support_L0","draws_Government_Support_L1","draws_Government_Support_L2","draws_Government_Support_L3","draws_Management_system_L0","draws_Management_system_L1","draws_Management_system_L2","draws_Management_system_L3","draws_Fingerling_cost","draws_asc"),
  intraDrawsType = "mlhs",
  intraNDraws    = 0,
  intraUnifDraws = c(),
  intraNormDraws = c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  # randcoeff[["b_Labor_intensity_L0"]] =mu_Labor_intensity_L0 + sigma_Labor_intensity_L0 * draws_Labor_intensity_L0 
  randcoeff[["b_Labor_intensity_L1"]] =mu_Labor_intensity_L1 + sigma_Labor_intensity_L1 * draws_Labor_intensity_L1 
  randcoeff[["b_Labor_intensity_L2"]] =mu_Labor_intensity_L2 + sigma_Labor_intensity_L2 * draws_Labor_intensity_L2 
  randcoeff[["b_Labor_intensity_L3"]] =mu_Labor_intensity_L3 + sigma_Labor_intensity_L3 * draws_Labor_intensity_L3 
  #  randcoeff[["b_Ecosystem_integration_L0"]] =  mu_Ecosystem_integration_L0 + sigma_Ecosystem_integration_L0 * draws_Ecosystem_integration_L0
  randcoeff[["b_Ecosystem_integration_L1"]] =  mu_Ecosystem_integration_L1 + sigma_Ecosystem_integration_L1 * draws_Ecosystem_integration_L1
  randcoeff[["b_Ecosystem_integration_L2"]] =  mu_Ecosystem_integration_L2 + sigma_Ecosystem_integration_L2 * draws_Ecosystem_integration_L2
  randcoeff[["b_Ecosystem_integration_L3"]] =  mu_Ecosystem_integration_L3 + sigma_Ecosystem_integration_L3 * draws_Ecosystem_integration_L3
  #  randcoeff[["b_Management_system_L0"]] =  mu_Management_system_L0 + sigma_Management_system_L0 * draws_Management_system_L0 
  randcoeff[["b_Management_system_L1"]] =  mu_Management_system_L1 + sigma_Management_system_L1 * draws_Management_system_L1 
  randcoeff[["b_Management_system_L2"]] =  mu_Management_system_L2 + sigma_Management_system_L2 * draws_Management_system_L2 
  randcoeff[["b_Management_system_L3"]] =  mu_Management_system_L3 + sigma_Management_system_L3 * draws_Management_system_L3 
  #  randcoeff[["b_Government_Support_L0"]] =  mu_Government_Support_L0 + sigma_Government_Support_L0 * draws_Government_Support_L0 
  randcoeff[["b_Government_Support_L1"]] =  mu_Government_Support_L1 + sigma_Government_Support_L1 * draws_Government_Support_L1 
  randcoeff[["b_Government_Support_L2"]] =  mu_Government_Support_L2 + sigma_Government_Support_L2 * draws_Government_Support_L2 
  randcoeff[["b_Government_Support_L3"]] =  mu_Government_Support_L3 + sigma_Government_Support_L3 * draws_Government_Support_L3 
  randcoeff[["b_Fingerling_cost"]] =  -exp(mu_log_Fingerling_cost + sigma_log_Fingerling_cost * draws_Fingerling_cost) 
  randcoeff[["ec"]] =  sigma_asc * draws_asc
  
  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 ange the following three commands
  ### AareaaFingerling_cost inputs and detaFingerling_cost 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_seareaings, order is irrelevant
  V = list()
  V[['NB']]  = asc + ec
  V[['alt1']]  = #b_Labor_intensity_L0 * (Labor_intensity1==0) + 
    b_Labor_intensity_L1 * (Labor_intensity1==1) + b_Labor_intensity_L2 * (Labor_intensity1==2) +b_Labor_intensity_L3 * (Labor_intensity1==3) + 
    #b_Ecosystem_integration_L0 * (Ecosystem_integration1==0)+
    b_Ecosystem_integration_L1 * (Ecosystem_integration1==1) +  b_Ecosystem_integration_L2 * (Ecosystem_integration1==2)+ b_Ecosystem_integration_L3 * (Ecosystem_integration1==3) + 
    #b_Management_system_L0 * (Management_system1==0) +
    b_Management_system_L1 * (Management_system1==1) + b_Management_system_L2 * (Management_system1==2)+b_Management_system_L3 * (Management_system1==3) + 
    #b_Government_Support_L0 * 
    (Government_Support1==0) +b_Government_Support_L1 * (Government_Support1==1) + b_Government_Support_L2 * (Government_Support1==2)+b_Government_Support_L3 * (Government_Support1==3) + 
    b_Fingerling_cost * Fingerling_cost1 
  V[['alt2']]  = #b_Labor_intensity_L0 * (Labor_intensity2==0) + 
    b_Labor_intensity_L1 * (Labor_intensity2==1) + b_Labor_intensity_L2 * (Labor_intensity2==2) +b_Labor_intensity_L3 * (Labor_intensity2==3) + 
    #b_Ecosystem_integration_L0 * (Ecosystem_integration2==0) + 
    b_Ecosystem_integration_L1 * (Ecosystem_integration2==1) +  b_Ecosystem_integration_L2 * (Ecosystem_integration2==2)+ b_Ecosystem_integration_L3 * (Ecosystem_integration2==3) + 
    #b_Management_system_L0 * (Management_system2==0) + 
    b_Management_system_L1 * (Management_system2==1) + b_Management_system_L2 * (Management_system2==2)+b_Management_system_L3 * (Management_system2==3) + 
    #b_Government_Support_L0 * (Government_Support2==0) +
    b_Government_Support_L1 * (Government_Support2==1) + b_Government_Support_L2 * (Government_Support2==2)+b_Government_Support_L3 * (Government_Support2==3) + 
    b_Fingerling_cost * Fingerling_cost2
  
  
  ### Define seareaings for MNL model component
  mnl_seareaings = list(
    alternatives  = c(alt1="1", alt2="2", NB="3"),
    avail         = list(alt1=1, alt2=1, NB=1),
    choiceVar     = ALT,
    utilities     = V
  )
  
  ### Compute probabilities using MNL model
  P[["model"]] = apollo_mnl(mnl_seareaings, 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, estimate_settings = list(maxIterations=500))#,bootstrapSE=30))

# ################################################################# #
#### MODEL OUTPUTS                                               ####
# ################################################################# #

# ----------------------------------------------------------------- #
#---- FORMAareaED OUTPUT (TO SCREEN)                               ----
# ----------------------------------------------------------------- #

#apollo_loadModel("Bangladesh CE")

apollo_modelOutput(model, list(printPVal=2))
Results:

Code: Select all

stimation method                           : bfgs
Model diagnosis                             : successful convergence
Optimisation diagnosis                      : Maximum found
     hessian properties                     : Negative definitive
     maximum eigenvalue                     : 0
Number of individuals                       : 252
Number of rows in database                  : 4032
Number of modelled outcomes                 : 4032

Number of cores used                        :  15 
Number of inter-individual draws            : 500 (mlhs)

LL(start)                                   : -4464.3
LL at equal shares, LL(0)                   : -4429.6
LL at observed shares, LL(C)                : -3790.01
LL(final)                                   : -3195.03
Rho-squared vs equal shares                  :  0.2787 
Adj.Rho-squared vs equal shares              :  0.2724 
Rho-squared vs observed shares               :  0.157 
Adj.Rho-squared vs observed shares           :  0.1501 
AIC                                         :  6446.06 
BIC                                         :  6622.52 

Estimated parameters                        : 28
Time taken (hh:mm:ss)                       :  00:10:25.61 
     pre-estimation                         :  00:02:22.44 
     estimation                             :  00:03:22.77 
          initial estimation                :  00:03:15.28 
          estimation after rescaling        :  00:00:7.5 
     post-estimation                        :  00:04:40.4 
Iterations                                  :  102  
     initial estimation                     :  101 
     estimation after rescaling             :  1 

Unconstrained optimisation.

Estimates:
                                  Estimate        s.e.   t.rat.(0)  p(2-sided)    Rob.s.e. Rob.t.rat.(0)  p(2-sided)
asc                               -3.75892     0.37678    -9.97637    0.000000     0.44821       -8.3865    0.000000
sigma_asc                         -3.66013     0.37002    -9.89167    0.000000     0.43789       -8.3587    0.000000
mu_Labor_intensity_L1             -0.02455     0.12378    -0.19836    0.842764     0.11631       -0.2111    0.832808
sigma_Labor_intensity_L1          -1.10012     0.12803    -8.59265    0.000000     0.13879       -7.9267   2.220e-15
mu_Labor_intensity_L2             -0.28015     0.12110    -2.31336    0.020703     0.13049       -2.1468    0.031807
sigma_Labor_intensity_L2          -0.85888     0.15869    -5.41230   6.222e-08     0.17201       -4.9931   5.941e-07
mu_Labor_intensity_L3             -0.27549     0.07694    -3.58033  3.4316e-04     0.06771       -4.0685   4.732e-05
sigma_Labor_intensity_L3           0.07842     0.12872     0.60924    0.542365     0.07750        1.0119    0.311601
mu_Ecosystem_integration_L1        0.44106     0.12424     3.55012  3.8506e-04     0.13163        3.3507  8.0610e-04
sigma_Ecosystem_integration_L1     0.97570     0.10982     8.88421    0.000000     0.12245        7.9680   1.554e-15
mu_Ecosystem_integration_L2        0.11479     0.09728     1.17992    0.238031     0.11329        1.0132    0.310965
sigma_Ecosystem_integration_L2    -0.07574     0.13440    -0.56357    0.573050     0.09016       -0.8401    0.400842
mu_Ecosystem_integration_L3        0.12817     0.09487     1.35098    0.176702     0.10943        1.1713    0.241488
sigma_Ecosystem_integration_L3     0.51177     0.12595     4.06317   4.841e-05     0.12536        4.0825   4.455e-05
mu_Government_Support_L1           0.48127     0.09919     4.85209   1.222e-06     0.09484        5.0744   3.887e-07
sigma_Government_Support_L1       -0.16938     0.17465    -0.96985    0.332122     0.12701       -1.3336    0.182329
mu_Government_Support_L2           0.68296     0.10267     6.65202   2.891e-11     0.09940        6.8707   6.391e-12
sigma_Government_Support_L2       -0.87376     0.11951    -7.31151   2.642e-13     0.11870       -7.3613   1.823e-13
mu_Government_Support_L3           0.23872     0.10403     2.29486    0.021741     0.10881        2.1940    0.028238
sigma_Government_Support_L3        0.62242     0.14810     4.20272   2.637e-05     0.15362        4.0517   5.085e-05
mu_Management_system_L1           -0.37760     0.12265    -3.07866    0.002079     0.12372       -3.0521    0.002273
sigma_Management_system_L1        -1.26939     0.12358   -10.27198    0.000000     0.11781      -10.7750    0.000000
mu_Management_system_L2           -0.30828     0.13379    -2.30426    0.021208     0.13928       -2.2135    0.026865
sigma_Management_system_L2        -1.16989     0.11492   -10.18017    0.000000     0.11515      -10.1593    0.000000
mu_Management_system_L3            0.09415     0.08764     1.07427    0.282702     0.08782        1.0720    0.283712
sigma_Management_system_L3        -0.07113     0.12247    -0.58076    0.561405     0.09719       -0.7319    0.464248
mu_log_Fingerling_cost         -1066.11898   1.484e+04    -0.07183    0.942736   205.05933       -5.1991   2.003e-07
sigma_log_Fingerling_cost        386.32047  5409.26766     0.07142    0.943065    74.73679        5.1691   2.352e-07
Top
stephanehess
Site Admin
Posts: 998
Joined: 24 Apr 2020, 16:29

Re: Unreasonable estimate for cost

Post by stephanehess »

Hi

this would imply that the median is very close to zero. Did you estimate an MNL model - what were the result like for that for this coefficient?

Stephane
--------------------------------
Stephane Hess
www.stephanehess.me.uk
Top
Nishihara
Posts: 9
Joined: 19 May 2023, 02:14

Re: Unreasonable estimate for cost

Post by Nishihara »

Dear Professor,

Thank you very much for your reply.

I tried to do the estimation through MNL but come out this error when I run the apollo_estimate:

"Error in apollo_probabilities(apollo_beta, apollo_inputs, functionality = "validate") :
object '"validate"' not found".

I attached my code below. Thank you so much for your kind help.

Code: Select all

### Clear memory
rm(list = ls())

### Load Apollo library
library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "Bangladesh CE",
  modelDescr      = "Mixed logit model on Bangladesh CE data",
  indivID         = "ID",  
  nCores          = 15,
  outputDirectory = "~"
)


des1 <- rotation.design(attribute.names = list(Labor_intensity = c ("L0_labor", "L1_labor", "L2_labor", "L3_labor"), Ecosystem_integration = c ("L0_ecosystem", "L1_ecosystem", "L2_ecosystem", "L3_ecosystem"), Government_Support = c ("L0_support", "L1_support", "L2_support", "L3_support"), Management_system = c ("L0_system", "L1_system", "L2_system", "L3_system"), Fingerling_cost = c ("500", "1000", "1500", "2000")), nalternatives = 2, nblocks = 1, row.renames = FALSE, randomize = TRUE, seed = 4592)

des2 <- rotation.design(attribute.names = list(Labor_intensity = c ("0", "1", "2", "3"), Ecosystem_integration = c ("0", "1", "2", "3"), Government_Support = c ("0", "1", "2", "3"), Management_system = c ("0", "1", "2", "3"), Fingerling_cost = c ("500", "1000", "1500", "2000")), nalternatives = 2, nblocks = 1, row.renames = FALSE, randomize = TRUE, seed =4592)

alt <- des2$alternatives  
alt1 <- alt$alt.1 
# Add "1" after variables' name
names(alt1)[names(alt1) == "Labor_intensity"] <- paste0("Labor_intensity1")
names(alt1)[names(alt1) == "Ecosystem_integration"] <- paste0("Ecosystem_integration1")
names(alt1)[names(alt1) == "Government_Support"] <- paste0("Government_Support1")
names(alt1)[names(alt1) == "Management_system"] <- paste0("Management_system1")
names(alt1)[names(alt1) == "Fingerling_cost"] <- paste0("Fingerling_cost1")

alt2 <- alt$alt.2
names(alt2)[names(alt2) == "Labor_intensity"] <- paste0("Labor_intensity2")
names(alt2)[names(alt2) == "Ecosystem_integration"] <- paste0("Ecosystem_integration2")
names(alt2)[names(alt2) == "Government_Support"] <- paste0("Government_Support2")
names(alt2)[names(alt2) == "Management_system"] <- paste0("Management_system2")
names(alt2)[names(alt2) == "Fingerling_cost"] <- paste0("Fingerling_cost2")

alt.all <- cbind(alt1,alt2)

alt.all <- alt.all[,-c(1,3,9:11)]

input1=read.csv("559.csv")

num <- length(unique(input1$ID))
alt.all.expand <- alt.all[rep(seq_len(nrow(alt.all)), times = num), ]

desmat1 <- make.design.matrix(choice.experiment.design = des1, optout = TRUE, categorical.attributes = c("Labor_intensity", "Ecosystem_integration", "Government_Support", "Management_system"), continuous.attributes = c("Fingerling_cost"), unlabeled = TRUE)

dataset <- make.dataset(respondent.dataset = input1, choice.indicators = c("q1", "q2", "q3", "q4", "q5", "q6", "q7", "q8", "q9", "q10", "q11", "q12", "q13", "q14", "q15", "q16"), design.matrix = desmat1)    
choice <- dataset[,c(1,3:5)]
choice <- subset(choice,RES=="TRUE")

database <- cbind(choice,alt.all.expand)
database$Fingerling_cost1 <- as.numeric(as.character(database$Fingerling_cost1))
database$Fingerling_cost2 <- as.numeric(as.character(database$Fingerling_cost2))


### Vector of parameters, including any that are kept fixed in estimation
apollo_beta = c(
                asc = 0,
                b_Labor_intensity_L1    = 0,
                b_Labor_intensity_L2    = 0,
                b_Labor_intensity_L3    = 0,
                b_Ecosystem_integration_L1    = 0,
                b_Ecosystem_integration_L2    = 0,
                b_Ecosystem_integration_L3    = 0,
                b_Government_Support_L1    = 0,
                b_Government_Support_L2    = 0,
                b_Government_Support_L3    = 0,
                b_Management_system_L1    = 0,
                b_Management_system_L2    = 0,
                b_Management_system_L3    = 0,
                b_Fingerling_cost    = 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()


# ################################################################# #
#### 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
  ### AareaaFingerling_cost inputs and detaFingerling_cost 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 bst use the same names as in mnl_seareaings, order is irrelevant
  V = list()
  V[['NB']]  = asc 
  V[['alt1']]  =  
    b_Labor_intensity_L1 * (Labor_intensity1==1) + b_Labor_intensity_L2 * (Labor_intensity1==2) +b_Labor_intensity_L3 * (Labor_intensity1==3) + 
    b_Ecosystem_integration_L1 * (Ecosystem_integration1==1) +  b_Ecosystem_integration_L2 * (Ecosystem_integration1==2)+ b_Ecosystem_integration_L3 * (Ecosystem_integration1==3) + 
    b_Management_system_L1 * (Management_system1==1) + b_Management_system_L2 * (Management_system1==2)+b_Management_system_L3 * (Management_system1==3) + 
    (Government_Support1==0) +b_Government_Support_L1 * (Government_Support1==1) + b_Government_Support_L2 * (Government_Support1==2)+b_Government_Support_L3 * (Government_Support1==3) + 
    b_Fingerling_cost * Fingerling_cost1 
  V[['alt2']]  = 
    b_Labor_intensity_L1 * (Labor_intensity2==1) + b_Labor_intensity_L2 * (Labor_intensity2==2) +b_Labor_intensity_L3 * (Labor_intensity2==3) + 
    b_Ecosystem_integration_L1 * (Ecosystem_integration2==1) +  b_Ecosystem_integration_L2 * (Ecosystem_integration2==2)+ b_Ecosystem_integration_L3 * (Ecosystem_integration2==3) + 
    b_Management_system_L1 * (Management_system2==1) + b_Management_system_L2 * (Management_system2==2)+b_Management_system_L3 * (Management_system2==3) + 
    b_Government_Support_L1 * (Government_Support2==1) + b_Government_Support_L2 * (Government_Support2==2)+b_Government_Support_L3 * (Government_Support2==3) + 
    b_Fingerling_cost * Fingerling_cost2
  
  
  ### Define seareaings for MNL model component
  mnl_seareaings = list(
    alternatives  = c(alt1=1, alt2=2, NB=3),
    avail         = list(alt1=1, alt2=1, NB=1),
    choiceVar     = ALT,
    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 ESTIMATION                                            ####
# ################################################################# #

model = apollo_estimate(apollo_beta, apollo_fixed,apollo_probabilities, apollo_inputs)#
Top
stephanehess
Site Admin
Posts: 998
Joined: 24 Apr 2020, 16:29

Re: Unreasonable estimate for cost

Post by stephanehess »

do you get any error messages before that? Also, what version of Apollo are you using?
--------------------------------
Stephane Hess
www.stephanehess.me.uk
Top
Nishihara
Posts: 9
Joined: 19 May 2023, 02:14

Re: Unreasonable estimate for cost

Post by Nishihara »

Dear Professor,

Thank you very much for your help.

No, I did not get any error before apollo_estimate, and I am using the latest version of 0.2.9.
Top
stephanehess
Site Admin
Posts: 998
Joined: 24 Apr 2020, 16:29

Re: Unreasonable estimate for cost

Post by stephanehess »

Hi

are you by any chance running R4.3.0 or another very recent version. The R team made a change that led to this rather unhelpful error message. If you upgrade to R4.3.1, you will get the right messages about your error

Stephane
--------------------------------
Stephane Hess
www.stephanehess.me.uk
Top
Nishihara
Posts: 9
Joined: 19 May 2023, 02:14

Re: Unreasonable estimate for cost

Post by Nishihara »

Dear Professor,

Thank you very much for your help.

It did work after I update to R4.3.1 and no error came out.

However, through MNL, the estimate for the cost is positive but the value is reasonable. I attached the code and results below.

Do you think there is a connection between the unreasonable estimate for cost from the Mixed logit model and the result from MNL?

Thank you very much in advance.

Code: Select all

### Clear memory
rm(list = ls())

### Load Apollo library
library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "Bangladesh CE",
  modelDescr      = "MNL on Bangladesh CE data",
  indivID         = "ID",  
  nCores          = 15,
  outputDirectory = "~"
)


des1 <- rotation.design(attribute.names = list(Labor_intensity = c ("L0_labor", "L1_labor", "L2_labor", "L3_labor"), Ecosystem_integration = c ("L0_ecosystem", "L1_ecosystem", "L2_ecosystem", "L3_ecosystem"), Government_Support = c ("L0_support", "L1_support", "L2_support", "L3_support"), Management_system = c ("L0_system", "L1_system", "L2_system", "L3_system"), Fingerling_cost = c ("500", "1000", "1500", "2000")), nalternatives = 2, nblocks = 1, row.renames = FALSE, randomize = TRUE, seed = 4592)



des2 <- rotation.design(attribute.names = list(Labor_intensity = c ("0", "1", "2", "3"), Ecosystem_integration = c ("0", "1", "2", "3"), Government_Support = c ("0", "1", "2", "3"), Management_system = c ("0", "1", "2", "3"), Fingerling_cost = c ("500", "1000", "1500", "2000")), nalternatives = 2, nblocks = 1, row.renames = FALSE, randomize = TRUE, seed =4592)

alt <- des2$alternatives  
alt1 <- alt$alt.1 
# Add "1" after variables' name
names(alt1)[names(alt1) == "Labor_intensity"] <- paste0("Labor_intensity1")
names(alt1)[names(alt1) == "Ecosystem_integration"] <- paste0("Ecosystem_integration1")
names(alt1)[names(alt1) == "Government_Support"] <- paste0("Government_Support1")
names(alt1)[names(alt1) == "Management_system"] <- paste0("Management_system1")
names(alt1)[names(alt1) == "Fingerling_cost"] <- paste0("Fingerling_cost1")

alt2 <- alt$alt.2
names(alt2)[names(alt2) == "Labor_intensity"] <- paste0("Labor_intensity2")
names(alt2)[names(alt2) == "Ecosystem_integration"] <- paste0("Ecosystem_integration2")
names(alt2)[names(alt2) == "Government_Support"] <- paste0("Government_Support2")
names(alt2)[names(alt2) == "Management_system"] <- paste0("Management_system2")
names(alt2)[names(alt2) == "Fingerling_cost"] <- paste0("Fingerling_cost2")

alt.all <- cbind(alt1,alt2)

alt.all <- alt.all[,-c(1,3,9:11)]


input1=read.csv("559.csv")

num <- length(unique(input1$ID))
alt.all.expand <- alt.all[rep(seq_len(nrow(alt.all)), times = num), ]

desmat1 <- make.design.matrix(choice.experiment.design = des1, optout = TRUE, categorical.attributes = c("Labor_intensity", "Ecosystem_integration", "Government_Support", "Management_system"), continuous.attributes = c("Fingerling_cost"), unlabeled = TRUE)

dataset <- make.dataset(respondent.dataset = input1, choice.indicators = c("q1", "q2", "q3", "q4", "q5", "q6", "q7", "q8", "q9", "q10", "q11", "q12", "q13", "q14", "q15", "q16"), design.matrix = desmat1)    
choice <- dataset[,c(1,3:5)]
choice <- subset(choice,RES=="TRUE")

database <- cbind(choice,alt.all.expand)
database$Fingerling_cost1 <- as.numeric(as.character(database$Fingerling_cost1))
database$Fingerling_cost2 <- as.numeric(as.character(database$Fingerling_cost2))


### Vector of parameters, including any that are kept fixed in estimation
apollo_beta = c(#asc_1 = 0,
                #asc_2 = 0,
                asc = 0,
                b_Labor_intensity_L1    = 0,
                b_Labor_intensity_L2    = 0,
                b_Labor_intensity_L3    = 0,
                b_Ecosystem_integration_L1    = 0,
                b_Ecosystem_integration_L2    = 0,
                b_Ecosystem_integration_L3    = 0,
                b_Government_Support_L1    = 0,
                b_Government_Support_L2    = 0,
                b_Government_Support_L3    = 0,
                b_Management_system_L1    = 0,
                b_Management_system_L2    = 0,
                b_Management_system_L3    = 0,
                b_Fingerling_cost    = 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()


# ################################################################# #
#### 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 ange the following three commands
  ### AareaaFingerling_cost inputs and detaFingerling_cost 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 bst use the same names as in mnl_seareaings, order is irrelevant
  V = list()
  V[['NB']]  = asc 
  V[['alt1']]  =  
    b_Labor_intensity_L1 * (Labor_intensity1==1) + b_Labor_intensity_L2 * (Labor_intensity1==2) +b_Labor_intensity_L3 * (Labor_intensity1==3) + 
    b_Ecosystem_integration_L1 * (Ecosystem_integration1==1) +  b_Ecosystem_integration_L2 * (Ecosystem_integration1==2)+ b_Ecosystem_integration_L3 * (Ecosystem_integration1==3) + 
    b_Management_system_L1 * (Management_system1==1) + b_Management_system_L2 * (Management_system1==2)+b_Management_system_L3 * (Management_system1==3) + 
    (Government_Support1==0) +b_Government_Support_L1 * (Government_Support1==1) + b_Government_Support_L2 * (Government_Support1==2)+b_Government_Support_L3 * (Government_Support1==3) + 
    b_Fingerling_cost * Fingerling_cost1 
  V[['alt2']]  = 
    b_Labor_intensity_L1 * (Labor_intensity2==1) + b_Labor_intensity_L2 * (Labor_intensity2==2) +b_Labor_intensity_L3 * (Labor_intensity2==3) + 
    b_Ecosystem_integration_L1 * (Ecosystem_integration2==1) +  b_Ecosystem_integration_L2 * (Ecosystem_integration2==2)+ b_Ecosystem_integration_L3 * (Ecosystem_integration2==3) + 
    b_Management_system_L1 * (Management_system2==1) + b_Management_system_L2 * (Management_system2==2)+b_Management_system_L3 * (Management_system2==3) + 
    b_Government_Support_L1 * (Government_Support2==1) + b_Government_Support_L2 * (Government_Support2==2)+b_Government_Support_L3 * (Government_Support2==3) + 
    b_Fingerling_cost * Fingerling_cost2
  
  
  ### Define seareaings for MNL model component
  mnl_settings = list(
    alternatives  = c(alt1=1, alt2=2, NB=3),
    avail         = list(alt1=1, alt2=1, NB=1),
    choiceVar     = ALT,
    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 ESTIMATION                                            ####
# ################################################################# #


model = apollo_estimate(apollo_beta, apollo_fixed,apollo_probabilities, apollo_inputs)#, estimate_settings = list(maxIterations=500))#,bootstrapSE=30))

# ################################################################# #
#### MODEL OUTPUTS                                               ####
# ################################################################# #

# ----------------------------------------------------------------- #
#---- FORMAareaED OUTPUT (TO SCREEN)                               ----
# ----------------------------------------------------------------- #

#apollo_loadModel("Bangladesh CE")

apollo_modelOutput(model, list(printPVal=2))

Results:

Code: Select all

Model run by Jie using Apollo 0.2.9 on R 4.3.1 for Darwin.
www.ApolloChoiceModelling.com

Model name                                  : Bangladesh CE
Model description                           : MNL on Bangladesh CE data
Model run at                                : 2023-06-18 15:30:06.251488
Estimation method                           : bfgs
Model diagnosis                             : successful convergence
Optimisation diagnosis                      : Maximum found
     hessian properties                     : Negative definitive
     maximum eigenvalue                     : -33.620427
Number of individuals                       : 252
Number of rows in database                  : 4032
Number of modelled outcomes                 : 4032

Number of cores used                        :  15 
Model without mixing

LL(start)                                   : -4399.06
LL at equal shares, LL(0)                   : -4429.6
LL at observed shares, LL(C)                : -3790.01
LL(final)                                   : -3746.09
Rho-squared vs equal shares                  :  0.1543 
Adj.Rho-squared vs equal shares              :  0.1511 
Rho-squared vs observed shares               :  0.0116 
Adj.Rho-squared vs observed shares           :  0.0084 
AIC                                         :  7520.18 
BIC                                         :  7608.4 

Estimated parameters                        : 14
Time taken (hh:mm:ss)                       :  00:00:33.63 
     pre-estimation                         :  00:00:30.38 
     estimation                             :  00:00:1.57 
          initial estimation                :  00:00:1.48 
          estimation after rescaling        :  00:00:0.1 
     post-estimation                        :  00:00:1.67 
Iterations                                  :  24  
     initial estimation                     :  23 
     estimation after rescaling             :  1 

Unconstrained optimisation.

Estimates:
                              Estimate        s.e.   t.rat.(0)  p(2-sided)    Rob.s.e. Rob.t.rat.(0)  p(2-sided)
asc                           -0.87549     0.13369   -6.548407   5.815e-11     0.17704     -4.945173   7.608e-07
b_Labor_intensity_L1           0.02344     0.08044    0.291366    0.770772     0.07137      0.328369    0.742632
b_Labor_intensity_L2          -0.24491     0.08528   -2.871874    0.004080     0.07777     -3.149327    0.001636
b_Labor_intensity_L3          -0.19265     0.06766   -2.847238    0.004410     0.05646     -3.411967  6.4496e-04
b_Ecosystem_integration_L1     0.32801     0.08401    3.904426   9.445e-05     0.08495      3.861238  1.1281e-04
b_Ecosystem_integration_L2  6.9484e-04     0.08003    0.008682    0.993073     0.08431      0.008242    0.993424
b_Ecosystem_integration_L3     0.05349     0.07349    0.727858    0.466700     0.07971      0.671022    0.502206
b_Government_Support_L1        0.55402     0.07804    7.099463   1.252e-12     0.06418      8.631818    0.000000
b_Government_Support_L2        0.73521     0.07047   10.433013    0.000000     0.06987     10.522770    0.000000
b_Government_Support_L3        0.35307     0.07868    4.487211   7.216e-06     0.07989      4.419425   9.896e-06
b_Management_system_L1        -0.24548     0.07947   -3.089073    0.002008     0.08585     -2.859529    0.004243
b_Management_system_L2        -0.23672     0.08649   -2.736802    0.006204     0.08010     -2.955154    0.003125
b_Management_system_L3         0.06865     0.07039    0.975334    0.329395     0.05913      1.160980    0.245650
b_Fingerling_cost           2.1246e-04   5.348e-05    3.972637   7.108e-05   4.337e-05      4.898439   9.660e-07
Top
stephanehess
Site Admin
Posts: 998
Joined: 24 Apr 2020, 16:29

Re: Unreasonable estimate for cost

Post by stephanehess »

Hi

you say the estimate is positive but the value is reasonable. How can a positive cost coefficient be reasoable?

If the estimate is positive in MNL, then there is clearly a problem, and if you then use a negative lognormal, it will collapse to close to zero. This is why should always start with MNL and try to understand what is happening with your data.

Also, you have 3 alternatives and should estimat 2 constants, not just 1 - there may be order effects for example

Stephane
--------------------------------
Stephane Hess
www.stephanehess.me.uk
Top
Nishihara
Posts: 9
Joined: 19 May 2023, 02:14

Re: Unreasonable estimate for cost

Post by Nishihara »

Dear Professor,

Thank you very much for your prompt reply.

Now I understand the problem.

May I ask if I have results with positive values of cost estimate from MNL, what should I do with the analysis?

Does it mean that the data is problematic and the analysis is undoable?

Thank you very much for your help.
Top
stephanehess
Site Admin
Posts: 998
Joined: 24 Apr 2020, 16:29

Re: Unreasonable estimate for cost

Post by stephanehess »

Difficult to say as it really depends on your data. I suggest you study your data before you do any modelling and understand whether people really choose the options more when they are more expensive as that's what your results seem to tell you
--------------------------------
Stephane Hess
www.stephanehess.me.uk
Top
Post Reply

Return to “Estimation results”