Merge branch 'rebase_sars_new' of https://github.com/mrc-ide/priority…

…-pathogens into rebase_sars_new
mrc-ide · Jun 27, 2024 · 17894a6 · 17894a6
2 parents 9d64986 + 03e8f1f
commit 17894a6
Show file tree

Hide file tree

Showing 3 changed files with 38 additions and 24 deletions.
diff --git a/sars_workflow.R b/sars_workflow.R
@@ -1,4 +1,4 @@
-## LASSA WORKFLOW
+## SARS WORKFLOW
 
 #install.packages("orderly2", repos = c("https://mrc-ide.r-universe.dev", "https://cloud.r-project.org"))
 #remotes::install_github("mrc-ide/orderly.sharedfile")

diff --git a/shared/lassa_functions.R b/shared/lassa_functions.R
@@ -36,12 +36,14 @@ data_curation <- function(articles, outbreaks, models, parameters, plotting) {
            parameter_uncertainty_lower_value = case_when(
              str_detect(str_to_lower(parameter_uncertainty_singe_type),"maximum") & no_unc ~ parameter_value,
              str_detect(str_to_lower(parameter_uncertainty_singe_type),"standard deviation") & no_unc ~ parameter_value-parameter_uncertainty_single_value,
+             str_detect(str_to_lower(parameter_uncertainty_singe_type),"variance") & no_unc ~ parameter_value-sqrt(parameter_uncertainty_single_value),
              str_detect(str_to_lower(parameter_uncertainty_singe_type),"standard error") & no_unc ~ parameter_value-parameter_uncertainty_single_value,
              str_detect(str_to_lower(distribution_type),"gamma") & no_unc ~ qgamma(0.05, shape = (distribution_par1_value/distribution_par2_value)^2, rate = distribution_par1_value/distribution_par2_value^2), 
              TRUE ~ parameter_uncertainty_lower_value),                                                 
            parameter_uncertainty_upper_value = case_when(
              str_detect(str_to_lower(parameter_uncertainty_singe_type),"maximum") & no_unc ~ parameter_uncertainty_single_value,
              str_detect(str_to_lower(parameter_uncertainty_singe_type),"standard deviation") & no_unc ~ parameter_value+parameter_uncertainty_single_value,
+             str_detect(str_to_lower(parameter_uncertainty_singe_type),"variance") & no_unc ~ parameter_value+sqrt(parameter_uncertainty_single_value),
              str_detect(str_to_lower(parameter_uncertainty_singe_type),"standard error") & no_unc ~ parameter_value+parameter_uncertainty_single_value,
              str_detect(str_to_lower(distribution_type),"gamma") & no_unc ~ qgamma(0.95, shape = (distribution_par1_value/distribution_par2_value)^2, rate = distribution_par1_value/distribution_par2_value^2), 
              TRUE ~ parameter_uncertainty_upper_value)) %>%
@@ -83,7 +85,7 @@ curation <- function(articles, outbreaks, models, parameters, plotting) {
 
 # function to produce forest plot for given dataframe
 
-forest_plot <- function(df, label, color_column, lims) {
+forest_plot <- function(df, label, color_column, lims, text_size = 11) {
 
   stopifnot(length(unique(df$parameter_unit[!is.na(df$parameter_unit)])) == 1)#values must have same units
 
@@ -106,12 +108,13 @@ forest_plot <- function(df, label, color_column, lims) {
   if (all(df$parameter_class=="Reproduction number")) {gg <- gg + geom_vline(xintercept = 1, linetype = "dashed", colour = "dark grey")}
 
   gg <- gg + scale_fill_lancet(palette = "lanonc") + scale_color_lancet(palette = "lanonc") +
-        scale_shape_manual(name = "Parameter Type",values = c(Mean = 21, Median = 22, Unspecified = 24),breaks = c("Mean", "Median", "Unspecified")) +
+        scale_shape_manual(name = "Parameter Type",values = c(Mean = 21, Median = 22, Unspecified = 24, Other = 23),breaks = c("Mean", "Median", "Unspecified", "Other")) +
         scale_x_continuous(limits = lims, expand = c(0, 0)) +
         scale_y_discrete(labels = setNames(df$refs, df$urefs)) +
         labs(x = label, y = NULL) +
         theme_minimal() + 
-        theme(panel.border = element_rect(color = "black", size = 1.25, fill = NA))
+        theme(panel.border = element_rect(color = "black", size = 1.25, fill = NA),
+              text = element_text(size = text_size))
   if (cats == 1) {
     gg <- gg + guides(fill = "none", color = FALSE, shape = guide_legend(title = NULL,order = 1))
   } else {

diff --git a/src/sars_delays/orderly.R b/src/sars_delays/orderly.R
@@ -22,11 +22,12 @@ orderly_dependency("db_compilation", "latest(parameter:pathogen == this:pathogen
                    c("articles.csv", "models.csv", "parameters.csv"))
 orderly_shared_resource("lassa_functions.R" = "lassa_functions.R")
 source("lassa_functions.R")
-orderly_artefact("sars-specific tables",c("figure_delays.png"))
+orderly_artefact("sars-specific tables",c("figure_delays.png","figure_delays.pdf"))
 
 ###################
 ## DATA CURATION ##
 ###################
+TEXT_SIZE <- 28
 
 articles   <- read_csv("articles.csv")
 models     <- read_csv("models.csv")
@@ -40,7 +41,9 @@ articles   <- epireview::assign_qa_score(articles = articles)$articles
 qa_scores  <- articles %>% dplyr::select(covidence_id,qa_score)
 
 models     <- dfs$models
-parameters <- dfs$parameters %>% left_join(qa_scores) 
+parameters <- dfs$parameters %>% left_join(qa_scores) %>%
+  mutate(article_label = make.unique(refs)) %>%
+  mutate(article_label = factor(article_label,levels=rev(unique(article_label))))
 
 # sample sd is the same as sd, so need to make this clear somehow / change param field to work for meta analysis
 
@@ -54,9 +57,14 @@ d3 <- parameters %>% filter(parameter_type == 'Human delay - time in care (lengt
 d4 <- parameters %>% filter(parameter_type == 'Human delay - symptom onset>discharge/recovery' |
                               parameter_type == 'Human delay - symptom onset>death')
 d5 <- parameters %>% filter(parameter_type == 'Human delay - serial interval')
+d6 <- parameters %>% filter(parameter_type == 'Human delay - infectious period')
+d7 <- parameters %>% filter(parameter_type == 'Human delay - generation time')
 
 # Serial interval sub-plot (both papers have qa above 0.5)
-SI_forest <- forest_plot(d5 %>% filter(qa_score>0.5),'Serial Interval (days)',"parameter_type",c(0,20))
+SI_forest <- forest_plot(d5 %>% filter(qa_score>0.5),'Serial Interval (days)',"parameter_type",c(0,20),text_size = TEXT_SIZE)
+IP_forest <- forest_plot(d6 %>% filter(qa_score>0.5),'Infectious Period (days)',"parameter_type",c(0,30),text_size = TEXT_SIZE)
+#IP_forest2 <- epireview::forest_plot_infectious_period(d6 %>% filter(qa_score>0.5), shape_by = 'parameter_value_type', ulim=30)
+GT_forest <- forest_plot(d7 %>% filter(qa_score>0.5),'Generation Time (days)',"parameter_type",c(0,20),text_size = TEXT_SIZE)
 
 # Incubation period
 #p1 <- forest_plot(d1,'Incubation Period (days)',"parameter_type",c(0,30))
@@ -65,7 +73,7 @@ d1_subgroups <- d1 %>% filter(population_group %in% c('General population','Mixe
   filter(qa_score>0.5) %>%
   mutate(parameter_value = coalesce(parameter_value,central))
 m1 <- metamean_wrap(dataframe = d1_subgroups, estmeansd_method = "Cai", 
-                    plot_study = FALSE, digits = 2, lims = c(0,10), colour = "dodgerblue3", label = "Incubation Period",
+                    plot_study = FALSE, digits = 2, lims = c(0,10), colour = "dodgerblue3", label = "Mean incubation period (days)",
                     width = 7500, height = 5750, resolution = 1000, subgroup = 'population_group', sort_by_subg = TRUE )
 
 #meta-analysis of onset-admission delay
@@ -78,24 +86,27 @@ m2 <- metamean_wrap(dataframe = d2, estmeansd_method = "Cai",
                     width = 9500, height = 4200, resolution = 1000)
 
 # admission to outcome... 
-d3_subgroup <- d3 %>% 
-  filter(!is.na(parameter_uncertainty_single_value)|!is.na(parameter_uncertainty_type) ) %>% 
-  mutate(group = stringr::str_to_title(str_replace(parameter_type,'Human delay - ','')),
-         parameter_type = 'Human delay') %>%  
-  filter(qa_score>0.5)
+outcome_forest <- forest_plot(d3 %>% filter(qa_score>0.5) %>%
+                                mutate(parameter_type = stringr::str_to_title(str_replace(parameter_type,'Human delay - ',''))),
+                              'Admission to outcome',"parameter_type",c(0,60),text_size = 22)
 
-outcome_forest <- forest_plot(d3 %>% filter(qa_score>0.5),'Admission to outcome',"parameter_type",c(0,50))
+layout <- "
+AABB#
+CCDD#
+EEFFF
+EEFFF
+"
+delays_plot <-  SI_forest + outcome_forest + GT_forest + IP_forest + m1$plot + m2$plot  + plot_layout(design = layout) + plot_annotation(tag_levels = 'A') 
 
+ggsave("figure_delays.png", plot = delays_plot, width = 39, height = 22)
+ggsave("figure_delays.pdf", plot = delays_plot, width = 39, height = 22)
 
-m3 <- metamean_wrap(dataframe = d3_subgroup, estmeansd_method = "Cai", 
-                    plot_study = FALSE, digits = 2, lims = c(8,25), colour = "dodgerblue3", label = "Admission to outcome",
-                    width = 9500, height = 4500, resolution = 1000, subgroup = 'group', sort_by_subg = TRUE )
 
-layout <- "
-AABBB
-CCDDD
-CCDDD
+layout2 <- "
+AAABBBDDD
+EEEEFFFFF
+EEEEFFFFF
 "
-delays_plot <-  SI_forest + outcome_forest + m1$plot + m2$plot  + plot_layout(design = layout) + plot_annotation(tag_levels = 'A')
-
-ggsave("figure_delays.png", plot = delays_plot, width = 28, height = 18)
+delays_plot2 <-  SI_forest + IP_forest + outcome_forest + m1$plot + theme(text = element_text(size = TEXT_SIZE)) + m2$plot + theme(text = element_text(size = TEXT_SIZE)) + 
+  plot_layout(design = layout2) + plot_annotation(tag_levels = 'A') 
+ggsave("figure_delays2.png", plot = delays_plot2, width = 39, height = 22)