在 Snakemake 中拆分文件

如何解决在 Snakemake 中拆分文件

我有一个简单的问题，但我自己就是想不通。

我有一个输入列表 (a,b,c)。对于每个输入，我需要提取一些数据（1 到 23）：

bcftools view -H a.vcf.gz -r 1 > a_chr1.txt
...
bcftools view -H a.vcf.gz -r 23 > a_chr23.txt

我可以使用 Snakemake 规则中的 FOR 循环来做到这一点：

IDS=['a','b','c']
chrs=range(1,23)
rule:
    input:
        expand("{id}.vcf.gz",id=IDS)
    output:
        expand("{id}_{chr}.txt",chr=chrs,id=IDS)
    run:
        for i in IDS:
            for c in chrs:
                shell("bcftools view -H {i}.vcf.gz -r {c} > {i}_chr{c}.txt")

，但 FOR 循环不会并行化它。我需要一个合适的 Snakemake 方式，如下所示，但它不起作用。

IDS=['a',23)

rule:
    input:
        expand("{id}.vcf.gz",id=IDS)
    params:
        c=expand("{chr}",chr=chrs)
    shell:
        "bcftools view -H {input} -r {params.c} > {output}"

你能帮忙吗？

解决方法

您没有利用这里的蛇形通配符。如果您在输入和输出中指定 expand，那么 snakemake 将只运行一次规则。它告诉snakemake运行规则需要所有的vcf文件，并且这个规则将产生所有的分割文件。您需要的是一个规则，它可以应用于任何 vcf 文件，并且只会产生一个分割（按 chr）的输出。

IDS=['a','b','c']
chrs=range(1,23)

rule all:
    input: expand("{id}_{chr}.txt",chr=chrs,id=IDS)

rule splitByChr:
    input:
        "{id}.vcf.gz"
    output:
        "{id}_{chr}.txt"
    shell:
        "bcftools view -H {input} -r {wildcards.chr} > {output}"

此处的 rule all 将根据需要多次触发 rule splitByChr。
另请注意，扩展函数中的 {id} 和 {chr} 不是通配符。它们是定义的扩展参数的占位符。