À propos du modèle SIR

Lorsque la déclaration d'urgence due à l'épidémie de COVID-19 a été déclarée le 7 avril 2020, le Dr Nishiura de l'Université de Kita a montré un graphique montrant l'estimation lorsque le taux de contact diminuait de 20%, 60% et 80%. Était là. En regardant ce qu'est le "modèle SIR" utilisé actuellement, il existe des exemples de calculs en Python et R, et la même expression peut être implémentée dans OpenModelica. Quand je l'ai essayé, j'ai pu le simuler, donc je vais l'écrire sous forme d'article.

Veuillez noter que les résultats obtenus à partir de ce calcul ne reflètent pas la situation réelle du COVID-19, et ne sont modélisés que par des ingénieurs qui ne sont pas des professionnels de la santé avec des paramètres temporaires. Mettre.

En modélisation, j'ai fait référence aux sites et articles suivants.

** Introduction d'un modèle de prédiction mathématique pour les maladies infectieuses (modèle SIR) ** https://qiita.com/kotai2003/items/3078f4095c3e94e5325c ** Équations différentielles et épidémies mathématiques de maladies infectieuses ** https://www.ms.u-tokyo.ac.jp/~inaba/inaba_science_2008.pdf

Implémenter le modèle SIR de base avec OpenModelica

Le modèle SIR est représenté par les trois équations différentielles suivantes.

// * Extrait de "Introduction of Mathematical Prediction Model for Infectious Diseases (SIR Model)" *

\begin{align}

\frac{dS}{dt} &= -\beta SI \\
\frac{dI}{dt} &=  \beta SI -\gamma I \\
\frac{dR}{dt} &=  \gamma I \\

\end{align}

\begin{align}
S &:Personne infectable\quad \text{(Susceptible)}　\\
I &:Personne infectée\quad \text{(Infectious)}　\\
R &:Ceux qui sont décédés après une infection ou qui ont acquis une immunité\quad \text{(Removed)}　\\
\beta &:Taux d'infection\quad \text{(The infectious rate)} \quad [1/day] \\
\gamma &:Taux d'exclusion\quad \text{(The Recovery rate)} \quad [1/day]　\\
\end{align}

Prérequis pour ce modèle SIR

Une personne qui a acquis l'immunité ne sera jamais infectée et ne perdra pas l'immunité. --Il n'y a pas d'entrée ou de sortie de l'extérieur dans la population totale. Personne n'est mort pour des causes autres que la naissance et l'infection.

// Citer jusqu'à ici

Utilisons cette formule pour implémenter le modèle SIR de base sur OpenModelica.

Les taux d'infection β et γ sont définis comme des paramètres, les valeurs initiales de S, I et R sont fournies en tant qu'entrées, et les valeurs de changement de S, I et R sont fournies en tant que sorties. Voici le code.

class cl_SIRsimple
  extends Modelica.Blocks.Icons.Block;
  parameter Real contact_rate = 0.5 / 1000 "Infectious rate (1/day)";
  parameter Real gamma_const = 0.2 "Recovery rate (1/day)";
  Modelica.Blocks.Interfaces.RealInput Si "Connector of initial S(Susceptible)" annotation(
    Placement(visible = true, transformation(origin = {-100, 60}, extent = {{-20, -20}, {20, 20}}, rotation = 0), iconTransformation(origin = {-120, 60}, extent = {{-20, -20}, {20, 20}}, rotation = 0)));
  Modelica.Blocks.Interfaces.RealInput Ii "Connector of initial I (Infectious)" annotation(
    Placement(visible = true, transformation(origin = {-100, 0}, extent = {{-20, -20}, {20, 20}}, rotation = 0), iconTransformation(origin = {-120, 0}, extent = {{-20, -20}, {20, 20}}, rotation = 0)));
  Modelica.Blocks.Interfaces.RealInput Ri "Connector of initial R (Removed)" annotation(
    Placement(visible = true, transformation(origin = {-100, -62}, extent = {{-20, -20}, {20, 20}}, rotation = 0), iconTransformation(origin = {-120, -60}, extent = {{-20, -20}, {20, 20}}, rotation = 0)));
  Modelica.Blocks.Interfaces.RealOutput Sy "Output connector of S(Susceptible)" annotation(
    Placement(visible = true, transformation(origin = {100, 60}, extent = {{-10, -10}, {10, 10}}, rotation = 0), iconTransformation(origin = {110, 60}, extent = {{-10, -10}, {10, 10}}, rotation = 0)));
  Modelica.Blocks.Interfaces.RealOutput Iy "Output connector of I (Infectious)" annotation(
    Placement(visible = true, transformation(origin = {100, 0}, extent = {{-10, -10}, {10, 10}}, rotation = 0), iconTransformation(origin = {110, 0}, extent = {{-10, -10}, {10, 10}}, rotation = 0)));
  Modelica.Blocks.Interfaces.RealOutput Ry "Output connector of R (Removed)" annotation(
    Placement(visible = true, transformation(origin = {100, -60}, extent = {{-10, -10}, {10, 10}}, rotation = 0), iconTransformation(origin = {110, -60}, extent = {{-10, -10}, {10, 10}}, rotation = 0)));
initial equation
  Sy = Si;
  Iy = Ii;
  Ry = Ri;
equation
  der(Sy) = -contact_rate * Sy * Iy;
  der(Iy) = contact_rate * Sy * Iy - Iy * gamma_const;
  der(Ry) = Iy * gamma_const;
  annotation(
    Icon(graphics = {Text(origin = {-14, 11}, extent = {{-56, 33}, {78, -51}}, textString = "SIR")}));
end cl_SIRsimple;

Simulation à l'aide de la classe cl_SIRsimple

Créons un modèle qui donne la valeur initiale sous forme de constante et déplacez-la. La valeur initiale est S = 999, I = 1, R = 0, en supposant que la population totale est de 1000. Les paramètres β = 0,5 / 1000 et γ = 0,2.

Le graphique suivant est obtenu, ce qui équivaut au modèle SIR calculé dans d'autres langues. Considérez les heures de l'unité comme remplaçant le jour.

Créez une classe cl_SIRbetaInput qui donne le taux d'infection en externe

Afin de simuler le modèle Nishiura, nous ajouterons une interface d'entrée pour changer la valeur β, et ajouterons également une sortie "Nombre de personnes nouvellement infectées (incrémentation par jour)".

Commencez par mettre en commentaire le contact_rate qui était un paramètre afin qu'il soit reçu de Input. Inew est la sortie du nombre de personnes nouvellement infectées.

//  parameter Real contact_rate = 0.2/1000 "The infectious rate (1/day)";
parameter Real gamma_const = 0.2 "The Recovery rate (1/day)";

Modelica.Blocks.Interfaces.RealInput contact_rate "Connector of infectious rate" annotation(
  Placement(visible = true, transformation(origin = {0, 100}, extent = {{-20, -20}, {20, 20}}, rotation = -90), iconTransformation(origin = {0, 120}, extent = {{-20, -20}, {20, 20}}, rotation = -90)));
Modelica.Blocks.Interfaces.RealOutput Inew "Output connector of New I(Infectious)" annotation(
  Placement(visible = true, transformation(origin = {0, -100}, extent = {{-10, -10}, {10, 10}}, rotation = -90), iconTransformation(origin = {0, -110}, extent = {{-10, -10}, {10, 10}}, rotation = -90)));

Ajout de Inew à la partie équation. Cela signifie "le nombre de personnes nouvellement infectées = réduction du nombre de personnes pouvant être infectées".

equation
  der(Sy) = -contact_rate * Sy * Iy;
  der(Iy) = contact_rate * Sy * Iy - Iy * gamma_const;
  der(Ry) = Iy * gamma_const;
  Inew = -der(Sy);

Le code entier ressemble à ceci.

class cl_SIRbetaInput
  extends Modelica.Blocks.Icons.Block;
    //  parameter Real contact_rate = 0.2/1000 "The infectious rate (1/day)";
    parameter Real gamma_const = 0.2 "The Recovery rate (1/day)";
    Modelica.Blocks.Interfaces.RealInput Si "Connector of initial S(Susceptible)" annotation(
      Placement(visible = true, transformation(origin = {-100, 60}, extent = {{-20, -20}, {20, 20}}, rotation = 0), iconTransformation(origin = {-120, 60}, extent = {{-20, -20}, {20, 20}}, rotation = 0)));
    Modelica.Blocks.Interfaces.RealInput Ii "Connector of initial I(Infectious)" annotation(
      Placement(visible = true, transformation(origin = {-100, 0}, extent = {{-20, -20}, {20, 20}}, rotation = 0), iconTransformation(origin = {-120, 0}, extent = {{-20, -20}, {20, 20}}, rotation = 0)));
    Modelica.Blocks.Interfaces.RealInput Ri "Connector of initial R(Removed)" annotation(
      Placement(visible = true, transformation(origin = {-100, -62}, extent = {{-20, -20}, {20, 20}}, rotation = 0), iconTransformation(origin = {-120, -60}, extent = {{-20, -20}, {20, 20}}, rotation = 0)));
    Modelica.Blocks.Interfaces.RealOutput Sy "Output connector of S(Susceptible)" annotation(
      Placement(visible = true, transformation(origin = {100, 60}, extent = {{-10, -10}, {10, 10}}, rotation = 0), iconTransformation(origin = {110, 60}, extent = {{-10, -10}, {10, 10}}, rotation = 0)));
    Modelica.Blocks.Interfaces.RealOutput Iy "Output connector of I(Infectious)" annotation(
      Placement(visible = true, transformation(origin = {100, 0}, extent = {{-10, -10}, {10, 10}}, rotation = 0), iconTransformation(origin = {110, 0}, extent = {{-10, -10}, {10, 10}}, rotation = 0)));
    Modelica.Blocks.Interfaces.RealOutput Ry "Output connector of R(Removed)" annotation(
      Placement(visible = true, transformation(origin = {100, -60}, extent = {{-10, -10}, {10, 10}}, rotation = 0), iconTransformation(origin = {110, -60}, extent = {{-10, -10}, {10, 10}}, rotation = 0)));
    Modelica.Blocks.Interfaces.RealInput contact_rate "Connector of infectious rate" annotation(
      Placement(visible = true, transformation(origin = {0, 100}, extent = {{-20, -20}, {20, 20}}, rotation = -90), iconTransformation(origin = {0, 120}, extent = {{-20, -20}, {20, 20}}, rotation = -90)));
    Modelica.Blocks.Interfaces.RealOutput Inew "Output connector of New I(Infectious)" annotation(
      Placement(visible = true, transformation(origin = {0, -100}, extent = {{-10, -10}, {10, 10}}, rotation = -90), iconTransformation(origin = {0, -110}, extent = {{-10, -10}, {10, 10}}, rotation = -90)));
  initial equation
    Sy = Si;
    Iy = Ii;
    Ry = Ri;
  equation
    der(Sy) = -contact_rate * Sy * Iy;
    der(Iy) = contact_rate * Sy * Iy - Iy * gamma_const;
    der(Ry) = Iy * gamma_const;
    Inew = -der(Sy);
    annotation(
      Icon(graphics = {Text(origin = {-14, 11}, extent = {{-56, 33}, {78, -51}}, textString = "SIR")}, coordinateSystem(initialScale = 0.1)));
end cl_SIRbetaInput;

Résultats de simulation avec des taux de contact variables

En utilisant cela, simulons en modifiant le taux de contact le 12e jour. Le triangle descendant de l'icône "SIR" est le connecteur de sortie pour le nombre de personnes nouvellement infectées.

Je ne pense pas que les paramètres aient été suffisamment ajustés, mais je pense que nous avons pu reproduire le graphique de l'évolution du nombre de personnes nouvellement infectées dans le «modèle Nishiura».

Résumé

OpenModelica a pu simuler un modèle SIR d'une maladie infectieuse. J'ai aussi fait un modèle qui montre l'effet du vaccin, un modèle qui donne un taux de mortalité au nombre de guérisseurs, et un modèle qui considère la réinfection, je voudrais donc l'introduire dans un autre article.