statsmodels 如何编码以字符串形式输入的 endog 变量?
How does statsmodels encode endog variables entered as strings?
我不熟悉使用 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
我可以 运行 'standard' 使用 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
==============================================================================
部门。输出中的可变单元格 table 识别但存在差异但结果相同。当 endog 变量作为字符串变量输入时,statsmodels 使用什么规则对其进行编码?
警告:此行为不是设计使然,而是通过 patsy 和 statsmodels 的交互产生的。
首先,patsy 对字符串公式和数据进行所有转换,以创建相应的设计矩阵,并可能对响应变量进行转换。
在响应变量 endog
或 y 是字符串的情况下,patsy 将其视为分类变量并使用分类变量的默认编码并构建相应的虚拟变量数组。此外,AFAIK patsy 按字母数字顺序对级别进行排序,这决定了列的顺序。
模型的主要部分,无论是 GLM 还是 Logit/Probit,只采用 patsy 提供的数组,并在可能的情况下以模型适当的方式对其进行解释,而无需太多特定的输入检查。
在示例中,patsy 创建了具有两列或三列的虚拟变量数组。 statsmodels 将其解释为 "success"、"failure" 计数。因此,按字母数字顺序排列的最低类别定义 "success"。行总和对应于观察中的试验次数,在本例中为 1。
如果它适用于三列,则必须缺少输入检查,这意味着它是第一个对其余的二进制响应。 (我猜这是实现细节的结果,而不是设计使然。)
离散模型 Logit 中有一个相关问题。 https://github.com/statsmodels/statsmodels/issues/2733
目前没有不需要大量二次猜测用户意图的明确解决方案。
所以现在最好对二元模型使用数值,尤其是因为定义 "success" 的内容和参数的符号取决于分类级别名称的字母数字顺序。
例如,尝试将 "success" 级别重命名为 "z" 而不是 "n"。
我不熟悉使用 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
我可以 运行 'standard' 使用 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
==============================================================================
部门。输出中的可变单元格 table 识别但存在差异但结果相同。当 endog 变量作为字符串变量输入时,statsmodels 使用什么规则对其进行编码?
警告:此行为不是设计使然,而是通过 patsy 和 statsmodels 的交互产生的。
首先,patsy 对字符串公式和数据进行所有转换,以创建相应的设计矩阵,并可能对响应变量进行转换。
在响应变量 endog
或 y 是字符串的情况下,patsy 将其视为分类变量并使用分类变量的默认编码并构建相应的虚拟变量数组。此外,AFAIK patsy 按字母数字顺序对级别进行排序,这决定了列的顺序。
模型的主要部分,无论是 GLM 还是 Logit/Probit,只采用 patsy 提供的数组,并在可能的情况下以模型适当的方式对其进行解释,而无需太多特定的输入检查。
在示例中,patsy 创建了具有两列或三列的虚拟变量数组。 statsmodels 将其解释为 "success"、"failure" 计数。因此,按字母数字顺序排列的最低类别定义 "success"。行总和对应于观察中的试验次数,在本例中为 1。 如果它适用于三列,则必须缺少输入检查,这意味着它是第一个对其余的二进制响应。 (我猜这是实现细节的结果,而不是设计使然。)
离散模型 Logit 中有一个相关问题。 https://github.com/statsmodels/statsmodels/issues/2733 目前没有不需要大量二次猜测用户意图的明确解决方案。
所以现在最好对二元模型使用数值,尤其是因为定义 "success" 的内容和参数的符号取决于分类级别名称的字母数字顺序。 例如,尝试将 "success" 级别重命名为 "z" 而不是 "n"。