fix no covs setting

Author: CoryMcCartan

R/ei_ridge.R

@@ -131,8 +131,7 @@ ei_ridge.ei_spec <- function(x, weights, penalty=NULL, scale=TRUE, ...) {
 
     form = as.formula(paste0(
         paste0(attr(spec, "ei_y"), collapse=" + "), " ~ ",
-        paste0(attr(spec, "ei_x"), collapse=" + "), " + ",
-        paste0(attr(spec, "ei_z"), collapse=" + ")
+        paste0(c(attr(spec, "ei_x"), attr(spec, "ei_z")), collapse=" + ")
     ))
 
     bp = hardhat::new_default_formula_blueprint(

R/ei_riesz.R

@@ -75,8 +75,7 @@ ei_riesz.ei_spec <- function(x, weights, penalty, scale=TRUE, ...) {
 
     form = as.formula(paste0(
         paste0(attr(spec, "ei_y"), collapse=" + "), " ~ ",
-        paste0(attr(spec, "ei_x"), collapse=" + "), " + ",
-        paste0(attr(spec, "ei_z"), collapse=" + ")
+        paste0(c(attr(spec, "ei_x"), attr(spec, "ei_z")), collapse=" + ")
     ))
 
     bp = hardhat::new_default_formula_blueprint(

R/rr_impl.R

@@ -79,15 +79,15 @@ ridge_hat_naive <- function(
 riesz_naive <- function(xz, p, total, weights, group=1, penalty=0) {
     n_x = ncol(xz) %/% (1L + p)
     use = c(group, n_x + p*(group-1) + seq_len(p))
-    Dz = crossprod(xz[, use], total) / mean(xz[, group] * total)
+    Dz = crossprod(xz[, use, drop=FALSE], total) / mean(xz[, group] * total)
     Lambda = diag(rep(penalty, ncol(xz)))
 
     XXinv = solve(crossprod(xz, weights*xz) + Lambda)
-    xzAinv = xz %*% XXinv[, use]
+    xzAinv = xz %*% XXinv[, use, drop=FALSE]
     alpha = c(xzAinv %*% Dz)
 
     h1m = 1 - ridge_hat_naive(xz, weights, XXinv, penalty)
-    xzi = xz[, use] * total / mean(xz[, group] * total)
+    xzi = xz[, use, drop=FALSE] * total / mean(xz[, group] * total)
     loo = c((alpha - rowSums(xzAinv * xzi)) / h1m)
 
     list(alpha = alpha, loo = loo, nu2 = NA)
@@ -98,18 +98,19 @@ riesz_naive <- function(xz, p, total, weights, group=1, penalty=0) {
 riesz_svd <- function(xz, udv, p, total, weights, sqrt_w, group=1, penalty=0) {
     n_x = ncol(xz) %/% (1L + p)
     use = c(group, n_x + p*(group-1) + seq_len(p))
-    Dz = colSums(xz[, use] * total) / mean(xz[, group] * total)
+    Dz = colSums(xz[, use, drop=FALSE] * total) / mean(xz[, group] * total)
     d_pen = c(udv$d / (udv$d^2 + penalty))
 
-    xzAinv = (udv$u / sqrt_w) %*% (d_pen * t(udv$v[use, ]))
+    xzAinv = (udv$u / sqrt_w) %*% (d_pen * t(udv$v[use, , drop=FALSE]))
     alpha = xzAinv %*% Dz
 
     h1m = 1 - ridge_hat_svd(udv, penalty)
-    xzi = xz[, use] * total / mean(xz[, group] * total)
+    xzi = xz[, use, drop=FALSE] * total / mean(xz[, group] * total)
     loo = rowSums(-xzAinv * shift_cols(xzi, Dz)) / h1m
 
     # Neyman-orthogonal estimate of criterion fn
-    nu2 = sum(crossprod(Dz, udv$v[use, ])^2 * (2/(udv$d^2 + penalty) - d_pen^2)) / nrow(xz)
+    nu2 = sum(crossprod(Dz, udv$v[use, , drop=FALSE])^2 *
+                  (2/(udv$d^2 + penalty) - d_pen^2)) / nrow(xz)
 
     list(alpha = alpha, loo = loo, nu2 = nu2)
 }

tests/testthat/test-ei_est.R

@@ -50,3 +50,16 @@ test_that("Shrinkage has the expected effect", {
     expect_lt(mean(est_w$std.error), mean(est_w0$std.error))
     # expect_lt(mean(est_d$std.error), mean(est_d0$std.error))
 })
+
+test_that("Estimation methods work with no predictors", {
+    spec = ei_spec(elec_1968, vap_white:vap_other, pres_dem_hum:pres_oth,
+                   total = pres_total)
+    m = ei_ridge(spec)
+    rr = ei_riesz(spec, penalty=m$penalty)
+
+    expect_no_error({
+        est_p = ei_est(m, data=spec)
+        est_w = ei_est(rr, data=spec)
+        est_d = ei_est(m, rr, data=spec)
+    })
+})