2

我的第一个命令行的"bcftools query -l {input.invcf} | head -n 1"输出打印了 vcf 文件的第一个个体的名称（即IND1）。我想selectvariants GATK在-sn IND1选项中使用该输出。如何将第一条命令行集成到snakemake 中以便在下一条中使用它的输出？

rule selectvar:
    input:
        invcf="{family}_my.vcf"
    params:
        ind= ???
        ref="ref.fasta"
    output:
        out="{family}.dn.vcf"
    shell:
        """
        bcftools query -l {input.invcf} | head -n 1 > {params.ind}
        gatk --java-options "-Xms2G -Xmx2g -XX:ParallelGCThreads=2" SelectVariants -R {params.ref} -V {input.invcf} -sn {params.ind} -O {output.out}
        """
4

3 回答 3

2

有几种选择，但最简单的一种是将结果存储到临时 bash 变量中：

rule selectvar:
   ...
   shell:
        """
        myparam=$(bcftools query -l {input.invcf} | head -n 1)
        gatk -sn "$myparam" ...
        """

正如@dariober 所指出的，如果修改pipefail行为，可能会出现意想不到的结果，请参阅此答案中的示例。

于 2022-02-04T12:55:20.020 回答
1

我想我找到了解决方案：

rule selectvar:
    input:
        invcf="{family}_my.vcf"
    params:
        ref="ref.fasta"
    output:
        out="{family}.dn.vcf"
    shell:
        """
        gatk --java-options "-Xms2G -Xmx2g -XX:ParallelGCThreads=2" SelectVariants -R {params.ref} -V {input.invcf} -sn `bcftools query -l {input.invcf} | head -n 1` -O {output.out}
        """
于 2022-02-04T12:52:37.637 回答
0

当我必须做这些事情时，我更喜欢使用 run 而不是 shell，然后只在最后使用 shell。这样做的原因是因为它可以让snakemake 对 run 语句进行 lint，并在出现问题时提前退出，而不是执行损坏的 shell 命令。

rule selectvar:
    input:
        invcf="{family}_my.vcf"
    params:
        ref="ref.fasta"
        gatk_opts='--java-options "-Xms2G -Xmx2g -XX:ParallelGCThreads=2" SelectVariants'
    output:
        out="{family}.dn.vcf"
    run:
        opts = "{params.gatk_opts} -R {params.ref} -V {input.invcf} -O {output.out}"
        sn_parameter = shell("bcftools query -l {input.invcf} | head -n 1")
        # we could add a sanity check here if necessary, before shelling out
        shell(f"gatk {options} {sn_parameter}")
        """
于 2022-02-07T04:30:22.560 回答