ActivitySim / activitysim

An Open Platform for Activity-Based Travel Modeling
https://activitysim.github.io
BSD 3-Clause "New" or "Revised" License
192 stars 99 forks source link

:exclamation: Possible math error in some scheduling models #605

Open jpn-- opened 2 years ago

jpn-- commented 2 years ago

There may be a hidden math error in some scheduling models: generally, any model where the is a shift term that is the square of a time period value (start time, end time, duration).

Trip and tour start and end times are stored as int8 in the various tables where they appear. By extension, durations also may be int8. This seems superficially fine, as the most ambitious model probably won't have less than 15 minute time periods, so the max number of periods in a day would by 96, less than the max int8 value of 127.

BUT: If you write an expression in a spec file doing math with time periods, AND if that expression might overflow the int8 e.g. by squaring, AND if you don't explicitly upcast the to a more appropriate dtype, THEN the math can overflow and wrap around. We observe this if we take the square of an int8 array of time period values:

>>> pd.Series(np.arange(24), dtype=np.int8) ** 2
0       0
1       1
2       4
3       9
4      16
5      25
6      36
7      49
8      64
9      81
10    100
11    121
12   -112
13    -87
14    -60
15    -31
16      0
17     33
18     68
19    105
20   -112
21    -71
22    -28
23     17
dtype: int8

The consequences of this math error may wash out if associated coefficients are zero, or if the squared values are always 11 or less by construction or coincidence, or potentially for other reasons. Some study will be required to determine if that is the case, but my preliminary guess is that it does not wash out at least sometimes, based on discovering this through anomalies between sharrow-with-time-windows and legacy code on the prototype MWCOG model.

Based on a cursory review of model spec files in our repository, I observe that scheduling for joint and non-mandatory tours may be affected in these models:

cc: @guyrousseau @JilanChen @mmoran3 @ray-ngo @jfdman @i-am-sijia

dhensle commented 2 years ago

Investigating this issue in the re-estimated SEMCOG configs (currently in the SEMCOG example pull request #603) shows that there are just a couple scheduling segments with non-zero squared coefficients. These segments and coefficients are: shopping tours (-0.0022), joint shopping tours, (-.0041), and maintenance tours (-0.0017).

The problem is perhaps even worse in SEMCOG for these cases because the models have 48 half hour time periods. The following table shows the worst case magnitude of the problem (joint shopping tours) for the implemented models. The utility difference between the actual and expected utility can be up to 4 utiles (!) for some time periods.

coef = -.004084
tt = pd.DataFrame()
tt['start'] = pd.Series(np.arange(48), dtype=np.int8)
tt['start^2'] = tt['start'] ** 2
tt['actual utility'] = (np.where((tt.start>17), (((17-tt.start)*(tt.start<=17) + (tt.start-17)*(tt.start>17)) ** 2), 0) * coef).round(2)
tt['expected utility'] = (np.where((tt.start>17), (((17-tt.start)*(tt.start<=17) + (tt.start-17)*(tt.start>17)).astype(np.float32) ** 2), 0) * coef).round(2)
tt['difference'] = (tt['expected utility'] - tt['actual utility']).round(2)
start start^2 actual utility expected utility difference
0 0 -0.00 -0.00 0.00
1 1 -0.00 -0.00 0.00
2 4 -0.00 -0.00 0.00
3 9 -0.00 -0.00 0.00
4 16 -0.00 -0.00 0.00
5 25 -0.00 -0.00 0.00
6 36 -0.00 -0.00 0.00
7 49 -0.00 -0.00 0.00
8 64 -0.00 -0.00 0.00
9 81 -0.00 -0.00 0.00
10 100 -0.00 -0.00 0.00
11 121 -0.00 -0.00 0.00
12 -112 -0.00 -0.00 0.00
13 -87 -0.00 -0.00 0.00
14 -60 -0.00 -0.00 0.00
15 -31 -0.00 -0.00 0.00
16 0 -0.00 -0.00 0.00
17 33 -0.00 -0.00 0.00
18 68 -0.00 -0.00 0.00
19 105 -0.02 -0.02 0.00
20 -112 -0.04 -0.04 0.00
21 -71 -0.07 -0.07 -0.00
22 -28 -0.10 -0.10 -0.00
23 17 -0.15 -0.15 -0.00
24 64 -0.20 -0.20 -0.00
25 113 -0.26 -0.26 0.00
26 -92 -0.33 -0.33 -0.00
27 -39 -0.41 -0.41 0.00
28 16 -0.49 -0.49 -0.00
29 73 0.46 -0.59 -1.05
30 -124 0.36 -0.69 -1.05
31 -63 0.25 -0.80 -1.05
32 0 0.13 -0.92 -1.05
33 65 -0.00 -1.05 -1.05
34 -124 -0.13 -1.18 -1.05
35 -55 -0.28 -1.32 -1.04
36 16 -0.43 -1.47 -1.04
37 89 0.46 -1.63 -2.09
38 -92 0.29 -1.80 -2.09
39 -15 0.11 -1.98 -2.09
40 64 -0.07 -2.16 -2.09
41 -111 -0.26 -2.35 -2.09
42 -28 -0.46 -2.55 -2.09
43 57 0.38 -2.76 -3.14
44 -112 0.16 -2.98 -3.14
45 -23 -0.07 -3.20 -3.13
46 68 -0.30 -3.43 -3.13
47 -95 0.51 -3.68 -4.19
ray-ngo commented 2 years ago

The problem would apply to MWCOG Model as it is also using 48 half-hour time periods.

guyrousseau commented 2 years ago

Indeed, the ARC ABM functions at a temporal resolution of 30 minutes, so in half-hour increments. Although activities are scheduled with a 30-minute resolution, level-of-service matrices are only created for five aggregate time-of-day periods: early A.M., A.M., Midday, P.M., and evening. Trips occurring in each time-of-day period reference the appropriate transportation network depending on their trip mode and the mid-point trip time. That said, we probably still need to look into this possible math error in our scheduling models, moving forward, so thanks for bringing this up to our attention.

dhensle commented 2 years ago

The affected SEMCOG tour scheduling models have been re-estimated after removing the squared terms discussed above. The changes can be seen in this commit and are incorporated in the updated SEMCOG example. By removing all squared terms from the scheduling configs, future re-estimation with these configs will not run into this issue unless they explicitly add some squared terms back in.