python - scipy.curve_fit 与 numpy.polyfit 不同的协方差矩阵

Question

我正在使用 Python 3.6 进行数据拟合。最近，我遇到了以下问题，我缺乏经验，因此我不确定如何处理。

如果我在同一组数据点上使用 numpy.polyfit(x, y, 1, cov=True) 和 scipy.curve_fit(lambda: x, a, b: a*x+b, x, y) ，我得到几乎相同的系数 a 和 b。但是 scipy.curve_fit 的协方差矩阵的值大约是 numpy.polyfit 的值的一半。

由于我想使用协方差矩阵的对角线来估计系数的不确定性(u = numpy.sqrt(numpy.diag(cov)))，所以我有三个问题:

1. 哪个协方差矩阵是正确的(我应该使用哪个)？
2. 为什么会有差异？
3. 需要什么才能使它们相等？

谢谢!

编辑:

import numpy as np
import scipy.optimize as sc

data = np.array([[1,2,3,4,5,6,7],[1.1,1.9,3.2,4.3,4.8,6.0,7.3]]).T

x=data[:,0]
y=data[:,1]

A=np.polyfit(x,y,1, cov=True)
print('Polyfit:', np.diag(A[1]))

B=sc.curve_fit(lambda x,a,b: a*x+b, x, y)
print('Curve_Fit:', np.diag(B[1]))

如果我使用 statsmodels.api，结果对应于 curve_fit 的结果。

Answer 1

我想这与 this 有关

593          # Some literature ignores the extra -2.0 factor in the denominator, but 
594          #  it is included here because the covariance of Multivariate Student-T 
595          #  (which is implied by a Bayesian uncertainty analysis) includes it. 
596          #  Plus, it gives a slightly more conservative estimate of uncertainty. 
597          if len(x) <= order + 2: 
598              raise ValueError("the number of data points must exceed order + 2 " 
599                               "for Bayesian estimate the covariance matrix") 
600          fac = resids / (len(x) - order - 2.0) 
601          if y.ndim == 1: 
602              return c, Vbase * fac 
603          else: 
604              return c, Vbase[:,:, NX.newaxis] * fac 

在这种情况下，len(x) - order 是 4 而 (len(x) - order - 2.0) 是 2，这可以解释为什么你的值是相差 2 倍。

这解释了问题 2。问题 3 的答案可能是“获取更多数据”。对于更大的 len(x)，差异可能可以忽略不计。

哪种表述正确(问题 1)可能是 Cross Validated 的问题，但我假设它是 curve_fit 因为它明确用于计算您所说的不确定性。来自documentation

pcov : 2d array

The estimated covariance of popt. The diagonals provide the variance of the parameter estimate. To compute one standard deviation errors on the parameters use perr = np.sqrt(np.diag(pcov)).

虽然上面的 polyfit 代码中的注释说它的意图更多是用于 Student-T 分析。