1.3 DCM模型的Python实践
我们在1.1节、1.2节系统地学习了选择行为的经济学理论、DCM的设计原理以及相关数学知识，本节将基于Python对旅行出行方式选择案例进行模型搭建及解读，希望读者能快速行动起来，将理论结合实际，尽快学以致用。
1.3.1 软件包和数据格式
DCM在Python中有可以直接调用的软件包，其中LR模型可以使用statsmodels软件包，MNL模型、NL模型可以使用pylogit软件包。这里推荐读者先安装好Anaconda（一款非常流行的数据分析平台），安装完毕后，再手动安装statsmodels、pylogit软件包，直接在终端输入代码清单1-1所示的命令即可。
代码清单1-1 安装相关软件包
pip install statsmodels pip install pylogit
有过建模经验的读者都知道，提供模型训练/预测的数据要严格符合模型要求，否则模型会运行失败。这一点对DCM尤为重要，这里先用一定的篇幅着重介绍DCM常用的两种数据格式：宽格式和长格式。
本节的案例数据来自William Greene《微观经济学建模及离散选择分析》课程中的旅行模式选择数据（Travel Mode Choice Data）。为了方便使用，我们先对数据进行一定的处理，处理后的数据形式如表1-2所示。
表1-2 旅行模式选择数据的字典
大多数读者应该对宽格式比较熟悉，因为常用的机器学习模型的输入数据大多是宽格式。如表1-3所示，一行数据代表一次选择过程，字段表示选择结果，其他与选择有关的信息被平行放置在一行中的各字段。LR模型的输入数据格式就是宽格式。
表1-3 宽格式示例
对于长格式，读者可能就比较陌生了，不过它是多数Logit模型使用的数据格式。如表1-4所示，一次选择行为包含多行数据，一行数据代表一次选择过程中的一个选项，其中会有一个字段表示该选项是否被选中（选中为1、未选中为0，一次选择只有一个选项被选中），对于备选项属性数据可以不同，决策者属性数据则相同。
表1-4 长格式示例
pylogit提供了长/宽数据格式相互转换的函数，操作实例如代码清单1-2所示。
代码清单1-2 使用pylogit进行长/宽数据格式相互转换
from collections import OrderedDict # OrderedDict用于记录模型的specification（声明）
import pylogit as pl # 引入Logit模型软件包pylogit
# 数据读入
long_data_path = u'long_data.csv'
long_df.to_csv(long_data_path, sep=',',index = False)
#---------------------#
# “长格式”转换为“宽格式”#
#---------------------#
# 指定决策者属性的列表
individual_specific_variables = ["HINC","PSIZE"]
# 指定备选项属性的列表
alternative_specific_variables = ['TTME', 'INVC', 'INVT', 'GC']
# 指定备选项属性的特殊说明列表
subset_specific_variables = {}
# “观测ID”，标识每次选择
observation_id_column = "OBS_ID"
# “备选项ID”，标识备选方案
alternative_id_column = "ALT_ID"
# “选择结果”，标识选择结果
choice_column = "MODE"
# 可选变量，记录与每个备选方案对应的名称，并允许在宽格式数据中创建有意义的列名
alternative_name_dict = {0: "AIR",
1: "TRAIN",
2: "BUS",
3: "CAR"}
wide_df = pl.convert_long_to_wide(long_df,
individual_specific_variables,
alternative_specific_variables,
subset_specific_variables,
observation_id_column,
alternative_id_column,
choice_column,
alternative_name_dict)
#---------------------#
# “宽格式”转换为“长格式”#
#---------------------#
# 创建决策者变量的列表
ind_variables =["HINC","PSIZE"]
# 指定每个备选项的属性所对应的字段，0/1/2/3代表备选项，后面的字段名为其属性在“宽格式数据”中的列名
alt_varying_variables = {"TTME": dict([(0, 'TTME_AIR'),
(1, 'TTME_TRAIN'),
(2, 'TTME_BUS'),
(3, 'TTME_CAR')]),
"INVC": dict([(0, 'INVC_AIR'),
(1, 'INVC_TRAIN'),
(2, 'INVC_BUS'),
(3, 'INVC_CAR')]),
"INVT": dict([(0, 'INVT_AIR'),
(1, 'INVT_TRAIN'),
(2, 'INVT_BUS'),
(3, 'INVT_CAR')]),
"GC": dict([(0, 'GC_AIR'),
(1, 'GC_TRAIN'),
(2, 'GC_BUS'),
(3, 'GC_CAR')])}
# 指定可用性变量，字典的键为可选项的ID，值为数据集中标记可用性的列
# 由于篇幅所限，前面的宽格式数据中省略了这3列，实际在做数据转化时需要标识每个选项的可用性
availability_variables = {0: 'availability_AIR',
1: 'availability_TRAIN',
2: 'availability_BUS',
3: 'availability_CAR'}
# “备选项ID”标识与每一行相关联的备选方案
custom_alt_id = "ALT_ID"
# “观测ID”标识每次选择，观测ID需要从1开始
obs_id_column = "OBS_ID"
wide_df[obs_id_column] = np.arange(wide_df.shape[0], dtype=int) + 1
# “选择结果”标识选择结果
choice_column = 'MODE'
# 执行长格式转换
long_df = pl.convert_wide_to_long(wide_df,
ind_variables,
alt_varying_variables,
availability_variables,
obs_id_column,
choice_column,
new_alt_id_name=custom_alt_id)
1.3.2 使用逻辑回归分析自驾选择问题
基于前文的介绍，相信读者已经迫不及待想使用MNL模型或NL模型进行建模分析了，这里先从LR模型的实操讲起。LR模型是目前应用最广泛的可解释二分类模型之一，深入了解LR模型对我们的日常工作有很大帮助。
通过对案例数据进行一定的处理，可以得到一份满足LR模型要求的宽格式数据。具体数据描述如表1-5所示，场景逻辑如图1-5所示。
表1-5 LR模型训练数据的字典
图1-5 LR模型的场景逻辑示意图
了解数据形式后，开始进行具体的模型搭建工作。
第1步：引入软件包，读取数据。重要的软件包在代码的备注中，如代码清单1-3所示。
代码清单1-3 引入软件包及读取数据
import numpy as np # 引入基础软件包numpy import pandas as pd # 引入基础软件包pandas import statsmodels.api as sm # 引入Logistic regression软件包statsmodels from sklearn.model_selection import train_test_split # 引入训练集/测试集构造工具包 from sklearn import metrics # 引入模型评价指标AUC计算工具包 import matplotlib.pyplot as plt # 引入绘图软件包 import scipy # 引入scipy软件包完成卡方检验 # 数据读入 data_path = 'wide_data.csv' raw_data = pd.read_table(data_path, sep=',', header=0)
第2步：数据预处理。数据预处理工作对于任何模型搭建都是必要的，这里结合LR模型及后续将介绍的MNL模型、NL模型的特点着重讲3个数据预处理的要点：①不要存在缺失值；②每一列数据均为数值型；③多枚举值离散变量输入模型前要进行哑变量处理，如代码清单1-4所示。
代码清单1-4 数据预处理
# 1. 缺失值探查&简单处理 model_data.info() # 查看每一列的数据类型和数值缺失情况 # | RangeIndex: 210 entries, 0 to 209 # | Data columns (total 9 columns): # | ... # | HINC 210 non-null int64 # | ... model_data = model_data.dropna() # 缺失值处理—删除 model_data = model_data.fillna(0) # 缺失值处理—填充（零、均值、中位数、预测值等） # 2. 数值型核查（连续变量应为int64或float数据类型） # 若上一步中存在应为连续数值变量的字段为object，则执行下列代码，这里假设'HINC'存在为字符串'null'的值 import re # 正则表达式工具包 float_patten = '^(-?\\d+)(\\.\\d+)?$' # 定义浮点数正则patten float_re = re.compile(float_patten) # 编译 model_data['HINC'][model_data['HINC'].apply(lambda x : 'not_float' if float_re.match(str(x)) == None else 'float') == 'not_float'] # 查看非浮点型数据 # | 2 null # | Name: distance, dtype: object model_data = model_data[model_data['HINC'] != 'null'] model_data['HINC'] = model_data['HINC'].astype(float)
第3步：单变量分析。在建模之前需要对每个自变量进行单变量分析，确定是否纳入模型。变量分为离散变量和连续变量两种，分析方式也有所不同。对于离散变量，我们使用k-1自由度的卡方检验，其中k为离散变量的值个数；对于连续变量，比较简单的分析方法是直接对单变量进行逻辑回归，查看回归系数的显著性，根据AUC分析自变量对y的解释能力。保留显著的自变量进入后续的操作，如代码清单1-5所示。
代码清单1-5 单变量分析
# 离散变量分析
crosstab = pd.crosstab( model_data['y'],model_data['PSIZE'])
p=scipy.stats.chi2_contingency(crosstab)[1]
print("PSIZE:",p)
# PSIZE: 0.0024577358937625327
# 连续变量分析
logistic = sm.Logit(model_data['y'],model_data['INVT_CAR']).fit()
p = logistic.pvalues['INVT_CAR']
y_predict = logistic.predict(model_data['INVT_CAR'])
AUC = metrics.roc_auc_score(model_data['y'],y_predict)
result = 'INVT_CAR:'+str(p)+' AUC:'+str(AUC)
print(result)
# INVT_CAR:2.971604856310474e-09 AUC:0.6242563699629587
第4步：共线性检验。由于LR模型是一种广义线性模型，变量间严重的共线性会对参数估计的准确性及泛化能力产生影响，因此需要对自变量间的共线性进行分析。若vif值大于10，可认为变量间具有很强的共线性，需要进行相应的处理，最简单的处理方式就是剔除自变量，保留单变量分析中AUC最大的变量。共线性检验示例如代码清单1-6所示。
代码清单1-6 共线性检验
from statsmodels.stats.outliers_influence import variance_inflation_factor
#共线性诊断包
X = raw_data[[ 'INVT_AIR', 'INVT_TRAIN','INVT_BUS', 'INVT_CAR']]
vif = pd.DataFrame()
vif['VIF Factor'] = [variance_inflation_factor(X.values, i) for i in range(X.shape[1])]
vif['features'] = X.columns
print('================多重共线性==============')
print(vif)
# | 0 14.229424 INVT_AIR
# | 1 72.782420 INVT_TRAIN
# | 2 80.279742 INVT_BUS
# | 3 35.003438 INVT_CAR
第5步：模型搭建。这里需要注意的是，对于3值及以上的离散变量要进行哑变量处理（需要记住去掉的枚举值），并且增加截距项Intercept，同时进行训练集和测试集的拆分（目的是防止模型过拟合，确定分析结论可以泛化），代码如清单1-7所示。
代码清单1-7 搭建LR模型
# 建模数据构造
X = model_data[[ 'HINC','PSIZE','TTME_TRAIN' , 'INVC_CAR']]
y = raw_data['y']
# 哑变量处理
dummies = pd.get_dummies(X['PSIZE'], drop_first=False)
dummies.columns = [ 'PSIZE'+'_'+str(x) for x in dummies.columns.values]
X = pd.concat([X, dummies], axis=1)
X = X.drop('PSIZE',axis=1) # 删去原离散变量
X = X.drop('PSIZE_4',axis=1) # 删去过于稀疏的字段
X = X.drop('PSIZE_5',axis=1) # 删去过于稀疏的字段
X = X.drop('PSIZE_6',axis=1) # 删去过于稀疏的字段
X['Intercept'] = 1 # 增加截距项
# 训练集与测试集的比例分别为80%和20%
X_train, X_test, y_train, y_test = train_test_split(X, y, train_size = 0.8, random_state=1234)
# 建模
logistic = sm.Logit(y_train,X_train).fit()
print(logistic.summary2())
# 重要返回信息
# | ------------------------------------------------------------------
# | Coef. Std.Err. z P>|z| [0.025 0.975]
# | ------------------------------------------------------------------
# | HINC 0.0264 0.0100 2.6477 0.0081 0.0068 0.0459
# | TTME_TRAIN 0.0389 0.0195 1.9916 0.0464 0.0006 0.0772
# | INVC_CAR -0.0512 0.0204 -2.5103 0.0121 -0.0913 -0.0112
# | PSIZE_1 -0.3077 0.7317 -0.4206 0.6741 -1.7419 1.1264
# | PSIZE_2 -1.0800 0.6417 -1.6829 0.0924 -2.3378 0.1778
# | PSIZE_3 -0.7585 0.7582 -1.0004 0.3171 -2.2444 0.7275
# | Intercept -1.8879 1.1138 -1.6951 0.0901 -4.0708 0.2950
# | =================================================================
# 模型评价
print("========训练集AUC========")
y_train_predict = logistic.predict(X_train)
print(metrics.roc_auc_score(y_train,y_train_predict))
print("========测试集AUC========")
y_test_predict = logistic.predict(X_test)
print(metrics.roc_auc_score(y_test,y_test_predict))
# | ========训练集AUC========
# | 0.7533854166666667
# | ========测试集AUC========
# | 0.6510263929618768
第6步：模型修正。可以看到，由于不显著变量的影响，模型的测试集AUC与训练集AUC存在较大差异，我们需要对不显著变量进行剔除。可以看到，新建模型的拟合优度尚可（AUC接近0.75），且自变量显著（p < 0.05），可以进行后续解读，如代码清单1-8所示。
代码清单1-8 修正LR模型
X = X.drop('PSIZE_1',axis=1)
X = X.drop('PSIZE_2',axis=1)
X = X.drop('PSIZE_3',axis=1)
# 训练集与测试集的比例分别为80%和20%
X_train, X_test, y_train, y_test = train_test_split(X, y, train_size = 0.8, random_state=1234)
# 建模
logistic = sm.Logit(y_train,X_train).fit()
print(logistic.summary2())
# 重要返回信息
# | ------------------------------------------------------------------
# | Coef. Std.Err. z P>|z| [0.025 0.975]
# | ------------------------------------------------------------------
# | HINC 0.0266 0.0096 2.7731 0.0056 0.0078 0.0454
# | TTME_TRAIN 0.0335 0.0161 2.0838 0.0372 0.0020 0.0650
# | INVC_CAR -0.0450 0.0168 -2.6805 0.0074 -0.0778 -0.0121
# | Intercept -2.3486 0.8275 -2.8384 0.0045 -3.9704 -0.7269
# | =================================================================
print("========训练集AUC========")
y_train_predict = logistic.predict(X_train)
print(metrics.roc_auc_score(y_train,y_train_predict))
print("========测试集AUC========")
y_test_predict = logistic.predict(X_test)
print(metrics.roc_auc_score(y_test,y_test_predict))
# | ========训练集AUC========
# | 0.7344618055555555
# | ========测试集AUC========
# | 0.7419354838709677
第7步：模型解读。DCM模型解读的对象可以分为概率（probability）和几率（odds）。在本例中，概率为“选择自驾的概率”，几率为“选择自驾的概率/不选择自驾的概率”。限于模型的数学性质，我们无法直接从模型参数中快速得到概率，而是需要经过一定的计算，这部分会在介绍复杂MNL模型及NL模型时展示。
得益于LR模型的数学性质，数据分析师可以基于模型参数直接对几率进行解读（这一点类似于线性回归）。模型解读的话术为“在其他条件保持不变的情况下，某因素增长一个单位（或属性a相对属性b），几率会变化（增长或降低）多少”，计算公式如下。
连续变量：odd(xi + 1)/odd(xi)-1=exp(βi)-1
离散变量：odd(xj = 1)/odd(xj=0)-1=exp(βj)-1
例如，根据模型可知：
在其他条件保持不变的情况下，家庭收入增长1个单位，选择自驾的odds会变化，exp(βHINC)-1=exp(0.0266)-1=0.027，即增加0.027倍。
在其他条件保持不变的情况下，自驾成本上升1个单位，选择自驾的odds会变化，exp(βINVC_CAR)-1=exp(-0.0450)-1=-0.044，即下降0.044倍。
1.3.3 使用多项Logit模型分析多种交通方式选择问题
1.3.2节我们使用逻辑回归模型分析了是否选择自驾的二项选择问题，如果要同时分析4种交通方式的选择问题，则需要使用MNL模型或NL模型。本节将介绍基于IIA假定的MNL模型，模型的问题场景映射如图1-6所示。
图1-6 MNL模型的场景逻辑示意图
需要注意的是，MNL模型的输入数据为长格式。不同于LR模型，MNL模型需要更加详细、复杂的初始化声明，以指定每种选项的效用函数形式。为了保证信息的完整性，尽量先保留自变量，定义如下模型。
Vair = ASCair + βttmeTTME + βinvcINVC + βinvtINVT + βhinc,airHINC + βpsize,airPSIZE
Vtrain = ASCtrain + βttmeTTME + βinvcINVC + βinvtINVT + βhinc,trainHINC + βpsize,trainPSIZE
Vbus = ASCbus + βttmeTTME + βinvcINVC + βinvtINVT + βhinc,busHINC + βpsize,busPSIZE
Vcar = βinvcINVC + βinvtINVT
根据设计好的模型结构搭建模型。1.3.2节已经介绍了数据清洗，这里不再赘述，直接进入模型搭建的代码学习，如代码清单1-9所示。
代码清单1-9 搭建MNL模型
# 第一步：模型初始化声明
basic_specification = OrderedDict()
basic_names = OrderedDict()
# 注意截距项包含选项个数减1
basic_specification["intercept"] = [0, 1, 2]
basic_names["intercept"] = ['ASC_air', 'ASC_train', 'ASC_bus']
# 可以灵活指定备选项属性的影响方式
basic_specification["TTME"] = [[0, 1, 2]]
basic_names["TTME"] = ['TTME']
basic_specification["INVC"] = [[0, 1, 2, 3]]
basic_names["INVC"] = ['INVC']
basic_specification["INVT"] = [[0, 1, 2, 3]]
basic_names["INVT"] = ['INVT']
# 也可以灵活指定决策者的影响方式，但需要注意的是，由于每个选项的决策者属性都一样，因此保证可估计性只对部分选项生效
basic_specification["HINC"] = [0, 1, 2]
basic_names["HINC"] = ['HINC_air', 'HINC_train', 'HINC_bus']
basic_specification["PSIZE"] = [0, 1, 2]
basic_names["PSIZE"] = ['PSIZE_air', 'PSIZE_train', 'PSIZE_bus']
# 第二步：创建模型
mnl = pl.create_choice_model(data = model_data,
alt_id_col="ALT_ID",
obs_id_col="OBS_ID",
choice_col="MODE",
specification=basic_specification,
model_type = "MNL",
names=basic_names)
# 第三步：模型估计和模型结果
mnl.fit_mle(np.zeros(12)) # 需要输入模型参数数量，根据之前的模型表达式即可得到
mnl.get_statsmodels_summary()
# | -------------------------------------------------------------
# | coef std.err z P>|z| [0.025 0.975]
# | -------------------------------------------------------------
# | ASC_air 6.0352 1.138 5.302 0.000 3.804 8.266
# | ASC_train 5.5735 0.711 7.836 0.000 4.179 6.968
# | ASC_bus 4.5047 0.796 5.661 0.000 2.945 6.064
# | TTME -0.1012 0.011 -9.081 0.000 -0.123 -0.079
# | INVC -0.0087 0.008 -1.101 0.271 -0.024 0.007
# | INVT -0.0041 0.001 -4.627 0.000 -0.006 -0.002
# | HINC_air 0.0075 0.013 0.567 0.571 -0.018 0.033
# | HINC_train -0.0592 0.015 -3.977 0.000 -0.088 -0.03
# | HINC_bus -0.0209 0.016 -1.278 0.201 -0.053 0.011
# | PSIZE_air -0.9224 0.259 -3.568 0.000 -1.429 -0.416
# | PSIZE_train 0.2163 0.234 0.926 0.355 -0.242 0.674
# | PSIZE_bus -0.1479 0.343 -0.432 0.666 -0.820 0.524
# |==============================================================
模型的搭建完成后，我们会发现有些变量不显著，此时需要进行模型的修正，如代码清单1-10所示。这里受篇幅限制，主要使用属性剔除及属性影响合并的方式进行修正，修正后的模型声明及模型效果如下。
Vair = ASCair + βttmeTTME + βinvtINVT + βpsize,airPSIZE
Vtrain = ASCtrain + βttmeTTME + βinvtINVT + βhinc,trainHINC
Vbus = ASCbus + βttmeTTME + βinvtINVT + βhinc,busHINC
Vcar = βinvtINVT
代码清单1-10 修正MNL模型
basic_specification = OrderedDict()
basic_names = OrderedDict()
basic_specification["intercept"] = [0, 1, 2]
basic_names["intercept"] = ['ASC_air', 'ASC_train', 'ASC_bus']
basic_specification["TTME"] = [[0, 1, 2]]
basic_names["TTME"] = ['TTME']
basic_specification["INVT"] = [[0, 1, 2, 3]]
basic_names["INVT"] = ['INVT']
basic_specification["HINC"] = [[1, 2]]
basic_names["HINC"] = [ 'HINC_train_bus']
basic_specification["PSIZE"] = [0]
basic_names["PSIZE"] = ['PSIZE_air']
mnl = pl.create_choice_model(data = model_data,
alt_id_col="ALT_ID",
obs_id_col="OBS_ID",
choice_col="MODE",
specification=basic_specification,
model_type = "MNL",
names=basic_names)
mnl.fit_mle(np.zeros(7))
mnl.get_statsmodels_summary()
# | -----------------------------------------------------------------
# | coef std.err z P>|z| [0.025 0.975]
# | -----------------------------------------------------------------
# | ASC_air 5.6860 0.937 6.068 0.000 3.849 7.523
# | ASC_train 5.4034 0.603 8.959 0.000 4.221 6.585
# | ASC_bus 5.0128 0.623 8.051 0.000 3.792 6.233
# | TTME -0.0992 0.011 -9.428 0.000 -0.12 -0.079
# | INVT -0.0039 0.001 -4.489 0.000 -0.006 -0.002
# | HINC_train_bus -0.0500 0.011 -4.484 0.000 -0.072 -0.028
# | PSIZE_air -0.8997 0.245 -3.680 0.000 -1.379 -0.420
# |==================================================================
根据模型系数可以初步判定模型的合理性，例如：TTME的系数为负，可以解释为当某个备选项站点等待时间延长，其被选择的概率会降低；HINC_train_bus的系数为负，可以解释为随着家庭收入增加，选择火车或长途汽车的概率会降低。这种定性的合理性判断有利于我们判断模型搭建是否合理。当然，如果想发挥模型真正的价值，还需要对模型进行量化解读。
对MNL模型的解读需要基于其预测功能，原理已经在1.2.3节进行了阐述，这里主要进行代码实现，如代码清单1-11所示。假设其他条件保持不变，因为火车提速，使得行程耗时降低20%，通过计算可知：
- 飞机的选择概率会由27.6%变为25.6%，降低了2.0%。
- 火车的选择概率会由30.0%变为36.2%，提升了6.2%。
- 长途汽车的选择概率会由14.3%变为12.7%，降低1.6%。
- 自驾的选择概率会由28.1%变为25.5%，降低2.6%。
代码清单1-11 解读MNL模型
# 创建用于预测的df
prediction_df = model_data[['OBS_ID', 'ALT_ID', 'MODE','TTME', 'INVT','HINC','PSIZE']]
choice_column = "MODE"
# 对火车耗时进行变化
def INVT(x,y):
if x == 1:
return y*0.8
else:
return y
prediction_df['INVT'] = prediction_df.apply(lambda x: INVT(x.ALT_ID, x.INVT), axis = 1)
# 默认情况下，predict()方法返回的结果是每个备选方案的选择概率
prediction_array = mnl.predict(prediction_df)
# 存储预测概率
prediction_df["MNL_Predictions"] = prediction_array
# 对比变化前后的概率
raw_probability = prediction_df.groupby(['ALT_ID'])['MODE'].mean()
new_probability = prediction_df.groupby(['ALT_ID'])['MNL_Predictions'].mean()
print("--------原概率--------")
print(raw_probability)
print("--------新概率--------")
print(new_probability)
# | --------原概率--------
# | ALT_ID
# | 0 0.276190
# | 1 0.300000
# | 2 0.142857
# | 3 0.280952
# | Name: MODE, dtype: float64
# | --------新概率--------
# | ALT_ID
# | 0 0.255643
# | 1 0.362788
# | 2 0.126937
# | 3 0.254632
1.3.4 使用嵌套Logit模型分析多种交通方式选择问题
本节我们要学习的是嵌套Logit模型，顾名思义，嵌套Logit模型就是人为设置选择的层次。我们假定选择的过程存在层次性，如图1-7所示，旅行家庭先进行飞机或陆地交通的选择，如果选择陆地交通，再决定到底是选火车、长途汽车还是自驾。该模型主要应用于不满足IIA假定的分析场景。
图1-7 NL模型的场景逻辑示意图
我们延续MNL模型的代码进行阐述，如代码清单1-12所示，NL模型与MNL模型的不同之处在于NL模型需要在声明每个备选项的效用表达式之前，单独声明NL模型的层次。
代码清单1-12 NL模型代码
# 声明嵌套形式
nest_membership = OrderedDict()
nest_membership["air_Modes"] = [0]
nest_membership["ground_Modes"] = [1, 2, 3]
# 声明备选项的效用函数
basic_specification = OrderedDict()
basic_names = OrderedDict()
basic_specification["intercept"] = [0, 1, 2]
basic_names["intercept"] = ['ASC_air', 'ASC_train', 'ASC_bus']
# 可以灵活指定备选项属性的影响方式
basic_specification["TTME"] = [[0, 1, 2]]
basic_names["TTME"] = ['TTME']
basic_specification["INVT"] = [[0, 1, 2, 3]]
basic_names["INVT"] = ['INVT']
# 也可以灵活指定决策者的影响方式，但需要注意的是，由于每个选项的决策者属性都一样，因此保证可估计性，只对部分选项生效
basic_specification["HINC"] = [[1, 2]]
basic_names["HINC"] = [ 'HINC_train_bus']
basic_specification["PSIZE"] = [0]
basic_names["PSIZE"] = ['PSIZE_air']
# 模型创建
nested_logit = pl.create_choice_model(data = model_data,
alt_id_col="ALT_ID",
obs_id_col="OBS_ID",
choice_col="MODE",
specification=basic_specification,
model_type = "Nested Logit",
names=basic_names,
nest_spec=nest_membership)
nested_logit.fit_mle(np.zeros(9))
nested_logit.summary
# | -----------------------------------------------------------------
# | parameters std_err t_stats p_values
# | -----------------------------------------------------------------
# | air_Modes 0.0000 NaN NaN NaN
# | ground_Modes 0.8187 0.668 1.225 0.220
# | ASC_air 4.0002 1.022 3.914 0.000
# | ASC_train 4.2224 0.744 5.672 0.000
# | ASC_bus 3.9471 0.747 5.285 0.000
# | TTME -0.0787 0.013 -6.180 0.000
# | INVT -0.0038 0.001 -4.718 0.000
# | HINC_train_bus -0.0364 0.010 -3.513 0.000
# | PSIZE_air -0.7535 0.244 -3.088 0.002
# |==================================================================
从模型结果来看，ground_Modes的系数并不显著，并且尝试用其他层次条件进行建模亦会发现，层次的系数也不显著，这从侧面表明了MNL模型在该问题上的应用是合理的。当然，多数情况下IIA假定并不一定被满足，NL模型会比MNL模型更科学和严谨。
NL模型的解读与MNL模型一致，这里不再赘述。