This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Snakemake 'Unexpected keyword ref in rule definition'

Hello everyone, I am trying to get more familiar with Snakemake. I have created this Snakefile:

rule alignment_haplotypecaller:
  output: 
      vcf = 'variants/{sample}.g.vcf.gz'
  input:  
      bam = 'alignment/{sample}.bam'
  reference:
      fasta = 'reference.fasta'
  shell:
    r"""
        samtools index {input.bam}
    ~/gatk-4.2.0.0/gatk HaplotypeCaller \
    --reference {reference.fasta} \
    --input  {input.bam} \
    --output {output.vcf} \
    --emit-ref-confidence BP_RESOLUTION 
    """

An I run Snakemake with the following command:

snakemake -j1 -F -p alignment_haplotypecaller

I get this error message:

SyntaxError in line 7 of /mnt/shared/usr/test_snakemake/Snakefile: Unexpected keyword ref in rule definition (Snakefile, line 7)

Would anyone know why the key line 7 ("fasta") is problematic?

snakemake python
• 6,761 views • link
updated by cpad0112 cpad0112 • written by Begonia_pavonina • Begonia_pavonina

Check if reference: can be used in rules.

• 0 views • link
written by cpad0112 cpad0112

2 answers

Before executing the snakemake is useful to run in dry mode (using -n ) to see if there's any problem in your syntax or input/output. You're using an expected keyword for the rule, reference . Replace it with params.

• 0 views • link
written by Shred Shred

Thanks you, the keyword was indeed the issue. I though they could take any name, but it is not the case. I didn't know about the dry mode. Thanks.

Now there is a wildcard error, I will try to manage them differently.

WorkflowError: Target rules may not contain wildcards. Please specify concrete files or a rule without wildcards at the command line, or have a rule without wildcards at the very top of your workflow (e.g. the typical "rule all" which just collects all results you want to generate in the end).

• 0 views • link
written by Begonia_pavonina • Begonia_pavonina

Snakemake wildcards are different from Bash ones. You need to inform Snakemake how to expand the wildcard in order to be consistent in input and output. In your case, for example, you could declare :

import glob
SAMPLES,= glob_wildcards('variants/{sample}.g.vcf.gz') # retrieve files via glob, instruct wildcards

Then, as suggested, declare a rule all, in which you declare how Snakemake rule's input will be expanded. Note that this is convenient because allows to use the same wildcards across multiple rule's input, just specifying how

rule all:
    input:
        expand('variants/{sample}.g.vcf.gz', sample=SAMPLES) # --> referred to alignment_haplotypecaller
        expand('alignment/{sample}.bam', sample=SAMPLES) # --> other rule example, which uses haplotypecaller input with same wildcard

Read more here (alternative to Snakemake's documentation )

• 0 views • link
written by Shred Shred

Thanks Shred! I have followed other documentation, but the link you provided will certainly help.

I have tried the following script, and got this error:

WorkflowError: Target rules may not contain wildcards. Please specify concrete files or a rule without wildcards at the command line, or have a rule without wildcards at the very top of your workflow (e.g. the typical "rule all" which just collects all results you want to generate in the end).

Would anyone know what is the cause of the error here?

# Input conditions and replicates to process

# Dictionnaries:
SAMPLES = glob_wildcards('alignment/Sample_{sample}.baits_conch.bam ') 

# Expand wildcards:    
    rule all:
      input:
          expand('alignment/Sample_{sample}.baits_conch.bam', sample=SAMPLES) 

    # First rule:
    rule alignment:
      output:
          snp = 'variants/Sample_{sample}.baits_conch.g.vcf.gz'
      input:
          ali = 'alignment/Sample_{sample}.baits_conch.bam'
      params:
          ref = 'path/to/reference.fasta'
      shell:
        r"""samtools index {input.ali}
        ~/gatk-4.2.0.0/gatk HaplotypeCaller \
        --reference {params.ref} \
        --input  {input.ali} \
        --output {output.snp} \
        --emit-ref-confidence BP_RESOLUTION
        """
• 0 views • link
written by Begonia_pavonina • Begonia_pavonina

why is your input in rule all is bam? Shouldn't be vcf.gz. (from output of alignment rule)? or is it a typo?

• 0 views • link
written by cpad0112 cpad0112

Thanks for the question.

It is not a typo. BAM files are the input in HaplotypeCaller, the shell command, which will output VCF files.

Do you mean that I should expand the output (VCF) rather than expanding the input (BAM)?

• 0 views • link
written by Begonia_pavonina • Begonia_pavonina

Yes. That is how snakemake builds dag. Here is an example snakemake workflow for variant calling using gatk. Look at the final snakemake file: https://github.com/snakemake-workflows/dna-seq-gatk-variant-calling/blob/main/workflow/Snakefile and rule for variant calling: https://github.com/snakemake-workflows/dna-seq-gatk-variant-calling/blob/main/workflow/rules/calling.smk

• 0 views • link
written by cpad0112 cpad0112

Thank you cpad0112, I have made a fork of this repository. I feel very confused about the output-focused approach of Snakemake, but working a bit on the gatk-variant-calling pipeline will be indeed helpful.

• 0 views • link
written by Begonia_pavonina • Begonia_pavonina

Each frame work comes with it's own quirks. Good luck with Snakemake. If you want to try another frame work, I suggest nextflow.

• 0 views • link
written by cpad0112 cpad0112

Thank Shred, the indent issue is a part of my not-well configured VIM. I will be careful to that. Do we need the tuple command you mentioned? There is one output for one input in this case.

• 0 views • link
written by Begonia_pavonina • Begonia_pavonina

Log in to answer this question.

More posts like this