本文介绍了Snakemake让检查点和聚合函数工作的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!
问题描述
我在执行Snakemake Aggregate命令时遇到问题。我希望获得一个给定的GTF文件,在GTF中查找单独的区域,如果找到,则将这些区域写入一个单独的文件。因此,我不确定每个输入GTF文件将创建的输出GTF文件的数量。为了解决此问题,我正在尝试使用蛇造检查站。
为此,我编写了一个名为collapse_gtf_file.py
的简短脚本,它只接受一个GTF文件,并根据找到的各个区域的数量生成N个文件。因此,如果给定的文件test_regions.gtf
有三个区域,它将分别生成test_regions_1.gtf, test_regions_2.gtf test_regions_3.gtf
。分离后,所有GTF文件都应转换为FASTA文件,并聚合。
但是,我无法使我的检查点命令工作。我可以让示例用例正常工作,但是当我尝试在这个检查点周围构建更大的管道时,它就会崩溃。到目前为止,我已经尝试按照这里的检查点教程https://snakemake.readthedocs.io/en/stable/snakefile/rules.html#dynamic-files
sample=config["samples"]
reference=config["reference"]
rule all:
input:
¦ expand("string_tie_assembly/{sample}.gtf", sample=sample),
¦ expand("string_tie_assembly_merged/merged_{sample}.gtf",sample=sample),
¦ #expand("split_gtf_file/{sample}", sample=sample),
¦ #expand("lncRNA_fasta/{sample}.fa", sample=sample),
¦ "combined_fasta/all_fastas_combined.fa",
¦ "run_CPC2/cpc_calc_output.txt"
rule samtools_sort:
input:
¦ "mapped_reads/{sample}.bam"
output:
¦ "sorted_reads/{sample}.sorted.bam"
shell:
¦ "samtools sort -T sorted_reads/{wildcards.sample} {input} > {output}"
rule samtools_index:
input:
"sorted_reads/{sample}.sorted.bam"
output:
"sorted_reads/{sample}.sorted.bam.bai"
shell:
"samtools index {input}"
rule generate_fastq:
input:
"sorted_reads/{sample}.sorted.bam"
output:
"fastq_reads/{sample}.fastq"
shell:
"samtools fastq {input} > {output}"
rule string_tie_assembly:
input:
"sorted_reads/{sample}.sorted.bam"
output:
"string_tie_assembly/{sample}.gtf"
shell:
"stringtie {input} -f 0.0 -a 0 -m 50 -c 3.0 -f 0.0 -o {output}"
rule merge_gtf_file_features:
input:
¦ "string_tie_assembly/{sample}.gtf"
output:
¦ "string_tie_assembly_merged/merged_{sample}.gtf"
shell:
¦ #prevents errors when there's no sequence
¦ """
¦ set +e
¦ stringtie --merge -o {output} -m 25 -c 3.0 {input}
¦ exitcode=$?
¦ if [ $exitcode -eq 1 ]
¦ then
¦ ¦ exit 0
¦ else
¦ ¦ exit 0
¦ fi
¦ """
#This is where the issue seems to arise from. Modeled after https://snakemake.readthedocs.io/en/stable/snakefile/rules.html#dynamic-files
checkpoint clustering:
input:
¦ "string_tie_assembly_merged/merged_{sample}.gtf"
output:
¦ clusters = directory("split_gtf_file/{sample}")
shell:
¦ """
¦ mkdir -p split_gtf_file/{wildcards.sample} ;
¦ python collapse_gtf_file.py -gtf {input} -o split_gtf_file/{wildcards.sample}/{wildcards.sample}
¦ """
rule gtf_to_fasta:
input:
¦ "split_gtf_file/{sample}/{sample}_{i}.gtf"
output:
¦ "lncRNA_fasta/{sample}/canidate_{sample}_{i}.fa"
wildcard_constraints:
¦ i="d+"
shell:
¦ "gffread -w {output} -g {reference} {input}"
def aggregate_input(wildcards):
checkpoint_output = checkpoints.clustering.get(**wildcards).output[0]
x = expand("lncRNA_fasta/{sample}/canidate_{sample}_{i}.fa",
¦ sample=wildcards.sample,
¦ i=glob_wildcards(os.path.join(checkpoint_output, "{i}.fa")).i)
return x
rule combine_fasta_file:
input:
¦ aggregate_input
output:
¦ "combined_fasta/all_fastas_combined.fa"
shell:
¦ "cat {input} > {output}"
Error Message:
InputFunctionException in combine_fasta_file:
WorkflowError: Missing wildcard values for sample
Wildcards:
在我看来,这表明我在Aggregate命令中调用上面的通配符的方式有问题,但我不知道是什么原因。如有任何提示,我们将不胜感激。
推荐答案
函数aggregate_input
需要变量wildcards.sample
,但您没有在rule combine_fasta_file
中指定的任何通配符。您可能希望在该规则中指定通配符,或者重构函数以使用全局变量(如果适用)。
这篇关于Snakemake让检查点和聚合函数工作的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持!