Explicit error terms production code

.github/workflows/core_tests.yml

@@ -141,6 +141,7 @@ jobs:
       matrix:
         region:
           - prototype_mtc
+          - prototype_arc
           - placeholder_psrc
           - prototype_marin
           - prototype_mtc_extended

activitysim/abm/models/disaggregate_accessibility.py

@@ -16,12 +16,13 @@
 from activitysim.abm.tables import shadow_pricing
 from activitysim.core import estimation, los, tracing, util, workflow
 from activitysim.core.configuration.base import (
+    ComputeSettings,
     PreprocessorSettings,
     PydanticReadable,
-    ComputeSettings,
 )
 from activitysim.core.configuration.logit import TourLocationComponentSettings
 from activitysim.core.expressions import assign_columns
+from activitysim.core.interaction_sample import _resolve_sample_method
 
 logger = logging.getLogger(__name__)
 
@@ -760,6 +761,21 @@ def get_disaggregate_logsums(
         state, "disaggregate_accessibility.yaml"
     )
 
+    if (
+        _resolve_sample_method(
+            state, getattr(disagg_model_settings, "compute_settings", None)
+        )
+        == "poisson"
+    ):
+        logger.warning(
+            "Using Poisson sampling method for disaggregate accessibility calculations. Currently the results will"
+            + " differ from those obtained with monte-carlo or eet sampling by a constant shift of"
+            + f" log({disagg_model_settings.DESTINATION_SAMPLE_SIZE}) if you are using the common correction factor"
+            + " `log(pick_count / prob)` in location choice specs. The results of the Poisson method are unbiased,"
+            + " i.e., they agree with the results obtained with a full destination sample, unlike those for"
+            + " monte-carlo or eet sampling."
+        )
+
     for model_name in [
         "workplace_location",
         "school_location",

activitysim/abm/models/joint_tour_participation.py

@@ -20,8 +20,8 @@
 )
 from activitysim.core.configuration.base import ComputeSettings, PreprocessorSettings
 from activitysim.core.configuration.logit import LogitComponentSettings
-from activitysim.core.util import assign_in_place, reindex
 from activitysim.core.exceptions import InvalidTravelError
+from activitysim.core.util import assign_in_place, reindex
 
 logger = logging.getLogger(__name__)
 
@@ -131,6 +131,7 @@ def participants_chooser(
     choosers: pd.DataFrame,
     spec: pd.DataFrame,
     trace_label: str,
+    nest_spec: Optional[dict, LogitNestSpec] = None,
 ) -> tuple[pd.Series, pd.Series]:
     """
     custom alternative to logit.make_choices for simulate.simple_simulate
@@ -250,14 +251,21 @@ def participants_chooser(
                     f"{num_tours_remaining} tours could not be satisfied after {iter} iterations"
                 )
 
-        choice_function = (
-            logit.make_choices_utility_based
-            if state.settings.use_explicit_error_terms
-            else logit.make_choices
-        )
-        choices, rands = choice_function(
-            state, probs_or_utils, trace_label=trace_label, trace_choosers=choosers
-        )
+        if state.settings.use_explicit_error_terms:
+            choices, rands = logit.make_choices_utility_based(
+                state,
+                probs_or_utils,
+                trace_label=trace_label,
+                trace_choosers=choosers,
+                nest_spec=nest_spec,
+            )
+        else:
+            choices, rands = logit.make_choices(
+                state,
+                probs_or_utils,
+                trace_label=trace_label,
+                trace_choosers=choosers,
+            )
         participate = choices == PARTICIPATE_CHOICE
 
         # satisfaction indexed by tour_id
@@ -429,14 +437,6 @@ def joint_tour_participation(
         if i not in model_settings.compute_settings.protect_columns:
             model_settings.compute_settings.protect_columns.append(i)
 
-    # TODO EET: this is related to the difference in nested logit and logit choice as per comment in
-    # make_choices_utility_based. As soon as alt_order_array is removed from arguments to
-    # make_choices_explicit_error_term_nl this guard can be removed
-    if state.settings.use_explicit_error_terms:
-        assert (
-            nest_spec is None
-        ), "Nested logit model custom chooser for EET requires name_mapping, currently not implemented in jtp"
-
     custom_chooser = participants_chooser
 
     choices = simulate.simple_simulate_by_chunk_id(

activitysim/abm/models/location_choice.py

@@ -15,10 +15,14 @@
     TourLocationComponentSettings,
     TourModeComponentSettings,
 )
-from activitysim.core.interaction_sample import interaction_sample
+from activitysim.core.exceptions import DuplicateWorkflowTableError
+from activitysim.core.interaction_sample import (
+    _resolve_sample_method,
+    interaction_sample,
+)
 from activitysim.core.interaction_sample_simulate import interaction_sample_simulate
+from activitysim.core.logit import AltsContext
 from activitysim.core.util import reindex
-from activitysim.core.exceptions import DuplicateWorkflowTableError
 
 """
 The school/workplace location model predicts the zones in which various people will
@@ -116,6 +120,8 @@ def _location_sample(
     chunk_tag,
     trace_label,
     zone_layer=None,
+    stable_alt_positions=None,
+    n_total_alts=None,
 ):
     """
     select a sample of alternative locations.
@@ -211,6 +217,8 @@ def _location_sample(
         chunk_tag=chunk_tag,
         trace_label=trace_label,
         zone_layer=zone_layer,
+        stable_alt_positions=stable_alt_positions,
+        n_total_alts=n_total_alts,
         explicit_chunk_size=model_settings.explicit_chunk,
         compute_settings=model_settings.compute_settings.subcomponent_settings(
             "sample"
@@ -226,6 +234,7 @@ def location_sample(
     persons_merged,
     network_los,
     dest_size_terms,
+    full_dest_size_terms,
     estimator,
     model_settings: TourLocationComponentSettings,
     chunk_size,
@@ -245,6 +254,16 @@ def location_sample(
 
     alt_dest_col_name = model_settings.ALT_DEST_COL_NAME
 
+    if state.settings.use_explicit_error_terms:
+        stable_alt_positions = full_dest_size_terms.index.get_indexer(
+            dest_size_terms.index
+        )
+        assert (stable_alt_positions >= 0).all()
+        n_total_alts = len(full_dest_size_terms)
+    else:
+        stable_alt_positions = None
+        n_total_alts = None
+
     choices = _location_sample(
         state,
         segment_name,
@@ -257,6 +276,8 @@ def location_sample(
         chunk_size,
         chunk_tag,
         trace_label,
+        stable_alt_positions=stable_alt_positions,
+        n_total_alts=n_total_alts,
     )
 
     return choices
@@ -366,6 +387,7 @@ def location_presample(
     chunk_size,
     chunk_tag,
     trace_label,
+    full_dest_size_terms=None,
 ):
     trace_label = tracing.extend_trace_label(trace_label, "presample")
 
@@ -378,6 +400,39 @@ def location_presample(
         state, dest_size_terms, network_los, model_settings
     )
 
+    if full_dest_size_terms is None:
+        full_dest_size_terms = dest_size_terms
+
+    if state.settings.use_explicit_error_terms:
+        full_taz_index = pd.Index(
+            network_los.map_maz_to_taz(full_dest_size_terms.index), name=DEST_TAZ
+        )
+        full_taz_index = full_taz_index[~full_taz_index.duplicated()]
+        stable_alt_positions = full_taz_index.get_indexer(TAZ_size_terms.index)
+        assert (stable_alt_positions >= 0).all()
+        n_total_alts = len(full_taz_index)
+
+        # The TAZ presample call below passes stable_alt_positions for both EET and Poisson sampling, so each TAZ is
+        # keyed to its position in the full TAZ universe. The MAZ-for-TAZ second stage only receives full_taz_index for
+        # Poisson: that stage uses one per-(chooser, TAZ) uniform to pick a MAZ within each sampled TAZ. Under Poisson
+        # each sampled TAZ appears at most once per chooser, so the per-TAZ uniform produces an independent MAZ choice.
+        # Under EET sampling (importance sampling with replacement) the same TAZ can appear multiple times in a
+        # chooser's sample and would all share one uniform, forcing every duplicate to pick the same MAZ. An EET-stable
+        # MAZ-for-TAZ would need a (TAZ, occurrence-rank)-keyed draw and many more random numbers per chooser; that's
+        # too expensive with the current RNG, revisit if a counter-based RNG is adapted.
+        sample_compute_settings = getattr(model_settings, "compute_settings", None)
+        if sample_compute_settings is not None:
+            sample_compute_settings = sample_compute_settings.subcomponent_settings(
+                "sample"
+            )
+        taz_sample_method = _resolve_sample_method(state, sample_compute_settings)
+        use_stable_taz_index = taz_sample_method == "poisson"
+    else:
+        full_taz_index = None
+        stable_alt_positions = None
+        n_total_alts = None
+        use_stable_taz_index = False
+
     # convert MAZ zone_id to 'TAZ' in choosers (persons_merged)
     # persons_merged[HOME_TAZ] = persons_merged[HOME_MAZ].map(maz_to_taz)
     assert HOME_MAZ in persons_merged
@@ -411,6 +466,8 @@ def location_presample(
         chunk_tag,
         trace_label,
         zone_layer="taz",
+        stable_alt_positions=stable_alt_positions,
+        n_total_alts=n_total_alts,
     )
 
     # print(f"taz_sample\n{taz_sample}")
@@ -423,7 +480,12 @@ def location_presample(
 
     # choose a MAZ for each DEST_TAZ choice, choice probability based on MAZ size_term fraction of TAZ total
     maz_choices = tour_destination.choose_MAZ_for_TAZ(
-        state, taz_sample, MAZ_size_terms, trace_label, model_settings
+        state,
+        taz_sample,
+        MAZ_size_terms,
+        trace_label,
+        model_settings,
+        full_taz_index=full_taz_index if use_stable_taz_index else None,
     )
 
     assert DEST_MAZ in maz_choices
@@ -462,6 +524,8 @@ def run_location_sample(
     23751,      14,           0.972732479292,  2
     """
 
+    full_dest_size_terms = dest_size_terms
+
     logger.debug(
         f"dropping {(~(dest_size_terms.size_term > 0)).sum()} "
         f"of {len(dest_size_terms)} rows where size_term is zero"
@@ -496,6 +560,7 @@ def run_location_sample(
             chunk_size,
             chunk_tag=f"{chunk_tag}.presample",
             trace_label=trace_label,
+            full_dest_size_terms=full_dest_size_terms,
         )
 
     else:
@@ -505,6 +570,7 @@ def run_location_sample(
             persons_merged,
             network_los,
             dest_size_terms,
+            full_dest_size_terms,
             estimator,
             model_settings,
             chunk_size,
@@ -603,6 +669,7 @@ def run_location_simulate(
     chunk_tag,
     trace_label,
     skip_choice=False,
+    alts_context: AltsContext | None = None,
 ):
     """
     run location model on location_sample annotated with mode_choice logsum
@@ -712,6 +779,7 @@ def run_location_simulate(
         compute_settings=model_settings.compute_settings.subcomponent_settings(
             "simulate"
         ),
+        alts_context=alts_context,
     )
 
     if not want_logsums:
@@ -756,6 +824,7 @@ def run_location_choice(
     model_settings : dict
     chunk_size : int
     trace_label : str
+    skip_choice : bool
 
     Returns
     -------
@@ -789,6 +858,13 @@ def run_location_choice(
             logger.info(f"{trace_label} skipping segment {segment_name}: no choosers")
             continue
 
+        if state.settings.use_explicit_error_terms:
+            # dest_size_terms contains 0-attraction zones so using this directly here, important for stable error terms
+            # when a zone goes from 0 base -> nonzero project
+            alts_context = AltsContext.from_series(dest_size_terms.index)
+        else:
+            alts_context = None
+
         # - location_sample
         location_sample_df = run_location_sample(
             state,
@@ -841,6 +917,7 @@ def run_location_choice(
                 trace_label, "simulate.%s" % segment_name
             ),
             skip_choice=skip_choice,
+            alts_context=alts_context,
         )
 
         if estimator:
@@ -1031,6 +1108,18 @@ def iterate_location_choice(
                 ]
                 persons_merged_df_ = persons_merged_df_.sort_index()
 
+        # reset rng offsets to identical state on each iteration. This ensures that the same set of random numbers is
+        # used on each iteration for the persons being re-simulated, so sampling and final choice draws are
+        # reproducible across shadow-pricing iterations.
+        # Scoped to the persons channel for these specific rows via reset_offsets_for_df so the dedicated
+        # shadow_pricing_persons channel (registered under EET) keeps its offset across iterations and advances
+        # naturally on each iteration's update_shadow_prices call.
+        if state.settings.use_explicit_error_terms and iteration > 1:
+            logger.debug(
+                f"{trace_label} resetting random number generator offsets for iteration {iteration}"
+            )
+            state.get_rn_generator().reset_offsets_for_df(persons_merged_df_)
+
         choices_df_, save_sample_df = run_location_choice(
             state,
             persons_merged_df_,
@@ -1092,6 +1181,11 @@ def iterate_location_choice(
             )
             break
 
+    # Drop the dedicated shadow_pricing RNG channel (registered lazily under EET by spc.update_shadow_prices) so it
+    # doesn't survive into the next location_choice model (e.g., school after work) — both models share the same
+    # channel name and would otherwise collide on the no-overlap assert in SimpleChannel.extend_domain. No-op for MC.
+    spc.cleanup_rng_channel(state)
+
     # - shadow price table
     if locutor:
         if spc.use_shadow_pricing and model_settings.SHADOW_PRICE_TABLE:

activitysim/abm/models/parking_location_choice.py

@@ -21,6 +21,7 @@
 from activitysim.core.configuration.base import PreprocessorSettings
 from activitysim.core.configuration.logit import LogitComponentSettings
 from activitysim.core.interaction_sample_simulate import interaction_sample_simulate
+from activitysim.core.logit import AltsContext
 from activitysim.core.tracing import print_elapsed_time
 from activitysim.core.util import assign_in_place, drop_unused_columns
 from activitysim.core.exceptions import DuplicateWorkflowTableError
@@ -112,6 +113,7 @@ def parking_destination_simulate(
     chunk_size,
     trace_hh_id,
     trace_label,
+    alts_context: AltsContext | None = None,
 ):
     """
     Chose destination from destination_sample (with od_logsum and dp_logsum columns added)
@@ -150,6 +152,7 @@ def parking_destination_simulate(
         trace_label=trace_label,
         trace_choice_name="parking_loc",
         explicit_chunk_size=model_settings.explicit_chunk,
+        alts_context=alts_context,
     )
 
     # drop any failed zero_prob destinations
@@ -211,6 +214,9 @@ def choose_parking_location(
     )
     destination_sample.index = np.repeat(trips.index.values, len(alternatives))
     destination_sample.index.name = trips.index.name
+    # use full land_use index to ensure AltsContext spans full range of potential zones
+    land_use = state.get_dataframe("land_use")
+    alts_context = AltsContext.from_series(land_use.index)
 
     destinations = parking_destination_simulate(
         state,
@@ -223,6 +229,7 @@ def choose_parking_location(
         chunk_size=chunk_size,
         trace_hh_id=trace_hh_id,
         trace_label=trace_label,
+        alts_context=alts_context,
     )
 
     if want_sample_table: