numpy.ma.polyfit#

ma.polyfit(x, y, deg, rcond=None, full=False, w=None, cov=False)[源代码]#

最小二乘多项式拟合.

备注

这是旧的多项式API的一部分.自版本1.4起,新的多项式API定义在 numpy.polynomial 中是首选.差异的总结可以在过渡指南中找到.

拟合一个多项式 p(x) = p[0] * x**deg + ... + p[deg] 的次数 deg 到点 (x, y).返回一个系数向量 p,该向量在顺序 deg, deg-1, … 0 中使平方误差最小化.

Polynomial.fit 类方法是新代码推荐使用的,因为其在数值上更稳定.有关该方法的更多信息,请参阅其文档.

参数:

xarray_like, 形状 (M,): M 个样本点 (x[i], y[i]) 的 x 坐标.
yarray_like, 形状 (M,) 或 (M, K): 样本点的 y 坐标.通过传入一个包含每列一个数据集的二维数组,可以同时拟合多个共享相同 x 坐标的样本点数据集.
degint: 拟合多项式的度数
rcond浮点型, 可选: 拟合的相对条件数.相对于最大奇异值,小于此值的奇异值将被忽略.默认值是 len(x)*eps,其中 eps 是浮点类型的相对精度,大多数情况下约为 2e-16.
fullbool, 可选: 确定返回值性质的开关.当它为 False（默认）时,仅返回系数,当为 True 时,还会返回来自奇异值分解的诊断信息.
warray_like, 形状 (M,), 可选: 权重.如果不是 None,权重 w[i] 应用于未平方的残差 y[i] - y_hat[i] 在 x[i] 处.理想情况下,权重的选择应使得产品 w[i]*y[i] 的误差都具有相同的方差.使用逆方差加权时,使用 w[i] = 1/sigma(y[i]).默认值是 None.
cov布尔值或字符串,可选: 如果给定且不为 False,不仅返回估计值,还返回其协方差矩阵.默认情况下,协方差按 chi2/dof 缩放,其中 dof = M - (deg + 1),即权重被认为在相对意义上是不可靠的,一切都缩放使得减少的 chi2 为单位.如果 cov='unscaled',则不进行此缩放,这在权重为 w = 1/sigma 且 sigma 已知为不确定性的可靠估计的情况下是相关的.

返回:

pndarray, 形状 (deg + 1,) 或 (deg + 1, K)

多项式系数,最高次幂在前.如果 y 是二维的,第 k 个数据集的系数在 p[:,k] 中.

残差, 秩, 奇异值, rcond

只有在 full == True 时,这些值才会被返回.

残差 – 最小二乘拟合的残差平方和
rank – 缩放范德蒙矩阵的有效秩
系数矩阵
singular_values – 缩放范德蒙矩阵的奇异值
系数矩阵
rcond – rcond 的值.

更多细节,请参见 numpy.linalg.lstsq.

Vndarray, 形状 (deg + 1, deg + 1) 或 (deg + 1, deg + 1, K)

仅在 full == False 且 cov == True 时呈现.这是多项式系数估计的协方差矩阵.该矩阵的对角线是每个系数的方差估计.如果 y 是一个 2-D 数组,那么第 k 个数据集的协方差矩阵在 V[:,:,k] 中.

警告:

RankWarning

在最小二乘拟合中,系数矩阵的秩不足.仅当 full == False 时,才会引发警告.

可以通过以下方式关闭警告:

>>> import warnings
>>> warnings.simplefilter('ignore', np.exceptions.RankWarning)

参见

polyval: 计算多项式值.
linalg.lstsq: 计算最小二乘拟合.
scipy.interpolate.UnivariateSpline: 计算样条拟合.

备注

在 x 中的任何掩码值都会在 y 中传播,反之亦然.

该解决方案最小化了平方误差

\[E = \sum_{j=0}^k |p(x_j) - y_j|^2\]

在方程中:

x[0]**n * p[0] + ... + x[0] * p[n-1] + p[n] = y[0]
x[1]**n * p[0] + ... + x[1] * p[n-1] + p[n] = y[1]
...
x[k]**n * p[0] + ... + x[k] * p[n-1] + p[n] = y[k]

系数 p 的系数矩阵是一个范德蒙矩阵.

当最小二乘拟合条件较差时,`polyfit` 会发出 RankWarning.这意味着由于数值误差,最佳拟合定义不明确.通过降低多项式次数或将 x 替换为 x - x.mean(),可以改善结果.`rcond` 参数也可以设置为小于其默认值,但由此产生的拟合可能是虚假的:包括来自小奇异值的贡献可能会增加数值噪声.

注意,当多项式的次数较高或样本点的区间中心不佳时,拟合多项式系数本质上是不稳定的.在这些情况下,应始终检查拟合质量.当多项式拟合不令人满意时,样条可能是很好的替代方法.

参考文献

[1]

Wikipedia, “曲线拟合”, https://en.wikipedia.org/wiki/Curve_fitting

[2]

Wikipedia, “多项式插值”, https://en.wikipedia.org/wiki/Polynomial_interpolation

示例

>>> import numpy as np
>>> import warnings
>>> x = np.array([0.0, 1.0, 2.0, 3.0,  4.0,  5.0])
>>> y = np.array([0.0, 0.8, 0.9, 0.1, -0.8, -1.0])
>>> z = np.polyfit(x, y, 3)
>>> z
array([ 0.08703704, -0.81349206,  1.69312169, -0.03968254]) # may vary

使用 poly1d 对象处理多项式非常方便:

>>> p = np.poly1d(z)
>>> p(0.5)
0.6143849206349179 # may vary
>>> p(3.5)
-0.34732142857143039 # may vary
>>> p(10)
22.579365079365115 # may vary

高阶多项式可能会剧烈振荡:

>>> with warnings.catch_warnings():
...     warnings.simplefilter('ignore', np.exceptions.RankWarning)
...     p30 = np.poly1d(np.polyfit(x, y, 30))
...
>>> p30(4)
-0.80000000000000204 # may vary
>>> p30(5)
-0.99999999999999445 # may vary
>>> p30(4.5)
-0.10547061179440398 # may vary

插图:

>>> import matplotlib.pyplot as plt
>>> xp = np.linspace(-2, 6, 100)
>>> _ = plt.plot(x, y, '.', xp, p(xp), '-', xp, p30(xp), '--')
>>> plt.ylim(-2,2)
(-2, 2)
>>> plt.show()