phylopomp
Phylodynamics for POMPs
All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
si2r.cc
Go to the documentation of this file.
1// SI2R: Two-deme model of superspreading (C++)
2#include "master.h"
3#include "popul_proc.h"
4#include "generics.h"
5#include "internal.h"
6
7static int ordinary = 0;
8static int superspreader = 1;
9
11typedef struct {
12 int S;
13 int I1;
14 int I2;
15 int R;
16 double N;
18
20typedef struct {
21 double Beta;
22 double mu;
23 double gamma;
24 double omega;
25 double psi1;
26 double psi2;
27 double sigma12;
28 double sigma21;
29 int S0;
30 int I0;
31 int R0;
33
36
37template<>
38std::string si2r_proc_t::yaml (std::string tab) const {
39 std::string t = tab + " ";
40 std::string p = tab + "parameter:\n"
42 + YAML_PARAM(mu)
47 + YAML_PARAM(sigma12)
48 + YAML_PARAM(sigma21)
49 + YAML_PARAM(S0)
50 + YAML_PARAM(I0)
51 + YAML_PARAM(R0);
52 std::string s = tab + "state:\n"
53 + YAML_STATE(S)
54 + YAML_STATE(I1)
55 + YAML_STATE(I2)
56 + YAML_STATE(R)
57 + YAML_STATE(N);
58 return p+s;
59}
60
61template<>
62void si2r_proc_t::update_params (double *p, int n) {
63 int m = 0;
70 PARAM_SET(sigma12);
71 PARAM_SET(sigma21);
72 if (m != n) err("wrong number of parameters!");
73}
74
75template<>
76void si2r_proc_t::update_IVPs (double *p, int n) {
77 int m = 0;
81 if (m != n) err("wrong number of initial-value parameters!");
82}
83
84template<>
85double si2r_proc_t::event_rates (double *rate, int n) const {
86 int m = 0;
87 double total = 0;
88 RATE_CALC(params.Beta * state.S * state.I1 / state.N);
89 RATE_CALC(params.Beta * state.S * state.I2 / state.N);
90 RATE_CALC(params.gamma * state.I1);
91 RATE_CALC(params.gamma * state.I2);
92 RATE_CALC(params.psi1 * state.I1);
93 RATE_CALC(params.psi2 * state.I2);
94 RATE_CALC(params.sigma12 * state.I1);
95 RATE_CALC(params.sigma21 * state.I2);
96 RATE_CALC(params.omega * state.R);
97 if (m != n) err("wrong number of events!");
98 return total;
99}
100
101template<>
103 state.S = params.S0;
104 state.I1 = params.I0;
105 state.I2 = 0;
106 state.R = params.R0;
107 state.N = double(params.S0+params.I0+params.R0);
109}
110
111template<>
113 switch (event) {
114 case 0:
115 state.S -= 1; state.I1 += 1; birth(ordinary,ordinary);
116 break;
117 case 1:
118 {
119 int n = 1+int(rgeom(1.0/params.mu));
120 if (state.S >= n) {
121 state.S -= n; state.I1 += n;
123 } else {
125 state.I1 += state.S; state.S = 0;
126 }
127 }
128 break;
129 case 2:
130 state.I1 -= 1; state.R += 1; death(ordinary);
131 break;
132 case 3:
133 state.I2 -= 1; state.R += 1; death(superspreader);
134 break;
135 case 4:
137 break;
138 case 5:
140 break;
141 case 6:
142 state.I1 -= 1; state.I2 += 1; migrate(ordinary,superspreader);
143 break;
144 case 7:
145 state.I1 += 1; state.I2 -= 1; migrate(superspreader,ordinary);
146 break;
147 case 8:
148 state.R -= 1; state.S += 1;
149 break;
150 default: // #nocov
151 assert(0); // #nocov
152 break; // #nocov
153 }
154}
155
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 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
double event_rates(double *rate, int n) const
Definition si2r.cc:85
std::string yaml(std::string tab) const
Definition si2r.cc:38
#define GENERICS(X, TYPE)
Definition generics.h:149
#define err(...)
Definition internal.h:18
#define mu
Definition lbdp_pomp.c:5
#define n
Definition lbdp_pomp.c:8
#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
#define N
Definition seirs_pomp.c:32
#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 Beta
Definition seirs_pomp.c:23
#define I0
Definition seirs_pomp.c:30
#define omega
Definition seirs_pomp.c:27
#define S
Definition seirs_pomp.c:33
static int ordinary
Definition si2r.cc:7
master_t< si2r_proc_t, 2 > si2r_genealogy_t
Definition si2r.cc:35
static int superspreader
Definition si2r.cc:8
popul_proc_t< si2r_state_t, si2r_parameters_t, 9 > si2r_proc_t
Definition si2r.cc:34
SI2R process parameters.
Definition si2r.cc:20
double gamma
Definition si2r.cc:23
double psi1
Definition si2r.cc:25
double omega
Definition si2r.cc:24
double psi2
Definition si2r.cc:26
double Beta
Definition si2r.cc:21
double sigma21
Definition si2r.cc:28
double sigma12
Definition si2r.cc:27
SI2R process state.
Definition si2r.cc:11
double N
Definition si2r.cc:16
#define psi1
#define I2
#define psi2
#define I1