Python - curve fitting of more complex function

Question

I wish to find the equation of the curve of best fit of the following graph:

Which has the equation in the form of:

I've attempted to find examples of curve fitting with numpy here and here, but they all only show how to plot only exponential or only sinusoidal, but I'd like to plot a graph combining the two functions.

How would I do this?

You can create your own function with your equation and then use `curve_fit` from SciPy. [Here](https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html) you can see some examples. — Sheldore, Oct 20 '18 at 12:27

score 1 · Accepted Answer · answered Oct 24 '18 at 02:58

Here's one approach you might find useful. This uses lmfit (http://lmfit.github.io/lmfit-py/), which provides a high-level approach to curve fitting:

import numpy as np
import matplotlib.pyplot as plt

from lmfit import Model

def decay_cosine(t, amp, beta, omega, phi):
    """model data as decaying cosine wave"""
    return amp * np.exp(-beta*t)* np.cos(omega*t + phi)

# create fake data to be fitted
t = np.linspace(0, 5, 101)
y = decay_cosine(t, 1.4, 0.9, 7.2, 0.23) + np.random.normal(size=len(t), scale=0.05)

# build model from decay_cosine
mod = Model(decay_cosine)

# create parameters, giving initial values
params = mod.make_params(amp=2.0, beta=0.5, omega=5, phi=0)

# you can place bounds on parameters:
params['phi'].max = np.pi/2
params['phi'].min = -np.pi/2
params['amp'].min = 0

# fit data to model

result = mod.fit(y, params, t=t)

# print out fit results
print(result.fit_report())

# plot data with best fit
plt.plot(t, y, 'bo', label='data')
plt.plot(t, result.best_fit, 'r')
plt.show()

This will print out a report like this:

[[Model]]
    Model(decay_cosine)
[[Fit Statistics]]
    # fitting method   = leastsq
    # function evals   = 46
    # data points      = 101
    # variables        = 4
    chi-square         = 0.25540159
    reduced chi-square = 0.00263301
    Akaike info crit   = -595.983903
    Bayesian info crit = -585.523421
[[Variables]]
    amp:    1.38812335 +/- 0.03034640 (2.19%) (init = 2)
    beta:   0.90760648 +/- 0.02820705 (3.11%) (init = 0.5)
    omega:  7.16579292 +/- 0.02891827 (0.40%) (init = 5)
    phi:    0.26249321 +/- 0.02225816 (8.48%) (init = 0)
[[Correlations]] (unreported correlations are < 0.100)
    C(omega, phi)  = -0.713
    C(amp, beta)   =  0.695
    C(amp, phi)    =  0.253
    C(amp, omega)  = -0.183
    C(beta, phi)   =  0.178
    C(beta, omega) = -0.128

and produce a plot like this:

Thank you, this is just what I was looking for. However, the error `ImportError: cannot import name 'Model' from 'lmfit'` popped up when I tried to run an exact copy of your code. Do you know what the issue could be? I've installed lmfit v. 0.9.11 via `pip` — George Tian, Oct 25 '18 at 09:17
Hm, sounds like it may not have been installed correctly. Can you do `import lmfit` in a python session? You may want to check where the module actually got installed. — M Newville, Oct 26 '18 at 00:19
Oh turns out I had named my file `lmfit.py` and it caused the error. Renaming the file worked. — George Tian, Oct 26 '18 at 02:58

score 0 · Answer 2 · answered Aug 18 '20 at 15:23

Here is a quite simple example using curve_fit and leastsq from scipy.optimize.

1. Setting parameter values, model and experimental data.

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

np.random.seed(0)  # choosing seed for reproducibility

# ==== theoretical parameter values ====
x0 = 1
beta = .5
omega = 2*np.pi
phi = 0
params = x0, beta, omega, phi

# ==== model ====
def decay_cosine(t, x0, beta, omega, phi):
    x = x0 * np.exp(-beta*t_data) * np.cos(omega*t_data + phi)
    return x

# ==== generating experimental data ====
t_data = np.linspace(0, 5, num=80)
noise = .05 * np.random.randn(t_data.size)
x_data = decay_cosine(t_data, *params) + noise

2. Fitting.

# ==== fitting using curve_fit ====
params_cf, _ = scipy.optimize.curve_fit(decay_cosine, t_data, x_data)

# ==== fitting using leastsq ====
def residuals(args, t, x):
    return x - decay_cosine(t, *args)

x0 = np.ones(len(params))  # initializing all params at one
params_lsq, _ = scipy.optimize.leastsq(residuals, x0, args=(t_data, x_data))

print(params_cf)
print(params_lsq)

array([ 1.04938794,  0.53877389,  6.30375113, -0.01850761])
array([ 1.04938796,  0.53877389,  6.30375103, -0.01850744])

3. Plotting.

plt.plot(t_data, x_data, '.', label='exp data')
plt.plot(t_data, decay_cosine(t_data, *params_cf), label='curve_fit')
plt.plot(t_data, decay_cosine(t_data, *params_lsq), '--', label='leastsq')
plt.legend()
plt.grid(True)
plt.show()

Python - curve fitting of more complex function

2 Answers2