Fitting of experimental data within two different regions

Question

I am fitting a set of experimental data (sample) within two different experimental regions and can be expressed with two mathematical functions as follows:

1st region:

y = m*x + c ( the slope can be constrained to zero)

2nd region:

y = d*exp(-k*x)

the experimental data is shown below and I coded it in python as follows:

def func(x, m, c, d, k):
   return m*x+ c + d*np.exp(-k*x) 
popt, pcov = curve_fit(func, t, y)

Unfortunately, my data is not fitting properly and fitted (returned) parameters do not make sense (see picture below).

Any assistance will be appreciated.

Answer 1

Very interesting question. As said by a_guest, you will have to fit to the two regions separately. However, I think you probably also want the two regions to connect smoothly at the point t0, the point where we switch from one model to the other. In order to do this, we need to add the constraint that y1 == y2 at the point t0.

In order to do this with scipy, look at scipy.optimize.minimize with the SLSQP method. However, I wrote a scipy wrapper to make this kind of thing easier, called symfit. I will show you how to do this with symfit, because I think it's better suited to the task, but with this example you should also be able to implement it with pure scipy if you prefer.

from symfit import parameters, variables, Fit, Piecewise, exp, Eq
import numpy as np
import matplotlib.pyplot as plt

t, y = variables('t, y')
m, c, d, k, t0 = parameters('m, c, d, k, t0')

# Help the fit by bounding the switchpoint between the models
t0.min = 0.6
t0.max = 0.9

# Make a piecewise model
y1 = m * t + c
y2 = d * exp(- k * t)
model = {y: Piecewise((y1, t <= t0), (y2, t > t0))}

# As a constraint, we demand equality between the two models at the point t0
# to do this, we substitute t -> t0 and demand equality using `Eq`
constraints = [Eq(y1.subs({t: t0}), y2.subs({t: t0}))]

# Read the data
tdata, ydata = np.genfromtxt('Experimental Data.csv', delimiter=',', skip_header=1).T

fit = Fit(model, t=tdata, y=ydata, constraints=constraints)
fit_result = fit.execute()
print(fit_result)

plt.scatter(tdata, ydata)
plt.plot(tdata, fit.model(t=tdata, **fit_result.params).y)
plt.show()