防止输入函数生成示例文件中不存在的文件

Question

我一直在研究一个我无法解决的 snakemake 问题。给定一个样本文件，例如：

tissue type replicate file
ear rep1 H3K4me3 00.data/chip_seq/H3K4me3/ear_H3K4me3_rep1.fastq
ear rep2 H3K4me3 00.data/chip_seq/H3K4me3/ear_H3K4me3_rep2.fastq
ear rep1 input 00.data/chip_seq/input/ear_input_rep1.fastq
ear rep2 input 00.data/chip_seq/input/ear_input_rep2.fastq
leaf rep1 H3K4me3 00.data/chip_seq/H3K4me3/ear_H3K4me3_rep1.fastq
leaf rep2 H3K4me3 00.data/chip_seq/H3K4me3/ear_H3K4me3_rep2.fastq
leaf rep1 input 00.data/chip_seq/input/ear_input_rep1.fastq
leaf rep2 input 00.data/chip_seq/input/ear_input_rep2.fastq
root rep1 input 00.data/chip_seq/input/ear_input_rep1.fastq
root rep2 input 00.data/chip_seq/input/ear_input_rep2.fastq

我用来输入此文件列表的 snakemake 函数 - 这里称为 get_chip_mods 生成实际上并不存在的通配符组合。因此，在这种情况下，get_chip_mods 会生成诸如 root_rep1_H3K4me3 之类的组合，即使示例中未指定所述文件。有没有办法阻止此函数生成示例文件中不存在的组合？

下面是我的管道的开始。

#Load Samples from the CSV file - index the important ones
samples = pd.read_csv(config["samples"], sep=' ').set_index(["tissue", "type", "replicate"], drop=False)
samples.index = samples.index.set_levels([i.astype(str) for i in samples.index.levels])  # enforce str in index



rule all:
    input:
    ¦   "00.data/reference/bowtie_idx.1.bt2",
    ¦   expand("00.data/trimmed_chip/{tissue}_{chip}_{replicate}_trimmed.fq" , tissue = samples["tissue"], 
chip = samples["type"], replicate = samples["replicate"]),


#This is where I believe I've been hitting issues. 
def get_chip_mods(wildcards):
    final_list = samples.loc[(wildcards.tissue, wildcards.type, wildcards.replicate), ["file"]].dropna()
    print(final_list)
    return final_list


rule trim_reads:
    input:
    ¦   get_chip_mods
    params:
    ¦   "00.data/trimmed_chip/log_files/{tissue}_{type}_{replicate}.log"
    output:
    ¦   "00.data/trimmed_chip/{tissue}_{type}_{replicate}_trimmed.fq"
    threads: 5
    message:"""Trimming"""
    shell:
    ¦   """
    ¦   java -jar /usr/local/apps/eb/Trimmomatic/0.36-Java-1.8.0_144/trimmomatic-0.36.jar \
    ¦   SE -threads {threads} -phred33 {input} {output} \
    ¦   ILLUMINACLIP:/scratch/jpm73279/04.lncRNA/02.Analysis/23.generate_all_metaplots/00.data/adapter.fa:2:30:10 \
    ¦   LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36
    ¦   """

我收到的错误如下

KeyError: 
Wildcards:
tissue=root
type=H3K4me3
replicate=rep1

Answer 1

解决歧义的唯一关键是包含可能组合的文件。这意味着您的脚本编写方式不应取决于可能的组合。

解决该问题的一种方法是将 all 规则中的三个通配符替换为一个通配符 {tissue_type_replicate}，并使用 python 函数生成可能的值。这将为 Snakemake 提供它需要生成哪些文件的信息。您也可以在其他规则中进行相同的更改（这是最简单可行的解决方案，只要您不需要 [=15] 中的 {tissue} {type} 和 {replicate} 的确切值=] 节）。无论如何，您仍然可以在其他规则中保留 {tissue} {type} 和 {replicate} 通配符：Snakemake 应该找到匹配项。

Answer 2

错误与规则 all 中的 expand 函数有关。默认情况下，该函数将使用 python itertools product 来生成所有可能的通配符组合。其中一些组合不存在于您的数据帧索引中，因此不存在错误。

expand 但是允许您自定义通配符的组合方式，因此您可以按以下方式重写函数来解决问题。

expand("00.data/trimmed_chip/{tissue}_{chip}_{replicate}_trimmed.fq".split(), zip, tissue = samples["tissue"], chip = samples["type"], replicate = samples["replicate"])

Source

防止输入函数生成示例文件中不存在的文件

Preventing input function from generating files not present in sample file

bioinformatics

snakemake