Failing a simple Cosine fit in Python

Question

Here's how I generate my data and the tried fit:

import matplotlib.pyplot as plt
from scipy import optimize
import numpy as np

def f(t,a,b):
    return a*np.cos(b*t)

v = 0 
x = 0.03

t = 0
dt = 0.001

time = []
pos = []

while t<3:
    a = (-5*x)/0.1
    v = v + a*dt
    x = x + v*dt
    time.append(t)
    pos.append(x)
    t = t+dt

pop, pcov = optimize.curve_fit(f,time,pos)
print(pop)

Even when I indicate initial values for the parameters (such as 0.03 for "a" and "7" for b), the resulting fit is still way off (see below, dashed line is the fit function).

Am I using the wrong library? or have I made an obvious blunder?

Thanks for any hints.

user_3014199 · Answer

As Tyberius noted, you need to provide better initial values. Why is that? optimize.curve_fit uses least_squares which finds a local minimum of the cost function. I believe in your case you are stuck in such a local minimum (that is not the global minimum). If you look at your diagram, your fit is approximately y=0. (It is a bit wavy because it is a cosine) If you were to increase a a bit the error would go up, so a stays close to zero. And if you were to increase b to fit the frequency of the data, the cost function would go up as well so that one stays low as well.

If you don't provide initial values, the parameters start at 1 each so it looks like this:

plt.plot(time, pos, 'black', label="data")
a,b = 1,1
init = [a*np.cos(b*t) for t in time]
plt.plot(time, init, 'b', label="a,b=1,1")
plt.legend()
plt.show()

a will go down and b will stay behind. I believe the scale is an additional problem. If you normalized your data to have an amplitude of 1 the humps might be more pronounced and easier to fit. If you start with a convenient value for a, b can find its way from an initial value as low as 5:

plt.plot(time, pos, 'black', label="data")
for i in [1, 4.8, 4.9, 5]:
    pop, pcov = optimize.curve_fit(f,time,pos, p0=(0.035,i))
    a,b = pop
    fit = [a*np.cos(b*t) for t in time]
    plt.plot(time, fit, label=f"$b_0 = {i}$")

plt.legend()
plt.show()

Failing a simple Cosine fit in Python

1 Answers