Handling dependencies, documentation, and modularity through an R package

base.r

@@ -7,8 +7,9 @@ library(tidyr)
 library(EnvStats)
 library(optparse)
 
-source('utils/read-data.r')
-source('utils/process-covariates.r')
+# Install R package
+remotes::install_local("rpackage/")
+library(covid19model)
 
 # Commandline options and parsing
 parser <- OptionParser()

docker/Dockerfile

@@ -16,24 +16,30 @@ RUN apt-get install -y r-cran-ggpubr
 RUN apt-get install -y r-cran-nortest
 RUN apt-get install -y r-cran-cowplot
 RUN apt-get install -y r-cran-bh
+RUN apt-get install -y r-cran-remotes
 
-RUN Rscript -e "update.packages(.libPaths()[1])"
+RUN Rscript -e "update.packages(.libPaths()[1]);remotes::install_local('rpackage/')"
 
 RUN install2.r --error --deps FALSE \
- EnvStats BH ggplot2 isoband \
+ EnvStats remotes BH ggplot2 isoband \
  && rm -rf /tmp/downloaded_packages/ /tmp/*.rds
 
 WORKDIR /var/model
 
 COPY data data/
 COPY stan-models stan-models/
 COPY utils utils/
+COPY rpackage rpackage/
 COPY *.r ./
 
+
 RUN mkdir results
 RUN mkdir figures
 RUN mkdir web
 RUN touch Rplots.pdf
 RUN chmod -R a+rwx data results figures web Rplots.pdf
+RUN touch Rplots.pdf
+RUN mkdir -p /usr/local/lib/R/site-library 
+RUN chmod 777 /usr/local/lib/R/site-library
 
 CMD Rscript base.r base

environment.yml

@@ -13,5 +13,6 @@ dependencies:
  - r-abind
  - r-svglite
  - r-jsonlite
+ - r-remotes
 channels:
  - conda-forge

plot-3-panel.r

@@ -13,8 +13,8 @@ library(bayesplot)
 library(cowplot)
 library(svglite)
 library(ggplot2)
+library(covid19model)
 
-source("utils/geom-stepribbon.r")
 #---------------------------------------------------------------------------
 make_three_pannel_plot <- function(){
 

plot-forecast.r

@@ -12,8 +12,8 @@ library(gridExtra)
 library(ggpubr)
 library(bayesplot)
 library(cowplot)
+library(covid19model)
 
-source("utils/geom-stepribbon.r")
 #---------------------------------------------------------------------------
 make_forecast_plot <- function(){
 

rpackage/.Rbuildignore

@@ -0,0 +1,3 @@
+^data-raw$
+^.*\.Rproj$
+^\.Rproj\.user$

rpackage/DESCRIPTION

@@ -0,0 +1,20 @@
+Package: covid19model
+Type: Package
+Title: Helper Functions for the covid19 Model from Imperial College London
+Version: 0.1
+Date: 2020-05-01
+Author: Mans Magnusson
+Maintainer: Mans Magnusson <[email protected]>
+Description: The package contain helper functions to run an analyze the covid19 Imperial College London model.
+License: MIT + file LICENSE
+RoxygenNote: 7.0.1
+Encoding: UTF-8
+Imports:
+ checkmate,
+ dplyr,
+ EnvStats,
+ ggplot2 (>= 3.3.0),
+ lubridate,
+ matrixStats
+Suggests: 
+ testthat (>= 2.1.0)

rpackage/LICENSE

@@ -0,0 +1,21 @@
+MIT License
+
+Copyright (c) 2020 Imperial College London
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.

rpackage/NAMESPACE

@@ -0,0 +1,12 @@
+# Generated by roxygen2: do not edit by hand
+
+export(GeomStepribbon)
+export(geom_stepribbon)
+export(process_covariates)
+export(read_ifr_data)
+export(read_interventions)
+export(read_obs_data)
+importFrom(ggplot2,GeomRibbon)
+importFrom(ggplot2,layer)
+importFrom(stats,ecdf)
+importFrom(utils,read.csv)

utils/geom-stepribbon.r → rpackage/R/geom-stepribbon.r

@@ -18,6 +18,7 @@
 #' last observation are missing. This processing is optimized for results
 #' from the survfit function.
 #' @examples
+#' library(ggplot2)
 #' huron <- data.frame(year = 1875:1972, level = as.vector(LakeHuron))
 #' h <- ggplot(huron, aes(year))
 #' h + geom_stepribbon(aes(ymin = level - 1, ymax = level + 1), fill = "grey70") +
@@ -50,7 +51,7 @@ geom_stepribbon <- function(
 #' @format NULL
 #' @usage NULL
 #' @export
-GeomStepribbon <- ggproto(
+GeomStepribbon <- ggplot2::ggproto(
  "GeomStepribbon", GeomRibbon, 
 
  extra_params = c("na.rm", "kmplot"),

utils/process-covariates.r → rpackage/R/process-covariates.r

@@ -1,15 +1,16 @@
-library(rstan)
-library(data.table)
-library(lubridate)
-library(gdata)
-library(dplyr)
-library(tidyr)
-library(EnvStats)
-library(scales)
-library(stringr)
-library(abind)
-
-process_covariates <- function(countries, interventions, d, ifr.by.country,N2){
+#' Process Covariates
+#' 
+#' @param countries Countries to use
+#' @param interventions Interventions object returned by [read_interventions()].
+#' @param d Observation data (see [read_obs_data()])
+#' @param ifr.by.country IFR by country (see [read_ifr_data()])
+#' @param N2 Increase if you need more forecast
+#' 
+#' @importFrom stats ecdf
+#' 
+#' 
+#' @export
+process_covariates <- function(countries, interventions, d, ifr.by.country, N2){
  serial.interval = read.csv("data/serial_interval.csv")
  # Pads serial interval with 0 if N2 is greater than the length of the serial
  # interval array
@@ -23,8 +24,8 @@ process_covariates <- function(countries, interventions, d, ifr.by.country,N2){
  # various distributions required for modeling
  mean1 <- 5.1; cv1 <- 0.86; # infection to onset
  mean2 <- 17.8; cv2 <- 0.45 # onset to death
- x1 <- rgammaAlt(1e6,mean1,cv1) # infection-to-onset distribution
- x2 <- rgammaAlt(1e6,mean2,cv2) # onset-to-death distribution
+ x1 <- EnvStats::rgammaAlt(1e6,mean1,cv1) # infection-to-onset distribution
+ x2 <- EnvStats::rgammaAlt(1e6,mean2,cv2) # onset-to-death distribution
 
  ecdf.saved <- ecdf(x1+x2)
  forecast <- 0
@@ -45,18 +46,18 @@ process_covariates <- function(countries, interventions, d, ifr.by.country,N2){
  region <-region[order(as.Date(region$DateRep)),] # ensure date ordering
 
  # padding in raw data backwards ex. portugal
- date_min <- dmy('31/12/2019') 
+ date_min <- lubridate::dmy('31/12/2019') 
  if (region$DateRep[1] > date_min){
  print(paste(Country,'In padding ECDC data'))
  pad_days <-region$DateRep[1] - date_min
- pad_dates <- date_min + days(1:pad_days[[1]]-1)
+ pad_dates <- date_min + lubridate::days(1:pad_days[[1]]-1)
  padded_data <- data.frame("Country" = rep(Country, pad_days),
  "DateRep" = pad_dates,
  "Cases" = as.integer(rep(0, pad_days)),
  "Deaths" = as.integer(rep(0, pad_days)),
  stringsAsFactors=F)
 
- region <- bind_rows(padded_data,region)
+ region <- dplyr::bind_rows(padded_data,region)
  }
  index = which(region$Cases>0)[1]
  index1 = which(cumsum(region$Deaths)>=10)[1] # also 5

utils/read-data.r → rpackage/R/read-data.r

@@ -1,8 +1,13 @@
-library(tidyr)
-library(lubridate)
-library(stringr)
-library(dplyr)
-
+#' Read Data
+#' 
+#' @details 
+#' The functions reads data from the data folder.
+#' 
+#' @param countries Countries to use include in data
+#' 
+#' @importFrom utils read.csv 
+#' 
+#' @export
 read_obs_data <- function(countries){
  # Read the deaths and cases data
  d <- readRDS('data/COVID-19-up-to-date.rds')
@@ -12,6 +17,8 @@ read_obs_data <- function(countries){
  return(d)
 }
 
+#' @rdname read_obs_data
+#' @export
 read_ifr_data <- function(){
  ifr.by.country <- read.csv("data/popt_ifr.csv")
  ifr.by.country$country <- as.character(ifr.by.country[,2])
@@ -20,12 +27,14 @@ read_ifr_data <- function(){
 
 }
 
+#' @rdname read_obs_data
+#' @export
 read_interventions <- function(countries){
  interventions = read.csv('data/interventions.csv', stringsAsFactors = FALSE)
  names_interventions = c('Schools + Universities','Self-isolating if ill', 'Public events', 'Lockdown', 'Social distancing encouraged')
  interventions <- interventions[interventions$Type %in% names_interventions,]
  interventions <- interventions[,c(1,2,4)]
- interventions <- spread(interventions, Type, Date.effective)
+ interventions <- eval(parse(text = "tidyr::spread(interventions, Type, Date.effective)"))
  names(interventions) <- c('Country','lockdown', 'public_events', 'schools_universities','self_isolating_if_ill', 'social_distancing_encouraged')
  interventions <- interventions[c('Country','schools_universities', 'self_isolating_if_ill', 'public_events', 'lockdown', 'social_distancing_encouraged')]
  interventions$schools_universities <- as.Date(interventions$schools_universities, format = "%d.%m.%Y")