Snakemake，RNA-seq :根据分析过的样本的特征，我如何执行管道的一个子部分或另一个子部分？

Question

我正在使用snakemake设计一个RNAseq-数据分析管道。虽然我设法做到了这一点，但我希望使我的管道尽可能具有适应性，并使它能够在相同的分析运行中处理单读( SE )数据或配对( PE )数据，而不是在一次运行中分析SE数据，而在另一次运行中处理PE数据。

我的管道应该是这样设计的：

提供一个文件(SE数据)或两个文件(PE数据)的

• 数据集

注意:A的所有规则都有一个输入和一个输出，B的所有规则都有两个输入和两个输出，它们各自的命令如下：

somecommand -i {input} -o {output}

• 2

• 1输入：输入：somecommand -i1 {input1} -i2 {input2} -o1 {output1} -o2 {output2}

注2:除输入/输出的差异外，A组和B组的所有规则都有相同的命令、参数/等。

换句话说，我希望我的管道能够根据示例在执行规则A或规则B之间切换，方法是在开始时在配置文件中向它提供有关示例的信息(示例1为SE，示例2为PE.或者要求snakemake对数据集下载后的文件数量进行计数，以便为每个示例选择合适的下一组规则。如果你发现了另一种方法，我们欢迎你告诉你。

我考虑过使用检查点、输入函数和if/ with语句，但我没有设法解决这些问题。

你有什么提示/建议/方法让“切换”发生吗？

Answer 1

如果事先知道布局，那么最简单的方法是将其存储在某个变量中，类似于以下内容(或者将其从配置文件中读取到字典中)：

layouts = {"sample1": "paired", "sample2": "single", ... etc}

然后，您可以这样“合并”您的规则(我猜您是在讨论修整和对齐，这就是我的例子)：

ruleorder: B > A

rule A:
    input:
        {sample}.fastq.gz
    output:
        trimmed_{sample}.fastq.gz
    shell:
        "somecommand -i {input} -o {output}"

rule B:
    input:
        input1={sample}_R1.fastq.gz,
        input2={sample}_R2.fastq.gz
    output:
        output1=trimmed_{sample}_R1.fastq.gz,
        output2=trimmed_{sample}_R2.fastq.gz
    shell:
        "somecommand -i1 {input.input1} -i2 {input.input2} -o1 {output.output1} -o2 {output.output2}"


def get_fastqs(wildcards):
    output = dict()
    if layouts[wildcards.sample] == "single":
        output["input"] = "trimmed_sample2.fastq.gz"
    elif layouts[wildcards.sample] == "paired":
        output["input1"] = "trimmed_sample1_R1.fastq.gz"
        output["input2"] = "trimmed_sample1_R2.fastq.gz"
    return output


rule alignment:  
    def input:
        unpack(get_fastqs)
    def output:
        somepath/{sample}.bam
    shell:
        ...

这里发生了很多事情。

• 首先您需要一个规则命令，因此snakemake知道如何处理模糊的情况，
• 规则A和B都必须存在(除非您对输出文件做了一些讨厌的事情)。对齐规则
• 需要一个输入函数来确定它需要哪个输入。

一些自我推销:我做了一个蛇形管道，它做了很多事情，包括RNA-seq和在线下载样本，并自动确定它们的布局(单端和成对端)。请看一看它是否解决了您的问题：https://vanheeringen-lab.github.io/seq2science/content/workflows/rna_seq.html

编辑

1. 当你说“合并”规则时，你是指规则A，B和对齐吗?

我的措辞不太清楚。通过合并，我的意思是“将单端逻辑、配对逻辑和配对逻辑合并在一起，这样您就可以继续使用单个规则(例如，计数表，您可以命名它)。

1. 规则顺序:为什么选择B>A？确保成对的样本不会在单端规则中运行？

一点儿没错!当规则需要trimmed_sample1_R1.fastq.gz时，Snakemake如何知道示例的名称？示例的名称是sample1，还是sample1_R1？这两者都有可能，这让snakemake抱怨说它不知道如何解决这个问题。当您添加规则命令时，您告诉Snakemake，当它不清楚时，按此顺序解决。

1. (对齐规则中的命令)需要1或2个输入。我打算在params指令中使用if/ the来选择输入。我这么想对吗？(我认为你在制作过程中也是这么做的)

是的，这就是我们解决问题的方法。我们是这样做的，因为我们希望每条规则都有自己的环境。如果您不使用单独的conda环境进行对齐，那么您可以使用更干净/更漂亮的环境，如下所示

rule alignment:  
    input:
        unpack(get_fastqs)
    output:
        somepath/{sample}.bam
    run:
        if layouts[wildcards.sample] == "single":
            shell("single-end command")
        if layouts[wildcards.sample] == "paired":
            shell("paired-end command")

我觉得这个选项比我们在seq2science管道中所做的要清晰得多。然而，在seq2science管道中，我们支持许多不同的对齐器，它们都有不同的conda环境，因此不能使用run指令。