4. MetaPhlAn 原理简介

聊完 Kraken2，咱们来看看另一个主流工具——MetaPhlAn。说实话，我最初接触宏基因组分析时，第一个用的就是它。那时候实验室老机器跑 Kraken2 特别慢，MetaPhlAn 反而挺快。后来才明白，它俩走的完全是两条技术路线。

Kraken2 拿整条序列的 k-mer 当武器，MetaPhlAn 则更“挑食”——只盯着每个物种里那些独一无二的标记基因 (marker genes)。说白了，它用最少的序列信息，干最精确的活。

核心思想：标记基因

什么是标记基因？就是基因组里一段保守又特异的序列。比如细菌的 16S rRNA 基因，但 MetaPhlAn 不止用这一种。

我记得 2019 年做过一个粪便样本的对比实验。同时跑了 Kraken2 和 MetaPhlAn，结果 Kraken2 检出了 300 多种“菌”，但很多丰度只有 0.01%。MetaPhlAn 只报出 80 多种，但后来培养验证，那 80 种基本都对得上。这就是 marker 策略的威力——宁可漏过，也不乱报。

MetaPhlAn 的开发者们干了一件极繁的事：从 NCBI 等公共库里，给每个已测序的细菌、古菌、真菌、病毒物种，手动筛选出上百个专有的标记基因。这些基因必须满足两个条件：

  • 物种内高度保守——不同菌株间差异极小。
  • 物种间绝对特异——在其他物种的基因组里找不到同源片段。

目前最新版数据库覆盖了超过 2 万个物种。你想想看，这个筛选工作量有多大。我尝试自己写过一个小脚本模拟这个过程，跑了三天三夜没跑完……

数据库构建

MetaPhlAn 的数据库不是让用户自己建的，官方直接提供打包好的数据库文件。但你得了解它长什么样，否则后续遇到问题都不知道去哪排查。

数据库本质上是一个 多序列 FASTA 文件，每条序列就是一段标记基因，头部标注了对应的物种 ID（比如 t__SGB1234）。此外还有一个 Bowtie2 索引，因为 MetaPhlAn 依赖 Bowtie2 做比对。

常见槽点：很多新手下载数据库时只拿了 FASTA，忘了建索引。直接跑 MetaPhlAn 就会报错：“No Bowtie2 index found”。解决办法：先 bowtie2-build 一把。

MetaPhlAn 的数据库下载命令如下（注意版本号，我用的是 v4.0 系列）：

metaphlan --install --bowtie2db /path/to/database
# 或者手动下载：
wget https://zenodo.org/record/5070819/files/metaphlan4_database_v0.1.tar.gz
tar -xzf metaphlan4_database_v0.1.tar.gz -C /path/to/database
bowtie2-build /path/to/database/metaphlan4_database_v0.1.fasta /path/to/database/metaphlan4_database_v0.1

比对的流程

MetaPhlAn 的分析流程其实很简洁，就三步：

  1. 比对：Bowtie2 把 reads 映射到标记基因数据库，生成 SAM 文件。
  2. 解析：MetaPhlAn 读 SAM，统计每个标记基因被多少 reads 覆盖。
  3. 丰度估计：用贝叶斯模型或最大期望算法，反推出每个物种的相对丰度。

为什么需要贝叶斯模型？因为多个标记基因的覆盖度可能不一致，偶尔还会出现跨物种的 cross-mapping。MetaPhlAn 内部有一套纠偏机制，碰到可疑的 read 会直接丢弃。我记得曾经处理一个口腔样本，有位点上的 marker 基因覆盖高达 500×，但其他 marker 只有 10×。算法自动判了低信心，最终给的丰度接近中间值——这个机制救了我一命。

运行命令非常简单：

metaphlan sample_R1.fastq,sample_R2.fastq --input_type fastq \
  --bowtie2db /path/to/database \
  --nproc 8 > profile.txt

结果文件 profile.txt 是一个制表符分隔的表，格式类似：

#clade_name      NCBI_tax_id      relative_abundance
k__Bacteria      2                95.12345
k__Bacteria|p__Firmicutes  1239    45.67890
...

每行的缩进层次就代表了分类等级（界、门、纲、目、科、属、种）。直接用 grep 或者 cut 就能提取你想要的级别。

Kraken2 vs MetaPhlAn：一张表看懂

对比维度Kraken2MetaPhlAn
比对数据库全基因组 k-mer标记基因集
假阳性风险高（容易匹配到保守区域）低（只认特有基因）
运行速度快（内存大则极快）中等（Bowtie2 比对耗时）
物种覆盖度广（含很多低质量物种）精（只含高质量参考基因组）
丰度准确性中等（依赖 k-mer 计数）高（贝叶斯纠偏）

MetaPhlAn 核心流程 SVG

MetaPhlAn 物种分类及丰度估计流程 输入 Reads (FASTQ) Bowtie2 比对 （标记基因数据库） 统计每个 marker 的 read 计数 贝叶斯/EM 丰度估计 输出丰度表格 关键优势 • 只比对特有基因 • 避免假阳性 • 丰度估计更准 • 数据库小、易维护
曾经踩过的坑：MetaPhlAn 默认使用 bowtie2 的 --very-sensitive 模式，这对短读长（比如 100bp）非常友好。但如果你是长读长数据（比如 PacBio CCS），建议加上 --local 参数，否则很多标记基因会被判定为不匹配，导致低估物种丰度。

核心参数解读

MetaPhlAn 的命令行参数很多，但日常分析我只关注这几个：

  • -t：指定输出分类等级，比如 -t rel_ab_w_read_stats 可以同时输出 reads 计数。
  • --unknown_estimation：开启后能估计未分类 reads 的比例，对于判断样本污染很有用。
  • --ignore_eukaryotes：如果只看细菌，加上它能让速度提升 30% 左右。

我习惯在正式分析前先跑一个小样本（比如随机抽 1 万条 reads），快速验证数据库路径和比对参数。用这个命令：

seqtk sample big_file.fastq 10000 | metaphlan - -o test_profile.txt

如果 5 分钟内能跑完，说明环境没问题。如果卡住不动，八成是 Bowtie2 索引路径写错了。

一句话总结

MetaPhlAn 是用最少的信息量，做最可靠的分类。它不像 Kraken2 那样“全基因组无差别射击”，而是像个狙击手——只打标记基因这个靶心。代价是覆盖度不如 Kraken2 广，但换来的是更高的精度和更低的假阳性。在临床样本、环境样本中，我一般优先用 MetaPhlAn 做第一轮分析，再用 Kraken2 做补盲。

← 上一章 📖 返回目录 下一章 →