Dear Stephan,
I have tried to run interactions with socio-demographic variables in my MNL model, but so far I have not been successful. My socio-demographic variables consist of age, use of pesticides or not, gender variable so far. I have more, but considering that I have not been successful in this, I wanted to consult with you if there is any mistake in the coding.
I have also tried to interact these conditions mainly with my production system due to my research question. Below you will find the code:
################################################################## #
#### LOAD LIBRARY AND DEFINE CORE SETTINGS ####
# ################################################################# #
### Clear memory
rm(list = ls())
### Load Apollo library
library(apollo)
# Load your dataset
database <- read.csv("Ecuador_database.csv")
# Ensure Provincia is a binary variable (0 or 1)
database$Provincia <- as.integer(database$Provincia)
# Split the dataset into two based on the Provincia variable
database_Guayas <- subset(database, Provincia == 1)
database_Napo <- subset(database, Provincia == 0)
#############################################################################################
########################## GUAYAS ################################################
#############################################################################################
# Set core controls for Provincia == 1
apollo_control = list(
modelName = "Model_Provincia_1",
modelDescr = "MNL model for Provincia == 1",
indivID = "code",
nCores = 1,
outputDirectory = "output"
)
# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS ####
# ################################################################# #
database = database_Guayas
# Ensure Province is a binary variable (0 or 1)
database$typeproduction <- as.factor(database$typeproduction)
################################################### #
#### DEFINE MODEL PARAMETERS ####
# ################################################################# #
apollo_beta = c(asc1 = 0,
asc2 = 0,
b_natural = 0,
b_organic = 0,
b_chemical = 0,
b_conventional = 0,
b_wood_cacao = 0,
b_woodfruit_cacao = 0,
b_fruit_cacao = 0,
b_cacao = 0,
b_government = 0,
b_farmers = 0,
b_scientifics = 0,
b_norecomm = 0,
b_pubicinstitu = 0,
b_intermediaries = 0,
b_privateenter = 0,
b_agriorganiza = 0,
b_compensation = 0,
#with women and type of production
#b_naturalXwomen = 0,
#b_organicXwomen = 0,
#b_chemicalXwomen = 0,
#with the pest and production system
#b_pestYesXnatural = 0,
#b_pestYesXorganic = 0,
#b_pestYesXchemical = 0,
#age with the production system
b_naturalXage = 0,
b_organicXage = 0,
b_chemicalXage = 0,
b_out = 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("b_conventional", "b_cacao", "b_norecomm", "b_intermediaries", "b_out")
# Validate inputs for Provincia == 1
apollo_inputs_1 = apollo_validateInputs()
# Define model and likelihood function for Provincia == 1
apollo_probabilities_1 = function(apollo_beta, apollo_inputs, functionality="estimate") {
# Function initialisation: do not change the following three commands
apollo_attach(apollo_beta, apollo_inputs)
on.exit(apollo_detach(apollo_beta, apollo_inputs))
# Create list of probabilities P
P = list()
#with women and production system
#b1_naturalXwomen = b_natural * (P1 == 1) + b_naturalXwomen * women_d
#b1_organicXwomen = b_organic * (P1 == 2) + b_organicXwomen * women_d
#b1_chemicalXwomen = b_chemical * (P1 == 3) + b_chemicalXwomen * women_d
#b2_naturalXwomen = b_natural * (P2 == 1) + b_naturalXwomen * women_d
#b2_organicXwomen = b_organic * (P2 == 2) + b_organicXwomen * women_d
#b2_chemicalXwomen = b_chemical * (P2 == 3) + b_chemicalXwomen * women_d
#with the pest and the organic production system
#pest1YesXnatural= b_natural * (P1 == 1) + b_pestYesXnatural * (pest==1)
#pest1YesXorganic = b_organic * (P1 == 2) + b_pestYesXorganic * (pest==1)
#pest1YesXchemical = b_chemical * (P1 == 3) + b_pestYesXchemical * (pest==1)
#pest2YesXnatural= b_natural * (P2 == 1) + b_pestYesXnatural * (pest==1)
#pest2YesXorganic = b_organic * (P2 == 2) + b_pestYesXorganic * (pest==1)
#pest2YesXchemical = b_chemical * (P2 == 3) + b_pestYesXchemical * (pest==1)
#production system and the age
b1_naturalXage = b_natural * (P1 == 1) + b_naturalXage * age
b1_organicXage = b_organic * (P1 == 2) + b_organicXage * age
b1_chemicalXage = b_chemical * (P1 == 3) + b_chemicalXage * age
b2_naturalXage = b_natural * (P2 == 1) + b_naturalXage * age
b2_organicXage = b_organic * (P2 == 2) + b_organicXage * age
b2_chemicalXage = b_chemical * (P2 == 3) + b_chemicalXage * age
# Create list of utilities
V = list()
V[["alt1"]] = asc1 + b_natural * (P1 == 1) + b_organic * (P1 == 2) + b_chemical * (P1 == 3) + b_conventional * (P1 == 4) +
#b1_naturalXwomen + b1_organicXwomen + b1_chemicalXwomen +
#pest1YesXnatural + pest1YesXorganic + pest1YesXchemical +
b1_naturalXage + b1_organicXage + b1_chemicalXage +
b_wood_cacao * (Micro1 == 1) + b_woodfruit_cacao * (Micro1 == 2) +
b_fruit_cacao * (Micro1 == 3) + b_cacao * (Micro1 == 4) +
b_government * (R1 == 1) + b_farmers * (R1 == 2) + b_scientifics * (R1 == 3) +
b_norecomm * (R1 == 4) + b_pubicinstitu * (CC1 == 1) + b_intermediaries * (CC1 == 2) +
b_privateenter * (CC1 == 3) + b_agriorganiza * (CC1 == 4) + b_compensation * Comp1
V[["alt2"]] = asc2 + b_natural * (P2 == 1) + b_organic * (P2 == 2) + b_chemical * (P2 == 3) + b_conventional * (P2 == 4) +
#b2_naturalXwomen + b2_organicXwomen + b2_chemicalXwomen +
#pest2YesXnatural + pest2YesXorganic + pest2YesXchemical +
b2_naturalXage + b2_organicXage + b2_chemicalXage +
b_wood_cacao * (Micro2 == 1) + b_woodfruit_cacao * (Micro2 == 2) +
b_fruit_cacao * (Micro2 == 3) + b_cacao * (Micro2 == 4) +
b_government * (R2 == 1) + b_farmers * (R2 == 2) + b_scientifics * (R2 == 3) +
b_norecomm * (R2 == 4) + b_pubicinstitu * (CC2 == 1) + b_intermediaries * (CC2 == 2) +
b_privateenter * (CC2 == 3) + b_agriorganiza * (CC2 == 4) + b_compensation * Comp2
V[["optout"]] = b_out
# Define settings for MNL model component
mnl_settings = list(
alternatives = c(alt1=1, alt2=2, optout=3),
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)
# Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}
# Model estimation for Provincia == 1
model_Guayas = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities_1, apollo_inputs_1)
# Output results for Provincia == 1
apollo_modelOutput(model_Guayas)
apollo_modelOutput(model_Guayas, list(printPVal = TRUE))
apollo_saveOutput(model_Guayas)
The results are:
Preparing user-defined functions.
Testing likelihood function...
Setting "avail" is missing, so full availability is assumed.
Overview of choices for MNL model component :
alt1 alt2 optout
Times available 15936.00 15936.00 15936.00
Times chosen 6144.00 5460.00 4332.00
Percentage chosen overall 38.55 34.26 27.18
Percentage chosen when available 38.55 34.26 27.18
Pre-processing likelihood function...
Testing influence of parameters
Starting main estimation
BGW using analytic model derivatives supplied by caller...
Iterates will be written to:
output/Model_Provincia_1_iterations.csv
it nf F RELDF PRELDF RELDX MODEL stppar
0 1 1.750748543e+04
1 4 1.698849281e+04 2.964e-02 2.814e-02 1.00e+00 G 1.11e+00
2 5 1.676987776e+04 1.287e-02 1.297e-02 6.11e-01 G 6.37e-02
3 6 1.674023903e+04 1.767e-03 2.125e-03 3.48e-01 G 1.63e-02
4 7 1.673198389e+04 4.931e-04 5.866e-04 2.01e-01 G 6.82e-16
5 8 1.673074272e+04 7.418e-05 7.247e-05 1.59e-01 S -6.35e-16
6 9 1.673073270e+04 5.989e-07 5.367e-07 2.31e-01 S 6.35e-16
7 10 1.673073188e+04 4.883e-08 5.039e-08 2.29e-01 S 1.37e-15
8 11 1.673073186e+04 1.158e-09 1.139e-09 2.45e-02 S -1.37e-15
9 12 1.673073186e+04 1.136e-11 1.126e-11 4.84e-02 S -1.37e-15
***** Singular convergence *****
Estimated parameters:
Estimate
asc1 -0.554140
asc2 -0.978841
b_natural 0.426446
b_organic 0.272560
b_chemical 0.159269
b_conventional 0.000000
b_wood_cacao 0.039191
b_woodfruit_cacao 0.356813
b_fruit_cacao 0.325766
b_cacao 0.000000
b_government -0.114919
b_farmers 0.197462
b_scientifics -0.168027
b_norecomm 0.000000
b_pubicinstitu -0.032341
b_intermediaries 0.000000
b_privateenter 0.405543
b_agriorganiza 0.229303
b_compensation 0.011912
b_naturalXage -0.009174
b_organicXage 0.006463
b_chemicalXage 0.007309
b_out 0.000000
Final LL: -16730.7319
WARNING: Estimation failed. No covariance matrix to compute.
Current process will resume in 3 seconds unless interrupted by the user...
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
Your model was estimated using the BGW algorithm. Please acknowledge this by citing Bunch et al. (1993) - DOI 10.1145/151271.151279
> apollo_modelOutput(model_Guayas)
Model run by satama using Apollo 0.3.2 on R 4.2.3 for Windows.
Please acknowledge the use of Apollo by citing Hess & Palma (2019)
DOI 10.1016/j.jocm.2019.100170
www.ApolloChoiceModelling.com
Model name : Model_Provincia_1
Model description : MNL model for Provincia == 1
Model run at : 2024-07-15 17:19:54
Estimation method : bgw
Model diagnosis : Singular convergence
Number of individuals : 114
Number of rows in database : 15936
Number of modelled outcomes : 15936
Number of cores used : 1
Model without mixing
LL(start) : -17507.49
LL at equal shares, LL(0) : -17507.49
LL at observed shares, LL(C) : -17346.91
LL(final) : -16730.73
Rho-squared vs equal shares : 0.0444
Adj.Rho-squared vs equal shares : 0.0433
Rho-squared vs observed shares : 0.0355
Adj.Rho-squared vs observed shares : 0.0346
AIC : 33497.46
BIC : 33635.64
Estimated parameters : 18
Time taken (hh:mm:ss) : 00:00:5.99
pre-estimation : 00:00:2.21
estimation : 00:00:3.73
post-estimation : 00:00:0.05
Iterations : 9 (Singular convergence)
Unconstrained optimisation.
Estimates:
Estimate s.e. t.rat.(0) Rob.s.e. Rob.t.rat.(0)
asc1 -0.554140 NA NA NA NA
asc2 -0.978841 NA NA NA NA
b_natural 0.426446 NA NA NA NA
b_organic 0.272560 NA NA NA NA
b_chemical 0.159269 NA NA NA NA
b_conventional 0.000000 NA NA NA NA
b_wood_cacao 0.039191 NA NA NA NA
b_woodfruit_cacao 0.356813 NA NA NA NA
b_fruit_cacao 0.325766 NA NA NA NA
b_cacao 0.000000 NA NA NA NA
b_government -0.114919 NA NA NA NA
b_farmers 0.197462 NA NA NA NA
b_scientifics -0.168027 NA NA NA NA
b_norecomm 0.000000 NA NA NA NA
b_pubicinstitu -0.032341 NA NA NA NA
b_intermediaries 0.000000 NA NA NA NA
b_privateenter 0.405543 NA NA NA NA
b_agriorganiza 0.229303 NA NA NA NA
b_compensation 0.011912 NA NA NA NA
b_naturalXage -0.009174 NA NA NA NA
b_organicXage 0.006463 NA NA NA NA
b_chemicalXage 0.007309 NA NA NA NA
b_out 0.000000 NA NA NA NA
Thank you very much for your response, and although it looks simple, I am quite confused that the interaction conditions are not working properly.
Important: Read this before posting to this forum
- 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.
- 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
- Before asking a question on the forum, users are kindly requested to follow these steps:
- Check that the same issue has not already been addressed in the forum - there is a search tool.
- 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
- Check the frequently asked questions section on the Apollo website, which discusses some common issues/failures. Please see http://www.apollochoicemodelling.com/faq.html
- Make sure that R is using the latest official release of Apollo.
- Users can check which version they are running by entering packageVersion("apollo").
- Then check what is the latest full release (not development version) at http://www.ApolloChoiceModelling.com/code.html.
- To update to the latest official version, just enter install.packages("apollo"). To update to a development version, download the appropriate binary file from http://www.ApolloChoiceModelling.com/code.html, and install the package from file
- If the above steps do not resolve the issue, then users should follow these steps when posting a question:
- provide full details on the issue, including the entire code and output, including any error messages
- posts will not immediately appear on the forum, but will be checked by a moderator first. We check the forum at least twice a week. It may thus take a couple of days for your post to appear and before we reply. There is no need to submit the post multiple times.
Singular convergence
-
stephanehess
- Site Admin
- Posts: 1348
- Joined: 24 Apr 2020, 16:29
Re: Singular convergence
Hi
I think your model is not quite doing what you intended. You have the following:
and then:
First, you are reusing b_natural in the interactions part. Second, you are adding and this will be added to both alternatives, as will etc. These terms will not be identified.
I think what you want is:
I think your model is not quite doing what you intended. You have the following:
Code: Select all
b1_naturalXage = b_natural * (P1 == 1) + b_naturalXage * age
b1_organicXage = b_organic * (P1 == 2) + b_organicXage * age
b1_chemicalXage = b_chemical * (P1 == 3) + b_chemicalXage * age
b2_naturalXage = b_natural * (P2 == 1) + b_naturalXage * age
b2_organicXage = b_organic * (P2 == 2) + b_organicXage * age
b2_chemicalXage = b_chemical * (P2 == 3) + b_chemicalXage * ageCode: Select all
V[["alt1"]] = asc1 + b_natural * (P1 == 1) + b_organic * (P1 == 2) + b_chemical * (P1 == 3) + b_conventional * (P1 == 4) +
#b1_naturalXwomen + b1_organicXwomen + b1_chemicalXwomen +
#pest1YesXnatural + pest1YesXorganic + pest1YesXchemical +
b1_naturalXage + b1_organicXage + b1_chemicalXage +
Code: Select all
b_naturalXage * ageCode: Select all
b_organicXage * ageI think what you want is:
Code: Select all
b1_naturalXage = b_naturalXage * (P1 == 1) * age