在python中将数据拟合到ode的数值解 [英] Fitting data to numerical solution of an ode in python

问题描述

我有一个由两个一阶ODE组成的系统，它们是非线性的，因此很难以封闭形式解析地求解.我想将这个ODE系统的数值解拟合到一个数据集.我的数据集仅适用于ODE系统中两个变量之一.我该怎么办?这并没有帮助，因为那里只有一个变量.

我当前导致错误的代码是:

 将numpy导入为np从scipy.integrate导入odeint从scipy.optimize导入curve_fitdef f(y，t，a，b，g):S，我= y#S，我应该是我的变量Sdot = -a * S * IIdot =(a-b)* S * I +(b-g-b * I)* Idydt = [Sdot，Idot]返还染料def y(t，a，b，g，y0):y = odeint(f，y0，t，args =(a，b，g))返回y.ravel()I_data = []#我只有我的数据，没有给S的数据文件=打开('./ratings_showdown.csv')对于file.read().split('\ r \ n')中的e_raw:尝试:e = float(e_raw);I_data.append(e)除了ValueError:继续data_t =范围(len(I_data))popt，cov = curve_fit(y，data_t，I_data，[.05，0.02，0.01，[0.99,0.01]])#想让我适合变量I的数据解决方案的一部分#ERROR在这里，ValueError:使用序列设置数组元素a_opt，b_opt，g_opt，y0_opt = poptprint("a =％g"％a_opt)print("b =％g"％b_opt)print("g =％g"％g_opt)打印("y0 =％g"％y0_opt)导入matplotlib.pyplot作为pltt = np.linspace(0，len(data_y)，2000)plt.plot(data_t，data_y，'.'，t，y(t，a_opt，b_opt，g_opt，y0_opt)，'-')plt.gcf().set_size_inches(6，4)#plt.savefig('out.png'，dpi = 96)#保存拟合结果plt.show() 

解决方案

在文档了解更多:)

I have a system of two first order ODEs, which are nonlinear, and hence difficult to solve analytically in a closed form. I want to fit the numerical solution to this system of ODEs to a data set. My data set is for only one of the two variables that are part of the ODE system. How do I go about this? This didn't help because there's only one variable there.

My code which is currently leading to an error is:

import numpy as np
from scipy.integrate import odeint
from scipy.optimize import curve_fit

def f(y, t, a, b, g):
    S, I = y # S, I are supposed to be my variables
    Sdot = -a * S * I
    Idot = (a - b) * S * I + (b - g - b * I) * I
    dydt = [Sdot, Idot]
    return dydt

def y(t, a, b, g, y0):
    y = odeint(f, y0, t, args=(a, b, g))
    return y.ravel()

I_data =[] # I have data only for I, not for S
file = open('./ratings_showdown.csv')
for e_raw in file.read().split('\r\n'):
    try:
        e=float(e_raw); I_data.append(e)
    except ValueError:
        continue

data_t = range(len(I_data))
popt, cov = curve_fit(y, data_t, I_data, [.05, 0.02, 0.01, [0.99,0.01]]) 
#want to fit I part of solution to data for variable I
#ERROR here, ValueError: setting an array element with a sequence
a_opt, b_opt, g_opt, y0_opt = popt

print("a = %g" % a_opt)
print("b = %g" % b_opt)
print("g = %g" % g_opt)
print("y0 = %g" % y0_opt)

import matplotlib.pyplot as plt
t = np.linspace(0, len(data_y), 2000)
plt.plot(data_t, data_y, '.',
         t, y(t, a_opt, b_opt, g_opt, y0_opt), '-')
plt.gcf().set_size_inches(6, 4)
#plt.savefig('out.png', dpi=96) #to save the fit result
plt.show()

解决方案

This type of ODE fitting becomes a lot easier in symfit, which I wrote specifically as a user friendly wrapper to scipy. I think it will be very useful for your situation because the decreased amount of boiler-plate code simplifies things a lot.

From the docs and applied roughly to your problem:

from symfit import variables, parameters, Fit, D, ODEModel

S, I, t = variables('S, I, t')
a, b, g = parameters('a, b, g')

model_dict = {
    D(S, t): -a * S * I,
    D(I, t): (a - b) * S * I + (b - g - b * I) * I
}

ode_model = ODEModel(model_dict, initial={t: 0.0, S: 0.99, I: 0.01})

fit = Fit(ode_model, t=tdata, I=I_data, S=None)
fit_result = fit.execute()

Check out the docs for more :)

这篇关于在python中将数据拟合到ode的数值解的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！