我使用snakemake管道运行GATK MarkDuplicate命令,其中包含来自不同读取组的多个输入bam文件。
rule mark_duplicates:
input:
get_dedup_input
output:
bam=temp("bams/{patient}.{sample_type}.markdups.bam"),
md5=temp("bams/{patient}.{sample_type}.markdups.bam.md5"),
metrics="qc/gatk/{patient}_{sample_type}_dup_metrics.txt"
conda:
"../envs/gatk.yml"
shell:
"""
gatk MarkDuplicates -I {input} -O {output.bam} -M {output.metrics} \
--CREATE_MD5_FILE true --ASSUME_SORT_ORDER "queryname"
"""get_dedup_input返回输入bam文件的列表。MarkDuplicates要求使用-I标志指定每个输入文件。如果我只有一个bam文件,我可以简单地编写-I {input},但是这失败了,因为这指定了-I file1.bam file2.bam,它需要是-I file1.bam -I file2.bam。将输入格式化以便将每个输入文件指定为-I [input file]的最佳方式是什么
下面是两个场景,说明如果我手动运行该命令,输入、输出和shell命令将是什么样子。为了简洁起见,我省略了一些非必要的MarkDuplicate标志:
一读组
Inputs: patient101.normal.rg1.bam
Output: patient101.normal.markdups.bam
Shell:
gatk MarkDuplicates -I patient101.normal.rg1.bam \
-O patient101.normal.markdups.bam \
-M metrics.txt两个读组
Inputs: patient101.normal.rg1.bam, patient101.normal.rg2.bam
Output: patient101.normal.markdups.bam
Shell:
gatk MarkDuplicates -I patient101.normal.rg1.bam \
-I patient101.normal.rg2.bam \
-O patient101.normal.markdups.bam \
-M metrics.txt发布于 2020-03-03 07:24:22
谢谢你提供最新的答复。
因此,最好的方法可能是使一个参数成为我们需要的字符串,如下所示:
rule mark_duplicates:
input:
get_dedup_input
output:
bam=temp("bams/{patient}.{sample_type}.markdups.bam"),
md5=temp("bams/{patient}.{sample_type}.markdups.bam.md5"),
metrics="qc/gatk/{patient}_{sample_type}_dup_metrics.txt"
conda:
"../envs/gatk.yml"
params:
input=lambda wildcards, input: " -I ".join(input)
shell:
"""
gatk MarkDuplicates -I {params.input} -O {output.bam} -M {output.metrics} \
--CREATE_MD5_FILE true --ASSUME_SORT_ORDER "queryname"
"""如果我们的输入是一个bam,那么params.input将是patient101.normal.rg1.bam,-I将如常地添加到前面。
如果我们有两个输入bam,那么lambda函数将-I放在它们之间,而shell命令在前面添加-I。
https://stackoverflow.com/questions/60498647
复制相似问题