scienceverse example for family-wise error control v2

sim_study_scienceverse_2.R

# Scienceverse Sim

# install scienceverse
# devtools::install_github("scienceverse/scienceverse")

library(scienceverse)
library(faux)

set.seed(2) # set.seed(2) is a random draw where H1 is corroborated.

nsim <- 1

res_sim_1 <-numeric(nsim) #set up empty container for all results
res_sim_2 <-numeric(nsim) #set up empty container for all results
res_sim_3 <-numeric(nsim) #set up empty container for all results

for (i in 1:nsim) {
  
  #Set up the study
  sim_study <- study("Simulating Null Effects",
                     author = c("Daniel Lakens", "Lisa DeBruine"))
  
  #Add a hypothesis
  sim_study <- add_hypothesis(
    study = sim_study,
    description = "dv1 will show an effect",
    id = "H1"
  )
  sim_study <- add_hypothesis(
    study = sim_study,
    description = "dv2 will show an effect",
    id = "H2"
  )
  sim_study <- add_hypothesis(
    study = sim_study,
    description = "dv3 will show an effect",
    id = "H3"
  )
  
  # Add an independent t-test as analysis
  # Note we know our dataframe is called dat, and has a column condition, with values control & treatment
  # We also have columns dv1 dv2 dv3 which are our dvs. 
  sim_study <- add_analysis(sim_study,
                            id = "ttest_1",
                            code = t.test(dat[which(dat$condition == "control"),]$dv1, 
                                          dat[which(dat$condition == "treatment"),]$dv1, 
                                          paired = FALSE, 
                                          conf.level = (1-alpha_level)), 
                            software = R.version.string)
  
  sim_study <- add_analysis(sim_study,
                            id = "ttest_2",
                            code = t.test(dat[which(dat$condition == "control"),]$dv2, 
                                          dat[which(dat$condition == "treatment"),]$dv2, 
                                          paired = FALSE, 
                                          conf.level = (1-alpha_level)), 
                            software = R.version.string)
  
  sim_study <- add_analysis(sim_study,
                            id = "ttest_3",
                            code = t.test(dat[which(dat$condition == "control"),]$dv3, 
                                          dat[which(dat$condition == "treatment"),]$dv3, 
                                          paired = FALSE, 
                                          conf.level = (1-alpha_level)), 
                            software = R.version.string)
  
  # Add criterion
  
  sim_study <- add_criterion(
    sim_study,
    id = "p_t_1",
    hypothesis_id = "H1",
    analysis_id = "ttest_1",
    result = "p.value",
    operator = "<",
    comparator = alpha_level) 
  
  sim_study <- add_criterion(
    sim_study,
    id = "p_t_2",
    hypothesis_id = "H2",
    analysis_id = "ttest_2",
    result = "p.value",
    operator = "<",
    comparator = alpha_level) 
  
  sim_study <- add_criterion(
    sim_study,
    id = "p_t_3",
    hypothesis_id = "H3",
    analysis_id = "ttest_3",
    result = "p.value",
    operator = "<",
    comparator = alpha_level) 
  
  
  # Add evaluation. 
  sim_study <- add_eval(sim_study,
                        hypothesis_id = "H1",
                        "corroboration",
                        description = "The hypothesis is corroborated if dv1 is significant.",
                        evaluation = "p_t_1")
  sim_study <- add_eval(sim_study, 
                        hypothesis_id = "H2",
                        "corroboration",
                        description = "The hypothesis is corroborated if dv2 is significant.",
                        evaluation = "p_t_2")
  sim_study <- add_eval(sim_study, 
                        hypothesis_id = "H3",
                        "corroboration",
                        description = "The hypothesis is corroborated if dv3 is significant.",
                        evaluation = "p_t_3")
  sim_study <- add_eval(sim_study,
                        hypothesis_id = "H1",
                        "falsification",
                        description = "The hypothesis is falsified if dv1 is not significant.",
                        evaluation = "!p_t_1")
  sim_study <- add_eval(sim_study, 
                        hypothesis_id = "H2",
                        "falsification",
                        description = "The hypothesis is falsified if dv2 is not significant.",
                        evaluation = "!p_t_2")
  sim_study <- add_eval(sim_study, 
                        hypothesis_id = "H3",
                        "falsification",
                        description = "The hypothesis is falsified if dv3 is not significant.",
                        evaluation = "!p_t_3")
  
  # Simulate some data and add it to the study
  sim_study <- add_sim_data(
    sim_study,
    data_id = "dat",
    within = list(dv = c("dv1", "dv2", "dv3")),
    between = list(condition = c("control", "treatment")),
    n = 100, 
    mu = c(100, 100, 100, 100, 100, 100), 
    sd = 10)
  
  # Take a look at the data
  # sim_study$data[[1]]$data
  
  # Analyze the results
  sim_study <- study_analyze(sim_study)
  
  res_sim_1[i] <- sim_study[["hypotheses"]][[1]][["conclusion"]] == "corroborate"
  res_sim_2[i] <- sim_study[["hypotheses"]][[2]][["conclusion"]] == "corroborate"
  res_sim_3[i] <- sim_study[["hypotheses"]][[3]][["conclusion"]] == "corroborate"
  
  if(nsim == 1){study_report(sim_study, template = "postreg", filename = "study_2.html")}
}

sum(res_sim_1/nsim)
sum(res_sim_2/nsim)
sum(res_sim_3/nsim)