diff --git a/cpp/CMakeLists.txt b/cpp/CMakeLists.txt index 3ce345b3fd..6af0d5d74e 100644 --- a/cpp/CMakeLists.txt +++ b/cpp/CMakeLists.txt @@ -183,6 +183,7 @@ if(MEMILIO_BUILD_MODELS) add_subdirectory(models/ide_secir) add_subdirectory(models/ide_seir) add_subdirectory(models/ode_seir) + add_subdirectory(models/ode_seirv) add_subdirectory(models/ode_seair) add_subdirectory(models/ode_sir) add_subdirectory(models/sde_sir) diff --git a/cpp/examples/CMakeLists.txt b/cpp/examples/CMakeLists.txt index 39a12539e1..7b87f40bd2 100644 --- a/cpp/examples/CMakeLists.txt +++ b/cpp/examples/CMakeLists.txt @@ -68,6 +68,10 @@ add_executable(seir_flows_example ode_seir_flows.cpp) target_link_libraries(seir_flows_example PRIVATE memilio ode_seir) target_compile_options(seir_flows_example PRIVATE ${MEMILIO_CXX_FLAGS_ENABLE_WARNING_ERRORS}) +add_executable(ode_seirv_example ode_seirv_example.cpp) +target_link_libraries(ode_seirv_example PRIVATE memilio ode_seirv) +target_compile_options(ode_seirv_example PRIVATE ${MEMILIO_CXX_FLAGS_ENABLE_WARNING_ERRORS}) + add_executable(ode_secir_example ode_secir.cpp) target_link_libraries(ode_secir_example PRIVATE memilio ode_secir) target_compile_options(ode_secir_example PRIVATE ${MEMILIO_CXX_FLAGS_ENABLE_WARNING_ERRORS}) diff --git a/cpp/examples/ode_seirv_example.cpp b/cpp/examples/ode_seirv_example.cpp new file mode 100644 index 0000000000..46b5928369 --- /dev/null +++ b/cpp/examples/ode_seirv_example.cpp @@ -0,0 +1,120 @@ +/* +* Copyright (C) 2020-2025 MEmilio +* +* Authors: Henrik Zunker +* +* Contact: Martin J. Kuehn +* +* Licensed under the Apache License, Version 2.0 (the "License"); +* you may not use this file except in compliance with the License. +* You may obtain a copy of the License at +* +* http://www.apache.org/licenses/LICENSE-2.0 +* +* Unless required by applicable law or agreed to in writing, software +* distributed under the License is distributed on an "AS IS" BASIS, +* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +* See the License for the specific language governing permissions and +* limitations under the License. +*/ + +#include "ode_seirv/model.h" +#include "memilio/compartments/simulation.h" +#include "memilio/utils/logging.h" + +/** + * @brief set_initial_population sets the initial population of the model + * + * We assume that no one is initially infected or exposed at the beginning of the season. + * Infections are seeded later via the OutsideFoI parameter. + * The destinction into the 2 layers of susceptibility is done via the parameters. + * @tparam FP floating point type + * @param[in, out] model SEIRV model + * @param[in] total_pop total population size + */ +template +void set_initial_population(mio::oseirv::Model& model, const FP total_pop) +{ + auto& params = model.parameters; + auto& pop = model.populations; + const size_t num_groups = (size_t)params.get_num_groups(); + + for (size_t i = 0; i < num_groups; ++i) { + pop.template set_difference_from_group_total( + {mio::AgeGroup(i), mio::oseirv::InfectionState::Susceptible}, total_pop / num_groups); + + // Total population N_i as currently stored in group totals + FP Ni = pop.get_group_total(mio::AgeGroup(i)); + + FP s_age = params.template get>()[mio::AgeGroup(i)]; + FP s_frac = params.template get>(); + FP vc = params.template get>()[mio::AgeGroup(i)]; + FP ve = params.template get>()[mio::AgeGroup(i)]; + + pop[{mio::AgeGroup(i), mio::oseirv::InfectionState::Exposed}] = 0; + pop[{mio::AgeGroup(i), mio::oseirv::InfectionState::ExposedVaccinated}] = 0; + pop[{mio::AgeGroup(i), mio::oseirv::InfectionState::Infected}] = 0; + pop[{mio::AgeGroup(i), mio::oseirv::InfectionState::InfectedVaccinated}] = 0; + + pop[{mio::AgeGroup(i), mio::oseirv::InfectionState::Susceptible}] = s_frac * s_age * (FP(1) - vc) * Ni; + pop[{mio::AgeGroup(i), mio::oseirv::InfectionState::SusceptibleVaccinated}] = + s_frac * s_age * (FP(1) - ve) * vc * Ni; + + pop[{mio::AgeGroup(i), mio::oseirv::InfectionState::Recovered}] = (FP(1) - s_frac * s_age) * (FP(1) - vc) * Ni; + pop[{mio::AgeGroup(i), mio::oseirv::InfectionState::RecoveredVaccinated}] = + (FP(1) - s_frac * s_age * (FP(1) - ve)) * vc * Ni; + } +} + +// Example usage of the SEIRV ODE model presented in https://doi.org/10.1186/s12879-017-2344-6 +int main() +{ + using FP = double; + mio::set_log_level(mio::LogLevel::debug); + + const FP t0 = 0., tmax = 42.; + const FP dt = 0.1; + + const size_t num_groups = 1; + const auto total_pop = 1e5; + mio::oseirv::Model model((int)num_groups); + auto& parameters = model.parameters; + + FP cont_freq = 10.0; + FP cont_freq_group = FP(1) / FP(num_groups); + + // Healthy contacts + mio::ContactMatrixGroup& contacts_healthy = parameters.get>(); + contacts_healthy[0] = mio::ContactMatrix( + Eigen::MatrixXd::Constant((Eigen::Index)num_groups, (Eigen::Index)num_groups, cont_freq_group * cont_freq)); + + // Sick contacts (here 20% reduced) + mio::ContactMatrixGroup& contacts_sick = parameters.get>(); + contacts_sick[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant( + (Eigen::Index)num_groups, (Eigen::Index)num_groups, cont_freq_group * cont_freq * 0.8)); + + // Parameters + parameters.get>() = 1.2; + parameters.get>() = 1.0 / 2.0; + parameters.get>() = 1.0; + parameters.get>() = 0.0; + parameters.get>() = 0.0; + parameters.get>() = 1e-6; + parameters.get>() = 0.9; + parameters.get>() = 2.0; + + for (size_t i = 0; i < num_groups; ++i) { + parameters.get>()[mio::AgeGroup(i)] = 1.0; + parameters.get>()[mio::AgeGroup(i)] = 0.3; + parameters.get>()[mio::AgeGroup(i)] = 0.5; + model.populations.set_difference_from_group_total( + {mio::AgeGroup(i), mio::oseirv::InfectionState::Susceptible}, 1e5 / num_groups); + } + + set_initial_population(model, total_pop); + + auto seirv = mio::simulate(t0, tmax, dt, model); + + std::vector vars = {"S", "SV", "E", "EV", "I", "IV", "R", "RV"}; + seirv.print_table(vars); +} diff --git a/cpp/models/ode_seirv/CMakeLists.txt b/cpp/models/ode_seirv/CMakeLists.txt new file mode 100644 index 0000000000..b6068e1b9d --- /dev/null +++ b/cpp/models/ode_seirv/CMakeLists.txt @@ -0,0 +1,12 @@ +add_library(ode_seirv + infection_state.h + model.h + model.cpp + parameters.h +) +target_link_libraries(ode_seirv PUBLIC memilio) +target_include_directories(ode_seirv PUBLIC + $ + $ +) +target_compile_options(ode_seirv PRIVATE ${MEMILIO_CXX_FLAGS_ENABLE_WARNING_ERRORS}) diff --git a/cpp/models/ode_seirv/README.md b/cpp/models/ode_seirv/README.md new file mode 100644 index 0000000000..0d2e2fc65e --- /dev/null +++ b/cpp/models/ode_seirv/README.md @@ -0,0 +1,23 @@ +# ODE SEIRV (seasonal influenza) + +A age‑structured **SEIRV** model implemented into the MEmilio ODE framework. +It extends the classic SEIR by a vaccinated layer (Sᵛ, Eᵛ, Iᵛ, Rᵛ) and follows the seasonal influenza +model structure from: + +> **Weidemann, F., Remschmidt, C., Buda, S. et al.** +> *Is the impact of childhood influenza vaccination less than expected: a transmission modelling study.* +> **BMC Infectious Diseases** 17, 258 (2017). +> https://doi.org/10.1186/s12879-017-2344-6 + +## Get started + +To get started, check out the +[official documentation](https://memilio.readthedocs.io/en/latest/cpp/models/oseirv.html) +or the [code example](../../examples/ode_seirv_example.cpp). + +## Reference + +Please cite the original study if you use this model: + +> Weidemann, F., Remschmidt, C., Buda, S. et al. *Is the impact of childhood influenza vaccination less than expected: a transmission modelling study.* **BMC Infect Dis** 17, 258 (2017). +> https://doi.org/10.1186/s12879-017-2344-6 diff --git a/cpp/models/ode_seirv/infection_state.h b/cpp/models/ode_seirv/infection_state.h new file mode 100644 index 0000000000..7eb0078eb9 --- /dev/null +++ b/cpp/models/ode_seirv/infection_state.h @@ -0,0 +1,42 @@ +/* +* Copyright (C) 2020-2025 MEmilio +* +* Authors: Henrik Zunker +* +* Contact: Martin J. Kuehn +* +* Licensed under the Apache License, Version 2.0 (the "License"); +* you may not use this file except in compliance with the License. +* You may obtain a copy of the License at +* +* http://www.apache.org/licenses/LICENSE-2.0 +* +* Unless required by applicable law or agreed to in writing, software +* distributed under the License is distributed on an "AS IS" BASIS, +* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +* See the License for the specific language governing permissions and +* limitations under the License. +*/ +#ifndef SEIRV_INFECTIONSTATE_H +#define SEIRV_INFECTIONSTATE_H +namespace mio +{ +namespace oseirv +{ + +enum class InfectionState +{ + Susceptible, + SusceptibleVaccinated, + Exposed, + ExposedVaccinated, + Infected, + InfectedVaccinated, + Recovered, + RecoveredVaccinated, + Count +}; + +} // namespace oseirv +} // namespace mio +#endif diff --git a/cpp/models/ode_seirv/model.cpp b/cpp/models/ode_seirv/model.cpp new file mode 100644 index 0000000000..5c62bf43cf --- /dev/null +++ b/cpp/models/ode_seirv/model.cpp @@ -0,0 +1,20 @@ +/* +* Copyright (C) 2020-2025 MEmilio +* +* Authors: Henrik Zunker +* +* Contact: Martin J. Kuehn +* +* Licensed under the Apache License, Version 2.0 (the "License"); +* you may not use this file except in compliance with the License. +* You may obtain a copy of the License at +* +* http://www.apache.org/licenses/LICENSE-2.0 +* +* Unless required by applicable law or agreed to in writing, software +* distributed under the License is distributed on an "AS IS" BASIS, +* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +* See the License for the specific language governing permissions and +* limitations under the License. +*/ +#include "ode_seirv/model.h" diff --git a/cpp/models/ode_seirv/model.h b/cpp/models/ode_seirv/model.h new file mode 100644 index 0000000000..54a0addc7c --- /dev/null +++ b/cpp/models/ode_seirv/model.h @@ -0,0 +1,172 @@ +/* +* Copyright (C) 2020-2025 MEmilio +* +* Authors: Henrik Zunker +* +* Contact: Martin J. Kuehn +* +* Licensed under the Apache License, Version 2.0 (the "License"); +* you may not use this file except in compliance with the License. +* You may obtain a copy of the License at +* +* http://www.apache.org/licenses/LICENSE-2.0 +* +* Unless required by applicable law or agreed to in writing, software +* distributed under the License is distributed on an "AS IS" BASIS, +* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +* See the License for the specific language governing permissions and +* limitations under the License. +*/ +#ifndef SEIRV_MODEL_H +#define SEIRV_MODEL_H + +#include "memilio/compartments/flow_model.h" +#include "memilio/config.h" +#include "memilio/epidemiology/age_group.h" +#include "memilio/epidemiology/populations.h" +#include "memilio/utils/time_series.h" +#include "ode_seirv/infection_state.h" +#include "ode_seirv/parameters.h" + +GCC_CLANG_DIAGNOSTIC(push) +GCC_CLANG_DIAGNOSTIC(ignored "-Wshadow") +#include +#include +GCC_CLANG_DIAGNOSTIC(pop) + +#include +#include + +namespace mio +{ +namespace oseirv +{ + +// clang-format off +using Flows = TypeList, + Flow, + Flow, + Flow, + Flow, + Flow>; +// clang-format on + +template +class Model + : public FlowModel, Parameters, Flows> +{ + using Base = FlowModel, Parameters, Flows>; + +public: + using ParameterSet = typename Base::ParameterSet; + using PopulationsType = typename Base::Populations; + + Model(const PopulationsType& pop, const ParameterSet& params) + : Base(pop, params) + { + } + Model(int num_agegroups) + : Base(PopulationsType({AgeGroup(num_agegroups), InfectionState::Count}), ParameterSet(AgeGroup(num_agegroups))) + { + } + + void get_flows(Eigen::Ref> pop, Eigen::Ref> y, FP t, + Eigen::Ref> flows) const override + { + const auto& params = this->parameters; + const size_t num_groups = (size_t)params.get_num_groups(); + + // Get effective contact matrix effective_contacts = contacts_healthy + m*contacts_sick + Eigen::MatrixX contacts_healthy = + params.template get>().get_cont_freq_mat().get_matrix_at(SimulationTime(t)); + Eigen::MatrixX contacts_sick = + params.template get>().get_cont_freq_mat().get_matrix_at(SimulationTime(t)); + const FP m = params.template get>(); + Eigen::MatrixX effective_contacts = contacts_healthy + m * contacts_sick; + + // Normalization ν(m): prepare susceptibility scaling and next-generation matrix (NG) + Eigen::VectorX susceptibility(num_groups); + for (size_t i = 0; i < num_groups; ++i) + susceptibility[(Eigen::Index)i] = params.template get>()[AgeGroup(i)] * + params.template get>(); + // Sigma is a diagonal matrix applying susceptibility on the receiving (row) side of contacts. + Eigen::MatrixX Sigma = susceptibility.asDiagonal(); + + // NG = Sigma * effective_contacts is the (unnormalized) next-generation matrix. Its spectral radius will be used + // below to scale contacts effective_contacts so that structure (mixing pattern) and transmissibility magnitude (BaselineTransmissibility) + // are cleanly separated. + Eigen::MatrixX NG = Sigma * effective_contacts; + + Eigen::ComplexEigenSolver> ces; + ces.compute(NG); + FP scale_transmissibility = ces.eigenvalues().cwiseAbs().maxCoeff(); + if (scale_transmissibility > FP(0.)) { + effective_contacts /= scale_transmissibility; + } + + const FP transmissibility_baseline = + params.template get>(); // baseline transmissibility scaling (R-like factor) + const FP recovery_rate = + params.template get>(); // progression/recovery rate (1 / mean duration) + const FP season_amplitude = params.template get>(); // amplitude of seasonal modulation + const FP season_shift_per_subtype = + params.template get>(); // base phase shift for seasonality + const FP season_shift_per_season = + params.template get>(); // additional phase shift (e.g. season-specific) + const FP outside_foi = params.template get>(); // constant external force of infection + const FP clustering_exp = params.template get>(); // clustering exponent; clustering_exp<1 dampens, clustering_exp>1 amplifies prevalence effect + + const FP season = + std::exp(season_amplitude * std::sin(FP(2) * std::numbers::pi_v * + (t / FP(52.0) - season_shift_per_subtype + season_shift_per_season))); + + for (auto i : make_index_range((mio::AgeGroup)num_groups)) { + // Flat indices + const size_t Si = this->populations.get_flat_index({i, InfectionState::Susceptible}); + const size_t SVi = this->populations.get_flat_index({i, InfectionState::SusceptibleVaccinated}); + const size_t Ei = this->populations.get_flat_index({i, InfectionState::Exposed}); + const size_t EVi = this->populations.get_flat_index({i, InfectionState::ExposedVaccinated}); + const size_t Ii = this->populations.get_flat_index({i, InfectionState::Infected}); + const size_t IVi = this->populations.get_flat_index({i, InfectionState::InfectedVaccinated}); + + FP sum = FP(0); + for (auto j : make_index_range((mio::AgeGroup)num_groups)) { + const size_t Sj = this->populations.get_flat_index({j, InfectionState::Susceptible}); + const size_t SVj = this->populations.get_flat_index({j, InfectionState::SusceptibleVaccinated}); + const size_t Ej = this->populations.get_flat_index({j, InfectionState::Exposed}); + const size_t EVj = this->populations.get_flat_index({j, InfectionState::ExposedVaccinated}); + const size_t Ij = this->populations.get_flat_index({j, InfectionState::Infected}); + const size_t IVj = this->populations.get_flat_index({j, InfectionState::InfectedVaccinated}); + const size_t Rj = this->populations.get_flat_index({j, InfectionState::Recovered}); + const size_t RVj = this->populations.get_flat_index({j, InfectionState::RecoveredVaccinated}); + + const FP Nj = pop[Sj] + pop[SVj] + pop[Ej] + pop[EVj] + pop[Ij] + pop[IVj] + pop[Rj] + pop[RVj]; + const FP inf_frac = (Nj > FP(0)) ? ((pop[Ij] + pop[IVj]) / Nj) : FP(0); + const FP transmission_rate = (clustering_exp == FP(1)) ? inf_frac : std::pow(inf_frac, clustering_exp); + sum += effective_contacts(i.get(), j.get()) * transmission_rate; + } + const FP foi_i = transmissibility_baseline * recovery_rate * season * sum + outside_foi; + + // Flows: S->E, S^V->E^V, E->I, E^V->I^V, I->R, I^V->R^V + flows[Base::template get_flat_flow_index(i)] = + foi_i * y[Si]; + flows[Base::template get_flat_flow_index(i)] = foi_i * y[SVi]; + + flows[Base::template get_flat_flow_index(i)] = + recovery_rate * y[Ei]; + flows[Base::template get_flat_flow_index(i)] = recovery_rate * y[EVi]; + + flows[Base::template get_flat_flow_index(i)] = + recovery_rate * y[Ii]; + flows[Base::template get_flat_flow_index(i)] = recovery_rate * y[IVi]; + } + } +}; + +} // namespace oseirv +} // namespace mio +#endif diff --git a/cpp/models/ode_seirv/parameters.h b/cpp/models/ode_seirv/parameters.h new file mode 100644 index 0000000000..e8149230fc --- /dev/null +++ b/cpp/models/ode_seirv/parameters.h @@ -0,0 +1,357 @@ +/* +* Copyright (C) 2020-2025 MEmilio +* +* Authors: Henrik Zunker +* +* Contact: Martin J. Kuehn +* +* Licensed under the Apache License, Version 2.0 (the "License"); +* you may not use this file except in compliance with the License. +* You may obtain a copy of the License at +* +* http://www.apache.org/licenses/LICENSE-2.0 +* +* Unless required by applicable law or agreed to in writing, software +* distributed under the License is distributed on an "AS IS" BASIS, +* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +* See the License for the specific language governing permissions and +* limitations under the License. +*/ +#ifndef SEIRV_PARAMETERS_H +#define SEIRV_PARAMETERS_H + +#include "memilio/epidemiology/age_group.h" +#include "memilio/epidemiology/uncertain_matrix.h" +#include "memilio/utils/custom_index_array.h" +#include "memilio/utils/uncertain_value.h" +#include "memilio/utils/parameter_set.h" + +namespace mio +{ +namespace oseirv +{ + +/** + * @brief Baseline transmissibility R_e (dimensionless). + * + * Multiplies the force of infection after contact-matrix normalization. Main scalar for transmission + * intensity. Default: 1.0. + */ +template +struct BaselineTransmissibility { + using Type = UncertainValue; + static Type get_default(AgeGroup) + { + return Type(1.0); + } + static std::string name() + { + return "BaselineTransmissibility"; + } +}; + +/** + * @brief Transition/recovery rate in 1/week. + * + * Controls E→I and I→R and is included in the force of infection to decouple the recovery rate from R_e. + * Use the model time unit (weeks). Default: 1/2.0 + */ +template +struct RecoveryRate { + using Type = UncertainValue; + static Type get_default(AgeGroup) + { + return Type(1.0 / 2.0); + } + static std::string name() + { + return "RecoveryRate"; + } +}; + +/** + * @brief Seasonality amplitude (dimensionless). + * + * Amplitude of the seasonal modulation. + * Set to 0.0 to disable seasonality. Default: 0.0. + */ +template +struct SeasonalityAmplitude { + using Type = UncertainValue; + static Type get_default(AgeGroup) + { + return Type(0.0); + } + static std::string name() + { + return "SeasonalityAmplitude"; + } +}; + +/** + * @brief Subtype-specific seasonal shift t_z. + * + * Phase shift for the sinusoid in the force of infection; may be adjusted by a seasonal-shift correction if implemented. Default: 0.0. + */ +template +struct SeasonalityShiftPerSubtype { + using Type = UncertainValue; + static Type get_default(AgeGroup) + { + return Type(0.0); + } + static std::string name() + { + return "SeasonalityShiftPerSubtype"; + } +}; + +/** + * @brief Season-specific fine shift t_s. + * + * Additional small seasonal timing correction per season. Default: 0.0. + */ +template +struct SeasonalityShiftPerSeason { + using Type = UncertainValue; + static Type get_default(AgeGroup) + { + return Type(0.0); + } + static std::string name() + { + return "SeasonalityShiftPerSeason"; + } +}; + +/** + * @brief Outside force of infection in 1/week. + * + * Additive external hazard that can seed infections when no infectives are present. + * Used in the beginning of the season (to import infections) or to model travel/imported cases. + * Default: 0.0. + */ +template +struct OutsideFoI { + using Type = UncertainValue; + static Type get_default(AgeGroup) + { + return Type(0.0); + } + static std::string name() + { + return "OutsideFoI"; + } +}; + +/** + * @brief Clustering/concavity parameter ρ (dimensionless). + * + * Exponent on the infectious fraction in force of infection. Must be > 0. Default: 1.0. + */ +template +struct ClusteringExponent { + using Type = UncertainValue; + static Type get_default(AgeGroup) + { + return Type(1.0); + } + static std::string name() + { + return "ClusteringExponent"; + } +}; + +/** + * @brief Mixing parameter m for “sick” contacts (dimensionless). + * + * Weight for the symptomatically sick contact matrix in β_eff = ν(m)^{-1}(β_healthy + m·β_sick). + * Combines share of symptomatic cases and their higher infectiousness. Default: 1.0. + */ +template +struct SickMixing { + using Type = UncertainValue; + static Type get_default(AgeGroup) + { + return Type(1.0); + } + static std::string name() + { + return "SickMixing"; + } +}; + +/** + * @brief Contact patterns of healthy people (age-structured contact frequencies). + * + * Used in β_eff, can be time-varying via damping. Unit: contacts per model time step (week). + */ +template +struct ContactPatternsHealthy { + using Type = UncertainContactMatrix; + static Type get_default(AgeGroup size) + { + return Type(1, static_cast((size_t)size)); + } + static std::string name() + { + return "ContactPatternsHealthy"; + } +}; + +/** + * @brief Contact patterns of symptomatically sick people (age-structured contact frequencies). + * + * Used in β_eff together with m, can be time-varying via damping. + * Unit: contacts per model time step (week). + */ +template +struct ContactPatternsSick { + using Type = UncertainContactMatrix; + static Type get_default(AgeGroup size) + { + return Type(1, static_cast((size_t)size)); + } + static std::string name() + { + return "ContactPatternsSick"; + } +}; + +/** + * @brief Age-specific baseline susceptibility σ_i. + * + * Represents pre-existing subtype/age-specific immunity. + * Default: 1.0 (fully susceptible). + */ +template +struct SusceptibilityByAge { + using Type = CustomIndexArray, AgeGroup>; + static Type get_default(AgeGroup size) + { + return Type(size, 1.0); + } + static std::string name() + { + return "SusceptibilityByAge"; + } +}; + +/** + * @brief Vaccination coverage VC_i at t0 (dimensionless, in [0,1]). + * + * Age-specific share vaccinated at season start; used in initial conditions. + * Default: 0.0. + */ +template +struct VaccineCoverage { + using Type = CustomIndexArray, AgeGroup>; + static Type get_default(AgeGroup size) + { + return Type(size, 0.0); + } + static std::string name() + { + return "VaccineCoverage"; + } +}; + +/** + * @brief Vaccine effectiveness VE_i (dimensionless, in [0,1]). + * + * Reduces effective susceptibility of vaccinated individuals in the initial state via (1 − VE_i). Default: 0.0. + */ +template +struct VaccineEffectiveness { + using Type = CustomIndexArray, AgeGroup>; + static Type get_default(AgeGroup size) + { + return Type(size, 0.0); + } + static std::string name() + { + return "VaccineEffectiveness"; + } +}; + +/** + * @brief Fraction of the population that remains susceptible at t0 phi (dimensionless, typically in [0,1]). + */ +template +struct SusceptibleFraction { + using Type = UncertainValue; + static Type get_default(AgeGroup) + { + return Type(1.0); + } + static std::string name() + { + return "SusceptibleFraction"; + } +}; + +template +using ParametersBase = + ParameterSet, RecoveryRate, SeasonalityAmplitude, + SeasonalityShiftPerSubtype, SeasonalityShiftPerSeason, OutsideFoI, ClusteringExponent, + SickMixing, ContactPatternsHealthy, ContactPatternsSick, SusceptibilityByAge, + VaccineCoverage, VaccineEffectiveness, SusceptibleFraction>; + +/** + * @brief Parameter set for the age-resolved SEIRV model (S,E,I,R plus vaccinated states) as per the appendix. + * + * Contains seasonality, contact-pattern and immunity/vaccination parameters. + * Model time unit is week; contact matrices may be time-dependent (damping). + */ +template +class Parameters : public ParametersBase +{ +public: + /** + * @brief Construct with the number of age groups. + * @param ng Number of age groups. + */ + Parameters(AgeGroup ng) + : ParametersBase(ng) + , m_num_groups{ng} + { + } + + /** + * @brief Returns the number of age groups. + */ + AgeGroup get_num_groups() const + { + return m_num_groups; + } + + /** + * @brief Apply minimal constraints to keep parameters valid. + * + * @return true if one or more parameters were corrected; false otherwise. + */ + bool apply_constraints() + { + bool corrected = false; + if (this->template get>() <= 0.0) { + this->template get>() = 1.0; + corrected = true; + } + if (this->template get>() < 0.0) { + this->template get>() = 0.0; + corrected = true; + } + if (this->template get>() < 0.0) { + this->template get>() = 0.0; + corrected = true; + } + return corrected; + } + +private: + AgeGroup m_num_groups; +}; + +} // namespace oseirv +} // namespace mio +#endif diff --git a/cpp/tests/CMakeLists.txt b/cpp/tests/CMakeLists.txt index 01db381a7e..f2933e0855 100644 --- a/cpp/tests/CMakeLists.txt +++ b/cpp/tests/CMakeLists.txt @@ -8,6 +8,7 @@ set(TESTSOURCES testmain.cpp test_populations.cpp test_odeseir.cpp + test_odeseirv.cpp test_odesir.cpp test_numericalIntegration.cpp test_smoother.cpp diff --git a/cpp/tests/test_odeseirv.cpp b/cpp/tests/test_odeseirv.cpp new file mode 100644 index 0000000000..5919da987b --- /dev/null +++ b/cpp/tests/test_odeseirv.cpp @@ -0,0 +1,286 @@ +/* +* Copyright (C) 2020-2025 MEmilio +* +* Authors: Henrik Zunker +* +* Contact: Martin J. Kuehn +* +* Licensed under the Apache License, Version 2.0 (the "License"); +* you may not use this file except in compliance with the License. +* You may obtain a copy of the License at +* +* http://www.apache.org/licenses/LICENSE-2.0 +* +* Unless required by applicable law or agreed to in writing, software +* distributed under the License is distributed on an "AS IS" BASIS, +* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +* See the License for the specific language governing permissions and +* limitations under the License. +*/ + +#include "memilio/config.h" +#include "memilio/math/integrator.h" +#include "memilio/math/euler.h" +#include "memilio/utils/time_series.h" +#include "memilio/compartments/simulation.h" +#include "memilio/compartments/flow_simulation.h" + +#include "ode_seirv/model.h" +#include "ode_seirv/infection_state.h" +#include "ode_seirv/parameters.h" + +#include + +#include +#include + +TEST(TestOdeSeirv, simulateDefault) +{ + double t0 = 0.0; + double tmax = 1.0; + double dt = 0.1; + + mio::oseirv::Model model(1); + auto result = simulate(t0, tmax, dt, model); // generic template simulate + + EXPECT_NEAR(result.get_last_time(), tmax, 1e-12); + EXPECT_EQ(result.get_num_elements(), (Eigen::Index)mio::oseirv::InfectionState::Count); +} + +TEST(TestOdeSeirv, populationConservation) +{ + mio::oseirv::Model model(1); + const double S = 300; + const double SV = 200; + const double E = 50; + const double EV = 25; + const double I = 40; + const double IV = 10; + const double R = 0; + const double RV = 0; + double total = S + SV + E + EV + I + IV + R + RV; + + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Susceptible}] = S; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::SusceptibleVaccinated}] = SV; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Exposed}] = E; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::ExposedVaccinated}] = EV; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Infected}] = I; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::InfectedVaccinated}] = IV; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Recovered}] = R; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::RecoveredVaccinated}] = RV; + + mio::ContactMatrixGroup& cm_h = model.parameters.get>(); + cm_h[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(1, 1, 1)); + mio::ContactMatrixGroup& cm_s = model.parameters.get>(); + cm_s[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(1, 1, 1)); + + model.parameters.set>(1.0); + model.parameters.set>(1.0); + + double t0 = 0.0, tmax = 5.0, dt = 0.5; + auto sim = simulate(t0, tmax, dt, model); + auto last = sim.get_last_value(); + EXPECT_NEAR(last.sum(), total, 1e-8); +} + +TEST(TestOdeSeirv, applyConstraints) +{ + mio::oseirv::Model model(1); + // First: defaults should need no correction + EXPECT_FALSE(model.parameters.apply_constraints()); + + model.parameters.set>(-0.5); // invalid, must become >0 (1.0) + model.parameters.set>(-2.); // invalid -> 0 + model.parameters.set>(-0.5); // invalid -> 0 + EXPECT_TRUE(model.parameters.apply_constraints()); + EXPECT_EQ(model.parameters.get>(), 1.0); + EXPECT_EQ(model.parameters.get>(), 0.0); + EXPECT_EQ(model.parameters.get>(), 0.0); +} + +TEST(TestOdeSeirv, flowsSingleAgeGroup) +{ + mio::oseirv::Model model(1); + + // Populations + const double S = 300, SV = 200, E = 50, EV = 25, I = 40, IV = 10; // R, RV = 0 + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Susceptible}] = S; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::SusceptibleVaccinated}] = SV; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Exposed}] = E; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::ExposedVaccinated}] = EV; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Infected}] = I; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::InfectedVaccinated}] = IV; + + // Parameters: identity contacts, RecoveryRate=1, ClusteringExponent=1, BaselineTransmissibility=1, others zero ⇒ λ = (I+IV)/N + mio::ContactMatrixGroup& cm_h = model.parameters.get>(); + cm_h[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(1, 1, 1)); + mio::ContactMatrixGroup& cm_s = model.parameters.get>(); + cm_s[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(1, 1, 1)); + model.parameters.set>(1.0); + model.parameters.set>(1.0); + model.parameters.set>(1.0); + model.parameters.set>(0.0); + model.parameters.set>(0.0); + + auto y0 = model.get_initial_values(); + auto pop = y0; // same vector for signature + Eigen::VectorXd flows(6); // number of flows for one age group + flows.setZero(); + model.get_flows(pop, y0, 0.0, flows); + + const auto recovery_rate = model.parameters.get>(); + const double N = S + SV + E + EV + I + IV; // (R,RV = 0) + const double lambda = (I + IV) / N; // expected force of infection + const double f_SE = S * lambda; + const double f_SV_EV = SV * lambda; + const double f_E_I = recovery_rate * E; + const double f_EV_IV = recovery_rate * EV; + const double f_I_R = recovery_rate * I; + const double f_IV_RV = recovery_rate * IV; + + auto idx_SE = model.template get_flat_flow_index(mio::AgeGroup(0)); + auto idx_SV = model.template get_flat_flow_index(mio::AgeGroup(0)); + auto idx_EI = + model.template get_flat_flow_index( + mio::AgeGroup(0)); + auto idx_EIV = + model.template get_flat_flow_index(mio::AgeGroup(0)); + auto idx_IR = model.template get_flat_flow_index(mio::AgeGroup(0)); + auto idx_IVR = + model.template get_flat_flow_index(mio::AgeGroup(0)); + + EXPECT_NEAR(flows[idx_SE], f_SE, 1e-12); + EXPECT_NEAR(flows[idx_SV], f_SV_EV, 1e-12); + EXPECT_NEAR(flows[idx_EI], f_E_I, 1e-12); + EXPECT_NEAR(flows[idx_EIV], f_EV_IV, 1e-12); + EXPECT_NEAR(flows[idx_IR], f_I_R, 1e-12); + EXPECT_NEAR(flows[idx_IVR], f_IV_RV, 1e-12); +} + +TEST(TestOdeSeirv, flowsTwoAgeGroupsIdentityContacts) +{ + mio::oseirv::Model model(2); + mio::ContactMatrixGroup& cm_h = model.parameters.get>(); + // Let each group have only contacts with itself and set other parameters, so that λ_i = I_i / N_i + Eigen::MatrixXd Id2 = Eigen::MatrixXd::Identity(2, 2); + cm_h[0] = mio::ContactMatrix(Id2); + mio::ContactMatrixGroup& cm_s = model.parameters.get>(); + cm_s[0] = mio::ContactMatrix(Eigen::MatrixXd::Zero(2, 2)); + model.parameters.set>(1.0); + model.parameters.set>(1.0); + model.parameters.set>(1.0); + model.parameters.set>(0.0); + model.parameters.set>(0.0); + + // Only population in the non-vaccinated susceptible and infected compartments + // Group 0 + double S0 = 100, I0 = 20; // others zero + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Susceptible}] = S0; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Infected}] = I0; + // Group 1 + double S1 = 80, I1 = 10; + model.populations[{mio::AgeGroup(1), mio::oseirv::InfectionState::Susceptible}] = S1; + model.populations[{mio::AgeGroup(1), mio::oseirv::InfectionState::Infected}] = I1; + + auto y0 = model.get_initial_values(); + auto pop = y0; + Eigen::VectorXd flows(12); // 6 flows * 2 age groups + flows.setZero(); + model.get_flows(pop, y0, 0.0, flows); + + double N0 = S0 + I0; + double N1 = S1 + I1; + double lambda0 = I0 / N0; // identity contacts => only own group contributes + double lambda1 = I1 / N1; + + auto idx_SE_0 = model.template get_flat_flow_index(mio::AgeGroup(0)); + auto idx_SE_1 = model.template get_flat_flow_index(mio::AgeGroup(1)); + auto idx_EI_0 = + model.template get_flat_flow_index( + mio::AgeGroup(0)); + auto idx_EI_1 = + model.template get_flat_flow_index( + mio::AgeGroup(1)); + auto idx_IR_0 = model.template get_flat_flow_index(mio::AgeGroup(0)); + auto idx_IR_1 = model.template get_flat_flow_index(mio::AgeGroup(1)); + + EXPECT_NEAR(flows[idx_SE_0], S0 * lambda0, 1e-12); + EXPECT_NEAR(flows[idx_SE_1], S1 * lambda1, 1e-12); + // Exposed are zero => progressions must be zero + EXPECT_NEAR(flows[idx_EI_0], 0.0, 1e-12); + EXPECT_NEAR(flows[idx_EI_1], 0.0, 1e-12); + EXPECT_NEAR(flows[idx_IR_0], I0, 1e-12); + EXPECT_NEAR(flows[idx_IR_1], I1, 1e-12); +} + +TEST(TestOdeSeirv, zeroPopulationNoNan) +{ + mio::oseirv::Model model(1); + model.populations.set_total(0.0); + auto y0 = model.get_initial_values(); + auto pop = y0; + Eigen::VectorXd flows(6); + flows.setZero(); + model.get_flows(pop, y0, 0.0, flows); + for (int i = 0; i < flows.size(); ++i) { + EXPECT_FALSE(std::isnan(flows[i])); + EXPECT_EQ(flows[i], 0.0); + } +} + +TEST(TestOdeSeirv, simulationEuler) +{ + mio::oseirv::Model model(1); + // Simple initial values + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Exposed}] = 10; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Infected}] = 5; + model.populations.set_difference_from_group_total( + {mio::AgeGroup(0), mio::oseirv::InfectionState::Susceptible}, 100.0); + + // Identity contacts + mio::ContactMatrixGroup& cm_h = model.parameters.get>(); + cm_h[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(1, 1, 1)); + mio::ContactMatrixGroup& cm_s = model.parameters.get>(); + cm_s[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(1, 1, 0.)); + model.parameters.set>(1.0); + model.parameters.set>(1.0); + + double t0 = 0.0, tmax = 1.0, dt = 0.5; + std::unique_ptr> integrator = std::make_unique>(); + auto sim = simulate(t0, tmax, dt, model, std::move(integrator)); + EXPECT_EQ(sim.get_num_time_points(), 3); // t=0,0.5,1.0 + // Sanity: all values non-negative + for (Eigen::Index i = 0; i < sim.get_last_value().size(); ++i) { + EXPECT_GE(sim.get_last_value()[i], 0.0); + } +} + +TEST(TestOdeSeirv, flowSimulationEuler) +{ + mio::oseirv::Model model(1); + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Infected}] = 5; + model.populations.set_difference_from_group_total( + {mio::AgeGroup(0), mio::oseirv::InfectionState::Susceptible}, 100.0); + + mio::ContactMatrixGroup& cm_h = model.parameters.get>(); + cm_h[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(1, 1, 1)); + model.parameters.set>(1.0); + model.parameters.set>(1.0); + + double t0 = 0.0, tmax = 0.5, dt = 0.5; + std::unique_ptr> integrator = std::make_unique>(); + auto sim = simulate_flows(t0, tmax, dt, model, std::move(integrator)); + EXPECT_EQ(sim[0].get_num_time_points(), 2); + EXPECT_EQ(sim[1].get_num_time_points(), 2); + // Flow vector size should be 6 for one age group + EXPECT_EQ(sim[1].get_last_value().size(), 6); +} diff --git a/docs/source/cpp/models/oseirv.rst b/docs/source/cpp/models/oseirv.rst new file mode 100644 index 0000000000..cd0591a013 --- /dev/null +++ b/docs/source/cpp/models/oseirv.rst @@ -0,0 +1,242 @@ +ODE-based SEIRV model +====================== + +The ODE-SEIRV module extends the classic SEIR model by an explicit vaccinated layer (``S^V, E^V, I^V, R^V``) and targets +seasonal influenza–type applications (cf. Weidemann et al. 2017). It is age-structured and uses a normalized contact +matrix formulation with explicit seasonality. Immunity after recovery is assumed to last for the considered season. Also, there is no transition between the unvaccinated and vaccinated compartments during the simulation horizon (i.e., no vaccinations after season start). + +The infection states and transitions are illustrated in the following figure. + +.. image:: https://martinkuehn.eu/research/images/seirv.png + :alt: SEIRV_model + + +Infection States +---------------- + +The model contains the following **InfectionState**\s: + +.. code-block:: RST + + `Susceptible` + `SusceptibleVaccinated` + `Exposed` + `ExposedVaccinated` + `Infected` + `InfectedVaccinated` + `Recovered` + `RecoveredVaccinated` + + +Infection State Transitions +--------------------------- + +The SEIRV model is implemented as a **FlowModel**. Thus, in each time step, the flows (new infections, progressions, +recoveries) are computed explicitly in addition to compartment values. The defined transitions `FromState, ToState` are: + +.. code-block:: RST + + `Susceptible, Exposed` + `SusceptibleVaccinated, ExposedVaccinated` + `Exposed, Infected` + `ExposedVaccinated, InfectedVaccinated` + `Infected, Recovered` + `InfectedVaccinated, RecoveredVaccinated` + + +Sociodemographic Stratification +-------------------------------- + +The population can be stratified by one sociodemographic dimension denoted **AgeGroup** (can be interpreted more +broadly). The number of age groups is specified in the constructor: + +.. code-block:: cpp + + mio::oseirv::Model model(num_agegroups); + +For stratifications with two or more dimensions, see :doc:`Model Creation <../ode_creation>`. + + +Parameters +---------- + +The model uses the following parameters (time unit: week): + +.. list-table:: + :header-rows: 1 + :widths: 20 25 55 + + * - Mathematical Symbol + - C++ Name / Type + - Description + * - :math:`R_e` + - ``BaselineTransmissibility`` + - Baseline transmissibility (dimensionless); scales the normalized force of infection. + * - :math:`\gamma` + - ``RecoveryRate`` + - Transition / recovery rate (1/week) for E -> I and I -> R. + * - :math:`\delta` + - ``SeasonalityAmplitude`` + - Amplitude of the seasonal modulation :math:`\exp(\delta\,\sin(2\pi(t/52 - t_z + t_s)))`. + * - :math:`t_z` + - ``SeasonalityShiftPerSubtype`` + - Coarse (subtype-specific) seasonal phase shift. + * - :math:`t_s` + - ``SeasonalityShiftPerSeason`` + - Fine seasonal phase adjustment per season. + * - :math:`\lambda_0` + - ``OutsideFoI`` + - External (additive) force of infection, can seed infections. + * - :math:`\rho` + - ``ClusteringExponent`` + - Clustering exponent on the infectious fraction. + * - :math:`m` + - ``SickMixing`` + - Mixing weight for symptomatic (“sick”) contacts in the blended contact matrix. + * - :math:`C^{H}` + - ``ContactPatternsHealthy`` + - Age-structured contact matrix (healthy). Can be time-dependent via damping. + * - :math:`C^{S}` + - ``ContactPatternsSick`` + - Age-structured contact matrix (symptomatic), combined using :math:`m`. + * - :math:`\sigma_i` + - ``CustomIndexArray`` + - Age-specific baseline susceptibility (pre-existing immunity modifier). + * - :math:`VC_i` + - ``VaccineCoverage`` + - Vaccination coverage per age group at season start (share vaccinated). + * - :math:`VE_i` + - ``VaccineEffectiveness`` + - Vaccine effectiveness (reducing effective susceptibility). + * - :math:`\phi_0` + - ``SusceptibleFraction`` + - Fraction of the total population forming the effectively susceptible pool at :math:`t_0`. + +Note: ``VaccineCoverage`` and ``VaccineEffectiveness`` are only used for initialization. Transitions presently +apply identical hazards to vaccinated and unvaccinated susceptible compartments. Future extensions may introduce +differential infection hazards. + + +Initial Conditions +------------------ + +Initial conditions are handled via the **Populations** class. Example for a single age group: + +.. code-block:: cpp + + mio::oseirv::Model model(1); + // Set total population in age group 0 + model.populations.set_total(total0); + + // Initialize vaccinated susceptibles (simple example) + double vc0 = 0.4; // vaccination coverage + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::SusceptibleVaccinated}] = vc0 * total0; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Infected}] = initial_infected; + model.populations[{mio::AgeGroup(0), mio::oseirv::InfectionState::Exposed}] = initial_exposed; + + // Other states (Recovered / RecoveredVaccinated) often 0 at season start + + // Set remaining susceptibles as difference + model.populations.set_difference_from_total( + {mio::AgeGroup(0), mio::oseirv::InfectionState::Susceptible}, total0); + +For age-resolved simulations, repeat for each age group; ``set_difference_from_group_total`` ensures correct +susceptible counts: + +.. code-block:: cpp + + for (auto a = mio::AgeGroup(0); a < num_agegroups; ++a) { + model.populations[{a, mio::oseirv::InfectionState::Exposed}] = exposed0 / num_agegroups; + model.populations[{a, mio::oseirv::InfectionState::Infected}] = infected0 / num_agegroups; + model.populations[{a, mio::oseirv::InfectionState::SusceptibleVaccinated}] = vc[a.get()] * group_size[a.get()]; + model.populations.set_difference_from_group_total( + {a, mio::oseirv::InfectionState::Susceptible}, group_size[a.get()]); + } + + +Simulation +---------- + +Like other ODE models in MEmilio, the SEIRV model can be simulated with standard compartment output or with explicit +flows. Once integrated with utility wrappers (analogous to ``oseir::simulate`` / ``simulate_flows``) usage follows the +same pattern. Example with a Runge–Kutta integrator: + +.. code-block:: cpp + + double t0 = 0.0; // start (weeks) + double tmax = 20.0; // end + double dt = 0.1; // initial step size + + auto integrator = std::make_unique(); + integrator->set_dt_min(0.01); + integrator->set_dt_max(0.5); + integrator->set_rel_tolerance(1e-4); + integrator->set_abs_tolerance(1e-6); + + auto sim = mio::simulate(t0, tmax, dt, model, std::move(integrator)); + +Flow simulation (when explicit flows are required): + +.. code-block:: cpp + + auto flowsim = mio::simulate_flows(t0, tmax, dt, model); + // flowsim[0] = compartment sizes, flowsim[1] = flows + + +Output +------ + +The result of a standard simulation is a ``mio::TimeSeries``: + +.. code-block:: cpp + + auto n_points = static_cast(sim.get_num_time_points()); + Eigen::VectorXd val_i = sim.get_value(i); + double time_i = sim.get_time(i); + auto last_val = sim.get_last_value(); + +Printing and CSV export: + +.. code-block:: cpp + + sim.print_table(); + std::vector labels = {"S","S_V","E","E_V","I","I_V","R","R_V"}; + sim.print_table(labels); + sim.export_csv("seirv_results.csv"); + + + +Contact Changes / Interventions +-------------------------------- + +Time-dependent changes of contact patterns (holidays, interventions) can be modeled via dampings (``add_damping``) on +``ContactPatternsHealthy`` and/or ``ContactPatternsSick``: + +.. code-block:: cpp + + mio::ContactMatrixGroup& cm_h = model.parameters.get>(); + mio::ContactMatrixGroup& cm_s = model.parameters.get>(); + cm_h[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(num_agegroups, num_agegroups, base_contacts)); + cm_s[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(num_agegroups, num_agegroups, base_contacts_sick)); + + // Reduce healthy contacts by 40% starting at week 5 + cm_h[0].add_damping(0.6, mio::SimulationTime(5.0)); + + +Visualization +------------- + +For visualization you can use the Python package :doc:`m-plot <../../python/m-plot>` as in the other models. + + +Literature +---------- + +* Weidemann, F., Remschmidt, C., Buda, S. et al. *Is the impact of childhood influenza vaccination less than expected: a transmission modelling study.* BMC Infectious Diseases 17, 258 (2017). https://doi.org/10.1186/s12879-017-2344-6 + + +Overview of the ``oseirv`` namespace: +------------------------------------- + +.. doxygennamespace:: mio::oseirv +