Snakemake规则（rule）详解：输入输出、shell命令、run块

我个人习惯把 Snakemake 的规则看作是流水线上的工位——你告诉它要什么原料，它取出原料，干完活，再吐出成品。而 rule 里最核心的三个东西就是 input、output，以及干活用的 shell 或 run。嗯，咱们一个一个拆开看。

1. 规则的基本骨架

一个最朴素的 rule 长什么样？你看下面这例子：

rule align_reads:
    input:
        "data/reads.fastq",
        "reference/genome.fa"
    output:
        "results/aligned.bam"
    shell:
        "bwa mem {input[1]} {input[0]} > {output}"

这里 input 是两个文件，output 一个文件，shell 里用 {input[1]} 和 {input[0]} 引用输入。注意顺序是从 0 开始，和 Python 列表一样。我在项目中见过有人把顺序搞反，结果出来的 BAM 文件全是乱的——因为参考基因组和测序数据放反了。这种错误很难排查，因为能跑完但结果不对。

2. input 和 output：文件的声明方式

input 说白了就是规则的前置依赖。Snakemake 会检查 output 文件是否存在，如果存在就跳过规则；如果不存在，它会回溯去查找生成这些 input 的其他规则。这就是所谓的“按需执行”。

output 也有讲究——它必须能够被通配符解析。举个例子：

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

这里的 {sample} 就是通配符。Snakemake 根据最终目标文件的 sample 值，自动把通配符替换成实际名字。要注意的是，output 里的通配符必须和 input 里的保持一致，否则 Snakemake 会报错无法匹配。我曾经花了两小时查一个 bug，最后发现是一个字母大小写不一致导致的——Windows 下不区分大小写，但在 Linux 上跑就崩了。

特性	input	output
通配符使用	通常与 output 通配符一致	定义通配符的地方
文件存在性	Snakemake 会查找是否有 rule 能生成	作为后续规则的依赖
多文件写法	可以用列表或字典	同左
自动创建目录	否	不，output 目录需要手动创建

3. shell 命令块：直接运行命令行

shell 块是最直接的方式——你平时怎么在终端敲命令，就怎么写。Snakemake 把大括号里的变量替换成实际值，然后交给系统的 shell 执行。比如：

rule fastqc:
    input:
        "data/{sample}.fastq"
    output:
        "results/fastqc/{sample}_fastqc.html",
        "results/fastqc/{sample}_fastqc.zip"
    shell:
        "fastqc -o results/fastqc {input}"

这里 output 有两个文件，Shell 命令会生成它们，Snakemake 会检查两者是否都存在。如果只有一个文件生成，规则会失败并报错。

嗯，这里要注意引号的嵌套。如果 shell 命令本身包含单引号或双引号，建议用双引号包裹 shell 块，命令内的引号用单引号或转义。我曾经因为没注意引号，导致命令里的 awk '{print $1}' 被 Python 字符串解析搞崩，调试了半天才发现。

4. run 块：用 Python 脚本执行

run 块允许你直接写 Python 代码。Snakemake 把 input、output、params 等变量作为 Python 对象传进来。比如：

rule merge_bams:
    input:
        bams=["results/{sample}_1.bam", "results/{sample}_2.bam"],
        bai=["results/{sample}_1.bam.bai", "results/{sample}_2.bam.bai"]
    output:
        "results/{sample}.merged.bam"
    run:
        shell("samtools merge {output} " + " ".join(input.bams))

run 块内你可以做复杂的事情：调用 Python 库，动态构建命令，甚至进行数据预处理。但要注意，run 块默认是在子进程执行的，如果抛出异常，Snakemake 会终止整个流程。

特性	shell	run
执行方式	通过系统 shell 执行命令	在 Python 环境中执行代码
变量引用	使用大括号 {input} 替换	直接使用 Python 变量
调试难度	低，命令本身可复制执行	中，需要理解 Python 上下文
性能开销	每次调用 shell 子进程	Python 进程内执行，略快

我说个真实的教训：有一次我需要在 run 块里调用外部程序，但忘记加 shell() 函数直接写了系统命令，结果报错说不认识那个命令。其实最简单的方法就是 run 里写 subprocess.run()，但 Snakemake 提供了更方便的 shell() 函数，它会自动处理变量替换和日志。

5. 实战避坑指南

最后汇总几个我亲自踩过的坑：

• 通配符不匹配：output 里的通配符必须在 input 或 params 里对应出现，否则 Snakemake 不知道从哪里取值。我曾经在 output 里写了 {sample}.bam，但 input 里写的是 {sample}.fastq，表面看没问题，但如果同时有多组样本，通配符解析会错乱。
• 目录未创建：前面说了，output 目录需要手动创建。可以在 shell 开头加 mkdir -p $(dirname {output})。
• 输出文件未全部生成：如果一个 rule 声明了多个 output，但只生成了部分文件，Snakemake 会认为这个规则失败。我在做质控流程时经常遇到——fastqc 有时只生成一个结果文件，需要检查是不是命令参数错了。
• shell 变量与环境变量冲突：在 shell 块里使用 $PATH 或者 $HOME 这类系统变量时，Snakemake 的大括号替换和 shell 的美元符号会打架。解决方案：在 shell 块内用 $...$ 转义，或者用 params 传入。

下面这张 SVG 图把 rule 的核心结构串起来了，你可以对照上面的例子理解每个部分的作用。

掌握 rule 的 input、output、shell/run 这三个基本块，你就有了构建流程的砖瓦。剩下的就是拼接和调试，这部分我们后面会深入讲。嗯，今天的先到这，你可以打开编辑器写一个最简单的 rule 试试。