Search code examples
pythonnumpycurve-fitting

numpy Polynomial.fit returns unrealistic values

I would like to fit a polynomial to a series of x/y-datapoints and evaluate it at arbitrary x-values. I am using the new numpy polynomial API:

import matplotlib.pyplot as plt
import numpy as np

x = [50.0, 150.0, 250.0, 400.0, 600.0, 850.0, 1000.0]
y = [3.2, 10.1, 16.3, 43.7, 69.1, 45.2, 10.8]

mypol = np.polynomial.polynomial.Polynomial.fit(
    x = x,
    y = y,
    deg = 3,
)
np.polyval(mypol.coef, x)

Unfortunately, this returns values which are far greater than the original y-values:

array([7.34024294e+06, 1.96122504e+08, 9.06027339e+08, 3.70658027e+09,
       1.25012342e+10, 3.55289563e+10, 5.78446450e+10])

I played around with the order of the polynomial, but wasn't able to get anything useful from the function. What am I doing wrong?

enter image description here

Solution

  • Numpy's polyval should be evaluated at the x-points, not the y-points. On top of that, the coefficients should be given in order of highest to lowest degree. That said, mypol.coef does not necessarily return the expanded polynomial coefficients, so mypol.coef may be useless for np.polyval.

    Instead, the numpy.polynomial.polynomial.Polynomial object returned by numpy.polynomial.Polynomial.fit can be evaluated at the x-points of interest, i.e. you can do mypol(x). Just note that x here should be a numpy array. Also note that the fit polynomial is a least-squares fit, meaning it is not an exact interpolation.

    import matplotlib.pyplot as plt
import numpy as np
from numpy.polynomial import Polynomial

x = np.array([50.0, 150.0, 250.0, 400.0, 600.0, 850.0, 1000.0])
y = np.array([3.2, 10.1, 16.3, 43.7, 69.1, 45.2, 10.8])

mypol = Polynomial.fit(x, y, 3)

x_fine = np.linspace(x.min(), x.max(), 100)
y_fine = mypol(x_fine)

fig, ax = plt.subplots()
ax.scatter(x, y, label="data")
ax.plot(x_fine, y_fine, "tab:orange", label="fit")
ax.set(xlabel="x", ylabel="y", title="NumPy Polynomial Fit")
ax.legend()