Add snakemake step to extract marker genes from the simulated data

Snakefile.py

@@ -31,11 +31,17 @@ res_dirs = expand(
   sim_par = sim_pars
 )
 
+mg_dirs = expand(
+  config["sim_data_folder"] + "{mg_type}-{sim_par}.rds",
+  sim_par = sim_pars, 
+  mg_type = config["mg_type"]
+)
+
 shell("module load r/4.0.0")
 shell("module load snakemake/5.18.1")
 
 rule all:
-  input: sim_dirs, res_dirs
+  input: sim_dirs, res_dirs, mg_dirs
 
 rule simulate_data:
   priority: 100
@@ -51,6 +57,16 @@ rule simulate_data:
 		sim_pars={input.sim_pars} wcs={wildcards}"\
 		{input.script} {log}'''
 
+rule extract_marker_genes:
+	input:  script = config["code_folder"] + "extract-marker-genes.R",
+	    sim = config["sim_data_folder"] + "{sim_id}.rds",
+	output: mgs = config["sim_data_folder"] + "{mg_type}-{sim_id}.rds"
+	log:	config["logs_folder"] + "extract_mgs-{mg_type}-{sim_id}.Rout"
+	shell:	'''R CMD BATCH --no-restore --no-save\
+		"--args sim={input.sim} mg_type={wildcards.mg_type}\
+		mgs={output.mgs}"\
+		{input.script} {log}'''
+
 rule run_method:
 	priority: 99
 	input:	script = config["code_folder"] + "run-method.R",
@@ -61,7 +77,7 @@ rule run_method:
 	log:	config["logs_folder"] + "run_method-{sim_id}-{method_id}.Rout"
 	resources: 
 	  mem_mb=lambda wildcards, attempt: (attempt ** 2) * 3000,
-	  time_min=lambda wildcards, attempt: [20, 60, 120][attempt - 1]
+	  time_min=lambda wildcards, attempt: [20, 60, 360][attempt - 1]
 	shell:	'''R CMD BATCH --no-restore --no-save\
 		"--args sim={input.sim} fun={input.fun}\
 		meth_pars={input.method_pars} res={output.res}"\

code/extract-marker-genes.R

+# Script to extract true marker genes from the simulation.
+
+# Get the passed arguments from the command line via Snakemake.
+args <- R.utils::commandArgs(trailingOnly = TRUE, asValues = TRUE)
+
+# Debugging.
+print(args)
+
+source(here::here("code", "find-marker-genes.R"))
+
+sim <- readRDS(args$sim)
+mg_type <- args$mg_type
+
+stopifnot(mg_type %in% c("umg", "sumg"))
+
+if (mg_type == "umg") {
+  mgs <- find_unique_marker_genes(sim)
+} else {
+  mgs <- find_semi_unique_marker_genes(sim)
+}
+
+saveRDS(mgs, args$mgs)
+

config.yaml

 R: "/opt/R/4.0.3/bin/R"
 
+mg_type: ["umg", "sumg"]
+
 code_folder: "code/"
 logs_folder: "logs/"