-
-
Save ncctiglao/d33feaa020702f804f99c386de7df150 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
# ################################################################# # | |
#### LOAD LIBRARY AND DEFINE CORE SETTINGS #### | |
# ################################################################# # | |
### Clear memory | |
rm(list = ls()) | |
### Load libraries | |
library(apollo) | |
library(tidyverse) | |
### Initialise code | |
apollo_initialise() | |
### Set core controls | |
apollo_control = list( | |
modelName = "Hybrid Model for Bicycle to Work", | |
modelDescr = "ICLV model", | |
indivID = "id", | |
mixing = TRUE, | |
nCores = 4, | |
outputDirectory = "output" | |
) | |
# ################################################################# # | |
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS #### | |
# ################################################################# # | |
database = read.csv("ils-biketoworkdata.csv", header=TRUE) | |
#database <- na.omit(database) | |
#View(database) | |
### Subtract mean of indicator variables to center them on zero | |
database$bikelane = database$bikelane - mean(database$bikelane) | |
database$pavement = database$pavement - mean(database$pavement) | |
database$safety = database$safety - mean(database$safety) | |
database$scenery = database$scenery - mean(database$scenery) | |
database$shade = database$shade - mean(database$shade) | |
database$shorttime = database$shorttime - mean(database$shorttime) | |
database$shortestdist = database$shortestdist - mean(database$shortestdist) | |
database$lesstraffic = database$lesstraffic - mean(database$lesstraffic) | |
database$slowvehspeed = database$slowvehspeed - mean(database$slowvehspeed) | |
database$resstreet = database$resstreet - mean(database$resstreet) | |
# ################################################################# # | |
#### DEFINE MODEL PARAMETERS #### | |
# ################################################################# # | |
### Vector of parameters, including any that are kept fixed in estimation | |
apollo_beta=c(asc_yes = 0, | |
b_age =0, | |
b_gender =0, | |
b_income = 0, | |
b_hhead =0, | |
b_educ =0, | |
b_bikeown = 0, | |
b_sector = 0, | |
lambda_safebiking = 0, | |
lambda_environment = 0, | |
lambda_easeoftravel = 0, | |
lambda_lowtraffic = 0, | |
gamma_bikelane = 1, | |
gamma_pavement = 1, | |
gamma_safety = 1, | |
gamma_scenery = 1, | |
gamma_shade = 1, | |
gamma_shorttime = 1, | |
gamma_shortestdist = 1, | |
gamma_lesstraffic = 1, | |
gamma_slowvehspeed = 1, | |
gamma_resstreet = 1, | |
sigma_bikelane = 1, | |
sigma_pavement = 1, | |
sigma_safety = 1, | |
sigma_scenery = 1, | |
sigma_shade = 1, | |
sigma_shorttime = 1, | |
sigma_shortestdist = 1, | |
sigma_lesstraffic = 1, | |
sigma_slowvehspeed = 1, | |
sigma_resstreet = 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="halton", | |
interNDraws=100, | |
interUnifDraws=c(), | |
interNormDraws=c("eta_safebiking", "eta_environment", "eta_easeoftravel", "eta_lowtraffic"), | |
intraDrawsType="", | |
intraNDraws=0, | |
intraUnifDraws=c(), | |
intraNormDraws=c() | |
) | |
### Create random parameters | |
apollo_randCoeff=function(apollo_beta, apollo_inputs){ | |
randcoeff = list() | |
randcoeff[["LV_safebiking"]] = eta_safebiking | |
randcoeff[["LV_environment"]] = eta_environment | |
randcoeff[["LV_easeoftravel"]] = eta_easeoftravel | |
randcoeff[["LV_lowtraffic"]] = eta_lowtraffic | |
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() | |
### Likelihood of indicators | |
normalDensity_settings1 = list(outcomeNormal=bikelane, | |
xNormal=gamma_bikelane*LV_safebiking, | |
mu=0, | |
sigma=sigma_bikelane, | |
componentName='indic_safebiking_1') | |
normalDensity_settings2 = list(outcomeNormal=pavement, | |
xNormal=gamma_pavement*LV_safebiking, | |
mu=0, | |
sigma=sigma_pavement, | |
componentName='indic_safebiking_2') | |
normalDensity_settings3 = list(outcomeNormal=safety, | |
xNormal=gamma_safety*LV_safebiking, | |
mu=0, | |
sigma=sigma_safety, | |
componentName='indic_safebiking_3') | |
normalDensity_settings4 = list(outcomeNormal=scenery, | |
xNormal=gamma_scenery*LV_environment, | |
mu=0, | |
sigma=sigma_scenery, | |
componentName='indic_environment_1') | |
normalDensity_settings5 = list(outcomeNormal=shade, | |
xNormal=gamma_shade*LV_environment, | |
mu=0, | |
sigma=sigma_shade, | |
componentName='indic_environment_2') | |
normalDensity_settings6 = list(outcomeNormal=shorttime, | |
xNormal=gamma_shorttime*LV_easeoftravel, | |
mu=0, | |
sigma=sigma_shorttime, | |
componentName='indic_easoftravel_1') | |
normalDensity_settings7 = list(outcomeNormal=shortestdist, | |
xNormal=gamma_shortestdist*LV_easeoftravel, | |
mu=0, | |
sigma=sigma_shortestdist, | |
componentName='indic_easoftravel_2') | |
normalDensity_settings8 = list(outcomeNormal=lesstraffic, | |
xNormal=gamma_lesstraffic*LV_lowtraffic, | |
mu=0, | |
sigma=sigma_lesstraffic, | |
componentName='indic_lowtraffic_1') | |
normalDensity_settings9 = list(outcomeNormal=slowvehspeed, | |
xNormal=gamma_slowvehspeed*LV_lowtraffic, | |
mu=0, | |
sigma=sigma_slowvehspeed, | |
componentName='indic_lowtraffic_2') | |
normalDensity_settings10= list(outcomeNormal=resstreet, | |
xNormal=gamma_resstreet*LV_lowtraffic, | |
mu=0, | |
sigma=sigma_resstreet, | |
componentName='indic_lowtraffic_3') | |
P[["indic_safebiking_1"]] = apollo_normalDensity(normalDensity_settings1, functionality) | |
P[["indic_safebiking_2"]] = apollo_normalDensity(normalDensity_settings2, functionality) | |
P[["indic_safebiking_3"]] = apollo_normalDensity(normalDensity_settings3, functionality) | |
P[["indic_environment_1"]] = apollo_normalDensity(normalDensity_settings4, functionality) | |
P[["indic_environment_2"]] = apollo_normalDensity(normalDensity_settings5, functionality) | |
P[["indic_easeoftravel_1"]] = apollo_normalDensity(normalDensity_settings6, functionality) | |
P[["indic_easeoftravel_2"]] = apollo_normalDensity(normalDensity_settings7, functionality) | |
P[["indic_lowtraffic_1"]] = apollo_normalDensity(normalDensity_settings8, functionality) | |
P[["indic_lowtraffic_2"]] = apollo_normalDensity(normalDensity_settings9, functionality) | |
P[["indic_lowtraffic_3"]] = apollo_normalDensity(normalDensity_settings10, functionality) | |
### Likelihood of choices | |
### List of utilities: these must use the same names as in mnl_settings, order is irrelevant | |
V = list() | |
V[['No']] = 0 | |
V[['Yes']] = asc_yes + b_age * age + b_gender * gender + b_income * income + + b_hhead * hhead + b_educ * educ + b_sector * sector + b_bikeown * bikeown + | |
lambda_safebiking * LV_safebiking + lambda_environment * LV_environment + lambda_easeoftravel * LV_easeoftravel + lambda_lowtraffic * LV_lowtraffic | |
### Define settings for MNL model component | |
mnl_settings = list( | |
alternatives = c(No=0,Yes=1), | |
avail = 1, | |
choiceVar = biketowork, | |
V = V, | |
componentName = "MNL_biketowork" | |
) | |
### Compute probabilities for MNL model component | |
P[["biketowork"]] = apollo_mnl(mnl_settings, functionality) | |
### Likelihood of the whole model | |
P = apollo_combineModels(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) | |
# ################################################################# # | |
#### MODEL OUTPUTS #### | |
# ################################################################# # | |
# ----------------------------------------------------------------- # | |
#---- FORMATTED OUTPUT (TO SCREEN) ---- | |
# ----------------------------------------------------------------- # | |
apollo_modelOutput(model) | |
# ----------------------------------------------------------------- # | |
#---- FORMATTED OUTPUT (TO FILE, using model name) ---- | |
# ----------------------------------------------------------------- # | |
apollo_saveOutput(model) | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment