python - 需要返回协方差的 Python 多项式拟合函数

Question

我想对数据集 (X,Y,Yerr) 进行最小二乘多项式拟合，并获得拟合参数的协方差矩阵。此外，由于我有很多数据集，CPU 时间是一个问题，所以我正在寻求一种分析(=快速)解决方案。我发现了以下(非理想)选项:

numpy.polyfit 拟合，但不考虑误差 Yerr，也不返回协方差；

numpy.polynomial.polynomial.polyfit 确实接受 Yerr 作为输入(以权重的形式)，但也不返回协方差；

scipy.optimize.curve_fit 和 scipy.optimize.leastsq 可以定制以拟合多项式并返回协方差矩阵，但是 - 作为迭代方法 - 这些方法要慢得多比 polyfit 例程(产生解析解)；

Python 是否提供返回拟合参数协方差的分析多项式拟合例程(或者我必须自己编写一个 :-)？

更新:看来在 Numpy 1.7.0 中，numpy.polyfit 现在不仅 接受权重，而且还返回系数的协方差矩阵……所以，问题解决了! :-)

Answer 1

您想要一个无需额外开销即可返回协方差矩阵的快速加权最小二乘模型吗？一般来说，正确的协方差矩阵将取决于数据生成过程 (DGP)，因为不同的 DGP(比如误差的异方差性)意味着参数估计的不同分布(想想怀特与 OLS 标准误差)。但是，如果您可以假设 WLS 是正确的方法，并且我相信您会使用 WLS 的 beta 的渐近方差估计，(1/n X'V^-1X)^-1，其中 V 是加权矩阵从 Yerrs 创建。如果 numpy.polynomial.polynomial.polyfit 适合您，那将是一个非常简单的公式。

我寻找在线引用，但找不到。但参见 Fumio Hayashi 的 Ecomometrics，2000 年，普林斯顿大学出版社，p. 133 - 137 进行推导和讨论。

2012 年 12 月 4 日更新:还有另一个堆栈溢出问题很接近: numpy.polyfit has no keyword 'cov'这对如何使用 scikits.statsmodels 做你想做的事情有一个很好的解释(带代码)。我相信您会想要替换该行:

result = sm.OLS(Y,reg_x_data).fit()

到

result = sm.WLS(Y,reg_x_data, weights).fit()

您可以像以前一样使用 numpy.polynomial.polynomial.polyfit 将权重定义为 Yerr 的函数。有关将 statsmodels 与 WLS 结合使用的更多详细信息，请访问statsmodels website .