首页 > 解决方案 > 使用snakemake通配符重命名文件

问题描述

我无法解决我在 snakemake 中遇到的问题。我的样本目前被命名为“1-2-Brain_187_006_S77_L002_R1_001.fastq.gz”。我想最终将它们重命名为更短的名称,例如“1_2_Brain_S77_L002_R1”,然后使用扩展名“_trim.fastq.gz”作为我的规则。我正在使用 bbduk 进行修剪。对于我的输入,我想调用我的字典列表 allSamples。然后我想访问每个字典中的值。具体来说,“shortName1”和“shortName2”值。我的问题是在我的干燥中,它显示整个列表作为单次运行的输入。我不确定如何让它注册每个元素都是它自己的运行。当我实际有 114 个文件名时,我将显示 3 个文件名作为示例。所以,我希望我的空运行有 114 个用于修剪作业的计数。

配置文件

{
   "allSamples" : ['1_2_Brain_S77_L002', '10_4_Kidney_S82_L002', '11_4_BB_S105_L002' ......],

   "1_2_Brain_S77_L002":{
        "sampleName1": "1-2-Brain_187_006_S77_L002_R1_001.fastq.gz",
        "sampleName2": "1-2-Brain_187_006_S77_L002_R2_001.fastq.gz",
        "shortName1": "1_2_Brain_S77_L002_R1",
        "shortName2": "1_2_Brain_S77_L002_R2",
        "stemName": "1_2_Brain_S77_L002"
        }, ....
}

我正在获取位于 rawReads/ 中的文件并将新修剪的文件存储在 trimmedReads/ 中。

蛇文件

configfile: "refs/config.json"

# variables
sampleDict = config["allSamples"]
sampleNames1 = [config[i]["sampleName1"] for i in sampleDict]
sampleNames2 = [config[i]["sampleName2"] for i in sampleDict]
shortNames1 = [config[i]["shortName1"] for i in sampleDict]
shortNames2 = [config[i]["shortName2"] for i in sampleDict]

rule all:
    input: 
        expand("trimmedReads/{trim1}_trim.fastq.gz", trim1 = shortNames1),
        expand("trimmedReads/{trim2}_trim.fastq.gz", trim2 = shortNames2)

rule trim:
    input:
        R1 = expand("rawReads/{sample1}", sample1 = sampleNames1),
        R2 = expand("rawReads/{sample2}", sample2 = sampleNames2)
    output:
        trim1 = expand("trimmedReads/{trim1}_trim.fastq.gz", trim1 = shortNames1),
        trim2 = expand("trimmedReads/{trim2}_trim.fastq.gz", trim2 = shortNames2)
    shell:
        """
        bbduk.sh in1={input.R1} in2={input.R2} out1={output.trim1} out2={output.trim2} ref=ref/adapters.fa ktrim=r k=23 mink=11 hdist=1 tpe tbo
        """

当我进行空跑时,我得到了这个。

Building DAG of jobs...
Job counts:
    count   jobs
    1   all
    1   trim
    2

[Mon May 24 22:42:36 2021]
rule trim:
    input: rawReads/1-2-Brain_187_006_S77_L002_R1_001.fastq.gz, rawReads/10-4-Kidney_127_066_S82_L002_R1_001.fastq.gz, rawReads/11-4_BB_041_152_S105_L002_R1_001.fastq.gz, ...
    output: trimmedReads/1_2_Brain_S77_L002_R1_trim.fastq.gz, trimmedReads/10_4_Kidney_S82_L002_R1_trim.fastq.gz, trimmedReads/11_4_BB_S105_L002_R1_trim.fastq.gz, ...
    jobid: 1


bbduk.sh in1=rawReads/1-2-Brain_187_006_S77_L002_R1_001.fastq.gz rawReads/10-4-Kidney_127_066_S82_L002_R1_001.fastq.gz rawReads/11-4_BB_041_152_S105_L002_R1_001.fastq.gz ... out1=trimmedReads/1_2_Brain_S77_L002_R1_trim.fastq.gz trimmedReads/10_4_Kidney_S82_L002_R1_trim.fastq.gz trimmedReads/11_4_BB_S105_L002_R1_trim.fastq.gz ... ref=ref/adapters.fa ktrim=r k=23 mink=11 hdist=1 tpe tbo
        

[Mon May 24 22:42:36 2021]
localrule all:
    input: trimmedReads/1_2_Brain_S77_L002_R1_trim.fastq.gz, trimmedReads/10_4_Kidney_S82_L002_R1_trim.fastq.gz, trimmedReads/11_4_BB_S105_L002_R1_trim.fastq.gz, ...
    jobid: 0

Job counts:
    count   jobs
    1   all
    1   trim
    2
This was a dry-run (flag -n). The order of jobs does not reflect the order of execution.

标签: pythonbioinformaticssnakemake

解决方案


的通配符rule test:是一个空字典。此规则中没有wildcards.sample指定通配符值。每个通配符都应在该output:部分中指定,并且该部分对于该规则是空的。实际上rule test:,除非您将其明确指定为目标,否则您绝对没有任何效果:如果没有指定任何输出,Snakemake 只会忽略这个无用的规则,不会产生任何结果。

我猜这些文件["rawReads/1_2_Brain_S77_L002", "rawReads/17_6_Brain_S83_L002"]已经存在,所以 Snakemake 发现目标存在于磁盘上,并且什么都不做,产生“无输出”。

我不明白您所说的“最终将它们重命名为更短的名称”是什么意思,但这里有一个如何复制文件的方法。将此作为“如何使用通配符访问我的示例名称”的模式:

rule all:
    input: ["path_to_target/foo_SampleName1_bar", "path_to_target/foo_SampleName2_bar"]
    # List the files you expect to get as a target

rule copy:
    input:
        "path_to_source/blablabla_{sample}_bazz"
    output:
        "path_to_target/foo_{sample}_bar"
    shell:
        "echo {input}; cp {input} {output}"

这个怎么运作:

  1. Snakemake 发现它需要生成一些文件(如果这些文件是“path_to_target/foo_SampleName1_bar”、“path_to_target/foo_SampleName2_bar”)。
  2. Snakemake 发现rule copy:声明输出(如果{sample}用 value替换"SampleName1")与文件名匹配"path_to_target/foo_SampleName1_bar"
  3. 如果文件"path_to_source/blablabla_SampleName1_bazz"存在,Snakemake 就满足了,并且知道如何生成文件"path_to_target/foo_SampleName1_bar"
  4. {sample}对值重复步骤 2、3 "SampleName2"
  5. 现在它知道rule copy:应该运行两次:每个文件运行一次。
  6. 所有的依赖都解决了,Snakemake 就可以启动管道了。

推荐阅读