Dev nf

bin/scflow_dge.r

@@ -146,6 +146,14 @@ required$add_argument(
   help = "p-value cutoff for DE [default %(default)s]"
 )
 
+required$add_argument(
+  "--n_label",
+  type = "integer",
+  default = 5,
+  metavar = "number",
+  help = "Number of genes to be highlighted on volcano plot"
+)
+
 required$add_argument(
   "--ensembl_mappings",
   help = "path to ensembl mappings file",
@@ -179,7 +187,9 @@ args$pseudobulk <- as.logical(args$pseudobulk)
 args$force_run <- as.logical(args$force_run)
 if (tolower(args$random_effects_var) == "null") args$random_effects_var <- NULL
 
-args$max_cores <- if (toupper(args$max_cores) == "NULL") NULL else {
+args$max_cores <- if (toupper(args$max_cores) == "NULL") {
+  NULL
+} else {
   as.numeric(as.character(args$max_cores))
 }
 
@@ -202,6 +212,7 @@ cli::cli_alert(sprintf(
   n_cores
 ))
 
+
 library(scFlow)
 
 #   ____________________________________________________________________________
@@ -221,8 +232,10 @@ if (args$pseudobulk) {
   sce_subset <- pseudobulk_sce(
     sce_subset,
     keep_vars = c(
-      args$dependent_var, args$confounding_vars, args$random_effects_var
-      ),
+      args$dependent_var,
+      args$confounding_vars,
+      args$random_effects_var
+    ),
     assay_name = "counts",
     celltype_var = args$celltype_var,
     sample_var = args$sample_var
@@ -249,26 +262,49 @@ de_results <- perform_de(
   species = getOption("scflow_species")
 )
 
-file_name <- paste0(args$celltype, "_",
-                    args$de_method, pb_str, "_")
+file_name <- paste0(
+  args$celltype, "_",
+  args$de_method, pb_str, "_"
+)
 
 for (result in names(de_results)) {
   if (dim(de_results[[result]])[[1]] > 0) {
     write.table(de_results[[result]],
-                file = file.path(getwd(),
-                                paste0(file_name, result, "_DE.tsv")),
-                quote = FALSE, sep = "\t", col.names = TRUE, row.names = FALSE)
+      file = file.path(
+        getwd(),
+        paste0(file_name, result, "_DE.tsv")
+      ),
+      quote = FALSE, sep = "\t", col.names = TRUE, row.names = FALSE
+    )
+
     report_de(de_results[[result]],
-              report_folder_path = file.path(getwd()),
-              report_file = paste0(file_name, result, "_scflow_de_report"))
+      fc_threshold = args$fc_threshold,
+      pval_cutoff = args$pval_cutoff,
+      n_label = args$n_label,
+      report_folder_path = file.path(getwd()),
+      report_file = paste0(file_name, result, "_scflow_de_report")
+    )
+
     print("report generated")
-    png(file.path(getwd(),
-                  paste0(file_name, result, "_volcano_plot.png")),
-        width = 247, height = 170, units = "mm", res = 600)
-    print(attr(de_results[[result]], "plot"))
-    dev.off()
 
+    p <- scFlow::volcano_plot(
+      dt = de_results[[result]],
+      fc_threshold = args$fc_threshold,
+      pval_cutoff = args$pval_cutoff,
+      n_label = args$n_label
+    )
+
+    ggplot2::ggsave(
+      filename = file.path(
+        getwd(),
+        paste0(file_name, result, "_volcano_plot.png")
+      ),
+      plot = p,
+      width = 7, height = 5, units = "in", dpi = 600
+    )
+
+    print("Volcano plot generated")
   } else {
     print(sprintf("No DE genes found for %s", result))
-    }
+  }
 }

bin/scflow_finalize_sce.r

@@ -5,12 +5,11 @@
 #   ____________________________________________________________________________
 #   Initialization                                                          ####
 
+options(mc.cores = future::availableCores())
+
 ##  ............................................................................
 ##  Load packages                                                           ####
 library(argparse)
-library(scFlow)
-library(magrittr)
-library(SingleCellExperiment)
 
 ##  ............................................................................
 ##  Parse command-line arguments                                            ####
@@ -106,6 +105,13 @@ required$add_argument(
   metavar = "N"
 )
 
+required$add_argument(
+  "--max_cores",
+  default = NULL,
+  help = "override for lower cpu core usage",
+  metavar = "N",
+  required = TRUE
+)
 
 ### . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ..
 ### Pre-process args                                                        ####
@@ -117,6 +123,33 @@ args$metric_vars <- strsplit(args$metric_vars, ",")[[1]]
 options("scflow_reddimplot_pointsize" = args$reddimplot_pointsize)
 options("scflow_reddimplot_alpha" = args$reddimplot_alpha)
 
+args$max_cores <- if (toupper(args$max_cores) == "NULL") {
+  NULL
+} else {
+  as.numeric(as.character(args$max_cores))
+}
+
+#   ____________________________________________________________________________
+#   Delay Package Loading for Optional Max Cores Override
+
+n_cores <- future::availableCores(methods = "mc.cores")
+
+if (is.null(args$max_cores)) {
+  options(mc.cores = n_cores)
+} else {
+  options(mc.cores = min(args$max_cores, n_cores))
+}
+
+cli::cli_alert(sprintf(
+  "Using %s cores on system with %s available cores.",
+  getOption("mc.cores"),
+  n_cores
+))
+
+library(scFlow)
+library(magrittr)
+library(SingleCellExperiment)
+
 ##  ............................................................................
 ##  Start                                                                   ####
 
@@ -163,26 +196,28 @@ colnames(celltypes) <- c("celltype", "n_cells")
 write.table(
   data.frame(celltypes),
   file = "celltypes.tsv",
-  row.names = FALSE, col.names = TRUE, quote = FALSE, sep = "\t")
+  row.names = FALSE, col.names = TRUE, quote = FALSE, sep = "\t"
+)
 
 ### Save Marker Gene Plots
 folder_path <- file.path(getwd(), "celltype_marker_plots")
 dir.create(folder_path)
 
 for (group in names(sce@metadata$markers)) {
-  pwidth <- max(10,
-                length(
-                  unique(sce@metadata$markers[[group]]$marker_plot$data$Group)
-                  )
+  pwidth <- max(
+    10,
+    length(unique(sce@metadata$markers[[group]]$marker_plot$data$Group))
   )
-  pheight <- length(
-    unique(sce@metadata$markers[[group]]$marker_plot$data$Gene)
-    )
+  pheight <- length(unique(sce@metadata$markers[[group]]$marker_plot$data$Gene))
+
   p <- sce@metadata$markers[[group]]$marker_plot
+
   plot_file_name <- paste0(group, "_markers")
+
   # save PNG
   png(file.path(folder_path, paste0(plot_file_name, ".png")),
-      width = pwidth * 12, height = pheight * 5, units = "mm", res = 600)
+    width = pwidth * 12, height = pheight * 5, units = "mm", res = 600
+  )
   print(p)
   dev.off()
 
@@ -195,14 +230,12 @@ for (group in names(sce@metadata$markers)) {
     units = "mm",
     scale = 1
   )
-
 }
 
 ### Save Marker Gene Tables
 folder_path <- file.path(getwd(), "celltype_marker_tables")
 dir.create(folder_path)
 for (group in names(sce@metadata$markers)) {
-
   marker_test_file_name <- paste0(group, "_markers_test.tsv")
   top_markers_file_name <- paste0(group, "_top_markers.tsv")
 
@@ -221,7 +254,6 @@ for (group in names(sce@metadata$markers)) {
     col.names = TRUE,
     sep = "\t"
   )
-
 }
 
 
@@ -231,5 +263,3 @@ write_sce(
   folder_path = file.path(getwd(), "final_sce")
 )
 
-##  ............................................................................
-##  Clean up                                                                ####