【问题标题】:Snakemake How to design downstream rules on empty output files?Snakemake 如何设计空输出文件的下游规则?
【发布时间】:2021-10-24 21:39:34
【问题描述】:

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

如果使用 PROKKA 对 0 字节文件运行基因组注释,则会返回错误:

 '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
        """

【问题讨论】:

  • 您能否说明您希望实现的目标?似乎您希望下游作业即使在它们的依赖项产生错误时也能运行?
  • 我想避免 snakemake 工作流在数百个作业中的一个产生错误时停止。这实际上是通过 -k 参数实现的

标签: snakemake fastq


【解决方案1】:

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

第一个基本上就是你所做的:使用 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)

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

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

【讨论】:

    猜你喜欢
    • 2020-03-21
    • 1970-01-01
    • 2021-12-29
    • 2021-09-29
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2022-08-18
    相关资源
    最近更新 更多