Блог

Я довольно новичок в программировании на python, поэтому буду признателен за любую помощь, и если потребуется разъяснение, пожалуйста, дайте мне знать.

Ниже приведен мой рабочий код. Это модель SIR, которая содержит систему из 7 дифференциальных уравнений первого порядка с некоторыми значениями параметров и начальными условиями.

 import scipy.integrate
import matplotlib.pyplot as plt
import numpy as np

def ode(t, y):
    sigma1 = 1/10
    sigma2 = 1/4
    alpha = 1/10.4
    gamma = 1/5
    f = 0.18122
    h = 0.24045
    dh = 0.06218
    d = 0.02591
    beta = 0.4
    N = 500000

    return np.array([-(beta*(y[2] y[3])/N)*y[0],
                     (beta*(y[2] y[3])/N)*y[0]-gamma*y[1],
                     f*gamma*y[1]-sigma1*y[2],
                     (1-f)*gamma*y[1]-sigma2*y[3],
                     h*sigma1*y[2]-alpha*y[4],
                     (1-h-d)*sigma1*y[2] sigma2*y[3] (1-dh)*alpha*y[4], 
                     d*sigma1*y[2] dh*alpha*y[4]]) 
t0 = 0 
t_bound = 100 
y0 = np.array([480000,0,10000,10000,0,0,0])
sol = scipy.integrate.RK45(ode, t0, y0, t_bound)
t = []
y = []
while sol.status == "running":
    t.append(sol.t)
    y.append(sol.y)
    sol.step()

plt.plot(np.array(t), np.array(y))
plt.legend(("susceptible", "exposed","infectious_s", "infectious_a", "hospitalized", "recovered", "deaths"))
plt.xlabel("time")
plt.ylabel("cases")
plt.show()
 

Я хотел бы изменить свой код так, чтобы я мог иметь массив значений для всех моих параметров вместо одного числового значения. Просто в качестве примера, что-то вроде этого:

     h = ([0.182,0.055,0.055,0.055,0.068,0.068,0.139,0.139,0.139,0.139,0.251,0.251,0.251,0.512,0.512,0.512,0.617])
    dh = ([0.002, 0, 0, 0, 0.002, 0.002, 0.009, 0.009, 0.009, 0.009, 0.036, 0.036, 0.036, 0.149, 0.149, 0.149, 0.328])
    d = ([0.001, 0, 0, 0, 0.001, 0.001, 0.004, 0.004, 0.004, 0.004, 0.014, 0.014, 0.014, 0.059, 0.059, 0.059, 0.129])
 

Я намерен векторизовать вычисления относительно параметров, поэтому вместо одного решения на каждом временном шаге я получаю набор решений. При этом мой вопрос заключается в том, как мне правильно перебирать значения параметров, чтобы получить такое решение?

Попробуйте создать новый массив для хранения каждого результата цикла и цикл for для получения каждого параметра.

 solutions_array = []

for index in range(len(h)):
    actual_h = h[index]
    actual_dh = dh[index]
    actual_d = d[index]
    # The rest of your calculations
    # End of calculation in variable result
    solutions_array.append(result)
 

Здесь в каждом цикле actual_x переменные будут содержать требуемый параметр.

Также вам нужно удалить () из ваших массивов, потому что там у вас есть один список (массив) внутри кортежа, и это избыточно. Вместо этого у вас должно быть:

 h = [0.182,0.055,0.055,0.055,0.068,0.068,0.139,0.139,0.139,0.139,0.251,0.251,0.251,0.512,0.512,0.512,0.617]
dh = [0.002, 0, 0, 0, 0.002, 0.002, 0.009, 0.009, 0.009, 0.009, 0.036, 0.036, 0.036, 0.149, 0.149, 0.149, 0.328]
d = [0.001, 0, 0, 0, 0.001, 0.001, 0.004, 0.004, 0.004, 0.004, 0.014, 0.014, 0.014, 0.059, 0.059, 0.059, 0.129]
 

Более элегантным способом было бы иметь простой массив, содержащий каждый параметр для каждого цикла, например:

 parameters = [[0.182, 0.02, 0.001],[0.055, 0, 0],[0.055, 0, 0],[0.068, 0.002, 0]] #[h, dh, d]

solutions_array = []
for h, dh, d in parameters:
    # The rest of your calculations
    # End of calculation in variable result
    solutions_array.append(result)