IDE-based SEIR-type model

The IDE-based SEIR type module models and simulates an epidemic using integro-differential equations. The model is particularly suited for pathogens with a latent infection state. The model assumes perfect immunity after recovery and is thus only suited for epidemic use cases. In the following, we present the model in detail. Note that this model is solved differently tha the IDE-SECIR model and the structure presented in IDE models applies only partially.

The simulation runs in discrete time steps using a trapezoidal rule. The model and the numerical scheme is based on the paper Keimer A, Pflug L. (2020). Modeling infectious diseases using integro-differential equations: Optimal control strategies for policy decisions and Applications in COVID-19. Preprint. http://dx.doi.org/10.13140/RG.2.2.10845.44000.

For a detailed description and application of the model, see: Plötzke L. (2020). Modellierung epidemischer Infektionskrankheiten auf der Basis von gewöhnlichen und Integro-Differentialgleichungen. Bachelor’s thesis, University of Cologne. https://elib.dlr.de/143504/.

Infection States

The model contains the following list of InfectionStates:

`Susceptible`
`Exposed`
`Infected`
`Recovered`

Parameters

The model implements the following parameters.

Mathematical variable	C++ variable name	Description
\(\phi\)	`ContactFrequency`	Average number of contacts of an individual per day.
\(\rho\)	`TransmissionRisk`	Transmission risk for people located in the susceptible compartments.
\(T_{E}\)	`LatencyTime`	Time in days an individual stays in the Exposed compartment.
\(T_{I}\)	`InfectiousTime`	Time in days an individual stays in the Infected compartment.

Initial conditions

To initialize the model, the following inputs need to be passed to the model constructor:

a TimeSeries containing the number of `Susceptible`s for a large enough number of time points before the start of the simulation,
the time step size \(dt\) used for numerical integration,
the size of the population of the considered region \(N\).

The initialization of the model can be done as follows where we set the Susceptible`s from :math:-15, dots, 0` based on the total population and the time of the previous time point.

using Vec = mio::TimeSeries<double>::Vector;


int N     = 810000;
double dt = 0.1;
mio::TimeSeries<double> init(1);

init.add_time_point<Eigen::VectorXd>(-15.0, Vec::Constant(1, N * 0.95));
while (init.get_last_time() < 0) {
    init.add_time_point(init.get_last_time() + dt,
                        Vec::Constant(1, (double)init.get_last_value()[0] + init.get_last_time()));
}

// Initialize model.
mio::iseir::Model<double> model(std::move(init), dt, N);

Nonpharmaceutical Interventions

In the SEIR model, nonpharmaceutical interventions (NPIs) are implemented through dampings in the contact matrix. These dampings reduce the contact rates between different groups to simulate interventions.

Basic dampings can be added to the contact matrix as follows:

// Create a contact matrix with constant contact rates between all groups.
ScalarType cont_freq = 10.;
mio::ContactMatrixGroup& contact_matrix = model.parameters.get<mio::osecir::ContactPatterns<ScalarType>>();
contact_matrix[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(1, 1, cont_freq));

// Add a damping that reduces contacts by 70% starting at day 30.
contact_matrix[0].add_damping(0.7, mio::SimulationTime(30.));

Simulation

After defining \(t_{\max}\), we can simulate, which means that we calculate the value for the compartment \(S\).

int tmax  = 15;
model.simulate(tmax);

The values of the remaining compartments \(E\), \(I\) and \(R\) are calculated using the parameters LatencyTime and InfectiousTime and obtain a time series containing the values of all compartments.

auto result = model.calculate_EIR();

Output

The output of the simulationis a TimeSeries object containing the size of the compartments at all time points. You can access the results as follows:

// Get the number of time points.
auto num_points = static_cast<size_t>(result.get_num_time_points());

// Access data at a specific time point.
Eigen::VectorX value_at_time_i = result.get_value(i);
ScalarType time_i = result.get_time(i);

// Access the last time point.
Eigen::VectorX last_value = result.get_last_value();
ScalarType last_time = result.get_last_time();

The order of the compartments follows the definition in the InfectionState enum.

You can print the simulation results as a formatted table:

// Print results to console with default formatting.
result.print_table();

// Print with custom column labels.
std::vector<std::string> labels = {"S", "E", "I", "R"};
result.print_table(labels);

Additionally, you can export the results to a CSV file:

// Export results to CSV with default settings.
result.export_csv("simulation_results.csv");

Visualization

To visualize the results of a simulation, you can use the Python package m-plot and its documentation.

You can export your simulation results to CSV format as described above.

Examples

An example can be found at:

examples/ide_seir.cpp

The code documentation for the model can be found at mio::iseir .