如何避免Snakemake在expand返回空文件列表时执行无通配符聚合规则？

我正在开发一个带有可选规则的Snakemake工作流，这些规则仅当制表符分隔的样本文件中存在可选元数据时才会执行。简单来说，原始样本数据可以先运行前十几条规则，后续步骤需要耗时收集的重复样本信息，而前十几条规则的输出能辅助这项收集工作。

我已经设置了分开的all规则来运行工作流的不同部分，但希望all规则能尽可能执行所有可运行的步骤：没有元数据时只运行前十几条规则，有元数据时则运行全部规则。

目前的实现大体可行，但存在一个问题：当元数据缺失时，一条无通配符的聚合规则会尝试在没有输入的情况下运行。该规则如下：

rule combine_idr_peaks:
    input:
        expand(
            "results/idr/{dataset}/merged.idr.narrowPeak",
            dataset=DATASETS,
        ),
    output:
        temp("results/idr/merged.idr.narrowPeak"),
    conda:
        "../envs/bedtools_coreutils.yml"
    shell:
        """
        cat {input:q} | \
        sortBed | \
        mergeBed > {output:q}
        """

当元数据缺失时，使用--printshellcmds参数查看命令，会得到如下内容：

cat  | sortBed | mergeBed > results/idr/merged.idr.narrowPeak

显然，工作流会挂起，因为cat命令在等待输入。
DATASETS变量由蛇形文件顶部的元数据填充，测试时故意缺失该可选元数据。

请问如何让Snakemake在输入为空时不执行该规则？

解决方案 #

方法1：动态控制all目标的输入列表 #

既然DATASETS是决定输入是否为空的核心变量，你可以在all规则中根据DATASETS是否为空，动态决定是否包含聚合规则的输出：

rule all:
    input:
        # 前十几条规则的固定输出
        "results/initial_step1.txt",
        "results/initial_step2.txt",
        # 仅当DATASETS非空时，才添加聚合规则的输出
        *([] if not DATASETS else ["results/idr/merged.idr.narrowPeak"])

当DATASETS为空时，all规则不会要求生成聚合文件，Snakemake也就不会触发combine_idr_peaks规则，从根源上避免了空输入问题。

方法2：在规则内添加输入有效性判断 #

如果需要保留规则定义，可以在shell脚本中添加条件判断，输入为空时直接创建空输出文件，避免cat命令挂起：

rule combine_idr_peaks:
    input:
        peaks=expand(
            "results/idr/{dataset}/merged.idr.narrowPeak",
            dataset=DATASETS,
        ),
    output:
        temp("results/idr/merged.idr.narrowPeak"),
    conda:
        "../envs/bedtools_coreutils.yml",
    params:
        has_input=lambda wildcards, input: len(input.peaks) > 0,
    shell:
        """
        if [[ {params.has_input} -eq 1 ]]; then
            cat {input.peaks:q} | sortBed | mergeBed > {output:q}
        else
            # 输入为空时创建空文件，满足Snakemake的输出要求
            touch {output:q}
        fi
        """

方法3：使用条件规则（Snakemake 7.0+） #

如果你的Snakemake版本在7.0及以上，可直接用条件语句包裹规则，仅当DATASETS非空时才定义该规则：

if DATASETS:
    rule combine_idr_peaks:
        input:
            expand(
                "results/idr/{dataset}/merged.idr.narrowPeak",
                dataset=DATASETS,
            ),
        output:
            temp("results/idr/merged.idr.narrowPeak"),
        conda:
            "../envs/bedtools_coreutils.yml",
        shell:
            """
            cat {input:q} | sortBed | mergeBed > {output:q}
            """

这种方式最直观，DATASETS为空时规则不会被加载，自然不会被执行。

内容的提问来源于stack exchange，提问作者hepcat72