revised bootstrapping method for when parallel == FALSE

harvard-ufds · Dec 18, 2023 · eca21cc · eca21cc
1 parent 4208a58
commit eca21cc
Show file tree

Hide file tree

Showing 3 changed files with 76 additions and 13 deletions.
diff --git a/R/unit_zi.R b/R/unit_zi.R
@@ -49,7 +49,7 @@ unit_zi <- function(samp_dat,
   }
 
   if (parallel && is(future::plan(), "sequential")) {
-    warning("In order for the internal processes to be run in parallel a `future::plan()` must be specified by the user")
+    message("In order for the internal processes to be run in parallel a `future::plan()` must be specified by the user")
     message("See <https://future.futureverse.org/reference/plan.html> for reference on how to use `future::plan()`")
   }
 
@@ -66,8 +66,6 @@ unit_zi <- function(samp_dat,
     "\\w+"
   ))
 
-
-
   original_out <- fit_zi(
     samp_dat,
     pop_dat,
@@ -76,8 +74,7 @@ unit_zi <- function(samp_dat,
     domain_level
   )
 
-
-  if (mse_est == T) {
+  if (mse_est) {
 
     zi_model_coefs <- mse_coefs(
       original_out$lmer,
@@ -232,10 +229,37 @@ unit_zi <- function(samp_dat,
       list(res_df, log_lst)
 
     }
-
-    with_progress({
-      boot_res <- boot_rep_with_progress_bar(1:B)
-    })
+
+    if (parallel) {
+
+      with_progress({
+        boot_res <- boot_rep_with_progress_bar(1:B)
+      }) 
+
+    } else {
+
+      res <- 
+        map(.x = 1:B,
+            .f = \(i) boot_rep(boot_pop_data,
+                               samp_dat,
+                               domain_level,
+                               num_plots,
+                               boot_lin_formula,
+                               boot_log_formula,
+                               boot_truth,
+                               by_domains))
+
+      res_lst <- res |>
+        map(.f = ~ .x$sqerr)
+
+      res_df <- do.call("rbind", res_lst)
+
+      log_lst <- res |>
+        map(.f = ~ .x$log)
+
+      boot_res <- list(res_df, log_lst)
+
+    }
 
     final_df <- setNames(
       aggregate(sq_error ~ domain,

diff --git a/tests/testthat/_snaps/unit_zi.md b/tests/testthat/_snaps/unit_zi.md
@@ -0,0 +1,29 @@
+# result is as expected
+
+    Code
+      result
+    Output
+      
+      Call:
+      unit_zi(samp_dat = samp, pop_dat = pop, lin_formula = lin_formula,      domain_level = "COUNTYFIPS", B = 5, mse_est = TRUE, parallel = FALSE)
+      
+      Linear Model: 
+      - Fixed effects: 
+       (Intercept)        tcc16         elev 
+      17.787586797  1.238840297 -0.004047397 
+      
+      - Random effects: 
+       Groups     Name        Std.Dev.
+       COUNTYFIPS (Intercept) 25.494  
+       Residual               68.230  
+      
+      Logistic Model: 
+      - Fixed effects: 
+       (Intercept)        tcc16         elev 
+      -4.381130205  0.104867536  0.001734101 
+      
+      - Random effects: 
+       Groups     Name        Std.Dev.
+       COUNTYFIPS (Intercept) 0.87583 
+      
+
diff --git a/tests/testthat/test-unit_zi.R b/tests/testthat/test-unit_zi.R
@@ -3,10 +3,20 @@ library(saeczi)
 data(pop)
 data(samp)
 
-lin_formula <- DRYBIO_AG_TPA_live_ADJ ~ tcc16 + tmean + tri
-
-result <- unit_zi(samp, pop, lin_formula, domain_level = "COUNTYFIPS", mse_est = TRUE, parallel = FALSE)
-
+lin_formula <- DRYBIO_AG_TPA_live_ADJ ~ tcc16 + elev
+
+set.seed(5)
+result <- unit_zi(samp,
+                  pop, 
+                  lin_formula,
+                  domain_level = "COUNTYFIPS",
+                  mse_est = TRUE,
+                  B = 5,
+                  parallel = FALSE)
+
+test_that("result is as expected", {
+  expect_snapshot(result)
+})
 
 test_that("result[[2]] is a df", {
   expect_s3_class(result[[2]], "data.frame")