Подгонка данных к численному решению оды в питоне

У меня есть система из двух ОДУ первого порядка, которые нелинейны и, следовательно, их трудно решить аналитически в замкнутой форме. Я хочу подогнать численное решение этой системы ОДУ к набору данных. Мой набор данных предназначен только для одной из двух переменных, которые являются частью системы ODE. Как мне это сделать? Это не помогло, потому что там только одна переменная.

Мой код, который в настоящее время приводит к ошибке:

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

python numpy scipy data-analysis curve-fitting
person Arkya    schedule 23.10.2016    source источник
comment
Добро пожаловать в Stackoverflow. Мы не будем писать код за вас, но если вы поделитесь с нами ошибкой, которую видите, мы сможем помочь.   -  person Tammo Heeren    schedule 23.10.2016

Ответы (2)


arrow_upward
2
arrow_downward

Этот тип подгонки ODE становится намного проще в symfit, который Я написал специально как удобную оболочку для 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()

Подробнее см. в документах :)

person tBuLi    schedule 23.10.2016

arrow_upward
0
arrow_downward

Итак, я понял проблему. Функция curve_fit(), по-видимому, возвращает список в качестве второго возвращаемого значения. Итак, вместо того, чтобы передавать начальные условия в виде списка [0,99,0,01], я передал их отдельно как 0,99 и 0,01.

person Arkya    schedule 23.10.2016