update

README.md

@@ -171,38 +171,62 @@ chip/
 
 ### 1. FastQC on raw reads
 
-`fastqc`
+```
+fastqc
+```
 
 ### 2. Alignment
 
-`bowtie2`
+```
+bowtie2-build
+
+```
 
 ### 3. Filtering
 
-`SAMtools`
+```
+SAMtools sort
+SAMtools view
 
-`picard`
+picard MarkDuplicates
+```
 
 ### 4. Alignment QC & Library Complexity
 
-`SAMtools`
+```
+SAMtool
 
-`picard`
+picard
 
-`preseq`
+preseq
+
+```
 
 ### 5. deepTools
 
-`deepTools`
+```
+deepTools
+```
 
 ### 7. phantomPeakQuals
 
-`spp`
+```
+spp
+```
 
 ### 8. Call narrow peaks (MACS2)
 
-`macs2 callpeaks`
+```
+macs2 callpeaks
+```
 
 ### 9. Create consensus peaksets
 ### 10. Annotate peaks relative to gene features (HOMER)
 ### 11. Present QC for raw read, alignment, peak-calling in MultiQC
+
+
+### 12. Plot DAG
+
+```
+snakemake --dag | dot -Tsvg > dag.svg
+```

Snakefile

@@ -41,7 +41,6 @@ CASE_BAM = expand("results/bowtie2/{sample}.sorted.bam", sample=CASES)
 
 ALL_SAMPLES = CASES + CONTROLS_UNIQUE
 ALL_BAM     = CONTROL_BAM + CASE_BAM
-ALL_FASTQ   = expand("results/fastqc/{sample}.fastq", sample = ALL_SAMPLES)
 ALL_FASTQC  = expand("results/fastqc/{sample}_fastqc.zip", sample = ALL_SAMPLES)
 ALL_INDEX = expand("results/bowtie2/{sample}_PPq20.sorted.dedup.bai", sample = ALL_SAMPLES)
 ALL_LIBCOMPLEXITY = expand("logs/{sample}.ccurve.txt", sample = ALL_SAMPLES)
@@ -54,7 +53,7 @@ ALL_PHANTOM1 = expand("results/phantomPeaks/{sample}.spp.png", sample = CASES)
 ALL_PHANTOM2 = expand("results/phantomPeaks/{sample}.spp.Rdata", sample = CASES)
 ALL_PHANTOM3 = expand("results/phantomPeaks/{sample}.spp.out", sample = CASES)
 ALL_BIGWIG = expand("results/deepTools/{sample}.SeqDepthNorm.bw", sample = ALL_SAMPLES)
-ALL_QC = ["multiQC/multiQC_log.html"]
+ALL_QC = ["results/multiqc/multiQC_log.html"]
 ALL_PEAKSxls = expand("results/macs2/{case}_vs_{control}_macs2_peaks.xls", zip, case=CASES, control=CONTROLS)
 ALL_PEAKSbed = expand("results/macs2/{case}_vs_{control}_macs2_summits.bed", zip, case=CASES, control=CONTROLS)
 ALL_PEAKSpeak = expand("results/macs2/{case}_vs_{control}_macs2_peaks.narrowPeak", zip, case=CASES, control=CONTROLS)
@@ -65,7 +64,7 @@ ALL_PEAKSpeak = expand("results/macs2/{case}_vs_{control}_macs2_peaks.narrowPeak
 
 localrules: all
 rule all:
-    input: ALL_FASTQC + ALL_BAM + ALL_INDEX + ALL_PHATOM + ALL_PEAKS + ALL_BIGWIG + ALL_inputSubtract_BIGWIG + ALL_FASTQ + ALL_FLAGSTAT + ALL_QC + ALL_SUPER
+    input: ALL_FASTQC + ALL_BAM + ALL_INDEX + ALL_LIBCOMPLEXITY + ALL_PICARDMETRIC + ALL_FINGERPRINT1 + ALL_FINGERPRINT2 + ALL_FINGERPRINT3 + ALL_FLAGSTAT + ALL_PHANTOM1 + ALL_PHANTOM2 + ALL_PHANTOM3 + ALL_BIGWIG + ALL_QC + ALL_PEAKSxls + ALL_PEAKSbed + ALL_PEAKSpeak
 
 # ===============================================================================================
 #  1. FASTQC
@@ -340,7 +339,7 @@ rule phantomPeakQuals:
     log:
         "logs/{sample}.phantomPeak"
     shell:
-        "Rscript run_spp.R -c={input} -savp={output.savp} -savd={output.savd} -out={output.out} &> {log}"
+        "Rscript scripts/run_spp.R -c={input} -savp={output.savp} -savd={output.savd} -out={output.out} &> {log}"
 
 # ===============================================================================================
 #  8. Call peaks (MACS2)
@@ -395,10 +394,10 @@ rule multiqc:
         expand("results/phantomPeaks/{sample}.spp.Rdata", sample = CASES),
         expand("results/phantomPeaks/{sample}.spp.out", sample = CASES),
         # macs2
-        expand("analysis/macs2/{case}_vs_{control}_model.r", zip, case=CASES, control=CONTROLS),
-        expand("analysis/macs2/{case}_vs_{control}_macs2_peaks.narrowPeak", zip, case=CASES, control=CONTROLS),
-        expand("analysis/macs2/{case}_vs_{control}_macs2_peaks.xls", zip, case=CASES, control=CONTROLS),
-        expand("analysis/macs2/{case}_vs_{control}_macs2_summits.bed", zip, case=CASES, control=CONTROLS)
+        expand("results/macs2/{case}_vs_{control}_model.r", zip, case=CASES, control=CONTROLS),
+        expand("results/macs2/{case}_vs_{control}_macs2_peaks.narrowPeak", zip, case=CASES, control=CONTROLS),
+        expand("results/macs2/{case}_vs_{control}_macs2_peaks.xls", zip, case=CASES, control=CONTROLS),
+        expand("results/macs2/{case}_vs_{control}_macs2_summits.bed", zip, case=CASES, control=CONTROLS)
     params:
         "results/fastqc/",
         "results/bowtie2/",
@@ -428,5 +427,13 @@ rule multiqc:
         -o results/multiqc/ \
         -n multiqc_report.html"
 
+rule plot_quals:
+    input:
+        "calls/all.vcf"
+    output:
+        "plots/quals.svg"
+    script:
+        "scripts/plot-quals.py"
+
 #Results generated by MultiQC collate pipeline QC from FastQC, TrimGalore, samtools flagstat, samtools idxstats, samtools stats,
 #picard CollectMultipleMetrics, picard MarkDuplicates, Preseq, deepTools plotProfile, deepTools plotFingerprint, phantompeakqualtools

envs/align.yaml

 channels:
     - bioconda
-    - conda-forge
 dependencies:
     - bowtie2=2.3.4.3
     - samtools=1.10-2

envs/qc.yaml

 channels:
     - bioconda
+    - conda-forge
 dependencies:
     - fastqc=0.11.8
     - samtools=1.10-2