further work on the dr bootstrap variance estimator

ncn-foreigners · Feb 25, 2025 · 7b3ed03 · 7b3ed03
1 parent 24264da
commit 7b3ed03
Show file tree

Hide file tree

Showing 3 changed files with 6 additions and 17 deletions.
diff --git a/R/boot_dr.R b/R/boot_dr.R
@@ -61,15 +61,9 @@ boot_dr <- function(selection,
         data_prob <- svydesign$variables[strap_rand_svy, ]
         data_prob$weight <- weights_rand_strap_svy[strap_rand_svy]
 
-        svyd_call <- as.list(svydesign$call)
-        svyd_call[[1]] <- NULL
-        svyd_call$ids <- as.formula(svyd_call$ids)
-        svyd_call$weights <- as.formula(svyd_call$weights)
-        svyd_call$strata <- as.formula(svyd_call$strata)
+        svyd_call <- svydesign$call
         svyd_call$data <- as.name("data_prob")
-
-        # Method 1: Using do.call
-        svydesign_b <- do.call(survey::svydesign, svyd_call)
+        svydesign_b <- eval(svyd_call)
 
         strap_nons <- sample.int(replace = TRUE, n = NROW(data), prob = 1 / weights)
 
@@ -222,15 +216,9 @@ boot_dr <- function(selection,
           data_prob <- svydesign$variables[strap_rand_svy, ]
           data_prob$weight <- weights_rand_strap_svy[strap_rand_svy]
 
-          svyd_call <- as.list(svydesign$call)
-          svyd_call[[1]] <- NULL
-          svyd_call$ids <- as.formula(svyd_call$ids)
-          svyd_call$weights <- as.formula(svyd_call$weights)
-          svyd_call$strata <- as.formula(svyd_call$strata)
+          svyd_call <- svydesign$call
           svyd_call$data <- as.name("data_prob")
-
-          # Method 1: Using do.call
-          svydesign_b <- do.call(survey::svydesign, svyd_call)
+          svydesign_b <- eval(svyd_call)
 
           strap_nons <- sample.int(replace = TRUE, n = NROW(data), prob = 1 / weights)
 

diff --git a/R/nonprob.R b/R/nonprob.R
@@ -28,6 +28,7 @@ nonprob <- function(data,
   method_selection <- match.arg(method_selection)
   family_outcome <- match.arg(family_outcome)
   method_outcome <- match.arg(method_outcome)
+  #stopifnot("We currently support `family_outcome` with a single entry." = NROW(family_outcome) == 1)
 
   data <- if (!is.data.frame(data)) data.frame(data) else data
   weights <- if (is.null(weights)) rep(1, nrow(data)) else weights

diff --git a/R/nonprob_dr.R b/R/nonprob_dr.R
@@ -365,7 +365,7 @@ nonprob_dr <- function(selection,
                                pop_size_fixed = pop_size_fixed)
 
     ## doubly robust estimator
-    mu_hat <- nonprobsvy:::mu_hatDR(y_hat = results_mi$output$mean,
+    mu_hat <- mu_hatDR(y_hat = results_mi$output$mean,
                        y_resid = do.call("cbind", results_mi$ys_resid),
                        weights = weights,
                        weights_nons = results_ipw$ipw_weights,