宏基因组分析：物种注释、功能注释与 MAG
宏基因组测序（shotgun metagenomics）把样本中的全部 DNA 随机打断后测序，不经过特定基因的扩增，因此可以同时回答「有哪些微生物」和「它们携带哪些功能基因」。主要有两条分析路线：基于读段的分析，把读段与参考数据库比对，用 Kraken2、MetaPhlAn 得到物种组成，用 HUMAnN 得到基因家族和代谢通路丰度；基于组装的分析，先把读段组装成 contig，再分箱得到宏基因组组装基因组（MAG），可以发现数据库中没有的新物种，并研究其完整的代谢能力。典型数据是粪便、土壤、水体等样本的双端短读长测序。与 16S 相比，宏基因组分辨率可到种甚至菌株，但成本更高，还需处理宿主 DNA 污染。
什么时候用
- 需要种或菌株水平的物种组成，而 16S 只能可靠到属
- 研究短链脂肪酸合成、耐药基因、毒力因子等功能潜力
- 从环境或人体样本中重建新物种的基因组（MAG）
- 同时研究细菌、古菌、病毒和真菌等不同类群
- 与代谢组数据联合，把功能通路与代谢物变化对应起来
不适合的情况
- 宿主 DNA 占绝对多数的组织或血液样本，微生物读段太少
- 样本量大、预算有限且只需属水平组成，16S 更经济
- 想评估基因的实际表达活性，需要宏转录组或宏蛋白组
需要准备的数据
- 双端 FASTQ
- 每个样本的原始测序数据，注明平台和读长；低生物量样本应同时测序阴性对照。
- 宿主参考基因组
- 人或小鼠等宿主基因组的比对索引，用于去除宿主读段；样本来源不同时需更换。
- 参考数据库
- Kraken2 数据库、MetaPhlAn 标记基因库、HUMAnN 所需的核酸与蛋白数据库，全项目版本统一。
- 样本元数据
- 分组、采样时间、个体 ID、测序批次和 DNA 提取批次，用于统计建模和批次检查。
分析步骤
- 1
质控与去宿主
用 fastp 或 Trimmomatic 去除接头和低质量碱基，再用 Bowtie2 比对宿主基因组去除宿主读段，KneadData 可串联这两步。记录每个样本的宿主读段比例，比例极高的样本需单独评估。
- 2
物种组成
用 Kraken2 加 Bracken 估计种水平读段丰度，或用 MetaPhlAn 基于标记基因估计相对丰度。查看未分类读段比例，过高提示数据库覆盖不足或仍有宿主残留。
- 3
功能谱
用 HUMAnN 得到 UniRef 基因家族和 MetaCyc 通路丰度，按需重组为 KO 或 EC，并标准化为 CPM 或相对丰度。检查未比对和未归类部分的比例是否在组间相近。
- 4
组装
单样本或多样本联合，用 MEGAHIT 或 metaSPAdes 组装。评估 contig 数量、N50 和读段回比率，回比率低说明大量读段没有被组装利用，低丰度物种可能缺失。
- 5
分箱与精炼
把读段回比到 contig 得到各样本覆盖度，用 MetaBAT2、MaxBin2、CONCOCT 等工具分箱，再用 DAS Tool 整合多个分箱结果，取各工具的优势。
- 6
MAG 评估注释
用 CheckM 或 CheckM2 评估完整度和污染度，用 GTDB-Tk 做分类，用 dRep 去冗余，再用 Prodigal 预测基因并做功能注释，得到每个 MAG 的代谢潜力。
- 7
专题功能注释
按研究问题选择数据库：耐药基因用 CARD，碳水化合物活性酶用 dbCAN，KEGG 通路可用 eggNOG-mapper 或 KofamScan。注释结果按基因长度和测序深度标准化后再比较。
关键参数与决策
Kraken2 置信度
confidence 参数默认为 0，只要有 k-mer 命中就分类，假阳性物种较多。提高到 0.1 左右（经验值）可减少低丰度假阳性，代价是未分类读段比例上升。
低丰度物种过滤
读段法常给出大量极低丰度物种，其中不少来自数据库误分类。做多样性和差异分析前，常按相对丰度或流行率过滤，阈值需事先确定并写入方法。
单样本与联合组装
单样本组装能保留菌株差异，计算量较小；联合组装提高低丰度物种的覆盖，但可能把不同菌株混成嵌合 contig。样本较多时可按分组或相似度分批联合组装。
MAG 质量标准
常参照 MIMAG 标准：高质量 MAG 要求完整度大于 90%、污染度小于 5%，并含 rRNA 基因和足够的 tRNA；中等质量要求完整度至少 50%、污染度小于 10%。
去冗余阈值
dRep 常以 95% 平均核苷酸一致性（ANI）作为物种水平的去冗余阈值，更高阈值（如 99%）用于菌株水平。阈值决定最终保留多少个代表基因组。
常见错误
- 不去宿主直接注释：宿主读段可能被误分为微生物，原始数据中的人源序列也涉及隐私，公开数据前必须去除。
- 混用不同版本的数据库：物种命名和标记基因随版本变化，跨批次比较会出现与生物学无关的虚假差异。
- 把 Kraken2 读段比例当作细胞丰度：基因组越大的物种贡献读段越多，读段比例与细胞比例不同，不能和 MetaPhlAn 结果互换。
- 只看 MAG 完整度不看污染度：污染度高的 MAG 可能混入多个物种的 contig，据此推断的代谢能力会张冠李戴。
- 测序深度不足却追求低丰度物种 MAG：覆盖度不够时组装高度碎片化，低丰度物种几乎无法成功分箱。
结果怎么看
物种组成堆叠图
展示种水平相对丰度。注意未分类比例在组间是否不同：某组未分类比例明显更高，可能是宿主残留或含较多未知物种，会扭曲已分类物种的相对丰度比较。
分层功能通路表
HUMAnN 输出的通路丰度可按物种分层，显示某通路主要由哪些物种贡献。若组间通路差异几乎全部来自单一物种，结论更接近物种变化，而不是功能重塑。
MAG 质量散点图
横轴为完整度、纵轴为污染度，每个点是一个 MAG。高质量 MAG 集中在完整度高、污染度低的区域；大量 MAG 污染度偏高，说明分箱参数或组装策略需要调整。
MAG 分类与注释表
GTDB-Tk 分类结果和各 MAG 的功能注释。无法归入已知种的 MAG 可能是新物种，可结合与参考基因组的 ANI 判断；中等质量 MAG 缺少某通路，可能只是基因组不完整。
耐药基因丰度表
列出每个样本中各耐药基因的标准化丰度。比较前按每百万读段或单拷贝标记基因校正，组间差异需排除测序深度和宿主读段比例不同造成的影响。
和相近方法怎么选
Kraken2 与 MetaPhlAn
Kraken2 基于 k-mer 对每条读段分类，速度快，可覆盖细菌、病毒和真菌，但需配合 Bracken 重估丰度且假阳性较多；MetaPhlAn 只用分支特异的标记基因，特异性高，丰度更接近细胞比例，但对数据库外物种不敏感。
基于读段与基于组装
读段法依赖参考数据库，速度快、样本间可比性好，适合大样本定量描述；组装法能发现新物种和完整基因簇，但计算量大，低丰度物种难以组装。两者常结合：读段法做定量，组装法做新物种和基因组层面的发现。
宏基因组与宏转录组
宏基因组反映群落拥有哪些基因，即功能潜力；宏转录组测量这些基因是否正在表达。宏转录组需要 RNA 保存和去除 rRNA，实验要求更高，常与宏基因组配对，计算表达量相对于基因拷贝数的变化。
常见问题
Kraken2 和 MetaPhlAn 结果为什么差别很大？
两者原理不同。Kraken2 统计读段被分到哪个物种，基因组越大的物种读段越多；MetaPhlAn 统计物种特异标记基因的覆盖，更接近细胞数比例。两者的数据库收录范围也不同。同一项目选定一种作为主分析，不要混用两种结果。
宏基因组每个样本测多少数据合适？
取决于研究目标。只做物种和功能谱描述，所需深度相对较低；要组装低丰度物种的 MAG 或分析菌株变异，需要明显更高的深度。可先在部分样本上评估物种累积曲线和组装质量，再决定其余样本的测序量。
HUMAnN 的结果怎么标准化？
HUMAnN 默认输出 RPK，即按基因长度校正后的读段数，不同测序深度的样本之间不能直接比较。用 humann_renorm_table 转为 CPM 或相对丰度后再做统计，并说明使用的是总表还是按物种分层的结果。
MAG 和参考基因组有什么区别？
MAG 是从混合样本中组装和分箱得到的基因组草图，通常由多个 contig 组成，可能不完整或带有少量污染，代表的是样本中一个群体的共识序列。分离培养后测序的参考基因组通常更完整，但很多微生物难以培养，MAG 正好弥补这一点。
在 Codar 里怎么做
在 Codar 中提供测序数据、宿主类型和研究问题后，会先确认数据库选择和是否需要组装，再进行宏基因组物种分类与丰度估计，得到种水平组成；进行宏基因组功能组成分析，得到基因家族和通路丰度；需要时进行宏基因组组装与 MAG 重建，评估完整度和污染度并完成分类。交付物种组成图、功能丰度表、MAG 质量表、分析报告和可复现脚本。
我有 40 份人粪便宏基因组双端测序数据，想去除宿主后用 MetaPhlAn 和 HUMAnN 比较两组的物种与通路，并尝试重建高质量 MAG。