Fitting SEIR and SIRS models to real data

Question

Fitting SEIR and SIRS models to real data

The study A SIR model assumption for the spread of COVID-19 in different communities uses a SIR model for prediction. I have taken the differential equations from the paper and transformed them to obtain the differential equations that define the SEIR and SIRS models.

SEIR model:

def deriv(r, beta, alpha, gamma):
    S, E, I, R = r
    dS = -beta*S*I
    dE = beta*S*I - alpha*E
    dI = alpha*E - gamma*I
    dR = gamma*I
    return np.array([dS, dE, dI, dR])

SIRS model:

def deriv(r, beta, gamma, rho):
    S, I, R = r
    dS = rho*R - beta*S*I
    dI = beta*S*I - gamma*I
    dR = gamma*I - rho*R
    return np.array([dS, dI, dR])

I want to compare the two models with each other and with the original model used in the paper. I want to do this by fitting all three models to the same coronal data and seeing which one fits the real data best (this of course requires a fit index). I solved the differential equations in all three cases using the Runge-Kutta method:

import numpy as np
import matplotlib.pyplot as plt

S0 = 1
I0 = 5*10**(-3)
R0 = 0
r = np.array([S0, I0, R0])

beta = 0.130
gamma = 0.048
rho = 0.005 # erre sincs egy konkrét, elméletileg 0 és 0,1 között vehet fel értékeket

def deriv(r, beta, gamma, rho):
    S, I, R = r
    dS = rho*R - beta*S*I
    dI = beta*S*I - gamma*I
    dR = gamma*I - rho*R
    return np.array([dS, dI, dR])

def RK4(r, h, beta, gamma, rho):
    k1 = h*deriv(r, beta, gamma, rho)
    k2 = h*deriv(r + 0.5*k1, beta, gamma, rho)
    k3 = h*deriv(r + 0.5*k2, beta, gamma, rho)
    k4 = h*deriv(r + k3, beta, gamma, rho)
    return r + (k1 + 2*k2 + 2*k3 + k4)/6

Ss = []
Is = []
Rs = []
dt = 1
T = 600
t = 0

while t < T:
    Ss.append(r[0])
    Is.append(r[1])
    Rs.append(r[2])
    r = RK4(r, dt, beta, gamma, rho)
    t += dt

plt.plot(range(T), Ss, label='Susceptible')
plt.plot(range(T), Is, label='Infected')
plt.plot(range(T), Rs, label='Recovered')
plt.xlabel('Idő (napokban mérve)')
plt.ylabel('A kategóriák megoszlása a társadalomban')
plt.title('SIRS-modell')
plt.legend()
plt.show()

This is an example of the SIRS model. I need a python code to compare the three models as described above.

I tried to fit it to the US data (several countries are included in the study).

The original study I used:

https://www.sciencedirect.com/science/article/pii/S0960077920304549