From 19980a77cc65572c27c892cd9c2d66e928de4cc8 Mon Sep 17 00:00:00 2001 From: "George G. Vega Yon" Date: Mon, 18 Nov 2024 17:40:18 -0700 Subject: [PATCH] Updating example notes --- epiworld.hpp | 10 +- examples/14-community-hosp/README.md | 99 ++++++++++++++++ examples/14-community-hosp/main.cpp | 110 +++++++++++++----- .../epiworld/agent-meat-virus-sampling.hpp | 6 +- include/epiworld/misc.hpp | 2 +- include/epiworld/model-meat.hpp | 2 +- 6 files changed, 193 insertions(+), 36 deletions(-) create mode 100644 examples/14-community-hosp/README.md diff --git a/epiworld.hpp b/epiworld.hpp index f6a3022b5..01073f310 100644 --- a/epiworld.hpp +++ b/epiworld.hpp @@ -595,7 +595,7 @@ inline bool IN(const Ta & a, const std::vector< Ta > & b) noexcept * @return int If -1 then it means that none got sampled, otherwise the index * of the entry that got drawn. */ -template +template inline int roulette( const std::vector< TDbl > & probs, Model * m @@ -7661,7 +7661,7 @@ template inline epiworld_double & Model::operator()(std::string pname) { if (parameters.find(pname) == parameters.end()) - throw std::range_error("The parameter "+ pname + "is not in the model."); + throw std::range_error("The parameter '"+ pname + "' is not in the model."); return parameters[pname]; @@ -13078,7 +13078,7 @@ namespace sampler { * @return Virus* of the selected virus. If none selected (or none * available,) returns a nullptr; */ -template +template inline std::function*,Model*)> make_update_susceptible( std::vector< epiworld_fast_uint > exclude = {} ) @@ -13237,7 +13237,7 @@ inline std::function*,Model*)> make_update_susceptible( * @return Virus* of the selected virus. If none selected (or none * available,) returns a nullptr; */ -template +template inline std::function*(Agent*,Model*)> make_sample_virus_neighbors( std::vector< epiworld_fast_uint > exclude = {} ) @@ -13405,7 +13405,7 @@ inline std::function*(Agent*,Model*)> make_sample_virus_ * @return Virus* of the selected virus. If none selected (or none * available,) returns a nullptr; */ -template +template inline Virus * sample_virus_single(Agent * p, Model * m) { diff --git a/examples/14-community-hosp/README.md b/examples/14-community-hosp/README.md new file mode 100644 index 000000000..0f4bae93d --- /dev/null +++ b/examples/14-community-hosp/README.md @@ -0,0 +1,99 @@ +# Community-hospital model + +In this model, we have three states: + +- Susceptible individuals (in the community), +- Infected individuals (in the community), and +- Infected hospitalized individuals. + +Susceptible individuals may become hospitalized or not (so they have two possible transitions), and infected individuals may recover or (if hospitalized) stay infected but be discharged (so they go back as infected but in the community). + +The model has the following parameters: + +- Prob hospitalization: 0.1 +- Prob recovery: 0.33 +- Discharge infected: 0.1 + +Here is an example of the run + +```bash +root ➜ /workspaces/epiworld/examples/14-community-hosp (example-karim) $ make +g++ -std=c++17 -O3 -g main.cpp -o main.o +root ➜ /workspaces/epiworld/examples/14-community-hosp (example-karim) $ ./main.o +_________________________________________________________________________ +Running the model... +||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| done. + done. +________________________________________________________________________________ +________________________________________________________________________________ +SIMULATION STUDY + +Name of the model : (none) +Population size : 1000 +Agents' data : (none) +Number of entities : 0 +Days (duration) : 100 (of 100) +Number of viruses : 1 +Last run elapsed t : 4.00ms +Last run speed : 20.81 million agents x day / second +Rewiring : off + +Global events: + (none) + +Virus(es): + - MRSA + +Tool(s): + (none) + +Model parameters: + - Discharge infected : 0.1000 + - Prob hospitalization : 0.1000 + - Prob recovery : 0.3300 + +Distribution of the population at time 100: + - (0) Susceptible : 990 -> 937 + - (1) Infected : 10 -> 49 + - (2) Infected (hospitalized) : 0 -> 14 + +Transition Probabilities: + - Susceptible 0.98 0.02 0.00 + - Infected 0.32 0.61 0.07 + - Infected (hospitalized) 0.35 0.07 0.58 +``` + +## Notes + +An few key observations from this example. + +1. **We have a sampling function that exclude population**. These two functions are used in the update functions so, when susceptible (in the community) sample, the sampling excludes individuals who are hospitalized. Likewise, hospitalized + + + ```cpp + // A sampler that excludes infected from the hospital + auto sampler_suscept = sampler::make_sample_virus_neighbors<>( + {States::Infected_Hospitalized} + ); + ``` + +2. **All update functions are built from scratch**. For instance, susceptible individuals are updated according to the following function: + + ```cpp + inline void update_susceptible(Agent * p, Model * m) + { + + auto virus = sampler_suscept(p, m); + if (virus != nullptr) + { + if (m->par("Prob hospitalization") > m->runif()) + p->set_virus(*virus, m, States::Infected_Hospitalized); + else + p->set_virus(*virus, m, States::Infected); + } + + + return; + + } + ``` diff --git a/examples/14-community-hosp/main.cpp b/examples/14-community-hosp/main.cpp index 7a020551b..de3bc1369 100644 --- a/examples/14-community-hosp/main.cpp +++ b/examples/14-community-hosp/main.cpp @@ -9,34 +9,88 @@ using namespace epiworld; -template -inline void my_suscept( - Agent * p, - Model * m - ) -{ +enum States : size_t { + Susceptible = 0u, + Infected, + Infected_Hospitalized +}; + +// A sampler that excludes infected from the hospital +auto sampler_suscept = sampler::make_sample_virus_neighbors<>( + {States::Infected_Hospitalized} + ); + +/** + * - Susceptibles only live in the community. + * - Infection from individuals in the community only. + * - Once they become infected, they may be hospitalized or not. + */ - // Counting how many entities the agent has - auto n_entities = p->get_n_entities(); +inline void update_susceptible(Agent * p, Model * m) +{ - if (n_entities == 0) + auto virus = sampler_suscept(p, m); + if (virus != nullptr) { - // You can move agents to an entity - // p->add_entity(m->get_entity("Hospital")); - return; - + if (m->par("Prob hospitalization") > m->runif()) + p->set_virus(*virus, m, States::Infected_Hospitalized); + else + p->set_virus(*virus, m, States::Infected); } - Virus * virus = sampler::sample_virus_single(p, m); - // You can remove them from the entity - // p->rm_entity(m->get_entity("Hospital")); - - if (virus == nullptr) - return; + return; + +} + +/** + * Infected individuals may: + * + * - Stay the same + * - Recover + * - Be hospitalized + * + * Notice that the roulette makes the probabilities to sum to 1. + */ +inline void update_infected(Agent * p, Model * m) +{ + + // Vector of probabilities + std::vector< epiworld_double > probs = { + m->par("Prob hospitalization"), + m->par("Prob recovery") + }; - p->set_virus(*virus, m); + // Sampling: + // - (-1) Nothing happens + // - (0) Hospitalization + // - (1) Recovery + int res = roulette<>(probs, m); + if (res == 0) + p->change_state(m, States::Infected_Hospitalized); + else if (res == 1) + p->rm_virus(m, States::Susceptible); + + return; + +} + +/** + * Infected individuals who are hospitalized may: + * - Stay infected. + * - Recover (and then be discharged) + * - Stay the same and be discharged. + */ +inline void update_infected_hospitalized(Agent * p, Model * m) +{ + + if (m->par("Prob recovery") > m->runif()) { + p->rm_virus(m, States::Susceptible); + } else if (m->par("Discharge infected") > m->runif()) { + p->change_state(m, States::Infected); + } + return; } @@ -47,12 +101,12 @@ int main() { Model<> model; // model.add_state("Susceptible", default_update_susceptible<>); // State 0 - model.add_state("Susceptible", my_suscept<>); // State 0 - model.add_state("Infected", default_update_exposed<>); // State 1 - - Entity<> hospital("Hospital"); - hospital.set_distribution(distribute_entity_randomly<>(.5, true, true)); - model.add_entity(hospital); + model.add_state("Susceptible", update_susceptible); // State 0 + model.add_state("Infected", update_infected); // State 1 + model.add_state( + "Infected (hospitalized)", + update_infected_hospitalized + ); // State 2 // Adding a new virus Virus<> mrsa("MRSA"); @@ -66,6 +120,10 @@ int main() { // Add a population model.agents_smallworld(1000, 4, 0.1, false); + model.add_param(0.1, "Prob hospitalization"); + model.add_param(0.33, "Prob recovery"); + model.add_param(0.1, "Discharge infected"); + // Adding a new population model.run(100, 1231); model.print(); diff --git a/include/epiworld/agent-meat-virus-sampling.hpp b/include/epiworld/agent-meat-virus-sampling.hpp index d72899377..bfe668eeb 100644 --- a/include/epiworld/agent-meat-virus-sampling.hpp +++ b/include/epiworld/agent-meat-virus-sampling.hpp @@ -20,7 +20,7 @@ namespace sampler { * @return Virus* of the selected virus. If none selected (or none * available,) returns a nullptr; */ -template +template inline std::function*,Model*)> make_update_susceptible( std::vector< epiworld_fast_uint > exclude = {} ) @@ -179,7 +179,7 @@ inline std::function*,Model*)> make_update_susceptible( * @return Virus* of the selected virus. If none selected (or none * available,) returns a nullptr; */ -template +template inline std::function*(Agent*,Model*)> make_sample_virus_neighbors( std::vector< epiworld_fast_uint > exclude = {} ) @@ -347,7 +347,7 @@ inline std::function*(Agent*,Model*)> make_sample_virus_ * @return Virus* of the selected virus. If none selected (or none * available,) returns a nullptr; */ -template +template inline Virus * sample_virus_single(Agent * p, Model * m) { diff --git a/include/epiworld/misc.hpp b/include/epiworld/misc.hpp index dbeb73d14..6e0376b1c 100644 --- a/include/epiworld/misc.hpp +++ b/include/epiworld/misc.hpp @@ -124,7 +124,7 @@ inline bool IN(const Ta & a, const std::vector< Ta > & b) noexcept * @return int If -1 then it means that none got sampled, otherwise the index * of the entry that got drawn. */ -template +template inline int roulette( const std::vector< TDbl > & probs, Model * m diff --git a/include/epiworld/model-meat.hpp b/include/epiworld/model-meat.hpp index 5df3befb5..49160165d 100644 --- a/include/epiworld/model-meat.hpp +++ b/include/epiworld/model-meat.hpp @@ -714,7 +714,7 @@ template inline epiworld_double & Model::operator()(std::string pname) { if (parameters.find(pname) == parameters.end()) - throw std::range_error("The parameter "+ pname + "is not in the model."); + throw std::range_error("The parameter '"+ pname + "' is not in the model."); return parameters[pname];