Plots credible intervals and median for the observed data under the posterior predictive distribution, and for a specific observation type. The user can control the interval levels (i.e. 30%, 50% etc.) and the plotted group(s). This is a generic function.
plot_obs(object, ...) # S3 method for epimodel plot_obs( object, type, groups = NULL, dates = NULL, date_breaks = "2 weeks", date_format = "%Y-%m-%d", cumulative = FALSE, by_100k = FALSE, bar = TRUE, levels = c(30, 60, 90), log = FALSE, ... ) spaghetti_obs( object, type, draws = min(500, posterior_sample_size(object)), alpha = 1/sqrt(draws), groups = NULL, dates = NULL, date_breaks = "2 weeks", date_format = "%Y-%m-%d", cumulative = FALSE, by_100k = FALSE, bar = TRUE, log = FALSE, smooth = 1, ... )
| object | A fitted model object returned by |
|---|---|
| ... | Additional arguments for
|
| type | A string specifying the name of the observations to plot. This should match one
of the names of the response variables in the |
| groups | Either |
| dates | A length 2 vector of |
| date_breaks | A string giving the distance between date tick labels.
Default is |
| date_format | This function attempts to coerce the |
| cumulative | If |
| by_100k | If |
| bar | If |
| levels | A numeric vector defining the levels of the plotted credible intervals. |
| log | If |
| draws | The number of sample paths to plot. |
| alpha | Sets transparency of sample paths. |
| smooth | An integer specifying the window used to smooth the reproduction rates. The
default is |
A ggplot object which can be further modified.
# \donttest{ data("EuropeCovid2") data <- EuropeCovid2$data data <- dplyr::filter(data, date > date[which(cumsum(deaths) > 10)[1] - 30]) data <- dplyr::filter(data, date < as.Date("2020-05-05")) rt <- epirt( formula = R(country, date) ~ 0 + (1 + public_events + schools_universities + self_isolating_if_ill + social_distancing_encouraged + lockdown || country) + public_events + schools_universities + self_isolating_if_ill + social_distancing_encouraged + lockdown, prior = shifted_gamma(shape=1/6, scale = 1, shift = log(1.05)/6), prior_covariance = rstanarm::decov(shape = c(2, rep(0.5, 5)),scale=0.25), link = scaled_logit(6.5) ) inf <- epiinf(gen = EuropeCovid$si, seed_days = 6) deaths <- epiobs( formula = deaths ~ 1, i2o = EuropeCovid2$inf2death, prior_intercept = rstanarm::normal(0,0.2), link = scaled_logit(0.02) ) args <- list(rt=rt, inf=inf, obs=deaths, data=data, seed=12345) args$group_subset <- c("Italy", "Austria", "Germany") args$algorithm <- "fullrank" args$iter <- 1e4 args$tol_rel_obj <- 1e-3 fm <- do.call(epim, args)#> Chain 1: ------------------------------------------------------------ #> Chain 1: EXPERIMENTAL ALGORITHM: #> Chain 1: This procedure has not been thoroughly tested and may be unstable #> Chain 1: or buggy. The interface is subject to change. #> Chain 1: ------------------------------------------------------------ #> Chain 1: #> Chain 1: #> Chain 1: #> Chain 1: Gradient evaluation took 0.000528 seconds #> Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 5.28 seconds. #> Chain 1: Adjust your expectations accordingly! #> Chain 1: #> Chain 1: #> Chain 1: Begin eta adaptation. #> Chain 1: Iteration: 1 / 250 [ 0%] (Adaptation) #> Chain 1: Iteration: 50 / 250 [ 20%] (Adaptation) #> Chain 1: Iteration: 100 / 250 [ 40%] (Adaptation) #> Chain 1: Iteration: 150 / 250 [ 60%] (Adaptation) #> Chain 1: Iteration: 200 / 250 [ 80%] (Adaptation) #> Chain 1: Iteration: 250 / 250 [100%] (Adaptation) #> Chain 1: Success! Found best value [eta = 0.1]. #> Chain 1: #> Chain 1: Begin stochastic gradient ascent. #> Chain 1: iter ELBO delta_ELBO_mean delta_ELBO_med notes #> Chain 1: 100 -13688.867 1.000 1.000 #> Chain 1: 200 -5751.465 1.190 1.380 #> Chain 1: 300 -4595.706 0.877 1.000 #> Chain 1: 400 -3865.385 0.705 1.000 #> Chain 1: 500 -3424.722 0.590 0.251 #> Chain 1: 600 -2755.889 0.532 0.251 #> Chain 1: 700 -2586.841 0.465 0.243 #> Chain 1: 800 -2307.851 0.422 0.243 #> Chain 1: 900 -2066.580 0.388 0.189 #> Chain 1: 1000 -2010.719 0.352 0.189 #> Chain 1: 1100 -1692.408 0.271 0.188 #> Chain 1: 1200 -1708.032 0.134 0.129 #> Chain 1: 1300 -1580.999 0.117 0.121 #> Chain 1: 1400 -1432.244 0.108 0.117 #> Chain 1: 1500 -1299.048 0.106 0.104 #> Chain 1: 1600 -1237.094 0.086 0.103 #> Chain 1: 1700 -1151.716 0.087 0.103 #> Chain 1: 1800 -1118.903 0.078 0.080 #> Chain 1: 1900 -1162.818 0.070 0.074 #> Chain 1: 2000 -1056.263 0.078 0.080 #> Chain 1: 2100 -1060.052 0.059 0.074 #> Chain 1: 2200 -1028.301 0.061 0.074 #> Chain 1: 2300 -1046.946 0.055 0.050 #> Chain 1: 2400 -1010.439 0.048 0.038 #> Chain 1: 2500 -994.414 0.040 0.036 #> Chain 1: 2600 -968.220 0.037 0.031 #> Chain 1: 2700 -985.944 0.032 0.029 #> Chain 1: 2800 -968.507 0.031 0.027 #> Chain 1: 2900 -962.022 0.028 0.018 #> Chain 1: 3000 -971.286 0.018 0.018 #> Chain 1: 3100 -955.881 0.020 0.018 #> Chain 1: 3200 -942.191 0.018 0.018 #> Chain 1: 3300 -945.570 0.017 0.016 #> Chain 1: 3400 -946.000 0.013 0.016 #> Chain 1: 3500 -939.521 0.012 0.015 #> Chain 1: 3600 -939.380 0.009 0.010 #> Chain 1: 3700 -934.574 0.008 0.007 #> Chain 1: 3800 -931.710 0.007 0.007 #> Chain 1: 3900 -931.896 0.006 0.005 #> Chain 1: 4000 -935.988 0.005 0.004 #> Chain 1: 4100 -931.628 0.004 0.004 #> Chain 1: 4200 -933.590 0.003 0.004 #> Chain 1: 4300 -929.327 0.003 0.004 #> Chain 1: 4400 -931.544 0.003 0.004 #> Chain 1: 4500 -927.411 0.003 0.004 #> Chain 1: 4600 -937.806 0.004 0.004 #> Chain 1: 4700 -930.599 0.004 0.004 #> Chain 1: 4800 -923.500 0.005 0.005 #> Chain 1: 4900 -922.723 0.005 0.005 #> Chain 1: 5000 -925.086 0.005 0.005 #> Chain 1: 5100 -919.915 0.005 0.005 #> Chain 1: 5200 -927.492 0.006 0.006 #> Chain 1: 5300 -918.643 0.006 0.008 #> Chain 1: 5400 -923.476 0.006 0.008 #> Chain 1: 5500 -913.747 0.007 0.008 #> Chain 1: 5600 -919.763 0.006 0.008 #> Chain 1: 5700 -925.350 0.006 0.007 #> Chain 1: 5800 -914.064 0.007 0.007 #> Chain 1: 5900 -918.877 0.007 0.007 #> Chain 1: 6000 -919.968 0.007 0.007 #> Chain 1: 6100 -924.044 0.007 0.007 #> Chain 1: 6200 -918.862 0.007 0.006 #> Chain 1: 6300 -913.823 0.006 0.006 #> Chain 1: 6400 -913.453 0.006 0.006 #> Chain 1: 6500 -911.831 0.005 0.006 #> Chain 1: 6600 -913.527 0.004 0.005 #> Chain 1: 6700 -913.369 0.004 0.004 #> Chain 1: 6800 -918.911 0.003 0.004 #> Chain 1: 6900 -919.172 0.003 0.002 #> Chain 1: 7000 -911.534 0.003 0.004 #> Chain 1: 7100 -907.337 0.003 0.005 #> Chain 1: 7200 -915.610 0.004 0.005 #> Chain 1: 7300 -906.436 0.004 0.005 #> Chain 1: 7400 -910.459 0.005 0.005 #> Chain 1: 7500 -910.673 0.005 0.005 #> Chain 1: 7600 -911.118 0.004 0.005 #> Chain 1: 7700 -909.753 0.005 0.005 #> Chain 1: 7800 -910.293 0.004 0.004 #> Chain 1: 7900 -908.857 0.004 0.004 #> Chain 1: 8000 -909.975 0.003 0.002 #> Chain 1: 8100 -907.878 0.003 0.002 #> Chain 1: 8200 -908.813 0.002 0.002 #> Chain 1: 8300 -906.957 0.002 0.002 #> Chain 1: 8400 -909.405 0.001 0.002 #> Chain 1: 8500 -906.035 0.002 0.002 #> Chain 1: 8600 -916.867 0.003 0.002 #> Chain 1: 8700 -908.479 0.004 0.002 #> Chain 1: 8800 -907.653 0.004 0.002 #> Chain 1: 8900 -903.013 0.004 0.003 #> Chain 1: 9000 -911.091 0.005 0.004 #> Chain 1: 9100 -904.717 0.005 0.005 #> Chain 1: 9200 -906.106 0.005 0.005 #> Chain 1: 9300 -903.341 0.005 0.005 #> Chain 1: 9400 -906.727 0.006 0.005 #> Chain 1: 9500 -906.461 0.005 0.005 #> Chain 1: 9600 -903.935 0.004 0.004 #> Chain 1: 9700 -900.322 0.004 0.004 #> Chain 1: 9800 -906.898 0.004 0.004 #> Chain 1: 9900 -902.999 0.004 0.004 #> Chain 1: 10000 -907.005 0.004 0.004 #> Chain 1: Informational Message: The maximum number of iterations is reached! The algorithm may not have converged. #> Chain 1: This variational approximation is not guaranteed to be meaningful. #> Chain 1: #> Chain 1: Drawing a sample of size 1000 from the approximate posterior... #> Chain 1: COMPLETED.#> Warning: Pareto k diagnostic value is 2.28. Resampling is disabled. Decreasing tol_rel_obj may help if variational algorithm has terminated prematurely. Otherwise consider using sampling instead.# different ways of using plot_rt p <- plot_rt(fm) # default, plots all groups and dates p <- plot_rt(fm, dates=c("2020-03-21", NA)) # plot 21 March 2020 onwards p <- plot_rt(fm, dates=c(NA, "2020-03-20")) # plot up to 20 March 2020 p <- plot_rt(fm, dates=c("2020-03-20", "2020-04-20")) p <- plot_rt(fm, dates=c("2020-03-20", "2020-04-20"), date_breaks="1 day") # ticks every day p <- plot_rt(fm, dates=c("2020-20-03", "2020-20-04"), date_format="%Y-%d-%m") # (different date format) # other plotting functions p <- plot_obs(fm, type = "deaths") p <- plot_infections(fm) p <- plot_infectious(fm) # }