From 93d4ec9331cb3f2987047155ff4c9a4efa5bf89f Mon Sep 17 00:00:00 2001
From: Beibei Chen <beibei.chen@utsouthwestern.edu>
Date: Tue, 27 Dec 2016 16:08:21 -0600
Subject: [PATCH] deeptools script works
---
workflow/scripts/runDeepTools.py | 24 ++++++++++++++----------
workflow/scripts/runDiffBind.R | 2 +-
2 files changed, 15 insertions(+), 11 deletions(-)
diff --git a/workflow/scripts/runDeepTools.py b/workflow/scripts/runDeepTools.py
index 091b59f..dad0cec 100644
--- a/workflow/scripts/runDeepTools.py
+++ b/workflow/scripts/runDeepTools.py
@@ -50,22 +50,26 @@ def run_signal(files, labels, genome):
gene_bed = "/project/BICF/BICF_Core/bchen4/chipseq_analysis/test/genome/"+genome+"/gene.bed"
bw_commands = []
for f in files:
- bw_commands.append("bamCoverage -b "+f+" -o "+f.replace("bam","bw"))
- work_pool = Pool(min(len(files), 12))
- work_pool.map(bam2bw_wrapper, bw_commands)
- work_pool.close()
- work_pool.join()
+ bw_commands.append("bamCoverage -bs 10 -b "+f+" -o "+f.replace("bam","bw"))
+ #work_pool = Pool(min(len(files), 12))
+ #work_pool.map(bam2bw_wrapper, bw_commands)
+ #work_pool.close()
+ #work_pool.join()
- #cm_command = "computeMatrix scale-regions "
+ cm_command = "computeMatrix scale-regions -R "+gene_bed+" -a 3000 -b 3000 --regionBodyLength 5000 --skipZeros -S *.bw -o samples.deeptools_generegionscalematrix.gz"
+ p = subprocess.Popen(cm_command, shell=True)
+ p.communicate()
+ hm_command = "plotHeatmap -m samples.deeptools_generetionscalematrix.gz -out samples.deeptools_readsHeatmap.png"
+ p = subprocess.Popen(hm_command, shell=True)
+ p.communicate()
def run(dfile,genome):
- #Get annotation address
- tss_bed = "/project/BICF/BICF_Core/bchen4/chipseq_analysis/test/genome/"+genome+"/tss.bed"
#parse dfile, suppose data files are the same folder as design file
dfile = pd.read_csv(dfile)
#QC: multiBamSummary and plotCorrelation
- run_qc(dfile['bamReads'], dfile['bamControl'], dfile['SampleID'])
-
+ #run_qc(dfile['bamReads'], dfile['bamControl'], dfile['SampleID'])
+ #signal plots
+ run_signal(dfile['bamReads'],dfile['SampleID'],genome)
def main():
argparser = prepare_argparser()
diff --git a/workflow/scripts/runDiffBind.R b/workflow/scripts/runDiffBind.R
index 2eb624b..c4cffa0 100644
--- a/workflow/scripts/runDiffBind.R
+++ b/workflow/scripts/runDiffBind.R
@@ -2,7 +2,7 @@ library("DiffBind")
#build dba object from sample sheet and do analysis
data <- dba(sampleSheet="samplesheet.csv")
-data <- dba.count(data, summits=250)
+data <- dba.count(data)
data <- dba.contrast(data, minMembers = 2, categories=DBA_CONDITION)
data <- dba.analyze(data)
--
GitLab