Use model housing eligibility for take-up (#996)

* Use model housing eligibility for takeup

* Add housing takeup changelog

* Use cloned counties for housing takeup eligibility

* Pin policyengine-us housing eligibility release

Author: MaxGhenis

changelog.d/housing_takeup.fixed.md

@@ -0,0 +1 @@
+Use model housing assistance eligibility for housing take-up draws.

policyengine_us_data/build_outputs/us_augmentations.py

@@ -6,6 +6,9 @@
 from typing import Any, Callable, Mapping
 
 import numpy as np
+from policyengine_us.variables.gov.hud.is_eligible_for_housing_assistance import (
+    housing_assistance_eligibility_from_income_limits,
+)
 
 from policyengine_us_data.calibration.block_assignment import (
     derive_geography_from_blocks,
@@ -424,6 +427,11 @@ def _build_takeup_results(
             "spm_unit": len(subentity_source_indices["spm_unit"]),
         }
         reported_anchors = _build_reported_takeup_anchors(data, time_period)
+        eligibility_masks = self._build_eligibility_masks(
+            data=data,
+            context=context,
+            subentity_source_indices=subentity_source_indices,
+        )
         voluntary_filing_inputs = self._build_voluntary_filing_inputs(
             context=context,
             tax_unit_source_indices=subentity_source_indices["tax_unit"],
@@ -446,6 +454,7 @@ def _build_takeup_results(
                 else None
             ),
             reported_anchors=reported_anchors,
+            eligibility_masks=eligibility_masks,
             voluntary_filing_inputs=voluntary_filing_inputs,
         )
 
@@ -495,6 +504,72 @@ def _build_voluntary_filing_inputs(
             )[tax_unit_source_indices],
         }
 
+    def _build_eligibility_masks(
+        self,
+        *,
+        data: PayloadData,
+        context: PayloadBuildContext,
+        subentity_source_indices: Mapping[str, np.ndarray],
+    ) -> dict[str, np.ndarray]:
+        time_period = context.time_period
+        spm_unit_source_indices = subentity_source_indices["spm_unit"]
+        spm_unit_household_indices = (
+            context.reindexed.subentity_household_clone_indices["spm_unit"].astype(
+                np.int64
+            )
+        )
+        household_county_fips = _required_period_array(
+            data,
+            "county_fips",
+            time_period,
+            "US take-up housing assistance eligibility requires county_fips "
+            "from USGeographyPostProcessor",
+        )
+        if (
+            "receives_housing_assistance" in data
+            and time_period in data["receives_housing_assistance"]
+        ):
+            receives_housing_assistance = data["receives_housing_assistance"][
+                time_period
+            ].astype(bool)
+        else:
+            receives_housing_assistance = calculate_variable_values(
+                context.simulation,
+                "receives_housing_assistance",
+                period=time_period,
+                map_to="spm_unit",
+            )[spm_unit_source_indices].astype(bool)
+
+        return {
+            "takes_up_housing_assistance_if_eligible": (
+                housing_assistance_eligibility_from_income_limits(
+                    county_fips=np.asarray(household_county_fips)[
+                        spm_unit_household_indices
+                    ],
+                    annual_income=calculate_variable_values(
+                        context.simulation,
+                        "hud_annual_income",
+                        period=time_period,
+                        map_to="spm_unit",
+                    )[spm_unit_source_indices],
+                    spm_unit_size=calculate_variable_values(
+                        context.simulation,
+                        "spm_unit_size",
+                        period=time_period,
+                        map_to="spm_unit",
+                    )[spm_unit_source_indices],
+                    spm_unit_tenure_type=calculate_variable_values(
+                        context.simulation,
+                        "spm_unit_tenure_type",
+                        period=time_period,
+                        map_to="spm_unit",
+                    )[spm_unit_source_indices],
+                    receives_housing_assistance=receives_housing_assistance,
+                    year=time_period,
+                ).astype(bool)
+            )
+        }
+
 
 def default_us_postprocessors() -> tuple[
     USEntityPostProcessor | USGeographyPostProcessor | USTakeupPostProcessor, ...

policyengine_us_data/calibration/entity_clone.py

@@ -16,6 +16,9 @@
 
 import h5py
 import numpy as np
+from policyengine_us.variables.gov.hud.is_eligible_for_housing_assistance import (
+    housing_assistance_eligibility_from_income_limits,
+)
 
 from policyengine_us_data.calibration.block_assignment import (
     derive_geography_from_blocks,
@@ -139,6 +142,13 @@ def _build_reported_takeup_anchors(data: dict, time_period: int) -> dict:
         reported_anchors["takes_up_medicaid_if_eligible"] = data[
             "has_medicaid_health_coverage_at_interview"
         ][time_period].astype(bool)
+    if (
+        "receives_housing_assistance" in data
+        and time_period in data["receives_housing_assistance"]
+    ):
+        reported_anchors["takes_up_housing_assistance_if_eligible"] = data[
+            "receives_housing_assistance"
+        ][time_period].astype(bool)
     return reported_anchors
 
 
@@ -423,6 +433,47 @@ def materialize_clone_household_chunk(
                 map_to="tax_unit",
             ).values[entity_clone_idx["tax_unit"]],
         }
+        if (
+            "receives_housing_assistance" in data
+            and time_period in data["receives_housing_assistance"]
+        ):
+            receives_housing_assistance = data["receives_housing_assistance"][
+                time_period
+            ].astype(bool)
+        else:
+            receives_housing_assistance = (
+                sim.calculate(
+                    "receives_housing_assistance",
+                    time_period,
+                    map_to="spm_unit",
+                )
+                .values[entity_clone_idx["spm_unit"]]
+                .astype(bool)
+            )
+        eligibility_masks = {
+            "takes_up_housing_assistance_if_eligible": (
+                housing_assistance_eligibility_from_income_limits(
+                    county_fips=clone_geo["county_fips"][entity_hh_indices["spm_unit"]],
+                    annual_income=sim.calculate(
+                        "hud_annual_income",
+                        time_period,
+                        map_to="spm_unit",
+                    ).values[entity_clone_idx["spm_unit"]],
+                    spm_unit_size=sim.calculate(
+                        "spm_unit_size",
+                        time_period,
+                        map_to="spm_unit",
+                    ).values[entity_clone_idx["spm_unit"]],
+                    spm_unit_tenure_type=sim.calculate(
+                        "spm_unit_tenure_type",
+                        time_period,
+                        map_to="spm_unit",
+                    ).values[entity_clone_idx["spm_unit"]],
+                    receives_housing_assistance=receives_housing_assistance,
+                    year=time_period,
+                ).astype(bool)
+            )
+        }
         takeup_results = apply_block_takeup_to_arrays(
             hh_blocks=active_blocks,
             hh_state_fips=clone_geo["state_fips"].astype(np.int32),
@@ -433,6 +484,7 @@ def materialize_clone_household_chunk(
             time_period=time_period,
             takeup_filter=takeup_filter,
             reported_anchors=reported_anchors,
+            eligibility_masks=eligibility_masks,
             voluntary_filing_inputs=voluntary_filing_inputs,
         )
         for variable, values in takeup_results.items():

policyengine_us_data/datasets/cps/cps.py

@@ -27,7 +27,6 @@
 from policyengine_us_data.datasets.cps.takeup import (
     align_reported_ssi_disability,
     prioritize_reported_recipients,
-    very_low_income_renter_mask,
 )
 from policyengine_us_data.datasets.org import (
     ORG_BOOL_VARIABLES,
@@ -579,15 +578,14 @@ def add_takeup(self):
         ),
         dtype=bool,
     )
-    very_low_income_renter = very_low_income_renter_mask(
-        baseline.calculate("hud_income_level").values,
-        baseline.calculate("spm_unit_tenure_type").values,
-    )
+    housing_assistance_eligible = baseline.calculate(
+        "is_eligible_for_housing_assistance"
+    ).values
     data["takes_up_housing_assistance_if_eligible"] = prioritize_reported_recipients(
         reported_housing_assistance,
         housing_assistance_rate,
         rng.random(n_spm_units),
-        eligible_mask=very_low_income_renter,
+        eligible_mask=housing_assistance_eligible,
     )
 
     # WIC: resolve draws to bools using category-specific rates

policyengine_us_data/datasets/cps/takeup.py

@@ -1,35 +1,12 @@
 import numpy as np
 
-VERY_LOW_INCOME_LEVELS = {"ESPECIALLY_LOW", "VERY_LOW"}
-
 
 def _validate_same_shape(*arrays: np.ndarray) -> None:
     shapes = {np.asarray(array).shape for array in arrays}
     if len(shapes) != 1:
         raise ValueError("All arrays must have the same shape")
 
 
-def _enum_names(values: np.ndarray) -> np.ndarray:
-    return np.asarray(
-        [
-            getattr(value, "name", str(value))
-            for value in np.asarray(values, dtype=object)
-        ]
-    )
-
-
-def very_low_income_renter_mask(
-    income_level: np.ndarray,
-    tenure_type: np.ndarray,
-) -> np.ndarray:
-    income_level = _enum_names(income_level)
-    tenure_type = _enum_names(tenure_type)
-    _validate_same_shape(income_level, tenure_type)
-    return np.isin(income_level, list(VERY_LOW_INCOME_LEVELS)) & (
-        tenure_type == "RENTER"
-    )
-
-
 def prioritize_reported_recipients(
     reported_receipt: np.ndarray,
     target_rate: float,

policyengine_us_data/utils/takeup.py

@@ -190,6 +190,7 @@ def assign_takeup_with_reported_anchors(
     rates,
     reported_mask: Optional[np.ndarray] = None,
     group_keys: Optional[np.ndarray] = None,
+    eligible_mask: Optional[np.ndarray] = None,
 ) -> np.ndarray:
     """Apply the SSI/SNAP-style reported-first takeup pattern.
 
@@ -206,22 +207,30 @@ def assign_takeup_with_reported_anchors(
         if len(rates_arr) != len(draws):
             raise ValueError("rates and draws must align")
 
+    if eligible_mask is None:
+        eligible_mask = np.ones(len(draws), dtype=bool)
+    else:
+        eligible_mask = np.asarray(eligible_mask, dtype=bool)
+        if len(eligible_mask) != len(draws):
+            raise ValueError("eligible_mask and draws must align")
+
     baseline = draws < rates_arr
     if reported_mask is None:
-        return baseline
+        return eligible_mask & baseline
 
     reported_mask = np.asarray(reported_mask, dtype=bool)
     if len(reported_mask) != len(draws):
         raise ValueError("reported_mask and draws must align")
 
+    eligible_mask = eligible_mask | reported_mask
     result = reported_mask.copy()
 
     if group_keys is None:
         unique_rates = np.unique(rates_arr)
         if len(unique_rates) != 1:
             raise ValueError("group_keys required when rates vary by entity")
-        target_count = int(unique_rates[0] * len(draws))
-        non_reporters = ~reported_mask
+        target_count = int(unique_rates[0] * int(eligible_mask.sum()))
+        non_reporters = eligible_mask & ~reported_mask
         remaining_needed = max(0, target_count - int(reported_mask.sum()))
         adjusted_rate = (
             remaining_needed / int(non_reporters.sum()) if non_reporters.any() else 0
@@ -238,10 +247,11 @@ def assign_takeup_with_reported_anchors(
         group_rates = np.unique(rates_arr[group_mask])
         if len(group_rates) != 1:
             raise ValueError("Each takeup group must have a single rate")
-        target_count = int(group_rates[0] * int(group_mask.sum()))
+        group_eligible = group_mask & eligible_mask
+        target_count = int(group_rates[0] * int(group_eligible.sum()))
         group_reported = reported_mask[group_mask]
         remaining_needed = max(0, target_count - int(group_reported.sum()))
-        group_non_reporters = group_mask & ~reported_mask
+        group_non_reporters = group_eligible & ~reported_mask
         adjusted_rate = (
             remaining_needed / int(group_non_reporters.sum())
             if group_non_reporters.any()
@@ -423,6 +433,7 @@ def compute_block_takeup_for_entities(
     entity_hh_ids: np.ndarray = None,
     entity_clone_ids: np.ndarray = None,
     reported_mask: Optional[np.ndarray] = None,
+    eligible_mask: Optional[np.ndarray] = None,
 ) -> np.ndarray:
     """Compute boolean takeup via block-level seeded draws."""
     draws = compute_block_takeup_draws_for_entities(
@@ -448,6 +459,7 @@ def compute_block_takeup_for_entities(
         rates,
         reported_mask=reported_mask,
         group_keys=group_keys,
+        eligible_mask=eligible_mask,
     )
 
 
@@ -660,6 +672,7 @@ def apply_block_takeup_to_arrays(
     takeup_filter: List[str] = None,
     precomputed_rates: Optional[Dict[str, Any]] = None,
     reported_anchors: Optional[Dict[str, np.ndarray]] = None,
+    eligibility_masks: Optional[Dict[str, np.ndarray]] = None,
     voluntary_filing_inputs: Optional[Dict[str, np.ndarray]] = None,
 ) -> Dict[str, np.ndarray]:
     """Compute takeup draws from raw arrays.
@@ -686,13 +699,18 @@ def apply_block_takeup_to_arrays(
         precomputed_rates: Optional {rate_key: rate_or_dict} cache.
             When provided, skips ``load_take_up_rate`` calls and
             uses cached values instead.
+        reported_anchors: Optional {takeup variable: bool array}; reported
+            recipients are always assigned take-up.
+        eligibility_masks: Optional {takeup variable: bool array}; non-reported
+            take-up is drawn only from the matching eligible entity rows.
 
     Returns:
         {variable_name: bool_array} for each takeup variable.
     """
     filter_set = set(takeup_filter) if takeup_filter is not None else None
     result = {}
     reported_anchors = reported_anchors or {}
+    eligibility_masks = eligibility_masks or {}
 
     for spec in SIMPLE_TAKEUP_VARS:
         var_name = spec["variable"]
@@ -716,6 +734,9 @@ def apply_block_takeup_to_arrays(
         reported_mask = reported_anchors.get(var_name)
         if reported_mask is not None and len(reported_mask) != n_ent:
             raise ValueError(f"reported anchor for {var_name} has wrong length")
+        eligible_mask = eligibility_masks.get(var_name)
+        if eligible_mask is not None and len(eligible_mask) != n_ent:
+            raise ValueError(f"eligibility mask for {var_name} has wrong length")
         if var_name == "would_file_taxes_voluntarily":
             if voluntary_filing_inputs is None:
                 raise ValueError(
@@ -739,6 +760,7 @@ def apply_block_takeup_to_arrays(
                 ent_hh_ids,
                 ent_clone_indices,
                 reported_mask=reported_mask,
+                eligible_mask=eligible_mask,
             )
         result[var_name] = bools
 

pyproject.toml

@@ -22,7 +22,7 @@ classifiers = [
     "Programming Language :: Python :: 3.14",
 ]
 dependencies = [
-    "policyengine-us==1.693.1",
+    "policyengine-us==1.693.2",
     # policyengine-core 3.26.1 is the current 3.26.x runtime and includes the fix for
     # PolicyEngine/policyengine-core#482 (user-set ETERNITY inputs lost
     # after _invalidate_all_caches) and is required by policyengine-us 1.682.1+.

tests/unit/build_outputs/test_us_augmentations.py

@@ -36,6 +36,10 @@ def calculate(self, variable, period=None, map_to=None):
             ("tax_unit_child_dependents", "tax_unit"): np.array([1, 2]),
             ("employment_income", "person"): np.array([10, 20, 30]),
             ("age_head", "tax_unit"): np.array([40, 50]),
+            ("hud_annual_income", "spm_unit"): np.array([200_000, 50_000]),
+            ("spm_unit_size", "spm_unit"): np.array([4, 4]),
+            ("spm_unit_tenure_type", "spm_unit"): np.array(["RENTER", "RENTER"]),
+            ("receives_housing_assistance", "spm_unit"): np.array([False, False]),
         }[(variable, map_to)]
         return _Calculation(values)
 
@@ -307,6 +311,10 @@ def fake_sum_person_values(person_values, person_tax_unit_ids, tax_unit_ids):
     )
     assert seen["entity_counts"] == {"person": 3, "tax_unit": 2, "spm_unit": 2}
     assert seen["takeup_filter"] == ["takes_up_snap_if_eligible"]
+    np.testing.assert_array_equal(
+        seen["eligibility_masks"]["takes_up_housing_assistance_if_eligible"],
+        np.array([True, False]),
+    )
     np.testing.assert_array_equal(
         seen["voluntary_filing_inputs"]["tax_unit_child_dependents"],
         np.array([2, 1]),