Snakemake从每个样本目录问题输入fastq文件，用于元基因组学分析
Snakemake从每个样本目录问题输入fastq文件，用于元基因组学分析
提问于 2021-05-20 01:55:52
我正在开发一个新的snakemake元基因组学管道，用于裁剪fastq文件，并通过kraken运行它们。每个示例都有一个包含正向和反向读取的目录。
Sample_1/r1_paired.fq.gz
Sample_1/r2_paired.fq.gz
Sample_2/r1_paired.fq.gz
Sample_2/r2_paired.fq.gz我提供了一个样本表，用户可以上传，其中包含样本名称和已读名称。我使用pandas来解析样例工作表，并提供snakefile所需的名称。这是我的蛇形文件。
#Extract sample names from CSV
import pandas as pd
import os
df = pd.read_csv("sample_table_test.csv")
print(df)
samples = df.library.to_list()
print("Samples being processed:", samples)
R1 = df.r1_file.to_list()
R2 = df.r2_file.to_list()
print(R1,R2)
rule all:
input:
expand("{sample}.bracken", sample=samples),
#Trimmomatic to trim paired end reads
rule trim_reads:
input:
"{sample}/{R1}",
"{sample}/{R2}",
output:
"{sample}/{R1}_1_trim_paired.fq.gz",
"{sample}/{R2}_2_trim_paired.fq.gz",
conda:
"env.yaml",
shell:
"trimmomatic PE -threads 8 {input} {input} {output} {output} SLIDINGWINDOW:4:30 LEADING:2 TRAILING:2 MINLEN:50"
#Kraken2 to bin reads and assign taxonomy
rule kraken2:
input:
"{sample}/{R1}_1_trim_paired.fq.gz",
"{sample}/{R2}_2_trim_paired.fq.gz",
output:
"{sample}_report.txt",
"{sample}_kraken_cseqs#.fq",
conda:
"env.yaml",
shell:
"kraken2 --gzip-compressed --paired --classified-out {output} {input} {input} --db database/minikraken2_v1_8GB/ --report {sample}_report.txt --threads 1"
#Bracken estimates abundance of a species within a sample
rule bracken:
input:
"{sample}_report.txt",
output:
"{sample}.bracken",
conda:
"env.yaml",
shell:
"bracken -d database/minikraken2_v1_8GB/ -i {input} -o {output} -r 150"我收到了下面的错误，并且一直在努力寻找一种更好的方法来编写我的snakefile来避免这个问题。这里的任何帮助都将不胜感激。
WildcardError in line 19 of /Metagenomics/Metagenomics/snakemake/Snakefile:
Wildcards in input files cannot be determined from output files:
'R1'谢谢!
回答 1
Stack Overflow用户
发布于 2021-05-20 06:39:07
问题出在您的rule kraken2中
rule kraken2:
input:
"{sample}/{R1}_1_trim_paired.fq.gz",
"{sample}/{R2}_2_trim_paired.fq.gz",
output:
"{sample}_report.txt",
"{sample}_kraken_cseqs#.fq",规则中的所有通配符应由output部分确定。每条规则的逻辑是它提供某些文件作为可能的输出。在您的示例中，规则提供了文件"{sample}_report.txt"和"{sample}_kraken_cseqs#.fq"，其中{sample}成为一种自由度，并被替换为将模式解析为文件名的特定值。现在，Snakemake可以确定此规则的输入，但前提是它拥有所有信息。好的，{sample}的值是根据输出定义的，但是{R1}和{R2}的值是什么呢
您有几个选项。第一种方法是在output:中的某个位置定义这些值。看起来那不是你的案子。第二个选项是全局定义这些值(您可能正在尝试这样做)：
R1 = df.r1_file.to_list()
R2 = df.r2_file.to_list()在这种情况下，{R1}和{R2}不应该是通配符，而是expand函数的参数：
rule kraken2:
input:
expand("{{sample}}/{R1}_1_trim_paired.fq.gz", R1=R1),
expand("{{sample}}/{R1}_1_trim_paired.fq.gz", R2=R2)
output:
"{sample}_report.txt",
"{sample}_kraken_cseqs#.fq",或者更好：
expand("{{sample}}/{R}_1_trim_paired.fq.gz", R=R1+R2)请注意，通配符{sample}现在必须放在双花括号中，以区别于expand函数的参数。
还有其他选项，如从lambda wildcards: ...等其他vildcards的值中解析出{R1}的值，但我猜这不是您所需要的。
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接：
https://stackoverflow.com/questions/67608482
复制相关文章
相似问题
腾讯云开发者
