Merge pull request #754 from maxplanck-ie/develop

Develop - for the 2.5 release

.azure-pipelines/testCreateEnv.py

@@ -0,0 +1,38 @@
+#!/usr/bin/env python
+
+import subprocess as sp
+import sys
+import argparse
+
+sys.path.append("..")
+from snakePipes import common_functions as cf
+
+
+def parse_args():
+    parser = argparse.ArgumentParser()
+    required = parser.add_argument_group('Required Arguments')
+    required.add_argument("-d",
+                          dest="envs_chunk",
+                          help="chunk of env to be used when running `snakePipes createEnvs`",
+                          required=True)
+    return parser
+
+
+parser = parse_args()
+args = parser.parse_args()
+
+keys = []
+length = len(cf.set_env_yamls().items())
+index = length / 3
+index = round(index)
+if args.envs_chunk == str(1):
+    envs = range(0, index)
+elif args.envs_chunk == str(2):
+    envs = range(index, 2 * index)
+elif args.envs_chunk == str(3):
+    envs = range(2 * index, length)
+
+for i in range(index):
+    keys.append(list(cf.set_env_yamls())[i])
+envs_to_test = " ".join(keys)
+sp.check_output("snakePipes createEnvs --only {} --force".format(envs_to_test), shell=True)

.ci_stuff/test_dag.sh

@@ -22,7 +22,7 @@ touch SE_input/sample1_R1.fastq.gz \
       SE_input/sample5_R1.fastq.gz \
       SE_input/sample6_R1.fastq.gz
 # Needed by ChIP and ATAC workflows
-mkdir -p BAM_input/deepTools_qc/bamPEFragmentSize BAM_input/filtered_bam BAM_input/Sambamba BAM_input/bamCoverage 
+mkdir -p BAM_input/deepTools_qc/bamPEFragmentSize BAM_input/filtered_bam BAM_input/Sambamba BAM_input/bamCoverage
 touch BAM_input/sample1.bam \
       BAM_input/sample2.bam \
       BAM_input/sample3.bam \
@@ -53,7 +53,7 @@ touch BAM_input/sample1.bam \
       BAM_input/bamCoverage/sample3.filtered.seq_depth_norm.bw \
       BAM_input/bamCoverage/sample4.filtered.seq_depth_norm.bw \
       BAM_input/bamCoverage/sample5.filtered.seq_depth_norm.bw \
-      BAM_input/bamCoverage/sample6.filtered.seq_depth_norm.bw 
+      BAM_input/bamCoverage/sample6.filtered.seq_depth_norm.bw
 
 mkdir -p allelic_BAM_input/allelic_bams allelic_BAM_input/filtered_bam  allelic_BAM_input/deepTools_qc/bamPEFragmentSize allelic_BAM_input/Sambamba allelic_BAM_input/bamCoverage/allele_specific
 touch allelic_BAM_input/allelic_bams/sample1.genome1.sorted.bam \
@@ -104,7 +104,7 @@ touch allelic_BAM_input/allelic_bams/sample1.genome1.sorted.bam \
       allelic_BAM_input/bamCoverage/allele_specific/sample3.genome1.seq_depth_norm.bw \
       allelic_BAM_input/bamCoverage/allele_specific/sample4.genome1.seq_depth_norm.bw \
       allelic_BAM_input/bamCoverage/allele_specific/sample5.genome1.seq_depth_norm.bw \
-      allelic_BAM_input/bamCoverage/allele_specific/sample6.genome1.seq_depth_norm.bw 
+      allelic_BAM_input/bamCoverage/allele_specific/sample6.genome1.seq_depth_norm.bw
 mkdir -p output
 touch /tmp/genes.gtf /tmp/genome.fa /tmp/genome.fa.fai /tmp/rmsk.txt /tmp/genes.bed /tmp/spikein_genes.gtf
 mkdir -p allelic_input
@@ -156,29 +156,29 @@ if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 1297 ]; then exit 1 ; fi
 
 # ChIP-seq
 WC=`ChIP-seq -d BAM_input --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp" .ci_stuff/organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
-if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 401 ]; then exit 1 ; fi
+if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 403 ]; then exit 1 ; fi
 WC=`ChIP-seq -d BAM_input --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp" --peakCaller Genrich .ci_stuff/organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
 if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 403 ]; then exit 1 ; fi
 WC=`ChIP-seq -d BAM_input --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp" --singleEnd .ci_stuff/organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
 if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 401 ]; then exit 1 ; fi
 WC=`ChIP-seq -d BAM_input --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp" --bigWigType log2ratio .ci_stuff/organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
-if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 363 ]; then exit 1 ; fi
+if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 365 ]; then exit 1 ; fi
 # fromBAM
 WC=`ChIP-seq -d outdir --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp" --fromBAM BAM_input/filtered_bam/ .ci_stuff/organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
-if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 715 ]; then exit 1 ; fi
+if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 717 ]; then exit 1 ; fi
 # spikein
 WC=`ChIP-seq -d BAM_input --useSpikeInForNorm --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp"  .ci_stuff/spikein_organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
-if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 433 ]; then exit 1 ; fi
+if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 435 ]; then exit 1 ; fi
 # fromBAM and spikein
 WC=`ChIP-seq -d outdir --useSpikeInForNorm --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp" --fromBAM BAM_input/filtered_bam/ .ci_stuff/spikein_organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
-if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 591 ]; then exit 1 ; fi
+if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 593 ]; then exit 1 ; fi
 WC=`ChIP-seq -d outdir --useSpikeInForNorm --getSizeFactorsFrom TSS --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp" --fromBAM BAM_input/filtered_bam/ .ci_stuff/spikein_organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
-if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 591 ]; then exit 1 ; fi
+if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 593 ]; then exit 1 ; fi
 WC=`ChIP-seq -d outdir --useSpikeInForNorm --getSizeFactorsFrom input --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp" --fromBAM BAM_input/filtered_bam/ .ci_stuff/spikein_organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
-if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 579 ]; then exit 1 ; fi
+if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 581 ]; then exit 1 ; fi
 # allelic
 WC=`ChIP-seq -d allelic_BAM_input --sampleSheet .ci_stuff/test_sampleSheet.tsv --snakemakeOptions " --dryrun --conda-prefix /tmp"  .ci_stuff/organism.yaml .ci_stuff/ChIP.sample_config.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`
-if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 331 ]; then exit 1 ; fi
+if [ ${PIPESTATUS[0]} -ne 0 ] || [ $WC -ne 333 ]; then exit 1 ; fi
 
 # mRNA-seq
 WC=`mRNA-seq -i PE_input -o output --snakemakeOptions " --dryrun --conda-prefix /tmp" .ci_stuff/organism.yaml | tee >(cat 1>&2) | grep -v "Conda environment" | sed '/^\s*$/d' | wc -l`

azure-pipelines.yml

@@ -17,7 +17,7 @@ jobs:
   - template: .azure-pipelines/setup.yml
   - bash: |
       source activate foo
-      flake8 --ignore=E501,E722 --exclude docs/conf.py .
+      flake8 --ignore=E501,E722,E402 --exclude docs/conf.py .
     displayName: flake8
 
 - job: 'CI_tests'
@@ -70,10 +70,23 @@ jobs:
   - bash: echo "##vso[task.prependpath]/usr/share/miniconda/bin"
     displayName: Add conda to PATH
   - template: .azure-pipelines/setup.yml
+  - task: UsePythonVersion@0
+    displayName: 'Use Python 3.7'
+    inputs:
+      versionSpec: 3.7
   - bash: |
       source activate foo
-      snakePipes createEnvs
-    displayName: "createEnvs"
+      python .azure-pipelines/testCreateEnv.py -d 1
+    displayName: "createEnvs1"
+  - bash: |
+      source activate foo
+      python .azure-pipelines/testCreateEnv.py -d 2
+    displayName: "createEnvs2"
+  - bash: |
+      source activate foo
+      python .azure-pipelines/testCreateEnv.py -d 3
+    displayName: "createEnvs3"
+
 
 - job: 'envs_OSX'
   pool:
@@ -95,8 +108,20 @@ jobs:
         sudo chown 501:20 /Users/vsts/.conda/pkgs/urls.txt
       fi
     displayName: Fix OSX permissions
-  - template: .azure-pipelines/macosx_setup.yml
+  - template: .azure-pipelines/setup.yml
+  - task: UsePythonVersion@0
+    displayName: 'Use Python 3.7'
+    inputs:
+      versionSpec: 3.7
+  - bash: |
+      source activate foo
+      sudo -E python .azure-pipelines/testCreateEnv.py -d 1
+    displayName: "createEnvs1"
+  - bash: |
+      source activate foo
+      sudo -E python .azure-pipelines/testCreateEnv.py -d 2
+    displayName: "createEnvs2"
   - bash: |
       source activate foo
-      sudo snakePipes createEnvs
-    displayName: "createEnvs"
+      sudo -E python .azure-pipelines/testCreateEnv.py -d 3
+    displayName: "createEnvs3"

conda-recipe/meta.yaml

@@ -1,6 +1,6 @@
 package:
   name: snakepipes
-  version: 2.4.3
+  version: 2.5.0
 
 source:
   path: ../

docs/content/News.rst

@@ -1,6 +1,14 @@
 snakePipes News
 ===============
 
+snakePipes 2.5.0
+----------------
+
+* Added tbb = 2020.2 to dna_mapping, creatIndices and rnae_seq yaml file
+* Added cut and tag parameters to macs2 and bowtei2. The parameters have been used in Kaya-Okur et al. 2019 and can be called by using --cut_n_tag.
+* Updated azure tests. python37 create envs constantly failing due its long run time. This test is now split into smaller chunks. set_macos is removed since it was completely redundant with the set.yaml
+
+
 snakePipes 2.4.3
 ----------------
 
@@ -89,7 +97,7 @@ snakePipes 2.2.0
 * Added STARsolo report for all samples to STARsolo output folder
 * FASTQ1 and FASTQ2 are not localrules anymore due to buggy logging
 * Included optional differential splicing analysis using rmats within mRNA-seq workflow
-* Symlinks in the output path are relative 
+* Symlinks in the output path are relative
 * Increased BBmap version
 * Increased STAR version to 2.7.4a in scRNAseq, noncoding-RNA-seq and mRNA-seq workflows
 * Fixed snakemake version at 5.18.0 due to a bug in DAG handling

snakePipes/__init__.py

@@ -1 +1 @@
-__version__ = '2.4.3'
+__version__ = '2.5.0'

snakePipes/shared/defaults.yaml

@@ -6,6 +6,7 @@
 # Note that due to limitations in yaml.dump, only very basic structures are
 # permitted here.
 ################################################################################
+#
 snakemakeOptions: ' --use-conda --conda-prefix /package/anaconda3/envs/ '
 organismsDir: 'shared/organisms'
 clusterConfig: 'shared/cluster.yaml'

snakePipes/shared/rules/Bowtie2.snakefile

@@ -10,7 +10,8 @@ if pairedEnd:
         log: "Bowtie2/logs/{sample}.sort.log"
         params:
             bowtie2_index=bowtie2_index,
-            alignerOpts = str(alignerOpts or ''),
+            alignerOpts = str(alignerOpts or ' ') if not cutntag else " --end-to-end --very-sensitive "\
+            "--no-mixed --no-discordant --phred33 -I 10 -X 700 ",
             mateOrientation = mateOrientation,
             insertSizeMax = insertSizeMax,
             tempDir = tempDir

snakePipes/shared/rules/ChIP_peak_calling.snakefile

@@ -16,22 +16,24 @@ if pairedEnd:
     rule MACS2:
         input:
             chip = "filtered_bam/{chip_sample}.filtered.bam",
-            control =
-                lambda wildcards: "filtered_bam/"+get_control(wildcards.chip_sample)+".filtered.bam" if get_control(wildcards.chip_sample)
-                else [],
-            insert_size_metrics = "deepTools_qc/bamPEFragmentSize/fragmentSize.metric.tsv"
+            frag_size = "deepTools_qc/bamPEFragmentSize/fragmentSize.metric.tsv"
         output:
             peaks = "MACS2/{chip_sample}.filtered.BAM_peaks.xls",
             peaksPE = "MACS2/{chip_sample}.filtered.BAMPE_peaks.xls"
         params:
-            genome_size = genome_size,
             broad_calling =
-                lambda wildcards: "--broad" if is_broad(wildcards.chip_sample) else "",
+                lambda wildcards: "--broad " if is_broad(wildcards.chip_sample) else "",
             control_param =
-                lambda wildcards: "-c filtered_bam/"+get_control(wildcards.chip_sample)+".filtered.bam" if get_control(wildcards.chip_sample)
+                lambda wildcards: " -c filtered_bam/"+get_control(wildcards.chip_sample)+".filtered.bam" if get_control(wildcards.chip_sample)
                 else "",
-            qval_cutoff=qval,
-            mfold=mfold
+            genome_size = str(genome_size),
+            ext_size =
+                lambda wildcards: " --nomodel --extsize "+get_pe_frag_length(wildcards.chip_sample,
+                                                                            "deepTools_qc/bamPEFragmentSize/fragmentSize.metric.tsv") \
+                                                                            if not cutntag else " ",
+            peakCaller_options = lambda wildcards: str(peakCallerOptions or '') if not cutntag else " -p 1e-5 ",
+            bampe_options = lambda wildcards: str(BAMPEPeaks or '')if not cutntag else " ",
+            bam_options = lambda wildcards: str(BAMPeaks or '') if not cutntag else " "
         log:
             out = "MACS2/logs/MACS2.{chip_sample}.filtered.out",
             err = "MACS2/logs/MACS2.{chip_sample}.filtered.err"
@@ -41,18 +43,20 @@ if pairedEnd:
         shell: """
             macs2 callpeak -t {input.chip} {params.control_param} \
                 -f BAM \
-                -g {params.genome_size} --qvalue {params.qval_cutoff}\
+                {params.bam_options} \
+                -g {params.genome_size} \
+                {params.ext_size} \
                 --keep-dup all \
                 --outdir MACS2 \
                 --name {wildcards.chip_sample}.filtered.BAM \
-                --nomodel \
-                --mfold {params.mfold}\
-                --extsize $(cat {input.insert_size_metrics} | grep filtered_bam/{wildcards.chip_sample}.filtered.bam | awk '{{printf("%i",$6)}}') \
+                {params.peakCaller_options} \
                 {params.broad_calling} > {log.out} 2> {log.err}
 
             # also run MACS2 in paired-end mode BAMPE for comparison with single-end mode
             macs2 callpeak -t {input.chip} \
-                {params.control_param} -f BAMPE --qvalue {params.qval_cutoff}\
+                {params.control_param} -f BAMPE \
+                {params.bampe_options} \
+                {params.peakCaller_options} \
                 -g {params.genome_size} --keep-dup all \
                 --outdir MACS2 --name {wildcards.chip_sample}.filtered.BAMPE \
                 {params.broad_calling} > {log.out}.BAMPE 2> {log.err}.BAMPE
@@ -67,26 +71,27 @@ else:
         output:
             peaks = "MACS2/{chip_sample}.filtered.BAM_peaks.xls",
         params:
-            genome_size = int(genome_size),
+            genome_size = str(genome_size),
             broad_calling =
                 lambda wildcards: "--broad" if is_broad(wildcards.chip_sample)
                 else "",
             control_param =
-                lambda wildcards: "-c filtered_bam/"+get_control(wildcards.chip_sample)+".filtered.bam" if get_control(wildcards.chip_sample)
+                lambda wildcards: " -c filtered_bam/"+get_control(wildcards.chip_sample)+".filtered.bam" if get_control(wildcards.chip_sample)
                 else "",
             frag_size=fragmentLength,
-            mfold=mfold,
-            qval_cutoff=qval
+            peakCaller_options = str(peakCallerOptions or ''),
+            bam_options = str(BAMPeaks or '')
         log:
             out = "MACS2/logs/MACS2.{chip_sample}.filtered.out",
             err = "MACS2/logs/MACS2.{chip_sample}.filtered.err"
         benchmark:
             "MACS2/.benchmark/MACS2.{chip_sample}.filtered.benchmark"
         conda: CONDA_CHIPSEQ_ENV
         shell: """
-            macs2 callpeak -t {input.chip} {params.control_param} -f BAM -g {params.genome_size} --qvalue {params.qval_cutoff} --keep-dup all --outdir MACS2 \
-                --name {wildcards.chip_sample}.filtered.BAM --mfold {params.mfold} --extsize {params.frag_size}\
-                {params.broad_calling} > {log.out} 2> {log.err}
+            macs2 callpeak -t {input.chip} {params.control_param} -f BAM -g {params.genome_size} \
+            {params.peakCaller_options} --keep-dup all --outdir MACS2 \
+            --name {wildcards.chip_sample}.filtered.BAM {params.bam_options} --extsize {params.frag_size} \
+            {params.broad_calling} > {log.out} 2> {log.err}
             """