使用 leavenberg marquardt 的 leastsq 方法查找参数
Finding params using leastsq method for levenberg marquardt
我正在尝试为我的 leastsq 模型和方程绘制数据,因为我不断收到形状错误。我只是想猜测适合数据的参数,但如果看不到图表,我就做不到。
这是我目前的情况:
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()
直接错误是调用diauxic_resid(p0,time,lsq_out[0])
。
- 您真的要在同一图中绘制原始数据和错误吗?
resid
和 run
的参数都是错误的。在要绘制改编结果的参数中有初始点 p0
已经很清楚了,有些事情是不对的。
因此替换为diauxic_run(lsq_out[0],time)
。或者更好的是,拆分绘图命令并增加曲线的样本密度
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))
来自
生成的测试数据
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
导致拟合参数
lsq_out[0]: [ 0.23199391 0.5998453 20.67961621 0.19636029 2.16841159 2.32688635]
我正在尝试为我的 leastsq 模型和方程绘制数据,因为我不断收到形状错误。我只是想猜测适合数据的参数,但如果看不到图表,我就做不到。
这是我目前的情况:
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()
直接错误是调用diauxic_resid(p0,time,lsq_out[0])
。
- 您真的要在同一图中绘制原始数据和错误吗?
resid
和run
的参数都是错误的。在要绘制改编结果的参数中有初始点p0
已经很清楚了,有些事情是不对的。
因此替换为diauxic_run(lsq_out[0],time)
。或者更好的是,拆分绘图命令并增加曲线的样本密度
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))
来自
生成的测试数据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
导致拟合参数
lsq_out[0]: [ 0.23199391 0.5998453 20.67961621 0.19636029 2.16841159 2.32688635]