snakemake 学习(2)

2019-05-30  本文已影响0人  BioLearner

举一个小例子:

通过在Snakefile中指定规则来定义Snakemake工作流。 规则通过指定如何从输入文件集创建输出文件集将工作流分解为小步骤(例如,单个工具的应用)。Snakemake 通过匹配文件名自动确定规则之间的依赖关系

1、比对reads

# 首先创建Snakefile文件,然后在文件中编写规则,写完后要记得给Snakefile加可操作权限-x
如果不是写在Snakefile文件,则snakemake --snakefile 文件名

rule bwa_map:
    input:
        "data/genome.fa",   #常见错误是忘记输入或输出项之间的逗号。由于Python连接后续字符串,这可能会导致意外行为。
        "data/samples/A.fastq"
    output:
        "mapped_reads/A.bam"
    shell:
        "bwa mem {input} | samtools view -Sb - > {output}"
#以下文字文字为上述规则解读: 
一个Snakemake规则有一个名称(此处为bwa_map)和一些指示(此处为input,output和shell)
在input和output指令之后是那些预计将在规则中使用或创建的文件列表在最简单的情况下,这些只是Python字符串
shell指令后跟一个包含要执行的shell命令,此处通过{input}和{output}来指定输入文件和输出文件
由于规则有多个输入文件,所以Snakemake将用空格将它们连接起来
换句话说,Snakemake在执行命令之前将用data/genome.fa data/samples/A.fastq替换{input}
bwa mem的输出通过管道成为samtools view的输入,最终输出到规则中的{output}中

2、bwa_map规则批处理

显然,该规则仅适用于具有文件读取的单个样本data/samples/A.fastq。但是,Snakemake允许使用命名通配符来使规则可批处理。只需用通配符{sample}替换第二个输入文件A和输出文件。

SAMPLES = ['A','B','C']

rule all:
    input:
        expand('mapped_reads/{sample}.bam', sample=SAMPLES)

rule bwa_map:
    input:
        "data/genome.fa",
        "data/samples/{sample}.fastq"
    output:
        "mapped_reads/{sample}.bam"
    shell:
        "bwa mem {input} | samtools view -Sb - > {output}"
#以下为注释
此前碰到这个问题:WorkflowError:
Target rules may not contain wildcards. Please specify concrete files or a rule without wildcards.
直译过来是:目标规则可能不包含通配符,请指定没有通配符的具体文件或规则
看似是通配符的使用问题,实则和 rule all 的 input 错误有关
Snakemake的rule all需要输入一个列表,这个列表包括所有的步骤会产生的文件(也就是output)
实际上,Snakemake的工作机制是先通过shell产生output文件,再检查这些output文件是否在指定的目录下存在,确认存在后再执行下一个rule 
这个检查过程就是通过寻找rule all下的input列表完成的
因此,在input列表没写清楚的情况下,如果在rule中还使用的通配符,就会发生相应的报错
改正方法就是修改你的input列表,确认每一个rule输出的文件和他们的目录都存在于input列表下

此外,请注意,如果rule有多个输出文件,Snakemake要求它们都具有完全相同的通配符。否则,可能会发生同一rule中的两个job写同一个文件。

snakemake对于没有关联的rule,其顺序按rule名字母顺序来

上一篇下一篇

猜你喜欢

热点阅读