numpy.polynomial.polynomial.polyfit#

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

对数据进行多项式的最小二乘拟合。

返回一个 deg 次多项式的系数,该多项式是给定点 x 处数据值 y 的最小二乘拟合。如果 y 是一维的,返回的系数也将是一维的。如果 y 是二维的,将进行多次拟合,对 y 的每一列进行一次拟合,结果系数存储在二维返回值的相应列中。拟合的多项式形式为

\[p(x) = c_0 + c_1 * x + ... + c_n * x^n,\]

其中 ndeg

参数:
xarray_like, 形状 (M,)

M 个样本(数据)点 (x[i], y[i]) 的 x 坐标。

yarray_like, 形状 (M,) 或 (M, K)

样本点的 y 坐标。通过为 y 传入一个每列包含一个数据集的二维数组,可以在一次 polyfit 调用中(独立地)拟合多个共享相同 x 坐标的样本点集。

deg整数或一维 array_like

拟合多项式的次数。如果 deg 是一个整数,则拟合中包含直至(包括)deg 次的所有项。对于 NumPy 版本 >= 1.11.0,可以使用一个整数列表来指定要包含的项的次数。

rcond浮点数, 可选

拟合的相对条件数。小于 rcond 的奇异值,相对于最大的奇异值,将被忽略。默认值是 len(x)*eps,其中 eps 是平台浮点类型的相对精度,在大多数情况下约为 2e-16。

full布尔值, 可选

确定返回值性质的开关。当为 False(默认值)时,仅返回系数;当为 True 时,还会返回奇异值分解(用于求解拟合的矩阵方程)的诊断信息。

warray_like, 形状 (M,), 可选

权重。如果不是 None,权重 w[i] 适用于 x[i] 处的未平方残差 y[i] - y_hat[i]。理想情况下,权重的选择应使乘积 w[i]*y[i] 的误差都具有相同的方差。使用逆方差加权时,使用 w[i] = 1/sigma(y[i])。默认值为 None。

返回:
coefndarray, 形状 (deg + 1,) 或 (deg + 1, K)

多项式系数,按从低到高的顺序排列。如果 y 是二维的,则 coef 中第 k 列的系数表示对 y 的第 k 列数据进行的拟合多项式。

[residuals, rank, singular_values, rcond]列表

这些值仅在 full == True 时返回。

  • residuals – 最小二乘拟合的平方残差之和

  • rank – 缩放的范德蒙德矩阵的数值秩

  • singular_values – 缩放的范德蒙德矩阵的奇异值

  • rcond – rcond 的值。

更多详情,请参见 numpy.linalg.lstsq

抛出:
RankWarning

如果最小二乘拟合中的矩阵是秩亏的,则会抛出此警告。此警告仅在 full == False 时抛出。可以通过以下方式关闭警告:

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

注意

解是使加权平方误差之和最小的多项式 p 的系数

\[E = \sum_j w_j^2 * |y_j - p(x_j)|^2,\]

其中 \(w_j\) 是权重。通过建立(通常是)超定矩阵方程来解决此问题

\[V(x) * c = w * y,\]

其中 Vx 的加权伪范德蒙德矩阵,c 是待求解的系数,w 是权重,y 是观测值。然后使用 V 的奇异值分解来求解此方程。

如果 V 的某些奇异值太小而被忽略(且 full == False),则会抛出 RankWarning。这意味着系数的值可能确定得不好。拟合较低次的多项式通常可以消除此警告(当然,这可能不是您想要的;如果您有独立的原因选择不起作用的次数,您可能需要:a) 重新考虑这些原因,和/或 b) 重新考虑数据的质量)。rcond 参数也可以设置为小于其默认值,但由此产生的拟合可能是伪造的,并且包含大量的舍入误差。

使用双精度进行的多项式拟合在大约(多项式)20 次时会“失败”。使用切比雪夫或勒让德级数进行的拟合通常条件更好,但很大程度上仍取决于样本点的分布和数据的平滑度。如果拟合质量不足,样条函数可能是一个很好的替代方案。

示例

>>> import numpy as np
>>> from numpy.polynomial import polynomial as P
>>> x = np.linspace(-1,1,51)  # x "data": [-1, -0.96, ..., 0.96, 1]
>>> rng = np.random.default_rng()
>>> err = rng.normal(size=len(x))
>>> y = x**3 - x + err  # x^3 - x + Gaussian noise
>>> c, stats = P.polyfit(x,y,3,full=True)
>>> c # c[0], c[1] approx. -1, c[2] should be approx. 0, c[3] approx. 1
array([ 0.23111996, -1.02785049, -0.2241444 ,  1.08405657]) # may vary
>>> stats # note the large SSR, explaining the rather poor results
[array([48.312088]),                                        # may vary
 4,
 array([1.38446749, 1.32119158, 0.50443316, 0.28853036]),
 1.1324274851176597e-14]

没有额外噪声的相同内容

>>> y = x**3 - x
>>> c, stats = P.polyfit(x,y,3,full=True)
>>> c # c[0], c[1] ~= -1, c[2] should be "very close to 0", c[3] ~= 1
array([-6.73496154e-17, -1.00000000e+00,  0.00000000e+00,  1.00000000e+00])
>>> stats # note the minuscule SSR
[array([8.79579319e-31]),
 np.int32(4),
 array([1.38446749, 1.32119158, 0.50443316, 0.28853036]),
 1.1324274851176597e-14]