Page 1 of 1

Correlation of MMNL random parameters calculated from unconditionals

Posted: 12 Feb 2025, 10:34
by MarianneLefebvre
Dear Stephane
In Apollo manual (p20), you indicate that correlations between MMNL random parameters can be calculated empirically from the draws produced by apollo_unconditionals. You you please provide the code to get that ? I could not find into the manual.
Thanks
Marianne

My code for the MMNL with 1 lognormal and 6 normal coefficients is:

Code: Select all

# ################################################################# #
#### LOAD LIBRARY AND DEFINE CORE SETTINGS                       ####
# ################################################################# #

### Clear memory

rm(list = ls())

### Load Apollo library
library(apollo)
#library(xtable)

### Initialise code
apollo_initialise()

### set working directory

setwd("C:/Users/marianne.lefebvre/Nextcloud/RECHERCHE/BEHAVE/BEHAVEPartagePauline/EAU/DCE_Eau/Data_analysis/Rapollo_RevisionERAE")


### Set core controls
apollo_control = list(
  modelName       = "MMNL_DCEEau_WTPspace",
  modelDescr      = "MMNL on Behave water DCE data in WTP space",
  indivID         = "idLS",
  nCores          = 3, 
  outputDirectory = "output"
)

# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS                     ####
# ################################################################# #

### Loading data from package
database = read.csv("Final_Result_Wide.csv",header=TRUE)
                     
settings<-list(
alternative_column="Quchoice",
alternative_specific_attributes=c("DeclaFarm","DeclaPlotSoil","DeclaAuto","ManualWeekly","Smart","Vol"),
choice_column="Choice",
ID_column="idLS",
observation_column="Choicenumber"
)


#database = apollo_longToWide(databaseLong2,longToWide_settings=settings)

# ################################################################# #
#### DEFINE MODEL PARAMETERS                                     ####
# ################################################################# #

#Estimates MNL
#                    Estimate    s.e.       t.rat.(0)
#ASC                 0.96904    0.242249      4.0002
#b_Vol               0.01871    0.002577      7.2585
#b_DeclaFarm         0.66822    0.138640      4.8198
#b_DeclaPlotSoil     0.85329    0.170235      5.0124
#b_DeclaAuto         0.31873    0.168379      1.8929
#b_ManualWeekly     -0.12685    0.230557     -0.5502
#b_Smart            -0.20943    0.220074     -0.9516

### Vector of parameters, including any that are kept fixed in estimation
apollo_beta = c(ASC_mu      =  13.4403, 
                ASC_sig      = 123.2308, 
                v_Vol_mu=  -2.6916 ,
                v_Vol_sig=0.7512,
                v_DeclaFarm_mu     = 47.3078,
                v_DeclaFarm_sig    = 29.6476, 
                v_DeclaPlotSoil_mu     = 47.7545,
                v_DeclaPlotSoil_sig    = 10.1956, 
                v_DeclaAuto_mu     = 26.4988,
                v_DeclaAuto_sig    = -30.1430, 
                v_ManualWeekly_mu     = -20.3407, 
                v_ManualWeekly_sig    = -53.2512, 
                v_Smart_mu     = -15.6283,	
                v_Smart_sig    = 31.1956)

### 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: inter because we have panel data with several choices per ind. we want heterogeneity at ind level (not choice level)
apollo_draws = list(
  interDrawsType = "halton",
  interNDraws    = 2000,
  interNormDraws = c("draws_DeclaFarm","draws_DeclaPlotSoil","draws_DeclaAuto","draws_ManualWeekly","draws_Smart","draws_Vol","draws_ASC" )
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  randcoeff[["b_Vol"]] = exp(v_Vol_mu + v_Vol_sig * draws_Vol )
  randcoeff[["v_DeclaFarm"]] = v_DeclaFarm_mu + v_DeclaFarm_sig * draws_DeclaFarm
  randcoeff[["v_DeclaPlotSoil"]] = v_DeclaPlotSoil_mu + v_DeclaPlotSoil_sig * draws_DeclaPlotSoil
  randcoeff[["v_DeclaAuto"]] = v_DeclaAuto_mu + v_DeclaAuto_sig * draws_DeclaAuto
  randcoeff[["v_ManualWeekly"]] = v_ManualWeekly_mu + v_ManualWeekly_sig * draws_ManualWeekly
  randcoeff[["v_Smart"]] = v_Smart_mu + v_Smart_sig * draws_Smart
  randcoeff[["ASC"]] = ASC_mu + ASC_sig * 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 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[["alt1"]] = (v_DeclaFarm * DeclaFarm1 + v_DeclaPlotSoil * DeclaPlotSoil1 + v_DeclaAuto * DeclaAuto1 + v_ManualWeekly * ManualWeekly1 + v_Smart * Smart1 + Vol1)*b_Vol
  V[["alt2"]] = (v_DeclaFarm * DeclaFarm2 + v_DeclaPlotSoil * DeclaPlotSoil2 + v_DeclaAuto * DeclaAuto2 + v_ManualWeekly * ManualWeekly2 + v_Smart * Smart2 + Vol2)*b_Vol
  V[["alt3"]] = (ASC + v_DeclaFarm * DeclaFarm3 + v_DeclaPlotSoil * DeclaPlotSoil3 + v_DeclaAuto * DeclaAuto3 + v_ManualWeekly * ManualWeekly3 + v_Smart * Smart3 + Vol3)*b_Vol
  
  ### Define settings for MNL model component
  mnl_settings = list(
    alternatives  = c(alt1=1, alt2=2,alt3=3),
    avail         = list(alt1=1, alt2=1,alt3=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                                            ####
# ################################################################# #

### Optional speedTest
#speedTest_settings=list(
#   nDrawsTry = c(100, 200, 500),
#   nCoresTry = c(1,3,5,7),
#   nRep      = 30
#)

#apollo_speedTest(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs, speedTest_settings)

model = apollo_estimate(apollo_beta, apollo_fixed,
                        apollo_probabilities, apollo_inputs)

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

# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO SCREEN)                               ----
# ----------------------------------------------------------------- #

apollo_modelOutput(model,modelOutput_settings = list(printPVal=1)) # for one sided t-test

# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO FILE, using model name)               ----
# ----------------------------------------------------------------- #

apollo_saveOutput(model,saveOutput_settings = list(printPVal=1)) # for one sided t-test


#As expected, the estimated coefficients for data sharing are positive (interpreted as a decrease in WTA: They ask for less compensation in terms of prioritized water volumes)
#The effect size is largely higher for water data reporting (farmers ask for more/less ? compensation in terms of )

#EXPORT TABLE TEX: commande à vérifier
#print(xtable(model,digits = 2,caption="RPL in WTA space",include.rownames = F)

# ################################################################# #
##### POST-PROCESSING                                            ####
# ################################################################# #

### Print outputs of additional diagnostics to new output file (remember to close file writing when complete)
apollo_sink()

beta= apollo_unconditionals(model,apollo_probabilities,apollo_inputs)

Re: Correlation of MMNL random parameters calculated from unconditionals

Posted: 27 Feb 2025, 16:06
by stephanehess
Hi

e.g.

Code: Select all

cor(as.vector(beta[[1]]),as.vector(beta[[2]]))
Stephane

Re: Correlation of MMNL random parameters calculated from unconditionals

Posted: 27 Feb 2025, 17:19
by MarianneLefebvre
Thanks for your reply. I now have the correlation matrix for all the coefficients in the MMNL model.
[,1] [,2] [,3] [,4] [,5] [,6] [,7]
[1,] 1.0e+00 -9.6e-05 -2.0e-04 4.8e-05 8.5e-05 -2.0e-05 -5.5e-05
[2,] -9.6e-05 1.0e+00 -4.7e-05 3.0e-05 3.2e-05 -5.0e-05 -2.2e-05
[3,] -2.0e-04 -4.7e-05 1.0e+00 5.2e-05 3.2e-05 -5.1e-05 -3.3e-05
[4,] 4.8e-05 3.0e-05 5.2e-05 1.0e+00 -5.1e-05 4.4e-05 3.9e-05
[5,] 8.5e-05 3.2e-05 3.2e-05 -5.1e-05 1.0e+00 4.5e-05 5.4e-05
[6,] -2.0e-05 -5.0e-05 -5.1e-05 4.4e-05 4.5e-05 1.0e+00 -3.4e-05
[7,] -5.5e-05 -2.2e-05 -3.3e-05 3.9e-05 5.4e-05 -3.4e-05 1.0e+00

Can you confirm this is sufficient to say we have estimated the model with complete covariance ?
Can we conclude that all the out-of-diagonal coefficients are small therefore not justifying a correlated MMNL ?

Thanks a lot for your help with interpretation!
Best
Marianne

Re: Correlation of MMNL random parameters calculated from unconditionals

Posted: 27 Feb 2025, 17:35
by stephanehess
Hi

no, you cannot reach that conclusion. You have not specified any correlation between the betas, so there can't be any. I assumed you had allowed for correlation, and wanted to see how strong it was. If you specify uncorrelated betas, they will be uncorrelated

Hope this helps

Stephane