Explicit error terms production code

[janzill]

This PR brings the explicit error term (EET) work to production standard. It contains the following major changes compared to the PoC implementation:
1. It adds testing and documentation for EET
2. It decouples the sampling method from the simulation method and adds Poisson sampling (based on @m-richards implementation in #1065), including tests and documentation. It also addresses runtime and memory issues with the sampling method "eet".
3. It removes several inconsistencies in the EET simulation branch where edge cases could have led to unexpected changes in choices for individual choosers due to non-alignment of error terms.
4. It consistently implements EET for nested logit models

Note that the default simulation method remains Monte Carlo, with all existing unit and integration tests unchanged and passing. Users should therefore not see any differences in their model runs unless explicitly opting in to EET. This is done by adding use_explicit_error_terms: True to the settings. This is a drop-in replacement, with default settings leading to a 3-10% runtime increase of a single demand model run for the models we have tested so far.

This is a large PR and I will add details in comments below so we can keep discussions focused.

[janzill]

@dhensle target is now main as discussed earlier. I cannot resolve the copilot suggestions, but they are all fixed.

[dhensle]

Looks good! Documentation is thorough. Just a couple of things for you to review.

[dhensle]

Suggest we parameterize this 3.0 value in logit.py (named like UTIL_LARGE_ENOUGH or something). Shows up in interaction_sample_simulate too in the same context.

[janzill]

Done. Replaced all three occurrences and added a comment on usage too.

[dhensle]

Is this fix something that should get just get implemented right now? Seem better than just leaving this fixme in here for a who-knows amount of time. (The proposed fix make sense to me)

[janzill]

Ideally yes, but it is not as simple as the idea suggests and therefore not within budget. Something to discuss in the next engineering meeting.

[dhensle]

Do we need to be careful about the sampling method for estimation? The actual choice method (eet vs monte-carlo) shouldn't matter because the choice just gets overwritten by the survey choice anyways. But is the sampling methodology important?

I don't think it should matter, but would like to double check.

[janzill]

Currently, destination choice logsums are different for Poisson because of the different meaning of prob (Poisson reports inclusion probabilities). That could be worked around, and in general I don't think the method matters theoretically, only in the sense that some additional things would need to be implemented. I do not see any advantage to this currently, and this is already a very large PR, so I would suggest just requiring MC for estimation. Looking at Jeff's update of the EET branch at the start of this phase, it looks like he added a "estimation needs MC sampling" explicitly (https://github.com/ActivitySim/activitysim/pull/1042/changes#diff-9ec5e393569d30759b9ced82065107dd053d6f1901a4ad20df23825f19fe5f22). @jpn-- anything I missed here?

[dhensle]

Combining this into a single resolve_sample_method seems cleaner. I don't think we need another function here.

[janzill]

Agreed and done.

[dhensle]

If this is true, then we should enforce it in the settings. But I am not sure why monte_carlo should be required for estimation.

[janzill]

Regarding if monte_carlo is required, see my answer in base.py. On the validation I agree in principle, but the sample_method can be set globally and per-component in compute_settings, and it looks like enabling estimation comes from a separate config, which is only read on first use in a model.

[dhensle]

Remove duplicate call.

[janzill]

Done

[dhensle]

I am not sure I really like how these retries are being applied. I see two issues:

First, if an alternative gets selected in a re-try, I think the inclusion probability correction np.log(pick_count / prob) is slightly off:
The current correction is
np.log(1/ $q_{i}$ )
where
$q_{i} = 1 - (1-p_{i})^{S}$ is the sampling probability for in the first attempt.
For n > 1, the prob should be corrected to include the probability of an empty sample
$P_{0} = \prod_{j} (1 - q_{j})$
as the probability of selection would be $q_i / (1.0 - P_0)$

I think this correction would be generally quite small, but the rows that actually require retries are the rows that would comparatively have a large $P_{0}$.

Second, setting the prob = 1 for the selected alternatives that pass the fallback method again isn't quite right as the correction factor therefore becomes np.log(1/prob) = log(1/1) = 0. But, the probability of being selected isn't 1 -- it's whatever the selection probability is in the fallback method.

Can we just avoid retrying and if we do not return anything in the first iteration, just select all alternatives that have a non-zero probability? Then their selection prob is 1 and the correction factor becomes moot and we don't waste time and complexity iterating?

I am also wondering whether we should be checking for a "zero probability" error here before we start this looping (or choosing) to catch cases where no alternative is available.

It's difficult for me to say how much any of this matters in practice....

[janzill]

Yes, these are good points and I agree on the math. The variable-number retries also potentially could mess up error term alignment. In practice, these issues did not show up for sample sizes of 30 or above and "peaky-enough" alternative sets, which is what is currently used in all models I know of, but I agree that it's good to fix this. I am a bit hesitant to add all alternatives here because this would slow down things considerably due to mode choice logsums being very slow. However, there is a solution: Choose the sample_size highest-utility alternatives as a fallback set when Poisson sampling came up empty. Those two events are mutually exclusive and the total inclusion probability becomes $prob_i = q_i + P_0 \times \delta_{i \in fallback}$. This solves the correction probability issue that was there before, as well as the EET error term alignment issue. I've added more documentation in both the code and the official documents on this. I've also added a runtime optimization where the second term is only calculated if $P_0$ is larger than $1e-12$ to avoid calculating if an alternative is in the top sample_size values, which requires a sort of the array, see comments in code.

[janzill]

Regarding checking for zero probabilities, that's handled in utils to probs

[janzill]

On another note, Poisson sampling could also be replaced with probabilistic top-k sampling, there is a known closed-form expression for inclusion probabilities and this would be easy to implement and some preliminary tests suggested that it has the same level of stability as Poisson sampling. However, with the above changes I do not think this is necessary.

[dhensle]

should this allow_zero_probs be True? Shouldn't this be passed in from our global/model-level settings?

[janzill]

It should not be, good catch. Fixed.

[dhensle]

if not active.any() (i.e. there are no available choices), then the argmax below just selects the first one? I am concerned that we aren't correctly catching the "zero probability" case in this function. We should be reporting zero probabilities and skipping their choices if that setting is enabled.

[janzill]

I thought we can just rely on validate_utils always being called upstream in simulate, sample_simulate, etc. On second thought, this is not a good pattern. For MC there is a utils_to_probs step outside make_sample_choices (the equivalent to validate_utils), and then there is an argument called allow_bad_probs in make_choices. I've added a corresponding argument to make_choices_explicit_error_term that mirrors the allow_bad_probs functionality. There's also more documentation now to make things clearer.

[dhensle]

Suggest we check the actual index (and columns?) here instead of just the shape.

[janzill]

Done

[janzill]

Thanks for the thorough review @dhensle, some good suggestions there. Your sample correction factor for Poisson retries led to another idea which does not require any retries and removes potential error-term misalignment between scenarios, let me know what you think.