我是使用蛇制作的新手,在使用PICARD MergeSamFiles将bam文件合并为一个bam文件时遇到问题。我想将1_sorted.bam 2_sorted.bam ... 10_sorted.bam合并到一个具有目录名称的b​​am文件中。

import snakemake.io
import os.path

PICARD="/data/src/picard.jar"
(SAMPLES,)=glob_wildcards("bam/{sample}_sorted.bam")
NAME=os.path.dirname

def bam_inputs(wildcards):
    files = expand("bam/{sample}_sorted.bam", sample=SAMPLES)
    INPUT = "I="+files
    return INPUT

rule all:
    input: "bam/{NAME}.bam"

rule merge_bams:
    input: bam_inputs
    output: "bam/{NAME}.bam"
    params: mrkdup_jar="/data/src/picard.jar"
    shell: "java -Xmx16G -jar {params.mrkdup_jar} MergeSamFiles \
    {input} \
    O={output} \
    SORT_ORDER=coordinate \
    ASSUME_SORTED=false \
    USE_THREADING=true"


错误:

Building DAG of jobs...
WildcardError in line 12 of /data/data/Samples/snakemake-example/WGS-test/step3.smk:
Wildcards in input files cannot be determined from output files:
'NAME'


我不知道如何将所有bam文件合并为一个,也不知道如何将目录名设置为最终bam文件的变量。请指教。

更新:

import snakemake.io

PICARD="/data/src/picard.jar"
(SAMPLES,)=glob_wildcards("bam/{sample}_sorted.bam")
#NAME=os.path.dirname
NAME="test"

rule all:
    input: "bam/{name}.bam".format(name=NAME)

rule merge_bams:
    input: expand("bam/{sample}_sorted.bam",sample=SAMPLES)
    output: "bam/{name}.bam".format(name=NAME)
    params: mrkdup_jar="/data/src/picard.jar"
    shell: """java -Xmx16G -jar {params.mrkdup_jar} MergeSamFiles \
    {"I=" + input} \
    O={output} \
    SORT_ORDER=coordinate \
    ASSUME_SORTED=false \
    USE_THREADING=true """

ERROR:

RuleException in line 11 of /data/data/Samples/snakemake-example/WGS-test/step3.smk:
NameError: The name '"I=" + input' is unknown in this context. Please make sure that you defined that variable. Also note that braces not used for variable access have to be escaped by repeating them, i.e. {{print $1}}

MergeSamFiles \
I= sub1_sorted.bam I=sub2_sorted.bam I=sub3_sorted.bam \
O= sub.bam \
SORT_ORDER=coordinate \
        ASSUME_SORTED=false \
        USE_THREADING=true

最佳答案

让我们来考虑rule all。您需要向snakemake展示您实际希望作为目标构建的文件。没有通配符:只有明确的内容。您说应该是带有目录名称的b​​am文件吗?

rule all:
    input: f"bam/{NAME}.bam"


请注意,使用f字符串,我将通配符中的{NAME}转换为来自变量NAME的确切字符串值。您可以选择其他任何方式,例如"bam/{name}.bam".format(name=NAME)

接下来,请记住,现在“ all”规则中的{NAME}和“ merge_bams”规则中的{NAME}是不同的实体,因此它们没有任何共同之处。此外,通配符不必等于您在第6行中定义的NAME变量。为避免误解,我将以某种方式对通配符进行命名。

还有一件事:我不确定您在bam_inputs函数中正在做什么:

INPUT = "I="+files


expand函数的结果应足以为merge_bams规则指定输入。如果您需要为列表中的每个文件添加“ I =“,请尝试在shell:部分中进行操作:

rule merge_bams:
    input: bam_inputs
    output: "bam/{NAME}.bam"
    params: mrkdup_jar="/data/src/picard.jar"
    shell: f"""java -Xmx16G -jar {{params.mrkdup_jar}} MergeSamFiles
        {" ".join(["I=" + s for s in input])}
        O={{output}}
        SORT_ORDER=coordinate
        ASSUME_SORTED=false
        USE_THREADING=true"""

关于python - snakemake PICARD合并bam文件,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/56288018/

10-10 05:57