Refactor concordance analysis

analysis/concordance.Rmd

@@ -18,44 +18,35 @@ library(tidyr)
 library(SingleCellExperiment)
 ```
 
-
 ```{r load-data}
 source(here::here("code", "top-genes.R"))
 
-results_folder <- here::here("results")
-file_names <- list.files(results_folder, full.names = TRUE)
-
-top_genes <- tibble(path = file_names) %>% 
-  # HACK: Only use one simulation for now.
-  filter(grepl("sim_1", path)) %>% 
-  expand_grid(n_genes = c(20, 100, 1000)) %>% 
-  rowwise() %>% 
-  mutate(
-    output = list(readRDS(path)), 
-    top_genes = list(get_top_genes(output, n_genes)),
-  ) 
-
-# FIXME: YUCK!!!!
-plot_data <- top_genes %>% 
-  select(path, n_genes, top_genes) %>% 
-  mutate(file_name = basename(path), .keep = "unused") %>% 
-  mutate(file_name = substr(file_name, 1,  nchar(file_name) - 4)) %>% 
-  separate(file_name, sep = "-", into = c("sim", "method")) %>% 
-  separate(method, sep = "_", into = c("method_name", "pars"),
-           extra = "merge") %>% 
-  separate(sim, sep = "_", into = c("sim", "sim_id")) %>% 
-  select(-c(sim))
+top_genes <- list()
+# First file is the countsimQC report.
+for (i in setdiff(seq_along(res_paths), c(206, 1))) {
+  print(i)
+  res <- readRDS(res_paths[[i]])
+  if (!(length(res$result) == 0)) {
+    top_genes[[i]] <- get_top_genes(res, 20)
+  }
+  rm(res)
+}
+
+res_names <- substr(basename(res_paths), 1, nchar(basename(res_paths)) - 4)
+top_genes <- tibble(res_name = res_names, top_genes = top_genes)
+# Remove the countsimQC report
+top_genes <- top_genes[-1, ]
+
+concordance_data <- retrive_simulation_parameters() %>% 
+  mutate(res_name = substr(file_name, 1, nchar(file_name) - 4)) %>% 
+  left_join(top_genes, by = "res_name")
 ```
 
-```{r plot-all-20}
-# Yuck, yuck, yuck!
-# Results are consistent between groups.
-# We can only show n_genes = 20 otherwise too messy.
-test <- plot_data %>% 
+```{r plot-all-10}
+test <- concordance_data %>% 
+  filter(sim_label == "standard_sim" & rep == 2) %>% 
   rowwise() %>% 
-  filter(n_genes == 20) %>% 
-  # Group 1!
-  mutate(top_genes = list(top_genes[[1]])) %>%
+  mutate(top_genes = list(top_genes[[1]])) %>% 
   unnest(top_genes) %>% 
   nest_by(top_genes) %>% 
   mutate(data = list(data[["pars"]])) %>%