Skip to content
Snippets Groups Projects
Commit 68fe00c7 authored by Laura Cook's avatar Laura Cook
Browse files

altered to run on 10M subsampled BAM files

parent d7d436c2
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
dunnart/.DS_Store
+ 0
− 0
View file @ 68fe00c7
No preview for this file type
dunnart/Snakefile
+ 45
− 47
View file @ 68fe00c7
......@@ -88,21 +88,20 @@ rule all:
# expand("results_10M/logs/{sample}_filt_15Mreads.SE.spp.log", sample=all_samples),
# expand("results_10M/qc/{sample}_filt_15Mreads.SE.cc.qc", sample=all_samples),
# expand("results_10M/qc/{sample}_filt_15Mreads.SE.cc.plot.pdf", sample=all_samples),
expand("results_10M/macs2/{case}_vs_{control}_macs2_peaks.narrowPeak", zip, case=IPS, control=INPUTS),
expand("results_10M/macs2/{case}_vs_{control}_macs2_peaks.xls", zip, case=IPS, control=INPUTS),
expand("results_10M/macs2/{case}_vs_{control}_macs2_summits.bed", zip, case=IPS, control=INPUTS),
expand("results_10M/macs2/{case}_vs_{control}_macs2_default_peaks.narrowPeak", zip, case=IPS, control=INPUTS),
expand("results_10M/macs2/{case}_vs_{control}_macs2_default_peaks.xls", zip, case=IPS, control=INPUTS),
expand("results_10M/macs2/{case}_vs_{control}_macs2_default_summits.bed", zip, case=IPS, control=INPUTS),
#expand("results_10M/qc/{case}-vs-{control}-narrowpeak-count_mqc.json", zip, case=IPS, control=INPUTS),
# expand("results_10M/bowtie2/{case}.bedpe", case=IPS),
# expand("results_10M/ggplot(promoter_annot, aes(x=distanceToTSS, y=width)) +
#expand("results_10M/logs/{case}.bamToBed", case=IPS),
# expand("results_10M/qc/{case}_vs_{control}.frip.txt", case=IPS, control=INPUTS),
# "results_10M/macs2/H3K4me3_pooled_macs2_peaks.narrowPeak",
# "results_10M/macs2/H3K27ac_pooled_macs2_peaks.narrowPeak",
# "results_10M/macs2/H3K4me3_overlap.narrowPeak",
# "results_10M/macs2/H3K27ac_overlap.narrowPeak",
# "results_10M/qc/H3K4me3_overlap.frip",
# "results_10M/qc/H3K27ac_overlap.frip"
# directory("results_10M/multiqc/multiqc_report_data/"),
expand("results_10M/bowtie2/{case}.bedpe", case=IPS),
expand("results_10M/logs/{case}.bamToBed", case=IPS),
expand("results_10M/qc/{case}_vs_{control}.frip_default.txt", case=IPS, control=INPUTS),
"results_10M/macs2/H3K4me3_pooled_macs2_default_peaks.narrowPeak",
"results_10M/macs2/H3K27ac_pooled_macs2_default_peaks.narrowPeak",
"results_10M/macs2/H3K4me3_overlap_default.narrowPeak",
"results_10M/macs2/H3K27ac_overlap_default.narrowPeak",
"results_10M/qc/H3K4me3_overlap_default.frip",
"results_10M/qc/H3K27ac_overlap_default.frip"
#directory("results_10M/multiqc/multiqc_report_data/"),
# "results_10M/multiqc/multiqc_report.html"
# ===============================================================================================
# 1. FASTQC
......@@ -518,17 +517,17 @@ rule call_peaks_macs2:
control = "results_10M/bowtie2/{control}_PPq30.sorted.dedup.bam",
case = "results_10M/bowtie2/{case}_PPq30.sorted.dedup.bam"
output:
"results_10M/macs2/{case}_vs_{control}_macs2_peaks.xls",
"results_10M/macs2/{case}_vs_{control}_macs2_summits.bed",
"results_10M/macs2/{case}_vs_{control}_macs2_peaks.narrowPeak",
"results_10M/macs2/{case}_vs_{control}_macs2_default_peaks.xls",
"results_10M/macs2/{case}_vs_{control}_macs2_default_summits.bed",
"results_10M/macs2/{case}_vs_{control}_macs2_default_peaks.narrowPeak",
log:
"results_10M/logs/{case}_vs_{control}_call_peaks_macs2.log"
"results_10M/logs/{case}_vs_{control}_call_peaks_macs2_default.log"
params:
name = "{case}_vs_{control}_macs2_P10-2",
name = "{case}_vs_{control}_macs2_default",
shell:
" macs2 callpeak -f BAMPE -t {input.case} \
-c {input.control} --keep-dup all \
--outdir results_10M/macs2/ -p 0.01 \
--outdir results_10M/macs2/ \
-n {params.name} \
-g 2740338543 2> {log} "
......@@ -538,22 +537,22 @@ rule call_peaks_macs2_pooled_replicates:
H3K27ac = "results_10M/bowtie2/H3K27ac_pooled_PPq30.sorted.dedup.bam",
input = "results_10M/bowtie2/input_pooled_PPq30.sorted.dedup.bam"
output:
"results_10M/macs2/H3K4me3_pooled_macs2_peaks.xls",
"results_10M/macs2/H3K4me3_pooled_macs2_summits.bed",
"results_10M/macs2/H3K4me3_pooled_macs2_peaks.narrowPeak",
"results_10M/macs2/H3K27ac_pooled_macs2_peaks.xls",
"results_10M/macs2/H3K27ac_pooled_macs2_summits.bed",
"results_10M/macs2/H3K27ac_pooled_macs2_peaks.narrowPeak"
"results_10M/macs2/H3K4me3_pooled_macs2_default_peaks.xls",
"results_10M/macs2/H3K4me3_pooled_macs2_default_summits.bed",
"results_10M/macs2/H3K4me3_pooled_macs2_default_peaks.narrowPeak",
"results_10M/macs2/H3K27ac_pooled_macs2_default_peaks.xls",
"results_10M/macs2/H3K27ac_pooled_macs2_default_summits.bed",
"results_10M/macs2/H3K27ac_pooled_macs2_default_peaks.narrowPeak"
log:
H3K4me3 ="results_10M/logs/H3K4me3_pooled_call_peaks_macs2.log",
H3K27ac ="results_10M/logs/H3K27ac_pooled_call_peaks_macs2.log"
H3K4me3 ="results_10M/logs/H3K4me3_pooled_call_peaks_macs2_default.log",
H3K27ac ="results_10M/logs/H3K27ac_pooled_call_peaks_macs2_default.log"
params:
H3K4me3 = "H3K4me3_pooled_macs2",
H3K27ac = "H3K27ac_pooled_macs2"
H3K4me3 = "H3K4me3_pooled_macs2_default",
H3K27ac = "H3K27ac_pooled_macs2_default"
run:
shell(" macs2 callpeak -f BAMPE -t {input.H3K4me3} \
-c {input.input} --keep-dup all \
--outdir results_10M/macs2/ -p 0.01 \
--outdir results_10M/macs2/ \
-n {params.H3K4me3} \
-g 2740338543 2> {log.H3K4me3} ")
shell("macs2 callpeak -f BAMPE -t {input.H3K27ac} \
......@@ -598,9 +597,9 @@ rule bamToBed:
rule frip:
input:
bed = "results_10M/bowtie2/{case}.bedpe",
peak = "results_10M/macs2/{case}_vs_{control}_macs2_peaks.narrowPeak"
peak = "results_10M/macs2/{case}_vs_{control}_macs2_default_peaks.narrowPeak"
output:
"results_10M/qc/{case}_vs_{control}.frip.txt"
"results_10M/qc/{case}_vs_{control}.frip_default.txt"
shell:
"python2.7 scripts/encode_frip.py {input.bed} {input.peak} > {output}"
......@@ -614,22 +613,22 @@ rule frip:
rule overlap_peaks_H3K4me3:
input:
peak1="results_10M/macs2/A-2_H3K4me3_vs_A-1_input_macs2_peaks.narrowPeak",
peak2="results_10M/macs2/B-2_H3K4me3_vs_B-1_input_macs2_peaks.narrowPeak",
pooled="results_10M/macs2/H3K4me3_pooled_macs2_peaks.narrowPeak"
peak1="results_10M/macs2/A-2_H3K4me3_vs_A-1_input_macs2_default_peaks.narrowPeak",
peak2="results_10M/macs2/B-2_H3K4me3_vs_B-1_input_macs2_default_peaks.narrowPeak",
pooled="results_10M/macs2/H3K4me3_pooled_macs2_default_peaks.narrowPeak"
output:
"results_10M/macs2/H3K4me3_overlap.narrowPeak"
"results_10M/macs2/H3K4me3_overlap_default.narrowPeak"
shell:
"python2.7 scripts/overlap_peaks.py {input.peak1} {input.peak2} {input.pooled} {output}"
rule overlap_peaks_H3K27ac:
input:
peak1="results_10M/macs2/A-3_H3K27ac_vs_A-1_input_macs2_peaks.narrowPeak",
peak2="results_10M/macs2/B-3_H3K27ac_vs_B-1_input_macs2_peaks.narrowPeak",
pooled="results_10M/macs2/H3K27ac_pooled_macs2_peaks.narrowPeak"
peak1="results_10M/macs2/A-3_H3K27ac_vs_A-1_input_macs2_default_peaks.narrowPeak",
peak2="results_10M/macs2/B-3_H3K27ac_vs_B-1_input_macs2_default_peaks.narrowPeak",
pooled="results_10M/macs2/H3K27ac_pooled_macs2_default_peaks.narrowPeak"
output:
"results_10M/macs2/H3K27ac_overlap.narrowPeak"
"results_10M/macs2/H3K27ac_overlap_default.narrowPeak"
shell:
"python2.7 scripts/overlap_peaks.py {input.peak1} {input.peak2} {input.pooled} {output}"
......@@ -638,17 +637,16 @@ rule overlap_frip:
input:
H3K4me3bam = "results_10M/bowtie2/H3K4me3_pooled_PPq30.sorted.dedup.bam",
H3K27acbam = "results_10M/bowtie2/H3K27ac_pooled_PPq30.sorted.dedup.bam",
H3K4me3bed = "results_10M/macs2/H3K4me3_overlap.narrowPeak",
H3K27acbed = "results_10M/macs2/H3K27ac_overlap.narrowPeak"
H3K4me3bed = "results_10M/macs2/H3K4me3_overlap_default.narrowPeak",
H3K27acbed = "results_10M/macs2/H3K27ac_overlap_default.narrowPeak"
output:
H3K4me3frip = "results_10M/qc/H3K4me3_overlap.frip",
H3K27acfrip = "results_10M/qc/H3K27ac_overlap.frip"
H3K4me3frip = "results_10M/qc/H3K4me3_overlap_default.frip",
H3K27acfrip = "results_10M/qc/H3K27ac_overlap_default.frip"
run:
shell("python2.7 scripts/encode_frip.py {input.H3K4me3bam} {input.H3K4me3bed} > {output.H3K4me3frip}")
shell("python2.7 scripts/encode_frip.py {input.H3K27acbam} {input.H3K27acbed} > {output.H3K27acfrip}")
# ===============================================================================================
# 10. Combine all QC into multiqc output
tyui# 10. Combine all QC into multiqc output
# ===============================================================================================
#
# rule multiqc:
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment