phylopomp
Phylodynamics for POMPs
All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
seir.cc
Go to the documentation of this file.
1 // SEIR: Classical susceptible-exposed-infected-recovered model (C++)
2 #include "master.h"
3 #include "popul_proc.h"
4 #include "generics.h"
5 #include "internal.h"
6 
7 static int Exposed = 0;
8 static int Infectious = 1;
9 
11 typedef struct {
12  int S;
13  int E;
14  int I;
15  int R;
16  double N;
17 } seir_state_t;
18 
20 typedef struct {
21  double Beta;
22  double sigma;
23  double gamma;
24  double psi;
25  double omega;
26  int S0;
27  int E0;
28  int I0;
29  int R0;
31 
34 
35 template<>
36 std::string seir_proc_t::yaml (std::string tab) const {
37  std::string t = tab + " ";
38  std::string p = tab + "parameter:\n"
39  + YAML_PARAM(Beta)
40  + YAML_PARAM(sigma)
41  + YAML_PARAM(gamma)
42  + YAML_PARAM(psi)
43  + YAML_PARAM(omega)
44  + YAML_PARAM(S0)
45  + YAML_PARAM(E0)
46  + YAML_PARAM(I0)
47  + YAML_PARAM(R0);
48  std::string s = tab + "state:\n"
49  + YAML_STATE(S)
50  + YAML_STATE(E)
51  + YAML_STATE(I)
52  + YAML_STATE(R)
53  + YAML_STATE(N);
54  return p+s;
55 }
56 
57 template<>
58 void seir_proc_t::update_params (double *p, int n) {
59  int m = 0;
60  PARAM_SET(Beta);
63  PARAM_SET(psi);
65  if (m != n) err("wrong number of parameters!");
66 }
67 
68 template<>
69 void seir_proc_t::update_IVPs (double *p, int n) {
70  int m = 0;
71  PARAM_SET(S0);
72  PARAM_SET(E0);
73  PARAM_SET(I0);
74  PARAM_SET(R0);
75  if (m != n) err("wrong number of initial-value parameters!");
76 }
77 
78 template<>
79 double seir_proc_t::event_rates (double *rate, int n) const {
80  int m = 0;
81  double total = 0;
82  RATE_CALC(params.Beta * state.S * state.I / state.N);
83  RATE_CALC(params.sigma * state.E);
84  RATE_CALC(params.gamma * state.I);
85  RATE_CALC(params.psi * state.I);
86  RATE_CALC(params.omega * state.R);
87  if (m != n) err("wrong number of events!");
88  return total;
89 }
90 
91 template<>
92 void seir_genealogy_t::rinit (void) {
93  state.S = params.S0;
94  state.E = params.E0;
95  state.I = params.I0;
96  state.R = params.R0;
97  state.N = double(params.S0+params.E0+params.I0+params.R0);
98  graft(0,params.E0);
99  graft(1,params.I0);
100 }
101 
102 template<>
103 void seir_genealogy_t::jump (int event) {
104  switch (event) {
105  case 0:
106  state.S -= 1; state.E += 1; birth(Infectious,Exposed);
107  break;
108  case 1:
109  state.E -= 1; state.I += 1; migrate(Exposed,Infectious);
110  break;
111  case 2:
112  state.I -= 1; state.R += 1; death(Infectious);
113  break;
114  case 3:
116  break;
117  case 4:
118  state.R -= 1; state.S += 1;
119  break;
120  default: // #nocov
121  assert(0); // #nocov
122  break; // #nocov
123  }
124 }
125 
Encodes the master process.
Definition: master.h:22
void birth(name_t i=0, name_t j=0, int n=1)
n births into deme j with parent in deme i
Definition: master.h:142
void sample(name_t i=0, int n=1)
sample in deme i
Definition: master.h:166
void death(name_t i=0)
death in deme i
Definition: master.h:153
void jump(int e)
makes a jump
Definition: lbdp.cc:72
void rinit(void)
initialize the state
Definition: lbdp.cc:66
void graft(name_t i=0, int m=1)
new root in deme i
Definition: master.h:159
void migrate(name_t i=0, name_t j=0)
migration from deme i to deme j
Definition: master.h:185
Population process class.
Definition: popul_proc.h:20
void update_IVPs(double *, int)
set initial-value parameters
Definition: lbdp.cc:48
void update_params(double *, int)
set parameters
Definition: lbdp.cc:39
std::string yaml(std::string tab) const
machine/human readable info
Definition: lbdp.cc:26
parameters_t params
Definition: popul_proc.h:35
double event_rates(double *rate, int n) const
compute event rates
Definition: lbdp.cc:55
state_t state
Definition: popul_proc.h:34
#define GENERICS(X, TYPE)
Definition: generics.h:149
#define err(...)
Definition: internal.h:18
#define n
Definition: lbdp_pomp.c:8
#define psi
Definition: lbdp_pomp.c:6
#define YAML_PARAM(X)
Definition: popul_proc.h:180
#define RATE_CALC(X)
Definition: popul_proc.h:179
#define YAML_STATE(X)
Definition: popul_proc.h:181
#define PARAM_SET(X)
Definition: popul_proc.h:177
static int Exposed
Definition: seir.cc:7
static int Infectious
Definition: seir.cc:8
#define N
Definition: seirs_pomp.c:32
#define E
Definition: seirs_pomp.c:34
#define R0
Definition: seirs_pomp.c:31
#define gamma
Definition: seirs_pomp.c:25
#define R
Definition: seirs_pomp.c:36
#define S0
Definition: seirs_pomp.c:28
#define I
Definition: seirs_pomp.c:35
#define sigma
Definition: seirs_pomp.c:24
#define Beta
Definition: seirs_pomp.c:23
#define I0
Definition: seirs_pomp.c:30
#define E0
Definition: seirs_pomp.c:29
#define omega
Definition: seirs_pomp.c:27
#define S
Definition: seirs_pomp.c:33
SEIR process parameters.
Definition: seir.cc:20
double psi
Definition: seir.cc:24
double omega
Definition: seir.cc:25
double Beta
Definition: seir.cc:21
double gamma
Definition: seir.cc:23
double sigma
Definition: seir.cc:22
SEIR process state.
Definition: seir.cc:11
int S
Definition: seir.cc:12
double N
Definition: seir.cc:16
int I
Definition: seir.cc:14
int R
Definition: seir.cc:15
int E
Definition: seir.cc:13