PICRUSt2 — 基于 16S rRNA 的微生物功能预测工具¶

一句话说明¶

PICRUSt2（Phylogenetic Investigation of Communities by Reconstruction of Unobserved States 2）利用 16S 测序数据和参考基因组数据库，预测微生物群落的功能基因丰度（MetaCyc/KO/EC 通路）。

安装与配置¶

# 创建专用环境
conda create -n picrust2 -c bioconda -c conda-forge picrust2 -y
# 当前最新版 v2.6.3（2025年），无需指定 Python 版本
# 等价分步安装：
# conda create -n picrust2 -y
# conda activate picrust2
# conda install -c bioconda -c conda-forge picrust2 -y

conda activate picrust2

# 验证安装
picrust2_pipeline.py --help      # 查看流程脚本帮助
place_seqs.py --help             # 查看序列定位帮助

# 下载参考数据（内置，安装时已包含，无需额外下载）
# 可选：下载扩展参考数据库（更大更全）

核心用法¶

一键运行完整流程（最推荐）¶

# PICRUSt2 完整流程：从 FASTA + 丰度表直接到功能预测
# -s：代表序列 FASTA（QIIME2 去噪后的 ASV/OTU 序列）
# -i：feature 表（ASV/OTU 丰度矩阵，BIOM 格式）
# -o：输出目录
# -p：并行进程数
picrust2_pipeline.py \
    -s repseqs.fasta \            # 代表序列（FASTA 格式）
    -i feature_table.biom \       # ASV/OTU 丰度表（BIOM 格式）
    -o picrust2_out/ \
    -p 16                         # 并行进程数

分步运行（便于排查问题）¶

# 步骤1：序列系统发育定位
# 将 16S 序列放到参考树中
place_seqs.py \
    -s repseqs.fasta \            # 查询序列
    -o placed_seqs.tre \          # 输出放置后的树
    -p 16 \
    --intermediate intermediate_place/ # 保存中间结果

# 步骤2：隐状态预测（Hidden State Prediction）
# 预测每个序列的基因家族拷贝数
hsp.py \
    -i 16S \                      # 预测 16S 拷贝数（标准化用）
    -t placed_seqs.tre \
    -o marker_predicted_and_nsti.tsv \
    -p 16 \
    -n                            # 输出 NSTI（近邻分类单元指数）

hsp.py \
    -i EC \                       # 预测 EC 酶类编号
    -t placed_seqs.tre \
    -o EC_predicted.tsv \
    -p 16

# 步骤3：样本功能丰度加权
metagenome_pipeline.py \
    -i feature_table.biom \       # ASV 丰度表
    -m marker_predicted_and_nsti.tsv \ # 标准化标记
    -f EC_predicted.tsv \         # 功能预测
    -o metagenome_out/            # 输出目录

# 步骤4：通路层次推断
# 从酶类 EC 推断 MetaCyc 通路丰度
pathway_pipeline.py \
    -i metagenome_out/pred_metagenome_unstrat.tsv.gz \
    -o pathways_out/ \
    -p 16

参数详解¶

参数	说明	默认值
`-s`	代表序列 FASTA	必填
`-i`	ASV/OTU 丰度表（BIOM）	必填
`-o`	输出目录	必填
`-p`	并行进程数	1
`--min_reads`	最少 reads 数过滤	1
`--max_nsti`	最大 NSTI 阈值（质控）	2.0
`--stratified`	输出贡献分层结果	关闭
`--coverage`	预测通路覆盖度	关闭
`--hsp_method`	HSP 方法（mp/emp/pic/scp/subtree_average）	mp

输出文件说明¶

ls picrust2_out/
# EC_metagenome_out/
#   pred_metagenome_unstrat.tsv.gz  — EC 酶类预测丰度（非分层，最常用）
#   pred_metagenome_strat.tsv.gz    — 按物种分层的预测丰度
# KO_metagenome_out/
#   pred_metagenome_unstrat.tsv.gz  — KEGG KO 预测丰度
# pathways_out/
#   path_abun_unstrat.tsv.gz        — MetaCyc 通路丰度（最常用于下游分析）
#   path_cov_unstrat.tsv.gz         — 通路覆盖度（完整性评估）
# placed_seqs.tre                   — 序列定位后的系统发育树
# marker_predicted_and_nsti.tsv     — 含 NSTI 质控指标

# NSTI 说明：值越小表示预测越可靠
# 通常 NSTI < 2 认为可靠；肠道样本通常很低（<0.1）

实战案例¶

# 场景：分析 T2D 患者和健康人肠道菌群 16S 数据的功能差异

# 前置：用 QIIME2 完成 16S 去噪，导出 BIOM 表和代表序列
# 假设已有：
# - repseqs.fasta（ASV 代表序列）
# - feature_table.biom（样本 × ASV 丰度表）

# 1. 运行 PICRUSt2 全流程
picrust2_pipeline.py \
    -s repseqs.fasta \
    -i feature_table.biom \
    -o picrust2_results/ \
    -p 16 \
    --max_nsti 2.0                # 过滤 NSTI > 2 的低质量预测

# 2. 解压并查看结果
gunzip picrust2_results/pathways_out/path_abun_unstrat.tsv.gz

# 3. 在 R 中分析通路差异（与 DESeq2 结合）
# 见 359_DESeq2 章节中功能预测数据分析

# 4. 查看 NSTI 质控指标
awk '$NF>0.2' picrust2_results/marker_predicted_and_nsti.tsv | wc -l
# 统计 NSTI > 0.2 的 ASV 数（这些预测精度较低）

常见报错与解决¶

报错信息	原因	解决方法
`BIOM format error`	输入表格式错误	检查 BIOM 文件格式，用 QIIME2 导出
`No sequences placed`	序列与参考数据库差异太大	检查序列是否真为 16S V4 区域
`epa-ng not found`	epa-ng 未安装	`conda install -c bioconda epa-ng`
预测结果全为 0	NSTI 过高，序列被过滤	放宽 `--max_nsti` 或检查序列质量
内存不足	样本量太大	减少并行进程 `-p`，分批处理

速查表¶

# 一键全流程（推荐）
picrust2_pipeline.py -s repseqs.fasta -i table.biom -o out/ -p 16

# 查看通路预测结果
zcat out/pathways_out/path_abun_unstrat.tsv.gz | head -5

# QIIME2 导出（PICRUSt2 的前置步骤）
# qiime tools export --input-path dada2-rep-seqs.qza --output-path exported/
# qiime tools export --input-path dada2-table.qza --output-path exported/
# 转换 BIOM 为 TSV
biom summarize-table -i exported/feature-table.biom  # 查看表摘要