Mortality and reinfection in class cl_SIRbetaInput created in the previous article, "Implementing a mathematical model of infectious diseases" SIR model "in OpenModelica" Let's create a full version with the rate added.
Prerequisites for this model ――A certain percentage of people die from infection to recovery ――Even those who have acquired immunity once lose immunity at a certain rate and return to infectious people over time. --Only survivors return from recoverers to infectious individuals
** Please note that all of the following simulations are only verification of the model with temporary parameters and do not reflect or anticipate the actual values. ** **
It is possible to set "dying_rate = mortality rate" as a parameter. The default value is 0.01 (1%).
parameter Real dying_rate = 0.01 "Dying rate";
The additional equation parts are as follows When transitioning from I (infected person) to R (recoverer), the number of recoverers becomes death in proportion to the mortality rate, and the rest become "survival recovery people".
equation
der(Rdead) = Iy * gamma_const * dying_rate;
// der(Rdead) = der(Ry) * dying_rate;
Rlive = Ry - Rdead;
For the formula considering reinfection, refer to "R and Pandemic Mathematical Model New Coronavirus (COVID-19) Research as an Example". I did. If ρ is the reinfection rate, the formula is as follows.
\begin{align}
\frac{dS}{dt} &= -\beta SI + ρ R\\
\frac{dI}{dt} &= \beta SI -\gamma I \\
\frac{dR}{dt} &= \gamma I - ρ R\\
\end{align}
Allows you to set "reinfect_rate = reinfection rate" as a parameter. The default value is 0.001 (0.1%).
parameter Real reinfect_rate = 0.001 "Reinfection rate";
The additional equation parts are as follows When transitioning from I (infected person) to R (recoverer), the number of recoverers becomes death in proportion to the mortality rate, and the rest become "survival recovery people". Among "survival recoverers", S (potentially infectious) shifts in proportion to the reinfection rate, and R (recovery) decreases accordingly.
equation
der(Sy) = (-contact_rate * Sy * Iy) + Rlive * reinfect_rate;
der(Iy) = contact_rate * Sy * Iy - Iy * gamma_const;
der(Ry) = Iy * gamma_const - Rlive * reinfect_rate;
Inew = contact_rate * Sy * Iy;
The whole code looks like this.
class cl_SIRreinfect
extends Modelica.Blocks.Icons.Block;
parameter Real gamma_const = 0.2 "Recovery rate (1/day)";
parameter Real reinfect_rate = 0.001 "Reinfection rate";
parameter Real dying_rate = 0.01 "Dying rate";
Modelica.Blocks.Interfaces.RealInput Si "Connector of initial S(Susceptible)" annotation(
Placement(visible = true, transformation(origin = {-120, 60}, extent = {{-20, -20}, {20, 20}}, rotation = 0), iconTransformation(origin = {-120, 62}, extent = {{-20, -20}, {20, 20}}, rotation = 0)));
Modelica.Blocks.Interfaces.RealInput Ii "Connector of initial I (Infectious)" annotation(
Placement(visible = true, transformation(origin = {-120, 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 = {-120, -60}, extent = {{-20, -20}, {20, 20}}, rotation = 0), iconTransformation(origin = {-120, -62}, extent = {{-20, -20}, {20, 20}}, rotation = 0)));
Modelica.Blocks.Interfaces.RealOutput Sy "Output connector of S(Susceptible)" annotation(
Placement(visible = true, transformation(origin = {110, 60}, extent = {{-10, -10}, {10, 10}}, rotation = 0), iconTransformation(origin = {110, 62}, extent = {{-10, -10}, {10, 10}}, rotation = 0)));
Modelica.Blocks.Interfaces.RealOutput Iy "Output connector of I (Infectious)" annotation(
Placement(visible = true, transformation(origin = {110, 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 = {110, -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, 120}, 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 = {-60, -110}, extent = {{-10, -10}, {10, 10}}, rotation = -90), iconTransformation(origin = {-60, -110}, extent = {{-10, -10}, {10, 10}}, rotation = -90)));
Modelica.Blocks.Interfaces.RealOutput Rlive "Output connector of R (Live)" annotation(
Placement(visible = true, transformation(origin = {60, -110}, extent = {{-10, -10}, {10, 10}}, rotation = -90), iconTransformation(origin = {60, -110}, extent = {{-10, -10}, {10, 10}}, rotation = -90)));
Modelica.Blocks.Interfaces.RealOutput Rdead "Output connector of R (Dead)" annotation(
Placement(visible = true, transformation(origin = {0, -110}, 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;
Rdead = 0;
equation
der(Sy) = (-contact_rate * Sy * Iy) + Rlive * reinfect_rate;
der(Iy) = contact_rate * Sy * Iy - Iy * gamma_const;
der(Ry) = Iy * gamma_const - Rlive * reinfect_rate;
Inew = contact_rate * Sy * Iy;
der(Rdead) = Iy * gamma_const * dying_rate;
// der(Rdead) = der(Ry) * dying_rate;
Rlive = Ry - Rdead;
annotation(
Diagram,
Icon(coordinateSystem(initialScale = 0.1), graphics = {Text(origin = {-62, 53}, extent = {{-52, 29}, {174, -133}}, textString = "SIR\nReI")}));
end cl_SIRreinfect;
Use the created cl_SIRreinfect to create and simulate the model mdl_SIRmodel_reinfect. Vaccines are not used here.
For clarity, we will simulate with a fairly high mortality rate of 20%.
As a result, the number of deaths has converged to about 18% of the total.
Next, let's simulate with a model in which only the reinfection rate is set. Try a higher 1% reinfection rate.
In the simulation up to 100 days, the number of people who can be infected continues to increase and it is not possible to see what will happen in the future. When simulating up to 500 days, S, I, and R all oscillate and gradually converge. Looking at I (number of infected people), we can see that there are 2nd and 3rd waves.
Try setting the above two parameters at the same time.
By modeling with OpenModelica, I think it was easy to add elements to the basic SIR model little by little and simulate it. Actually, it seems that both the mortality rate and the reinfection rate are more than an order of magnitude lower than this simulation, but if more accurate actual data can be obtained, it may be possible to identify the parameters.
Implementing the mathematical model "SIR model" of infectious diseases with OpenModelica Implementing a mathematical model "SIR model" of infectious diseases with OpenModelica (see the effect of the vaccine) Introduction of mathematical prediction model for infectious diseases (SIR model)
Recommended Posts