1
$\begingroup$

I'm trying to automate some of my previous bash scripts with snakemake. Below is one of the recipe I created:

import os
shell.executable("/bin/bash")

workdir: "/proj/uppstore2018136/publication/jupyterNotebooks"

ENSEMBL_RELEASE='94'
ENSEMBL_GRCm38_BASE='ftp://ftp.ensembl.org/pub/release-'+ENSEMBL_RELEASE+'/fasta/mus_musculus/dna'
ENSEMBL_GRCm38_GTF_BASE='ftp://ftp.ensembl.org/pub/release-'+ENSEMBL_RELEASE+'/gtf/mus_musculus'
GTF_FILE='Mus_musculus.GRCm38.'+ENSEMBL_RELEASE+'.gtf'

rule all:
    input:
        "genome/genome.fa",
        "genome/"+GTF_FILE,
        "indexes/bowtie2/reference.1.bt2",
        "indexes/bowtie2/reference.2.bt2",
        "indexes/bowtie2/reference.3.bt2",
        "indexes/bowtie2/reference.4.bt2",
        "indexes/bowtie2/reference.rev1.bt2",
        "indexes/bowtie2/reference.rev2.bt2",
        "indexes/hisat2/reference.1.ht2",
        "indexes/hisat2/reference.2.ht2",
        "indexes/hisat2/reference.3.ht2",
        "indexes/hisat2/reference.4.ht2",
        "indexes/hisat2/reference.5.ht2",
        "indexes/hisat2/reference.6.ht2",
        "indexes/hisat2/reference.7.ht2",
        "indexes/hisat2/reference.8.ht2"

rule getGenome:
    output:
        output1="genome/genome.fa",
        output2="genome/"+GTF_FILE
    shell:
           '''
           curl -o genome/Mus_musculus.GRCm38.dna.primary_assembly.fa.gz {ENSEMBL_GRCm38_BASE}/Mus_musculus.GRCm38.dna.primary_assembly.fa.gz
           gunzip genome/Mus_musculus.GRCm38.dna.primary_assembly.fa.gz
           mv genome/Mus_musculus.GRCm38.dna.primary_assembly.fa genome/genome.fa
           curl -o genome/{GTF_FILE}.gz {ENSEMBL_GRCm38_GTF_BASE}/{GTF_FILE}.gz
           gunzip genome/{GTF_FILE}.gz
           '''

rule bowtie2Build:
    input: "genome/genome.fa"
    params:
        basename="indexes/bowtie2/reference"
    output:
        output1="indexes/bowtie2/reference.1.bt2",
        output2="indexes/bowtie2/reference.2.bt2",
        output3="indexes/bowtie2/reference.3.bt2",
        output4="indexes/bowtie2/reference.4.bt2",
        outputrev1="indexes/bowtie2/reference.rev1.bt2",
        outputrev2="indexes/bowtie2/reference.rev2.bt2"
    threads: 20
    shell: 
        '''
        module load bowtie2
        bowtie2-build --threads {threads} {input} {params.basename}
        '''

rule hisat2Build:
    input: "genome/genome.fa",
        "genome/"+GTF_FILE 
    params:
        basename="indexes/hisat2/reference"
    output:
        output1="indexes/hisat2/reference.1.ht2",
        output2="indexes/hisat2/reference.2.ht2",
        output3="indexes/hisat2/reference.3.ht2",
        output4="indexes/hisat2/reference.4.ht2",
        output5="indexes/hisat2/reference.5.ht2",
        output6="indexes/hisat2/reference.6.ht2",
        output7="indexes/hisat2/reference.7.ht2",
        output8="indexes/hisat2/reference.8.ht2",
    threads: 20
    shell: 
        '''
        module load HISAT2
        hisat2_extract_splice_sites.py genome/{GTF_FILE} > genome/genome.ss
        hisat2_extract_exons.py genome/{GTF_FILE} > genome/genome.exon
        hisat2-build -f -p {threads} {input} --ss genome/genome.ss --exon genome/genome.exon {params.basename}
        rm genome/genome.ss
        rm genome/genome.exon
        '''

rule clean:
    shell:
        '''
        rm genome/*
        rm indexes/bowtie2/*
        rm indexes/hisat2/*
        '''

And this is my driver script which I'm submitting as a job itself:

#!/bin/bash

module load snakemake

snakemake --jobs 1 --local-cores 20 --cluster-config /crex/proj/uppstore2018136/publication/jupyterNotebooks/recipes/cluster.json --latency-wait 120 --max-jobs-per-second 1 -c "sbatch -A {cluster.account} -p {cluster.partition} -n {-cluster.n} -t {cluster.time} -C {cluster.C}" -s indexMake

But I'm getting this error:

Building DAG of jobs...
Using shell: /bin/bash
Provided cluster nodes: 1
Job counts:
    count   jobs
    1   all
    1   bowtie2Build
    1   getGenome
    1   hisat2Build
    4

rule getGenome:
    output: genome/genome.fa, genome/Mus_musculus.GRCm38.94.gtf
    jobid: 2

RuleException in line 30 of /crex/proj/uppstore2018136/publication/jupyterNotebooks/recipes/indexMake:
NameError: The name '-cluster' is unknown in this context. Please make sure that you defined that variable. Also note that braces not used for variable access have to be escaped by repeating them, i.e. {{print $1}

}

I really am not sure what I'm doing wrong here ? Is there a conflict with variables or there is something wrong with the snakemake version I'm using (4.5.0).

$\endgroup$
  • $\begingroup$ (vote to leave open) I know that Snakemake can be used for non-bioinformatics problems. But in practice it is used solely by bioinfo community therefore I think debugging snakefiles is appropriate for this Q&A site. $\endgroup$ – Kamil S Jaron Dec 10 '18 at 14:32
  • $\begingroup$ (vote to close) The official Snakemake docs point folks to Stack Overflow. I suggest we use StackExchange's migrate mechanism to send generic Snakemake questions to SO to keep it all in one place. $\endgroup$ – Daniel Standage Dec 26 '18 at 18:10
  • $\begingroup$ I'm voting to close this question as off-topic because it belongs on StackOverflow (although I lack the reputation to initiate the migration). $\endgroup$ – Daniel Standage Dec 26 '18 at 18:12
  • $\begingroup$ @DanielStandage it does seem to be aimed primarily at bioinformaticians, and is indeed a tool commonly used by the community, so it seems reasonable that we should support it. $\endgroup$ – terdon Jan 6 '19 at 1:13
5
$\begingroup$

In your submission script you have -n {-cluster.n}, hence the error from snakemake I guess. Are you sure it is not a typo and it should be -n {cluster.n} instead?

| improve this answer | |
$\endgroup$
  • $\begingroup$ Yes, that was indeed the error. I'm not sure how it crept in - working on Sundays is not that productive - $\endgroup$ – Siddharth Dec 9 '18 at 19:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.