Plot regional or national Rt for COVID-19 based on FHM data to reproduce FHM reports

restim.R

# Plot regional or national Rt for COVID-19 based on FHM data, to reproduce the
# results published on https://www.folkhalsomyndigheten.se/smittskydd-beredskap/utbrott/aktuella-utbrott/covid-19/statistik-och-analyser/analys-och-prognoser/
# Adapted from:
# Churches (2020, Feb. 18). Tim Churches Health Data Science Blog: Analysing COVID-19 (2019-nCoV) outbreak data with R - part 1. Retrieved from https://timchurches.github.io/blog/posts/2020-02-18-analysing-covid-19-2019-ncov-outbreak-data-with-r-part-1/, CC BY-SA

library(readxl)
library(EpiEstim)
# path to latest FHM Excel file
secov <- read_excel("data/FHM/FHM_latest.xlsx")
secov$Statistikdatum <- as.Date(secov$Statistikdatum)
regframe <- function(reg, sdate, edate) {
    df <- data.frame(secov$Statistikdatum, secov[[reg]])
    colnames(df) <- c("dates", "I")
    return(df[df$dates>=sdate & df$dates<=edate, ])
}
secov_par_si <- function(rf) {
    return(estimate_R(rf, method = "parametric_si",
        config = make_config(list(mean_si = 4.8, std_si = 2.3))))
}
secov_par_si_u <- function(rf) {
    return(estimate_R(rf, method = "uncertain_si",
        config = make_config(list(mean_si = 4.8, std_mean_si = 2, min_mean_si = 1.2,
            max_mean_si = 8.4, std_si = 2.3, std_std_si = 1,
            min_std_si = 0.6, max_std_si = 4, n1 = 100, n2 = 100))))
}

# plot_Ri(secov_par_si(regframe("Uppsala", "2020-09-01", "2021-02-25")))
plot_Ri <- function(estimate_R_obj) {
    p_I <- plot(estimate_R_obj, "incid", add_imported_cases = TRUE)  # plots the incidence
    p_SI <- plot(estimate_R_obj, "SI")  # plots the serial interval distribution
    p_Ri <- plot(estimate_R_obj, "R")
    return(gridExtra::grid.arrange(p_I, p_SI, p_Ri, ncol = 1))
}