Line data    Source code

       1             : /* 
       2             : * Copyright (C) 2020-2024 MEmilio
       3             : *
       4             : * Authors: Wadim Koslow, Daniel Abele, Martin J. Kühn
       5             : *
       6             : * Contact: Martin J. Kuehn <Martin.Kuehn@DLR.de>
       7             : *
       8             : * Licensed under the Apache License, Version 2.0 (the "License");
       9             : * you may not use this file except in compliance with the License.
      10             : * You may obtain a copy of the License at
      11             : *
      12             : *     http://www.apache.org/licenses/LICENSE-2.0
      13             : *
      14             : * Unless required by applicable law or agreed to in writing, software
      15             : * distributed under the License is distributed on an "AS IS" BASIS,
      16             : * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
      17             : * See the License for the specific language governing permissions and
      18             : * limitations under the License.
      19             : */
      20             : #ifndef ODESECIRVVS_PARAMETER_SPACE_H
      21             : #define ODESECIRVVS_PARAMETER_SPACE_H
      22             : 
      23             : #include "memilio/mobility/metapopulation_mobility_instant.h"
      24             : #include "memilio/utils/logging.h"
      25             : #include "ode_secirvvs/model.h"
      26             : #include "ode_secirvvs/infection_state.h"
      27             : 
      28             : #include <assert.h>
      29             : 
      30             : namespace mio
      31             : {
      32             : namespace osecirvvs
      33             : {
      34             : /**
      35             :      * draws a sample from the specified distributions for all parameters related to the demographics, e.g. population.
      36             :      * @tparam FP floating point type, e.g., double
      37             :      * @param[inout] model Model including contact patterns for alle age groups
      38             :      */
      39             : template <typename FP = double>
      40         108 : void draw_sample_demographics(Model<FP>& model)
      41             : {
      42         108 :     model.parameters.template get<ICUCapacity<FP>>().draw_sample();
      43         108 :     model.parameters.template get<TestAndTraceCapacity<FP>>().draw_sample();
      44             : 
      45         666 :     for (auto i = AgeGroup(0); i < model.parameters.get_num_groups(); i++) {
      46         558 :         double group_total = model.populations.get_group_total(i);
      47             : 
      48             :         //sample initial compartments (with exceptions)
      49       15624 :         for (auto inf_state = Index<InfectionState>(0); inf_state < InfectionState::Count; ++inf_state) {
      50       30132 :             if (inf_state != InfectionState::SusceptibleNaive && //not sampled, fixed after sampling everything else
      51       29574 :                 inf_state != InfectionState::DeadNaive && //not sampled, fixed from data
      52       73098 :                 inf_state != InfectionState::DeadPartialImmunity && //not sampled, fixed from data
      53       28458 :                 inf_state != InfectionState::DeadImprovedImmunity) { //not sampled, fixed from data
      54       12834 :                 model.populations[{i, inf_state}].draw_sample();
      55             :             }
      56             :         }
      57             : 
      58             :         //set susceptibles so the total number stays the same as before sampling.
      59             :         //if the new total without susceptibles is already bigger than the previous total
      60             :         //subtract the overflow from SusceptibleImprovedImmunity, susceptibles will then be approximately zero.
      61         558 :         model.populations[{i, InfectionState::SusceptibleNaive}] = 0;
      62         558 :         double diff                                              = model.populations.get_group_total(i) - group_total;
      63         558 :         if (diff > 0) {
      64           0 :             model.populations[{i, InfectionState::SusceptibleImprovedImmunity}] -= diff;
      65           0 :             if (model.populations[{i, InfectionState::SusceptibleImprovedImmunity}] < 0.0) {
      66           0 :                 log_error("Negative Compartment after sampling.");
      67             :             }
      68           0 :             assert(std::abs(group_total - model.populations.get_group_total(i)) < 1e-10 && "Sanity check.");
      69             :         }
      70         558 :         model.populations.template set_difference_from_group_total<AgeGroup>({i, InfectionState::SusceptibleNaive},
      71             :                                                                              group_total);
      72             :     }
      73         108 : }
      74             : 
      75             : /**
      76             :      * draws a sample from the specified distributions for all parameters related to the infection.
      77             :      * @tparam FP floating point type, e.g., double
      78             :      * @param[inout] model Model including contact patterns for alle age groups
      79             :      */
      80             : template <typename FP = double>
      81          99 : void draw_sample_infection(Model<FP>& model)
      82             : {
      83          99 :     model.parameters.template get<Seasonality<FP>>().draw_sample();
      84             : 
      85             :     //not age dependent
      86          99 :     model.parameters.template get<TimeExposed<FP>>()[AgeGroup(0)].draw_sample();
      87          99 :     model.parameters.template get<TimeInfectedNoSymptoms<FP>>()[AgeGroup(0)].draw_sample();
      88          99 :     model.parameters.template get<RelativeTransmissionNoSymptoms<FP>>()[AgeGroup(0)].draw_sample();
      89          99 :     model.parameters.template get<RiskOfInfectionFromSymptomatic<FP>>()[AgeGroup(0)].draw_sample();
      90          99 :     model.parameters.template get<MaxRiskOfInfectionFromSymptomatic<FP>>()[AgeGroup(0)].draw_sample();
      91             : 
      92          99 :     model.parameters.template get<ReducExposedPartialImmunity<FP>>()[AgeGroup(0)].draw_sample();
      93          99 :     model.parameters.template get<ReducExposedImprovedImmunity<FP>>()[AgeGroup(0)].draw_sample();
      94          99 :     model.parameters.template get<ReducInfectedSymptomsPartialImmunity<FP>>()[AgeGroup(0)].draw_sample();
      95          99 :     model.parameters.template get<ReducInfectedSymptomsImprovedImmunity<FP>>()[AgeGroup(0)].draw_sample();
      96          99 :     model.parameters.template get<ReducInfectedSevereCriticalDeadPartialImmunity<FP>>()[AgeGroup(0)].draw_sample();
      97          99 :     model.parameters.template get<ReducInfectedSevereCriticalDeadImprovedImmunity<FP>>()[AgeGroup(0)].draw_sample();
      98          99 :     model.parameters.template get<ReducTimeInfectedMild<FP>>()[AgeGroup(0)].draw_sample();
      99             : 
     100         603 :     for (auto i = AgeGroup(0); i < model.parameters.get_num_groups(); i++) {
     101             :         //not age dependent
     102         504 :         model.parameters.template get<TimeExposed<FP>>()[i] =
     103        1008 :             model.parameters.template get<TimeExposed<FP>>()[AgeGroup(0)];
     104         504 :         model.parameters.template get<TimeInfectedNoSymptoms<FP>>()[i] =
     105        1008 :             model.parameters.template get<TimeInfectedNoSymptoms<FP>>()[AgeGroup(0)];
     106         504 :         model.parameters.template get<RelativeTransmissionNoSymptoms<FP>>()[i] =
     107        1008 :             model.parameters.template get<RelativeTransmissionNoSymptoms<FP>>()[AgeGroup(0)];
     108         504 :         model.parameters.template get<RiskOfInfectionFromSymptomatic<FP>>()[i] =
     109        1008 :             model.parameters.template get<RiskOfInfectionFromSymptomatic<FP>>()[AgeGroup(0)];
     110         504 :         model.parameters.template get<MaxRiskOfInfectionFromSymptomatic<FP>>()[i] =
     111        1008 :             model.parameters.template get<MaxRiskOfInfectionFromSymptomatic<FP>>()[AgeGroup(0)];
     112             : 
     113         504 :         model.parameters.template get<ReducExposedPartialImmunity<FP>>()[i] =
     114        1008 :             model.parameters.template get<ReducExposedPartialImmunity<FP>>()[AgeGroup(0)];
     115         504 :         model.parameters.template get<ReducExposedImprovedImmunity<FP>>()[i] =
     116        1008 :             model.parameters.template get<ReducExposedImprovedImmunity<FP>>()[AgeGroup(0)];
     117         504 :         model.parameters.template get<ReducInfectedSymptomsPartialImmunity<FP>>()[i] =
     118        1008 :             model.parameters.template get<ReducInfectedSymptomsPartialImmunity<FP>>()[AgeGroup(0)];
     119         504 :         model.parameters.template get<ReducInfectedSymptomsImprovedImmunity<FP>>()[i] =
     120        1008 :             model.parameters.template get<ReducInfectedSymptomsImprovedImmunity<FP>>()[AgeGroup(0)];
     121         504 :         model.parameters.template get<ReducInfectedSevereCriticalDeadPartialImmunity<FP>>()[i] =
     122        1008 :             model.parameters.template get<ReducInfectedSevereCriticalDeadPartialImmunity<FP>>()[AgeGroup(0)];
     123         504 :         model.parameters.template get<ReducInfectedSevereCriticalDeadImprovedImmunity<FP>>()[i] =
     124        1008 :             model.parameters.template get<ReducInfectedSevereCriticalDeadImprovedImmunity<FP>>()[AgeGroup(0)];
     125         504 :         model.parameters.template get<ReducTimeInfectedMild<FP>>()[i] =
     126        1008 :             model.parameters.template get<ReducTimeInfectedMild<FP>>()[AgeGroup(0)];
     127             : 
     128             :         //age dependent
     129         504 :         model.parameters.template get<TimeInfectedSymptoms<FP>>()[i].draw_sample();
     130         504 :         model.parameters.template get<TimeInfectedSevere<FP>>()[i].draw_sample();
     131         504 :         model.parameters.template get<TimeInfectedCritical<FP>>()[i].draw_sample();
     132             : 
     133         504 :         model.parameters.template get<TransmissionProbabilityOnContact<FP>>()[i].draw_sample();
     134         504 :         model.parameters.template get<RecoveredPerInfectedNoSymptoms<FP>>()[i].draw_sample();
     135         504 :         model.parameters.template get<DeathsPerCritical<FP>>()[i].draw_sample();
     136         504 :         model.parameters.template get<SeverePerInfectedSymptoms<FP>>()[i].draw_sample();
     137         504 :         model.parameters.template get<CriticalPerSevere<FP>>()[i].draw_sample();
     138             :     }
     139          99 : }
     140             : 
     141             : /** Draws a sample from Model parameter distributions and stores sample values
     142             :     * as Parameters parameter values (cf. UncertainValue and Parameters classes)
     143             :     * @tparam FP floating point type, e.g., double
     144             :     * @param[inout] model Model including contact patterns for alle age groups
     145             :     */
     146             : template <typename FP = double>
     147             : void draw_sample(Model<FP>& model)
     148             : {
     149             :     draw_sample_infection(model);
     150             :     draw_sample_demographics(model);
     151             :     model.parameters.template get<ContactPatterns<FP>>().draw_sample();
     152             :     model.apply_constraints();
     153             : }
     154             : 
     155             : /**
     156             :     * Draws samples for each model node in a graph.
     157             :     * Some parameters are shared between nodes and only sampled once.
     158             :     * @tparam FP floating point type, e.g., double
     159             :     * @param graph Graph to be sampled.
     160             :     * @param variant_high If true, use high value for infectiousness of variant.
     161             :     * @return Graph with nodes and edges from the input graph sampled.
     162             :     */
     163             : template <typename FP = double>
     164          99 : Graph<Model<FP>, MobilityParameters<FP>> draw_sample(Graph<Model<FP>, MobilityParameters<FP>>& graph, bool variant_high)
     165             : {
     166          99 :     Graph<Model<FP>, MobilityParameters<FP>> sampled_graph;
     167             : 
     168             :     //sample global parameters
     169          99 :     auto& shared_params_model = graph.nodes()[0].property;
     170          99 :     draw_sample_infection(shared_params_model);
     171          99 :     auto& shared_contacts = shared_params_model.parameters.template get<ContactPatterns<FP>>();
     172          99 :     shared_contacts.draw_sample_dampings();
     173          99 :     auto& shared_dynamic_npis = shared_params_model.parameters.template get<DynamicNPIsInfectedSymptoms<FP>>();
     174          99 :     shared_dynamic_npis.draw_sample();
     175          99 :     auto& shared_dynamic_npis_delay = shared_params_model.parameters.template get<DynamicNPIsImplementationDelay<FP>>();
     176          99 :     shared_dynamic_npis_delay.draw_sample();
     177             : 
     178             :     double delta_fac;
     179          99 :     if (variant_high) {
     180          99 :         delta_fac = 1.6;
     181             :     }
     182             :     else {
     183           0 :         delta_fac = 1.4;
     184             :     }
     185             : 
     186             :     //infectiousness of virus variants is not sampled independently but depend on base infectiousness
     187         603 :     for (auto i = AgeGroup(0); i < shared_params_model.parameters.get_num_groups(); ++i) {
     188         504 :         shared_params_model.parameters.template get<InfectiousnessNewVariant<FP>>()[i] = delta_fac;
     189             :     }
     190             : 
     191         315 :     for (auto& params_node : graph.nodes()) {
     192         108 :         auto& node_model = params_node.property;
     193             : 
     194             :         //sample local parameters
     195         108 :         draw_sample_demographics(params_node.property);
     196             : 
     197             :         //copy global parameters
     198             :         //save demographic parameters so they aren't overwritten
     199         108 :         auto local_icu_capacity = node_model.parameters.template get<ICUCapacity<FP>>();
     200         108 :         auto local_tnt_capacity = node_model.parameters.template get<TestAndTraceCapacity<FP>>();
     201         108 :         auto local_holidays     = node_model.parameters.template get<ContactPatterns<FP>>().get_school_holidays();
     202         108 :         auto local_daily_v1     = node_model.parameters.template get<DailyFirstVaccination<FP>>();
     203         108 :         auto local_daily_v2     = node_model.parameters.template get<DailyFullVaccination<FP>>();
     204         108 :         node_model.parameters   = shared_params_model.parameters;
     205         108 :         node_model.parameters.template get<ICUCapacity<FP>>()                           = local_icu_capacity;
     206         108 :         node_model.parameters.template get<TestAndTraceCapacity<FP>>()                  = local_tnt_capacity;
     207         108 :         node_model.parameters.template get<ContactPatterns<FP>>().get_school_holidays() = local_holidays;
     208         108 :         node_model.parameters.template get<DailyFirstVaccination<FP>>()                 = local_daily_v1;
     209         108 :         node_model.parameters.template get<DailyFullVaccination<FP>>()                  = local_daily_v2;
     210             : 
     211         108 :         node_model.parameters.template get<ContactPatterns<FP>>().make_matrix();
     212         108 :         node_model.apply_constraints();
     213             : 
     214         108 :         sampled_graph.add_node(params_node.id, node_model);
     215             :     }
     216             : 
     217         117 :     for (auto& edge : graph.edges()) {
     218           9 :         auto edge_params = edge.property;
     219             :         //no dynamic NPIs
     220             :         //TODO: add switch to optionally enable dynamic NPIs to edges
     221           9 :         sampled_graph.add_edge(edge.start_node_idx, edge.end_node_idx, edge_params);
     222             :     }
     223             : 
     224         198 :     return sampled_graph;
     225          99 : }
     226             : } // namespace osecirvvs
     227             : } // namespace mio
     228             : 
     229             : #endif // ODESECIRVVS_PARAMETER_SPACE_H