Snakemake可变数量的文件

Question

我处于一种情况，我想将我的工作流程分散到可变数量的块中，我事先并不知道。也许通过具体来解释问题是最容易的：

有人递给我使用 bcl2fastq 和 no-lane-splitting 选项解复用的 FASTQ 文件。我想根据车道拆分这些文件，分别映射每个车道，然后最后再次收集所有内容。但是，我事先不知道车道数。

理想情况下，我想要这样的解决方案，

rule split_fastq_file: (...)  # results in N FASTQ files
rule map_fastq_file: (...)  # do this N times
rule merge_bam_files: (...)  # merge the N BAM files

但我不确定这是否可能。 expand 函数要求我知道车道数，也不知道如何使用通配符。

我应该说我是 Snakemake 的新手，我可能完全误解了 Snakemake 的工作原理。我花了一些时间来习惯通过专注于输出文件然后向后工作来思考事情 "upside-down"。

Answer 1

一个选项是在拆分 fastq 时使用 checkpoint，这样您就可以在以后动态地重新评估 DAG 以获得结果通道。

这是 MWE 的一步一步：

设置并制作一个示例 fastq 文件。

# Requires Python 3.6+ for f-strings, Snakemake 5.4+ for checkpoints
import pathlib
import random

random.seed(1)

rule make_fastq:
    output:
        fastq = touch("input/{sample}.fastq")

创建一个介于 1 和 9 之间的随机车道数，每个车道具有从 1 到 9 的随机标识符。请注意，我们将此声明为 checkpoint，而不是规则，以便我们稍后可以访问结果.此外，我们将此处的输出声明为特定于示例的目录，以便我们稍后可以在其中进行 glob 以获取创建的通道。

checkpoint split_fastq:
    input:
        fastq = rules.make_fastq.output.fastq
    output:
        lane_dir = directory("temp/split_fastq/{sample}")
    run:
        pathlib.Path(output.lane_dir).mkdir(exist_ok=True)
        n_lanes = random.randrange(1, 10)-
        lane_numbers = random.sample(range(1, 10), k = n_lanes)
        for lane_number in lane_numbers:
            path = pathlib.Path(output.lane_dir) / f"L00{lane_number}.fastq"
            path.touch()

做一些中间处理。

rule map_fastq:
    input:
        fastq = "temp/split_fastq/{sample}/L00{lane_number}.fastq"
    output:
        bam = "temp/map_fastq/{sample}/L00{lane_number}.bam"
    run:
        bam = pathlib.Path(output.bam)
        bam.parent.mkdir(exist_ok=True)
        bam.touch()

为了合并所有已处理的文件，我们使用输入函数来访问在 split_fastq 中创建的通道，以便我们可以对这些进行动态 expand。我们对中间处理步骤链中的最后一条规则执行 expand，在本例中为 map_fastq，以便我们请求正确的输入。

def get_bams(wildcards):
    lane_dir = checkpoints.split_fastq.get(**wildcards).output[0]
    lane_numbers = glob_wildcards(f"{lane_dir}/L00{{lane_number}}.fastq").lane_number
    bams = expand(rules.map_fastq.output.bam, **wildcards, lane_number=lane_numbers)
    return bams

这个输入函数现在让我们可以轻松访问我们希望合并的 bam 文件，无论它们有多少，也不管它们叫什么。

rule merge_bam:
    input:
        get_bams
    output:
        bam = "temp/merge_bam/{sample}.bam"
    shell:
        "cat {input} > {output.bam}"

此示例运行，random.seed(1) 恰好创建了三个通道（l001、l002 和 l005）。

如果您不想使用 checkpoint，我认为您可以通过为 merge_bam 创建一个输入函数来实现类似的功能，该函数打开原始输入 fastq，扫描读取名称车道信息，并预测输入文件应该是什么。然而，这似乎不太稳健。

Snakemake可变数量的文件

Snakemake variable number of files

snakemake