2

我正在使用 Genrich 开发一个 ATACseq 管道以与 Snakemake 一起运行。

事实上，Genrich 允许在同一步骤中调用多个重复中的峰，从而避免了额外的步骤（即 IDR）。

在 Snakemake 中，我找到了同时返回所有我想要的样本（即从一个条件复制）的方法，但是如果每个文件都被引用，Genrich 要求逗号分隔的文件作为输入或空格分隔的文件。

通常，输入返回一个以空格分隔的文件列表（即file1 file2 file3），由于我不知道如何使它返回逗号分隔的文件，所以我尝试引用它们。

理论上，在 Snakemake 5.8.0 版本之后，您可以参考{input:q}规则的 shell 命令中的输入来返回引用的输入，如此处所述。

但是，在我的情况下，返回的输入没有被引用。

我创建了一个测试规则来查看输入是如何返回的：

rule genrich_merge_test:
    input:
        lambda w: expand("{condition}.sorted.bam", condition = SAMPLES.loc[SAMPLES["CONDITION"] == w.condition].NAME),
    output:
        "{condition}_peaks.narrowPeak",
    shell:
        """
        echo {input:q} > {output}     
        """

存储在输出文件中的返回输入是：

rep1.sorted.bam rep2.sorted.bam

有人知道如何解决这个问题并返回引用的输入或返回逗号分隔的文件列表而不是空格分隔的文件吗？

谢谢你。

4

2 回答 2

0

假设您的输入文件名不包含空格（如果包含空格，我强烈建议避免使用它们），您可以简单地将文件列表放在引号中，您不需要引用列表中的每个文件：

rule genrich:
    input:
        t= ['a.bam', 'b.bam'],
    ...
    shell:
        r"""
        Genrich -t '{input.t}' ...
        """

（注意周围的单引号'{input.t}'）

于 2020-10-27T09:00:00.297 回答
0

我在想 echo 和 shell 可能会在管道输出之前剥离引号，但是检查snakemake -p正在执行的命令表明它们不存在。当存在空格时，引号似乎只显示单个文件名。

Dariober 的答案应该可以引用该列表，但为了完整起见，如果您想要一个以逗号分隔的文件列表，请在 params 指令中使用 lambda 函数：

rule genrich_merge_test:
    input:
        lambda w: expand("{condition}.sorted.bam", 
                         condition=SAMPLES.loc[SAMPLES["CONDITION"] == w.condition].NAME),
    params:
        files=lambda wildcards, input: ','.join(input)
    output:
        "{condition}_peaks.narrowPeak",
    shell:
        """
        echo {params.files} > {output}     
        """

编辑

这是一个玩具示例，演示了在输入中使用参数：

# snakefile
inputs = expand('{wc}.out', wc=range(4))

rule all:
    input: "test_peaks.narrowPeak"

rule genrich:
    input:
        inputs
    params:
        files=lambda wildcards, input: ','.join(input)
    output:
        "test_peaks.narrowPeak",
    shell:
        """
        echo {params.files} > {output}     
        """

rule generator:
    output: touch('{file}.out')

$ snakemake -np
...
rule genrich:
    input: 0.out, 1.out, 2.out, 3.out
    output: test_peaks.narrowPeak
    jobid: 1


        echo 0.out,1.out,2.out,3.out > test_peaks.narrowPeak 
...

也如此处所示

请注意，与 input 指令相比，params 指令可以选择采用比通配符更多的参数，即输入、输出、线程和资源。

于 2020-10-27T12:48:52.513 回答