CODEWITHSUNDEEP
pythonarraysvectorizationode
J Wright
J Wright

Reputation: 21

Lopping over an array of parameter values

I am fairly new to coding in python, so any help would be appreciated and if clarification is needed please let me know.

Below is my working code. It is an SIR model that contains a system of 7 first order differential equations with some parameter values and initial conditions.

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

I would like to modify my code so that I can have an array of values for all of my parameters instead of one numeric value. Just as an example, something like this:

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

I intend to vectorize the calculations with respect to the parameters, so instead of a single solution at each time step, I get a collection of solutions. That being said, my question is how do I properly loop over the parameter values to get such a solution?

Upvotes: 0

Views: 56

Answers (1)

colmo007
colmo007

Reputation: 38

Try having a new array to store every loop result, and a for loop to get every parameter.

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)

Here in each cycle the actual_x variables will contain the desired parameter.

Also yo need to remove the () from your arrays, because there you have one list(array) inside a tuple and that's redundant. Insted you should have:

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]

A more elegant way would be to have a simple array containing each parameter for each cycle like:

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)

Upvotes: 1

Related Questions