当两个不同的规则路径可以生成一个给定的输出时，snakemake能避免歧义吗？

Question

初始工作流程

我有一个snakefile，可以从成对的结束数据生成一些输出。在

在这个snakefile中，我有一个规则“安装”存储在配置文件（get_raw_data）中的给定信息的数据。在

然后我有一个规则，它使用这些数据生成工作流其余部分依赖的中间文件（run_tophat）。在

以下是这些规则的输入和输出（OPJ代表os.path.join）：

rule get_raw_data:
    output:
        OPJ(raw_data_dir, "{lib}_1.fastq.gz"),
        OPJ(raw_data_dir, "{lib}_2.fastq.gz"),

（稍后将详细介绍该规则的实施情况）

^{pr2}$

我的主要原则是：

rule all:
    input:
        expand(OPJ(output_dir, "{lib}", "junctions.bed"), lib=LIBS),

将工作流扩展到单端数据

我现在必须在单端数据上运行我的工作流。在

我希望避免最终输出具有不同的名称模式，这取决于数据是单端还是成对端。在

我可以很容易地将上述两个规则的变体用于单端数据（get_raw_data_single_end和{}），它们的输入和输出如下：

rule get_raw_data_single_end:
    output:
        OPJ(raw_data_dir, "{lib}.fastq.gz")

rule run_tophat_single_end:
    input:
        transcriptome = OPJ(annot_dir, "dmel-all-r5.9.gff"),
        fq = OPJ(raw_data_dir, "{lib}.fastq.gz"),
    output:
        junctions = OPJ(output_dir, "{lib}", "junctions.bed"),
        bam = OPJ(output_dir, "{lib}", "accepted_hits.bam"),

如何为snakemake提供足够的信息来选择正确的规则路径？在

配置文件包含有关lib通配符是否与单端数据或成对结尾数据关联的信息：库名称是lib2raw或lib2raw_single_end字典中的键（这两个字典都是从配置文件中读取的）。在

我不希望同一个库名在两个字典中都是一个键。因此，在某种意义上，我希望执行工作流的单端分支还是成对端分支并不含糊。在

get_raw_data和{}都使用函数lib2data（使用这些字典）来确定运行哪个shell命令来“安装”数据。在

以下是此函数的简化版本（实际的函数包含一个额外的分支，用于从SRR标识符生成数据命令）：

def lib2data(wildcards):
    lib = wildcards.lib
    if lib in lib2raw:
        raw = lib2raw[lib]
        link_1 = "ln -s %s %s_1.fastq.gz" % (raw.format(mate="1"), lib)
        link_2 = "ln -s %s %s_2.fastq.gz" % (raw.format(mate="2"), lib)
        return "%s\n%s\n" % (link_1, link_2)
    elif lib in lib2raw_single_end:
        raw = lib2raw_single_end[lib]
        return "ln -s %s %s.fastq.gz\n" % (raw, lib)
    else:
        raise ValueError("Procedure to get raw data for %s unknown." % lib)

除了它们的输出之外，这两个get_raw_data*规则是相同的，其工作方式如下：

params:
    shell_command = lib2data,
shell:
    """
    (
    cd {raw_data_dir}
    {params.shell_command}
    )
    """

snakemake是否能够确定正确的规则路径？给定的信息不是在规则输入和输出中编码的，而是只在配置文件和函数中编码的？在

看来情况并非如此。实际上，我正在尝试测试我的新snakefile（添加了*_single_end规则），，但是在get_raw_data规则的执行过程中，KeyError发生了，而为其执行规则的库与单端数据关联。在

如何实现所需的行为（一个双分支工作流能够使用配置中的信息来选择正确的分支）？在

编辑：KeyError是由于lib2data

中的错误造成的

在使用正确的字典获取与库名称关联的数据后，我出现了以下错误：

AmbiguousRuleException:
Rules run_tophat and run_tophat_single_end are ambiguous for the file tophat_junction_discovery_revision_supplement/HWT3/junctions.bed.
Expected input files:
        run_tophat: ./HWT3_1.fastq.gz ./HWT3_2.fastq.gz Annotations/dmel-all-r5.9.gff
        run_tophat_single_end: ./HWT3.fastq.gz Annotations/dmel-all-r5.9.gff

编辑2：向get_raw_data*规则添加输入

在阅读了this post on the snakemake mailing list之后，我试图在规则中添加一些输入以避免歧义。在

def lib2data_input(wildcards):
    lib = wildcards.lib
    if lib in lib2raw:
        raw = lib2raw[lib]
        return [raw.format(mate="1"), raw.format(mate="2")]
    elif lib in lib2raw_single_end:
        raw = lib2raw_single_end[lib]
        return [raw]
    else:
        raise ValueError("Procedure to get raw data for %s unknown." % lib)

rule get_raw_data:
    input:
        lib2data_input
# [same output, params and shell as before]
# [same modification for the single-end case]

这将产生一个MissingInputException。奇怪的是，据说丢失的文件确实存在。这个把戏有用吗？（无法复制此项，现在将导致：）

^{9}$

我为“数据安装”规则指定输入的方法显然不足以引导snakemake找到正确的规则。在

Answer 1

作为user1829905的suggested，我曾试图使get_raw_data*规则成为一个规则，但由于该规则的输出是可变的，所以失败了。在

但是，我可以将run_tophat*规则合并为一个：它们具有相同的输出。在

rule run_tophat:
    input:
        transcriptome = OPJ(annot_dir, "dmel-all-r5.9.gff"),
        fq = lib2fq,
    output:
        junctions = OPJ(output_dir, "{lib}", "junctions.bed"),
        bam = OPJ(output_dir, "{lib}", "accepted_hits.bam"),

我尝试使用以下函数生成此融合规则的输入：

^{pr2}$

但此尝试失败，返回InputFunctionException：

^{3}$

然而，让第二条规则的输入明确地定义为第一条规则的输出，就解决了这个问题。在

def lib2fq(wildcards):
    lib = wildcards.lib
    if lib in lib2raw_single_end:
        return rules.get_raw_data_single_end.output
    else:
        return rules.get_raw_data.output

我不完全明白为什么会有这种不同。在

Answer 2

我不知道这是否有用，但您可以使用函数来定义规则的输入。这样，您就可以使用相同的规则来处理单端或成对端数据，前提是规则的输出是相同的。。。在

def my_inputs(wildcards):
    data_type = config["data_type"]
    if (data_type == "pe"):
        input = ...
    elif (data_type == "se"):
        input = ...
    return input

rule my_rule:
    input: my_inputs
    ...