python - 在 Scipy 中，curve_fit 如何以及为什么计算参数估计的协方差

Question

我一直在使用 scipy.optimize.leastsq 来拟合一些数据。我想获得这些估计值的一些置信区间，因此我查看了 cov_x 输出，但文档非常不清楚这是什么以及如何从中获取我的参数的协方差矩阵。

首先它说它是雅可比行列式，但在 notes它还说“cov_x 是 Hessian 的 Jacobian 近似”，因此它实际上不是 Jacobian，而是使用 Jacobian 的某种近似的 Hessian。这些说法中哪一个是正确的？

其次，这句话让我很困惑:

This matrix must be multiplied by the residual variance to get the covariance of the parameter estimates – see curve_fit.

我确实去看看他们所做的 curve_fit 的源代码:

s_sq = (func(popt, *args)**2).sum()/(len(ydata)-len(p0))
pcov = pcov * s_sq

对应于将 cov_x 乘以 s_sq 但我在任何引用资料中都找不到这个等式。有人可以解释为什么这个等式是正确的吗？我的直觉告诉我它应该是相反的，因为 cov_x 应该是一个衍生物(雅可比或黑森)所以我在想:cov_x * covariance(parameters) = sum of errors(residuals) 其中 sigma(parameters) 是我想要的。

我如何将curve_fit 正在做的事情与我在例如看到的事情联系起来。维基百科: http://en.wikipedia.org/wiki/Propagation_of_uncertainty#Non-linear_combinations

Answer 1

好的，我想我找到了答案。首先解决方案:cov_x*s_sq 只是您想要的参数的协方差。取对角元素的 sqrt 会给你标准偏差(但要小心协方差!)。

残差 = 减少卡方 = s_sq = sum[(f(x)-y)^2]/(N-n)，其中 N 是数据点的数量，n 是拟合参数的数量。 Reduced chi square .

我感到困惑的原因是，由 minimumsq 给出的 cov_x 实际上并不是其他地方所谓的 cov(x)，而是简化的 cov(x) 或分数 cov(x)。它没有出现在任何其他引用资料中的原因是它是一种简单的重新缩放，在数值计算中很有用，但与教科书无关。

关于 Hessian 与 Jacobian，文档措辞不佳。在这两种情况下计算的都是 Hessian 矩阵，这是显而易见的，因为 Jacobian 矩阵至少为零。他们的意思是他们正在使用雅可比矩阵的近似值来找到 Hessian 矩阵。

进一步说明。似乎curve_fit结果实际上并没有考虑到错误的绝对大小，而只是考虑了所提供的sigma的相对大小。这意味着即使错误栏改变了一百万倍，返回的 pcov 也不会改变。这当然是不对的，但似乎是标准做法，即。 Matlab 在使用他们的曲线拟合工具箱时做同样的事情。此处描述了正确的程序:https://en.wikipedia.org/wiki/Linear_least_squares_(mathematics)#Parameter_errors_and_correlation

一旦找到最优值，至少对于线性最小二乘而言，这样做似乎相当简单。