如何使用 Snakemake 处理特定的文件对？

Question

我是 Snakemake 的新手，目前正在绞尽脑汁研究如何使用 Snakemake 管道处理特定的文件对（如果问题很蹩脚，我深表歉意）。
我有一组来自肿瘤和匹配的正常样本的 fastq 文件。在执行变异调用等操作时，每个肿瘤样本都必须与其匹配的正常样本一起处理。因此，我准备了一个配置文件，其中列出了以下样本：

sample_list:
  - sample: 1
    tumor: AO1_04_RN_1_T_4_S4
    control: AO2_07_C007558T1Wa_S37
  - sample: 2
    tumor: AO2_01_C007589T1FTa_S2
    control: AO2_07_C007589T1Wa_S34
  - sample: 3
    tumor: AO9_09_FM_1_T_7_S13 
    control: AO2_07_C007558T1Wa_S37

然而，我遇到的第一个问题是，这些样本每个都有 4 个 fastq.gz 文件，我需要将其连接成 2 个 - 例如，肿瘤样本 AO1_04_RN_1_T_4_S4 的文件是：

- AO1_04_RN_1_T_4_S4_L001_R1_001.fastq.gz
- AO1_04_RN_1_T_4_S4_L002_R1_001.fastq.gz
- AO1_04_RN_1_T_4_S4_L001_R2_001.fastq.gz
- AO1_04_RN_1_T_4_S4_L002_R2_001.fastq.gz

并且我需要第一个规则（例如

rule concat_fastq

）来简单地将每个样本的R1和R2文件连接成两个单独的fastq

AO1_04_RN_1_T_4_S4_R1.fastq.gz

和

AO1_04_RN_1_T_4_S4_R2.fastq.gz

，因为我需要提供这些作为输入对齐。

目前，我发现指定这样的输入：

r1 = expand("{path}/{sample}_L{num}_R1_001.fastq.gz",
        path = config["input_path"],
        num = ["001","002"],
    sample = [sample["tumor"] for sample in config["sample_list"]] + [sample["control"] for sample in config["sample_list"]]),
r2 = expand("{path}/{sample}_L{num}_R2_001.fastq.gz",
        path = config["input_path"],
        num = ["001","002"],
    sample = [sample["tumor"] for sample in config["sample_list"]] + [sample["control"] for sample in config["sample_list"]])

显然每次都会连接太多文件，而像这样的东西，不使用

expand

：

r1 = ["{path}/{sample}_L001_R1_001.fastq.gz", "{path}/{sample}_L002_R1_001.fastq.gz"],
r2 = ["{path}/{sample}_L001_R2_001.fastq.gz", "{path}/{sample}_L002_R2_001.fastq.gz"]

导致以下错误：

Wildcards in input files cannot be determined from output files:
'path'

因为我无法指定 {path} 是什么。

有人可以提供一些关于如何使用 Snakemake 处理配对文件的建议吗？我还希望我的配置文件结构正确，以便处理肿瘤-正常对。

非常感谢您的帮助

Answer 1

使用规则 all 来设置示例通配符可以稍微简化这一点。然后在

concat_fastq

中，我们再次使用expand来获取每个读取对的两个通道。我还可能建议使用 CSV 来指定样本和元数据，因为 CSV 更容易解析 (imo)。

configfile: "config.yaml"


TUMORS = [d["tumor"] for d in config["sample_list"]]
CONTROLS = [d["control"] for d in config["sample_list"]]
SAMPLES = TUMORS + CONTROLS
READS_PATH = "data/reads"


rule all:
    input:
        reads=expand(
            "results/concat_fastq/{sample}_R{num}.fastq.gz", sample=SAMPLES, num=[1, 2]
        ),


rule concat_fastq:
    input:
        # Use {{ }} to escape sample wildcard in expand.
        r1=expand(READS_PATH + "{{sample}}_L00{lane}_R1.fastq.gz", lane=[1, 2]),
        r2=expand(READS_PATH + "{{sample}}_L00{lane}_R2.fastq.gz", lane=[1, 2]),
    output:
        r1="results/concat_fastq/{sample}_R1.fastq.gz",
        r2="results/concat_fastq/{sample}_R2.fastq.gz",
    shell:
        """
        cat {input.r1} > {output.r1}
        cat {input.r2} > {output.r2}
        """