update scripts

Author: kcotto

scripts/compare_junctions_hist_v2.R

@@ -5,27 +5,27 @@
 # load libraries
 library(data.table)
 library(plyr)
+library(tidyverse)
 
 debug = F
 
-system.time({
-if (debug){
-  tag = paste("_", "default", sep="")
-} else {
-  # get options tag
-  args = commandArgs(trailingOnly = TRUE)
-  tag = args[1]
-  input_file = args[2]
-  if ( substr(tag, 2, 3) == "--"){
-    stop("Please specify an option tag (e.g. \"default\", \"i20e5\")")
-  }
-}
- 
-# tag = 'I'
-# input_file = '~/Desktop/CHOL/all_splicing_variants_I.bed'
+# system.time({
+# if (debug){
+#   tag = paste("_", "default", sep="")
+# } else {
+#   # get options tag
+#   args = commandArgs(trailingOnly = TRUE)
+#   tag = args[1]
+#   input_file = args[2]
+#   if ( substr(tag, 2, 3) == "--"){
+#     stop("Please specify an option tag (e.g. \"default\", \"i20e5\")")
+#   }
+# }
+
+tag = 'I'
+input_file = '~/Desktop/CHOL/all_splicing_variants_I.bed'
 
 # All splicing relevant variants (union of rows from variants.bed files; add column with comma-separated list of sample names)
-input_file = args[2]
 all_splicing_variants = unique(data.table::fread(input_file), sep = '\t', header = T, stringsAsFactors = FALSE)
 colnames(all_splicing_variants) <- c("chrom", "start", "end", "samples")
 
@@ -214,6 +214,6 @@ all_splicing_variants <- as.data.table(all_splicing_variants)
 regtools_data = regtools_data %>% distinct()
 
 
-write.table(regtools_data, file=paste(input_file, "_out.tsv", sep=""), quote=FALSE, sep='\t', row.names = F)
+write.table(regtools_data, file=paste(input_file, "_out_test.tsv", sep=""), quote=FALSE, sep='\t', row.names = F)
 
-})
+# })

scripts/stats_wrapper.py

@@ -17,22 +17,28 @@
 tag = args.tag
 cwd = os.getcwd()
 
-lines_per_file = 50000
+lines_per_file = 1000
 smallfile = None
+num_small_file = 0
 with open(f'all_splicing_variants_{tag}.bed', 'r') as bigfile:
+    num_small_file +=1
     for lineno, line in enumerate(bigfile):
         if lineno % lines_per_file == 0:
             if smallfile:
                 smallfile.close()
             small_filename = 'small_file_{}.txt'.format(lineno + lines_per_file)
             smallfile = open(small_filename, "w")
+            if num_small_file > 1:
+                smallfile.write(f'chrom\tstart\tend\tsamples')
         smallfile.write(line)
     if smallfile:
+        num_small_file += 1
         smallfile.close()
 #get chunks
 files = glob.glob('small_file_*')
+files.sort()
 for file in files:
-    subprocess.run(f'Rscript --vanilla /home/ec2-user/workspace/regtools/scripts/compare_junctions_hist_v2.R {tag} ~{cwd}/{file}', shell=True, check=True)
+    subprocess.run(f'Rscript --vanilla ~/Git/regtools/scripts/compare_junctions_hist_v2.R {tag} {file}', shell=True, check=False)
 output_files = glob.glob("*_out.tsv")
 output_files.sort()  # glob lacks reliable ordering, so impose your own if output order matters
 with open(f'junction_pvalues_{tag}.tsv', 'wb') as outfile: