inst/doc/branch.R

In multiverse: 'Explorable Multiverse' Data Analysis and Reports

params <-
list(EVAL = TRUE)

## ----setup, message=FALSE, warning=FALSE--------------------------------------
library(knitr)
library(dplyr)
library(tidyr)
library(purrr)
library(ggplot2)
library(multiverse)

## ---- include=FALSE-----------------------------------------------------------
M = multiverse()

## ---- chunk-setup, include=FALSE----------------------------------------------
opts_chunk$set(
  echo = TRUE,
  fig.width = 6, 
  fig.height = 4,
  eval = if (isTRUE(exists("params"))) params$EVAL else FALSE
)

## -----------------------------------------------------------------------------
M = multiverse()

## -----------------------------------------------------------------------------
x = rnorm(100, 30, 10)

## ----h-0, inside = M, echo = FALSE, engine="multiverse"-----------------------
y = ifelse(x < 30, "low", "high")

## ----h-0, inside = M, echo = FALSE, engine="multiverse"-----------------------
y = ifelse(x < 20, "low", ifelse(x > 40, "high", "medium"))

## -----------------------------------------------------------------------------
expand(M)

## -----------------------------------------------------------------------------
data("userlogs")
data.userlogs.raw = userlogs %>%
  mutate( modality = factor(modality) ) %>%
  arrange( modality )

## -----------------------------------------------------------------------------
M = multiverse()

## ----h-1, inside = M, echo = FALSE, engine="multiverse"-----------------------
df <- data.userlogs.raw %>%
    select(modality, duration) %>%
    mutate(duration = duration)

## ----h-1, inside = M, echo = FALSE, engine="multiverse"-----------------------
df <- data.userlogs.raw %>%
    select(modality, duration) %>%
    mutate(duration = log(duration))

## -----------------------------------------------------------------------------
expand(M)

## -----------------------------------------------------------------------------
expand(M)$.code

## -----------------------------------------------------------------------------
M = multiverse()

## ----h-2, inside = M, echo = FALSE, engine="multiverse"-----------------------
data("hurricane")
hurricane_data = hurricane %>%
    rename(year = Year, name = Name, dam = NDAM, death = alldeaths, female = Gender_MF,
        masfem = MasFem, category = Category, pressure = Minpressure_Updated_2014,
        wind = HighestWindSpeed)
df.filtered = hurricane_data %>%
    filter(TRUE)

## ----h-2, inside = M, echo = FALSE, engine="multiverse"-----------------------
data("hurricane")
hurricane_data = hurricane %>%
    rename(year = Year, name = Name, dam = NDAM, death = alldeaths, female = Gender_MF,
        masfem = MasFem, category = Category, pressure = Minpressure_Updated_2014,
        wind = HighestWindSpeed)
df.filtered = hurricane_data %>%
    filter(name != "Katrina")

## ----h-2, inside = M, echo = FALSE, engine="multiverse"-----------------------
data("hurricane")
hurricane_data = hurricane %>%
    rename(year = Year, name = Name, dam = NDAM, death = alldeaths, female = Gender_MF,
        masfem = MasFem, category = Category, pressure = Minpressure_Updated_2014,
        wind = HighestWindSpeed)
df.filtered = hurricane_data %>%
    filter(!(name %in% c("Katrina", "Audrey")))

## -----------------------------------------------------------------------------
expand(M)

## ----h-3, inside = M, echo = FALSE, engine="multiverse"-----------------------
df.filtered = hurricane_data %>%
    filter(TRUE)

## ----h-3, inside = M, echo = FALSE, engine="multiverse"-----------------------
df.filtered = hurricane_data %>%
    filter(name != "Katrina")

## ----h-3, inside = M, echo = FALSE, engine="multiverse"-----------------------
df.filtered = hurricane_data %>%
    filter(!(name %in% c("Katrina", "Audrey")))

## -----------------------------------------------------------------------------
M = multiverse()

## ----default-h-6, inside = M, engine="multiverse"-----------------------------
{
    x = 1
}

## ----default-h-6, inside = M, engine="multiverse"-----------------------------
{
    x = 2
}

## -----------------------------------------------------------------------------
M = multiverse()

## ----h-7, inside = M, echo = FALSE, engine="multiverse"-----------------------
duration = log(data.userlogs.raw$duration)

## ----h-7, inside = M, echo = FALSE, engine="multiverse"-----------------------
duration = identity(data.userlogs.raw$duration)

## -----------------------------------------------------------------------------
execute_multiverse(M)

expand(M) %>%
  mutate(transformed_duration = map(.results, "duration"))

## -----------------------------------------------------------------------------
M = multiverse()

## ----h-8, inside = M, echo = FALSE, engine="multiverse"-----------------------
duration_transform = log
duration = duration_transform(data.userlogs.raw$duration)

## ----h-8, inside = M, echo = FALSE, engine="multiverse"-----------------------
duration_transform = identity
duration = duration_transform(data.userlogs.raw$duration)

## -----------------------------------------------------------------------------
execute_multiverse(M)

expand(M) %>%
  mutate( transformed_duration = map(.results, "duration" ) ) %>% 
  unnest(c(transformed_duration)) %>%
  head()

## -----------------------------------------------------------------------------
library(boot)
library(broom)
samplemean <- function(x, d) {return(mean(x[d]))}

M = multiverse()

t_test_ci <- function(x, y) {
    broom::tidy(t.test(x, conf.level = y))
}

bootstrapped_ci <- function(x, y) {
    boot(data = x, statistic = samplemean, R = 5000) %>%
        broom::tidy( conf.int = TRUE, conf.method = "bca", conf.level =  y) %>%
        rename(estimate = statistic)
}

## ----h-9, inside = M, echo = FALSE, engine="multiverse"-----------------------
data_trans <- log
fit <- t_test_ci
ci <- 0.95
df <- data.userlogs.raw %>%
    mutate(duration = data_trans(duration)) %>%
    mutate(modality.f = factor(modality)) %>%
    group_by(subject, modality.f, modalityname) %>%
    summarise(duration = mean(duration), .groups = "drop") %>%
    group_by(modality.f) %>%
    rename(value = duration) %>%
    summarise(data = list(value), .groups = "keep")
df <- df %>%
    mutate(fit = map(data, ~fit(.x, ci))) %>%
    select(-data) %>%
    unnest(cols = c(fit))

## ----h-9, inside = M, echo = FALSE, engine="multiverse"-----------------------
data_trans <- log
fit <- bootstrapped_ci
ci <- 0.95
df <- data.userlogs.raw %>%
    mutate(duration = data_trans(duration)) %>%
    mutate(modality.f = factor(modality)) %>%
    group_by(subject, modality.f, modalityname) %>%
    summarise(duration = mean(duration), .groups = "drop") %>%
    group_by(modality.f) %>%
    rename(value = duration) %>%
    summarise(data = list(value), .groups = "keep")
df <- df %>%
    mutate(fit = map(data, ~fit(.x, ci))) %>%
    select(-data) %>%
    unnest(cols = c(fit))

## ----h-9, inside = M, echo = FALSE, engine="multiverse"-----------------------
data_trans <- identity
fit <- t_test_ci
ci <- 0.95
df <- data.userlogs.raw %>%
    mutate(duration = data_trans(duration)) %>%
    mutate(modality.f = factor(modality)) %>%
    group_by(subject, modality.f, modalityname) %>%
    summarise(duration = mean(duration), .groups = "drop") %>%
    group_by(modality.f) %>%
    rename(value = duration) %>%
    summarise(data = list(value), .groups = "keep")
df <- df %>%
    mutate(fit = map(data, ~fit(.x, ci))) %>%
    select(-data) %>%
    unnest(cols = c(fit))

## ----h-9, inside = M, echo = FALSE, engine="multiverse"-----------------------
data_trans <- identity
fit <- bootstrapped_ci
ci <- 0.95
df <- data.userlogs.raw %>%
    mutate(duration = data_trans(duration)) %>%
    mutate(modality.f = factor(modality)) %>%
    group_by(subject, modality.f, modalityname) %>%
    summarise(duration = mean(duration), .groups = "drop") %>%
    group_by(modality.f) %>%
    rename(value = duration) %>%
    summarise(data = list(value), .groups = "keep")
df <- df %>%
    mutate(fit = map(data, ~fit(.x, ci))) %>%
    select(-data) %>%
    unnest(cols = c(fit))

## -----------------------------------------------------------------------------
execute_multiverse(M)

## ---- fig.width = 8, fig.height = 8-------------------------------------------
expand(M) %>%
  mutate(df = map(.results, "df" ) ) %>%
  unnest( cols = c(df) )  %>%
  mutate( modality.f = recode( modality.f, 
                 "1" = "physical-touch",
                 "2" = "no-touch",
                 "3" = "virtual-prop",
                 "4" = "virtual-mouse"
  ) ) %>%
  mutate( 
    estimate = ifelse(data_trans == "log", exp(estimate), estimate),
    conf.low = ifelse(data_trans == "log", exp(conf.low), conf.low),
    conf.high = ifelse(data_trans == "log", exp(conf.high), conf.high)
  ) %>%
  arrange(ci, desc(data_trans), desc(bootstrap)) %>%
  filter( bootstrap == "bca" ) %>%
  ggplot() + 
  #geom_vline( xintercept = 0,  colour = '#979797' ) +
  geom_point( aes(x = estimate, y = modality.f)) +
  geom_errorbarh( aes(xmin = conf.low, xmax = conf.high, y = modality.f), height = 0) +
  facet_grid(. ~ data_trans, scales = "free")

## -----------------------------------------------------------------------------
M.1 = multiverse()

## -----------------------------------------------------------------------------
M.2 = multiverse()

## -----------------------------------------------------------------------------
M = multiverse()

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(1L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(2L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(3L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(4L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(5L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(6L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(7L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(8L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(9L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(10L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(11L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(12L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(13L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(14L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(15L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(16L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(17L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(18L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(19L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(20L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(21L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(22L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(23L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(24L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## ----h-10, inside = M, echo = FALSE, engine="multiverse"----------------------
set.seed(25L)
x1 <- rnorm(25)
x2 <- rnorm(25)
y <- x1 + x2 + runif(25)
m <- lm(y ~ x1 + x2)
intervals <- broom::tidy(m, conf.int = TRUE)

## -----------------------------------------------------------------------------
execute_multiverse(M)

## ---- fig.width = 8, fig.height = 8-------------------------------------------
expand(M) %>%
  mutate(df = map(.results, "intervals" ) ) %>%
  unnest( cols = c(df) ) %>%
  unnest(seed) %>%  # won't be necessary once issue #34 is fixed
  ggplot(aes(x = .universe, y = estimate, ymin = conf.low, ymax = conf.high), alpha = 0.5) +
  geom_pointrange() +
  facet_grid( . ~ term )  +
  scale_x_continuous( breaks = seq(0, 25, by = 5) ) +
  labs(x = "Universe (i.e. the seed parameter)", y = "Model Parameter Values", title = "Estimates with 95% confidence intervals")