In nikosbosse/epipredictr: What the Package Does (One Line, Title Case)

#knitr::opts_knit$set(root.dir = 'epipredictr/')
#knitr::opts_knit$set(fig.path = 'inst/report/figure/')

library(epipredictr)
library(rstan)
options(mc.cores = 4)
rstan_options(auto_write = TRUE)
library(scoringutils)
library(dplyr)
library(cowplot)
# source("./R/utilities.R")
# source("./R/models.R")

Authors

Summary

Aim

Forecast the time-varying R_t estimates into the future to get a sense of the trajectory of the COVID-19 epidemic in different countries.

Abstract

Methods

We used three different models

1. windowed regression
2. BSTS with a global trend
3. BSTS with a local trend

to forecast transmission rates. Predictions were made using R and Stan and results were scored using CRPS and LogS using the scoringRules package.

Key Results

[ ] graph with predictions by every country for best model
[ ] table with 
    - CI for forecasts
    - Past forecast score?

Details

d <- readRDS("./data/time_varying_params.rds")[[1]]
y_true <- d$mean

## Windowed Regression 
# model_lin_reg <- stan_model(file = "./inst/stan/linear_regression.stan")
# res_lin <- fit_iteratively(incidences = y_true, model = model_lin_reg, 
#                          n_pred = 7, 
#                          max_n_past_obs = 7, vb = FALSE)

res_lin <- fit_iteratively(incidences = y_true, model = "lin_reg", n_pred = 7,                       max_n_past_obs = 7, vb = FALSE)



## BSTS global trend
# model_bsts <- stan_model(file = "./inst/stan/bsts.stan")
# res_bsts <- fit_iteratively(incidences = y_true, 
#                           model = model_bsts, n_pred = 7, 
#                           max_n_past_obs = Inf, vb = FALSE)
res_bsts <- fit_iteratively(incidences = y_true, 
                            model = "bsts", n_pred = 7, 
                            max_n_past_obs = Inf, vb = FALSE)



## BSTS local trend
# model_bsts_local <- stan_model(file = "./inst/stan/bsts_local_trend.stan")
# res_bsts_local <- fit_iteratively(incidences = y_true, 
#                                 model = model_bsts_local, 
#                                 n_pred = 7, 
#                                 max_n_past_obs = Inf, vb = FALSE)

res_bsts_local <- fit_iteratively(incidences = y_true, 
                            model = "bsts_local_trend", n_pred = 7, 
                            max_n_past_obs = Inf, vb = FALSE, 
                            length_local_trend = 7)

Model Fit

Prior vs. Posterior Distributions

BSTS model global trend

prior_post <- plot_prior_vs_posterior(res_bsts$stanfitobjects)
prior_post$plot

BSTS model local trend

prior_post <- plot_prior_vs_posterior(res_bsts_local$stanfitobjects)
prior_post$plot

Forecasts

Table with forecasts

Plots

p_reg <- plot_pred_vs_true(y_pred_samples = res_lin$predictive_samples, 
                        y_true = res_lin$y, 
                        forecast_run = res_lin$forecast_run)


p_bsts <- plot_pred_vs_true(y_pred_samples = res_bsts$predictive_samples, 
                        y_true = res_bsts$y, 
                        forecast_run = res_bsts$forecast_run)

p_bsts_local <- plot_pred_vs_true(y_pred_samples = res_bsts_local$predictive_samples, 
                        y_true = res_bsts_local$y, 
                        forecast_run = res_bsts_local$forecast_run)


plot_grid(p_reg, p_bsts, p_bsts_local, labels = "AUTO", ncol = 1)

Scoring of Forecasts

scoringutils::eval_forecasts(true_values = y_true[15:76], 
                             predictions = res_lin$predictive_samples[15:76, ], 
                             outcome_type = "continuous")

scoringutils::eval_forecasts(true_values = y_true[15:76], 
                             predictions = res_bsts$predictive_samples[15:76, ], 
                             outcome_type = "continuous")

scoringutils::eval_forecasts(true_values = y_true[15:76], 
                             predictions = res_bsts_local$predictive_samples[15:76, ], 
                             outcome_type = "continuous")

nikosbosse/epipredictr documentation built on April 7, 2020, 9:51 p.m.