statsmodels 如何编码以字符串形式输入的 endog 变量？

Question

我是使用 statsmodels 进行统计分析的新手。我大多数时候都会得到预期的答案，但有些事情我不太明白 statsmodels 在以字符串形式输入时定义逻辑回归的 endog（依赖）变量的方式。

可以定义如下所示的示例 Pandas 数据框来说明该问题。yN、yA 和 yA2 列代表定义 endog 变量的不同方式： yN 是编码为 0、1 的二进制变量；yA 是一个二进制变量，编码为“y”、“n”；yA2 是编码为“x”、“y”和“w”的变量：

import pandas as pd

df = pd.DataFrame({'yN':[0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1],
                   'yA':['y','y','y','y','y','y','y','n','n','n','n','n','n','n','n','n','n','n','n','n',],
                   'yA2':['y','y','y','w','y','w','y','n','n','n','n','n','n','n','n','n','n','n','n','n',],
                   'xA':['a','a','b','b','b','c','c','c','c','c','a','a','a','a','b','b','b','b','c','c']})

数据框看起来像：

   xA yA yA2  yN
0   a  y   y   0
1   a  y   y   0
2   b  y   y   0
3   b  y   w   0
4   b  y   y   0
5   c  y   w   0
6   c  y   y   0
7   c  n   n   1
8   c  n   n   1
9   c  n   n   1
10  a  n   n   1
11  a  n   n   1
12  a  n   n   1
13  a  n   n   1
14  b  n   n   1
15  b  n   n   1
16  b  n   n   1
17  b  n   n   1
18  c  n   n   1
19  c  n   n   1

我可以使用 0/1 编码的 endog 变量和分类 exog 变量 (xA) 运行“标准”逻辑回归，如下所示：

import statsmodels.formula.api as smf
import statsmodels.api as sm

phjLogisticRegressionResults = smf.glm(formula='yN ~ C(xA)',
                                       data=df,
                                       family = sm.families.Binomial(link = sm.genmod.families.links.logit)).fit()

print('\nResults of logistic regression model')
print(phjLogisticRegressionResults.summary())

这会产生以下结果，这与我的预期完全一致：

                 Generalized Linear Model Regression Results                  
==============================================================================
Dep. Variable:                     yN   No. Observations:                   20
Model:                            GLM   Df Residuals:                       17
Model Family:                Binomial   Df Model:                            2
Link Function:                  logit   Scale:                             1.0
Method:                          IRLS   Log-Likelihood:                -12.787
Date:                Thu, 18 Jan 2018   Deviance:                       25.575
Time:                        02:19:45   Pearson chi2:                     20.0
No. Iterations:                     4                                         
==============================================================================
                 coef    std err          z      P>|z|      [0.025      0.975]
------------------------------------------------------------------------------
Intercept      0.6931      0.866      0.800      0.423      -1.004       2.391
C(xA)[T.b]    -0.4055      1.155     -0.351      0.725      -2.669       1.858
C(xA)[T.c]     0.2231      1.204      0.185      0.853      -2.137       2.583
==============================================================================

但是，如果我使用二进制 endog 变量编码“y”和“n”运行相同的模型（但与前面示例中直观的 0/1 编码完全相反），或者如果我包含一个变量，其中某些“y”编码已被替换为“w”（即现在有 3 个结果），它仍然产生相同的结果，如下所示：

phjLogisticRegressionResults = smf.glm(formula='yA ~ C(xA)',
                                       data=df,
                                       family = sm.families.Binomial(link = sm.genmod.families.links.logit)).fit()


                 Generalized Linear Model Regression Results                  
==============================================================================
Dep. Variable:     ['yA[n]', 'yA[y]']   No. Observations:                   20
Model:                            GLM   Df Residuals:                       17
Model Family:                Binomial   Df Model:                            2
Link Function:                  logit   Scale:                             1.0
Method:                          IRLS   Log-Likelihood:                -12.787
Date:                Thu, 18 Jan 2018   Deviance:                       25.575
Time:                        02:29:06   Pearson chi2:                     20.0
No. Iterations:                     4                                         
==============================================================================
                 coef    std err          z      P>|z|      [0.025      0.975]
------------------------------------------------------------------------------
Intercept      0.6931      0.866      0.800      0.423      -1.004       2.391
C(xA)[T.b]    -0.4055      1.155     -0.351      0.725      -2.669       1.858
C(xA)[T.c]     0.2231      1.204      0.185      0.853      -2.137       2.583
==============================================================================

...和...

phjLogisticRegressionResults = smf.glm(formula='yA2 ~ C(xA)',
                                       data=df,
                                       family = sm.families.Binomial(link = sm.genmod.families.links.logit)).fit()


                       Generalized Linear Model Regression Results                        
==========================================================================================
Dep. Variable:     ['yA2[n]', 'yA2[w]', 'yA2[y]']   No. Observations:                   20
Model:                                        GLM   Df Residuals:                       17
Model Family:                            Binomial   Df Model:                            2
Link Function:                              logit   Scale:                             1.0
Method:                                      IRLS   Log-Likelihood:                -12.787
Date:                            Thu, 18 Jan 2018   Deviance:                       25.575
Time:                                    02:29:06   Pearson chi2:                     20.0
No. Iterations:                                 4                                         
==============================================================================
                 coef    std err          z      P>|z|      [0.025      0.975]
------------------------------------------------------------------------------
Intercept      0.6931      0.866      0.800      0.423      -1.004       2.391
C(xA)[T.b]    -0.4055      1.155     -0.351      0.725      -2.669       1.858
C(xA)[T.c]     0.2231      1.204      0.185      0.853      -2.137       2.583
==============================================================================

该部。输出表中的变量单元格可识别，但存在差异，但结果相同。当将 endog 变量作为字符串变量输入时，statsmodels 使用什么规则对它进行编码？

Answer 1

警告：此行为不是有意设计的，而是通过 patsy 和 statsmodels 的交互产生的。

首先，patsy 对字符串公式和数据进行所有转换以创建相应的设计矩阵，并可能对响应变量进行转换。如果响应变量endogy 是字符串，则 patsy 会将其视为分类变量，并使用分类变量的默认编码并构建相应的虚拟变量数组。另外，据我所知，patsy 按字母数字对级别进行排序，这决定了列的顺序。

模型的主要部分（GLM 或 Logit/Probit）仅采用 patsy 提供的数组，并在可能的情况下以模型适当的方式对其进行解释，并且无需进行太多具体的输入检查。

在示例中，patsy 创建具有两列或三列的虚拟变量数组。statsmodels 将其解释为“成功”、“失败”计数。因此，按字母数字顺序排列的最低类别定义为“成功”。行总和对应于观察中的试验次数，在本例中为 1。如果或它适用于三列，则必须缺乏输入检查，这意味着它是先对其余的二进制响应。（我猜这是实现细节的结果，而不是设计使然。）

一个相关的问题是离散模型 Logit。https://github.com/statsmodels/statsmodels/issues/2733 目前没有明确的解决方案不需要大量猜测用户的意图。

因此，目前最好对二进制模型使用数值，特别是因为定义“成功”的内容和参数的符号取决于分类级别名称的字母数字顺序。例如，尝试将“成功”级别重命名为“z”而不是“n”。