<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Choice Models
Step1:  In this guide, we'll take a look at building a discrete choice model using Larch.
We assume you have a decent grasp of the fundamentals of choice modeling -- if not,
we suggest reading the 
Discrete Choice Modeling
section of the Python for Transportation Modeling course.
Some addition advanced or detailed topics are broken out into seperate guides
Step2:  The basic structure of a choice model in Larch is contained in the 
Model object.
Step3:  Choices
The dependent variable for a discrete choice model is an array that describes the choices.
In Larch, there are three different ways to indicate choices, by assigning to different
attributes
Step4:  After setting the choice definition, the loaded or computed choice array
should be available as the 'ch' DataArray in the model's 
dataset.
Step5:  Availability
In addition to the choices, we can also define an array that describes
the availability of the various alternatives. Unlike the choices, for
the availability factors we expect that we'll need to toggle the 
availability on or off for potentially every alternative in each case,
so there's only two ways to define availability, by assigning to
attributes
Step6:  After setting the availability definition, the loaded or computed availability array
should be available as the 'av' DataArray in the model's 
dataset.
Step7:  Utility Functions
Choice models in Larch rely on utility expressions that are 
linear-in-parameters functions, which combine 
parameters P and data 
X. You can attach these function
expressions to the model in two ways
Step8:  The computed values for the utility function can be accessed using 
the utility method, which also permits
the user to set new values for various model parameters.
Step9:  Data Preparation
Larch works with two "tiers" of data
Step10:  Nesting Structures
By default, a model in Larch is assumed to be a simple multinomial 
logit model, unless a nesting structure is defined. That structure
is defined in a model's graph.
Step11:  Adding a nest can be accomplished the the new_node method,
which allows you to give a nesting node's child codes, a name, and attach a logsum parameter.
Step12:  The return value of new_node
is the code number of the new nest. This is assigned automatically so 
as to not overlap with any other alternatives or nests.  We can use this
to develop multi-level nesting structures, by putting that new code 
number as the child for yet another new nest.
Step13:  Nothing in Larch prevents you from overloading the nesting structure with
degenerate nests, as shown above.  You may have difficult with estimating
parameters if you are not careful with such complex structures.  If you
need to remove_node, you 
can do so by giving its code--but you'll likely find you'll be much better off
just fixing your code and starting over, as node removal can have some odd
side effects for complex structures.
Step14:  Parameter Estimation
Larch can automatically find all the model parameters contained in the model
specification, so we don't need to address them separately unless we want
to modify any defaults.
We can review the parameters Larch has found, as well as the current values
set for them, in the parameter frame, or pf.
Step15:  If we want to set certain parameters to be constrained to be certain values,
that can be accomplished with the lock_value method.
Because our sample data has so few observations, it won't be possible to estimate
values for all four parameters, so we can assert values for two of them.
Step16:  The default infinite bounds on the remaining parameters can be problematic
for some optimization algorithms, so it's usually good practice to set large
but finite limits for those values.  The set_cap method
can do just that, setting a minimum and maximum value for all the parameters
that otherwise have bounds outside the cap.
Step17:  To actually develop maximum likelihood estimates for the remaining
unconstrained parameters, use the 
maximize_loglike method.
Step18:  In a Jupyter notebook, this method displays a live-updating view of the
progress of the optmization algorithm, so that the analyst can interrupt
if something looks wrong.
The maximize_loglike method does 
not include the calculation of parameter covariance matrixes, standard
errors, or t-statistics.  For large models, this can be a computationally 
expensive process, and it is often but not always necessary. Those
computatations are made in the
calculate_parameter_covariance
method instead.  Once completed, things like t-statistics and standard
errors are available in the parameter frame.
Step19:  Reporting
Larch includes a variety of pre-packaged and a la carte reporting options.
Commonly used report tables are available directly in a Jupyter notebook
through a selection of reporting functions.
Step20:  To save a model report to an Excel file, use the to_xlsx method.
<ASSISTANT_TASK:>
Python Code:
# TEST
import numpy as np
import pandas as pd
import larch.numba as lx
from pytest import approx
import larch.numba as lx
Explanation: Choice Models
End of explanation
# HIDDEN
df_ca = pd.read_csv("example-data/tiny_idca.csv")
cats = df_ca['altid'].astype(pd.CategoricalDtype(['Car', 'Bus', 'Walk'])).cat
df_ca['altnum'] = cats.codes + 1
df_ca = df_ca.set_index(['caseid', 'altnum'])
data = lx.Dataset.construct.from_idca(df_ca.sort_index(), fill_missing=0)
data = data.drop_vars("_avail_")
data['ChosenCode'] = (data['Chosen'] * data['Chosen'].altnum).sum('altnum')
data.coords['alt_names'] = lx.DataArray(cats.categories, dims=('altnum'), coords={'altnum': data.altnum})
alts = dict(zip(data['altnum'].values, data['alt_names'].values))
for i, k in alts.items():
    data[f'{k}Time'] = data['Time'].sel(altnum=i)
data
Explanation: In this guide, we'll take a look at building a discrete choice model using Larch.
We assume you have a decent grasp of the fundamentals of choice modeling -- if not,
we suggest reading the 
Discrete Choice Modeling
section of the Python for Transportation Modeling course.
Some addition advanced or detailed topics are broken out into seperate guides:
Data Fundamentals for Choice Models
Working with Linear Functions
Machine Learning Integration
The examples below work with the tiny dataset introduced in the 
Data Fundamentals section.
End of explanation
m = lx.Model(data)
Explanation: The basic structure of a choice model in Larch is contained in the 
Model object.
End of explanation
m.choice_co_code = 'ChosenCode'
# TEST
ch = np.array([ 
    [ 1.,  0.,  0.,  1.],
    [ 0.,  1.,  0.,  1.],
    [ 0.,  0.,  1.,  1.],
    [ 0.,  0.,  1.,  1.],
])
assert m.loglike() == approx(np.log(1/3) * 4) 
np.testing.assert_array_almost_equal(m._data_arrays.ch, ch)
m.choice_co_vars = {
    1: 'ChosenCode == 1',
    2: 'ChosenCode == 2',
    3: 'ChosenCode == 3',
}
# TEST
assert m.loglike() == approx(np.log(1/3) * 4) 
np.testing.assert_array_almost_equal(m._data_arrays.ch, ch)
m.choice_ca_var = 'Chosen'
# TEST
assert m.loglike() == approx(np.log(1/3) * 4) 
np.testing.assert_array_almost_equal(m._data_arrays.ch, ch)
Explanation: Choices
The dependent variable for a discrete choice model is an array that describes the choices.
In Larch, there are three different ways to indicate choices, by assigning to different
attributes:
m.choice_ca_var
: The choices are given by indicator values (typically but not 
neccessarily dummy variables) in an idca variable.
m.choice_co_code
: The choices are given by altid values in an idco variable. 
These choice codes are then converted to binary indicators
by Larch.
m.choice_co_vars
: The choices are given by indicator values (typically but not 
neccessarily dummy variables) in multiple idco variables,
one for each alternative.
Given the dataset (which has all these formats defined), all 
the following choice definitions result in the same choice
representation:
End of explanation
m.dataset['ch']
Explanation: After setting the choice definition, the loaded or computed choice array
should be available as the 'ch' DataArray in the model's 
dataset.
End of explanation
m.availability_ca_var = "Time > 0"
# TEST
assert m.loglike_casewise() == approx([-1.098612, -0.693147, -1.098612, -1.098612])
m.availability_co_vars = {
    1: True,
    2: 'BusTime > 0',
    3: 'WalkTime > 0',
}
# TEST
assert m.availability_var is None
assert m.loglike_casewise() == approx([-1.098612, -0.693147, -1.098612, -1.098612])
Explanation: Availability
In addition to the choices, we can also define an array that describes
the availability of the various alternatives. Unlike the choices, for
the availability factors we expect that we'll need to toggle the 
availability on or off for potentially every alternative in each case,
so there's only two ways to define availability, by assigning to
attributes:
m.availability_ca_var
: The availability of alternatives is given by binary values 
(booleans, or equivalent integers) in an idca variable.
m.availability_co_vars
: The availability of alternatives is given by binary values 
(booleans, or equivalent integers) in multiple idco variables,
one for each alternative.
Given the dataset, both of the following availability definitions 
result in the same availability representation:
End of explanation
m.dataset['av']
Explanation: After setting the availability definition, the loaded or computed availability array
should be available as the 'av' DataArray in the model's 
dataset.
End of explanation
from larch import P, X
m.utility_ca = P.Time * X.Time + P.Cost * X.Cost
# TEST
assert m.loglike({'Time': -0.01, 'Cost': -0.02}) == approx(-4.040671968970389)
assert m.utility() == approx(np.array([
    [-3.3     , -2.4     , -0.2     , -0.055162],
    [-2.75    , -2.35    ,   -np.inf, -1.836985],
    [-2.9     , -2.      , -0.3     , -0.071306],
    [-4.65    , -3.2     , -0.1     , -0.045875],
]), rel=1e-5)
m.utility_co = {
    1: P.Income_Car * X.Income / 1000,
    2: P.Income_Bus * X.Income / 1000,
}
# TEST
assert m.loglike({'Time': -0.01, 'Cost': -0.02, 'Income_Car': 0.1}) == approx(-6.089136524348358)
assert m.utility({
    'Time': -0.01, 'Cost': -0.02, 'Income_Car': 0.1
}) == approx(np.array([
    [-0.3     , -2.4     , -0.2     ,  0.500937],
    [ 0.25    , -2.35    , -np.inf  ,  0.321645],
    [ 1.1     , -2.      , -0.3     ,  1.355918],
    [ 0.35    , -3.2     , -0.1     ,  0.860637],
]))
Explanation: Utility Functions
Choice models in Larch rely on utility expressions that are 
linear-in-parameters functions, which combine 
parameters P and data 
X. You can attach these function
expressions to the model in two ways:
m.utility_ca
: A linear function containing generic expressions 
that are evaluated against the idca portion of the dataset. 
These expression can technically also reference idco variables,
but to define a well-specified choice model with identifiable
parameters, each data term will need at least one idca
component.
m.utility_co
: A mapping of alternative-specific expressions that are evaluated 
against only the idco portion of the dataset. Each key gives
an alternative id, and the values are linear functions. 
These two utility function definitions are not mutually
exclusive, and you can mix both types of functions in the same
model.
End of explanation
m.utility(
    {'Time': -0.01, 'Cost': -0.02, 'Income_Car': 0.1},
    return_format='dataarray',
)
Explanation: The computed values for the utility function can be accessed using 
the utility method, which also permits
the user to set new values for various model parameters.
End of explanation
m.datatree
m.dataset
Explanation: Data Preparation
Larch works with two "tiers" of data:
m.datatree 
: A DataTree that holds the raw data used for the model. This can
consist of just a single Dataset, (which is automatically converted
into a one-node tree when you assign it to this attribute) or multiple
related datasets linked together using the sharrow library. 
m.dataset 
: The assembled arrays actually used in calculation, stored
as a Dataset that has variables for various required data elements
and dimensions structured to support the model design.
The dataset is wiped out when any aspect of the model structure is
changed, and rebuilt as needed for computation. For 
particular applications that require specialized
optimization, the dataset can be provided exogenously after the
model stucture is finalized, but generally 
it will be convenient for users to let Larch build the dataset 
automatically from a datatree.
End of explanation
m.graph
# TEST
assert sorted(m.graph.nodes) == [0, 1, 2, 3]
assert sorted(m.graph.edges) == [(0, 1), (0, 2), (0, 3)]
assert m.graph.standard_sort_names == ['Car', 'Bus', 'Walk', '_root_']
assert m.graph.standard_sort == (1, 2, 3, 0)
Explanation: Nesting Structures
By default, a model in Larch is assumed to be a simple multinomial 
logit model, unless a nesting structure is defined. That structure
is defined in a model's graph.
End of explanation
z = m.graph.new_node(parameter='Mu_Motorized', children=[1,2], name='Motorized')
m.graph
# TEST
assert sorted(m.graph.nodes) == [0, 1, 2, 3, 4]
assert sorted(m.graph.edges) == [(0, 3), (0, 4), (4, 1), (4, 2)]
assert m.graph.standard_sort_names == ['Car', 'Bus', 'Walk', 'Motorized', '_root_']
assert m.graph.standard_sort == (1, 2, 3, 4, 0)
Explanation: Adding a nest can be accomplished the the new_node method,
which allows you to give a nesting node's child codes, a name, and attach a logsum parameter.
End of explanation
m.graph.new_node(parameter='Mu_Omni', children=[z, 3], name='Omni')
m.graph
# TEST
assert sorted(m.graph.nodes) == [0, 1, 2, 3, 4, 5]
assert sorted(m.graph.edges) == [(0, 5), (4, 1), (4, 2), (5, 3), (5, 4)]
assert m.graph.standard_sort_names == [
    'Car', 'Bus', 'Walk', 'Motorized', 'Omni', '_root_'
]
assert m.graph.standard_sort == (1, 2, 3, 4, 5, 0)
Explanation: The return value of new_node
is the code number of the new nest. This is assigned automatically so 
as to not overlap with any other alternatives or nests.  We can use this
to develop multi-level nesting structures, by putting that new code 
number as the child for yet another new nest.
End of explanation
m.graph.remove_node(5)
m.graph
# TEST
assert sorted(m.graph.nodes) == [0, 1, 2, 3, 4]
assert sorted(m.graph.edges) == [(0, 3), (0, 4), (4, 1), (4, 2)]
assert m.graph.standard_sort_names == ['Car', 'Bus', 'Walk', 'Motorized', '_root_']
assert m.graph.standard_sort == (1, 2, 3, 4, 0)
Explanation: Nothing in Larch prevents you from overloading the nesting structure with
degenerate nests, as shown above.  You may have difficult with estimating
parameters if you are not careful with such complex structures.  If you
need to remove_node, you 
can do so by giving its code--but you'll likely find you'll be much better off
just fixing your code and starting over, as node removal can have some odd
side effects for complex structures.
End of explanation
m.pf
Explanation: Parameter Estimation
Larch can automatically find all the model parameters contained in the model
specification, so we don't need to address them separately unless we want
to modify any defaults.
We can review the parameters Larch has found, as well as the current values
set for them, in the parameter frame, or pf.
End of explanation
m.lock_value('Time', -0.01)
m.lock_value('Cost', -0.02)
m.pf
Explanation: If we want to set certain parameters to be constrained to be certain values,
that can be accomplished with the lock_value method.
Because our sample data has so few observations, it won't be possible to estimate
values for all four parameters, so we can assert values for two of them.
End of explanation
m.set_cap(100)
m.pf
Explanation: The default infinite bounds on the remaining parameters can be problematic
for some optimization algorithms, so it's usually good practice to set large
but finite limits for those values.  The set_cap method
can do just that, setting a minimum and maximum value for all the parameters
that otherwise have bounds outside the cap.
End of explanation
m.maximize_loglike()
Explanation: To actually develop maximum likelihood estimates for the remaining
unconstrained parameters, use the 
maximize_loglike method.
End of explanation
m.calculate_parameter_covariance()
m.pf
Explanation: In a Jupyter notebook, this method displays a live-updating view of the
progress of the optmization algorithm, so that the analyst can interrupt
if something looks wrong.
The maximize_loglike method does 
not include the calculation of parameter covariance matrixes, standard
errors, or t-statistics.  For large models, this can be a computationally 
expensive process, and it is often but not always necessary. Those
computatations are made in the
calculate_parameter_covariance
method instead.  Once completed, things like t-statistics and standard
errors are available in the parameter frame.
End of explanation
m.parameter_summary()
m.estimation_statistics()
m.most_recent_estimation_result
Explanation: Reporting
Larch includes a variety of pre-packaged and a la carte reporting options.
Commonly used report tables are available directly in a Jupyter notebook
through a selection of reporting functions.
End of explanation
m.to_xlsx("/tmp/larch-demo.xlsx")
Explanation: To save a model report to an Excel file, use the to_xlsx method.
End of explanation
<END_TASK>