Search code examples
pythonnumpyscipycurve-fittingdata-fitting

Exponential decay fitting in Scipy


I have two NumPy arrays x_data and y_data. When I try to fit my data using exponential decay function and scipy.optimize.curve_fit with this simple code:

import numpy as np
from scipy.optimize import curve_fit
import plotly.graph_objects as go

# Define the exponential decay function
def exp_decay(x, a, b, c):
    return a * np.exp(-b * x) + c

# Example data points
x_data = np.array([50.6733, 50.6736, 50.6739, 50.6742, 50.6744, 50.6747, 50.675,
       50.6753, 50.6756, 50.6758, 50.6761, 50.6764, 50.6767, 50.6769,
       50.6772, 50.6775, 50.6778, 50.6781, 50.6783, 50.6786, 50.6789,
       50.6792, 50.6794, 50.6797, 50.68  , 50.6803, 50.6808, 50.6808,
       50.6811, 50.6814, 50.6817, 50.6819, 50.6822, 50.6825, 50.6828,
       50.6831, 50.6833])
y_data = np.array([3.039, 3.035, 3.031, 3.025, 3.021, 3.017, 3.013, 3.008, 3.003,
       3.001, 2.998, 2.996, 2.994, 2.992, 2.99 , 2.988, 2.986, 2.985,
       2.984, 2.983, 2.982, 2.981, 2.98 , 2.979, 2.978, 2.978, 2.978,
       2.978, 2.977, 2.977, 2.976, 2.976, 2.976, 2.975, 2.975, 2.975,
       2.975])

# Fit the exponential decay model to the data
initial_guess = [1, 0.000001, 1]  # Initial guess for the parameters [a, b, c]
params, covariance = curve_fit(exp_decay, x_data, y_data, p0=initial_guess)

# Extract the fitted parameters
a_fitted, b_fitted, c_fitted = params
print("Fitted parameters: a =", a_fitted, ", b =", b_fitted, ", c =", c_fitted)


###Plotting data vs fitting
# Prepare the fitted curve data
x_fit = np.linspace(min(x_data), max(x_data), 100)
y_fit = exp_decay(x_fit, *params)

# Create Plotly plot
fig = go.Figure()

# Add scatter plot for data points
fig.add_trace(go.Scatter(x=x_data, y=y_data, mode='markers', name='Data Points'))

# Add line plot for the fitted model
fig.add_trace(go.Scatter(x=x_fit, y=y_fit, mode='lines', name='Fitted Model'))

# Update plot layout
fig.update_layout(title='Fit of Exponential Decay Model to Data',
                  xaxis_title='x',
                  yaxis_title='y',
                  legend_title='Legend')

# Show the figure
fig.show()

I get a fit the following fit enter image description here

and nonsense fitting parameters

Fitted parameters: a = 796.8340789795466 , b = 0.02189572018876207 , c = -259.70667935123606

What is the problem?


Solution

  • It's easier to fit the data when there isn't such a large disparity in the values. What you can do is offset the x_data so that it starts at 0 and fit that. It will be easier to get a good starting guess at that point (I did it through a couple of tries, but you can create a widget with sliders so you can find a good guess more easily). Then you can generate y_fit with the shifted x_data but plot it with the actual x_data.

    import numpy as np
    from scipy.optimize import curve_fit
    import matplotlib.pyplot as plt
    
    plt.close("all")
    
    
    def exp_decay(x, a, b, c):
        return a * np.exp(-b * x) + c
    
    
    x_data = np.array([50.6733, 50.6736, 50.6739, 50.6742, 50.6744, 50.6747, 50.675,
                       50.6753, 50.6756, 50.6758, 50.6761, 50.6764, 50.6767, 50.6769,
                       50.6772, 50.6775, 50.6778, 50.6781, 50.6783, 50.6786, 50.6789,
                       50.6792, 50.6794, 50.6797, 50.68, 50.6803, 50.6808, 50.6808,
                       50.6811, 50.6814, 50.6817, 50.6819, 50.6822, 50.6825, 50.6828,
                       50.6831, 50.6833])
    
    y_data = np.array([3.039, 3.035, 3.031, 3.025, 3.021, 3.017, 3.013, 3.008, 3.003,
                       3.001, 2.998, 2.996, 2.994, 2.992, 2.99, 2.988, 2.986, 2.985,
                       2.984, 2.983, 2.982, 2.981, 2.98, 2.979, 2.978, 2.978, 2.978,
                       2.978, 2.977, 2.977, 2.976, 2.976, 2.976, 2.975, 2.975, 2.975,
                       2.975])
    
    x_offset = x_data.min()
    x_data_offset = x_data - x_offset
    
    initial_guess = [0.07, 300, 3]
    params, covariance = curve_fit(
        exp_decay, x_data_offset, y_data, p0=initial_guess)
    
    a_fitted, b_fitted, c_fitted = params
    print(f"Fitted parameters: a = {a_fitted}, b = {b_fitted}, c = {c_fitted}")
    
    
    x_fit = np.linspace(min(x_data), max(x_data), 100)
    # generate y_fit with shifted x_data
    y_fit = exp_decay(x_fit - x_offset, *params)
    
    fig, ax = plt.subplots()
    ax.plot(x_data, y_data, ".")
    ax.plot(x_fit, y_fit)
    
    Fitted parameters: a = 0.07022304888360707, b = 344.53511176973166, c = 2.9720506329733154