# HG changeset patch # User jeremyjliu # Date 1415381898 18000 # Node ID 00627aa25b6c0804eb6aee6b4e4a700385443be9 # Parent ef85d5e787d7881a7b47ec705c8e043a51e07506 Update relative directories within intersect.r diff -r ef85d5e787d7 -r 00627aa25b6c region_motif_intersect.r --- a/region_motif_intersect.r Fri Nov 07 12:37:44 2014 -0500 +++ b/region_motif_intersect.r Fri Nov 07 12:38:18 2014 -0500 @@ -19,18 +19,18 @@ # Set common and data directories args <- commandArgs() workingDir = args[7] -commonDir = concat(workingDir, "/tools/my_tools") +dbDir = concat(workingDir, "/region_motif_db") dbCode = args[8] if (dbCode == "t") { - motifDB = concat(commonDir, "/region_motif_db/pouya_test_motifs.bed.bgz") + motifDB = concat(dbDir, "/pouya_test_motifs.bed.bgz") } else if (dbCode == "p") { - motifDB = concat(commonDir, "/region_motif_db/pouya_motifs.bed.bgz") + motifDB = concat(dbDir, "/pouya_motifs.bed.bgz") } else if (dbCode == "j") { - motifDB = concat(commonDir, "/region_motif_db/jaspar_jolma_motifs.bed.bgz") + motifDB = concat(dbDir, "/jaspar_jolma_motifs.bed.bgz") } else if (dbCode == "m") { - motifDB = concat(commonDir, "/region_motif_db/mm9_motifs_split.bed.bgz") + motifDB = concat(dbDir, "/mm9_motifs.bed.bgz") } else { - motifDB = concat(commonDir, "/region_motif_db/pouya_motifs.bed.bgz") + motifDB = concat(dbDir, "/pouya_motifs.bed.bgz") } # Set input and reference files, comment to toggle commmand line arguments @@ -46,16 +46,16 @@ cat("Running ... Started at:", format(startTime, "%a %b %d %X %Y"), "...\n") # Load dependencies -#cat("Loading dependencies...\n") +cat("Loading dependencies...\n") suppressPackageStartupMessages(library(Rsamtools, quietly=TRUE)) # Initializing hash table (as env) with motif names and loading tabix file -#cat("Loading motif database and initializing hash table...\n") +cat("Loading motif database and initializing hash table...\n") motifTable = new.env() motifTbx <- TabixFile(motifDB) # Loading input bed file, convert integer columns to numeric, name columns -#cat("Loading region file...\n") +cat("Loading region file...\n") regionsDF = read_bed(inBed) dfTemp = sapply(regionsDF, is.integer) regionsDF[dfTemp] = lapply(regionsDF[dfTemp], as.numeric) @@ -64,7 +64,7 @@ names(regionsDF)[names(regionsDF) == "V3"] = "end" # Filtering regions to exclude chromosomes not in motif database -#cat("Determining intersection counts...\n") +cat("Determining intersection counts...\n") motifTbxChrs = seqnamesTabix(motifTbx) regionsDFFilter = subset(regionsDF, chr %in% motifTbxChrs) @@ -76,7 +76,7 @@ regionsIntersects = scanTabix(motifTbx, param=param) # Parsing result list and updating motif count hash table -#cat("Parsing result list...\n") +cat("Parsing result list...\n") for(regionIntersects in regionsIntersects) { for(regionIntersect in strsplit(regionIntersects, " ")) { intersectMotif = strsplit(regionIntersect, "\t")[[1]][4] @@ -96,7 +96,7 @@ } # Outputting intersection counts to tab delineated file -#cat("Outputting to file...\n") +cat("Outputting to file...\n") write.table(counts, outTab, quote=FALSE, sep="\t", row.names=TRUE, col.names=FALSE) cat("Done. Job started at:", format(startTime, "%a %b %d %X %Y."), "Job ended at:", format(Sys.time(), "%a %b %d %X %Y."), "\n")