09离散因变量模型
离散因变量模型
要考察人们做出某种具体选择的情况及其影响因素时，可把这些离散的定性变量作为因变量进行分析，把影响因素作为自变量，这样建立的模型称之为离散选择模型。如出行交通工具选择的情况。
还有一种是因变量是以离散计数的方式描述的，分析自变量对计数因变量的影响所建立的模型，称之为计数模型。如发生交通事故的次数。
线性概率模型
离散选择模型在广义线性模型（generalized linear model）的框架下展开，并依赖结果是两个或多个选择将模型分位二项选择、多项选择模型和受限因变量模型
离散选择模型主要研究选择结果的概率与影响因素之间的关系，即
\[Prob(事件i发生) = Prob(Y=i)=F(影响因素)
\]
其中，影响因素可能包含做出选择的主体属性和选择方案属性。如选择何种交通工具出行，既受到选择主体收入程度、生活习惯等属性的影响，也收到交通工具的价格、便捷性等属性的影响。
示例：对影响手机购买意向的因素进行分析
购买意向为定性变量，有两种选择：0表示不购买，1表示购买。其影响因素可能有性别、年龄、收入、职位、行业等诸多因素。
设因变量 \(y\) 表示是否购买手机，则有
\[y= \begin{cases}
0 & 不购买 \\
1 & 购买
\end{cases}
\]
影响 \(y\) 的因素记为 \(x=(x_1,x_2,\cdots, x_n)\)，根据多元回归的思想，可得
\[y = \beta_0 + \beta_1 x_1+\beta_2 x_2+\cdots +\beta_n x_n + \varepsilon
\]
其中，\((\beta_1,\beta_1,\cdots, \beta_n)^T=\beta\) 表示回归模型中的参数即回归系数，则简化为
\[y = \beta_0 + \beta x + \varepsilon
\]
在因变量是离散变量的情况下，不能把 \(\beta_i(i=1,2,\cdots,n)\) 理解为保持其他因素不变的情况下对 \(y\) 的边际影响，因为 \(y\) 的取值为1或0。
回忆经典回归分析中 \(E(\varepsilon|x)=0\) 的假定，对回归模型左右两边求条件期望，则有
\[E(y|x) = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \cdots +\beta_n x_n
\]
因为 \(y\) 的取值为1或0，记 \(y=1\) 的概率为 \(p\)，依数学期望的定义有
\[E(y) = 0 \times (1-p) + 1 \times p = p = Prob(y=1|x)
\]
因此，做出购买决策的概率等于 \(y\) 的期望值，联立上式有
\[Prob(y=1|x) = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \cdots +\beta_n x_n
\]
做出购买觉得的概率 \(p(y=1|x)\) 是 \(x_i\) 的一个线性函数，该回归方程亦即线性概率模型（LPM, linear probability model），\(p(y=1|x)\) 可称之为响应概率。用普通最小二乘法对模型进行参数估计，可得如下回归方程
\[Prob(y=1|x) = \hat{y} =\hat{\beta}_0 + \hat{\beta}_1 x_1 + \hat{\beta}_2 x_2 + \cdots +
\hat{\beta}_n x_n = \hat{\beta}_0 + \hat{\beta} x
\]
\(\hat{y}\)则是 \(y=1\) 时的概率，在保持其他因素不变的情况下，\(\beta_1\) 表示因 \(x_i\) 的变化导致购买(\(y=1\))决策概率的变化，即
\[\Delta P(y=1|x)=\beta_i \Delta x_i
\]
二元选择模型
二元选择模型(binary choice model)具体是指离散因变量具有两个选项或两种属性。二元因变量的属性往往是对立或互斥的。一般可用1和0来描述。可用回归模型来描述这些变量作为因变量时与其他自变量之间的关系，因此二元选择模型又可称之为0-1因变量回归模型。
在应用二元选择模型时需注意，因变量中的0和1只是对应属性的标注或符号，不具备任何数值上的意义，不能进行直接数学运算。
线性概率模型的缺陷与改进
线性概率模型虽然能够测度出因变量事件发生的响应概率，可以考察出因自变量变动导致因变量决策的概率变化。但该模型存在如下几个比较严重的缺陷：
1.解释变量的合理变化会导致概率预测值溢出在[0,1]区间之外
2.随机误差项的分布未知
3.模型误差项具有异方差性，异方差性使参数估计不具有有效性
4.即使用加权最小二乘法修正异方差性也无法保证概率预测值在[0,1]区间之间
随机误差项分布的缺陷较容易克服，当样本量充分大时，其普通最小二乘参数估计量的结果近似服从正态分布。但是其他的问题需要对变量之间的线性关系进行变换，使得自变量 \(x\) 对应的所有概率值都在[0,1]区间范围内，且自变量 \(x\) 变化时，\(y\) 单调变化。
依据上述要求，最为常见且符合要求的形式便是利用分布函数 \(F(.)\) 进行变换。因为分布函数对于所有的自变量而言，其值均在区间[0,1]范围之内，且当自变量发生变化时，分布函数值单调变化。
二元选择模型的基本原理
依据分布函数的基本特征，二元选择模型可从一个满足经典线性模型假定的隐变量模型中推导出来。所谓隐变量便是不能够直接进行观测，但可以通过其他直接观测得到的变量（显变量）进行描述和反映的变量。
- 模型构建和参数估计过程
设 \(y^*\) 是一个由
\[y^* = \beta_0 + x\beta + \varepsilon^*\\
y = \begin{cases}
1 ,& y^*>0 \\
0 ,& y^* \leq 0
\end{cases}
\]
所决定的不能够直接观测得到的隐变量，隐变量通常是解释变量或影响因素的线性函数；且假定 \(\varepsilon^*,x\)相互独立，\(\varepsilon^*\) 服从某种对称于0的分布。则 \(y^* = \beta_0 + x\beta + \varepsilon^*\) 也可称之为隐变量模型。
记 \(F\) 为分布函数或链接函数(link fucntion)，因为 \(\varepsilon^*\) 假定对称于0，则有 \(F(z)=1-F(-z)\)，其中 \(z\) 是随机变量。
依据 \(y^*\) 与 \(y\) 的对应关系和对称分布函数性质，则 \(y\) 的响应概率为
\[\begin{align*}
P(y=1|x)
&= P(y^*>0|x) \\
&= P((\beta_0 + x\beta) + \varepsilon^*>0|x) \\
&= P(\varepsilon^* >-(\beta_0 + x\beta)|x) \\
&= P(\varepsilon^* < (\beta_0 + x\beta)|x) \\
&= F(\beta_0+x\beta) \\
&= 1 - F(-(\beta_0+x\beta))
\end{align*}
\]
如果已知 \(\varepsilon^*\) 分布函数 \(F(z)\) 的具体形式，那么决定 \(y\) 响应概率的模型便确定了。
在二元选择模型中，通常\(\varepsilon^*\) 服从的分布如下所示的两种最为常见的形式
| \(\varepsilon^*\)的分布 | 分布函数 \(F(z)\) | 二元选择模型 | 模型具体形式 |
|---|---|---|---|
| 标准正态分布 | \(\Phi(z)\) | Probit模型 | \(P(y=1|x)=F(\beta_0+x\beta)=\Phi(\beta_0+x\beta)\) |
| 逻辑分布 | \(\frac{e^z}{1+e^z}\) | Logit模型 | \(P(y=1|x)=F(\beta_0+x\beta)=\frac{e^{\beta_0+x\beta}}{1+ e^{\beta_0+x\beta}}\) |
\[\because E(y)=P(y=1|x)=F(\beta_0+x\beta) \\
\therefore y=E(y)+(y-E(y)) = E(y)+\varepsilon=F(\beta_0+x\beta)+\varepsilon
\]
二元选择模型中一般选用极大似然法进行参数估计，其似然函数为
\[L = \prod_{i=1}^{n}{\Big[F(\beta_0+x\beta)\Big]^{y_i}\Big[1-F(\beta_0+x\beta)\Big]^{1-y_i}}
\]
其对数似然函数为
\[\ln L = \sum_{i=1}^{n}\Big\{y_i\ln F(\beta_0+x\beta)+(1-y_i)\ln[1-F(\beta_0+x\beta)]\Big\}
\]
求对数似然函数对每个参数的偏导数，使之均为0，再求解方程组即可解得模型中的参数估计值。
对模型进行解读时注意：二元选择模型的参数估计结果不能理解为自变量变动对因变量的边际影响，而应当理解为自变量的变动，对因变量取1的概率的影响有多大，即
\[\frac{\partial P(y=1|x)}{\partial x_i}=\frac{\partial F(\beta_0+x\beta)}{\partial (\beta_0+x\beta)}\beta_i=f(\beta_0+x\beta)\beta_i
\]
其中，\(f(.)\)是分布函数 \(F(.)\) 的密度函数。
- 模型检验
二元选择模型与经典回归模型一样，同样可以用 \(Z\) 统计量对回归系数做显著性检验，该 \(Z\) 统计量由极大似然法给出。对于多个系数的约束显著性检验，还可以计算Wald统计量，利用Wald统计量近似服从 \(\chi^2\) 分布进行 \(\chi^2\) 检验。而对于模型的拟合优度检验，则可以利用 \(LR\) 似然比进行 \(\chi^2\) 检验。此外，还可以计算模型的AIC、BIC等信息指数在多个模型之间进行评价和模型选择。
Binary Probbit模型
二元（binary）Probbit模型对隐变量随机误差项假定服从标准正态分布，其模型具有如下形式
\[P(y=1|x)=F(\beta_0+x\beta)=\Phi(\beta_0+x\beta)
\]
示例
import pandas as pd
import statsmodels.api as sm
import scipy.stats as stats
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
# 对某产品做满意度调查，考察消费者的满意状况、对产品的抱怨状况、对产品购买使用的忠诚度状况，
# 以此对消费者是否会持续购买使用该产品进行分析
product = pd.read_csv('./data/product_usage.csv')
print(product.head(5))
# CSI Complaint Loyalty Attitude
# 0 8 3 8 1
# 1 9 4 10 1
# 2 6 1 10 1
# 3 8 4 7 1
# 4 8 7 7 1
# 本例中Attitude是调查研究对象，该变量有两个离散且互斥的值0和1，可以建立二元Probit模型进行分析
# 构建二元Probit模型
from statsmodels.formula.api import glm
formula = 'Attitude~CSI+Complaint+Loyalty'
product_m = glm(formula, data=product, family=sm.families.Binomial(sm.families.links.probit)).fit()
# glm建模过程中需要指定随机误差项的分布簇及其链接函数
# 常见分布簇有：Binomial,Gamma,Gaussian,InverseGaussian,NegativeBinomial,Possion等
# 其链接函数有：Probit,Logit,Log等15种转换形式
res = product_m.summary()
print(res)
# Generalized Linear Model Regression Results
# ==============================================================================
# Dep. Variable: Attitude No. Observations: 337
# Model: GLM Df Residuals: 333
# Model Family: Binomial Df Model: 3
# Link Function: probit Scale: 1.0000
# Method: IRLS Log-Likelihood: -134.37
# Date: Wed, 10 Jun 2020 Deviance: 268.74
# Time: 09:24:26 Pearson chi2: 395.
# No. Iterations: 7
# Covariance Type: nonrobust
# ==============================================================================
# coef std err z P>|z| [0.025 0.975]
# ------------------------------------------------------------------------------
# Intercept -1.7923 0.550 -3.260 0.001 -2.870 -0.715
# CSI 0.2479 0.065 3.792 0.000 0.120 0.376
# Complaint -0.1329 0.046 -2.890 0.004 -0.223 -0.043
# Loyalty 0.1984 0.044 4.512 0.000 0.112 0.285
# ==============================================================================
# 拟合优度检验：对数似然比(Log-Likelihood)和卡方统计量(Pearson chi2)比较法，任务拟合效果好
# 模型系数检验：z统计量进行检验，p值非常小，在\alpha=0.05的条件下，非常显著
# 解释性：CSI、Loyalty的系数为正，表明正向影响购买,Complaint的系数为负，表明反向影响购买
# 由于模型是通过链接函数建立的，所以估计出来的系数并不代表随着自变量变动，因变量概率概率变动的绝对数值。
# 但是可以通过参数估计之后建立的模型对现有的或新的观测样本进行考察或预测。
# 得出模型方程：P(y=1|x)=F(\beta_0+x\beta)=\Phi(\beta_0+x\beta)=
# \Phi(-1.7923 + 0.2479*CSI - 0.1329*Complaint + 0.1984*Loyalty)
# 预测
# CST=8, Complaint=4, Loyalty=7
res = product_m.predict(pd.DataFrame({'CSI': [8], 'Complaint': [4], 'Loyalty': [7]}))
print(res) # 0 0.8528
示例2
import pandas as pd
import statsmodels.api as sm
import scipy.stats as stats
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
# 针对不同性别、年龄的客户进行购买意向的因素分析，其中有定性变量作为自变量
G3 = pd.read_csv('./data/threeg.csv')
print(G3.head(5))
# Gender Age Purchase
# 0 1 35 0
# 1 1 59 1
# 2 1 56 1
# 3 1 46 1
# 4 1 39 0
# 本例中Purchase是调查研究对象，该变量有两个离散且互斥的值0和1，可以建立二元Probit模型进行分析
# 数据标签
G3['Gender'] = G3['Gender'].astype('category')
G3['Gender'].cat.categories = ['female', 'male']
G3['Gender'].cat.set_categories = ['female', 'male']
# 构建二元Probit模型
# 方法一：glm
from statsmodels.formula.api import glm
formula = 'Purchase~C(Gender)+Age'
G3_m = glm(formula, data=G3, family=sm.families.Binomial(sm.families.links.probit)).fit()
res = G3_m.summary()
print(res)
# Generalized Linear Model Regression Results
# ==============================================================================
# Dep. Variable: Purchase No. Observations: 65
# Model: GLM Df Residuals: 62
# Model Family: Binomial Df Model: 2
# Link Function: probit Scale: 1.0000
# Method: IRLS Log-Likelihood: -33.267
# Date: Wed, 10 Jun 2020 Deviance: 66.534
# Time: 10:18:10 Pearson chi2: 72.8
# No. Iterations: 5
# Covariance Type: nonrobust
# =====================================================================================
# coef std err z P>|z| [0.025 0.975]
# -------------------------------------------------------------------------------------
# Intercept -4.6363 1.186 -3.908 0.000 -6.961 -2.311
# C(Gender)[T.male] 0.7986 0.358 2.234 0.026 0.098 1.499
# Age 0.0962 0.025 3.789 0.000 0.046 0.146
# =====================================================================================
# Gender和Age做自变量对消费者购买的意愿影响显著；
# Intercept项估计值表示定性变量Gender为女时对响应概率对影响(Gender取女时，对应的参数估计值为0，故未显示)
# 从估计值符号得：男性比女性的购买概率要大
# 模型方程：P(y=1|x) = \Phi(-4.6363 + 0.7986*Gender_{male} + 0.0962*Age)
# 预测
# 男、女年龄均45
res = G3_m.predict(pd.DataFrame({'Gender': ['male', 'female'], 'Age': [45, 45]}))
print(res)
# 0 0.688636
# 1 0.379574
# 方法二：probit
from statsmodels.formula.api import probit
G3_m2 = probit(formula, data=G3).fit()
res = G3_m2.summary()
print(res)
# Optimization terminated successfully.
# Current function value: 0.511798
# Iterations 6
# Probit Regression Results
# ==============================================================================
# Dep. Variable: Purchase No. Observations: 65
# Model: Probit Df Residuals: 62
# Method: MLE Df Model: 2
# Date: Wed, 10 Jun 2020 Pseudo R-squ.: 0.2605
# Time: 10:28:05 Log-Likelihood: -33.267
# converged: True LL-Null: -44.985
# Covariance Type: nonrobust LLR p-value: 8.142e-06
# =====================================================================================
# coef std err z P>|z| [0.025 0.975]
# -------------------------------------------------------------------------------------
# Intercept -4.6363 1.168 -3.969 0.000 -6.926 -2.347
# C(Gender)[T.male] 0.7986 0.356 2.242 0.025 0.101 1.497
# Age 0.0962 0.025 3.810 0.000 0.047 0.146
# =====================================================================================
# 模型自变量对因变量的边际影响
res = G3_m2.get_margeff().summary()
print(res)
# Probit Marginal Effects
# =====================================
# Dep. Variable: Purchase
# Method: dydx
# At: overall
# =====================================================================================
# dy/dx std err z P>|z| [0.025 0.975]
# -------------------------------------------------------------------------------------
# C(Gender)[T.male] 0.2326 0.093 2.489 0.013 0.049 0.416
# Age 0.0280 0.005 6.027 0.000 0.019 0.037
# =====================================================================================
Binary Logit模型
二元（binary）Logit模型对隐变量随机误差项假定服从逻辑分布，其模型具有如下形式
\[P(y=1|x)=F(\beta_0+x\beta)=\frac{e^{\beta_0+x\beta}}{1+ e^{\beta_0+x\beta}}
\]
示例
import pandas as pd
import statsmodels.api as sm
import scipy.stats as stats
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
# 针对不同性别、年龄的客户进行购买意向的因素分析，其中有定性变量作为自变量
G3 = pd.read_csv('./data/threeg.csv')
print(G3.head(5))
# Gender Age Purchase
# 0 1 35 0
# 1 1 59 1
# 2 1 56 1
# 3 1 46 1
# 4 1 39 0
# 本例中Purchase是调查研究对象，该变量有两个离散且互斥的值0和1，可以建立二元Logit模型进行分析
# 数据标签
G3['Gender'] = G3['Gender'].astype('category')
G3['Gender'].cat.categories = ['female', 'male']
G3['Gender'].cat.set_categories = ['female', 'male']
# 构建二元Logit模型
# 方法一：glm
from statsmodels.formula.api import glm
formula = 'Purchase~C(Gender)+Age'
G3_m = glm(formula, data=G3, family=sm.families.Binomial(sm.families.links.logit)).fit()
res = G3_m.summary()
print(res)
# Generalized Linear Model Regression Results
# ==============================================================================
# Dep. Variable: Purchase No. Observations: 65
# Model: GLM Df Residuals: 62
# Model Family: Binomial Df Model: 2
# Link Function: logit Scale: 1.0000
# Method: IRLS Log-Likelihood: -33.110
# Date: Wed, 10 Jun 2020 Deviance: 66.219
# Time: 10:33:20 Pearson chi2: 71.6
# No. Iterations: 5
# Covariance Type: nonrobust
# =====================================================================================
# coef std err z P>|z| [0.025 0.975]
# -------------------------------------------------------------------------------------
# Intercept -8.0925 2.223 -3.640 0.000 -12.450 -3.735
# C(Gender)[T.male] 1.4311 0.628 2.279 0.023 0.200 2.662
# Age 0.1671 0.047 3.550 0.000 0.075 0.259
# =====================================================================================
# Gender和Age做自变量对消费者购买的意愿影响显著；
# Intercept项估计值表示定性变量Gender为女时对响应概率对影响(Gender取女时，对应的参数估计值为0，故未显示)
# 从估计值符号得：男性比女性的购买概率要大
# 模型方程：P(y=1|x) = F(\beta_0 + x\beta) = \frac{e^{\beta_0 + x\beta}}{1+e^{\beta_0 + x\beta}}
#
# 预测
# 男、女年龄均45
res = G3_m.predict(pd.DataFrame({'Gender': ['male', 'female'], 'Age': [45, 45]}))
print(res)
# 0 0.702411
# 1 0.360706
#
# 方法二：logit
from statsmodels.formula.api import logit
G3_m2 = logit(formula, data=G3).fit()
res = G3_m2.summary()
print(res)
# Optimization terminated successfully.
# Current function value: 0.509380
# Iterations 6
# Logit Regression Results
# ==============================================================================
# Dep. Variable: Purchase No. Observations: 65
# Model: Logit Df Residuals: 62
# Method: MLE Df Model: 2
# Date: Wed, 10 Jun 2020 Pseudo R-squ.: 0.2640
# Time: 11:48:51 Log-Likelihood: -33.110
# converged: True LL-Null: -44.985
# Covariance Type: nonrobust LLR p-value: 6.958e-06
# =====================================================================================
# coef std err z P>|z| [0.025 0.975]
# -------------------------------------------------------------------------------------
# Intercept -8.0925 2.223 -3.640 0.000 -12.450 -3.735
# C(Gender)[T.male] 1.4311 0.628 2.279 0.023 0.200 2.662
# Age 0.1671 0.047 3.550 0.000 0.075 0.259
# =====================================================================================
#
# 模型自变量对因变量的边际影响
res = G3_m2.get_margeff().summary()
print(res)
# Logit Marginal Effects
# =====================================
# Dep. Variable: Purchase
# Method: dydx
# At: overall
# =====================================================================================
# dy/dx std err z P>|z| [0.025 0.975]
# -------------------------------------------------------------------------------------
# C(Gender)[T.male] 0.2427 0.091 2.658 0.008 0.064 0.422
# Age 0.0283 0.005 6.153 0.000 0.019 0.037
# =====================================================================================
多重选择模型
对于因变量有多个选项是，可以用多重选择离散变量进行描述。对于分析具有多个属性的离散因变量与自变量之间的关系，可以利用多重选择模型（multiple choice）进行分析。
在多重选择模型中，不同选项或者选择结果之间的关系有两种情况。一种情况是各项选择是有顺序之分的，另一种情况是选项之间没有次序或顺序之分。
依据离散因变量选项的含义和次序不同，又可分为：顺序选择（ordinal）模型和无序选择（multinomial）模型。
对于顺序选择模型，有0、1、2、M种选择，则有
仍然设 \(y^*\) 是一个由
\[y^* = \beta_0 + x\beta + \varepsilon^*\\
y = \begin{cases}
0 & y^* < c_1 \\
1 & c_1 \leq y^* < c_2 \\
2 & c_2 \leq y^* < c_3 \\
\vdots \\
M & c_M \leq y^*
\end{cases}
\]
所决定的不能够直接观测得到的隐变量，\(c_1，c_2,\cdots,c_M\) 可称之为临界值或统称为截距项；\(\varepsilon_i^*\) 是独立同分布的随机变量，其分布函数为 \(F(z)\) ，则选择 \(j（j=0，1，2，\cdots,M)\) 方案的概率为
\[\begin{align*}
& P(y=0) = F(c_1-\alpha-x\beta)\\
& P(y=1) = F(c_2-\alpha-x\beta)-F(c_1-\alpha-x\beta)\\
& P(y=2) = F(c_3-\alpha-x\beta)-F(c_2-\alpha-x\beta)\\
& \cdots \\
& P(y=M) = 1 - F(c_M-\alpha-x\beta)
\end{align*}
\]
当分布函数 \(F(z)\) 为标准正态分布时，上述模型为ordinal Probit模型；当分布函数 \(F(z)\) 为逻辑分布时，称之为ordinal logit模型。
同理，无序选择模型也可以细分为 multinomial probit模型和multinomial logit模型。估计muttinomial probit模型的参数非常困难，在分析时较少采用该类模型。
对于无序选择模型，选择选项 \(j\) 的模型基本形式如 \(y_i^*=x_j\beta+\varepsilon^*\)
在影响因素 \(x_j\) 的作用下，选择选项 \(j\) 表示 \(y_j^*\) 在所欲选项中的效用时最大的。因此，选择选项 \(j\) 的概率模型为
\[P(y=j) = P(y_j^* > y_k^*), \forall k \neq j
\]
\(\varepsilon^*\) 是独立同分布的随机变量，一般假定其服从韦伯（weibull）分布，具体的分布函数形式为 \(F(z)=e^{-e^{-z}}\)，则有
\[P(y=j)=\frac{e^{x_j \beta}}{\sum_{k \neq j}{e^{x_k \beta}}}
\]
此即条件Logit模型。
对于多重选择模型，其估计和检验过程与二元选择模型类似。
此外，顺序选择模型中离散因变量的各种选择通常也是用阿拉伯数字来表示的。如对某个产品的喜好程度：1-喜欢，2-无所谓，3-不喜欢。这些数字之间只有顺序意义，没有数值上的意义。即这些所选择的数字仅能代表偏好的程度，所选择的数字代表之间不能够进行数值运算。
示例
import pandas as pd
import statsmodels.api as sm
import scipy.stats as stats
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
# 针对不同性别、年龄的客户进行购买意向的因素分析，将购买的意向细分为1-不购买，2-无所谓，3-购买
# 现假定购买意向Purchase的选项之间没有顺序关系，2中选择只反映购买态度，没有任何程度的偏好上的差异或顺序之分，
# 则该问题转换为一个无序选择问题
G3_m = pd.read_csv('./data/threeg_multi.csv')
print(G3_m.iloc[[0, 1, 28, 38, 58], :])
# gender age purchase
# 0 2 54 1
# 1 1 59 1
# 28 2 31 2
# 38 1 39 2
# 58 2 42 3
# 数据预处理，加值标签
G3_m['gender'] = G3_m['gender'].astype('category')
G3_m['gender'].cat.categories = ['female', 'male']
G3_m['gender'].cat.set_categories = ['female', 'male']
# G3_m['purchase'] = G3_m['purchase'].astype('category') # 编码后出错，暂未找到原因
# G3_m['purchase'].cat.categories = ['no buy', 'no sense', 'buy']
# G3_m['purchase'].cat.set_categories = ['no buy', 'no sense', 'buy']
# print(G3_m.head())
# 无序选择模型
from statsmodels.formula.api import mnlogit
formula = 'purchase~C(gender)+age'
G3_m_model = mnlogit(formula, data=G3_m).fit()
res = G3_m_model.summary()
print(res)
# Optimization terminated successfully.
# Current function value: 0.747045
# Iterations 7
# MNLogit Regression Results
# ==============================================================================
# Dep. Variable: purchase No. Observations: 65
# Model: MNLogit Df Residuals: 59
# Method: MLE Df Model: 4
# Date: Wed, 10 Jun 2020 Pseudo R-squ.: 0.2667
# Time: 14:20:49 Log-Likelihood: -48.558
# converged: True LL-Null: -66.217
# Covariance Type: nonrobust LLR p-value: 3.997e-07
# =====================================================================================
# purchase=2 coef std err z P>|z| [0.025 0.975]
# -------------------------------------------------------------------------------------
# Intercept 11.2468 3.207 3.507 0.000 4.961 17.533
# C(gender)[T.male] 1.5088 0.909 1.659 0.097 -0.274 3.291
# age -0.2343 0.067 -3.499 0.000 -0.366 -0.103
# -------------------------------------------------------------------------------------
# purchase=3 coef std err z P>|z| [0.025 0.975]
# -------------------------------------------------------------------------------------
# Intercept 10.3972 3.611 2.879 0.004 3.320 17.474
# C(gender)[T.male] 3.5062 1.107 3.168 0.002 1.337 5.675
# age -0.2634 0.078 -3.372 0.001 -0.416 -0.110
# =====================================================================================
# 如果因变量有k种选择，则在结果中只能得到k-1个响应概率，默认次序最小的选项对应的响应概率不显示。
# 分析过程与前面的类似；
# 预测
res = G3_m_model.predict(pd.DataFrame({'gender': ['male', 'female'], 'age': [45, 45]}))
print(res)
# 0 1 2
# 0 0.055831 0.509616 0.434553
# 1 0.307467 0.620711 0.071822
# 行分别与自变量gender的属性对应，列分别与因变量的3个属性对应
计数模型
当因变量表示时间发生的次数时，则因变量是一个离散形式的整数计数变量。这些变量都是以整数计数为表现形式的，分析这种离散计数因变量的模型即为计数模型（count model）
计数模型中的离散因变量的数字是有数值含义的，即次数之间是可以进行运算的。因此对于只能处理无数值意义的二元选择模型或多重选择模型而言，计数数据不适用。
对于某个时间、空间等范围之内发生次数的计数数据，一般都认为其近似服从泊松（Poisson）分布。
对于因变量和自变量之间的关系可以考虑建立回归模型。Poisson回归模型即考虑变量服从 Poisson 分布而建立一种回归模型，其假定因变量 \(y\) 服从参数为 \(\lambda\) 的泊松分布，则该模型的初始方程为
\[Prob(y=k)=\frac{e^{-\lambda}\lambda^k}{k!}, k=0,1,2,...
\]
其中，\(\lambda\) 表示所考察的事件在一定范围之内平均发生的次数；\(k\) 为整数，表示事件实际观测到的发生次数。
Poisson回归模型的研究可以从初始方程的唯一参数 \(\lambda\) 入手。所考察的事件在一定范围内平均发生的次数 \(\lambda\) 收到各种条件或影响因素 \(x\) 的影响，因此参数是变化的。经过对数变换，可以得到Poisson回归模型
\[\ln \lambda = \alpha + \beta_1 x_1 + \beta_2 x_2+\cdots + \beta_n x_n + \varepsilon = \alpha + x\beta + \varepsilon
\]
对于Poisson回归模型，可以使用极大似然法进行参数估计，其模型及参数估计值的检验问题类似前面。
示例
import pandas as pd
import statsmodels.api as sm
import scipy.stats as stats
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
# 对打印机发生故障的情况进行分析，发生故障次数可能受到：打印纸张数量、打印机使用时长、硒鼓（原装/兼容）等因素的影响。
printer = pd.read_csv('./data/printer.csv')
print(printer.head(5))
# cartridge Counts Pages Length
# 0 2 5 87.8 44.194
# 1 1 1 52.2 3.663
# 2 2 0 0.7 0.331
# 3 2 1 81.7 18.422
# 4 2 4 89.9 45.003
# 数据预处理：数值标签
printer['cartridge'] = printer['cartridge'].astype('category')
printer['cartridge'].cat.categories = ['原装', '兼容']
printer['cartridge'].cat.set_categories = ['原装', '兼容']
# 泊松回归
from statsmodels.formula.api import poisson
formula = 'Counts~C(cartridge)+Pages+Length'
printer_model = poisson(formula, data=printer).fit()
res = printer_model.summary()
print(res)
# Optimization terminated successfully.
# Current function value: 1.326793
# Iterations 6
# Poisson Regression Results
# ==============================================================================
# Dep. Variable: Counts No. Observations: 30
# Model: Poisson Df Residuals: 26
# Method: MLE Df Model: 3
# Date: Wed, 10 Jun 2020 Pseudo R-squ.: 0.3598
# Time: 15:03:38 Log-Likelihood: -39.804
# converged: True LL-Null: -62.175
# Covariance Type: nonrobust LLR p-value: 1.049e-09
# ======================================================================================
# coef std err z P>|z| [0.025 0.975]
# --------------------------------------------------------------------------------------
# Intercept -1.8102 0.524 -3.452 0.001 -2.838 -0.783
# C(cartridge)[T.兼容] 0.9585 0.359 2.668 0.008 0.254 1.663
# Pages 0.0167 0.005 3.267 0.001 0.007 0.027
# Length 0.0240 0.010 2.392 0.017 0.004 0.044
# ======================================================================================
# 在显著水平\alpha=0.05时，三个自变量影响显著性检验的P值比较小，均非常显著；
# 由系数的正负号可知，打印的页数和时长增加，发生故障的次数增加，使用兼容耗材比原装耗材发生的故障要多；
# 预测
res = printer_model.predict(pd.DataFrame({'cartridge':['原装'], 'Pages': [92], 'Length':[46.228]}))
print(res)
# 0 2.299743
浙公网安备 33010602011771号