Snakemake 如何设计空输出文件的下游规则?

Snakemake How to design downstream rules on empty output files?

我想知道如何处理下游规则中的空输出文件。 使用 SHOVILL 组装短读 fastq 数据可能会失败并产生 0 字节 contigs.fa

如果 PROKKA 的基因组注释是 运行 在 0 字节文件上,它 returns 一个错误:

 'input.fasta' is not a readable non-empty FASTA file

Snakemake 的组装和注释规则:

rule shovill:
    input:
        fw="input/{fastq}_R1.fastq.gz",
        rv="input/{fastq}_R2.fastq.gz"
    output:
        directory("results/shovill/{sample}")
    threads: 16
    shell:
        "shovill --cpus {threads} --R1 {input.fw} --R2 {input.rv} --outdir {output}\
            2> {log}"

我目前的解决方案是使用|| true 并初步生成结果目录。

rule prokka:
    input:
        "results/shovill/{sample}/contigs.fa"
    output:
        directory("results/prokka/{sample}")
    threads: 8
    params:
        prefix="{sample}",
        gcode=config["genetic_code"],
        outdir="results/prokka/{sample}"
    shell:
        """
        mkdir -p {params.outdir}
        prokka --cpus {threads}  --force {input} &> {log} || true
        """

我可以想到两种方法,都不是完美的。

第一个基本上就是您所做的:使用 bash 来变通。但是,我建议使用 -s 文件测试运算符。这样,您仍然会收到来自 prokka 的真正错误通知:

shell:
    """
    if [ -s {input} ]; then
       prokka --cpus {threads}  --force {input} > {log} 
    else;
       mkdir -p {params.outdir}
    """

或者,使用 checkpoints。这样就把所有的逻辑都放在了snakemake中,但是有点繁琐。类似这样的内容,但我不确定这是 100% 正确还是最佳版本:

checkpoint gather_outputs:
    input: expand("results/shovill/{sample}/contigs.fa", sample=samples)
    output: 'results/shovill/non_empty.list'
    shell: 
        with open(str(output), 'wt') as output_handle:
            for input_file in input:
                if os.path.exists(input_file) and os.path.getsize(input_file):
                    sample = os.path.basename(os.path.dirname(input_file))
                    output_handle.write(f"{sample}\n")

def get_prokka_outputs(wildcards):
    # this makes the rest of the workflow wait for the checkpoint
    chk_out = checkpoints.gather_outputs.get().output
    with open(chk_out) as sample_lines:
        samples = [line.strip() for line in sample_lines]

    return expand("results/prokka/{sample}", sample=samples)

您也许可以对每个样本进行检查,但我从未尝试过。

我想得越多,检查点对您的情况来说开销太大了。