ncn-foreigners / nonprobsvy

An R package for modern methods for non-probability surveys
https://ncn-foreigners.github.io/nonprobsvy/
Other
29 stars 4 forks source link

error when running "bootstrap" variance and pps defined as poisson_sampling in svydesign function #55

Open LukaszChrostowski opened 1 month ago

LukaszChrostowski commented 1 month ago
seed_for_sim` <- 2024-8-27
set.seed(seed_for_sim)
N <- 10000
n_A <- 500
p <- 50
alpha_vec1 <- c(-2, 1, 1, 1,1, rep(0, p-5))
alpha_vec2 <- c(0,0,0,3,3,3,3, rep(0, p-7))
beta_vec <- c(1,0,0,1,1,1,1, rep(0, p-7))
X <- cbind("(Intercept)"=1, matrix(rnorm(N*(p-1)), nrow=N, byrow=T, dimnames = list(NULL, paste0("X",1:(p-1)))))
Y <- 1 + as.numeric(X %*% beta_vec) +   rnorm(N) ## linear model
#Y <- 1 + exp(3*sin(as.numeric(X %*% beta_vec))) + X[, "X5"] + X[, "X6"] + rnorm(N) ## nonlinear model
pi_B <- plogis(as.numeric(X %*% alpha_vec1)) ## linear probability
#pi_B <- plogis(3.5 + as.numeric(log(X^2) %*% alpha_vec2) - sin(X[, "X3"] + X[, "X4"]) - X[,"X5"] + X[, "X6"]) ## nonlinear probability
flag_B <- rbinom(N, 1, prob = pi_B)
pi_A <- inclusionprobabilities(0.25 + abs(X[, "X1"]) + 0.03*abs(Y), n_A)
flag_A <- UPpoisson(pik = pi_A)
pop_data <- data.frame(pi_A, pi_B, flag_A, flag_B, Y, X[, 2:p])

X_totals <- colSums(X)
X_means <- colMeans(X[,-1])
sample_A_svy <- svydesign(ids = ~ 1, 
                          probs = ~ pi_A, 
                          pps = poisson_sampling(pop_data$pi_A[pop_data$flag_A == 1]), 
                          data = pop_data[pop_data$flag_A == 1, ])

sample_B <- pop_data[pop_data$flag_B == 1, ]

y_true <- mean(Y)
y_true

ipw_logit_sample_gee_boot <- nonprob(selection = ~ X1 + X2 + X3 + X4,
                                     target = ~ Y,
                                     data = sample_B,
                                     svydesign = sample_A_svy,
                                     method_selection = "logit",
                                     control_selection = controlSel(est_method_sel = "gee", h = 1),
                                     control_inference = controlInf(var_method = "bootstrap", num_boot = 100),
                                     verbose = TRUE)

generates an error

Error in design$fpc[, 2] : incorrect number of dimensions