2

I am trying to use a polynomial expression that would fit my function (signal). I am using numpy.polynomial.polynomial.Polynomial.fit function to fit my function(signal) using the coefficients. Now, after generating the coefficients, I want to put those coefficients back into the polynomial equation - get the corresponding y-values - and plot them on the graph. But I am not getting what I want (orange line) . What am I doing wrong here? enter image description here

Thanks.

import math
def getYValueFromCoeff(f,coeff_list): # low to high order
    y_plot_values=[]

    for j in range(len(f)):
        item_list= []
        for i in range(len(coeff_list)):
            item= (coeff_list[i])*((f[j])**i)
            item_list.append(item)
        y_plot_values.append(sum(item_list))
    print(len(y_plot_values))
    return y_plot_values

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

no_of_coef= 10

#original signal
x = np.linspace(0, 0.01, 10)
period = 0.01
y = np.sin(np.pi * x / period)

#poly fit
test1= poly.fit(x,y,no_of_coef)
coeffs= test1.coef
#print(test1.coef)
coef_y= getYValueFromCoeff(x, test1.coef)
#print(coef_y)

plt.plot(x,y)
plt.plot(x, coef_y)
MRR
  • 45
  • 1
  • 6

3 Answers3

0

If you check out the documentation, consider the two properties: poly.domain and poly.window. To avoid numerical issues, the range poly.domain = [x.min(), x.max()] of independent variable (x) that we pass to the fit() is being normalized to poly.window = [-1, 1]. This means the coefficients you get from poly.coef apply to this normalized range. But you can adjust this behaviour (sacrificing numerical stability) accordingly, that is, adjustig the poly.window will make your curves match:

...
test1 = poly.fit(x, y, deg=no_of_coef, window=[x.min(), x.max()])
...

But unless you have a good reason to do that, I'd stick to the default behaviour of fit().

As a side note: Evaluating polynomials or lists of coefficients is already implemented in numpy, e.g. using directly

coef_y = test1(x)

or alternatively using np.polyval.

flawr
  • 10,814
  • 3
  • 41
  • 71
0

I always like to see original solutions to problems. I urge you to continue to pursue that as that is the best way to learn how to fit functions programmatically. I also wanted to provide the solution that is much more tailored towards a standard numpy implementation. As for your custom function, you did really well. The only issue is that the coefficients are from high to low order, while you were counting up in powers from 0 to highest power. Simply counting down from highest power to 0, allows your function to give the correct result. Notice how your function overlays perfectly with the numpy polyval.

import numpy as np
import matplotlib.pyplot as plt     


def getYValueFromCoeff(f,coeff_list): # low to high order
    y_plot_values=[]
    for j in range(len(f)):
        item_list= []
        for i in range(len(coeff_list)):
            item= (coeff_list[i])*((f[j])**(len(coeff_list)-i-1))
            item_list.append(item)
        y_plot_values.append(sum(item_list))
    print(len(y_plot_values))
    return y_plot_values


no_of_coef = 10

#original signal
x = np.linspace(0, 0.01, 10)
period = 0.01
y = np.sin(np.pi * x / period)

#poly fit
coeffs = np.polyfit(x,y,no_of_coef)
coef_y = np.polyval(coeffs,x)

COEF_Y = getYValueFromCoeff(x,coeffs)

plt.figure()
plt.plot(x,y)
plt.plot(x, coef_y)
plt.plot(x, COEF_Y) 
plt.legend(['Original Function', 'Fitted Function', 'Custom Fitting'])
plt.show()

Output

enter image description here

  • thank you so much for appreciating my approach. and thanks for pointing out my mistake. I could only accept one answer, and the solution from @flawr points to a fundamental aspect of the poly.fit() function. Anyways, a massive thanks again :) – MRR Oct 27 '22 at 21:03
  • where are the coefficients printed from high to low? i printed test1.coef and this is what i get 0.9999939276745801 + 5.430942487918026e-15 x**1 - 1.2331182561455556 x**2 - 1.2411844215849982e-13 x**3 + 0.24595992358187177 x**4 + 5.754500392963442e-13 x**5 + 0.010857316629965433 x**6 - 9.144291551983254e-13 x**7 - 0.047377830552342154 x**8 + 4.577772389431185e-13 x**9 + 0.023684918811482034 x**10 [ 9.99993928e-01 5.43094249e-15 -1.23311826e+00 -1.24118442e-13 2.45959924e-01 5.75450039e-13 1.08573166e-02 -9.14429155e-13 -4.73778306e-02 4.57777239e-13 2.368e-02] – MRR Oct 27 '22 at 22:09
  • I'm referencing the numpy documentation for poly: https://numpy.org/doc/stable/reference/generated/numpy.poly.html. It returns the coeffs highest to lowest power in poly as well as polyfit. What your numbers mean is .9999939276745801*x^10 + 5.430942487918026e-15*x^9 + ... + 4.577772389431185e-13*x + 0.023684918811482034. Notice the coeff order is the same; however, the powers are opposite of what you have. This is consistent with the numpy documentation; moreover, the documentation explains that the highest power will always be 1 (or .999999... in this case). – greenerpastures Oct 28 '22 at 12:27
  • i am sorry but i still dont understand. 0.9999939276745801 x0+ 5.430942487918026e-15 x1 - 1.2331182561455556 x2... this equation is going from lowest to highest degree (power of x). then why did you write like this 9999939276745801*x^10 + 5.430942487918026e-15*x^9 + ... ? i think i am following this documentation https://numpy.org/doc/stable/reference/generated/numpy.polynomial.polynomial.Polynomial.fit.html – MRR Oct 28 '22 at 14:58
  • I see what has happened. I was confused by your import of Polynomial as 'poly'. There is a pre-existing function in numpy called 'poly' which is different. It seems they present the coeff vector opposite of one another. Your original function works with Polynomial. The function modification I suggested works with the built in 'poly' function. Sorry for the confusion. – greenerpastures Oct 28 '22 at 18:05
0

Here's the simple way of doing it if you didn't know that already...

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

no_of_coef= 10

#original signal
x = np.linspace(0, 0.01, 10)
period = 0.01
y = np.sin(np.pi * x / period)

#poly fit
test1= poly.fit(x,y,no_of_coef)

plt.plot(x, y, 'r', label='original y')
x = np.linspace(0, 0.01, 1000)
plt.plot(x, test1(x), 'b', label='y_fit')
plt.legend()
plt.show()

plot containing original signal and poly fit

hknjj
  • 319
  • 1
  • 8