test_that("OLS matches sens", { n <- 1000 p <- runif(n)^2 group <- sample(1:3, n, replace = TRUE) obs <- runif(n) <= p + (group == 2) * 0.1 df <- dplyr::tibble(p, group, obs) %>% dplyr::mutate(group = factor(group, levels = 1:3)) lm_fit <- lm(obs ~ p + group, data = df) %>% broom::tidy() %>% dplyr::filter(stringr::str_detect(term, "group")) %>% dplyr::pull(estimate) %>% rep(2) sens_fit <- sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200) %>% dplyr::filter(epsilon == 0) %>% dplyr::select(tidyselect::starts_with("beta")) %>% tidyr::pivot_longer(cols = tidyselect::starts_with("beta"), names_to = "term", values_to = "estimate") %>% dplyr::pull(estimate) expect_equal(lm_fit, sens_fit) }) test_that("Error if columns not in dataframe", { n <- 1000 p <- runif(n)^2 group <- sample(1:3, n, replace = TRUE) obs <- runif(n) <= p + (group == 2) * 0.1 df <- dplyr::tibble(p, group, obs) %>% dplyr::mutate(group = factor(group, levels = 1:3)) expect_error( sens(df, group_, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "must be present" ) expect_error( sens(df, group, obs_, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "must be present" ) expect_error( sens(df, group, obs, p_, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "must be present" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, lwr_col = lwr_, n_threads = 1, chunk_size = 200), "must be present" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, upr_col = upr_, n_threads = 1, chunk_size = 200), "must be present" ) # Pass non-coercible column names expect_error( sens(df, 2.0, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "Invalid column name" ) expect_error( sens(df, group, 2.0, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "Invalid column name" ) expect_error( sens(df, group, obs, 2.0, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "Invalid column name" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, lwr_col = 2.0, n_threads = 1, chunk_size = 200), "Invalid column name" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, upr_col = 2.0, n_threads = 1, chunk_size = 200), "Invalid column name" ) }) test_that("Arguments have the right type", { n <- 1000 p <- runif(n)^2 group <- sample(1:3, n, replace = TRUE) obs <- runif(n) <= p + (group == 2) * 0.1 df <- dplyr::tibble(p, group, obs) %>% dplyr::mutate(group = factor(group, levels = 1:3)) df0 <- as.list(df) expect_error( sens(df0, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "data frame" ) # Change group to a list column df1 <- df %>% dplyr::rowwise() %>% dplyr::mutate(group = list(group)) %>% dplyr::ungroup() expect_error( sens(df1, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "must be a factor" ) # Obs is not logical or coercible to logical df3 <- df %>% dplyr::mutate(obs = as.character(obs)) expect_error( sens(df3, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "must be a logical" ) # risk is not numeric df4 <- df %>% dplyr::mutate(p = as.character(p)) expect_error( sens(df4, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "must be a numeric" ) # Lwr and upr are not numeric df6 <- df %>% dplyr::mutate(lwr_ = as.character(0)) expect_error( sens(df6, group, obs, p, 1, 0.5, eta = 0.01, lwr_col = lwr_, n_threads = 1, chunk_size = 200), "must be a numeric" ) df7 <- df %>% dplyr::mutate(upr_ = as.character(1)) expect_error( sens(df7, group, obs, p, 1, 0.5, eta = 0.01, upr_col = upr_, n_threads = 1, chunk_size = 200), "must be a numeric" ) # Check that epsilon is numeric expect_error( sens(df, group, obs, p, 1, "a", eta = 0.01, n_threads = 1, chunk_size = 200), "must be a single, positive numeric" ) # Check that eta is numeric expect_error( sens(df, group, obs, p, 1, 0.5, eta = "a", n_threads = 1, chunk_size = 200), "must be a single, positive numeric" ) # Check that m is numeric expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, m = "a", n_threads = 1, chunk_size = 200), "must be a single positive integer" ) # Check that n_threads is numeric expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = "a", chunk_size = 200), "single positive integer" ) # Check that chunk size is numeric expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = "a"), "single positive integer" ) # Check that N is numeric expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, N = "a"), "single non-negative integer" ) # Check that alpha is numeric expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, alpha = "a"), "must be a single numeric" ) }) test_that("Errors when argument integrity violated", { n <- 1000 p <- runif(n)^2 group <- sample(1:3, n, replace = TRUE) obs <- runif(n) <= p + (group == 2) * 0.1 df <- dplyr::tibble(p, group, obs) %>% dplyr::mutate(group = factor(group, levels = 1:3)) # Use a basegroup that is not a vector expect_error( sens(df, group, obs, p, list(), 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "must occur in the group column" ) # Use a base group that is not in the data expect_error( sens(df, group, obs, p, 4, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "Base group" ) # Base group is not a single value expect_error( sens(df, group, obs, p, c(1, 2), 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "Base group" ) # Only one group df2 <- df %>% dplyr::filter(group == 1) expect_error( sens(df2, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "two distinct groups" ) # Risk is not between 0 and 1 df5 <- df %>% dplyr::mutate(p = p + 2) expect_error( sens(df5, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "between 0 and 1" ) # Lwr and upr are not between 0 and 1 df8 <- df %>% dplyr::mutate(lwr_ = -1) expect_error( sens(df8, group, obs, p, 1, 0.5, eta = 0.01, lwr_col = lwr_, n_threads = 1, chunk_size = 200), "Invalid lower bound" ) df9 <- df %>% dplyr::mutate(upr_ = 2) expect_error( sens(df9, group, obs, p, 1, 0.5, eta = 0.01, upr_col = upr_, n_threads = 1, chunk_size = 200), "Invalid upper bound" ) # Check that epsilon is positive, length one, and less than one. expect_error( sens(df, group, obs, p, 1, -1, eta = 0.01, n_threads = 1, chunk_size = 200), "must be a single, positive numeric" ) expect_error( sens(df, group, obs, p, 1, c(1, 2), eta = 0.01, n_threads = 1, chunk_size = 200), "must be a single, positive numeric" ) expect_warning( sens(df, group, obs, p, 1, 1.1, eta = 0.01, n_threads = 1, chunk_size = 200), "greater than 1" ) # Check that eta is strictly positive, length one, and not too large expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0, n_threads = 1, chunk_size = 200), "must be a single, positive numeric" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = c(1, 2), n_threads = 1, chunk_size = 200), "must be a single, positive numeric" ) expect_warning( sens(df, group, obs, p, 1, 0.5, eta = 2, n_threads = 1, chunk_size = 200), "large relative to", ) # Check that m is strictly positive and length one expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, m = 0, n_threads = 1, chunk_size = 200), "must be a single positive integer" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, m = c(1, 2), n_threads = 1, chunk_size = 200), "must be a single positive integer" ) # Check that n_threads is strictly positive, length one, and roundable to # integer, and not larger than the number of cores expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 0, chunk_size = 200), "single positive integer" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = c(1, 2), chunk_size = 200), "single positive integer" ) expect_warning( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1.1, chunk_size = 200), "not a whole number" ) # NOTE: To avoid using too many cores on CRAN, we only check that the number # of threads is not larger than the number of cores when not on CRAN. if (!testthat:::on_cran()) { expect_warning( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = parallel::detectCores() + 1, chunk_size = 200 ), "greater than the number of cores" ) } # Check that chunk size is strictly positive, length one, and roundable to # integer expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 0), "single positive integer" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = c(1, 2)), "single positive integer" ) expect_warning( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 1.1), "not a whole number" ) # Check that `N` is non-negative, length one, and roundable to integer expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, N = -1), "single non-negative integer" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, N = c(1, 2)), "single non-negative integer" ) expect_warning( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, N = 0.1), "not a whole number" ) # Check that alpha is between 0 and 1, length one, greater than one when multiplied by N, and # greater than 5 when multiplied by N expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, alpha = -1), "must be a single numeric" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, alpha = c(1, 2)), "must be a single numeric" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, alpha = 1.1), "between zero and one" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, alpha = 0), "between zero and one" ) expect_error( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, N = 1, alpha = 1 / 2), "too small" ) expect_warning( sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200, N = 4, alpha = 1 / 2), "small relative to" ) }) test_that("Errors when data integrity violated", { n <- 1000 p <- runif(n)^2 group <- sample(1:3, n, replace = TRUE) obs <- runif(n) <= p + (group == 2) * 0.1 # Upr and lwr formed from logit of p logit <- function(x) log(x / (1 - x)) inv_logit <- function(x) exp(x) / (1 + exp(x)) lwr <- inv_logit(logit(p) - 0.1) upr <- inv_logit(logit(p) + 0.1) df <- dplyr::tibble(p, group, obs, lwr, upr) %>% dplyr::mutate(group = factor(group, levels = 1:3)) # Error if less than two groups df2 <- df %>% dplyr::mutate(group = 1) expect_error( sens(df2, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "two distinct groups" ) # Warn if unused levels in group df3 <- df %>% dplyr::mutate(group = factor(group, levels = 1:4)) expect_warning( sens(df3, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "Some groups had no observations" ) # Error if no observations in a group or observation status df3 <- df %>% dplyr::filter(group != 1 | obs) expect_error( sens(df3, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "at least one observed and one unobserved" ) df4 <- df %>% dplyr::filter(obs == FALSE) expect_error( sens(df4, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "at least one observed and one unobserved" ) # Error if there are any missing values df5 <- df %>% dplyr::mutate(p = ifelse(runif(n) < 0.1, NA, p)) expect_error( sens(df5, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 200), "can be missing" ) # Error if not lwr <= p <= upr df6 <- df %>% dplyr::mutate(p = dplyr::if_else(runif(n) < 0.1, inv_logit(logit(lwr) - 0.1), p)) expect_error( sens(df6, group, obs, p, 1, 0.5, eta = 0.01, lwr_col = lwr, upr_col = upr, n_threads = 1, chunk_size = 200), "risks .* are invalid" ) df7 <- df %>% dplyr::mutate(p = dplyr::if_else(runif(n) < 0.1, inv_logit(logit(upr) + 0.1), p)) expect_error( sens(df7, group, obs, p, 1, 0.5, eta = 0.01, lwr_col = lwr, upr_col = upr, n_threads = 1, chunk_size = 200), "risks .* are invalid" ) # Error if the data are not sortable df8 <- df %>% dplyr::mutate( lwr = pmin(p, 1 - p), upr = pmax(p, 1 - p) ) expect_error( sens(df8, group, obs, p, 1, 0.5, eta = 0.01, lwr_col = lwr, upr_col = upr, n_threads = 1, chunk_size = 200), "sorted" ) }) test_that("Bootstrapped confidence intervals are valid", { # NOTE: This test can be moderately fragile, so we set the seed to ensure # reproucibility. set.seed(1) # This test is slow, so don't run it on CRAN skip_on_cran() # It's also probabilistic, and can occasionally fail, so don't run it on CI skip_on_ci() n <- 1000 p <- runif(n)^2 group <- sample(1:2, n, replace = TRUE) obs <- runif(n) <= p + (group == 2) * 0.1 # Upr and lwr formed from logit of p logit <- function(x) log(x / (1 - x)) inv_logit <- function(x) exp(x) / (1 + exp(x)) lwr <- inv_logit(logit(p) - 0.1) upr <- inv_logit(logit(p) + 0.1) df <- dplyr::tibble(p, group, obs, lwr, upr) %>% dplyr::mutate(group = factor(group, levels = 1:2)) # Check that the confidence intervals are valid sens_fit <- sens(df, group, obs, p, 1, 0.5, eta = 0.01, n_threads = 1, chunk_size = 100, N = 1000) expect_true(all(sens_fit$beta_min_2_02.5 <= sens_fit$beta_min_2)) expect_true(all(sens_fit$beta_min_2 <= sens_fit$beta_min_2_97.5)) expect_true(all(sens_fit$beta_max_2_02.5 <= sens_fit$beta_max_2)) expect_true(all(sens_fit$beta_max_2 <= sens_fit$beta_max_2_97.5)) })