diff --git a/.github/actions/build-py/action.yml b/.github/actions/build-py/action.yml index 418b3d0380..49cf79b69c 100644 --- a/.github/actions/build-py/action.yml +++ b/.github/actions/build-py/action.yml @@ -28,16 +28,16 @@ runs: cd pycode/memilio-${{ inputs.package }}/ # else stay in root directory fi - /opt/python/cp38-cp38/bin/python -m pip install --upgrade pip setuptools wheel + /opt/python/cp39-cp39/bin/python -m pip install --upgrade pip setuptools wheel /opt/python/cp312-cp312/bin/python -m pip install --upgrade pip setuptools wheel - /opt/python/cp38-cp38/bin/python -m pip install scikit-build scikit-build-core + /opt/python/cp39-cp39/bin/python -m pip install scikit-build scikit-build-core /opt/python/cp312-cp312/bin/python -m pip install scikit-build scikit-build-core # Install setuptools-scm only for memilio-simulation if [ "${{ inputs.package }}" == "simulation" ]; then - /opt/python/cp38-cp38/bin/python -m pip install setuptools-scm + /opt/python/cp39-cp39/bin/python -m pip install setuptools-scm /opt/python/cp312-cp312/bin/python -m pip install setuptools-scm fi - /opt/python/cp38-cp38/bin/python -m build --no-isolation --wheel + /opt/python/cp39-cp39/bin/python -m build --no-isolation --wheel /opt/python/cp312-cp312/bin/python -m build --no-isolation --wheel # Exclude memilio-generation, because its a pure python package, cmake is only used in the build process to retrieve data from cpp if [[ -f "CMakeLists.txt" ]] && [ "${{ inputs.package }}" != "generation" ]; then diff --git a/.github/actions/test-py/action.yml b/.github/actions/test-py/action.yml index 5eeb4d14a0..5d766fc7ba 100644 --- a/.github/actions/test-py/action.yml +++ b/.github/actions/test-py/action.yml @@ -7,7 +7,7 @@ inputs: version: description: "Python version to use for the test." required: false - default: "3.8" + default: "3.9" coverage: description: "Generate coverage report from running the tests, ON or OFF, default OFF." required: false diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index c914695f5f..cbcb97a7a1 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -339,7 +339,7 @@ jobs: if: github.event.pull_request.draft == false strategy: matrix: - version: ["3.8", "3.12"] + version: ["3.9", "3.12"] needs: build-py-generation runs-on: ubuntu-latest steps: @@ -353,7 +353,7 @@ jobs: needs: build-py-simulation strategy: matrix: - version: ["3.8", "3.12"] + version: ["3.9", "3.12"] runs-on: ubuntu-latest steps: - uses: actions/checkout@v6 @@ -366,7 +366,7 @@ jobs: needs: [build-py-surrogatemodel, build-py-simulation] strategy: matrix: - version: ["3.8", "3.12"] + version: ["3.9", "3.12"] runs-on: ubuntu-latest steps: - uses: actions/checkout@v6 diff --git a/cpp/examples/glct_secir.cpp b/cpp/examples/glct_secir.cpp index 1b8243fd8d..5fd7d4ac51 100644 --- a/cpp/examples/glct_secir.cpp +++ b/cpp/examples/glct_secir.cpp @@ -40,7 +40,7 @@ int main() // We need to double the choices of the number of subcompartments for some compartments, // as we define different strains for different transition possibilities. For both strains, the same number of // subcompartments is chosen. The transition probabilities are defined in the StartingProbabilities. - constexpr size_t NumExposed = 2, NumInfectedNoSymptoms = 6, NumInfectedSymptoms = 2, NumInfectedSevere = 2, + constexpr size_t NumExposed = 2, NumInfectedNoSymptoms = 6, NumInfectedSymptoms = 2, NumInfectedSevere = 3, NumInfectedCritical = 10; using Model = mio::glsecir::Model; @@ -62,6 +62,7 @@ int main() const ScalarType recoveredPerInfectedNoSymptoms = 0.09; const ScalarType severePerInfectedSymptoms = 0.2; const ScalarType criticalPerSevere = 0.25; + const ScalarType deathsPerSevere = 0.; const ScalarType deathsPerCritical = 0.3; // Define the initial values with the distribution of the population into subcompartments. @@ -77,7 +78,7 @@ int main() 10 * (1 - recoveredPerInfectedNoSymptoms), 20 * recoveredPerInfectedNoSymptoms, 10 * recoveredPerInfectedNoSymptoms, 10 * recoveredPerInfectedNoSymptoms}, {50 * severePerInfectedSymptoms, 50 * (1 - severePerInfectedSymptoms)}, - {50 * criticalPerSevere, 50 * (1 - criticalPerSevere)}, + {50 * criticalPerSevere, 50 * deathsPerSevere, 50 * (1 - criticalPerSevere - deathsPerSevere)}, {10 * deathsPerCritical, 10 * deathsPerCritical, 5 * deathsPerCritical, 3 * deathsPerCritical, 2 * deathsPerCritical, 10 * (1 - deathsPerCritical), 10 * (1 - deathsPerCritical), 5 * (1 - deathsPerCritical), 3 * (1 - deathsPerCritical), 2 * (1 - deathsPerCritical)}, @@ -136,9 +137,8 @@ int main() Eigen::VectorX StartingProbabilitiesInfectedNoSymptoms = Eigen::VectorX::Zero(LctState::get_num_subcompartments()); StartingProbabilitiesInfectedNoSymptoms[0] = 1 - recoveredPerInfectedNoSymptoms; - StartingProbabilitiesInfectedNoSymptoms[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = - recoveredPerInfectedNoSymptoms; + StartingProbabilitiesInfectedNoSymptoms[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = recoveredPerInfectedNoSymptoms; model.parameters.get>() = StartingProbabilitiesInfectedNoSymptoms; // Define equal TransitionMatrices for the strains. @@ -155,10 +155,9 @@ int main() // InfectedSymptoms. Eigen::VectorX StartingProbabilitiesInfectedSymptoms = Eigen::VectorX::Zero(LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedSymptoms[0] = severePerInfectedSymptoms; - StartingProbabilitiesInfectedSymptoms[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = - 1 - severePerInfectedSymptoms; + StartingProbabilitiesInfectedSymptoms[0] = severePerInfectedSymptoms; + StartingProbabilitiesInfectedSymptoms[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = 1 - severePerInfectedSymptoms; model.parameters.get>() = StartingProbabilitiesInfectedSymptoms; model.parameters.get>() = @@ -170,25 +169,29 @@ int main() // InfectedSevere. Eigen::VectorX StartingProbabilitiesInfectedSevere = Eigen::VectorX::Zero(LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedSevere[0] = criticalPerSevere; - StartingProbabilitiesInfectedSevere[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = - 1 - criticalPerSevere; + StartingProbabilitiesInfectedSevere[0] = criticalPerSevere; + StartingProbabilitiesInfectedSevere[(Eigen::Index)( + LctState::get_num_subcompartments() / 3.)] = deathsPerSevere; + StartingProbabilitiesInfectedSevere[2 * (Eigen::Index)( + LctState::get_num_subcompartments() / + 3.)] = 1 - criticalPerSevere - deathsPerSevere; model.parameters.get>() = StartingProbabilitiesInfectedSevere; model.parameters.get>() = mio::glsecir::TransitionMatrixInfectedSevereToInfectedCritical().get_default( - (size_t)(LctState::get_num_subcompartments() / 2.), timeInfectedSevere); + (size_t)(LctState::get_num_subcompartments() / 3.), timeInfectedSevere); + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedSevereToDead().get_default( + (size_t)(LctState::get_num_subcompartments() / 3.), timeInfectedSevere); model.parameters.get>() = mio::glsecir::TransitionMatrixInfectedSevereToRecovered().get_default( - (size_t)(LctState::get_num_subcompartments() / 2.), timeInfectedSevere); + (size_t)(LctState::get_num_subcompartments() / 3.), timeInfectedSevere); // InfectedCritical. Eigen::VectorX StartingProbabilitiesInfectedCritical = Eigen::VectorX::Zero(LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedCritical[0] = deathsPerCritical; - StartingProbabilitiesInfectedCritical[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = - 1 - deathsPerCritical; + StartingProbabilitiesInfectedCritical[0] = deathsPerCritical; + StartingProbabilitiesInfectedCritical[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = 1 - deathsPerCritical; model.parameters.get>() = StartingProbabilitiesInfectedCritical; model.parameters.get>() = diff --git a/cpp/models/glct_secir/model.h b/cpp/models/glct_secir/model.h index b27679c836..a83e6446eb 100644 --- a/cpp/models/glct_secir/model.h +++ b/cpp/models/glct_secir/model.h @@ -268,15 +268,24 @@ class Model params.template get>().transpose() * y.segment(LctState::template get_first_index(), dimensionInfectedSevereToInfectedCritical); + // Flow from InfectedSevere To Dead. + size_t dimensionInfectedSevereToDead = params.template get>().rows(); + dydt.segment(LctState::template get_first_index() + + dimensionInfectedSevereToInfectedCritical, + dimensionInfectedSevereToDead) += + params.template get>().transpose() * + y.segment(LctState::template get_first_index() + + dimensionInfectedSevereToInfectedCritical, + dimensionInfectedSevereToDead); // Flow from InfectedSevere To Recovered. size_t dimensionInfectedSevereToRecovered = params.template get>().rows(); dydt.segment(LctState::template get_first_index() + - dimensionInfectedSevereToInfectedCritical, + dimensionInfectedSevereToInfectedCritical + dimensionInfectedSevereToDead, dimensionInfectedSevereToRecovered) += params.template get>().transpose() * y.segment(LctState::template get_first_index() + - dimensionInfectedSevereToInfectedCritical, + dimensionInfectedSevereToInfectedCritical + dimensionInfectedSevereToDead, dimensionInfectedSevereToRecovered); // Add flow directly to Recovered compartment. dydt[LctState::template get_first_index()] += @@ -284,7 +293,7 @@ class Model Eigen::VectorX::Ones(dimensionInfectedSevereToRecovered)) .transpose() * y.segment(LctState::template get_first_index() + - dimensionInfectedSevereToInfectedCritical, + dimensionInfectedSevereToInfectedCritical + dimensionInfectedSevereToDead, dimensionInfectedSevereToRecovered); // --- InfectedCritical. --- @@ -325,7 +334,14 @@ class Model dimensionInfectedCriticalToRecovered); // --- Dead. --- - dydt[LctState::template get_first_index()] = + dydt[LctState::template get_first_index()] += + -(params.template get>() * + Eigen::VectorX::Ones(dimensionInfectedSevereToDead)) + .transpose() * + y.segment(LctState::template get_first_index() + + dimensionInfectedSevereToInfectedCritical, + dimensionInfectedSevereToDead); + dydt[LctState::template get_first_index()] += -(params.template get>() * Eigen::VectorX::Ones(dimensionInfectedCriticalToDead)) .transpose() * diff --git a/cpp/models/glct_secir/parameters.h b/cpp/models/glct_secir/parameters.h index d46c25bc05..8c72df0f8e 100644 --- a/cpp/models/glct_secir/parameters.h +++ b/cpp/models/glct_secir/parameters.h @@ -263,6 +263,31 @@ struct TransitionMatrixInfectedSevereToInfectedCritical { } }; +/** + * @brief Transition matrix of the phase-type distribution describing the stay time in the InfectedSevere + * compartment before dying. + */ +template +struct TransitionMatrixInfectedSevereToDead { + using Type = Eigen::MatrixX; + /** + * @brief Default parameters can be used to get an Erlang distributed stay time in InfectedSevere compartment + * before dying. + * @param[in] dimension Number of rows/columns of the transition matrix. + * @param[in] time Average time spent in InfectedSevere before dying in day unit. + */ + static Type get_default(size_t dimension, FP time = 1.) + { + Eigen::MatrixX def = Eigen::VectorX::Constant(dimension, -(FP)dimension / time).asDiagonal(); + def.diagonal(1).setConstant((FP)dimension / time); + return def; + } + static std::string name() + { + return "TransitionMatrixInfectedSevereToDead"; + } +}; + /** * @brief Transition matrix of the phase-type distribution describing the stay time in the InfectedSevere * compartment before recovery. @@ -461,10 +486,11 @@ using ParametersBase = TransitionMatrixInfectedNoSymptomsToRecovered, StartingProbabilitiesInfectedSymptoms, TransitionMatrixInfectedSymptomsToInfectedSevere, TransitionMatrixInfectedSymptomsToRecovered, StartingProbabilitiesInfectedSevere, TransitionMatrixInfectedSevereToInfectedCritical, - TransitionMatrixInfectedSevereToRecovered, StartingProbabilitiesInfectedCritical, - TransitionMatrixInfectedCriticalToDead, TransitionMatrixInfectedCriticalToRecovered, - TransmissionProbabilityOnContact, ContactPatterns, RelativeTransmissionNoSymptoms, - RiskOfInfectionFromSymptomatic, StartDay, Seasonality>; + TransitionMatrixInfectedSevereToDead, TransitionMatrixInfectedSevereToRecovered, + StartingProbabilitiesInfectedCritical, TransitionMatrixInfectedCriticalToDead, + TransitionMatrixInfectedCriticalToRecovered, TransmissionProbabilityOnContact, + ContactPatterns, RelativeTransmissionNoSymptoms, RiskOfInfectionFromSymptomatic, + StartDay, Seasonality>; /// @brief Parameters of an GLCT-SECIR model. template @@ -523,6 +549,8 @@ class Parameters : public ParametersBase this->template get>().rows()) || (this->template get>().cols() != this->template get>().rows()) || + (this->template get>().cols() != + this->template get>().rows()) || (this->template get>().cols() != this->template get>().rows()) || (this->template get>().cols() != @@ -559,6 +587,7 @@ class Parameters : public ParametersBase if (this->template get>().rows() != this->template get>().rows() + + this->template get>().rows() + this->template get>().rows()) { log_error("Constraint check: Dimensions of StartingProbabilitiesInfectedSevere and " "TransitionMatrices of InfectedSevere compartment are not matching."); @@ -653,6 +682,14 @@ class Parameters : public ParametersBase "flow InfectedSevereToInfectedCritical."); return true; } + if (((this->template get>() * + Eigen::VectorX::Ones(this->template get>().rows())) + .array() > 1e-10) + .any()) { + log_warning("Constraint check: The entries of TransitionMatrixInfectedSevereToDead lead to a negative " + "flow InfectedSevereToDead."); + return true; + } if (((this->template get>() * Eigen::VectorX::Ones(this->template get>().rows())) .array() > 1e-10) diff --git a/cpp/models/ide_secir/parameters_io.h b/cpp/models/ide_secir/parameters_io.h index c346553d8f..f48b91c06c 100644 --- a/cpp/models/ide_secir/parameters_io.h +++ b/cpp/models/ide_secir/parameters_io.h @@ -56,9 +56,6 @@ namespace isecir * using the means of the respective TransitionDistribution. * The flow InfectedNoSymptomsToInfectedSymptoms is calculated with the standard formula from the IDE model * using the results for ExposedToInfectedNoSymptoms. -* Throughout this RKI-based initialization, we assume that individuals can only die from InfectedCritical and not from -* InfectedSevere, i.e. this initialization routine assumes that the probability from transitioning from InfectedSevere -* to Dead is 0. * * The number of deaths used in the model is computed by shifting the reported RKI data according to the mean values of the transitions * InfectedSymptomsToInfectedSevere, InfectedSevereToInfectedCritical and InfectedCriticalToDead. diff --git a/cpp/models/lct_secir/infection_state.h b/cpp/models/lct_secir/infection_state.h index 342b891970..198e6f4e2c 100644 --- a/cpp/models/lct_secir/infection_state.h +++ b/cpp/models/lct_secir/infection_state.h @@ -56,10 +56,11 @@ enum class InfectionTransition InfectedSymptomsToInfectedSevere = 4, InfectedSymptomsToRecovered = 5, InfectedSevereToInfectedCritical = 6, - InfectedSevereToRecovered = 7, - InfectedCriticalToDead = 8, - InfectedCriticalToRecovered = 9, - Count = 10 + InfectedSevereToDead = 7, + InfectedSevereToRecovered = 8, + InfectedCriticalToDead = 9, + InfectedCriticalToRecovered = 10, + Count = 11 }; } // namespace lsecir diff --git a/cpp/models/lct_secir/model.h b/cpp/models/lct_secir/model.h index 925da29fbf..8816230fa3 100644 --- a/cpp/models/lct_secir/model.h +++ b/cpp/models/lct_secir/model.h @@ -260,7 +260,9 @@ class Model : public CompartmentalModel>()[Group]); + dydt[Ri] += dydt[ICri_first_index] * (1 - params.template get>()[Group] - + params.template get>()[Group]); + dydt[Di] = dydt[ICri_first_index] * params.template get>()[Group]; dydt[ICri_first_index] = dydt[ICri_first_index] * params.template get>()[Group]; for (size_t subcomp = 0; subcomp < LctStateGroup::template get_num_subcompartments() - 1; @@ -276,7 +278,7 @@ class Model : public CompartmentalModel>()[Group]) * y[Ri - 1]; dydt[Ri - 1] -= flow; dydt[Ri] = dydt[Ri] + (1 - params.template get>()[Group]) * flow; - dydt[Di] = params.template get>()[Group] * flow; + dydt[Di] = dydt[Di] + params.template get>()[Group] * flow; // Function call for next group if applicable. if constexpr (Group + 1 < num_groups) { diff --git a/cpp/models/lct_secir/parameters.h b/cpp/models/lct_secir/parameters.h index e8ab766a5e..34d722186b 100644 --- a/cpp/models/lct_secir/parameters.h +++ b/cpp/models/lct_secir/parameters.h @@ -236,6 +236,22 @@ struct CriticalPerSevere { } }; +/** + * @brief The percentage of dead patients per hospitalized patients for each group in the SECIR model. + */ +template +struct DeathsPerSevere { + using Type = Eigen::VectorX>; + static Type get_default(size_t size) + { + return Type::Constant(size, 1, 0.); + } + static std::string name() + { + return "DeathsPerSevere"; + } +}; + /** * @brief The percentage of dead patients per ICU patients for each group in the SECIR model. */ @@ -295,7 +311,7 @@ using ParametersBase = TimeInfectedCritical, TransmissionProbabilityOnContact, ContactPatterns, RelativeTransmissionNoSymptoms, RiskOfInfectionFromSymptomatic, RecoveredPerInfectedNoSymptoms, SeverePerInfectedSymptoms, CriticalPerSevere, - DeathsPerCritical, StartDay, Seasonality>; + DeathsPerSevere, DeathsPerCritical, StartDay, Seasonality>; /** * @brief Parameters of an LCT-SECIR model. @@ -392,6 +408,18 @@ class Parameters : public ParametersBase return true; } + if (this->template get>()[i] < 0.0 || + this->template get>()[i] > 1.0) { + log_error("Constraint check: Parameter DeathsPerSevere smaller {} or larger {}", 0, 1); + return true; + } + + if (this->template get>()[i] + this->template get>()[i] > 1.0) { + log_error("Constraint check: CriticalPerSevere + DeathsPerSevere exceed 1.0 for age group {}.", + static_cast(i)); + return true; + } + if (this->template get>()[i] < 0.0 || this->template get>()[i] > 1.0) { log_error("Constraint check: Parameter DeathsPerCritical smaller {} or larger {}", 0, 1); diff --git a/cpp/models/lct_secir_2_diseases/model.h b/cpp/models/lct_secir_2_diseases/model.h index a869b00988..b759d9253c 100644 --- a/cpp/models/lct_secir_2_diseases/model.h +++ b/cpp/models/lct_secir_2_diseases/model.h @@ -354,8 +354,10 @@ class Model : public CompartmentalModel>()[Group]); + // Split flow from last subcompartment of InfectedSevere_1a between Recovered_1a, Dead1a and InfectedCritical_1a. + dydt[Ri_a] += dydt[ICri_1a_first_index] * (1 - params.template get>()[Group] - + params.template get>()[Group]); + dydt[Di_a] = dydt[ICri_1a_first_index] * params.template get>()[Group]; dydt[ICri_1a_first_index] = dydt[ICri_1a_first_index] * params.template get>()[Group]; for (size_t subcomp = 0; subcomp < LctStateGroup::template get_num_subcompartments() - 1; @@ -366,8 +368,10 @@ class Model : public CompartmentalModel>()[Group]); + // Split flow from last subcompartment of InfectedSevere_1b between Recovered_1b, Dead_1b and InfectedCritical_1b. + dydt[Ri_b] += dydt[ICri_1b_first_index] * (1 - params.template get>()[Group] - + params.template get>()[Group]); + dydt[Di_b] = dydt[ICri_1b_first_index] * params.template get>()[Group]; dydt[ICri_1b_first_index] = dydt[ICri_1b_first_index] * params.template get>()[Group]; for (size_t subcomp = 0; subcomp < LctStateGroup::template get_num_subcompartments() - 1; @@ -512,8 +516,10 @@ class Model : public CompartmentalModel>()[Group]); + // Split flow from last subcompartment of InfectedSevere_2a between Recovered_2ab, Dead_2ab and InfectedCritical_2a. + dydt[Ri_ab] += dydt[ICri_2a_first_index] * (1 - params.template get>()[Group] - + params.template get>()[Group]); + dydt[Di_a] += dydt[ICri_2a_first_index] * params.template get>()[Group]; dydt[ICri_2a_first_index] = dydt[ICri_2a_first_index] * params.template get>()[Group]; for (size_t subcomp = 0; subcomp < LctStateGroup::template get_num_subcompartments() - 1; @@ -524,8 +530,10 @@ class Model : public CompartmentalModel>()[Group]); + // Split flow from last subcompartment of InfectedSevere_2b between Recovered_2ab and InfectedCritical_2b. + dydt[Ri_ab] += dydt[ICri_2b_first_index] * (1 - params.template get>()[Group] - + params.template get>()[Group]); + dydt[Di_b] += dydt[ICri_2b_first_index] * params.template get>()[Group]; dydt[ICri_2b_first_index] = dydt[ICri_2b_first_index] * params.template get>()[Group]; for (size_t subcomp = 0; subcomp < LctStateGroup::template get_num_subcompartments() - 1; diff --git a/cpp/models/lct_secir_2_diseases/parameters.h b/cpp/models/lct_secir_2_diseases/parameters.h index ab83447613..0a845895ed 100644 --- a/cpp/models/lct_secir_2_diseases/parameters.h +++ b/cpp/models/lct_secir_2_diseases/parameters.h @@ -369,6 +369,22 @@ struct CriticalPerSevere_a { } }; +/** + * @brief The percentage of dead patients per hospitalized patients for each group in the SECIR model. + */ +template +struct DeathsPerSevere_a { + using Type = Eigen::VectorX>; + static Type get_default(size_t size) + { + return Type::Constant(size, 1, 0.); + } + static std::string name() + { + return "DeathsPerSevere_a"; + } +}; + /** * @brief The percentage of dead patients per ICU patients for disease a for each group. */ @@ -433,6 +449,22 @@ struct CriticalPerSevere_b { } }; +/** + * @brief The percentage of dead patients per hospitalized patients for each group in the SECIR model. + */ +template +struct DeathsPerSevere_b { + using Type = Eigen::VectorX>; + static Type get_default(size_t size) + { + return Type::Constant(size, 1, 0.); + } + static std::string name() + { + return "DeathsPerSevere_b"; + } +}; + /** * @brief The percentage of dead patients per ICU patients for disease b for each group. */ @@ -487,16 +519,15 @@ struct Seasonality { }; template -using ParametersBase = - ParameterSet, TimeInfectedNoSymptoms_a, TimeInfectedSymptoms_a, TimeInfectedSevere_a, - TimeInfectedCritical_a, TimeExposed_b, TimeInfectedNoSymptoms_b, - TimeInfectedSymptoms_b, TimeInfectedSevere_b, TimeInfectedCritical_b, - TransmissionProbabilityOnContact_a, TransmissionProbabilityOnContact_b, ContactPatterns, - RelativeTransmissionNoSymptoms_a, RiskOfInfectionFromSymptomatic_a, - RecoveredPerInfectedNoSymptoms_a, SeverePerInfectedSymptoms_a, CriticalPerSevere_a, - DeathsPerCritical_a, RelativeTransmissionNoSymptoms_b, RiskOfInfectionFromSymptomatic_b, - RecoveredPerInfectedNoSymptoms_b, SeverePerInfectedSymptoms_b, CriticalPerSevere_b, - DeathsPerCritical_b, StartDay, Seasonality>; +using ParametersBase = ParameterSet< + TimeExposed_a, TimeInfectedNoSymptoms_a, TimeInfectedSymptoms_a, TimeInfectedSevere_a, + TimeInfectedCritical_a, TimeExposed_b, TimeInfectedNoSymptoms_b, TimeInfectedSymptoms_b, + TimeInfectedSevere_b, TimeInfectedCritical_b, TransmissionProbabilityOnContact_a, + TransmissionProbabilityOnContact_b, ContactPatterns, RelativeTransmissionNoSymptoms_a, + RiskOfInfectionFromSymptomatic_a, RecoveredPerInfectedNoSymptoms_a, SeverePerInfectedSymptoms_a, + CriticalPerSevere_a, DeathsPerSevere_a, DeathsPerCritical_a, RelativeTransmissionNoSymptoms_b, + RiskOfInfectionFromSymptomatic_b, RecoveredPerInfectedNoSymptoms_b, SeverePerInfectedSymptoms_b, + CriticalPerSevere_b, DeathsPerSevere_b, DeathsPerCritical_b, StartDay, Seasonality>; /** * @brief Parameters of an LCT-SECIR-2-DISEASES model. @@ -662,6 +693,32 @@ class Parameters : public ParametersBase return true; } + if (this->template get>()[i] < 0.0 || + this->template get>()[i] > 1.0) { + log_error("Constraint check: Parameter DeathsPerSevere_a smaller {:d} or larger {:d}", 0, 1); + return true; + } + + if (this->template get>()[i] < 0.0 || + this->template get>()[i] > 1.0) { + log_error("Constraint check: Parameter DeathsPerSevere_b smaller {:d} or larger {:d}", 0, 1); + return true; + } + + if (this->template get>()[i] + this->template get>()[i] > + 1.0) { + log_error("Constraint check: CriticalPerSevere_a + DeathsPerSevere_a exceed 1.0 for age group {}.", + static_cast(i)); + return true; + } + + if (this->template get>()[i] + this->template get>()[i] > + 1.0) { + log_error("Constraint check: CriticalPerSevere_b + DeathsPerSevere_b exceed 1.0 for age group {}.", + static_cast(i)); + return true; + } + if (this->template get>()[i] < 0.0 || this->template get>()[i] > 1.0) { log_error("Constraint check: Parameter DeathsPerCritical_a smaller {:d} or larger {:d}", 0, 1); diff --git a/cpp/tests/test_glct_secir.cpp b/cpp/tests/test_glct_secir.cpp index 64614cae67..267f90c173 100755 --- a/cpp/tests/test_glct_secir.cpp +++ b/cpp/tests/test_glct_secir.cpp @@ -1,4 +1,4 @@ -/* +/* * Copyright (C) 2020-2026 MEmilio * * Authors: Lena Ploetzke @@ -18,6 +18,7 @@ * limitations under the License. */ +#include "glct_secir/infection_state.h" #include "glct_secir/model.h" #include "glct_secir/parameters.h" #include "memilio/config.h" @@ -37,13 +38,18 @@ TEST(TestGLCTSecir, testEvalRightHandSide) // Define initial values, parameters and numbers of subcompartments according to the choices of the // testEvalRightHandSide of the LCT testing suite. For more details, // we refer to the example glct_secir.cpp. - using Model = mio::glsecir::Model; + using Model = mio::glsecir::Model; using LctState = Model::LctState; using InfectionState = LctState::InfectionState; Model model; // Set parameters such that the stay times are Erlang-distributed as in the corresponding LCT model. + ScalarType recoveredPerInfectedNoSymptoms = 0.09; + ScalarType severePerInfectedSymptoms = 0.2; + ScalarType criticalPerSevere = 0.25; + ScalarType deathsPerSevere = 0.0; + ScalarType deathsPerCritical = 0.3; // Exposed. // Default functions are used to set the parameters but the corresponding dimensions have to be set manually. model.parameters.get>() = @@ -55,9 +61,9 @@ TEST(TestGLCTSecir, testEvalRightHandSide) // InfectedNoSymptoms. Eigen::VectorX StartingProbabilitiesInfectedNoSymptoms = Eigen::VectorX::Zero( (Eigen::Index)LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedNoSymptoms[0] = 1 - 0.09; - StartingProbabilitiesInfectedNoSymptoms[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = 0.09; + StartingProbabilitiesInfectedNoSymptoms[0] = 1 - recoveredPerInfectedNoSymptoms; + StartingProbabilitiesInfectedNoSymptoms[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = recoveredPerInfectedNoSymptoms; model.parameters.get>() = StartingProbabilitiesInfectedNoSymptoms; model.parameters.get>() = @@ -69,9 +75,9 @@ TEST(TestGLCTSecir, testEvalRightHandSide) // InfectedSymptoms. Eigen::VectorX StartingProbabilitiesInfectedSymptoms = Eigen::VectorX::Zero( (Eigen::Index)LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedSymptoms[0] = 0.2; - StartingProbabilitiesInfectedSymptoms[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = 1 - 0.2; + StartingProbabilitiesInfectedSymptoms[0] = severePerInfectedSymptoms; + StartingProbabilitiesInfectedSymptoms[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = 1 - severePerInfectedSymptoms; model.parameters.get>() = StartingProbabilitiesInfectedSymptoms; model.parameters.get>() = @@ -83,23 +89,29 @@ TEST(TestGLCTSecir, testEvalRightHandSide) // InfectedSevere. Eigen::VectorX StartingProbabilitiesInfectedSevere = Eigen::VectorX::Zero( (Eigen::Index)LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedSevere[0] = 0.25; - StartingProbabilitiesInfectedSevere[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = 1 - 0.25; + StartingProbabilitiesInfectedSevere[0] = criticalPerSevere; + StartingProbabilitiesInfectedSevere[(Eigen::Index)( + LctState::get_num_subcompartments() / 3.)] = deathsPerSevere; + StartingProbabilitiesInfectedSevere[2 * (Eigen::Index)( + LctState::get_num_subcompartments() / + 3.)] = 1. - criticalPerSevere - deathsPerSevere; model.parameters.get>() = StartingProbabilitiesInfectedSevere; model.parameters.get>() = mio::glsecir::TransitionMatrixInfectedSevereToInfectedCritical().get_default( - (size_t)(LctState::get_num_subcompartments() / 2.), 9.5); + (size_t)(LctState::get_num_subcompartments() / 3.), 9.5); + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedSevereToDead().get_default( + (size_t)(LctState::get_num_subcompartments() / 3.), 9.5); model.parameters.get>() = mio::glsecir::TransitionMatrixInfectedSevereToRecovered().get_default( - (size_t)(LctState::get_num_subcompartments() / 2.), 9.5); + (size_t)(LctState::get_num_subcompartments() / 3.), 9.5); // InfectedCritical. Eigen::VectorX StartingProbabilitiesInfectedCritical = Eigen::VectorX::Zero( (Eigen::Index)LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedCritical[0] = 0.3; - StartingProbabilitiesInfectedCritical[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = 1 - 0.3; + StartingProbabilitiesInfectedCritical[0] = deathsPerCritical; + StartingProbabilitiesInfectedCritical[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = 1 - deathsPerCritical; model.parameters.get>() = StartingProbabilitiesInfectedCritical; model.parameters.get>() = @@ -130,7 +142,14 @@ TEST(TestGLCTSecir, testEvalRightHandSide) {30 * StartingProbabilitiesInfectedSymptoms[0], 20 * StartingProbabilitiesInfectedSymptoms[0], 30 * (1 - StartingProbabilitiesInfectedSymptoms[0]), 20 * (1 - StartingProbabilitiesInfectedSymptoms[0])}, {40 * StartingProbabilitiesInfectedSevere[0], 10 * StartingProbabilitiesInfectedSevere[0], - 40 * (1 - StartingProbabilitiesInfectedSevere[0]), 10 * (1 - StartingProbabilitiesInfectedSevere[0])}, + 40 * StartingProbabilitiesInfectedSevere[(Eigen::Index)( + LctState::get_num_subcompartments() / 3.)], + 10 * StartingProbabilitiesInfectedSevere[(Eigen::Index)( + LctState::get_num_subcompartments() / 3.)], + 40 * StartingProbabilitiesInfectedSevere + [2 * (Eigen::Index)(LctState::get_num_subcompartments() / 3.)], + 10 * StartingProbabilitiesInfectedSevere + [2 * (Eigen::Index)(LctState::get_num_subcompartments() / 3.)]}, {10 * StartingProbabilitiesInfectedCritical[0], 20 * StartingProbabilitiesInfectedCritical[0], 10 * (1 - StartingProbabilitiesInfectedCritical[0]), 20 * (1 - StartingProbabilitiesInfectedCritical[0])}, {20}, @@ -149,9 +168,15 @@ TEST(TestGLCTSecir, testEvalRightHandSide) model.get_derivatives(pop, pop, 0, dydt); // This vector is the equivalent of the result defined in the test suite testEvalRightHandSide of the LCT model. Eigen::VectorX compare(LctState::Count); - compare << -15.3409, -3.4091, 6.25, -17.5 * 0.91, 15 * 0.91, 0 * 0.91, -17.5 * 0.09, 15 * 0.09, 0 * 0.09, - 3.3052 * 0.2, 3.4483 * 0.2, 3.3052 * 0.8, 3.4483 * 0.8, -7.0417 * 0.25, 6.3158 * 0.25, -7.0417 * 0.75, - 6.3158 * 0.75, -2.2906 * 0.3, -2.8169 * 0.3, -2.2906 * 0.7, -2.8169 * 0.7, 12.3899, 1.6901; + compare << -15.3409, -3.4091, 6.25, -17.5 * (1. - recoveredPerInfectedNoSymptoms), + 15 * (1. - recoveredPerInfectedNoSymptoms), 0 * (1. - recoveredPerInfectedNoSymptoms), + -17.5 * recoveredPerInfectedNoSymptoms, 15 * recoveredPerInfectedNoSymptoms, 0 * recoveredPerInfectedNoSymptoms, + 3.3052 * severePerInfectedSymptoms, 3.4483 * severePerInfectedSymptoms, + 3.3052 * (1. - severePerInfectedSymptoms), 3.4483 * (1 - severePerInfectedSymptoms), + -7.0417 * criticalPerSevere, 6.3158 * criticalPerSevere, -7.0417 * deathsPerSevere, 6.3158 * deathsPerSevere, + -7.0417 * (1. - criticalPerSevere - deathsPerSevere), 6.3158 * (1. - criticalPerSevere - deathsPerSevere), + -2.2906 * deathsPerCritical, -2.8169 * deathsPerCritical, -2.2906 * (1. - deathsPerCritical), + -2.8169 * (1. - deathsPerCritical), 12.3899, 1.6901; for (size_t i = 0; i < LctState::Count; i++) { EXPECT_NEAR(compare[i], dydt[i], 1e-3) << "Condition failed at index: " << i; @@ -162,7 +187,7 @@ TEST(TestGLCTSecir, testEvalRightHandSide) class ModelTestGLCTSecir : public testing::Test { public: - using Model = mio::glsecir::Model; + using Model = mio::glsecir::Model; using LctState = Model::LctState; using InfectionState = LctState::InfectionState; @@ -171,6 +196,11 @@ class ModelTestGLCTSecir : public testing::Test { model = new Model(); // --- Set parameters. --- + ScalarType recoveredPerInfectedNoSymptoms = 0.09; + ScalarType severePerInfectedSymptoms = 0.2; + ScalarType criticalPerSevere = 0.25; + ScalarType deathsPerSevere = 0.0; + ScalarType deathsPerCritical = 0.3; // Exposed. model->parameters.get>() = mio::glsecir::StartingProbabilitiesExposed().get_default( @@ -181,9 +211,10 @@ class ModelTestGLCTSecir : public testing::Test // InfectedNoSymptoms. Eigen::VectorX StartingProbabilitiesInfectedNoSymptoms = Eigen::VectorX::Zero( (Eigen::Index)LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedNoSymptoms[0] = 1 - 0.09; - StartingProbabilitiesInfectedNoSymptoms[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = 0.09; + StartingProbabilitiesInfectedNoSymptoms[0] = 1 - recoveredPerInfectedNoSymptoms; + StartingProbabilitiesInfectedNoSymptoms[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = + recoveredPerInfectedNoSymptoms; model->parameters.get>() = StartingProbabilitiesInfectedNoSymptoms; model->parameters.get>() = @@ -195,9 +226,10 @@ class ModelTestGLCTSecir : public testing::Test // InfectedSymptoms. Eigen::VectorX StartingProbabilitiesInfectedSymptoms = Eigen::VectorX::Zero( (Eigen::Index)LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedSymptoms[0] = 0.2; - StartingProbabilitiesInfectedSymptoms[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = 1 - 0.2; + StartingProbabilitiesInfectedSymptoms[0] = severePerInfectedSymptoms; + StartingProbabilitiesInfectedSymptoms[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = + 1 - severePerInfectedSymptoms; model->parameters.get>() = StartingProbabilitiesInfectedSymptoms; model->parameters.get>() = @@ -209,23 +241,29 @@ class ModelTestGLCTSecir : public testing::Test // InfectedSevere. Eigen::VectorX StartingProbabilitiesInfectedSevere = Eigen::VectorX::Zero( (Eigen::Index)LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedSevere[0] = 0.25; - StartingProbabilitiesInfectedSevere[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = 1 - 0.25; + StartingProbabilitiesInfectedSevere[0] = criticalPerSevere; + StartingProbabilitiesInfectedSevere[(Eigen::Index)( + LctState::get_num_subcompartments() / 3.)] = deathsPerSevere; + StartingProbabilitiesInfectedSevere + [2 * (Eigen::Index)(LctState::get_num_subcompartments() / 3.)] = + 1. - criticalPerSevere - deathsPerSevere; model->parameters.get>() = StartingProbabilitiesInfectedSevere; model->parameters.get>() = mio::glsecir::TransitionMatrixInfectedSevereToInfectedCritical().get_default( - (size_t)(LctState::get_num_subcompartments() / 2.), 9.5); + (size_t)(LctState::get_num_subcompartments() / 3.), 9.5); + model->parameters.get>() = + mio::glsecir::TransitionMatrixInfectedSevereToDead().get_default( + (size_t)(LctState::get_num_subcompartments() / 3.), 9.5); model->parameters.get>() = mio::glsecir::TransitionMatrixInfectedSevereToRecovered().get_default( - (size_t)(LctState::get_num_subcompartments() / 2.), 9.5); + (size_t)(LctState::get_num_subcompartments() / 3.), 9.5); // InfectedCritical. Eigen::VectorX StartingProbabilitiesInfectedCritical = Eigen::VectorX::Zero( (Eigen::Index)LctState::get_num_subcompartments()); - StartingProbabilitiesInfectedCritical[0] = 0.3; - StartingProbabilitiesInfectedCritical[( - Eigen::Index)(LctState::get_num_subcompartments() / 2.)] = 1 - 0.3; + StartingProbabilitiesInfectedCritical[0] = deathsPerCritical; + StartingProbabilitiesInfectedCritical[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = 1 - deathsPerCritical; model->parameters.get>() = StartingProbabilitiesInfectedCritical; model->parameters.get>() = @@ -256,7 +294,11 @@ class ModelTestGLCTSecir : public testing::Test 10 * (1 - StartingProbabilitiesInfectedNoSymptoms[0]), 10 * (1 - StartingProbabilitiesInfectedNoSymptoms[0])}, {50 * StartingProbabilitiesInfectedSymptoms[0], 50 * (1 - StartingProbabilitiesInfectedSymptoms[0])}, - {50 * StartingProbabilitiesInfectedSevere[0], 50 * (1 - StartingProbabilitiesInfectedSevere[0])}, + {50 * StartingProbabilitiesInfectedSevere[0], + 50 * StartingProbabilitiesInfectedSevere[(Eigen::Index)( + LctState::get_num_subcompartments() / 3.)], + 50 * StartingProbabilitiesInfectedSevere + [2 * (Eigen::Index)(LctState::get_num_subcompartments() / 3.)]}, {10 * StartingProbabilitiesInfectedCritical[0], 10 * StartingProbabilitiesInfectedCritical[0], 5 * StartingProbabilitiesInfectedCritical[0], 3 * StartingProbabilitiesInfectedCritical[0], 2 * StartingProbabilitiesInfectedCritical[0], 10 * (1 - StartingProbabilitiesInfectedCritical[0]), @@ -319,7 +361,7 @@ TEST_F(ModelTestGLCTSecir, testConstraintsModel) EXPECT_FALSE(constraint_check); // Check if the number of subcompartments does not match the dimension of the vector with StartingProbabilities. - Eigen::VectorX wrong_size = Eigen::VectorX::Zero(3); + Eigen::VectorX wrong_size = Eigen::VectorX::Zero(4); wrong_size[0] = 1; // Exposed. model->parameters.get>() = wrong_size; @@ -404,7 +446,7 @@ TEST_F(ModelTestGLCTSecir, testConstraintsParameters) (size_t)(LctState::get_num_subcompartments() / 2.), 7.1); // Check non matching dimensions of TransitionMatrix and vector with StartingProbabilities. - Eigen::VectorX wrong_size = Eigen::VectorX::Zero(3); + Eigen::VectorX wrong_size = Eigen::VectorX::Zero(4); wrong_size[0] = 1; // Exposed. model->parameters.get>() = wrong_size; @@ -450,8 +492,14 @@ TEST_F(ModelTestGLCTSecir, testConstraintsParameters) model->parameters.get>()[1] = -0.1; constraint_check = model->parameters.check_constraints(); EXPECT_TRUE(constraint_check); + model->parameters.get>()[0] = 0.9; + model->parameters.get>()[1] = 0.1; + model->parameters.get>()[2] = 0.1; + constraint_check = model->parameters.check_constraints(); + EXPECT_TRUE(constraint_check); model->parameters.get>()[0] = 1.; model->parameters.get>()[1] = 0.; + model->parameters.get>()[2] = 0.; // --- Check with invalid transition matrices. --- // ExposedToInfectedNoSymptoms. @@ -531,3 +579,151 @@ TEST_F(ModelTestGLCTSecir, testCalculatePopWrongSize) // Reactive log output. mio::set_log_level(mio::LogLevel::warn); } + +TEST(TestGLCTSecir, deathsPerSevere_flows) +{ + // Test that DeathsPerSevere causes ISev->Dead flow. + // CriticalPerSevere=0 blocks the ISev->ICr->Dead path entirely, so + // DeathsPerCritical has no influence and is left at its default values. + + using Model = mio::glsecir::Model; + using LctState = Model::LctState; + using InfectionState = LctState::InfectionState; + + // Initialize a model with one group. + Model model; + + // Set parameters. + ScalarType recoveredPerInfectedNoSymptoms = 0.0; + ScalarType severePerInfectedSymptoms = 1.0; + ScalarType criticalPerSevere = 0.0; // block ISev->ICr->Dead path + ScalarType deathsPerSevere = 0.1; + ScalarType deathsPerCritical = 0.0; + // Exposed. + // Default functions are used to set the parameters but the corresponding dimensions have to be set manually. + model.parameters.get>() = + mio::glsecir::StartingProbabilitiesExposed().get_default( + LctState::get_num_subcompartments()); + model.parameters.get>() = + mio::glsecir::TransitionMatrixExposedToInfectedNoSymptoms().get_default( + LctState::get_num_subcompartments(), 3.2); + // InfectedNoSymptoms. + Eigen::VectorX StartingProbabilitiesInfectedNoSymptoms = Eigen::VectorX::Zero( + (Eigen::Index)LctState::get_num_subcompartments()); + StartingProbabilitiesInfectedNoSymptoms[0] = 1 - recoveredPerInfectedNoSymptoms; + StartingProbabilitiesInfectedNoSymptoms[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = recoveredPerInfectedNoSymptoms; + model.parameters.get>() = + StartingProbabilitiesInfectedNoSymptoms; + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedNoSymptomsToInfectedSymptoms().get_default( + (size_t)(LctState::get_num_subcompartments() / 2.), 2.); + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedNoSymptomsToRecovered().get_default( + (size_t)(LctState::get_num_subcompartments() / 2.), 2.); + // InfectedSymptoms. + Eigen::VectorX StartingProbabilitiesInfectedSymptoms = Eigen::VectorX::Zero( + (Eigen::Index)LctState::get_num_subcompartments()); + StartingProbabilitiesInfectedSymptoms[0] = severePerInfectedSymptoms; + StartingProbabilitiesInfectedSymptoms[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = 1 - severePerInfectedSymptoms; + model.parameters.get>() = + StartingProbabilitiesInfectedSymptoms; + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedSymptomsToInfectedSevere().get_default( + (size_t)(LctState::get_num_subcompartments() / 2.), 5.8); + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedSymptomsToRecovered().get_default( + (size_t)(LctState::get_num_subcompartments() / 2.), 5.8); + // InfectedSevere. + Eigen::VectorX StartingProbabilitiesInfectedSevere = Eigen::VectorX::Zero( + (Eigen::Index)LctState::get_num_subcompartments()); + StartingProbabilitiesInfectedSevere[0] = criticalPerSevere; + StartingProbabilitiesInfectedSevere[(Eigen::Index)( + LctState::get_num_subcompartments() / 3.)] = deathsPerSevere; + StartingProbabilitiesInfectedSevere[2 * (Eigen::Index)( + LctState::get_num_subcompartments() / + 3.)] = 1. - criticalPerSevere - deathsPerSevere; + model.parameters.get>() = + StartingProbabilitiesInfectedSevere; + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedSevereToInfectedCritical().get_default( + (size_t)(LctState::get_num_subcompartments() / 3.), 9.5); + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedSevereToDead().get_default( + (size_t)(LctState::get_num_subcompartments() / 3.), 9.5); + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedSevereToRecovered().get_default( + (size_t)(LctState::get_num_subcompartments() / 3.), 9.5); + // InfectedCritical. + Eigen::VectorX StartingProbabilitiesInfectedCritical = Eigen::VectorX::Zero( + (Eigen::Index)LctState::get_num_subcompartments()); + StartingProbabilitiesInfectedCritical[0] = deathsPerCritical; + StartingProbabilitiesInfectedCritical[(Eigen::Index)( + LctState::get_num_subcompartments() / 2.)] = 1 - deathsPerCritical; + model.parameters.get>() = + StartingProbabilitiesInfectedCritical; + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedCriticalToDead().get_default( + (size_t)(LctState::get_num_subcompartments() / 2.), 7.1); + model.parameters.get>() = + mio::glsecir::TransitionMatrixInfectedCriticalToRecovered().get_default( + (size_t)(LctState::get_num_subcompartments() / 2.), 7.1); + + // Define initial population distribution in infection states, one entry per subcompartment. + // We start with 1000 individuals in the first subcompartment of the InfectedSymtpoms compartment which refers to + // the strain of individuals that are transitioning from InfectedSymptoms to InfectedSevere. + const auto initial_pop = 1000.; + const auto dead_initial = 0.; + std::vector> initial_populations = {{0.}, {0., 0.}, {0., 0.}, {initial_pop, 0.}, + {0., 0., 0.}, {0., 0.}, {0.}, {dead_initial}}; + + // Transfer the initial values in initial_populations to the model. + std::vector flat_initial_populations; + for (auto&& vec : initial_populations) { + flat_initial_populations.insert(flat_initial_populations.end(), vec.begin(), vec.end()); + } + for (size_t i = 0; i < LctState::Count; i++) { + model.populations[mio::Index(i)] = flat_initial_populations[i]; + } + + // Simulate a small time step. + ScalarType t0 = 0, tmax = 2., dt = 0.1; + mio::TimeSeries result = mio::simulate(t0, tmax, dt, model); + + // With DeathsPerSevere=0.1, deaths from ISev must be > 0. + mio::TimeSeries population = model.calculate_compartments(result); + auto dead_final = population.get_last_value()[(size_t)mio::glsecir::InfectionState::Dead]; + EXPECT_GT(dead_final, dead_initial); + + // Now run again with DeathsPerSevere=0. + deathsPerSevere = 0.0; + StartingProbabilitiesInfectedSevere[(Eigen::Index)( + LctState::get_num_subcompartments() / 3.)] = deathsPerSevere; + StartingProbabilitiesInfectedSevere[2 * (Eigen::Index)( + LctState::get_num_subcompartments() / + 3.)] = 1. - criticalPerSevere - deathsPerSevere; + model.parameters.get>() = + StartingProbabilitiesInfectedSevere; + + mio::TimeSeries result_no_severe_deaths = mio::simulate(t0, tmax, dt, model); + mio::TimeSeries population_no_severe_deaths = model.calculate_compartments(result_no_severe_deaths); + auto dead_no_severe_deaths = + population_no_severe_deaths.get_last_value()[(size_t)mio::glsecir::InfectionState::Dead]; + + // With DeathsPerSevere=0 there should be no deaths from ISev, + // so dead compartment stays 0 (ICr is also 0 initially). + EXPECT_DOUBLE_EQ(dead_no_severe_deaths, dead_initial); + + // Population must be conserved in both runs. + ScalarType total = population.get_last_value().sum(); + ScalarType total_no_severe_deaths = population_no_severe_deaths.get_last_value().sum(); + EXPECT_NEAR(total, initial_pop, 1e-10); + EXPECT_NEAR(total_no_severe_deaths, initial_pop, 1e-10); + + // Recovery flow must be higher when DeathsPerSevere=0. + auto rec = population.get_last_value()[(size_t)mio::glsecir::InfectionState::Recovered]; + auto rec_no_severe_deaths = + population_no_severe_deaths.get_last_value()[(size_t)mio::glsecir::InfectionState::Recovered]; + EXPECT_GT(rec_no_severe_deaths, rec); +} diff --git a/cpp/tests/test_lct_parameters_io.cpp b/cpp/tests/test_lct_parameters_io.cpp index 1acc1c4ca4..f6d09f3a39 100644 --- a/cpp/tests/test_lct_parameters_io.cpp +++ b/cpp/tests/test_lct_parameters_io.cpp @@ -112,6 +112,7 @@ TEST(TestLCTParametersIo, ReadPopulationDataRKI) model.parameters.template get>()[0] = 0.2; model.parameters.template get>()[0] = 0.1; model.parameters.template get>()[0] = 0.3; + model.parameters.template get>()[0] = 0.0; model.parameters.template get>()[0] = 0.2; // Calculate initial value vector for subcompartments with RKI data. @@ -158,6 +159,7 @@ TEST(TestLCTParametersIo, ReadPopulationDataRKIAgeres) model.parameters.template get>()[i] = 0.2; model.parameters.template get>()[i] = 0.1; model.parameters.template get>()[i] = 0.3; + model.parameters.template get>()[i] = 0.0; model.parameters.template get>()[i] = 0.2; } @@ -235,6 +237,7 @@ TEST(TestLCTParametersIo, CheckScalingDIVI) model.parameters.template get>()[i] = 0.2; model.parameters.template get>()[i] = 0.1; model.parameters.template get>()[i] = 0.3; + model.parameters.template get>()[i] = 0.0; model.parameters.template get>()[i] = 0.2; } @@ -369,6 +372,7 @@ TEST(TestLCTParametersIo, ReadPopulationDataRKIFailure) model.parameters.template get>()[0] = 0.2; model.parameters.template get>()[0] = 0.1; model.parameters.template get>()[0] = 0.3; + model.parameters.template get>()[0] = 0.0; model.parameters.template get>()[0] = 0.2; // Deactivate temporarily log output for next tests. @@ -435,6 +439,7 @@ TEST(TestLCTParametersIo, ReadPopulationDataRKIFailureAgeres) model.parameters.template get>()[i] = 0.2; model.parameters.template get>()[i] = 0.1; model.parameters.template get>()[i] = 0.3; + model.parameters.template get>()[i] = 0.0; model.parameters.template get>()[i] = 0.2; } diff --git a/cpp/tests/test_lct_secir.cpp b/cpp/tests/test_lct_secir.cpp index e7cbc428ee..4728450116 100755 --- a/cpp/tests/test_lct_secir.cpp +++ b/cpp/tests/test_lct_secir.cpp @@ -361,8 +361,8 @@ TEST(TestLCTSecir, testEvalRightHandSideGroups) } } -/* Test if the function get_derivatives() is working for more than one group. - * The parameters and LctStates are chosen, such that the sum of the groups should produce the +/* Test if the function get_derivatives() is working for more than one group. + * The parameters and LctStates are chosen, such that the sum of the groups should produce the * same output as the function with one group. This is tested. */ TEST(TestLCTSecir, testEvalRightHandSideThreeGroupsEqual) { @@ -521,7 +521,7 @@ TEST_F(ModelTestLCTSecir, compareWithPreviousRun) } /* Test compares a simulation with the result of a previous run stored in a .csv file. - * Here, three groups with the same parametrization as the one used for the previous result is used. + * Here, three groups with the same parametrization as the one used for the previous result is used. * It is tested, if the sum of the results of the groups is equal to the result with one group. */ TEST(TestLCTSecir, compareWithPreviousRunThreeGroups) @@ -652,6 +652,7 @@ TEST(TestLCTSecir, testConstraintsParameters) parameters_lct.get>()[0] = 0.1; parameters_lct.get>()[0] = 0.1; parameters_lct.get>()[0] = 0.1; + parameters_lct.get>()[0] = 0.1; parameters_lct.get>()[0] = 0.1; // Check improper TimeExposed. @@ -719,6 +720,17 @@ TEST(TestLCTSecir, testConstraintsParameters) EXPECT_TRUE(constraint_check); parameters_lct.get>()[0] = 0.1; + // Check DeathsPerSevere. + parameters_lct.get>()[0] = -1; + constraint_check = parameters_lct.check_constraints(); + EXPECT_TRUE(constraint_check); + + // Check CriticalPerSevere + DeathsPerSevere. + parameters_lct.get>()[0] = 1.0; + constraint_check = parameters_lct.check_constraints(); + EXPECT_TRUE(constraint_check); + parameters_lct.get>()[0] = 0.1; + // Check DeathsPerCritical. parameters_lct.get>()[0] = -1; constraint_check = parameters_lct.check_constraints(); @@ -785,3 +797,80 @@ TEST(TestLCTSecir, testConstraintsModel) constraint_check = model.check_constraints(); EXPECT_FALSE(constraint_check); } + +TEST(TestLCTSecir, deathsPerSevere_flows) +{ + // Test that DeathsPerSevere causes ISev->Dead flow. + // CriticalPerSevere=0 blocks the ISev->ICr->Dead path entirely, so + // DeathsPerCritical has no influence and is left at its default values. + + using InfState = mio::lsecir::InfectionState; + using LctState = mio::LctInfectionState; + using Model = mio::lsecir::Model; + + // Initialize a model with one group. + Model model; + auto& params = model.parameters; + + // params.get>()[0] = 5.0; + // params.get>()[0] = 5.0; + params.get>()[0] = 0.0; // block ISev->ICr->Dead path + params.get>()[0] = 0.1; + + // Define initial population distribution in infection states, one entry per subcompartment. + // We start with 500 individuals in each InfectedSevere subcompartment, no individuals in other compartments. + const auto severe_initial = 1000.; + const auto dead_initial = 0.; + + std::vector> initial_populations = { + {0.}, + {0., 0.}, + {0., 0.}, + {0., 0.}, + {severe_initial / LctState::template get_num_subcompartments(), + severe_initial / LctState::template get_num_subcompartments()}, + {0., 0.}, + {0.}, + {dead_initial}}; + + // Transfer the initial values in initial_populations to the model. + std::vector flat_initial_populations; + for (auto&& vec : initial_populations) { + flat_initial_populations.insert(flat_initial_populations.end(), vec.begin(), vec.end()); + } + for (size_t i = 0; i < LctState::Count; i++) { + model.populations[i] = flat_initial_populations[i]; + } + + // Simulate a small time step. + ScalarType t0 = 0, tmax = 1., dt = 0.1; + mio::TimeSeries result = mio::simulate(t0, tmax, dt, model); + + // With DeathsPerSevere=0.1, deaths from ISev must be > 0. + mio::TimeSeries population = model.calculate_compartments(result); + auto dead_final = population.get_last_value()[(size_t)mio::lsecir::InfectionState::Dead]; + EXPECT_GT(dead_final, dead_initial); + + // Now run again with DeathsPerSevere=0. + params.get>()[0] = 0.0; + mio::TimeSeries result_no_severe_deaths = mio::simulate(t0, tmax, dt, model); + mio::TimeSeries population_no_severe_deaths = model.calculate_compartments(result_no_severe_deaths); + auto dead_no_severe_deaths = + population_no_severe_deaths.get_last_value()[(size_t)mio::lsecir::InfectionState::Dead]; + + // With DeathsPerSevere=0 there should be no deaths from ISev, + // so dead compartment stays 0 (ICr is also 0 initially). + EXPECT_DOUBLE_EQ(dead_no_severe_deaths, dead_initial); + + // Population must be conserved in both runs. + ScalarType total = population.get_last_value().sum(); + ScalarType total_no_severe_deaths = population_no_severe_deaths.get_last_value().sum(); + EXPECT_NEAR(total, severe_initial, 1e-10); + EXPECT_NEAR(total_no_severe_deaths, severe_initial, 1e-10); + + // Recovery flow must be higher when DeathsPerSevere=0. + auto rec = population.get_last_value()[(size_t)mio::lsecir::InfectionState::Recovered]; + auto rec_no_severe_deaths = + population_no_severe_deaths.get_last_value()[(size_t)mio::lsecir::InfectionState::Recovered]; + EXPECT_GT(rec_no_severe_deaths, rec); +} diff --git a/cpp/tests/test_lct_secir_2_diseases.cpp b/cpp/tests/test_lct_secir_2_diseases.cpp index 38cc691fdb..5bbf5da1cb 100644 --- a/cpp/tests/test_lct_secir_2_diseases.cpp +++ b/cpp/tests/test_lct_secir_2_diseases.cpp @@ -750,12 +750,14 @@ TEST(TestLCTSecir2d, testConstraintsParameters) parameters_lct2d.get>()[0] = 0.1; parameters_lct2d.get>()[0] = 0.1; parameters_lct2d.get>()[0] = 0.1; + parameters_lct2d.get>()[0] = 0.1; parameters_lct2d.get>()[0] = 0.1; parameters_lct2d.get>()[0] = 1; parameters_lct2d.get>()[0] = 1; parameters_lct2d.get>()[0] = 0.1; parameters_lct2d.get>()[0] = 0.1; parameters_lct2d.get>()[0] = 0.1; + parameters_lct2d.get>()[0] = 0.1; parameters_lct2d.get>()[0] = 0.1; // Check improper TimeExposed. @@ -878,6 +880,23 @@ TEST(TestLCTSecir2d, testConstraintsParameters) EXPECT_TRUE(constraint_check); parameters_lct2d.get>()[0] = 0.1; + // Check DeathsPerSevere and CriticalPerSevere + DeathsPerSevere. + parameters_lct2d.get>()[0] = -1; + constraint_check = parameters_lct2d.check_constraints(); + EXPECT_TRUE(constraint_check); + parameters_lct2d.get>()[0] = 1; + constraint_check = parameters_lct2d.check_constraints(); + EXPECT_TRUE(constraint_check); + parameters_lct2d.get>()[0] = 0.1; + + parameters_lct2d.get>()[0] = -1; + constraint_check = parameters_lct2d.check_constraints(); + EXPECT_TRUE(constraint_check); + parameters_lct2d.get>()[0] = 1; + constraint_check = parameters_lct2d.check_constraints(); + EXPECT_TRUE(constraint_check); + parameters_lct2d.get>()[0] = 0.1; + // Check DeathsPerCritical. parameters_lct2d.get>()[0] = -1; constraint_check = parameters_lct2d.check_constraints(); @@ -973,3 +992,115 @@ TEST(TestLCTSecir2d, testConstraintsModel) constraint_check = model.check_constraints(); EXPECT_FALSE(constraint_check); } + +TEST(TestLCTSecir2d, deathsPerSevere_flows) +{ + // Test that DeathsPerSevere causes ISev->Dead flow. + // CriticalPerSevere=0 blocks the ISev->ICr->Dead path entirely, so + // DeathsPerCritical has no influence and is left at its default values. + + using InfState = mio::lsecir2d::InfectionState; + using LctState = mio::LctInfectionState; + using Model = mio::lsecir2d::Model; + + // Initialize a model with one group. + Model model; + auto& params = model.parameters; + + params.get>()[0] = 0.0; // block ISev->ICr->Dead path for disease a + params.get>()[0] = 0.1; + params.get>()[0] = 0.0; // block ISev->ICr->Dead path for disease b + params.get>()[0] = 0.1; + + // Define initial population distribution in infection states, one entry per subcompartment. + // We start with 1000 individuals in each InfectedSevere compartment, no indiviuals in other compartments. + const auto severe_initial_1a = 1000., severe_initial_1b = 1000., severe_initial_2a = 1000., + severe_initial_2b = 1000.; + const auto severe_initial_total = severe_initial_1a + severe_initial_1b + severe_initial_2a + severe_initial_2b; + const auto dead_initial_a = 0., dead_initial_b = 0.; + + std::vector> initial_populations = { + {0.}, // S + {0.}, // E 1a + {0.}, // INS 1a + {0.}, // ISy 1a + {severe_initial_1a}, // ISev 1a + {0.}, // ICri 1a + {0.}, // R 1a + {dead_initial_a}, // D a + {0.}, // E 2a + {0.}, // INS 2a + {0.}, // ISy 2a + {severe_initial_2a}, // ISev 2a + {0.}, // ICri 2a + {0.}, // E 1b + {0.}, // INS 1b + {0.}, // ISy 1b + {severe_initial_1b}, // ISev 1b + {0.}, // ICri 1b + {0.}, // R 1b + {dead_initial_b}, // D b + {0.}, // E 2b + {0.}, // INS 2b + {0.}, // ISy 2b + {severe_initial_2b}, // ISev 2b + {0.}, // ICri 2b + {0.}, // R 2ab + }; + + // Transfer the initial values in initial_populations to the model. + std::vector flat_initial_populations; + for (auto&& vec : initial_populations) { + flat_initial_populations.insert(flat_initial_populations.end(), vec.begin(), vec.end()); + } + for (size_t i = 0; i < LctState::Count; i++) { + model.populations[i] = flat_initial_populations[i]; + } + + // Simulate a small time step. + ScalarType t0 = 0, tmax = 2., dt = 0.1; + mio::TimeSeries result = mio::simulate(t0, tmax, dt, model); + + // With DeathsPerSevere=0.1, deaths from ISev must be > 0 for both diseases. + mio::TimeSeries population = model.calculate_compartments(result); + auto dead_final_a = population.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Dead_a]; + auto dead_final_b = population.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Dead_b]; + EXPECT_GT(dead_final_a, dead_initial_a); + EXPECT_GT(dead_final_b, dead_initial_b); + + // Now run again with DeathsPerSevere=0. + params.get>()[0] = 0.0; + params.get>()[0] = 0.0; + mio::TimeSeries result_no_severe_deaths = mio::simulate(t0, tmax, dt, model); + mio::TimeSeries population_no_severe_deaths = model.calculate_compartments(result_no_severe_deaths); + auto dead_no_severe_deaths_a = + population_no_severe_deaths.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Dead_a]; + auto dead_no_severe_deaths_b = + population_no_severe_deaths.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Dead_b]; + + // With DeathsPerSevere=0 there should be no deaths from ISev, + // so dead compartment stays 0 (ICr is also 0 initially) for both diseases. + EXPECT_DOUBLE_EQ(dead_no_severe_deaths_a, dead_initial_a); + EXPECT_DOUBLE_EQ(dead_no_severe_deaths_b, dead_initial_b); + + // Population must be conserved in both runs. + ScalarType total = population.get_last_value().sum(); + ScalarType total_no_severe_deaths = population_no_severe_deaths.get_last_value().sum(); + EXPECT_NEAR(total, severe_initial_total, 1e-10); + EXPECT_NEAR(total_no_severe_deaths, severe_initial_total, 1e-10); + + // Recovery flow must be higher when DeathsPerSevere=0. + auto rec_1a = population.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Recovered_1a]; + auto rec_1b = population.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Recovered_1b]; + auto rec_2ab = population.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Recovered_2ab]; + auto rec_no_severe_deaths_1a = + population_no_severe_deaths.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Recovered_1a]; + auto rec_no_severe_deaths_1b = + population_no_severe_deaths.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Recovered_1b]; + auto rec_no_severe_deaths_2ab = + population_no_severe_deaths.get_last_value()[(size_t)mio::lsecir2d::InfectionState::Recovered_2ab]; + EXPECT_GT(rec_no_severe_deaths_1a, rec_1a); + EXPECT_GT(rec_no_severe_deaths_1b, rec_1b); + EXPECT_GT(rec_no_severe_deaths_2ab, rec_2ab); +} diff --git a/docs/source/cpp/models/glsecir.rst b/docs/source/cpp/models/glsecir.rst index 7a18be9e6b..1efd3adfbc 100644 --- a/docs/source/cpp/models/glsecir.rst +++ b/docs/source/cpp/models/glsecir.rst @@ -68,7 +68,7 @@ Below is an overview of the model variables: The model equations are given below. For a simpler description let :math:`\mathcal{Z}=\{E,I_{NS},I_{Sy},I_{Sev},I_{Cr}\}` be the set of the compartments that can be divided into subcompartments. -.. image:: https://github.com/SciCompMod/memilio/assets/70579874/e1da5e1d-e719-4c16-9f14-45374be7c353 +.. image:: http://martinkuehn.eu/research/images/glct_equations.png :alt: equations Note that the bold notation :math:`\mathbf{z}(t)` for :math:`z \in \mathcal{Z}` stands for a vector. If several transitions are possible from a compartment, the vector is split in order to be able to select the stay times until the transitions individually. For example, the order @@ -80,9 +80,21 @@ Note that the bold notation :math:`\mathbf{z}(t)` for :math:`z \in \mathcal{Z}` \mathbf{I_{\text{NS}}^{\text{R}}}(t) \end{bmatrix} -is used. Similar holds true for the other compartments :math:`z \in \mathcal{Z}`. +is used. Similar holds true for the other compartments :math:`z \in \mathcal{Z}`. In particular, we have three transitions +originating from the InfectedSevere compartment and the vector is given by -Implicitly, the matrices :math:`\mathbf{A_{z}^{*}}` for one :math:`z \in \mathcal{Z}` are a block of a matrix :math:`\mathbf{A_{z}}` corresponding to the whole vector :math:`\mathbf{z}(t)`. As we have no transitions in between the strains defined for different transition probabilities, we would have many zeros in the matrix. The matrix can be defined as +.. math:: + + \mathbf{I_{\text{NS}}}(t) = \begin{bmatrix} + \mathbf{I_{\text{NS}}^{\text{Sy}}}(t) \\ + \mathbf{I_{\text{NS}}^{\text{D}}}(t) \\ + \mathbf{I_{\text{NS}}^{\text{R}}}(t) + \end{bmatrix}. + +Implicitly, the matrices :math:`\mathbf{A_{z}^{*}}` for one :math:`z \in \mathcal{Z}` are a block of a matrix +:math:`\mathbf{A_{z}}` corresponding to the whole vector :math:`\mathbf{z}(t)`. As we have no transitions in between +the strains defined for different transition probabilities, we would have many zeros in the matrix. The matrix can be +defined as .. math:: @@ -92,7 +104,11 @@ Implicitly, the matrices :math:`\mathbf{A_{z}^{*}}` for one :math:`z \in \mathca \mathbf{0} & \mathbf{A_{z}^{*_2}} \end{bmatrix}, -where :math:`{*}_{1}` is the compartment of the first transition, e.g., :math:`I_{\text{Sy}}` for :math:`z=I_{\text{NS}}` and :math:`*_{2}` the compartment of the second possible transition, e.g., :math:`R`. Therefore, we just store the non-zero blocks of the matrix. Using these parameters, the phase-type distribution that defines the stay time in compartment :math:`z \in \mathcal{Z}` has the probability density function +where :math:`{*}_{1}` is the compartment of the first transition, e.g., :math:`I_{\text{Sy}}` for +:math:`z=I_{\text{NS}}` and :math:`*_{2}` the compartment of the second possible transition, e.g., :math:`R`. +Therefore, we just store the non-zero blocks of the matrix. The number of non-zero blocks corresponds to the number of +strains originating from the considered compartment. Using these parameters, the phase-type distribution that defines +the stay time in compartment :math:`z \in \mathcal{Z}` has the probability density function .. math:: @@ -145,6 +161,7 @@ We continue by defining some epidemiological parameters needed throughout the mo const ScalarType recoveredPerInfectedNoSymptoms = 0.09; const ScalarType severePerInfectedSymptoms = 0.2; const ScalarType criticalPerSevere = 0.25; + const ScalarType deathsPerSevere = 0.; const ScalarType deathsPerCritical = 0.3; Now, we define the initial values with the distribution of the population into subcompartments. Note that this method of defining the initial values using a vector of vectors is not necessary, but should show how the entries of the initial value vector relate to the defined template parameters of the model or the number of subcompartments. It is also possible to define the initial values directly. diff --git a/docs/source/cpp/models/lsecir.rst b/docs/source/cpp/models/lsecir.rst index 3dafb697b4..afca3cf21c 100644 --- a/docs/source/cpp/models/lsecir.rst +++ b/docs/source/cpp/models/lsecir.rst @@ -13,7 +13,7 @@ The model is particularly suited for pathogens with pre- or asymptomatic infecti Below is a visualization of the infection states and transitions without a stratification according to sociodemographic groups. -.. image:: https://github.com/SciCompMod/memilio/assets/70579874/6a5d5a95-20f9-4176-8894-c091bd48bfb7 +.. image:: http://martinkuehn.eu/research/images/lct.png :alt: tikzLCTSECIR @@ -118,6 +118,9 @@ The model implements the following parameters: * - :math:`\mu_{I_{Sev}}^{I_{Cr}}` - ``CriticalPerSevere`` - Probability of transition from compartment InfectedSevere to InfectedCritical. + * - :math:`\mu_{I_{Sev}}^{D}` + - ``DeathsPerSevere`` + - Probability of dying when in compartment InfectedSevere. * - :math:`\mu_{I_{Cr}}^{D}` - ``DeathsPerCritical`` - Probability of dying when in compartment InfectedCritical. diff --git a/docs/source/cpp/models/lsecir2d.rst b/docs/source/cpp/models/lsecir2d.rst index 02b7a5a6a6..d84bf495f3 100644 --- a/docs/source/cpp/models/lsecir2d.rst +++ b/docs/source/cpp/models/lsecir2d.rst @@ -259,6 +259,12 @@ The model implements the following parameters: * - :math:`\mu_{I_{Sev},b}^{I_{Cr}}` - ``CriticalPerSevere_b`` - Probability of transition from compartment InfectedSevere_1b to InfectedCritical_1b or from InfectedSevere_2b to InfectedCritical_2b. + * - :math:`\mu_{I_{Sev},a}^{D}` + - ``DeathsPerSevere_a`` + - Probability of dying when in compartment InfectedSevere_1a or InfectedSevere_2a. + * - :math:`\mu_{I_{Sev},b}^{D}` + - ``DeathsPerSevere_b`` + - Probability of dying when in compartment InfectedSevere_1b or InfectedSevere_2b. * - :math:`\mu_{I_{Cr},a}^{D}` - ``DeathsPerCritical_a`` - Probability of dying when in compartment InfectedCritical_1a or InfectedCritical_2a.