ApolloChoiceModelling forum

Code: Select all

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

### Clear memory
rm(list = ls())
s12 <- read_csv("discretechoice.csv")


### Load Apollo library
library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName       = "MOLnautos",
  modelDescr      = "Mixed ordered logit model fitted to MOCOLODO questionnaire",
  indivID         = "uuid2",
  mixing          = TRUE,
  nCores          = 4, 
  outputDirectory = "output"
)

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

### Loading data from package
### if data is to be loaded from a file (e.g. called data.csv),
### the code would be: database = read.csv("data.csv",header=TRUE)
database = s12
### for data dictionary, use ?apollo_drugChoiceData

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

#ref cat strategy: don't use one that is too small and think of what you want to explore. Otherwise choose the middle or extremes.

### Vector of parameters, including any that are kept fixed in estimation
apollo_beta = c(
  sigma_error                                       = 0,
  #beta_randstad_change_randstad                    = 0,
  beta_randstad_change_notrandstad                    = 0,
  beta_randstad_change_torandstad                    = 0,
  beta_randstad_change_fromrandstad                    = 0,
  # beta_leeftijd_group_10_19               = 0,
  # beta_leeftijd_group_20_29               = 0,
  # #beta_leeftijd_group_30_39               = 0,
  # beta_leeftijd_group_40_49               = 0,
  # beta_leeftijd_group_50_59               = 0,
  # beta_leeftijd_group_60_69               = 0,
  # beta_leeftijd_group_70_79               = 0,
  # beta_leeftijd_group_80_                 = 0,
  # beta_leeftijd_group_10_29               = 0,
  # #beta_leeftijd_group_30_49               = 0,
  # beta_leeftijd_group_50_69               = 0,
  # beta_leeftijd_group_70_                 = 0,
  # #beta_huishoudensinkomen__980            = 0,
  # beta_huishoudensinkomen_980_1870        = 0,
  # beta_huishoudensinkomen_1870_2680       = 0,
  # beta_huishoudensinkomen_2680_3800       = 0,
  # beta_huishoudensinkomen_3800_5460       = 0,
  # beta_huishoudensinkomen_5460_           = 0,
  beta_huishoudensinkomen__1870           = 0,
  #beta_huishoudensinkomen_1870_3800       = 0,
  beta_huishoudensinkomen_3800_           = 0,
  #beta_nelectrischefietsen_change_2       = 0,
  beta_nelectrischefietsen_change_3       = 0,
  beta_nelectrischefietsen_change_1      = 0,
  #beta_nnormalefietsen_change_2           = 0,
  beta_nnormalefietsen_change_3           = 0,
  beta_nnormalefietsen_change_1          = 0,
  # #beta_sted_woon_nu_1                     = 0,
  # beta_sted_woon_nu_2                     = 0,
  # beta_sted_woon_nu_3                     = 0,
  # beta_sted_woon_nu_4                     = 0,
  # beta_sted_woon_nu_5                     = 0,
  # beta_sted_woon_nu_NA                    = 0,
  beta_huishoudenstype_anders             = 0,
  beta_huishoudenstype_eengezinshuishouden= 0,
  beta_huishoudenstype_eengezinshuishoudenmetkinderen    = 0,
  beta_huishoudenstype_meergezinshuishoudenmetkinderen   = 0,
  #beta_huishoudenstype_meergezinshuishoudenzonderkinderen= 0,
  beta_weekdist_woonwerk_changecat_decreasedcommute    = 0,
  beta_weekdist_woonwerk_changecat_increasedcommute    = 0,
  beta_weekdist_woonwerk_changecat_samecommute         = 0,
  #beta_weekdist_woonwerk_changecat_nojob              = 0,
  beta_weekdist_woonwerk_changecat_lostjob             = 0,
  beta_weekdist_woonwerk_changecat_gotjob              = 0,
  beta_weekdist_woonwerk_changecat_alwaysworkfromhome  = 0,
  beta_weekdist_woonwerk_changecat_workfromhometocommuting = 0,
  beta_weekdist_woonwerk_changecat_commutingtoworkfromhome = 0,
  beta_weekdist_woonwerk_changecat_differingaddresses  = 0,
  tau_nautos_change_1                = 0,
  tau_nautos_change_2                = 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 = "sobol",
  interNDraws    = 10,
  interUnifDraws = c(),
  interNormDraws = paste0("draws_", c("error")),
  intraDrawsType = "",
  intraNDraws    = 0,
  intraUnifDraws = c(),
  intraNormDraws = c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["error"]]  = sigma_error*draws_error
  
  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()
  
  ### Calculate probabilities using Ordered Logit model
  ol_settings = list(outcomeOrdered = nautos_change,
                     utility        =
                       error +
                       beta_randstad_change_notrandstad*randstad_change_notrandstad +
                       beta_randstad_change_torandstad*randstad_change_torandstad +
                       beta_randstad_change_fromrandstad*randstad_change_fromrandstad +
                       #beta_randstad_randstad_randstad_randstad +
                       # beta_leeftijd_group_10_19*leeftijd_group_10_19 +
                       # beta_leeftijd_group_20_29*leeftijd_group_20_29 +
                       # #beta_leeftijd_group_30_39*leeftijd_group_30_39 +
                       # beta_leeftijd_group_40_49*leeftijd_group_40_49 +
                       # beta_leeftijd_group_50_59*leeftijd_group_50_59 +
                       # beta_leeftijd_group_60_69*leeftijd_group_60_69 +
                       # beta_leeftijd_group_70_79*leeftijd_group_70_79 +
                       # beta_leeftijd_group_80_*leeftijd_group_80_ +
                       # beta_leeftijd_group_10_29*leeftijd_group_10_29 +
                       # #beta_leeftijd_group_30_49*leeftijd_group_30_49 +
                     # beta_leeftijd_group_50_69*leeftijd_group_50_69 +
                     # beta_leeftijd_group_70_*leeftijd_group_70_ +
                     # #beta_huishoudensinkomen__980*huishoudensinkomen__980 +
                     # beta_huishoudensinkomen_980_1870*huishoudensinkomen_980_1870 +
                     # beta_huishoudensinkomen_1870_2680*huishoudensinkomen_1870_2680 +
                     # beta_huishoudensinkomen_2680_3800*huishoudensinkomen_2680_3800 +
                     # beta_huishoudensinkomen_3800_5460*huishoudensinkomen_3800_5460 +
                     # beta_huishoudensinkomen_5460_*huishoudensinkomen_5460_ +
                     beta_huishoudensinkomen__1870*huishoudensinkomen__1870 +
                       #beta_huishoudensinkomen_1870_3800*huishoudensinkomen_1870_3800 +
                       beta_huishoudensinkomen_3800_*huishoudensinkomen_3800_ +
                       #beta_nelectrischefietsen_change_2*nelectrischefietsen_change_2 +
                       beta_nelectrischefietsen_change_3*nelectrischefietsen_change_3 +
                       beta_nelectrischefietsen_change_1*nelectrischefietsen_change_1 +
                       #beta_nnormalefietsen_change_2*nnormalefietsen_change_2 +
                       beta_nnormalefietsen_change_3*nnormalefietsen_change_3 +
                       beta_nnormalefietsen_change_1*nnormalefietsen_change_1 +
                       # #beta_sted_woon_nu_1*sted_woon_nu_1 +
                       # beta_sted_woon_nu_2*sted_woon_nu_2 +
                       # beta_sted_woon_nu_3*sted_woon_nu_3 +
                       # beta_sted_woon_nu_4*sted_woon_nu_4 +
                       # beta_sted_woon_nu_5*sted_woon_nu_5 +
                       # beta_sted_woon_nu_NA*sted_woon_nu_NA +
                       beta_huishoudenstype_anders*huishoudenstype_anders +
                       beta_huishoudenstype_eengezinshuishouden*huishoudenstype_eengezinshuishouden +
                       beta_huishoudenstype_eengezinshuishoudenmetkinderen*huishoudenstype_eengezinshuishoudenmetkinderen +
                       beta_huishoudenstype_meergezinshuishoudenmetkinderen*huishoudenstype_meergezinshuishoudenmetkinderen +
                       #beta_huishoudenstype_meergezinshuishoudenzonderkinderen*huishoudenstype_meergezinshuishoudenzonderkinderen +
                       beta_weekdist_woonwerk_changecat_decreasedcommute*weekdist_woonwerk_changecat_decreasedcommute +
                       beta_weekdist_woonwerk_changecat_increasedcommute*weekdist_woonwerk_changecat_increasedcommute +
                       beta_weekdist_woonwerk_changecat_samecommute*weekdist_woonwerk_changecat_samecommute +
                       #beta_weekdist_woonwerk_changecat_nojob*weekdist_woonwerk_changecat_nojob +
                       beta_weekdist_woonwerk_changecat_lostjob*weekdist_woonwerk_changecat_lostjob +
                       beta_weekdist_woonwerk_changecat_gotjob*weekdist_woonwerk_changecat_gotjob +
                       beta_weekdist_woonwerk_changecat_alwaysworkfromhome*weekdist_woonwerk_changecat_alwaysworkfromhome +
                       beta_weekdist_woonwerk_changecat_workfromhometocommuting*weekdist_woonwerk_changecat_workfromhometocommuting +
                       beta_weekdist_woonwerk_changecat_commutingtoworkfromhome*weekdist_woonwerk_changecat_commutingtoworkfromhome +
                       beta_weekdist_woonwerk_changecat_differingaddresses*weekdist_woonwerk_changecat_differingaddresses,
                     tau             = list(tau_nautos_change_1, tau_nautos_change_2))
  #rows           = (task==1))
  
  P[["model"]] = apollo_ol(ol_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)

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

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

apollo_modelOutput(model)

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

apollo_saveOutput(model)

Code: Select all

Pre-processing likelihood function...
Preparing workers for multithreading...

Testing influence of parametersError in checkForRemoteErrors(lapply(cl, recvResult)) : 
  4 nodes produced errors; first error: Need to average over inter-individual draws! (see ?apollo_avgInterDraws)

Code: Select all

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "sobol",
  interNDraws    = 10,
  interUnifDraws = c(),
  interNormDraws = paste0("draws_", c("error")),
  intraDrawsType = "",
  intraNDraws    = 0,
  intraUnifDraws = c(),
  intraNormDraws = c()
)

Code: Select all

All checks on database completed.
Error in apollo_makeDraws(apollo_inputs, silent = silent) : 
Intra-person draws are used without a panel structure. This is not allowed!

Code: Select all

Need to average over inter-individual draws!

Code: Select all

apollo_avgInterDraws

Code: Select all

apollo_initialise()

apollo_control = list(
  modelName  ="OP_model",
  indivID    ="id_pess", 
  mixing=T,
  nCores = 9)

### Set up starting values for parameters to estimate
apollo_beta <- c(
  # parameters for predictive elements
  sd_multi_1 = 0,
  gamma_multi_1 = 0,
  sd_multi_2 = 3,
  gamma_multi_2 = 3
)
# Fixed parameters, to force identification.

apollo_fixed <- c()

# Set up the form of the maximum simulated likelihood
apollo_draws <- list(
  interDrawsType = "halton",
  interNDraws = 100,
  interNormDraws = c("eta_multi_1","eta_multi_2")
)

# Set up random coefficients definitions
apollo_randCoeff <- function(apollo_beta, apollo_inputs) {
  randcoef <- list()
  randcoef[['tau_multi_1']] = gamma_multi_1 + eta_multi_1*sd_multi_1
  randcoef[['tau_multi_2']] = gamma_multi_2 + eta_multi_2*sd_multi_2
  return(randcoef)
}

### Group and validate input
apollo_inputs <- apollo_validateInputs()

### Set up the probabilities for the choices.
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 for probabilities
  P <- list()

  op_settings_multi <- list(outcomeOrdered = outcome1,
                            V = 0,
                            tau = list(tau_multi_1, tau_multi_2),
                            coding = 1:3)
  #
  P[['model']] <- apollo_op(op_settings_multi, 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)

apollo_modelOutput(model)