Merge pull request #289 from DeclareDesign/na.rm-diagnosands

na.rm

R/declare_diagnosands.R

@@ -4,6 +4,7 @@
 #' @param keep_defaults A flag for whether to report the default diagnosands. Defaults to \code{TRUE}.
 #' @param subset A subset of the simulations data frame within which to calculate diagnosands e.g. \code{subset = p.value < .05}.
 #' @param alpha Alpha significance level. Defaults to \code{.05}.
+#' @param na.rm logical. Should simulations with missing estimates be dropped? Defaults to TRUE.
 #' @param label Label for the set of diagnosands.
 #' @param data a data.frame
 #'
@@ -13,14 +14,15 @@
 #' will contain the step label. This can be used as an additional dimension for use in diagnosis.
 #'
 #'
-#' @importFrom rlang eval_tidy quos  is_quosure quo_is_call %||%
+#' @importFrom rlang eval_tidy quos is_quosure quo_is_call %||%
 #' @rdname declare_diagnosands
 diagnosand_handler <- function(data, ..., 
                                select,
                                subtract,
                                keep_defaults = TRUE,
                                subset = NULL, 
                                alpha = 0.05,
+                               na.rm = TRUE,
                                label) {
 
   options <- quos(...)
@@ -67,16 +69,34 @@ diagnosand_handler <- function(data, ...,
   }
 
   ret <- vector("list", length(options))
+
+  if(na.rm){
+    n_row_original <- nrow(data)
+    data <- na.omit(data)
+    n_deleted <- n_row_original - nrow(data)
+
+    na_ret <- data.frame(diagnosand_label = "n_deleted",
+                         diagnosand = n_deleted,
+                         stringsAsFactors = FALSE)
+  }
+
   for (i in seq_along(options)) {
     ret[i] <- eval_tidy(options[[i]], data = data)
   }
+
   ret <- simplify2array(ret)
 
+  ret <- 
   data.frame(
     diagnosand_label = names(options),
     diagnosand = ret,
     stringsAsFactors = FALSE
   )
+
+  if(na.rm){
+    ret <- rbind_disjoint(list(ret, na_ret))
+  }
+  ret
 }
 
 
@@ -216,7 +236,17 @@ validation_fn(diagnosand_handler) <- function(ret, dots, label){
 #'
 declare_diagnosands <- make_declarations(diagnosand_handler, "diagnosand", "diagnosands")
 
-default_diagnosands <- function(data, alpha = .05){
+default_diagnosands <- function(data, alpha = .05, na.rm = TRUE){
+
+  if(na.rm){
+    n_row_original <- nrow(data)
+    data <- na.omit(data)
+    n_deleted <- n_row_original - nrow(data)
+
+    na_ret <- data.frame(diagnosand_label = "n_deleted",
+                         diagnosand = n_deleted,
+                         stringsAsFactors = FALSE)
+  }
 
   estimate <- data$estimate %||% NA
   estimand <- data$estimand %||% NA
@@ -235,6 +265,7 @@ default_diagnosands <- function(data, alpha = .05){
   type_s_rate <- mean((sign(estimate) != sign(estimand))[p.value < alpha])
   mean_estimand <- mean(estimand)
 
+  ret <- 
   data.frame(
     diagnosand_label = c(
       "bias",
@@ -257,8 +288,15 @@ default_diagnosands <- function(data, alpha = .05){
       mean_se,
       type_s_rate,
       mean_estimand
-    )
+    ),
+    stringsAsFactors = FALSE
   )
+
+  if(na.rm){
+    ret <- rbind_disjoint(list(ret, na_ret))
+  }
+
+  ret
 }
 
 

tests/testthat/test-bootstrap-diagnosands.R

@@ -30,7 +30,7 @@ test_that("test diagnosands", {
 
   diagnosis <- diagnose_design(my_design, sims = 2, bootstrap_sims = 2)
 
-  expect_equal(dim(diagnosis$diagnosands_df), c(2,23))
+  expect_equal(dim(diagnosis$diagnosands_df), c(2,25))
 
   expect_equal(dim(diagnosis$simulations_df), c(4,14))
 

tests/testthat/test-diagnosands.R

@@ -94,7 +94,7 @@ test_that("test diagnosands without estimands", {
   diagnosis <- diagnose_design(my_design2, sims = 2, diagnosands = my_dig, bootstrap_sims = FALSE)
 
 
-  expect_equal(dim(diagnosis$diagnosands_df), c(1,6))
+  expect_equal(dim(diagnosis$diagnosands_df), c(1,7))
 
 })
 
@@ -167,12 +167,30 @@ test_that("diagnosis, no estimator", {
     declare_estimand(foo = 2, bar = 3)
 
   diagnosand <- declare_diagnosands(z = mean(estimand > 0), keep_defaults = FALSE)
-
-  expect_equivalent(diagnose_design(d, diagnosands = diagnosand, sims = 5, bootstrap_sims = 5)$diagnosands_df,
-                    structure(list(design_label = structure(c(1L, 1L), .Label = "d", class = "factor"), 
-                                   estimand_label = c("bar", "foo"), z = c(1, 1), `se(z)` = c(0, 
-                                                                                              0), n_sims = c(5L, 5L)), row.names = c(NA, -2L), class = "data.frame")
-  )
+
+
+  expect_equivalent(
+    diagnose_design(
+      d,
+      diagnosands = diagnosand,
+      sims = 5,
+      bootstrap_sims = 5
+    )$diagnosands_df,
+    structure(
+      list(
+        design_label = structure(c(1L, 1L), .Label = "d", class = "factor"),
+        estimand_label = c("bar", "foo"),
+        z = c(1, 1),
+        `se(z)` = c(0,
+                    0),
+        n_deleted = c(0, 0),
+        `se(n_deleted)` = c(0, 0),
+        n_sims = c(5L,
+                   5L)
+      ),
+      class = "data.frame",
+      row.names = c(NA,-2L)
+    ))
 })
 
 test_that("Overriding join conditions",{
@@ -211,16 +229,18 @@ test_that("diagnosis, NAs if no estimand", {
   ols <- declare_estimator(extra ~ group)
   d <- declare_population(sleep) + ols
 
-sleep_ols <- structure(list(design_label = structure(1L, .Label = "d", class = "factor"), 
-                           estimator_label = "estimator", term = "group2", 
-                           bias = NA_real_, `se(bias)` = NA_real_, rmse = NA_real_, 
-                           `se(rmse)` = NA_real_, power = 0, `se(power)` = 0, coverage = NA_real_, 
-                           `se(coverage)` = NA_real_, mean_estimate = 1.58, `se(mean_estimate)` = 0, 
-                           sd_estimate = 0, `se(sd_estimate)` = 0, mean_se = 0.849091017238762, 
-                           `se(mean_se)` = 0, type_s_rate = NaN, `se(type_s_rate)` = NA_real_, 
-                           mean_estimand = NA_real_, `se(mean_estimand)` = NA_real_, 
-                           n_sims = 4L), row.names = c(NA, -1L), class = "data.frame")
-
+  sleep_ols <-
+    structure(list(design_label = structure(1L, .Label = "d", class = "factor"), 
+                   estimator_label = "estimator", term = "group2", bias = NA_real_, 
+                   `se(bias)` = NA_real_, rmse = NA_real_, `se(rmse)` = NA_real_, 
+                   power = 0, `se(power)` = 0, coverage = NA_real_, `se(coverage)` = NA_real_, 
+                   mean_estimate = 1.58, `se(mean_estimate)` = 0, sd_estimate = 0, 
+                   `se(sd_estimate)` = 0, mean_se = 0.849091017238762, `se(mean_se)` = 0, 
+                   type_s_rate = NaN, `se(type_s_rate)` = NA_real_, mean_estimand = NA_real_, 
+                   `se(mean_estimand)` = NA_real_, n_deleted = 0, `se(n_deleted)` = 0, 
+                   n_sims = 4L), class = "data.frame", row.names = c(NA, -1L
+                   ))
+
 expect_equivalent(diagnose_design(d, sims = 4, bootstrap_sims = 5)$diagnosands_df, sleep_ols)
 
 })
@@ -230,14 +250,15 @@ test_that("diagnosis, NAs if no estimand", {
   d <- declare_population(sleep) + mu
 
   sleep_ols <- structure(list(design_label = structure(1L, .Label = "d", class = "factor"), 
-                              estimand_label = "estimand", bias = NA_real_, `se(bias)` = NA_real_, 
-                              rmse = NA_real_, `se(rmse)` = NA_real_, power = NA_real_, 
-                              `se(power)` = NA_real_, coverage = NA_real_, `se(coverage)` = NA_real_, 
-                              mean_estimate = NA_real_, `se(mean_estimate)` = NA_real_, 
-                              sd_estimate = NA_real_, `se(sd_estimate)` = NA_real_, mean_se = NA_real_, 
-                              `se(mean_se)` = NA_real_, type_s_rate = NA_real_, `se(type_s_rate)` = NA_real_, 
-                              mean_estimand = 1.54, `se(mean_estimand)` = 0, n_sims = 4L), row.names = c(NA, 
-                                                                                                         -1L), class = "data.frame")
+                                    estimand_label = "estimand", bias = NA_real_, `se(bias)` = NA_real_, 
+                                    rmse = NA_real_, `se(rmse)` = NA_real_, power = NA_real_, 
+                                    `se(power)` = NA_real_, coverage = NA_real_, `se(coverage)` = NA_real_, 
+                                    mean_estimate = NA_real_, `se(mean_estimate)` = NA_real_, 
+                                    sd_estimate = NA_real_, `se(sd_estimate)` = NA_real_, mean_se = NA_real_, 
+                                    `se(mean_se)` = NA_real_, type_s_rate = NA_real_, `se(type_s_rate)` = NA_real_, 
+                                    mean_estimand = 1.54, `se(mean_estimand)` = 0, n_deleted = 0, 
+                                    `se(n_deleted)` = 0, n_sims = 4L), class = "data.frame", row.names = c(NA, 
+                                                                                                           -1L))
 
   expect_equivalent(diagnose_design(d, sims = 4)$diagnosands_df, sleep_ols)
 })
@@ -254,21 +275,21 @@ test_that("select, subtract, add diagnosands",{
   dx <- diagnose_design(my_design, diagnosands = my_diags, sims = 4, bootstrap_sims = FALSE)
 
   expect_true(all(dx$diagnosand_names %in% c("new_diag", "bias", "rmse", "power", "coverage", "mean_estimate", 
-                                 "sd_estimate", "mean_se", "type_s_rate", "mean_estimand")))
+                                 "sd_estimate", "mean_se", "type_s_rate", "mean_estimand", "n_deleted")))
 
   # add a diagnosand, dont keep
   my_diags <- declare_diagnosands(new_diag = mean(estimate - estimand), keep_defaults = FALSE)
   dx <- diagnose_design(my_design, diagnosands = my_diags, sims = 4, bootstrap_sims = FALSE)
-  expect_true(all(dx$diagnosand_names %in% c("new_diag")))
+  expect_true(all(dx$diagnosand_names %in% c("new_diag", "n_deleted")))
 
   # select
   my_diags <- declare_diagnosands(select = bias, keep_defaults = TRUE)
   dx <- diagnose_design(my_design, diagnosands = my_diags, sims = 4, bootstrap_sims = FALSE)
-  expect_true(all(dx$diagnosand_names %in% c("bias")))
+  expect_true(all(dx$diagnosand_names %in% c("bias", "n_deleted")))
 
   my_diags <- declare_diagnosands(select = c(bias, rmse))
   dx <- diagnose_design(my_design, diagnosands = my_diags, sims = 4, bootstrap_sims = FALSE)
-  expect_true(all(dx$diagnosand_names %in% c("bias", "rmse")))
+  expect_true(all(dx$diagnosand_names %in% c("bias", "rmse", "n_deleted")))
 
 
   # subtract
@@ -309,3 +330,54 @@ test_that("declare time errors",{
 
 })
 
+
+
+test_that("missingness",{
+
+  my_population <- declare_population(N = 50, noise = rnorm(N))
+  fixed_pop <-  my_population()
+  my_pop <- declare_population(fixed_pop)
+
+  my_odd_estimator <- function(data) {
+    estimate = rnorm(1)
+    if(estimate > 0){
+      estimate <- NA
+    }
+    data.frame(estimate = estimate)
+  }
+  # my_odd_estimator(my_pop)
+  estimator <- declare_estimator(handler = my_odd_estimator)
+  des <- my_pop + estimator
+
+  dx <- diagnose_design(des, sims = 50, bootstrap_sims = FALSE)
+  expect_equal(
+    names(dx$diagnosands_df),
+    c(
+      "design_label",
+      "bias",
+      "rmse",
+      "power",
+      "coverage",
+      "mean_estimate",
+      "sd_estimate",
+      "mean_se",
+      "type_s_rate",
+      "mean_estimand",
+      "n_deleted",
+      "n_sims"
+    )
+  )
+
+  diags <- declare_diagnosands(select = c(mean_estimate), na.rm = TRUE)
+  dx <- diagnose_design(des, sims = 50, diagnosands = diags, bootstrap_sims = FALSE)
+  expect_equal(names(dx$diagnosands_df), c("design_label", "mean_estimate", "n_deleted", "n_sims"))
+
+  diags <- declare_diagnosands(select = c(mean_estimate), na.rm = FALSE)
+  dx <- diagnose_design(des, sims = 50, diagnosands = diags, bootstrap_sims = FALSE)
+  expect_equal(names(dx$diagnosands_df), c("design_label", "mean_estimate", "n_sims"))
+
+
+})
+
+
+

tests/testthat/test-diagnose-design.R

@@ -114,8 +114,10 @@ test_that("default diagnosands work", {
     design_1 = design_1,
     sims = 2
   )
-  expect_equal(names(diag$diagnosands_df), c("design_label", "estimand_label", "estimator_label", "term", 
-                                               "med_bias", "se(med_bias)", "n_sims"))
+  expect_equal(names(diag$diagnosands_df), 
+               c("design_label", "estimand_label", "estimator_label", "term", 
+                 "med_bias", "se(med_bias)", "n_deleted", "se(n_deleted)", "n_sims"
+               ))
 
   # w/ set diagnosands each manually
 
@@ -133,7 +135,8 @@ test_that("default diagnosands work", {
   )
 
   expect_equal(names(diag$diagnosands_df), c("design_label", "estimand_label", "estimator_label", "term", 
-                                            "my_bias", "se(my_bias)", "my_power", "se(my_power)", "n_sims"))
+                                             "my_bias", "se(my_bias)", "n_deleted", "se(n_deleted)", "my_power", 
+                                             "se(my_power)", "n_sims"))
 
   # w/ none set
 
@@ -143,12 +146,13 @@ test_that("default diagnosands work", {
     sims = 2
   )
 
-  expect_equal(names(diag$diagnosands_df), c("design_label", "estimand_label", "estimator_label", "term", 
-                                             "bias", "se(bias)", "rmse", "se(rmse)", "power", "se(power)", 
-                                             "coverage", "se(coverage)", "mean_estimate", "se(mean_estimate)", 
-                                             "sd_estimate", "se(sd_estimate)", "mean_se", "se(mean_se)", "type_s_rate", 
-                                             "se(type_s_rate)", "mean_estimand", "se(mean_estimand)", "n_sims"
-  ))
+  expect_equal(names(diag$diagnosands_df), 
+               c("design_label", "estimand_label", "estimator_label", "term", 
+                 "bias", "se(bias)", "rmse", "se(rmse)", "power", "se(power)", 
+                 "coverage", "se(coverage)", "mean_estimate", "se(mean_estimate)", 
+                 "sd_estimate", "se(sd_estimate)", "mean_se", "se(mean_se)", "type_s_rate", 
+                 "se(type_s_rate)", "mean_estimand", "se(mean_estimand)", "n_deleted", 
+                 "se(n_deleted)", "n_sims"))
 
   # w/ none set and override
 
@@ -159,15 +163,14 @@ test_that("default diagnosands work", {
     sims = 2
   )
 
-  expect_equal(names(diag$diagnosands_df), c("design_label", "estimand_label", "estimator_label", "term", 
-                                             "med_bias", "se(med_bias)", "n_sims"))
-
+  expect_equal(names(diag$diagnosands_df), 
+               c("design_label", "estimand_label", "estimator_label", "term", 
+                 "med_bias", "se(med_bias)", "n_deleted", "se(n_deleted)", "n_sims"
+               ))
 
 
   # w/ mix of set and unset
 
-
-
   # // expand_designs list
 
   # w/ diagnosands set
@@ -176,16 +179,17 @@ test_that("default diagnosands work", {
 
   diag <- diagnose_design(designs, sims = 5, bootstrap_sims = 0)
 
-  expect_equal(names(diag$diagnosands_df), c("design_label", "N", "estimand_label", "estimator_label", "term", 
-                                             "med_bias", "n_sims"))
+  expect_equal(names(diag$diagnosands_df), 
+               c("design_label", "N", "estimand_label", "estimator_label", "term", 
+                 "med_bias", "n_deleted", "n_sims"))
 
   # w mix of diagnosands set
 
   attr(designs[[1]], "diagnosands") <- NULL
 
   diag <- diagnose_design(designs, sims = 5, bootstrap_sims = 0)
 
-  expect_equal(ncol(diag$diagnosands_df), 16)
+  expect_equal(ncol(diag$diagnosands_df), 17)
 
   # // simulation df
   sims <- set_diagnosands(simulate_design(designs, sims = 5), declare_diagnosands(med_bias = median(estimate - estimand)))

tests/testthat/test-factorial.R

@@ -37,6 +37,6 @@ test_that("Factorial", {
 
   expect_equal(diagnosis %>% get_simulations %>% dim, c(2, 14))
 
-  expect_equal(diagnosis %>%  get_diagnosands %>% dim, c(1,14))
+  expect_equal(diagnosis %>%  get_diagnosands %>% dim, c(1,15))
 
 })

tests/testthat/test-multiple-coefficients.R

@@ -27,6 +27,6 @@ test_that("Multiple Coefficients", {
 
   expect_equal(diagnosis %>% get_simulations %>% dim, c(4, 12))
 
-  expect_equal(diagnosis %>% get_diagnosands %>% dim, c(2,14))
+  expect_equal(diagnosis %>% get_diagnosands %>% dim, c(2,15))
 
 })