Skip to content

Commit

Permalink
fix: review suggestions
Browse files Browse the repository at this point in the history
  • Loading branch information
@rroutsong
rroutsong committed Dec 30, 2024
1 parent 4b32239 commit c8e2b06
Show file tree
Hide file tree
Showing 4 changed files with 58 additions and 61 deletions.
54 changes: 26 additions & 28 deletions workflow/rules/paired/peakcall.smk
View file
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,6 @@ tmpdir = config['options']['tmp_dir']
rule MACS2_broad:
input:
chip = join(bam_dir, "{name}.Q5DD.bam"),
txt = join(ppqt_dir, "{name}.Q5DD.ppqt.txt"),
c_option = lambda w: join(bam_dir, f"{chip2input[w.name]}.Q5DD.bam")
if chip2input[w.name] else [],
output:
Expand All @@ -37,6 +36,32 @@ rule MACS2_broad:
"""


rule MACS2_narrow:
input:
chip = join(bam_dir, "{name}.Q5DD.bam"),
c_option = lambda w: join(bam_dir, f"{chip2input[w.name]}.Q5DD.bam")
if chip2input[w.name] else [],
output:
join(macsN_dir, "{name}", "{name}_peaks.narrowPeak"),
params:
rname = 'MACS2_narrow',
gsize = config['references'][genome]['EFFECTIVEGENOMESIZE'],
macsver = config['tools']['MACSVER'],
flag = lambda w: "-c" if chip2input[w.name] else "",
shell:
"""
module load {params.macsver};
macs2 callpeak \\
-t {input.chip} {params.flag} {input.c_option} \\
-g {params.gsize} \\
-n {wildcards.name} \\
--outdir {macsN_dir}/{wildcards.name} \\
-q 0.01 \\
--keep-dup="all" \\
-f "BAMPE"
"""


rule SICER:
input:
chip = lambda w: join(bam_dir, w.name + ".Q5DD.bam"),
Expand Down Expand Up @@ -111,31 +136,4 @@ rule SICER:
mv ${{tmp}}/{params.name}.Q5DD-W300-G600.scoreisland {params.this_sicer_dir}
fi
"""


rule MACS2_narrow:
input:
chip = join(bam_dir, "{name}.Q5DD.bam"),
txt = join(ppqt_dir, "{name}.Q5DD.ppqt.txt"),
c_option = lambda w: join(bam_dir, f"{chip2input[w.name]}.Q5DD.bam")
if chip2input[w.name] else [],
output:
join(macsN_dir, "{name}", "{name}_peaks.narrowPeak"),
params:
rname = 'MACS2_narrow',
gsize = config['references'][genome]['EFFECTIVEGENOMESIZE'],
macsver = config['tools']['MACSVER'],
flag = lambda w: "-c" if chip2input[w.name] else "",
shell:
"""
module load {params.macsver};
macs2 callpeak \\
-t {input.chip} {params.flag} {input.c_option} \\
-g {params.gsize} \\
-n {wildcards.name} \\
--outdir {macsN_dir}/{wildcards.name} \\
-q 0.01 \\
--keep-dup="all" \\
-f "BAMPE"
"""
56 changes: 29 additions & 27 deletions workflow/rules/single/peakcall.smk
View file
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,35 @@ rule MACS2_broad:
"""


rule MACS2_narrow:
input:
chip = join(bam_dir, "{name}.Q5DD.bam"),
txt = join(ppqt_dir, "{name}.Q5DD.ppqt.txt"),
c_option = lambda w: join(bam_dir, f"{chip2input[w.name]}.Q5DD_tagAlign.gz")
if w.name and chip2input[w.name] else [],
output:
join(macsN_dir, "{name}", "{name}_peaks.narrowPeak"),
params:
rname = 'MACS2_narrow',
gsize = config['references'][genome]['EFFECTIVEGENOMESIZE'],
macsver = config['tools']['MACSVER'],
flag = lambda w: "-c" if chip2input[w.name] else "",
shell:
"""
module load {params.macsver};
ppqt_len=$(awk '{{print $1}}' {input.txt})
macs2 callpeak \\
-t {input.chip} {params.flag} {input.c_option} \\
-g {params.gsize} \\
-n {wildcards.name} \\
--outdir {macsN_dir}/{wildcards.name} \\
-q 0.01 \\
--keep-dup="all" \\
--nomodel \\
--extsize $ppqt_len
"""


rule SICER:
input:
chip = lambda w: join(bam_dir, w.name + ".Q5DD_tagAlign.gz"),
Expand Down Expand Up @@ -111,30 +140,3 @@ rule SICER:
"""


rule MACS2_narrow:
input:
chip = join(bam_dir, "{name}.Q5DD.bam"),
txt = join(ppqt_dir, "{name}.Q5DD.ppqt.txt"),
c_option = lambda w: join(bam_dir, f"{chip2input[w.name]}.Q5DD_tagAlign.gz")
if w.name and chip2input[w.name] else [],
output:
join(macsN_dir, "{name}", "{name}_peaks.narrowPeak"),
params:
rname = 'MACS2_narrow',
gsize = config['references'][genome]['EFFECTIVEGENOMESIZE'],
macsver = config['tools']['MACSVER'],
flag = lambda w: "-c" if chip2input[w.name] else "",
shell:
"""
module load {params.macsver};
ppqt_len=$(awk '{{print $1}}' {input.txt})
macs2 callpeak \\
-t {input.chip} {params.flag} {input.c_option} \\
-g {params.gsize} \\
-n {wildcards.name} \\
--outdir {macsN_dir}/{wildcards.name} \\
-q 0.01 \\
--keep-dup="all" \\
--nomodel \\
--extsize $ppqt_len
"""
1 change: 0 additions & 1 deletion workflow/rules/single/qc.smk
View file
Original file line number Diff line number Diff line change
Expand Up @@ -109,7 +109,6 @@ rule kraken:
--output {output.krakenout} \\
--gzip-compressed \\
{input.fq1}
fi
# Generate Krona Report
cut -f2,3 {output.krakenout} | \\
Expand Down
8 changes: 3 additions & 5 deletions workflow/rules/single/trim_align_dedup.smk
View file
Original file line number Diff line number Diff line change
Expand Up @@ -171,7 +171,6 @@ rule dedup:
out_bam = join(bam_dir, "{name}.Q5DD.bam"),
flagstat = join(bam_dir, "{name}.Q5DD.bam.flagstat"),
idxstat = join(bam_dir, "{name}.Q5DD.bam.idxstat"),
bwadups = join(bam_dir, "{name}.bwa.Q5.duplic"),
tagalign = join(bam_dir, "{name}.Q5DD_tagAlign.gz")
params:
rname = 'dedup',
Expand Down Expand Up @@ -281,7 +280,7 @@ rule bam2bw:
"""
input:
bam = lambda w: join(bam_dir, f"{w.name}.{w.ext}.bam"),
ppqt = lambda w: join(ppqt_dir, f"{w.name}.{w.ext}.ppqt.txt")
ppqt = lambda w: join(ppqt_dir, f"{w.name}.{w.ext}.ppqt.txt"),
output:
outbw = join(bw_dir, "{name}.{ext}.RPGC.bw"),
params:
Expand All @@ -299,8 +298,7 @@ rule bam2bw:
tmp=$(mktemp -d -p "{params.tmpdir}")
trap 'rm -rf "${{tmp}}"' EXIT
bam_cov_option=--centerReads
echo "printing out value of bam-cov-option $bam_cov_option"
ppqt_len=$(awk '{{print $1}}' {input.ppqt})
bamCoverage \\
--bam {input.bam} \\
Expand All @@ -310,5 +308,5 @@ rule bam2bw:
--numberOfProcessors {threads} \\
--normalizeUsing RPGC \\
--effectiveGenomeSize {params.effectivegenomesize} \\
${{bam_cov_option}};
-e ${{ppqt_len}}
"""

0 comments on commit c8e2b06

Please sign in to comment.