Skip to content
Snippets Groups Projects
Commit 05d9bfb4 authored by nc71qaxa's avatar nc71qaxa
Browse files

Update models

parent 7b4a4b0a
No related branches found
No related tags found
No related merge requests found
Expand all Show whitespace changes
Inline Side-by-side
Showing
with 3338 additions and 171 deletions
Scripts/MAKE_FILE.R
+ 1
− 0
View file @ 05d9bfb4
......@@ -67,6 +67,7 @@ mxl_wtp_case_b <- apollo_loadModel("Estimation_results/mxl/MXL_wtp Case B")
mxl_wtp_case_b_NR <- apollo_loadModel("Estimation_results/mxl/MXL_wtp NR B")
mxl_wtp_case_c <- apollo_loadModel("Estimation_results/mxl/MXL_wtp_Case_C")
mxl_wtp_case_c_NR <- apollo_loadModel("Estimation_results/mxl/MXL_wtp_NR_Case_C")
mxl_wtp_case_d <- apollo_loadModel("Estimation_results/mxl/MXL_wtp_Case_D base")
# rent interactions models
mxl_wtp_case_a_rentINT <- apollo_loadModel("Estimation_results/mxl/MXL_wtp Case A Rent Int")
......
Scripts/create_tables.R
+ 12
− 11
View file @ 05d9bfb4
#library(choiceTools, lib.loc = "/home/nc71qaxa/r-packages")
library(choiceTools)
library(choiceTools, lib.loc = "/home/nc71qaxa/r-packages")
#library(choiceTools)
dir.create("Tables/mxl")
dir.create("Tables/logit")
......@@ -131,7 +131,7 @@ texreg(l=list(logit_choice_treat_uni), stars = c(0.01, 0.05, 0.1), float.pos="tb
##### MXL #######
### Baseline case A
case_A <- quicktexregapollo(mxl_wtp_case_a_rentINT)
case_A <- quicktexregapollo(mxl_wtp_case_a)
coef_names <- case_A@coef.names
coef_names <- sub("^(mu_)(.*)(_T|_VT)$", "\\2\\3", coef_names)
......@@ -141,14 +141,14 @@ case_A@coef.names <- coef_names
case_A_cols <- map(c("^mu_", "^sig_", "_T$", "_VT$"), subcoef, case_A)
texreg(c(case_A_cols[1], remGOF(case_A_cols[2:4])),
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "rent" = "Rent",
"ASC_sq" = "ASC SQ", "_natural" = "Naturalness", "nat" = "Naturalness",
custom.coef.map = list("ASC_sq" = "ASC SQ", "natural" = "Naturalness", "walking" = "Walking Distance", "rent" = "Rent",
"_natural" = "Naturalness", "nat" = "Naturalness",
"wd" = "Walking Distance", "asc" = "ASC SQ"),
custom.model.names = c("Mean", "SD", "Treated", "Optional Treatment"), custom.note = "%stars (one-sided). Robust standard errors in parentheses.",
stars = c(0.01, 0.05, 0.1), float.pos="tb",
label = "tab:mxl_A",
caption = "Results of mixed logit model with treatment interactions for Case A.",
file="Tables/mxl/case_A_rent_INT.tex")
file="Tables/mxl/case_A.tex")
### Baseline case C
case_C <- quicktexregapollo(mxl_wtp_case_c_rentINT)
......@@ -184,8 +184,8 @@ case_C_NR@coef.names <- coef_names
case_C_cols_NR <- map(c("^mu_", "^sig_", "_vid1$", "_vid2$", "_nv1$", "_nv2$", "_no_info$", "_NR$"), subcoef, case_C_NR)
texreg(c(case_C_cols_NR[1], remGOF(case_C_cols_NR[2:8])),
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "rent" = "Rent",
"ASC_sq" = "ASC SQ", "_natural" = "Naturalness", "nat" = "Naturalness",
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "ASC_sq" = "ASC SQ", "rent" = "Rent",
"_natural" = "Naturalness", "nat" = "Naturalness",
"wd" = "Walking Distance", "asc" = "ASC SQ",
"ASC_sq_info" = "ASC SQ", "rent_info" = "Rent", "nat_info" = "Naturalness", "walking_info" = "Walking Distance"),
custom.model.names = c("Mean", "SD", "Video 1", "Video 2", "Text 1", "Text 2", "No Info", "NR"), custom.note = "%stars. Robust standard errors in parentheses.",
......@@ -196,7 +196,7 @@ texreg(c(case_C_cols_NR[1], remGOF(case_C_cols_NR[2:8])),
### New Case B
case_B <- quicktexregapollo(mxl_wtp_case_d_rentINT)
case_B <- quicktexregapollo(mxl_wtp_case_d)
coef_names <- case_B@coef.names
coef_names <- sub("^(mu_)(.*)(vol_treated|_treated|_no_info)$", "\\2\\3", coef_names)
......@@ -207,8 +207,8 @@ case_B@coef.names <- coef_names
case_B_cols <- map(c("^mu_", "^sig_", "_treated$", "_vol_treated$","_no_info$"), subcoef, case_B)
texreg(c(case_B_cols[1], remGOF(case_B_cols[2:5])),
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "rent" = "Rent",
"ASC_sq" = "ASC SQ", "_natural" = "Naturalness", "nat" = "Naturalness",
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "ASC_sq" = "ASC SQ", "rent" = "Rent",
"_natural" = "Naturalness", "nat" = "Naturalness",
"wd" = "Walking Distance", "asc" = "ASC SQ",
"ASC_sq_info" = "ASC SQ", "rent_info" = "Rent", "nat_info" = "Naturalness", "walking_info" = "Walking Distance"),
custom.model.names = c("Mean", "SD", "Treated", "Vol. Treated", "No Info"), custom.note = "%stars (one-sided). Robust standard errors in parentheses.",
......@@ -251,3 +251,4 @@ texreg(c(case_B_cols_NR[1], remGOF(case_B_cols_NR[2:6])),
# "mu_asc_T" = "ASC X Treated", "mu_nat_VT" = "Naturalness X Vol. Treated", "mu_wd_VT" = "Walking Distance X Vol. Treated",
# "mu_rent_VT" = "Rent X Vol. Treated", "mu_asc_VT" = "ASC X Vol. Treated"),
# stars = c(0.01, 0.05, 0.1), override.se = mxl_wtp_case_a_rentINT$robse, file="Tables/mxl/case_A_rent_INT.tex")
Scripts/data_prep.R
+ 136
− 134
View file @ 05d9bfb4
......@@ -3,7 +3,7 @@
database_full <- database_full %>% rename(Gender = "Q03W123", Education = "Q06W123", HHSize = "Q41W123",
WorkingTime = "Q44W123", Birthyear = "Q01W123", Rent_net = "Q07W123",
Number_Kids = "Q42W123", Employment_type = "Q43W123", Conseq_UGS = "Q28W3",
Conseq_Money = "Q29W3")
Conseq_Money = "Q29W3", UGS_visits = "Q23W123")
database_full <- database_full %>% mutate(Gender = dplyr::recode(Gender, "A1" = 1, "A2" = 2, "A3"=3),
......@@ -11,7 +11,9 @@ database_full <- database_full %>% mutate(Gender = dplyr::recode(Gender, "A1" =
Employment_type = dplyr::recode(Employment_type, "A1" = 1, "A2" = 2, "A3"=3, "A4" = 4,
"A5" = 5, "A6" = 6),
Conseq_UGS = dplyr::recode(Conseq_UGS, "A1" = 5, "A2" = 4, "A3"=3, "A4" = 2, "A5" = 1, "A6" = NA_real_),
Conseq_Money = dplyr::recode(Conseq_Money, "A1" = 5, "A2" = 4, "A3"=3, "A4" = 2, "A5" = 1, "A6" = NA_real_))
Conseq_Money = dplyr::recode(Conseq_Money, "A1" = 5, "A2" = 4, "A3"=3, "A4" = 2, "A5" = 1, "A6" = NA_real_),
UGS_visits = dplyr::recode(UGS_visits, "A1" = 1, "A2" = 2, "A3"=3, "A4" = 4, "A5" = 5, "A6" = 6,
"A7" = 7, "A8" = 8, "A9" = 9, "A10"= 10))
database_full <- database_full %>% mutate(Gender_female = case_when(Gender == 2 ~1, TRUE~0),
Age = 2023-Birthyear,
......
Scripts/logit/chr_vol_treat.R
+ 43
− 6
View file @ 05d9bfb4
library(margins)
library(pROC)
# Test treatment effect
......@@ -8,7 +9,9 @@ database_full <- database_full %>%
Dummy_Video_2 = case_when(Treatment_new == 5 ~ 1, TRUE ~ 0),
Dummy_no_info = case_when(Treatment_new == 3 ~ 1, TRUE~0),
Dummy_Info_nv1 = case_when(Treatment_new == 2 ~1, TRUE~0),
Dummy_Info_nv2 = case_when(Treatment_new == 4 ~1 , TRUE~0))
Dummy_Info_nv2 = case_when(Treatment_new == 4 ~1 , TRUE~0)) %>%
mutate(Engagement_ugs = groupTime1731 + groupTime1726 + groupTime1769 + groupTime1727 + groupTime1770 + groupTime1771 + groupTime1768 +
groupTime1738)
......@@ -33,9 +36,12 @@ summary(logit_choice_treat)
logit_choice_treat_uni<-glm(Choice_Treat ~ as.factor(Gender)+Z_Mean_NR+Age_mean + QFIncome +
Uni_degree + Kids_Dummy+ as.factor(Q27W123), data, family=binomial)
Uni_degree + Kids_Dummy + Engagement_ugs + UGS_visits, data, family=binomial)
summary(logit_choice_treat_uni)
logit_sig <- glm(Choice_Treat ~ Z_Mean_NR + Age_mean + QFIncome + UGS_visits, data, family=binomial)
# Calculate marginal effects
marginal_effects <- margins(logit_choice_treat_uni)
......@@ -43,14 +49,45 @@ marginal_effects <- margins(logit_choice_treat_uni)
summary(marginal_effects)
cut_off <- 0.45
# Make predictions
data$probabilities <- logit_choice_treat_uni %>% predict(data, type = "response")
data$predicted.classes <- ifelse(data$probabilities >= 0.5, 1, 0)
data$probabilities <- logit_sig %>% predict(data, type = "response")
data$predicted.classes <- ifelse(data$probabilities >= cut_off, 1, 0)
# Model accuracy
mean(data$predicted.classes == data$Choice_Treat, na.rm = T)
table(predicted.classes ,data$Choice_Treat)
table(data$predicted.classes ,data$Choice_Treat)
calculate_metrics <- function(confusion_matrix) {
# Extract values from the confusion matrix
TN <- confusion_matrix[1, 1]
FP <- confusion_matrix[2, 1]
FN <- confusion_matrix[1, 2]
TP <- confusion_matrix[2, 2]
# Calculate sensitivity
sensitivity <- round(TP / (TP + FN), 2)
# Calculate specificity
specificity <- round(TN / (TN + FP), 2)
# Return the results as a list
return(list(sensitivity = sensitivity, specificity = specificity))
}
# Calculate metrics
metrics <- calculate_metrics(table(data$predicted.classes ,data$Choice_Treat))
# Print results
print(paste("Sensitivity:", metrics$sensitivity))
print(paste("Specificity:", metrics$specificity))
roc_obj <- roc(data$Choice_Treat, data$probabilities)
plot(roc_obj, main = "ROC Curve", col = "blue", lwd = 2)
auc_value <- auc(roc_obj)
best_coords <- coords(roc_obj, "best", best.method="youden")
cut_off <- best_coords$threshold
Scripts/mxl/mxl_wtp_space_caseE.R 0 → 100644
+ 175
− 0
View file @ 05d9bfb4
#### Apollo standard script #####
library(apollo) # Load apollo package
# Test treatment effect
database_full$probabilities <- logit_sig %>% predict(database_full, type = "response")
database_full$predicted.classes <- ifelse(database_full$probabilities >= cut_off, 1, 0)
database <- database_full %>%
filter(!is.na(Treatment_new)) %>%
mutate(Dummy_Treated = case_when(Treatment_new == 1|Treatment_new == 2 ~ 1, TRUE ~ 0),
Dummy_Vol_Treated = case_when(Treatment_new == 5 |Treatment_new == 4 ~ 1, TRUE ~ 0),
Dummy_no_info = case_when(Treatment_new == 3 ~ 1, TRUE~0),
Dummy_pred_Treated = case_when(Treatment_new == 6 & predicted.classes == 1 ~ 1, TRUE~0))
table(database$Dummy_Treated)
table(database$Dummy_Vol_Treated)
table(database$Dummy_no_info)
#initialize model
apollo_initialise()
### Set core controls
apollo_control = list(
modelName = "MXL_wtp_Case_E base",
modelDescr = "MXL wtp space Case E base",
indivID ="id",
mixing = TRUE,
HB= FALSE,
nCores = n_cores,
outputDirectory = "Estimation_results/mxl"
)
##### Define model parameters depending on your attributes and model specification! ####
# set values to 0 for conditional logit model
apollo_beta=c(mu_natural = 15,
mu_walking = -1,
mu_rent = -2,
ASC_sq = 0,
mu_ASC_sq_treated = 0,
mu_ASC_sq_vol_treated = 0,
mu_ASC_sq_no_info = 0,
mu_ASC_sq_pred_treat = 0,
mu_nat_treated =0,
mu_nat_vol_treated = 0,
mu_nat_no_info = 0,
mu_nat_pred_treat = 0,
mu_walking_treated =0,
mu_walking_vol_treated = 0,
mu_walking_no_info = 0,
mu_walking_pred_treat = 0,
sig_natural = 15,
sig_walking = 2,
sig_rent = 2,
sig_ASC_sq = 2)
### specify parameters that should be kept fixed, here = none
apollo_fixed = c()
### Set parameters for generating draws, use 2000 sobol draws
apollo_draws = list(
interDrawsType = "sobol",
interNDraws = n_draws,
interUnifDraws = c(),
interNormDraws = c("draws_natural", "draws_walking", "draws_rent", "draws_asc"),
intraDrawsType = "halton",
intraNDraws = 0,
intraUnifDraws = c(),
intraNormDraws = c()
)
### Create random parameters, define distribution of the parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
randcoeff = list()
randcoeff[["b_mu_natural"]] = mu_natural + sig_natural * draws_natural
randcoeff[["b_mu_walking"]] = mu_walking + sig_walking * draws_walking
randcoeff[["b_mu_rent"]] = -exp(mu_rent + sig_rent * draws_rent)
randcoeff[["b_ASC_sq"]] = ASC_sq + sig_ASC_sq * draws_asc
return(randcoeff)
}
### validate
apollo_inputs = apollo_validateInputs()
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 (later integrated in mnl_settings below) ####
# Define utility functions here:
V = list()
V[['alt1']] = -(b_mu_rent)*(b_mu_natural * Naturalness_1 + b_mu_walking * WalkingDistance_1 +
+ mu_nat_treated * Naturalness_1 *Dummy_Treated + mu_nat_no_info * Naturalness_1 * Dummy_no_info
+ mu_nat_vol_treated * Naturalness_1 * Dummy_Vol_Treated
+ mu_walking_treated * WalkingDistance_1 *Dummy_Treated + mu_walking_no_info * WalkingDistance_1 * Dummy_no_info
+ mu_walking_vol_treated * WalkingDistance_1 * Dummy_Vol_Treated - Rent_1
+ mu_nat_pred_treat * Naturalness_1 * Dummy_pred_Treated + mu_walking_pred_treat * WalkingDistance_1 *Dummy_pred_Treated)
V[['alt2']] = -(b_mu_rent)*(b_mu_natural * Naturalness_2 + b_mu_walking * WalkingDistance_2
+ mu_nat_treated * Naturalness_2 *Dummy_Treated + mu_nat_no_info * Naturalness_2 * Dummy_no_info
+ mu_nat_vol_treated * Naturalness_2 * Dummy_Vol_Treated
+ mu_walking_treated * WalkingDistance_2 *Dummy_Treated + mu_walking_no_info * WalkingDistance_2 * Dummy_no_info
+ mu_walking_vol_treated * WalkingDistance_2 * Dummy_Vol_Treated- Rent_2
+ mu_nat_pred_treat * Naturalness_2 * Dummy_pred_Treated + mu_walking_pred_treat * WalkingDistance_2 *Dummy_pred_Treated)
V[['alt3']] = -(b_mu_rent)*(b_ASC_sq + b_mu_natural * Naturalness_3 + b_mu_walking * WalkingDistance_3
+ mu_nat_treated * Naturalness_3 *Dummy_Treated + mu_nat_no_info * Naturalness_3 * Dummy_no_info
+ mu_nat_vol_treated * Naturalness_3 * Dummy_Vol_Treated
+ mu_walking_treated * WalkingDistance_3 *Dummy_Treated + mu_walking_no_info * WalkingDistance_3 * Dummy_no_info
+ mu_walking_vol_treated * WalkingDistance_3 * Dummy_Vol_Treated
+ mu_ASC_sq_treated * Dummy_Treated + mu_ASC_sq_vol_treated * Dummy_Vol_Treated
+ mu_ASC_sq_no_info * Dummy_no_info - Rent_3
+ mu_nat_pred_treat * Naturalness_3 * Dummy_pred_Treated + mu_walking_pred_treat * WalkingDistance_3 *Dummy_pred_Treated
+ mu_ASC_sq_pred_treat * Dummy_pred_Treated)
### Define settings for MNL model component
mnl_settings = list(
alternatives = c(alt1=1, alt2=2, alt3=3),
avail = 1, # all alternatives are available in every choice
choiceVar = choice,
V = V#, # tell function to use list vector defined above
)
### 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 - nur bei Mixed Logit!
P = apollo_avgInterDraws(P, apollo_inputs, functionality)
### Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}
# ################################################################# #
#### MODEL ESTIMATION ##
# ################################################################# #
# estimate model with bfgs algorithm
mxl_wtp_case_e = apollo_estimate(apollo_beta, apollo_fixed,
apollo_probabilities, apollo_inputs,
estimate_settings=list(maxIterations=400,
estimationRoutine="bfgs",
hessianRoutine="analytic"))
# ################################################################# #
#### MODEL OUTPUTS ##
# ################################################################# #
apollo_saveOutput(mxl_wtp_case_e)
Scripts/mxl/mxl_wtp_space_caseF.R 0 → 100644
+ 183
− 0
View file @ 05d9bfb4
#### Apollo standard script #####
library(apollo) # Load apollo package
# Test treatment effect
database_full$probabilities <- logit_sig %>% predict(database_full, type = "response")
database_full$predicted.classes <- ifelse(database_full$probabilities >= cut_off, 1, 0)
database <- database_full %>%
filter(!is.na(Treatment_new)) %>%
mutate(Dummy_Treated = case_when(Treatment_new == 1|Treatment_new == 2 ~ 1, TRUE ~ 0),
Dummy_Vol_Treated = case_when(Treatment_new == 5 |Treatment_new == 4 ~ 1, TRUE ~ 0),
Dummy_no_info = case_when(Treatment_new == 3 ~ 1, TRUE~0),
Dummy_pred_Treated = case_when(Treatment_new == 6 & predicted.classes == 1 ~ 1, TRUE~0),
Dummy_Treated_predicted = case_when(Treatment_new == 1 & predicted.classes == 1|Treatment_new == 2 & predicted.classes == 1 ~ 1, TRUE ~ 0),
Dummy_Treated_notpredicted = case_when(Treatment_new == 1 & predicted.classes == 0|Treatment_new == 2 & predicted.classes == 0 ~ 1, TRUE ~ 0))
table(database$Dummy_Treated)
table(database$Dummy_Vol_Treated)
table(database$Dummy_no_info)
#initialize model
apollo_initialise()
### Set core controls
apollo_control = list(
modelName = "MXL_wtp_Case_F base",
modelDescr = "MXL wtp space Case F base",
indivID ="id",
mixing = TRUE,
HB= FALSE,
nCores = n_cores,
outputDirectory = "Estimation_results/mxl"
)
##### Define model parameters depending on your attributes and model specification! ####
# set values to 0 for conditional logit model
apollo_beta=c(mu_natural = 15,
mu_walking = -1,
mu_rent = -2,
ASC_sq = 0,
mu_ASC_sq_treated_predicted = 0,
mu_ASC_sq_treated_notpredicted = 0,
mu_ASC_sq_vol_treated = 0,
mu_ASC_sq_no_info = 0,
mu_ASC_sq_pred_treat = 0,
mu_nat_treated_predicted =0,
mu_nat_treated_notpredicted = 0,
mu_nat_vol_treated = 0,
mu_nat_no_info = 0,
mu_nat_pred_treat = 0,
mu_walking_treated_predicted =0,
mu_walking_treated_notpredicted =0,
mu_walking_vol_treated = 0,
mu_walking_no_info = 0,
mu_walking_pred_treat = 0,
sig_natural = 15,
sig_walking = 2,
sig_rent = 2,
sig_ASC_sq = 2)
### specify parameters that should be kept fixed, here = none
apollo_fixed = c()
### Set parameters for generating draws, use 2000 sobol draws
apollo_draws = list(
interDrawsType = "sobol",
interNDraws = n_draws,
interUnifDraws = c(),
interNormDraws = c("draws_natural", "draws_walking", "draws_rent", "draws_asc"),
intraDrawsType = "halton",
intraNDraws = 0,
intraUnifDraws = c(),
intraNormDraws = c()
)
### Create random parameters, define distribution of the parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
randcoeff = list()
randcoeff[["b_mu_natural"]] = mu_natural + sig_natural * draws_natural
randcoeff[["b_mu_walking"]] = mu_walking + sig_walking * draws_walking
randcoeff[["b_mu_rent"]] = -exp(mu_rent + sig_rent * draws_rent)
randcoeff[["b_ASC_sq"]] = ASC_sq + sig_ASC_sq * draws_asc
return(randcoeff)
}
### validate
apollo_inputs = apollo_validateInputs()
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 (later integrated in mnl_settings below) ####
# Define utility functions here:
V = list()
V[['alt1']] = -(b_mu_rent)*(b_mu_natural * Naturalness_1 + b_mu_walking * WalkingDistance_1 +
+ mu_nat_treated_predicted * Naturalness_1 *Dummy_Treated_predicted + mu_nat_no_info * Naturalness_1 * Dummy_no_info
+ mu_nat_vol_treated * Naturalness_1 * Dummy_Vol_Treated + mu_nat_treated_notpredicted * Naturalness_1 *Dummy_Treated_notpredicted
+ mu_walking_treated_predicted * WalkingDistance_1 *Dummy_Treated_predicted + mu_walking_treated_notpredicted * WalkingDistance_1 *Dummy_Treated_notpredicted
+ mu_walking_no_info * WalkingDistance_1 * Dummy_no_info
+ mu_walking_vol_treated * WalkingDistance_1 * Dummy_Vol_Treated - Rent_1
+ mu_nat_pred_treat * Naturalness_1 * Dummy_pred_Treated + mu_walking_pred_treat * WalkingDistance_1 * Dummy_pred_Treated)
V[['alt2']] = -(b_mu_rent)*(b_mu_natural * Naturalness_2 + b_mu_walking * WalkingDistance_2
+ mu_nat_treated_predicted * Naturalness_2 *Dummy_Treated_predicted + mu_nat_no_info * Naturalness_2 * Dummy_no_info
+ mu_nat_vol_treated * Naturalness_2 * Dummy_Vol_Treated + mu_nat_treated_notpredicted * Naturalness_2 *Dummy_Treated_notpredicted
+ mu_walking_treated_predicted * WalkingDistance_2 *Dummy_Treated_predicted + mu_walking_treated_notpredicted * WalkingDistance_2 *Dummy_Treated_notpredicted
+ mu_walking_no_info * WalkingDistance_2 * Dummy_no_info
+ mu_walking_vol_treated * WalkingDistance_2 * Dummy_Vol_Treated- Rent_2
+ mu_nat_pred_treat * Naturalness_2 * Dummy_pred_Treated + mu_walking_pred_treat * WalkingDistance_2 * Dummy_pred_Treated)
V[['alt3']] = -(b_mu_rent)*(b_ASC_sq + b_mu_natural * Naturalness_3 + b_mu_walking * WalkingDistance_3
+ mu_nat_treated_predicted * Naturalness_3 *Dummy_Treated_predicted + mu_nat_no_info * Naturalness_3 * Dummy_no_info
+ mu_nat_vol_treated * Naturalness_3 * Dummy_Vol_Treated + mu_nat_treated_notpredicted * Naturalness_3 *Dummy_Treated_notpredicted
+ mu_walking_treated_predicted * WalkingDistance_3 *Dummy_Treated_predicted + mu_walking_treated_notpredicted * WalkingDistance_3 *Dummy_Treated_notpredicted
+ mu_walking_no_info * WalkingDistance_3 * Dummy_no_info
+ mu_walking_vol_treated * WalkingDistance_3 * Dummy_Vol_Treated
+ mu_ASC_sq_treated_predicted * Dummy_Treated_predicted + mu_ASC_sq_vol_treated * Dummy_Vol_Treated
+ mu_ASC_sq_no_info * Dummy_no_info - Rent_3 + mu_ASC_sq_treated_notpredicted * Dummy_Treated_notpredicted
+ mu_nat_pred_treat * Naturalness_3 * Dummy_pred_Treated + mu_walking_pred_treat * WalkingDistance_3 * Dummy_pred_Treated
+ mu_ASC_sq_pred_treat * Dummy_pred_Treated)
### Define settings for MNL model component
mnl_settings = list(
alternatives = c(alt1=1, alt2=2, alt3=3),
avail = 1, # all alternatives are available in every choice
choiceVar = choice,
V = V#, # tell function to use list vector defined above
)
### 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 - nur bei Mixed Logit!
P = apollo_avgInterDraws(P, apollo_inputs, functionality)
### Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}
# ################################################################# #
#### MODEL ESTIMATION ##
# ################################################################# #
# estimate model with bfgs algorithm
mxl_wtp_case_f = apollo_estimate(apollo_beta, apollo_fixed,
apollo_probabilities, apollo_inputs,
estimate_settings=list(maxIterations=400,
estimationRoutine="bfgs",
hessianRoutine="analytic"))
# ################################################################# #
#### MODEL OUTPUTS ##
# ################################################################# #
apollo_saveOutput(mxl_wtp_case_f)
Scripts/mxl/mxl_wtp_space_caseG.R 0 → 100644
+ 176
− 0
View file @ 05d9bfb4
#### Apollo standard script #####
library(apollo) # Load apollo package
# Test treatment effect
database_full$probabilities <- logit_sig %>% predict(database_full, type = "response")
database_full$predicted.classes <- ifelse(database_full$probabilities >= cut_off, 1, 0)
database <- database_full %>%
filter(!is.na(Treatment_new)) %>%
mutate(Dummy_Treated = case_when(Treatment_new == 1|Treatment_new == 2 ~ 1, TRUE ~ 0),
Dummy_Vol_Treated = case_when(Treatment_new == 5 |Treatment_new == 4 ~ 1, TRUE ~ 0),
Dummy_no_info = case_when(Treatment_new == 3 ~ 1, TRUE~0),
Dummy_pred_Treated = case_when(Treatment_new == 6 & predicted.classes == 1 ~ 1, TRUE~0),
Dummy_Treated_predicted = case_when(Treatment_new == 1 & predicted.classes == 1|Treatment_new == 2 & predicted.classes == 1 ~ 1, TRUE ~ 0),
Dummy_Treated_notpredicted = case_when(Treatment_new == 1 & predicted.classes == 0|Treatment_new == 2 & predicted.classes == 0 ~ 1, TRUE ~ 0))
table(database$Dummy_Treated)
table(database$Dummy_Vol_Treated)
table(database$Dummy_no_info)
#initialize model
apollo_initialise()
### Set core controls
apollo_control = list(
modelName = "MXL_wtp_Case_G base",
modelDescr = "MXL wtp space Case G base",
indivID ="id",
mixing = TRUE,
HB= FALSE,
nCores = n_cores,
outputDirectory = "Estimation_results/mxl"
)
##### Define model parameters depending on your attributes and model specification! ####
# set values to 0 for conditional logit model
apollo_beta=c(mu_natural = 15,
mu_walking = -1,
mu_rent = -2,
ASC_sq = 0,
mu_ASC_sq_treated_predicted = 0,
mu_ASC_sq_treated_notpredicted = 0,
mu_ASC_sq_vol_treated = 0,
mu_ASC_sq_no_info = 0,
mu_nat_treated_predicted =0,
mu_nat_treated_notpredicted = 0,
mu_nat_vol_treated = 0,
mu_nat_no_info = 0,
mu_walking_treated_predicted =0,
mu_walking_treated_notpredicted =0,
mu_walking_vol_treated = 0,
mu_walking_no_info = 0,
sig_natural = 15,
sig_walking = 2,
sig_rent = 2,
sig_ASC_sq = 2)
### specify parameters that should be kept fixed, here = none
apollo_fixed = c()
### Set parameters for generating draws, use 2000 sobol draws
apollo_draws = list(
interDrawsType = "sobol",
interNDraws = n_draws,
interUnifDraws = c(),
interNormDraws = c("draws_natural", "draws_walking", "draws_rent", "draws_asc"),
intraDrawsType = "halton",
intraNDraws = 0,
intraUnifDraws = c(),
intraNormDraws = c()
)
### Create random parameters, define distribution of the parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
randcoeff = list()
randcoeff[["b_mu_natural"]] = mu_natural + sig_natural * draws_natural
randcoeff[["b_mu_walking"]] = mu_walking + sig_walking * draws_walking
randcoeff[["b_mu_rent"]] = -exp(mu_rent + sig_rent * draws_rent)
randcoeff[["b_ASC_sq"]] = ASC_sq + sig_ASC_sq * draws_asc
return(randcoeff)
}
### validate
apollo_inputs = apollo_validateInputs()
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 (later integrated in mnl_settings below) ####
# Define utility functions here:
V = list()
V[['alt1']] = -(b_mu_rent)*(b_mu_natural * Naturalness_1 + b_mu_walking * WalkingDistance_1 +
+ mu_nat_treated_predicted * Naturalness_1 *Dummy_Treated_predicted + mu_nat_no_info * Naturalness_1 * Dummy_no_info
+ mu_nat_vol_treated * Naturalness_1 * Dummy_Vol_Treated + mu_nat_treated_notpredicted * Naturalness_1 *Dummy_Treated_notpredicted
+ mu_walking_treated_predicted * WalkingDistance_1 *Dummy_Treated_predicted + mu_walking_treated_notpredicted * WalkingDistance_1 *Dummy_Treated_notpredicted
+ mu_walking_no_info * WalkingDistance_1 * Dummy_no_info
+ mu_walking_vol_treated * WalkingDistance_1 * Dummy_Vol_Treated - Rent_1)
V[['alt2']] = -(b_mu_rent)*(b_mu_natural * Naturalness_2 + b_mu_walking * WalkingDistance_2
+ mu_nat_treated_predicted * Naturalness_2 *Dummy_Treated_predicted + mu_nat_no_info * Naturalness_2 * Dummy_no_info
+ mu_nat_vol_treated * Naturalness_2 * Dummy_Vol_Treated + mu_nat_treated_notpredicted * Naturalness_2 *Dummy_Treated_notpredicted
+ mu_walking_treated_predicted * WalkingDistance_2 *Dummy_Treated_predicted + mu_walking_treated_notpredicted * WalkingDistance_2 *Dummy_Treated_notpredicted
+ mu_walking_no_info * WalkingDistance_2 * Dummy_no_info
+ mu_walking_vol_treated * WalkingDistance_2 * Dummy_Vol_Treated- Rent_2)
V[['alt3']] = -(b_mu_rent)*(b_ASC_sq + b_mu_natural * Naturalness_3 + b_mu_walking * WalkingDistance_3
+ mu_nat_treated_predicted * Naturalness_3 *Dummy_Treated_predicted + mu_nat_no_info * Naturalness_3 * Dummy_no_info
+ mu_nat_vol_treated * Naturalness_3 * Dummy_Vol_Treated + mu_nat_treated_notpredicted * Naturalness_3 *Dummy_Treated_notpredicted
+ mu_walking_treated_predicted * WalkingDistance_3 *Dummy_Treated_predicted + mu_walking_treated_notpredicted * WalkingDistance_3 *Dummy_Treated_notpredicted
+ mu_walking_no_info * WalkingDistance_3 * Dummy_no_info
+ mu_walking_vol_treated * WalkingDistance_3 * Dummy_Vol_Treated
+ mu_ASC_sq_treated_predicted * Dummy_Treated_predicted + mu_ASC_sq_vol_treated * Dummy_Vol_Treated
+ mu_ASC_sq_no_info * Dummy_no_info - Rent_3 + mu_ASC_sq_treated_notpredicted * Dummy_Treated_notpredicted)
### Define settings for MNL model component
mnl_settings = list(
alternatives = c(alt1=1, alt2=2, alt3=3),
avail = 1, # all alternatives are available in every choice
choiceVar = choice,
V = V#, # tell function to use list vector defined above
)
### 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 - nur bei Mixed Logit!
P = apollo_avgInterDraws(P, apollo_inputs, functionality)
### Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}
# ################################################################# #
#### MODEL ESTIMATION ##
# ################################################################# #
# estimate model with bfgs algorithm
mxl_wtp_case_g = apollo_estimate(apollo_beta, apollo_fixed,
apollo_probabilities, apollo_inputs,
estimate_settings=list(maxIterations=400,
estimationRoutine="bfgs",
hessianRoutine="analytic"))
# ################################################################# #
#### MODEL OUTPUTS ##
# ################################################################# #
apollo_saveOutput(mxl_wtp_case_g)
custom_app.lua 0 → 100644
+ 16
− 0
View file @ 05d9bfb4
local in_appendix = false
Header = function(h)
if h.level == 1 then
if h.classes:includes("appendix") then
in_appendix = true
h.attributes["visibility"] = "uncounted"
else
in_appendix = false
end
end
if h.level == 2 and in_appendix then
h.attributes["visibility"] = "uncounted"
end
return h
end
project_start.qmd
+ 17
− 20
View file @ 05d9bfb4
......@@ -20,7 +20,7 @@ format:
slide-number: true
smaller: true
logo: Grafics/iDiv_logo_item.PNG
footer: "WONV 2024: Hot Cities, Cool Choices"
footer: "DCE Network meeting"
scrollable: true
embed-resources: true
---
......@@ -155,7 +155,7 @@ list_ols <- list("(Intercept)" = "Intercept", "as.factor(Treatment_A)Treated" =
+ **Socio-demographics**: Age, Gender, Income, Education
+ **Attitudinal variable**: Measure derived from 21 items on **nature relatedness** (NR-Index) [@nisbet2009nature]
+ **Attitudinal variable**: Measure derived from 21 items on **nature relatedness** [@nisbet2009nature]
:::
......@@ -237,11 +237,11 @@ list_ols <- list("(Intercept)" = "Intercept", "as.factor(Treatment_A)Treated" =
<!-- ::: -->
## OLS: Engagement
::: panel-tabset
### OLS Specification
## OLS Specification
```{=tex}
\begin{equation}
......@@ -277,13 +277,12 @@ with
```
### Results
## OLS: Engagement
```{r}
ggpubr::ggarrange(plot_interview_A, plot_cc_A)
```
:::
## OLS: Information Recall & Consequentiality
......@@ -293,13 +292,11 @@ ggpubr::ggarrange(plot_mani_A, plot_cons_A)
## MXL: Effects on Stated Preferences
::: panel-tabset
### MXL Specification
```{=tex}
```{=latex}
\begin{equation}
U_i = -(\beta_{C_i} + \beta_{TreatC_i} \cdot v_{Treat}) \cdot (\beta_{X_i} \cdot v_{X_i} + \beta_{TreatX_i} \cdot v_{X_i} \cdot v_{Treat} - C_i) + \epsilon_i
U_i = -\beta_{C_i} \cdot (\beta_{X_i} \cdot v_{X_i} + \beta_{TreatX_i} \cdot v_{X_i} \cdot v_{Treat} - C_i) + \epsilon_i
\label{mxl_base}
\end{equation}
```
......@@ -322,20 +319,20 @@ with
\end{equation}
```
### Results
## MXL: Effects on Stated Preferences
::: {style="font-size: 90%;"}
::: {style="font-size: 80%;"}
```{r, results='asis'}
htmlreg(c(case_A_cols[1], remGOF(case_A_cols[2:4])), single.row = TRUE,
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "rent" = "Rent",
"ASC_sq" = "ASC SQ", "_natural" = "Naturalness", "nat" = "Naturalness",
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "ASC_sq" = "ASC SQ", "rent" = "Rent",
"_natural" = "Naturalness", "nat" = "Naturalness",
"wd" = "Walking Distance", "asc" = "ASC SQ"),
custom.model.names = c("Mean", "SD", "Treated", "Optional Treatment"), custom.note = "%stars. Standard errors in parentheses.",
stars = c(0.01, 0.05, 0.1), float.pos="tb",
label = "tab:mxl_A",
caption = "Results of mixed logit model with treatment interactions for Case A.")
```
:::
:::
......@@ -398,8 +395,8 @@ ggpubr::ggarrange(plot_mani_D, plot_cons_D)
::: {style="font-size: 80%;"}
```{r, results='asis'}
htmlreg(c(case_B_cols[1], remGOF(case_B_cols[2:5])),
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "rent" = "Rent",
"ASC_sq" = "ASC SQ", "_natural" = "Naturalness", "nat" = "Naturalness",
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "ASC_sq" = "ASC SQ", "rent" = "Rent",
"_natural" = "Naturalness", "nat" = "Naturalness",
"wd" = "Walking Distance", "asc" = "ASC SQ",
"ASC_sq_info" = "ASC SQ", "rent_info" = "Rent", "nat_info" = "Naturalness", "walking_info" = "Walking Distance"),
custom.model.names = c("Mean", "SD", "Treated", "Vol. Treated", "No Info"), custom.note = "%stars (one-sided). Robust standard errors in parentheses.",
......@@ -414,8 +411,8 @@ htmlreg(c(case_B_cols[1], remGOF(case_B_cols[2:5])),
::: {style="font-size: 80%;"}
```{r, results='asis'}
htmlreg(c(case_B_cols_NR[1], remGOF(case_B_cols_NR[2:6])),
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "rent" = "Rent",
"ASC_sq" = "ASC SQ", "_natural" = "Naturalness", "nat" = "Naturalness",
custom.coef.map = list("natural" = "Naturalness", "walking" = "Walking Distance", "ASC_sq" = "ASC SQ", "rent" = "Rent",
"_natural" = "Naturalness", "nat" = "Naturalness",
"wd" = "Walking Distance", "asc" = "ASC SQ", "ASC" ="ASC SQ",
"ASC_sq_info" = "ASC SQ", "rent_info" = "Rent", "nat_info" = "Naturalness", "walking_info" = "Walking Distance"),
custom.model.names = c("Mean", "SD", "Treated", "Vol. Treated", "No Info", "NR-Index"), custom.note = "%stars (one-sided). Robust standard errors in parentheses.",
......
project_start_wonv.qmd 0 → 100644
+ 598
− 0
View file @ 05d9bfb4
This diff is collapsed.
references.bib 0 → 100644
+ 1981
− 0
View file @ 05d9bfb4
This diff is collapsed.
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment