fix(linalg): cluster-count check precedes saturated-design NaN guard (CI review P1)

The saturated-design guard added earlier (n_eff <= k -> NaN vcov) was placed
before the cluster-robust "need >= 2 clusters" validation, so a 1-cluster
*saturated* cluster-robust request returned all-NaN instead of raising the
documented ValueError (verified: compute_robust_vcov(eye(2), 0, zeros(2)) now
raises; a non-saturated 1-cluster fit always raised). The cluster-count check
(with the zero-total-weight effective-cluster exclusion) now runs first on the
clustered path; the saturated NaN return only applies once that assumption
holds. Multi-cluster saturated designs still return NaN (no ZeroDivisionError),
and full-rank fits are unchanged.

Tests: 1-cluster saturated -> ValueError; 2-cluster saturated -> NaN vcov
(tests/test_linalg.py::TestComputeRobustVcov).

Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>

Author: igerber

diff_diff/linalg.py

@@ -2623,6 +2623,22 @@ def _compute_robust_vcov_numpy(
     # ------------------------------------------------------------------
     assert vcov_type == "hc1"
 
+    # Cluster-robust validity check FIRST: a cluster-robust request with fewer
+    # than 2 clusters is invalid and must raise the documented error — this has
+    # to precede the saturated-design guard below so a 1-cluster *saturated* fit
+    # still raises rather than being masked by the NaN return. (Mirrors the
+    # effective-cluster count computed in the cluster branch, including the
+    # zero-total-weight exclusion.)
+    if cluster_ids is not None:
+        n_clusters_check = len(np.unique(cluster_ids))
+        if weights is not None and weight_type != "fweight" and np.any(weights == 0):
+            cluster_weight_sums = pd.Series(weights).groupby(cluster_ids).sum()
+            n_clusters_check = int((cluster_weight_sums > 0).sum())
+        if n_clusters_check < 2:
+            raise ValueError(
+                f"Need at least 2 clusters for cluster-robust SEs, got {n_clusters_check}"
+            )
+
     # Saturated design (no residual degrees of freedom): both the HC1 adjustment
     # n_eff/(n_eff-k) and the CR1 adjustment (n_eff-1)/(n_eff-k) divide by
     # (n_eff - k), which is zero when the design exactly determines y. Return a

tests/test_linalg.py

@@ -531,6 +531,22 @@ def test_cluster_robust_symmetric(self, ols_data):
 
         np.testing.assert_array_almost_equal(vcov, vcov.T)
 
+    def test_cluster_count_check_precedes_saturated_guard(self):
+        """A 1-cluster cluster-robust request raises 'need >= 2 clusters' even on
+        a saturated design — the cluster-count validation must precede the
+        saturated (no residual DOF) NaN guard, not be masked by it."""
+        # 1 cluster AND saturated (n == k): must still raise, not return NaN.
+        with pytest.raises(ValueError, match="at least 2 clusters"):
+            compute_robust_vcov(np.eye(2), np.zeros(2), np.zeros(2))
+
+    def test_saturated_multi_cluster_returns_nan(self):
+        """A saturated design (no residual DOF) with >= 2 clusters returns a NaN
+        vcov rather than raising ZeroDivisionError from the CR1 adjustment."""
+        # 2 clusters, n == k == 4 (saturated 2x2 with one obs per cluster-period).
+        X = np.column_stack([np.ones(4), [0, 0, 1, 1], [0.0, 1.0, 0.0, 1.0], [0, 0, 0, 1.0]])
+        vcov = compute_robust_vcov(X, np.zeros(4), np.array([0, 0, 1, 1]))
+        assert np.all(np.isnan(vcov))
+
     def test_numerical_instability_fallback_warns(self, ols_data):
         """Test that numerical instability in Rust backend triggers warning and fallback."""
         from unittest.mock import patch
@@ -2140,9 +2156,7 @@ def test_tracker_appended_once_on_truncation_only(self):
         rng = np.random.RandomState(3)
         n = 50
         # Deficient: constant column collinear with intercept.
-        A_def = self._gram(
-            np.column_stack([np.ones(n), np.full(n, 2.0), rng.standard_normal(n)])
-        )
+        A_def = self._gram(np.column_stack([np.ones(n), np.full(n, 2.0), rng.standard_normal(n)]))
         tracker = []
         _rank_guarded_inv(A_def, tracker=tracker)
         assert len(tracker) == 1  # exactly one condition-number sample