Finding params using leastsq method for levenberg marquardt

Question

I am trying to plot the data for my leastsq model and equation as I keep getting shape errors. I'm just trying to guess the parameters that fit the data but I can't do that if I can't see the graphs.

Here is what I have so far:

from scipy.integrate import odeint
Import numpy as np
Import matplotlib.pyplot as plt

#data
didata1 = np.loadtxt("diauxic1.txt")
time = didata1[:,0]
pop2 = didata1[:,1]

# model equations
def diauxic_ode(x,t,params):
    r1,r2,k = params  
    y,S1,S2 = x
    derivs = [r1*S1*y+(k/(k+S1))*r2*S2*y, -r1*S1*y, -(k/(k+S1))*r2*S2*y]
    return derivs

# runs a simulation and returns the population size
def diauxic_run(pars,t):
    r1,r2,k,y0,S10,S20 = pars
    ode_params=[r1,r2,k]
    ode_starts=[y0,S10,S20]
    out = odeint(diauxic_ode, ode_starts, t, args=(ode_params,))
    return out[:,0] 

# residual function
def diauxic_resid(pars,t,data):
    r1,r2,k,y0,S10,S20 = pars
    ode_params=[r1,r2,k]
    ode_starts=[y0,S10,S20]
    out = odeint(diauxic_ode, ode_starts, t, args=(ode_params,))
    return diauxic_run(pars,t)-data


p0 =[1,1,1,1,1,1]
lsq_out = leastsq(diauxic_resid, p0, args=(time,pop2))

plt.plot(time,pop2,'o',time,diauxic_resid(p0,time,lsq_out[0]))

plt.show()

Answer 1

The immediate error is the call diauxic_resid(p0,time,lsq_out[0]).

• Do you really want to plot the original data and the errors in the same plot?
• The arguments are wrong, both for resid and run. That something is not right is already clear in having the initial point p0 in the arguments where you want to plot the adapted result.

Thus replace with diauxic_run(lsq_out[0],time). Or even better, split the plot commands and increase the sample density for the curve

plt.plot(time,pop2,'o');
time = np.linspace(time[0], time[-1], len(time)*10-9)
plt.plot(time,diauxic_run(lsq_out[0], time))

from test data generated by

ptest = [0.5, 0.2, 20, 0.2,3.,1.5]
time = np.linspace(0,10,21);
pop2 = diauxic_run(ptest, time)+ np.random.randn(len(time))*0.01

leading to fitted parameters

lsq_out[0]: [ 0.23199391  0.5998453  20.67961621  0.19636029  2.16841159  2.32688635]