通过 scipy.stats.fit() 从 MLE 中检索优化结果?

Retrieve optimization results from MLE by scipy.stats.fit()?

我正在尝试通过 scipy.stats."some_dist".fit() 来估计不同分布的参数,但我在检索有关正在使用的优化过程的任何信息时遇到了极大的困难。具体来说,我正在寻找大多数实现算法输出的 Hessian,如记录 here.

对于它的所有优化,SciPy returns 一个名为 OptimizeResult 的对象,其中包含有用的信息,例如 Hessian,所以在 .fit() 之后如何获得它'ing 一些数据,调用相同的优化?

看起来需要一点 source-diving;幸运的是,它不需要复制大量源代码。这就是基本合身的工作原理:

from scipy.stats import cauchy
data = [0, 3, 4, 4, 5, 9]
res = cauchy.fit(data)   # (3.9798237305661255, 0.9205374643383732)

这就是它被修改为 return OptimizeResult:

的方式
from scipy.optimize import minimize
args = cauchy._fitstart(data)
x0, func, restore, args = cauchy._reduce_func(args, {})
res = minimize(func, x0, args=(data,), method='BFGS')

现在 res

      fun: 14.337039523098689
 hess_inv: array([[ 0.23321703, -0.0117229 ],
       [-0.0117229 ,  0.36807373]])
      jac: array([ 0.0000000e+00, -1.1920929e-07])
  message: 'Optimization terminated successfully.'
     nfev: 32
      nit: 5
     njev: 8
   status: 0
  success: True
        x: array([3.9798262 , 0.92055376])

您可能会认为这些参数是 res.x。被最小化的函数是 "penalized NNLF"(非负似然函数)。


顺便说一句,

For all its optimizations, SciPy returns an object called OptimizeResult

是一个 over-generalization。 minimize 方法也是如此。 scipy.stats.fit 默认使用 fmin,return 不是这样的。

因此,如果想要与 fmin 相同的输出,可以安排,但那里不会有额外的信息。

from scipy.optimize import fmin
args = cauchy._fitstart(data)
x0, func, restore, args = cauchy._reduce_func(args, {})
res = fmin(func, x0, args=(data,))
print(tuple(res))   #(3.9798237305661255, 0.9205374643383732)

使用 minimize 和 method='Nelder-Mead' 具有基本相同的效果。您确实获得了一些额外的信息,但鉴于该方法是 simplex-based(未计算导数),此信息的用途有限。

res = minimize(func, x0, args=(data,), method='Nelder-Mead')
print(res)
print(tuple(res.x))

打印

 final_simplex: (array([[3.97982373, 0.92053746],
       [3.97983057, 0.92060317],
       [3.97977536, 0.92059568]]), array([14.33703952, 14.33703953, 14.33703953]))
           fun: 14.337039523477827
       message: 'Optimization terminated successfully.'
          nfev: 58
           nit: 31
        status: 0
       success: True
             x: array([3.97982373, 0.92053746])
(3.9798237305661255, 0.9205374643383732)