diff --git a/Projects/Rbook/parameters_Rbook.R b/Projects/Rbook/parameters_Rbook.R
index 0f32bf6a247bd9a5c7bc005abda94cbb8a298354..aa27915a492b2957b9f4bb385540ebae9384d6a6 100644
--- a/Projects/Rbook/parameters_Rbook.R
+++ b/Projects/Rbook/parameters_Rbook.R
@@ -1,5 +1,7 @@
+rm(list=ls())
+ source("functions.R")
designpath<- "Projects/Rbook/Designs/"
#notes <- "This design consists of different heuristics. One group did not attend the methan attribute, another group only decided based on the payment"
@@ -23,13 +25,6 @@ bheight3=0.50
destype <- "spdesign"
-
-
-
-
-
-
-
#place your utility functions here
u<- list(u1= list(
v1 =V.1 ~ bsq * alt1.sq,
@@ -40,7 +35,7 @@ u<- list(u1= list(
-## logBonus
+sim_all()
diff --git a/functions.R b/functions.R
index 2092551891864a1d1e8065e959551eddba42d963..007dc77625855d76b838049c55dd131ab961f140 100644
--- a/functions.R
+++ b/functions.R
@@ -78,7 +78,113 @@ all_designs[["time"]]=time
return(all_designs)
}
-##test
+
+
+
+readdesign <- function(design = designfile, designtype=destype) {
+
+
+ design <- switch(designtype,
+
+ "ngene" =readr::read_delim(design,delim = "\t",
+ escape_double = FALSE,
+ trim_ws = TRUE ,
+ col_select = c(-Design, -starts_with("...")) ,
+ name_repair = "universal" , show_col_types = FALSE) %>%
+ filter(!is.na(Choice.situation)),
+
+ "spdesign" =as.data.frame(readRDS(design)) %>%
+ mutate(Choice.situation=1:n()) %>%
+ rename_with(~ stringr::str_replace(.,pattern = "_","\\."), everything()),
+
+
+ stop("Invalid value for design. Please provide either 'ngene' or 'spdesign'.")
+ )
+
+}
+
+
+simulate_choices <- function(data=datadet, utility =utils, setspp) { #the part in dataset that needs to be repeated in each run
+
+
+
+ by_formula <- function(equation){ #used to take formulas as inputs in simulation utility function
+ # //! cur_data_all may get deprecated in favor of pick
+ # pick(everything()) %>%
+ cur_data_all() %>%
+ transmute(!!lhs(equation) := !!rhs(equation) )
+ }
+
+ cat(" \n does sou_gis exist: ", exists("sou_gis"), "\n")
+
+ if (exists("sou_gis") && is.function(sou_gis)) {
+ sou_gis()
+
+ cat("\n source of gis has been done \n")
+ }
+
+ # source("Projects/ValuGaps/code/GIS_data.R")
+
+ if(!exists("manipulations")) manipulations=list() ## If no user input on further data manipulations
+
+ n=seq_along(1:length(utility)) # number of utility functions
+
+
+
+
+
+ cat("\n dataset final_set exists: ",exists("final_set"), "\n")
+
+ if(exists("final_set")) data = left_join(data,final_set, by = "ID")
+
+ cat("\n decisiongroups exists: " ,exists("decisiongroups"))
+
+ if(exists("decisiongroups")) { ### create a new variable to classify decision groups.
+
+ data = mutate(data,group = as.numeric(cut(row_number(),
+ breaks = decisiongroups * n(),
+ labels = seq_along(decisiongroups[-length(decisiongroups)]),
+ include.lowest = TRUE)))
+
+ print(table(data$group))
+ } else {
+
+ data$group=1
+ }
+
+
+ data<- data %>%
+ group_by(ID) %>%
+ mutate(!!! manipulations)
+
+
+
+ subsets<- split(data,data$group)
+
+ subsets <- map2(.x = seq_along(u),.y = subsets,
+ ~ mutate(.y,map_dfc(u[[.x]],by_formula)))
+
+ data <-bind_rows(subsets)
+
+ data<- data %>%
+ rename_with(~ stringr::str_replace(.,pattern = "\\.","_"), everything()) %>%
+ mutate(across(.cols=n,.fns = ~ rgumbel(setspp,loc=0, scale=1), .names = "{'e'}_{n}" ),
+ across(starts_with("V_"), .names = "{'U'}_{n}") + across(starts_with("e_")) ) %>% ungroup() %>%
+ mutate(CHOICE=max.col(.[,grep("U_",names(.))])
+ ) %>%
+ as.data.frame()
+
+
+
+
+ cat("\n data has been made \n")
+
+ cat("\n First few observations \n ")
+ print(head(data))
+ cat("\n \n ")
+ return(data)
+
+}
sim_choice <- function(designfile, no_sim=10, respondents=330, mnl_U,utils=u[[1]] ,destype) {
@@ -101,34 +207,13 @@ sim_choice <- function(designfile, no_sim=10, respondents=330, mnl_U,utils=u[[1]
require("rlang")
-
- readdesign <- function(design = designfile, designtype=destype) {
-
-
- design <- switch(designtype,
-
- "ngene" =readr::read_delim(design,delim = "\t",
- escape_double = FALSE,
- trim_ws = TRUE ,
- col_select = c(-Design, -starts_with("...")) ,
- name_repair = "universal" , show_col_types = FALSE) %>%
- filter(!is.na(Choice.situation)),
-
- "spdesign" =as.data.frame(readRDS(design)) %>%
- mutate(Choice.situation=1:n()) %>%
- rename_with(~ stringr::str_replace(.,pattern = "_","\\."), everything()),
-
-
- stop("Invalid value for design. Please provide either 'ngene' or 'spdesign'.")
- )
-
- }
+
mnl_U <-paste(map_chr(utils,as.character,keep.source.attr = TRUE),collapse = "",";") %>%
str_replace_all( c( "priors\\[\"" = "" , "\"\\]" = "" , "~" = "=", "\\." = "_" , " b" = " @b" , "V_"="U_", " alt"="$alt"))
cat("mixl \n")
- print(mnl_U)
+ cat(mnl_U)
cat("\n Simulation \n")
@@ -136,93 +221,19 @@ sim_choice <- function(designfile, no_sim=10, respondents=330, mnl_U,utils=u[[1]
- by_formula <- function(equation){ #used to take formulas as inputs in simulation utility function
- # //! cur_data_all may get deprecated in favor of pick
- # pick(everything()) %>%
- cur_data_all() %>%
- transmute(!!lhs(equation) := !!rhs(equation) )
- }
+
- simulate_choices <- function(data=datadet) { #the part in dataset that needs to be repeated in each run
-cat(" \n does sou_gis exist: ", exists("sou_gis"), "\n")
-
- if (exists("sou_gis") && is.function(sou_gis)) {
- sou_gis()
-
- cat("\n source of gis has been done \n")
- }
-
- # source("Projects/ValuGaps/code/GIS_data.R")
-
- if(!exists("manipulations")) manipulations=list() ## If no user input on further data manipulations
-
- n=seq_along(1:length(utils)) # number of utility functions
-
-#browser()
-
-
- cat("\n dataset final_set exists: ",exists("final_set"), "\n")
-
- if(exists("final_set")) data = left_join(data,final_set, by = "ID")
-
- cat("\n decisiongroups exists: " ,exists("decisiongroups"))
-
- if(exists("decisiongroups")) { ### create a new variable to classify decision groups.
- data = mutate(data,group = as.numeric(cut(row_number(),
- breaks = decisiongroups * n(),
- labels = seq_along(decisiongroups[-length(decisiongroups)]),
- include.lowest = TRUE)))
-
- print(table(data$group))
- } else {
-
- data$group=1
- }
-
-
- data<- data %>%
- group_by(ID) %>%
- mutate(!!! manipulations)
-
-
-
-subsets<- split(data,data$group)
-
-subsets <- map2(.x = seq_along(u),.y = subsets,
- ~ mutate(.y,map_dfc(u[[.x]],by_formula)))
-
-data <-bind_rows(subsets)
-
- data<- data %>%
- rename_with(~ stringr::str_replace(.,pattern = "\\.","_"), everything()) %>%
- mutate(across(.cols=n,.fns = ~ rgumbel(setpp,loc=0, scale=1), .names = "{'e'}_{n}" ),
- across(starts_with("V_"), .names = "{'U'}_{n}") + across(starts_with("e_")) ) %>% ungroup() %>%
- mutate(CHOICE=max.col(.[,grep("U_",names(.))])
- ) %>%
- as.data.frame()
-
-
-
-
-print("\n data has been made \n")
-
-cat("\n First few observations \n ")
-print(head(data))
- cat("\n \n ")
- return(data)
-
- }
estimate_sim <- function(run=1) { #start loop
cat("This is Run number ", run)
- database <- simulate_choices(datadet)
+ database <- simulate_choices(datadet, utility = utils, setspp=setpp )
cat("This is the utility functions \n" , mnl_U)
@@ -236,15 +247,15 @@ print(head(data))
designs_all <- list()
- design<- readdesign()
+ design<- readdesign(design = designfile)
if (!exists("design$Block")) design$Block=1
nsets<-nrow(design)
nblocks<-max(design$Block)
setpp <- nsets/nblocks # Choice Sets per respondent; in this 'no blocks' design everyone sees all 24 sets
- #respondents <- replications*nblocks
+
replications <- respondents/nblocks
- #browser()
+
@@ -259,24 +270,23 @@ as.data.frame()
- database <- simulate_choices()
+database <- simulate_choices(data=datadet, utility = utils, setspp = setpp)
- model_spec <- mixl::specify_model(mnl_U, database, disable_multicore=F)
+model_spec <- mixl::specify_model(mnl_U, database, disable_multicore=F)
- est=setNames(rep(0,length(model_spec$beta_names)), model_spec$beta_names)
+est=setNames(rep(0,length(model_spec$beta_names)), model_spec$beta_names)
- availabilities <- mixl::generate_default_availabilities(
+availabilities <- mixl::generate_default_availabilities(
database, model_spec$num_utility_functions)
- plan(multisession, workers = 8)
-
- output<- 1:no_sim %>% map(estimate_sim)
+
+output<- 1:no_sim %>% map(estimate_sim)
- coefs<-map(1:length(output),~summary(output[[.]])[["coefTable"]][c(1,8)] %>%
+coefs<-map(1:length(output),~summary(output[[.]])[["coefTable"]][c(1,8)] %>%
tibble::rownames_to_column() %>%
pivot_wider(names_from = rowname, values_from = c(est, rob_pval0)) ) %>%
bind_rows(.id = "run")
diff --git a/generatemd.R b/generatemd.R
index d3d9b3a1244fdbcbb542bbd291175d17ca1628c3..5185cbb6dc750f5d65ef6f506aecbf965100c55e 100644
--- a/generatemd.R
+++ b/generatemd.R
@@ -3,9 +3,9 @@ rm(list=ls())
source("functions.R")
-file <- "Projects/feedadditives/parameters_feedadd.R"
+file <- "Projects/Rbook/parameters_Rbook.R"
# file <- "Projects/CSA/parameters_csa.R"
-make_md()
\ No newline at end of file
+make_md()