-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathREADME.Rmd
175 lines (109 loc) · 4.8 KB
/
README.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
[![CRAN_Status_Badge](http://www.r-pkg.org/badges/version/simulacr)](https://cran.r-project.org/package=simulacr)
# Welcome to the *simulacr* package!
This package is under development. Do not use it without contacting the authors first.
## Installing the package
To install the current stable, CRAN version of the package, type:
```{r install, eval = FALSE}
install.packages("simulacr")
```
To benefit from the latest features and bug fixes, install the development, *github* version of the package using:
```{r install2, eval = FALSE}
devtools::install_github("reconhub/simulacr")
```
Note that this requires the package *devtools* installed.
# What does it do?
*simulacr* implements outbreak simulation using branching processes. The main
features of the package include:
* `simulate_outbreak`: the main function, allowing to simulate outbreaks for
specified values of reproduction number, incubation period, duration of
infectiousness, and optionally reporting delays; it outputs a linelist stored
as a `data.frame` with the class `outbreak`, including information on
transmission chains; the output can be converted to `epicontacts` objects for
visualisation, and plotted using `plot(...)`
* `as_epicontacts`: a function to convert `outbreak` object to `epicontact`;
called implicitely when plotting `outbreak` objects
* `make_disc_gamma`: wrapper function to build discretised Gamma distributions
easily
* `draw_labels`: wrapper function to generate random alphanumeric labels
<br>
# Worked example
This brief example illustrates how we can simulate an outbreak, with different
ways of specifying input distributions.
## Specifying delay distributions
Distributions can be specified in different ways:
* as a discretised distribution stored as `distcrete` object, e.g. generated via
`make_disc_gamma`
* as a function with a single argument (the delay) computing probability mass
function (pmf) of different delays
* as a vector of $n$ positive numbers, taken to be the pmf of delays of
0:($n$-1)
We illustrate these different options below, with:
* a flat distribution of the incubation time 1-4 days (all other values have a
null probability)
* an infectious period defined as a discretised Gamma distribution stored as a
`distcrete` object
* a reporting delay following a Poisson distribution, specified as a function
```{r distributions}
library("simulacr")
incubation <- c(0, 1, 1, 1, 1) # numbers = unscaled PMF
infectious_period <- make_disc_gamma(10, 7) # distcrete object
reporting <- function(x) dpois(x, 5) # PMF function
```
## Simulating an outbreak
```{r simulations}
set.seed(1)
x <- simulate_outbreak(R = runif(100, 1, 3), # random values on [1;3]
dist_incubation= incubation,
dist_infectious_period = infectious_period,
dist_reporting = reporting)
```
The output is a `data.frame` with the class `outbreak`, which contains a
linelist of cases:
```{r view_x}
class(x)
dim(x)
head(x)
tail(x)
```
This object can be plotted using (this will open an interactive graph:
```{r plot_x}
plot(x)
```
For any work relying on transmission trees, it may be easiest to convert the
`outbreak` object to an `epicontacts`:
```{r convert_x}
net <- as_epicontacts(x)
net
```
## Simulating contacts
Contacts can be simulated and added to a simulated outbreak, using similar
procedures to the one used in `simulate_outbreak`, with a few differences:
* `n_contacts`: the numbers of contacts per index case is specified in the same
way as the reproduction number `R` in `simulate_outbreaks`
* the two distributions relate to the time between onset of the index case and
the beginning of secondar exposures (`dist_time_to_contact`) and the duration
of the exposure window (`duration`) after the first day of contact
```{r sim_contacts}
## exposure starts 0-2 days post onset
time_to_contact = c(1, 1, 1)
## geom dist for duration of exposure
duration <- function(x) dgeom(x, prob = .9)
x_with_contacts <- simulate_contacts(
x[1:10, ],
n_contacts = 1:10, # 1 to 10 contacts
dist_time_to_contact = time_to_contact,
dist_duration = duration)
## check output
class(x_with_contacts)
dim(x_with_contacts)
head(x_with_contacts)
plot(x_with_contacts)
```
# Resources
## Vignettes
No vignette currently available.
## Getting help online
Bug reports and feature requests should be posted on *github* using the [*issue*](http://github.com/reconhub/simulacr/issues) system. All other questions should be posted on the **RECON forum**: <br>
[http://www.repidemicsconsortium.org/forum/](http://www.repidemicsconsortium.org/forum/)
Contributions are welcome via **pull requests**.
Please note that this project is released with a [Contributor Code of Conduct](CONDUCT.md). By participating in this project you agree to abide by its terms.