numpy polynomial.Polynomial.fit() gives different coefficients than polynomial.polyfit()

Question

I do not understand why polynomial.Polynomial.fit() gives coefficients very different from the expected coefficients :

import numpy as np

x = np.linspace(0, 10, 50)
y = x**2 + 5 * x + 10

print(np.polyfit(x, y, 2))
print(np.polynomial.polynomial.polyfit(x, y, 2))
print(np.polynomial.polynomial.Polynomial.fit(x, y, 2))

Gives :

[ 1.  5. 10.]
[10.  5.  1.]
poly([60. 75. 25.])

The two first results are OK, and thanks to this answer I understand why the two arrays are in reversed order.

However, I do not understand the signification of the third result. The coefficients looks wrong, though the polynomial that I got this way seems to give correct predicted values.

Limpid

Answer 1

The answer is slightly hidden in the docs, of course. Looking at the class numpy.polynomial.polynomial.Polynomial(coef, domain=None, window=None) It is clear that in general the coefficients [a, b, c, ...] are for the polynomial a + b * x + c * x**2 + .... However, there are the keyword parameters domain and window both with default [-1,1]. I am not into that class, so I am not sure about the purpose, but it is clear that a remapping takes place. Now in the case of polynomial.Polynomial.fit() one has a class method that automatically takes the x data as domain, but still makes the mapping to the window. Hence, in the OP [0-10] is mapped onto [-1,1]. This is done by x = x' / 5 - 1 or x' -> 5 * x + 5. Putting the latter in the OP polynomial we get

( 5 x' + 5 )**2 + 5 * ( 5 * x' + 5 ) + 10 = 25 * x'**2 + 75 * x' + 60

Voila.

To get the expected result one has to put

print(np.polynomial.polynomial.Polynomial.fit(x, y, 2, window=[0, 10] ) )

wich gives

poly([10.  5.  1.])

Unsatisfactory

I had finished by getting the expected result by using domain=[] but I still didn't understand the default result. Your explanation is very clear. – Limpid

In the documentation, it is said : "The default value was the class domain in NumPy 1.4 and None in later versions". I have no clue why they changed it so that default value gives now correct but unexpected results. – Limpid

Answer 2

Buried in the docs:

Note that the coefficients are given in the scaled domain defined by the linear mapping between the window and domain. convert can be used to get the coefficients in the unscaled data domain.

So use:

poly.convert()    

This will rescale your coefficients to what you are probably expecting.

Example for data generated from 1 + 2x + 3x^2:

from numpy.polynomial import Polynomial

test_poly = Polynomial.fit([0, 1, 2, 3, 4, 5],
                           [1, 6, 17, 34, 57, 86],
                           2)

print(test_poly)
print(test_poly.convert())

Output:

poly([24.75 42.5  18.75])
poly([1. 2. 3.])

Domicile