This function implements a sequential Bayesian estimation method of R0 due to +Bettencourt and Riberio (PloS One, 2008). See details for important +implementation notes.
+Arguments
+ + +- cases +
Vector of case counts. The vector must only contain non-negative +integers, and have at least two positive integers.
+
+
+- mu +
Mean of the serial distribution. This must be a positive number. +The value should match the case counts in time units. For example, if case +counts are weekly and the serial distribution has a mean of seven days, +then
mu
should be set to1
. If case counts are daily and the serial +distribution has a mean of seven days, thenmu
should be set to7
.
+
+
+- kappa +
Largest possible value of the uniform prior (defaults to
20
). +This must be a number greater than or equal to1
. It describes the prior +belief on the ranges of R0, and should be set to a higher value if R0 is +believed to be larger.
+
+
+- post +
Whether to return the posterior distribution of R0 instead of the +estimate of R0 (defaults to
FALSE
). This must be a value identical to +TRUE
orFALSE
.
+
+
Value
+If post
is identical to TRUE
, a list containing the following
+components is returned:
supp
- the support of the posterior distribution of R0
+pmf
- the probability mass function of the posterior distribution of R0
+
Otherwise, if post
is identical to FALSE
, only the estimate of R0 is
+returned. Note that the estimate is equal to sum(supp * pmf)
(i.e., the
+posterior mean).
Details
+The method sets a uniform prior distribution on R0 with possible values
+between 0
and kappa
, discretized to a fine grid. The distribution of R0
+is then updated sequentially, with one update for each new case count
+observation. The final estimate of R0 is the mean of the (last) posterior
+distribution. The prior distribution is the initial belief of the
+distribution of R0, which is the uninformative uniform distribution with
+values between 0
and kappa
. Users can change the value of kappa
only
+(i.e., the prior distribution cannot be changed from the uniform). As more
+case counts are observed, the influence of the prior distribution should
+lessen on the final estimate.
This method is based on an approximation of the SIR model, which is most +valid at the beginning of an epidemic. The method assumes that the mean of +the serial distribution (sometimes called the serial interval) is known. The +final estimate can be quite sensitive to this value, so sensitivity testing +is strongly recommended. Users should be careful about units of time (e.g., +are counts observed daily or weekly?) when implementing.
+Our code has been modified to provide an estimate even if case counts equal +to zero are present in some time intervals. This is done by grouping the +counts over such periods of time. Without grouping, and in the presence of +zero counts, no estimate can be provided.
+See also
+vignette("seq_bayes_post", package = "Rnaught")
for examples of
+using the posterior distribution.
Examples
+# Weekly data.
+cases <- c(1, 4, 10, 5, 3, 4, 19, 3, 3, 14, 4)
+
+# Obtain R0 when the serial distribution has a mean of five days.
+seq_bayes(cases, mu = 5 / 7)
+#> [1] 1.026563
+
+# Obtain R0 when the serial distribution has a mean of three days.
+seq_bayes(cases, mu = 3 / 7)
+#> [1] 1.015938
+
+# Obtain R0 when the serial distribution has a mean of seven days, and R0 is
+# believed to be at most 4.
+estimate <- seq_bayes(cases, mu = 1, kappa = 4)
+
+# Same as above, but return the posterior distribution of R0 instead of the
+# estimate.
+posterior <- seq_bayes(cases, mu = 1, kappa = 4, post = TRUE)
+
+# Display the support and probability mass function of the posterior.
+posterior$supp
+#> [1] 0.00 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.10 0.11 0.12 0.13 0.14
+#> [16] 0.15 0.16 0.17 0.18 0.19 0.20 0.21 0.22 0.23 0.24 0.25 0.26 0.27 0.28 0.29
+#> [31] 0.30 0.31 0.32 0.33 0.34 0.35 0.36 0.37 0.38 0.39 0.40 0.41 0.42 0.43 0.44
+#> [46] 0.45 0.46 0.47 0.48 0.49 0.50 0.51 0.52 0.53 0.54 0.55 0.56 0.57 0.58 0.59
+#> [61] 0.60 0.61 0.62 0.63 0.64 0.65 0.66 0.67 0.68 0.69 0.70 0.71 0.72 0.73 0.74
+#> [76] 0.75 0.76 0.77 0.78 0.79 0.80 0.81 0.82 0.83 0.84 0.85 0.86 0.87 0.88 0.89
+#> [91] 0.90 0.91 0.92 0.93 0.94 0.95 0.96 0.97 0.98 0.99 1.00 1.01 1.02 1.03 1.04
+#> [106] 1.05 1.06 1.07 1.08 1.09 1.10 1.11 1.12 1.13 1.14 1.15 1.16 1.17 1.18 1.19
+#> [121] 1.20 1.21 1.22 1.23 1.24 1.25 1.26 1.27 1.28 1.29 1.30 1.31 1.32 1.33 1.34
+#> [136] 1.35 1.36 1.37 1.38 1.39 1.40 1.41 1.42 1.43 1.44 1.45 1.46 1.47 1.48 1.49
+#> [151] 1.50 1.51 1.52 1.53 1.54 1.55 1.56 1.57 1.58 1.59 1.60 1.61 1.62 1.63 1.64
+#> [166] 1.65 1.66 1.67 1.68 1.69 1.70 1.71 1.72 1.73 1.74 1.75 1.76 1.77 1.78 1.79
+#> [181] 1.80 1.81 1.82 1.83 1.84 1.85 1.86 1.87 1.88 1.89 1.90 1.91 1.92 1.93 1.94
+#> [196] 1.95 1.96 1.97 1.98 1.99 2.00 2.01 2.02 2.03 2.04 2.05 2.06 2.07 2.08 2.09
+#> [211] 2.10 2.11 2.12 2.13 2.14 2.15 2.16 2.17 2.18 2.19 2.20 2.21 2.22 2.23 2.24
+#> [226] 2.25 2.26 2.27 2.28 2.29 2.30 2.31 2.32 2.33 2.34 2.35 2.36 2.37 2.38 2.39
+#> [241] 2.40 2.41 2.42 2.43 2.44 2.45 2.46 2.47 2.48 2.49 2.50 2.51 2.52 2.53 2.54
+#> [256] 2.55 2.56 2.57 2.58 2.59 2.60 2.61 2.62 2.63 2.64 2.65 2.66 2.67 2.68 2.69
+#> [271] 2.70 2.71 2.72 2.73 2.74 2.75 2.76 2.77 2.78 2.79 2.80 2.81 2.82 2.83 2.84
+#> [286] 2.85 2.86 2.87 2.88 2.89 2.90 2.91 2.92 2.93 2.94 2.95 2.96 2.97 2.98 2.99
+#> [301] 3.00 3.01 3.02 3.03 3.04 3.05 3.06 3.07 3.08 3.09 3.10 3.11 3.12 3.13 3.14
+#> [316] 3.15 3.16 3.17 3.18 3.19 3.20 3.21 3.22 3.23 3.24 3.25 3.26 3.27 3.28 3.29
+#> [331] 3.30 3.31 3.32 3.33 3.34 3.35 3.36 3.37 3.38 3.39 3.40 3.41 3.42 3.43 3.44
+#> [346] 3.45 3.46 3.47 3.48 3.49 3.50 3.51 3.52 3.53 3.54 3.55 3.56 3.57 3.58 3.59
+#> [361] 3.60 3.61 3.62 3.63 3.64 3.65 3.66 3.67 3.68 3.69 3.70 3.71 3.72 3.73 3.74
+#> [376] 3.75 3.76 3.77 3.78 3.79 3.80 3.81 3.82 3.83 3.84 3.85 3.86 3.87 3.88 3.89
+#> [391] 3.90 3.91 3.92 3.93 3.94 3.95 3.96 3.97 3.98 3.99 4.00
+posterior$pmf
+#> [1] 4.396081e-14 6.866777e-14 1.069979e-13 1.663113e-13 2.578585e-13
+#> [6] 3.987896e-13 6.151736e-13 9.465244e-13 1.452563e-12 2.223289e-12
+#> [11] 3.393931e-12 5.167074e-12 7.845296e-12 1.187914e-11 1.793742e-11
+#> [16] 2.700983e-11 4.055633e-11 6.072366e-11 9.065827e-11 1.349567e-10
+#> [21] 2.003117e-10 2.964355e-10 4.373744e-10 6.433725e-10 9.435054e-10
+#> [26] 1.379386e-09 2.010358e-09 2.920742e-09 4.229916e-09 6.106257e-09
+#> [31] 8.786364e-09 1.260143e-08 1.801328e-08 2.566338e-08 3.643912e-08
+#> [36] 5.156328e-08 7.271379e-08 1.021837e-07 1.430935e-07 1.996715e-07
+#> [41] 2.776227e-07 3.846102e-07 5.308815e-07 7.300785e-07 1.000278e-06
+#> [46] 1.365319e-06 1.856496e-06 2.514686e-06 3.393017e-06 4.560205e-06
+#> [51] 6.104659e-06 8.139541e-06 1.080892e-05 1.429523e-05 1.882819e-05
+#> [56] 2.469542e-05 3.225496e-05 4.194988e-05 5.432506e-05 7.004650e-05
+#> [61] 8.992291e-05 1.149299e-04 1.462361e-04 1.852322e-04 2.335598e-04
+#> [66] 2.931434e-04 3.662201e-04 4.553696e-04 5.635406e-04 6.940727e-04
+#> [71] 8.507116e-04 1.037615e-03 1.259347e-03 1.520857e-03 1.827440e-03
+#> [76] 2.184679e-03 2.598363e-03 3.074378e-03 3.618573e-03 4.236597e-03
+#> [81] 4.933708e-03 5.714562e-03 6.582973e-03 7.541662e-03 8.591990e-03
+#> [86] 9.733699e-03 1.096466e-02 1.228062e-02 1.367503e-02 1.513890e-02
+#> [91] 1.666067e-02 1.822624e-02 1.981901e-02 2.142007e-02 2.300844e-02
+#> [96] 2.456146e-02 2.605520e-02 2.746508e-02 2.876640e-02 2.993509e-02
+#> [101] 3.094837e-02 3.178548e-02 3.242832e-02 3.286214e-02 3.307605e-02
+#> [106] 3.306344e-02 3.282231e-02 3.235539e-02 3.167013e-02 3.077854e-02
+#> [111] 2.969682e-02 2.844487e-02 2.704571e-02 2.552471e-02 2.390888e-02
+#> [116] 2.222595e-02 2.050365e-02 1.876888e-02 1.704700e-02 1.536123e-02
+#> [121] 1.373213e-02 1.217724e-02 1.071084e-02 9.343864e-03 8.083917e-03
+#> [126] 6.935429e-03 5.899896e-03 4.976201e-03 4.160989e-03 3.449076e-03
+#> [131] 2.833858e-03 2.307727e-03 1.862441e-03 1.489475e-03 1.180313e-03
+#> [136] 9.266886e-04 7.207803e-04 5.553459e-04 4.238132e-04 3.203273e-04
+#> [141] 2.397617e-04 1.777009e-04 1.304008e-04 9.473459e-05 6.812875e-05
+#> [146] 4.849552e-05 3.416469e-05 2.381841e-05 1.643092e-05 1.121447e-05
+#> [151] 7.572111e-06 5.057422e-06 3.340936e-06 2.182658e-06 1.410047e-06
+#> [156] 9.006636e-07 5.687529e-07 3.550312e-07 2.190488e-07 1.335660e-07
+#> [161] 8.047851e-08 4.791165e-08 2.817909e-08 1.637135e-08 9.394205e-09
+#> [166] 5.323528e-09 2.978849e-09 1.645703e-09 8.975385e-10 4.831664e-10
+#> [171] 2.566999e-10 1.345806e-10 6.961600e-11 3.552605e-11 1.788285e-11
+#> [176] 8.878079e-12 4.346434e-12 2.098061e-12 9.984180e-13 4.683320e-13
+#> [181] 2.165108e-13 9.863385e-14 4.427200e-14 1.957601e-14 8.526007e-15
+#> [186] 3.657023e-15 1.544555e-15 6.422519e-16 2.628849e-16 1.059047e-16
+#> [191] 4.198412e-17 1.637588e-17 6.283525e-18 2.371425e-18 8.801378e-19
+#> [196] 3.211841e-19 1.152248e-19 4.063045e-20 1.407979e-20 4.794057e-21
+#> [201] 1.603601e-21 5.268629e-22 1.699923e-22 5.385323e-23 1.674809e-23
+#> [206] 5.112219e-24 1.531308e-24 4.500307e-25 1.297374e-25 3.668155e-26
+#> [211] 1.016962e-26 2.764082e-27 7.363758e-28 1.922485e-28 4.917601e-29
+#> [216] 1.232200e-29 3.023823e-30 7.265876e-31 1.709166e-31 3.935069e-32
+#> [221] 8.865403e-33 1.954015e-33 4.212549e-34 8.880821e-35 1.830437e-35
+#> [226] 3.687665e-36 7.260114e-37 1.396466e-37 2.623677e-38 4.813725e-39
+#> [231] 8.622619e-40 1.507577e-40 2.572153e-41 4.281381e-42 6.950773e-43
+#> [236] 1.100362e-43 1.698170e-44 2.554224e-45 3.743317e-46 5.343929e-47
+#> [241] 7.429450e-48 1.005609e-48 1.324832e-49 1.698377e-50 2.118019e-51
+#> [246] 2.568782e-52 3.029042e-53 3.471688e-54 3.866425e-55 4.182989e-56
+#> [251] 4.394865e-57 4.482891e-58 4.438084e-59 4.263117e-60 3.972112e-61
+#> [256] 3.588768e-62 3.143145e-63 2.667736e-64 2.193531e-65 1.746742e-66
+#> [261] 1.346659e-67 1.004822e-68 7.254056e-70 5.065108e-71 3.419533e-72
+#> [266] 2.231354e-73 1.406840e-74 8.567343e-76 5.037567e-77 2.859011e-78
+#> [271] 1.565585e-79 8.268889e-81 4.210849e-82 2.066735e-83 9.773085e-85
+#> [276] 4.450901e-86 1.951504e-87 8.234368e-89 3.342436e-90 1.304664e-91
+#> [281] 4.895156e-93 1.764793e-94 6.110914e-96 2.031551e-97 6.481579e-99
+#> [286] 1.983744e-100 5.821848e-102 1.637658e-103 4.413549e-105 1.139116e-106
+#> [291] 2.814324e-108 6.652979e-110 1.504185e-111 3.251137e-113 6.714617e-115
+#> [296] 1.324527e-116 2.494329e-118 4.482263e-120 7.682227e-122 1.255211e-123
+#> [301] 1.954236e-125 2.897728e-127 4.090211e-129 5.493210e-131 7.015875e-133
+#> [306] 8.517133e-135 9.822891e-137 1.075711e-138 1.117987e-140 1.102133e-142
+#> [311] 1.030047e-144 9.121635e-147 7.649726e-149 6.072101e-151 4.559445e-153
+#> [316] 3.236849e-155 2.171329e-157 1.375543e-159 8.224662e-162 4.638796e-164
+#> [321] 2.466499e-166 1.235627e-168 5.828612e-171 2.587326e-173 1.080141e-175
+#> [326] 4.238233e-178 1.562047e-180 5.404238e-183 1.754003e-185 5.337064e-188
+#> [331] 1.521491e-190 4.061129e-193 1.014255e-195 2.368537e-198 5.168362e-201
+#> [336] 1.053099e-203 2.002304e-206 3.550040e-209 5.865079e-212 9.022847e-215
+#> [341] 1.291610e-217 1.719185e-220 2.126178e-223 2.441418e-226 2.600913e-229
+#> [346] 2.568760e-232 2.350200e-235 1.990379e-238 1.559119e-241 1.128742e-244
+#> [351] 7.546371e-248 4.655454e-251 2.647985e-254 1.387535e-257 6.692519e-261
+#> [356] 2.968866e-264 1.210268e-267 4.529956e-271 1.555446e-274 4.895396e-278
+#> [361] 1.410954e-281 3.720881e-285 8.970133e-289 1.975060e-292 3.968208e-296
+#> [366] 7.268446e-300 1.212601e-303 1.840840e-307 2.540527e-311 3.184384e-315
+#> [371] 3.621946e-319 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
+#> [376] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
+#> [381] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
+#> [386] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
+#> [391] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
+#> [396] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
+#> [401] 0.000000e+00
+
+# Note that the following always holds:
+estimate == sum(posterior$supp * posterior$pmf)
+#> [1] TRUE
+