Snakemake：将多个输入用于具有多个子组的一个输出的规则

我有一个工作管道，用于下载、比对和对公共测序数据执行变体调用。问题是它目前只能在每个样本的基础上工作（i.e作为每个单独测序实验的样本）。如果我想对一组实验（例如样本的生物和/或技术复制）执行变体调用，则它不起作用。我试图解决它，但我无法让它发挥作用。

这是对齐规则的简化：

rule alignment:
    input:
        rules.download.output.fastq
    output:
        '{group}/alignment/{sample}.bam'
    shell:
        "bash scripts/02_alignment.sh {wildcards.group} {wildcards.samples}"

对于变体调用也是如此：

rule variant_calling:
    input:
        rules.alignment.output
    output:
        '{group}/variants/{sample}.vcf.gz'
    shell:
        "bash scripts/03_variant_calling.sh {wildcards.sample} {wildcards.group}"

这工作得很好，因为有一个.vcf为每个对齐生成的文件.bam文件。但我想做的是生成一个.vcf任意数量的文件.bam文件。我有一个pandas包含所有的数据框sample名称及其对应的group。我基本上想改变output第二条规则为'{group}/variants/{group}.vcf'，但我所做的一切都在某种程度上失败了。

我的想法是提供与所有每组对齐的规则.bam文件作为输入，然后只需为它运行的脚本提供它们所在的目录。问题是我找不到一种方法以这种按组的方式提供输入：要么我给出每个样本（作为工作管道），要么我给出所有.bam每个组变体调用的文件，无论它们实际属于哪个组。我不能只使用通配符，因为{sample}最后的输出中不存在通配符。我也尝试使用函数作为输入，但这会导致与上面相同的问题。

问题的关键似乎是分组的层次：如果我想对所有对齐的执行变体调用.bam数据集中的文件作为一个整体，可能会很好地工作，给出我上面提到的问题。问题在于整个数据集的子组：

  sample1      sample2             sample1      sample2      sample3
     |            |                   |            |            |
     |            |                   |            |            |
     --------------                   ---------------------------
            |                                      |
            |                                      |
          group1                                 group2

关于如何解决这个问题有什么想法吗？

您必须使用某种结构来将样本分组：

GROUPS = {
    "group1":["sample1","sample2"],
    "group2":["sample1","sample2","sample3"]
}

然后是这样的：

rule all:
    input:
         expand("{group}/variants/{group}.vcf.gz", group=list(GROUPS.keys()))

rule alignment:
    input:
        rules.download.output.fastq
    output:
        '{group}/alignment/{sample}.bam'
    shell:
        "bash scripts/02_alignment.sh {wildcards.group} {wildcards.samples}"

rule variant_calling:
    input:
        lambda wildcards: expand("{group}/alignment/{sample}.bam", group=wildcards.group, sample=GROUPS[wildcards.group]
    output:
        '{group}/variants/{group}.vcf.gz'
    shell:
        "bash scripts/03_variant_calling.sh {input} {output}"

当然，有些规则是您没有展示的，但我想您会明白的！

规则variant_calling中的shell命令可能很难处理，但您始终可以将目录定义为参数，例如：

params: groupAlignDir = "{group}/alignment"

并在 shell 中使用它：

"bash scripts/03_variant_calling.sh {params.groupAlignDir} {output}"

并获取脚本“variant_calling.sh”目录中的所有 bam 文件