statsmodel 中的 MNLogit 返回 nan

Question

我正在尝试在著名的 iris 数据集上使用 statsmodels 的 MNLogit 函数。当我尝试拟合模型时，我得到：“当前函数值：nan”。这是我正在使用的代码：

import statsmodels.api as st
iris = st.datasets.get_rdataset('iris','datasets')
y = iris.data.Species
x = iris.data.ix[:, 0:4]
x = st.add_constant(x, prepend = False)
mdl = st.MNLogit(y, x)
mdl_fit = mdl.fit()
print (mdl_fit.summary())

Answer 1

在鸢尾花的例子中，我们可以完美地预测 Setosa。这会导致 Logit 和 MNLogit 中（部分）完美分离的问题。

完美的分离有利于预测，但 logit 的参数会趋于无穷大。在这种情况下，我得到一个奇异矩阵错误，而不是使用相对较新版本的 statsmodels master（在 Windows 上）的 Nans 错误。

离散模型的默认优化器是 Newton，当 Hessian 矩阵变得奇异时，它会失败。其他不使用 Hessian 信息的优化器也能够完成优化。例如使用“bfgs”，我得到

>>> mdl_fit = mdl.fit(method='bfgs')
Warning: Maximum number of iterations has been exceeded.
         Current function value: 0.057112
         Iterations: 35
         Function evaluations: 37
         Gradient evaluations: 37
e:\josef\eclipsegworkspace\statsmodels-git\statsmodels-all-new2_py27\statsmodels\statsmodels\base\model.py:471: ConvergenceWarning: Maximum Likelihood optimization failed to converge. Check mle_retvals
  "Check mle_retvals", ConvergenceWarning)

Setosa 的预测概率本质上是 (1, 0, 0)，即它们被完美预测

>>> fitted = mdl_fit.predict()
>>> fitted[y=='setosa'].min(0)
array([  9.99497636e-01,   2.07389867e-11,   1.71740822e-38])
>>> fitted[y=='setosa'].max(0)
array([  1.00000000e+00,   5.02363854e-04,   1.05778255e-20])

然而，由于完全分离，参数无法识别，这些值主要由优化器的停止标准确定，并且标准误差非常大。

>>> print(mdl_fit.summary())
                          MNLogit Regression Results                          
==============================================================================
Dep. Variable:                Species   No. Observations:                  150
Model:                        MNLogit   Df Residuals:                      140
Method:                           MLE   Df Model:                            8
Date:                Mon, 20 Jul 2015   Pseudo R-squ.:                  0.9480
Time:                        04:08:04   Log-Likelihood:                -8.5668
converged:                      False   LL-Null:                       -164.79
                                        LLR p-value:                 9.200e-63
=====================================================================================
Species=versicolor       coef    std err          z      P>|z|      [95.0% Conf. Int.]
--------------------------------------------------------------------------------------
Sepal.Length          -1.4959    444.817     -0.003      0.997      -873.321   870.330
Sepal.Width           -8.0560    282.766     -0.028      0.977      -562.267   546.155
Petal.Length          11.9301    374.116      0.032      0.975      -721.323   745.184
Petal.Width            1.7039    759.366      0.002      0.998     -1486.627  1490.035
const                  1.6444   1550.515      0.001      0.999     -3037.309  3040.597
--------------------------------------------------------------------------------------
Species=virginica       coef    std err          z      P>|z|      [95.0% Conf. Int.]
-------------------------------------------------------------------------------------
Sepal.Length         -8.0348    444.835     -0.018      0.986      -879.896   863.827
Sepal.Width         -15.8195    282.793     -0.056      0.955      -570.083   538.444
Petal.Length         22.1797    374.155      0.059      0.953      -711.152   755.511
Petal.Width          14.0603    759.384      0.019      0.985     -1474.304  1502.425
const                -6.5053   1550.533     -0.004      0.997     -3045.494  3032.483
=====================================================================================

关于statsmodels中的实现

Logit 专门检查完美分离并引发异常，可以选择将其弱化为警告。对于像 MNLogit 这样的其他模型，还没有对完美分离的明确检查，主要是因为缺乏良好的测试用例和易于识别的一般条件。（像https://github.com/statsmodels/statsmodels/issues/516这样的几个问题仍然悬而未决）

我的总体策略：

当收敛失败时，尝试不同的优化器和不同的起始值（start_params）。如果某些优化器成功，那么这可能是一个困难的优化问题，要么是目标函数的曲率，要么是尺度不当的解释变量或类似问题。一个有用的检查是使用稳健优化器的参数估计（例如nm或 ）作为更严格的优化器（例如或 ）powell的起始值。newtonbfgs

如果某些优化器收敛后结果仍然不好，那么可能是数据本身的问题，例如 Logit、Probit 和其他几种模型中的完美分离或奇异或接近奇异的设计矩阵。在这种情况下，必须更改模型。通过互联网搜索可以找到完美分离的建议。