Pythonのodeの数値解にデータを当てはめる

Question

私は2次のODEの系を持ちます。これは非線形であり、閉じた形で解析的に解くのは難しいです。私は数値解をこのODEシステムにデータセットに適合させたいと考えています。私のデータセットは、ODEシステムの一部である2つの変数のうちの1つのみです。これについてどうすればいいですか？ Thisは変数が1つしかないので役に立たなかった。現在、エラーにつながっている

私のコードは次のとおりです。

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() 

Answer 1

だから私はこの問題を考え出しました。 curve_fit（）関数は明らかに2番目の戻り値としてリストを返します。したがって、最初の条件をリスト[0.99,0.01]として渡す代わりに、私は0.99と0.01として別々に渡しました。

Answer 2

ではこのタイプのODEフィッティングが非常に簡単になりました。具体的には、ユーザーフレンドリーなラッパーとしてscipyに記述しました。私はあなたの状況に非常に便利だと思います。なぜなら、ボイラープレートコードの量が減ることで物事が簡単になるからです。ドキュメントから

し、あなたの問題にほぼ適用：

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() 

にStackOverflowより:)