用 Snakemake 构建 GLORI-seq m6A 检测流程

表观转录组学——研究 RNA 上的化学修饰——是分子生物学中最令人兴奋的前沿之一。在已知的数百种 RNA 修饰中，N6-甲基腺苷（m6A）是信使 RNA 上最丰富的内部修饰。它影响着从 mRNA 稳定性、剪接到翻译效率的方方面面，而其调控异常与癌症、神经退行性疾病和发育缺陷都有关系。

但要在单碱基分辨率下检测 m6A 历来很困难。早期的大多数技术（如 MeRIP-seq）只能将 m6A 定位到约 100-200 个核苷酸的区域内。GLORI 的出现改变了这一切——这是一种化学生物学方法，能够以单碱基分辨率定量检测 m6A。在读博期间，我需要一个可靠、可复现的流程来处理 GLORI-seq 数据，于是我构建了 snakemake-epitranscriptome。

为什么选择 Snakemake？

如果你曾经用一堆零散的 shell 脚本来处理二代测序数据，你一定体会过那种痛苦：脚本换台机器就跑不通，中间文件被覆盖，复现一篇已发表的分析几乎不可能。Snakemake 通过将有向无环图（DAG）来表达分析流程来解决这个问题——每条规则定义了输入、输出以及如何从前者产生后者。Snakemake 自动判断哪些步骤需要运行、处理并行化，并确保可复现性。

GLORI 方法简述

GLORI 的原理很巧妙：化学处理会将未修饰的腺苷（A）转化为肌苷（I），在测序时被读作鸟苷（G）。而 m6A 修饰的腺苷则被保护，不会发生这种转化。因此，通过比较处理组和对照组在每个位置的 A-to-G 转化率，就可以定量计算单碱基分辨率的 m6A 比例。

这意味着生物信息学流程需要做一些不寻常的事：构建一个修饰过的参考基因组，将其中的 A 全部替换为 G，模拟化学处理在序列层面的效果。

流程概述

整个工作流分为三个主要阶段：

1. 数据准备

原始 FASTQ 文件经过质量控制、接头去除和去重。由于 GLORI-seq 常使用唯一分子标识符（UMI）来区分真正的生物重复和 PCR 扩增假象，流程包含了完整的 UMI 处理和提取步骤，避免过高估计读数。

2. 构建索引

这是最有趣的部分。流程不是将 reads 比对到标准参考基因组，而是构建 AG 转换参考——将基因组中的 A 替换为 G 的修改版本。这模拟了 GLORI 处理过程中的化学转化，使比对器能够正确映射转换后的 reads。处理组和对照组分别构建索引。

3. m6A 检测

将 reads 比对到相应的转换参考后，流程使用 GLORI-tools 进行 m6A 位点检测。在每个腺苷位置，比较处理组和对照组之间的转化率，进行统计检验，输出高置信度的 m6A 位点及其修饰水平。

发表在 Snakemake Workflow Catalog

我的目标之一是让这个流程方便他人使用。很高兴 snakemake-epitranscriptome 已经被收录在 Snakemake Workflow Catalog 中，这意味着任何人都可以找到它、克隆它，只需最少的配置就能在自己的数据上运行。被收录需要遵循 Snakemake 在流程结构、文档和测试方面的最佳实践——这本身就是一次有价值的锻炼。

可复现性为什么重要

在表观转录组学这样一个年轻的领域，可复现性就是一切。不同实验室在修剪、比对和统计阈值上使用的参数略有不同。没有标准化的流程，跨研究的比较就是鸡同鸭讲。通过将整个分析打包成带有版本锁定依赖的 Snakemake 流程，我希望让任何人都能自信地复现和扩展分析。

学到了什么

构建这个流程让我在多个方面收获颇丰：

GLORI 方法原理： 理解 A-to-G 转化背后的化学原理，以及它如何影响生物信息学决策（比如构建转换参考基因组），让我更加深刻地体会到实验方法如何约束计算分析。
修饰参考基因组： 构建 AG 转换参考不是大多数标准流程会做的事。它需要仔细处理基因组索引、链信息和坐标映射回原始参考等问题。
Snakemake 最佳实践： 为 Snakemake Catalog 组织流程结构，让我学会了规则模块化、配置 schema、每条规则独立的 conda 环境管理，以及编写有意义的集成测试。
发布到流程目录： 准备公开分发的流程——文档、示例数据、持续集成测试——本身就是一次有价值的科学软件工程实践。

如果你正在处理 GLORI-seq 数据，或者只是对 m6A 检测感兴趣，欢迎访问 snakemake-epitranscriptome。问题和 Pull Request 随时欢迎。