3. 常用生物信息学流程：QIIME2、MOTHUR等

宏基因组分析，说白了就是一场从海量读段里“淘金”的活儿。工具选得好，事半功倍；选不好，光折腾格式就能让你怀疑人生。我自己刚入行那会儿，QIIME和MOTHUR都试过，踩了不少坑，今天把心得都倒给你们。

这两个流程覆盖了从原始数据质控、OTU/ASV聚类、物种注释到多样性分析的全链路。你不用全都精通，但至少得知道它们各自擅长什么、瓶颈在哪。

3.1 两大主流流程的定位

流程	核心哲学	适合场景	我常用的版本
QIIME 2	插件化、可追溯、强依赖Python生态系统	批量项目、需要Conda管理的复杂环境	2024.10 稳定版
MOTHUR	单体二进制、C++底层、命令链式调用	服务器无Root权限、快速尝鲜	1.48.0

我个人习惯：如果团队里有生物信息学基础参差不齐的人，我会强推QIIME 2，因为它的artifact系统能保证你每一步的“输入-输出”都不乱。但如果只是自己快速跑个16S看看群落结构，MOTHUR的批处理模式反而更省心。

3.2 QIIME 2 的核心工作流

QIIME 2 把分析拆成了一个个独立的“插件”（qiime2.plugins）。每个插件处理一个步骤，比如去噪、分类、建树。它的数据格式不是普通的FASTQ，而是.qza（QIIME 2 Artifact）和.qzv（可视化结果）。

为什么要这么设计？

你想，普通流程里你跑完一步，生成一个FASTA，下一步还得自己写脚本判断版本。QIIME 2 的.qza文件里不光有数据，还嵌入了处理元数据、参数、引用。说白了，这是“数据结构化”策略。我在一个跨国合作项目里用过这个特性——对方只需要解压.qzv就能直接看交互式图表，根本不用装任何生物信息学软件。

3.2.1 安装与环境配置

我个人建议用Miniforge3 + Conda环境，别用系统自带的pip。QIIME 2的依赖非常复杂，搞不好就出现插件冲突。

# 创建一个干净的环境，别混装
conda create -n qiime2-2024.10 \
    --override-channels \
    --strict-channel-priority \
    -c conda-forge \
    -c bioconda \
    -c qiime2/label/r2024.10 \
    qiime2

conda activate qiime2-2024.10

嗯，这里要注意：--strict-channel-priority 这行不能省。我见过有人因为没加这个参数，装出来的Cutadapt版本跟QIIME 2不兼容，报错报了整整一个下午。

3.2.2 数据导入与去噪（DADA2）

导入这一步，很多人直接拿原始FASTQ就开干，结果运行时报错“不支持的格式类型”。别急，先确认你的测序数据是单端还是双端，以及质量分数编码是Phred33还是Phred64。

# 导入双端数据（最常见的格式）
qiime tools import \
  --type 'SampleData[PairedEndSequencesWithQuality]' \
  --input-path /path/to/your/seqs \
  --input-format CasavaOneEightSingleLanePerSampleDirFmt \
  --output-path demux.qza

然后就是核心的去噪和ASV聚类。我建议用DADA2而不是Deblur，因为DADA2对低质量区的纠错能力更强，适合真实的环境样本。

qiime dada2 denoise-paired \
  --i-demultiplexed-seqs demux.qza \
  --p-trunc-len-f 250 \
  --p-trunc-len-r 200 \
  --p-n-threads 0 \
  --o-table table.qza \
  --o-representative-sequences rep-seqs.qza \
  --o-denoising-stats stats.qza

3.2.3 物种注释与多样性分析

分类器我习惯用预训练的Naive Bayes，别自己从头训练——除非你做的样本特别偏，比如深海热泉或极端酸性矿排水。

# 使用Silva 138.1数据库的预训练分类器
qiime feature-classifier classify-sklearn \
  --i-classifier silva-138.1-20240916-classifier.qza \
  --i-reads rep-seqs.qza \
  --o-classification taxonomy.qza

多样性分析就更机械化了，直接套用官方的core-metrics-phylogenetic插件。

qiime diversity core-metrics-phylogenetic \
  --i-table table.qza \
  --i-phylogeny rooted-tree.qza \
  --m-metadata-file sample-metadata.tsv \
  --p-sampling-depth 10000 \
  --output-dir core-metrics

3.3 MOTHUR 的“硬核”模式

MOTHUR 跟QIIME 2完全不同。它就一个可执行文件，用命令行参数把整个分析链串起来。它的优势是：对系统依赖极少，哪怕你用的是一台只有g++运行时库的旧服务器，也能跑。

但有个坑：MOTHUR的自定义脚本是它自己的一套内部指令（batch模式），不支持Python或R的灵活扩展。一旦遇到复杂统计模型，你只能算完结果导出再另行处理。

3.3.1 经典分析脚本骨架

我一般把MOTHUR的指令写在一个文本文件（比如 mothur.batch）里，然后用 mothur mothur.batch 跑。这样方便修改参数，也方便别人复现。

# mothur.batch 示例
set.dir(input=./data, output=./results)

make.contigs(file=stability.files, processors=8)

screen.seqs(fasta=current, group=current, maxambig=0, maxhomop=8)

unique.seqs(fasta=current)

count.seqs(name=current, group=current)

align.seqs(fasta=current, reference=silva.v4.fasta)

screen.seqs(fasta=current, count=current, start=1968, end=11554)

filter.seqs(fasta=current, vertical=T, trump=.)

pre.cluster(fasta=current, count=current, diffs=2)

chimera.vsearch(fasta=current, count=current, dereplicate=T)

remove.seqs(fasta=current, accnos=current)

classify.seqs(fasta=current, count=current, reference=trainset16_022016.rdp.fasta, taxonomy=trainset16_022016.rdp.tax, cutoff=80)

dist.seqs(fasta=current, cutoff=0.03)

cluster(column=current, count=current)

make.shared(list=current, count=current, label=0.03)

classify.otu(list=current, count=current, taxonomy=current, label=0.03)

summary.tax(list=current, count=current, taxonomy=current, label=0.03)

quit()

3.3.2 输出结果解读

MOTHUR的关键输出有：

• .shared 文件：OTU表，可以直接读入R的vegan或phyloseq包。
• .cons.taxonomy 文件：每个OTU的共识分类。
• .alpha/diversity 目录：Shannon、Chao1等多样性指数。

我个人觉得MOTHUR最大的优点是 表格干净。没有一堆.qza的神秘包装，文本文件打开就能看。对于需要跟实验组同事频繁对齐数据格式的场景，特别友好。

3.4 流程选择决策图

3.5 避坑指南与个人经验集

我曾经遇到过一个很典型的问题：用QIIME 2的DADA2跑完，feature-table里居然有超过一半的序列被过滤掉了。后来排查发现，是--p-trunc-len-f设得太大，导致大部分读段在修剪后长度不足被丢弃。记住：宁可在参数上宽松一点（比如设220而不是250），留出更多序列用于下游分析，也不要为了追求短读长的完美质量而丢失太多真实信号。

另一个坑是关于MOTHUR的数据库路径：如果你在MOTHUR里用set.dir(reference=...)设了相对路径，那么MOTHUR生成的中间文件里所有引用都会跟着相对路径跑。一旦你移动了batch文件或者数据库位置，所有输出都会报“找不到文件”。我个人习惯：要么用绝对路径，要么在batch脚本开头显式cd到固定工作目录。