📚 课程目录第4章规则进阶

4、规则进阶：通配符的使用、多输出规则、日志文件配置

Snakemake 的规则一旦碰上多样本、多文件，就会变得特别有意思。说白了，通配符就是让规则能批量干活的核心。我刚开始学的时候，总觉得通配符像个黑魔法，后来踩过几次坑，才明白它其实挺朴素的——就是把样本名字当作变量，在输入输出里来回传递。

这一章我们聊三个进阶点：通配符怎么玩、多输出规则怎么写、日志怎么配。嗯，这三个东西在实际项目中几乎天天用，你掌握好了，写 Snakefile 就能省下一大半体力。

通配符：规则的灵魂

先看一个最简单的例子。假设你有一堆样本，每个样本都有一个原始测序文件 data/sampleA.fastq，你想把它们分别比对到参考基因组上。用通配符，规则可以这样写：

rule align:
    input:
        "data/{sample}.fastq",
        "reference/genome.fa"
    output:
        "results/{sample}.bam"
    shell:
        "bwa mem {input[1]} {input[0]} | samtools sort -o {output}"

看到 {sample} 了吗？它就是通配符。Snakemake 会从目标文件或 expand() 函数中自动推导出通配符的具体值。比如你执行 snakemake results/sampleA.bam，它就知道 sample=A，然后自动把 input 里的 {sample} 替换成 sampleA。

我个人习惯给通配符起有语义的名字，比如 {sample}、{chr}、{tool}，而不是 {x}、{y}。这样你看到输出路径，一眼就知道是哪个样本、哪个步骤。我曾经在项目里用 {a}、{b} 命名，结果两个月后自己都看不懂哪个对哪个，改起来特别痛苦。

通配符的约束：让规则更可控

如果通配符的值范围很大（比如所有样本名字），你可以在规则里用 wildcard_constraints 来限定它的格式：

rule align:
    input:
        "data/{sample}.fastq"
    output:
        "results/{sample}.bam"
    wildcard_constraints:
        sample = "[A-Za-z0-9_]+"
    shell:
        "bwa mem {input} > {output}"

这样做的好处是：

避免通配符匹配到不该匹配的路径（比如含斜杠、点等特殊字符）
让规则更清晰，别人看代码就知道样本名的合法字符集
有时能防止意外的文件覆盖

遇到跨路径的通配符（比如样本名里包含 sample 子目录），我建议你直接用 wildcard_constraints 把路径分隔符排除掉，否则 Snakemake 可能会把目录当成通配符的一部分，造成匹配混乱。

多输出规则：一个规则生成多个文件

很多生物信息学工具会同时输出多个文件。比如 GATK 的 HaplotypeCaller 会生成一个 GVCF 文件和一个索引文件；又或者一个质控步骤会产生 report.html 和 report.json。Snakemake 支持一个规则写多个输出文件，并且输出文件之间可以共享通配符。

写法很简单，把 output 写成一个列表：

rule quality_control:
    input:
        "raw/{sample}.fastq"
    output:
        "qc/{sample}/report.html",
        "qc/{sample}/report.json"
    shell:
        "fastqc -o qc/{sample} {input}"

注意，这里的两个输出文件都使用了同一个通配符 {sample}。当 Snakemake 需要生成 qc/sampleA/report.html 时，它自动推导出 sample = sampleA，然后检查或生成第二个输出文件 qc/sampleA/report.json。

重要提示： 如果多个输出文件中的通配符不一致，Snakemake 会报错。比如 output: "a/{x}.txt", "b/{y}.txt"，它无法保证 {x} 和 {y} 的值有唯一关联。所以多输出规则中，所有输出文件最好共享同一套通配符。

我在项目中遇到过一种情况：工具输出文件名字固定，但样本名不同。比如 coverage_metrics.txt 和 variant.vcf，没有通配符可以串联。这时候我一般用 output: "results/{sample}/coverage_metrics.txt", "results/{sample}/variant.vcf" 来强制保持通配符一致。

多输出与 aggregation

多输出规则还经常用于聚合步骤。比如你想把所有样本的质控报告收集到一个表格里：

rule collect_qc:
    input:
        expand("qc/{sample}/report.json", sample=SAMPLES)
    output:
        "results/all_qc.json"
    run:
        # 读取所有 JSON，合并输出
        ...

这里 expand() 会生成所有样本 report.json 的列表，Snakemake 会自动推导需要先运行 quality_control 规则来生成这些文件。多输出与聚合搭配，能组成高效的 DAG 图。

日志文件配置：别让错误信息淹没在终端

你想想看，一个流水线跑几百个样本，如果某个样本报错，你总不能去终端翻几千行输出吧？配置日志文件，就是把每个规则的执行日志单独存到文件里，方便事后排查。

Snakemake 在规则里提供了 log 指令，用法和 input、output 类似：

rule align:
    input:
        "data/{sample}.fastq"
    output:
        "results/{sample}.bam"
    log:
        "logs/align/{sample}.log"
    shell:
        "bwa mem ref.fa {input} 2> {log}"

注意：log 文件不会被 Snakemake 视为输出文件，所以不会触发下游规则的输入依赖。它只是纯粹用来存日志的路径。如果某个规则失败了，你就可以直接去看对应的日志文件，而不是满屏刷屏。

我曾经犯过一个错误：把日志文件写到了 output 里。结果下游规则想依赖这个日志文件，导致 Snakemake 误以为日志是重要的中间产物，还试图用它来触发重新运行。嗯，所以记住：日志只走 log，别走 output。

日志文件的路径规范

我个人的习惯是：把日志放在一个独立的 logs/ 目录下，按规则名和样本名组织。比如：

logs/align/sampleA.log
logs/qc/sampleA.log
logs/variant_calling/sampleA.log

这样当你需要定位某个样本的某个步骤时，路径清晰可辨。如果项目规模很大，还可以按日期或批次再分一层：logs/2025-04/align/sampleA.log。

日志文件与重跑逻辑

Snakemake 默认不会因为日志文件的变化而重新运行规则——日志文件只写，不读取。但是，如果你在 shell 命令里把日志写在标准输出（> {log}），并且日志文件提前存在，Snakemake 会不会认为规则已经完成？答案是不会。`log` 文件不参与判据，`output` 文件才是核心。

避坑提醒： 别把 log 路径跟 output 路径写重叠。Snakemake 会优先把匹配 output 的通配符值传给 log，但如果某个 log 文件恰好跟某个 output 文件同名，Snakemake 的行为可能让你意外。我建议日志文件名加个 .log 后缀，跟输出文件区分开。

本章知识结构图

下面这张图帮你梳理了通配符、多输出、日志之间的关系，以及它们如何融入 Snakemake 的数据流：

我的小技巧： 刚开始练习时，你可以先写一个单样本的硬编码规则（不包含通配符），跑通后再把样本名替换成 {sample}。这样能确保命令本身没问题，避免通配符和命令错误混在一起排查。毕竟，调试 Snakemake 最烦的就是“不知道是通配符出事还是工具出事”。

三者搭配：实际案例

最后给一个完整的代码片段，把通配符、多输出、日志都串起来：

# 定义样本列表
SAMPLES = ["A", "B", "C"]

# 比对规则
rule bwa_align:
    input:
        fq = "data/{sample}.fastq",
        ref = "reference/genome.fa"
    output:
        bam = "bam/{sample}.bam",
        bai = "bam/{sample}.bam.bai"
    log:
        "logs/bwa_align/{sample}.log"
    shell:
        "(bwa mem {input.ref} {input.fq} "
        "| samtools sort -o {output.bam} && "
        "samtools index {output.bam}) 2> {log}"

# 聚合质控
rule collect_stats:
    input:
        expand("bam/{sample}.bam.bai", sample=SAMPLES)
    output:
        "stats/coverage_summary.txt"
    log:
        "logs/collect_stats/coverage_summary.log"
    shell:
        "echo 'Sample\tReads' > {output} && "
        "for bam in {input}; do "
        "  samtools stats $bam >> {output}; done 2> {log}"

这个例子中，bwa_align 规则输出了 BAM 和 BAI 两个文件，并且把 stderr 都定向到了日志。之后 collect_stats 规则通过 expand() 聚合所有样本的 BAI 文件，生成汇总统计。你运行 snakemake stats/coverage_summary.txt，Snakemake 会自动推导出所有中间文件，包括日志。

如果某个样本的比对失败了，你直接去 logs/bwa_align/ 目录下找对应的 .log 文件，就能看到错误原因，而不需要重跑整个流水线。我个人觉得，日志配置这个习惯，一旦养成了，就再也回不去了——太省心了。