1: MTC MNL Mode Choice#

import pandas as pd

import larch as lx

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

m = lx.Model(d)

Then we can build up the utility function. We’ll use some :ref:idco data first, using the Model.utility.co attribute. This attribute is a dict-like object, to which we can assign :class:LinearFunction objects for each alternative code.

from larch import PX, P, X

m.utility_co[2] = P("ASC_SR2") + P("hhinc#2") * X("hhinc")
m.utility_co[3] = P("ASC_SR3P") + P("hhinc#3") * X("hhinc")
m.utility_co[4] = P("ASC_TRAN") + P("hhinc#4") * X("hhinc")
m.utility_co[5] = P("ASC_BIKE") + P("hhinc#5") * X("hhinc")
m.utility_co[6] = P("ASC_WALK") + P("hhinc#6") * X("hhinc")

Next we’ll use some idca data, with the utility_ca attribute. This attribute is only a single :class:LinearFunction that is applied across all alternatives using :ref:idca data. Because the data is structured to vary across alternatives, the parameters (and thus the structure of the :class:LinearFunction) does not need to vary across alternatives.

m.utility_ca = PX("tottime") + PX("totcost")

Lastly, we need to identify :ref:idca data that gives the availability for each alternative, as well as the number of times each alternative is chosen. (In traditional discrete choice analysis, this is often 0 or 1, but it need not be binary, or even integral.)

m.availability_ca_var = "avail"
m.choice_ca_var = "chose"

And let’s give our model a descriptive title.

m.title = "MTC Example 1 (Simple MNL)"

We can view a summary of the choices and alternative availabilities to make sure the model is set up correctly.

m.choice_avail_summary()

	name	chosen	available
1	DA	3637	4755
2	SR2	517	5029
3	SR3+	161	5029
4	Transit	498	4003
5	Bike	50	1738
6	Walk	166	1479
< Total All Alternatives >		5029	<NA>

We’ll set a parameter cap (bound) at +/- 20, which helps improve the numerical stability of the optimization algorithm used in estimation.

m.set_cap(20)

Having created this model, we can then estimate it:

result = m.maximize_loglike(stderr=True)

┏━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━┓
┃ Larch Model Dashboard                                                                        ┃
┡━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━┩
│ optimization complete                                                                        │
│ Log Likelihood Current =      -3626.185791 Best =      -3626.185791                          │
│ ┏━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━┳━━━━━━━━━━━━━━━━┳━━━━━━━━━━━━━━━━┳━━━━━━━━━━━━┳━━━━━━━━━━┓ │
│ ┃ Parameter                      ┃ Estimate       ┃ Std. Error     ┃ t-Stat     ┃ Null Val ┃ │
│ ┡━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━╇━━━━━━━━━━━━━━━━╇━━━━━━━━━━━━━━━━╇━━━━━━━━━━━━╇━━━━━━━━━━┩ │
│ │ ASC_BIKE                       │ -2.376399      │  0.30451715    │ -7.804     │  0       │ │
│ │ ASC_SR2                        │ -2.1779427     │  0.10463986    │ -20.81     │  0       │ │
│ │ ASC_SR3P                       │ -3.7254057     │  0.17768404    │ -20.97     │  0       │ │
│ │ ASC_TRAN                       │ -0.67102525    │  0.13259123    │ -5.061     │  0       │ │
│ │ ASC_WALK                       │ -0.20671203    │  0.19410205    │ -1.065     │  0       │ │
│ │ hhinc#2                        │ -0.0021717871  │  0.0015533642  │ -1.398     │  0       │ │
│ │ hhinc#3                        │  0.00036546088 │  0.0025373173  │  0.144     │  0       │ │
│ │ hhinc#4                        │ -0.0052828484  │  0.0018287523  │ -2.889     │  0       │ │
│ │ hhinc#5                        │ -0.012807539   │  0.0053243754  │ -2.405     │  0       │ │
│ │ hhinc#6                        │ -0.0096854051  │  0.0030330766  │ -3.193     │  0       │ │
│ │ totcost                        │ -0.0049202335  │  0.00023889064 │ -20.6      │  0       │ │
│ │ tottime                        │ -0.051346785   │  0.0030995551  │ -16.57     │  0       │ │
│ └────────────────────────────────┴────────────────┴────────────────┴────────────┴──────────┘ │
└──────────────────────────────────────────────────────────────────────────────────────────────┘

m.parameter_summary()

	Value	Std Err	t Stat	Signif	Null Value
Parameter
ASC_BIKE	-2.38	0.305	-7.80	***	0.00
ASC_SR2	-2.18	0.105	-20.81	***	0.00
ASC_SR3P	-3.73	0.178	-20.97	***	0.00
ASC_TRAN	-0.671	0.133	-5.06	***	0.00
ASC_WALK	-0.207	0.194	-1.06		0.00
hhinc#2	-0.00217	0.00155	-1.40		0.00
hhinc#3	0.000365	0.00254	0.14		0.00
hhinc#4	-0.00528	0.00183	-2.89	**	0.00
hhinc#5	-0.0128	0.00532	-2.41	*	0.00
hhinc#6	-0.00969	0.00303	-3.19	**	0.00
totcost	-0.00492	0.000239	-20.60	***	0.00
tottime	-0.0513	0.00310	-16.57	***	0.00

It is a little tough to read this report because the parameters show up in alphabetical order. We can use the reorder method to fix this and group them systematically:

m.ordering = (
    (
        "LOS",
        "totcost",
        "tottime",
    ),
    (
        "ASCs",
        "ASC.*",
    ),
    (
        "Income",
        "hhinc.*",
    ),
)

m.parameter_summary()

		Value	Std Err	t Stat	Signif	Null Value
Category	Parameter
LOS	totcost	-0.00492	0.000239	-20.60	***	0.00
LOS	tottime	-0.0513	0.00310	-16.57	***	0.00
ASCs	ASC_BIKE	-2.38	0.305	-7.80	***	0.00
	ASC_SR2	-2.18	0.105	-20.81	***	0.00
	ASC_SR3P	-3.73	0.178	-20.97	***	0.00
	ASC_TRAN	-0.671	0.133	-5.06	***	0.00
	ASC_WALK	-0.207	0.194	-1.06		0.00
Income	hhinc#2	-0.00217	0.00155	-1.40		0.00
	hhinc#3	0.000365	0.00254	0.14		0.00
	hhinc#4	-0.00528	0.00183	-2.89	**	0.00
	hhinc#5	-0.0128	0.00532	-2.41	*	0.00
	hhinc#6	-0.00969	0.00303	-3.19	**	0.00

m.estimation_statistics()

Statistic	Aggregate	Per Case
Number of Cases	5029
Log Likelihood at Convergence	-3626.19	-0.72
Log Likelihood at Null Parameters	-7309.60	-1.45
Rho Squared w.r.t. Null Parameters	0.504

1: MTC MNL Mode Choice#

This Page