Python - Scipy and Numpy don't get along - problem using numpy arrays in scipy.optimize.curve

I am trying to fit some experimental data (x and y) with a custom function (Srt) and using scipy.optimize.curve_fit():

Reading the data and defining the function, using dummy values (10,10) for Km and Vmax (which are to be determined using the curve fit) works fine, as long as I use np.asarray():

from scipy.special import lambertw
from scipy.optimize import curve_fit
import numpy as np
import scipy

def Srt(t,s,Km,Vmax):
    print("t",type(t))
    print("t",t)
    print("last element of t:",t[-1])
    print("s",type(s))
    print("s",s)
    print("last element of s:",s[-1])
    Smax = s[-1] # Substrate concentration at end of reaction 
    t0 = t[0]    # time=0 (beginning of reaction)
    s0 = s[0]    # Substrate concentration at time = 0 (beginning of reaction)
    E = np.exp(((Smax - s0) - Vmax*(t t0))/Km)
    L = lambertw(((Smax - s0)/Km)*E)
    y = Smax - Km*L
    return y

x=[2.780000e-03,2.778000e-02,5.278000e-02,7.778000e-02,1.027800e-01
,1.277800e-01,1.527800e-01,1.777800e-01,2.027800e-01,2.277800e-01
,2.527800e-01,2.777800e-01,3.027800e-01,3.277800e-01,3.527800e-01]

y=[0.44236,0.4308,0.42299,0.41427,0.40548,0.39908,0.39039,0.3845,0.37882
,0.37411,0.36759,0.36434,0.35864,0.35508,0.35138]

xdata = np.asarray(x)
ydata = np.asarray(y)
Srt(xdata, ydata,10,10)

If I do not use np.asarray, I get a "Type Error":

Srt(x, y,10,10)

When I continue to use curve_fit to make the fit for Vmax and Km with:

parameters, covariance = scipy.optimize.curve_fit(Srt, xdata, ydata)

I get into trouble:

If I understand the error message correctly, for some reason, the array ydata is not an array anymore when it is read in as s ?!?

What do I have to change in my code so that I can work with my function Srt and curve_fit ?

CodePudding user response：

Please have a closer look at the documentation of the curve_fit function. Where it states that ydata must nominaly be the result of func(xdata... ). So the ydata that you hand to curve_fit is never passed as argument of the call of Srt as you indicated in the manual call.

Furthermore, the parameters to be estimated must have the same shape, which means that you have to define Smax and s0 as float input. I modified your example such that it actually runs:

from scipy.special import lambertw
from scipy.optimize import curve_fit
import numpy as np
import scipy

def Srt(t, Smax: float, s0: float, Km: float, Vmax: float):
    t0 = t[0]    # time=0 (beginning of reaction)
    E = np.exp(((Smax - s0) - Vmax*(t t0))/Km)
    # L = lambertw(((Smax - s0)/Km)*E)  # this apparently can be complex which causes another Error
    L = np.abs(lambertw(((Smax - s0)/Km)*E))
    y = Smax - Km*L
    return y

x=[2.780000e-03,2.778000e-02,5.278000e-02,7.778000e-02,1.027800e-01
,1.277800e-01,1.527800e-01,1.777800e-01,2.027800e-01,2.277800e-01
,2.527800e-01,2.777800e-01,3.027800e-01,3.277800e-01,3.527800e-01]

y=[0.44236,0.4308,0.42299,0.41427,0.40548,0.39908,0.39039,0.3845,0.37882
,0.37411,0.36759,0.36434,0.35864,0.35508,0.35138]

xdata = np.array(x)
ydata = np.array(y)

parameters, covariance = scipy.optimize.curve_fit(Srt, xdata, ydata)

NOTE:
The np.abs inside the function does not make sense, but the complex result of lambertw apparently can be complex. In this case an error is raised as there is no safe casting rule, causing curvefit to abort.

CodePudding user response：

Your first error is produced by the t t0 expression. It t is a list x, that's a list "concatenate" expression, which is fine for [1,2,3] [4,5] but not [1,2,3] 5. That's why x and y have to arrays.

In the second error, what did the

print("s",type(s))
print("s",s)

show? Apparently s is not an array, or even a list.