运行 snakemake 迭代规则

Run snakemake rules iteratively

所以我以为我终于掌握了 snakemake,但是当尝试 运行 几个不同的数据文件时,我意识到它并不像我想的那样工作。这是 Snakefile:

import pandas as pd

configfile: "config.json"
experiments = pd.read_csv(config["experiments"], sep = '\t')
experiments['Name'] = [filename.split('/')[-1].split('.fa')[0] for filename in experiments['Files']]

rule all:
    input:
        expand("{output}/Preprocess/Trimmomatic/quality_trimmed_{name}{fr}.fq", output = config["output"],
            fr = (['_forward_paired', '_reverse_paired'] if experiments["Files"].str.contains(',').tolist() else ''),
               name = experiments['Name'])

rule preprocess:
    input:
        experiments["Files"].str.split(',')
    output:
        expand("{output}/Preprocess/Trimmomatic/quality_trimmed_{name}{fr}.fq", output = config["output"],
            fr = (['_forward_paired', '_reverse_paired'] if experiments["Files"].str.contains(',').tolist() else ''),
               name = experiments['Name'])
    threads:
        config["threads"]
    run:
        shell("python preprocess.py -i {reads} -t {threads} -o {output} -adaptdir MOSCA/Databases/illumina_adapters -rrnadbs MOSCA/Databases/rRNA_databases -d {data_type}",
            output = config["output"], data_type = experiments["Data type"].tolist(), reads = ",".join(input))

这是配置文件:

{
  "output": "test_snakemake",
  "threads": 14,
  "experiments": "experiments.tsv"
}

这是实验文件

Files   Sample  Data type   Condition
path/to/mg_R1.fastq,path/to/mg_R2.fastq Sample  dna
path/to/a/0.01/mt_0.01a_R1.fastq,path/to/a/0.01/mt_0.01a_R2.fastq   Sample  rna c1
path/to/b/0.01/mt_0.01b_R1.fastq,path/to/b/0.01/mt_0.01b_R2.fastq   Sample  rna c1
path/to/c/0.01/mt_0.01c_R1.fastq,path/to/c/0.01/mt_0.01c_R2.fastq   Sample  rna c1
path/to/a/1/mt_1a_R1.fastq,path/to/a/1/mt_1a_R2.fastq   Sample  rna c2
path/to/b/1/mt_1b_R1.fastq,path/to/b/1/mt_1b_R2.fastq   Sample  rna c2
path/to/c/1/mt_1c_R1.fastq,path/to/c/1/mt_1c_R2.fastq   Sample  rna c2
path/to/a/100/mt_100a_R1.fastq,path/to/a/100/mt_100a_R2.fastq   Sample  rna c3
path/to/b/100/mt_100b_R1.fastq,path/to/b/100/mt_100b_R2.fastq   Sample  rna c3
path/to/c/100/mt_100c_R1.fastq,path/to/c/100/mt_100c_R2.fastq   Sample  rna c3

我想要做的是让预处理规则分别处理每一行。我认为那是 shell 解释命令的方式,它会 运行 命令 python preprocess.py -i path/to/mg_R1.fastq,path/to/mg_R2.fastq -t 14 -o test_snakemake -adaptdir MOSCA/Databases/illumina_adapters -rrnadbs MOSCA/Databases/rRNA_databases -d dna,而不是尝试将所有行和 运行 同时连接到所有样本 python preprocess.py -i path/to/mg_R1.fastq,path/to/mg_R2.fastq,path/to/a/0.01/mt_0.01a_R1.fastq,path/to/a/0.01/mt_0.01a_R2.fastq,path/to/b/0.01/mt_0.01b_R1.fastq,path/to/b/0.01/mt_0.01b_R2.fastq,... -t 14 -o test_snakemake -adaptdir MOSCA/Databases/illumina_adapters -rrnadbs MOSCA/Databases/rRNA_databases -d dna rna rna rna rna rna rna rna rna rna.

如何让 snakemake 分别考虑每一行?

这是一个很常见的错误。要记住的是,规则应该适用于单个样本。 Snakemake 将采用您的路径(使用通配符)并根据规则生成特定的作业。你已经写了一些接受所有输入和所有输出的东西,然后我推测,preprocess.py 期望一个 input/output.

相反,一次考虑一个文件。对于输出 "{output}/Preprocess/Trimmomatic/quality_trimmed_{name}{fr}.fq",您如何生成该文件?您必须使用名称作为键来匹配实验数据框中的输入文件。

def preprocess_input(wildcards):
    # get files with matching names
    df = experiments.loc[experiments['Name'] == wildcards.name, 'Files']
    # get first value (in case multiple) and split on commas
    return df.iloc[0].split(',')

rule preprocess:
    input:
        preprocess_input
    output:
        "{output}/Preprocess/Trimmomatic/quality_trimmed_{name}{fr}.fq"
    threads:
        config["threads"]
    shell:
        'python preprocess.py -i {reads} -t {threads} -o {config[output]} ...'

它使用输入函数从输出文件中找到正确的输入文件。它并不完美,但应该能让您朝着正确的方向前进。