--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/signac-subset.R	Tue Aug 02 19:11:27 2022 +0000
@@ -0,0 +1,104 @@
+#!/usr/bin/env Rscript
+
+# Load optparse we need to check inputs
+
+suppressPackageStartupMessages(require(optparse))
+
+# Load common functions
+
+suppressPackageStartupMessages(require(workflowscriptscommon))
+
+# parse options
+
+option_list = list(
+  make_option(
+    c("--signac-object"),
+    action = "store",
+    default = NA,
+    type = 'character',
+    help = ""
+  ),
+  make_option(
+    c("--peak-region-fragments-min"),
+    action = "store",
+    default = NA,
+    type = 'character',
+    help = "."
+  ),
+  make_option(
+    c("--peak-region-fragments-max"),
+    action = "store",
+    default = NA,
+    type = 'character',
+    help = "."
+  ),
+  make_option(
+    c("--pct-reads-in-peaks"),
+    action = "store",
+    default = NA,
+    type = 'character',
+    help = "."
+  ),
+  make_option(
+    c("--blacklist-ratio"),
+    action = "store",
+    default = NA,
+    type = 'character',
+    help = "."
+  ),
+  make_option(
+    c("--nucleosome-signal"),
+    action = "store",
+    default = NA,
+    type = 'character',
+    help = "."
+  ),
+  make_option(
+    c("--tss-enrichment"),
+    action = "store",
+    default = NA,
+    type = 'character',
+    help = "."
+  ),
+  make_option(
+    c("--output-object-file"),
+    action = "store",
+    default = NA,
+    type = 'character',
+    help = "File name in which to store serialized R matrix object."
+  ) 
+)
+
+opt <- wsc_parse_args(option_list)
+
+suppressPackageStartupMessages(require(Signac))
+
+# extract gene annotations from EnsDb
+signac_object <- readRDS(file = opt$signac_object)
+
+
+## transform input parameters to numeric
+peak_region_fragments_min <- as.numeric(opt$peak_region_fragments_min)
+peak_region_fragments_max <- as.numeric(opt$peak_region_fragments_max)
+pct_reads_in_peaks_var <- as.numeric(opt$pct_reads_in_peaks)
+blacklist_ratio_var <- as.numeric(opt$blacklist_ratio)
+nucleosome_signal_var <- as.numeric(opt$nucleosome_signal)
+tss_enrichment_var <- as.numeric(opt$tss_enrichment)
+
+print("Signac object before filtering:")
+
+signac_object
+
+signac_object <- subset(signac_object, peak_region_fragments > peak_region_fragments_min)
+signac_object <- subset(signac_object, peak_region_fragments < peak_region_fragments_max)
+signac_object <- subset(signac_object, pct_reads_in_peaks > pct_reads_in_peaks_var)
+signac_object <- subset(signac_object, blacklist_ratio < blacklist_ratio_var)
+signac_object <- subset(signac_object, nucleosome_signal < nucleosome_signal_var)
+signac_object <- subset(signac_object, TSS.enrichment > tss_enrichment_var)
+
+print("Signac object after filtering:")
+
+signac_object
+
+# Output to a serialized R object
+saveRDS(signac_object, file = opt$output_object_file)