30: MTC MNL Constrained Mode Choice#

import pandas as pd

import larch as lx
from larch import P, X

This example is a mode choice model built using the MTC example dataset. First we create the Dataset and Model objects:

m = lx.example(1)
m.title = "MTC Example 30 (Constrained Simple MNL)"
m.compute_engine = "numba"

Suppose we want to ensure that the implied value of time must be at least $5 per hour. Our data expresses cost in cents and time in minutes. Given our simple utility formulation, the target value of time is achieved if the ratio of the cost parameter to the time parameter exceeds 3.0.

We could explicitly set the ratio at exactly 3.0 by writing the model with only one parameter.

m_explicit = m.copy()
m_explicit.utility_ca = P.tottime * X.tottime + P.tottime * 3 * X.totcost
m_explicit.remove_unused_parameters()
dropped 1 parameters: totcost
result_explicit = m_explicit.maximize_loglike(stderr=True)

Iteration 049 [Optimization terminated successfully]

Best LL = -3788.5731875906613

value best initvalue minimum maximum nullvalue holdfast
param_name
ASC_BIKE -3.160e+00 -3.160e+00 0.0 -20.0 20.0 0.0 0
ASC_SR2 -2.455e+00 -2.455e+00 0.0 -20.0 20.0 0.0 0
ASC_SR3P -4.080e+00 -4.080e+00 0.0 -20.0 20.0 0.0 0
ASC_TRAN -1.893e+00 -1.893e+00 0.0 -20.0 20.0 0.0 0
ASC_WALK -1.987e+00 -1.987e+00 0.0 -20.0 20.0 0.0 0
hhinc#2 -1.987e-03 -1.987e-03 0.0 -20.0 20.0 0.0 0
hhinc#3 4.573e-04 4.573e-04 0.0 -20.0 20.0 0.0 0
hhinc#4 -6.162e-03 -6.162e-03 0.0 -20.0 20.0 0.0 0
hhinc#5 -1.409e-02 -1.409e-02 0.0 -20.0 20.0 0.0 0
hhinc#6 -9.379e-03 -9.379e-03 0.0 -20.0 20.0 0.0 0
tottime -1.739e-03 -1.739e-03 0.0 -20.0 20.0 0.0 0
from larch.model.constraints import RatioBound

m.pmaximum = {"totcost": 0, "tottime": 0}

m.constraints = [
    RatioBound("totcost", "tottime", min_ratio=3.0, max_ratio=999.0, scale=100),
]

Having created this model, we can then estimate it:

result = m.maximize_loglike(stderr=True)

Iteration 052 [Optimization terminated successfully]

Best LL = -3788.5735425913663

value best initvalue minimum maximum nullvalue holdfast
param_name
ASC_BIKE -3.161e+00 -3.161e+00 0.0 -20.0 20.0 0.0 0
ASC_SR2 -2.455e+00 -2.455e+00 0.0 -20.0 20.0 0.0 0
ASC_SR3P -4.079e+00 -4.079e+00 0.0 -20.0 20.0 0.0 0
ASC_TRAN -1.894e+00 -1.894e+00 0.0 -20.0 20.0 0.0 0
ASC_WALK -1.987e+00 -1.987e+00 0.0 -20.0 20.0 0.0 0
hhinc#2 -1.986e-03 -1.986e-03 0.0 -20.0 20.0 0.0 0
hhinc#3 4.526e-04 4.526e-04 0.0 -20.0 20.0 0.0 0
hhinc#4 -6.150e-03 -6.150e-03 0.0 -20.0 20.0 0.0 0
hhinc#5 -1.417e-02 -1.417e-02 0.0 -20.0 20.0 0.0 0
hhinc#6 -9.372e-03 -9.372e-03 0.0 -20.0 20.0 0.0 0
totcost -5.217e-03 -5.217e-03 0.0 -20.0 0.0 0.0 0
tottime -1.739e-03 -1.739e-03 0.0 -20.0 0.0 0.0 0
/opt/hostedtoolcache/Python/3.10.17/x64/lib/python3.10/site-packages/larch/model/constraints.py:119: FutureWarning: Series.__getitem__ treating keys as positions is deprecated. In a future version, integer keys will always be treated as labels (consistent with DataFrame behavior). To access a value by position, use `ser.iloc[pos]`
  return x[self.i_num] * self.cmin_num + x[self.i_den] * self.cmin_den
/opt/hostedtoolcache/Python/3.10.17/x64/lib/python3.10/site-packages/larch/model/constraints.py:122: FutureWarning: Series.__getitem__ treating keys as positions is deprecated. In a future version, integer keys will always be treated as labels (consistent with DataFrame behavior). To access a value by position, use `ser.iloc[pos]`
  return x[self.i_num] * self.cmax_num + x[self.i_den] * self.cmax_den
result
keyvalue
x
0
ASC_BIKE -3.161e+00
ASC_SR2 -2.455e+00
ASC_SR3P -4.079e+00
ASC_TRAN -1.894e+00
ASC_WALK -1.987e+00
hhinc#2 -1.986e-03
hhinc#3 4.526e-04
hhinc#4 -6.150e-03
hhinc#5 -1.417e-02
hhinc#6 -9.372e-03
totcost -5.217e-03
tottime -1.739e-03
loglossnp.float64(0.7533453057449525)
d_logloss
0
ASC_BIKE 3.803e-05
ASC_SR2 -8.158e-06
ASC_SR3P 6.016e-06
ASC_TRAN -2.386e-06
ASC_WALK -2.598e-05
hhinc#2 -4.142e-04
hhinc#3 -6.993e-06
hhinc#4 -7.637e-04
hhinc#5 2.198e-03
hhinc#6 -5.378e-04
totcost 4.809e-01
tottime -1.442e+00
nit52
nfev71
njev52
status0
message'Optimization terminated successfully'
successnp.True_
elapsed_time0:00:00.757537
method'slsqp'
n_cases5029
iteration_number52
loglikenp.float64(-3788.5735425913663)
m.parameter_summary()
  Value Std Err t Stat Signif Null Value Constrained
Parameter            
ASC_BIKE -3.16  0.309 -10.21 *** 0.00
ASC_SR2 -2.46  0.103 -23.95 *** 0.00
ASC_SR3P -4.08  0.175 -23.26 *** 0.00
ASC_TRAN -1.89  0.112 -16.87 *** 0.00
ASC_WALK -1.99  0.169 -11.78 *** 0.00
hhinc#2 -0.00199  0.00154 -1.29 0.00
hhinc#3  0.000453  0.00251  0.18 0.00
hhinc#4 -0.00615  0.00180 -3.41 *** 0.00
hhinc#5 -0.0142  0.00551 -2.57 * 0.00
hhinc#6 -0.00937  0.00306 -3.07 ** 0.00
totcost -0.00522  0.000243 -21.49 *** 0.00 totcost / tottime ≥ 3.0
tottime -0.00174  8.09e-05 -21.49 *** 0.00 totcost / tottime ≥ 3.0

It is a little tough to read this report because the parameters can show up in pretty much any order, as they are not sorted when they are automatically discovered by Larch. We can use the reorder method to fix this:

m.ordering = (
    (
        "LOS",
        "totcost",
        "tottime",
    ),
    (
        "ASCs",
        "ASC.*",
    ),
    (
        "Income",
        "hhinc.*",
    ),
)
m.parameter_summary()
    Value Std Err t Stat Signif Null Value Constrained
Category Parameter            
LOS totcost -0.00522  0.000243 -21.49 *** 0.00 totcost / tottime ≥ 3.0
tottime -0.00174  8.09e-05 -21.49 *** 0.00 totcost / tottime ≥ 3.0
ASCs ASC_BIKE -3.16  0.309 -10.21 *** 0.00
ASC_SR2 -2.46  0.103 -23.95 *** 0.00
ASC_SR3P -4.08  0.175 -23.26 *** 0.00
ASC_TRAN -1.89  0.112 -16.87 *** 0.00
ASC_WALK -1.99  0.169 -11.78 *** 0.00
Income hhinc#2 -0.00199  0.00154 -1.29 0.00
hhinc#3  0.000453  0.00251  0.18 0.00
hhinc#4 -0.00615  0.00180 -3.41 *** 0.00
hhinc#5 -0.0142  0.00551 -2.57 * 0.00
hhinc#6 -0.00937  0.00306 -3.07 ** 0.00