Modified SEIR Python Model

Modified-SEIR-Model.py

from random import uniform

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

# Define constants
n_countries = 100
t_start = 0
t_end = 150
beta = 0.5  # Transmission rate
sigma = 0.04  # Basic reproduction number
mu = 0.04  # Recovery rate
_lambda = np.array([0] * n_countries)
s0 = 0.99
e0 = 0.01
i0 = 0.0
r0 = 0.0


def seir_function(_, y):
    s, e, i, r = y
    return np.array([-beta * i * s,
                     beta * i * s - sigma * e,
                     sigma * e - mu * i + np.sum(_lambda[:n_countries - 1]),
                     mu * i])


def solve():
    return solve_ivp(seir_function, [t_start, t_end], [s0, e0, i0, r0], rtol=1e-6)

# Solve
reference = solve()
sols = []
for _ in range(100):
    _lambda = np.array([uniform(0, 0.00005) for _ in range(n_countries)])
    sols.append(solve())

# Create figure
fig = plt.figure(figsize=(10, 10))
ax = plt.axes()

# Plot lines
for sol in sols:
    ax.plot(sol.t, sol.y[2], color="yellow")
ax.plot(reference.t, reference.y[2], color="black")

# Show figure
plt.show()