文章目录

• conda安装qiime2

• 导入数据

  • 制作Manifest和Metadata表

  • Import数据

  • 查看原始数据质量

• DADA2去噪、去嵌合体和生成OTU

• 构建进化树

• 绘制稀释曲线

• 计算物种多样性

现如今，测序已经将近实现“千元基因组“，近百元就能测一个扩增子样本，产出数据能广泛揭示不同环境样本的物种组成和丰度等众多信息。扩增子测序就像下饭菜，大米饭里拌几勺吃起来就美味而可口，随手“一片”SCI，Qiime2扩增子处理流程确定不了解一下？

conda安装qiime2

本教程介绍在Linux上的qiime2 命令行版q2cli(qiime2 command line interface)，最新版为qiime2-2020.2。安装方式采用最便捷的conda安装方法：

# 下载qiime2的conda安装配置文件

导入数据

理论上qiime2可以从任何一步分析的数据导入进行重新分析，这里以最为常见的公司给的拆分过barcode的pair-end带质量信息的fastq序列示例。(目前，qiime仅能拆分单端数据，对于双端数据的拆分需要自己撸代码或者借助其他工具，暂且不表。)

制作Manifest和Metadata表

在导入数据之前需要制作两个文件：
manifest表用于显示样本名称和对应的双端序列路径
metadata表用于记录样本对应的元数据信息(例如分组、年龄、性别等等)
manifest示例manifest.tsv：
第一列样本名称，第二列为forward序列路径，第三列为reverse序列路径
metadata示例metadata.tsv：
第二列重复样本名称是便于后面显示每个样本的alpha多样性曲线。

• 注意，这两个表都是制表符分割，制作好后可以先在google浏览器插件Keemei上验证文件格式是否标准：
  Add ones-Keemei-Validate qiime2 metadata file
  验证通过示例：

Import数据

qiime2要求规范数据格式，主要是.qza的数据文件和.qzv的可视化文件，.qza的文件可以通过unzip(推荐)和export方式解压查看，.qzv的文件可通过View qiime2查看。

import的目的就是生成一个标准的包含原始测序数据的.qza文件(官网称之为Artifact)

# import

一旦导入数据生成.qza的文件，你就可以开始qiime2流程的各种分析，首先查看 原始数据测序质量。

查看原始数据质量

qiime demux summarize \    --i-data paired-end-demux.qza \ #上一步结果    --o-visualization demux-raw-summary.qzv

在Overview可以看到所有样本reads数的基本统计信息：最小、最大值，中位数、平均数。本次分析样本reads数多数在77551附近。
在Interactivate quality plot可以看到双端序列的碱基质量boxplot图，黑框代表50%的数据情况，在3端末尾质量有所下降，后续可视情况去除。

DADA2去噪、去嵌合体和生成OTU

qiime dada2 denoise-paired \    --i-demultiplexed-seqs paired-end-demux.qza \    --p-trunc-len-f 0 \ #截掉3端forward序列X位置后的序列    --p-trunc-len-r 0 \ #截掉3端reverse序列X位置后的序列    --p-n-threads 30 \ #线程    --o-table table.qza \ #FeatureTable[Frequency]    --o-representative-sequences rep-seqs.qza \ #FeatureData[Sequence],每一个feature对应的一个joined的双端代表序列    --o-denoising-stats denoising-stats.qza \#嵌合体等统计信息

# 可视化统计信息qiime metadata tabulate \    --m-input-file denoising-stats.qza \    --o-visualization denoising-stats.qzv

现在最重要的OTU特征表和代表序列已经生成，可以通过unzip解压到指定目录unzip rep-seqs.qza -d OTU/ && unzip rep-seqs.qza -d OTU/，解压后feature-table.biom和dna-sequences.fasta即是特征表和代表序列。
dna-sequences.fasta:eature-table.biom需要经过biom转换成文本查看：biom convert -i feature-table.biom -o feature-table.tsv --to-tsv
feature-table.ts,第一行为注释，第二行为样本名称，第三行起为OTU的丰度，每一个OTU ID都可以在代表性序列中找到。
打开denoising-stats.qzv可看到不同步骤处理后剩余的reads，并且可以排序，这里看到的最少的非嵌合体序列也在35000以上，一般满足分析了。

构建进化树

接下来利用代表序列构建进化树，用于计算物种多样性

#使用mafft比对的fasttree构建进化树

绘制稀释曲线

为了检查样本测序量是否足够反映真实的物种情况,需要抽样检测：原理是按照一定的梯度抽取样本reads，检测OTU数量，当OTU数量随着梯度增加不再增加时，则说明测序深度足够。

qiime diversity alpha-rarefaction \    --i-table table.qza \    --i-phylogeny rooted-tree.qza \#可选，有的话可以在计算observed_otus, shannon之外，计算faith_pd    --p-max-depth 30000 \#最大采样深度，低于该阈值的样本将被丢弃    --m-metadata-file metadata.tsv \#结果可以按照元数据信息分别显示    --o-visualization alpha-rarefaction.qzv

可选择不同的评估方法(observed_otus, shannon和faith_pd)或者查看在不同分组之间的情况。要注意的是系统变量LANG需要设置为en_US.utf-8，否则容易出现字符编码错误,报错：
需要export LANG=en_US.utf-8，或者写到环境变量

计算物种多样性

基于feature table，通过生成一系列的系统发育和非系统发育多样性度量对物种数量和丰度二维信息进行描述，从而评估物种多样性。

qiime diversity core-metrics-phylogenetic \    --i-phylogeny rooted-tree.qza \    --i-table table.qza \    --p-sampling-depth 10000 \#为了比较序列深度不均匀的样本，样本中的最小序列数可用作采样深度，也可以稍微降低一点，可根据feature表的统计结果中Frequency per sample中minimum frequency来确定；需要注意的是，低于这个阈值的样本将会被丢弃。    --m-metadata-file ${metadata} \    --output-dir core-metrics-results

core-metrics-results目录同时包含qza和qzv数据：

• α多样性：
  香农(Shannon’s)多样性指数(群落丰富度的定量度量，即包括丰富度richness和均匀度evenness两个层面)
  Observed OTUs(群落丰富度的定性度量，只包括丰富度)
  Faith’s系统发育多样性(包含特征之间的系统发育关系的群落丰富度的定性度量)
  均匀度(或 Pielou’s均匀度；群落均匀度的度量)

• β多样性
  Jaccard距离(群落差异的定性度量，即只考虑种类，不考虑丰度)
  Bray-Curtis距离(群落差异的定量度量)
  非加权UniFrac距离(包含特征之间的系统发育关系的群落差异定性度量)
  加权UniFrac距离(包含特征之间的系统发育关系的群落差异定量度量)

同样打开qzv文件可查看结果，以Bray_curties_emperor.qzv为例展示基于Bray_curties距离的PCoA：
最为常用的操作包括：
依据分组信息显示不同的颜色和形状(第二个特征)，将背景设置成白色，保存为矢量图等，设置好以后如下图：

如果觉得qiime在线绘图效果不满意，可以解压相应qza文件获得画图数据，再通过R自行绘制。

至此，我们已经得到特征表、稀释曲线和物种多样性的结果，后续将介绍物种注释和qiime2常用操作，欢迎围观。