QIIME 2用户文档.4人体微生物组分析Moving Pictures

前情提要

QIIME 2用户文档. 4人体各部位微生物组

“Moving Pictures” tutorial

https://docs.qiime2.org/2019.7/tutorials/moving-pictures/

本节1.6万字，14张图。阅读时间大约40分钟。

在本教程中，你将使用QIIME 2在五个时间点对来自两个人四个身体部位的微生物组样本进行分析，第一个时间点紧接着是抗生素的使用。基于这些样本的研究文章《Moving pictures of the human microbiome》在2011年发表于Genome Biology。本教程中使用的数据基于Illumina HiSeq产出，使用地球微生物组计划扩增16S rRNA基因高变区4（V4）测序的方法。

对于熟悉QIIME 1的用户，本数据也出现在QIIME的教程中。

在开始本教程前，我们需要进入工作环境创建新目录并进入

本节视频视频教程

https://v.qq.com/x/page/w0918ebti6m.html

a文件准备和样本拆分

https://v.qq.com/x/page/l0918vwb1no.html

b结果查看、质控方法dada2/deblur并生成特征表

https://v.qq.com/x/page/c09194lgqb5.html

c进化树构建，多样性分析统计和可视化，物种注释和柱状图展示，差异比较

查看更多视频和相应专辑，访问下方链接至作者个人频道，持续更新ing

http://v.qq.com/vplus/22b577627f014f0ca25e9827b38c171e

启动QIIME2运行环境

对于上文提到了两种常用安装方法，我们每次在分析数据前，需要打开工作环境，根据情况选择对应的打开方式。

比如我的工作目录为，这是与Github中同步的目录，方便同行下载测试数据。用户可以随便定义你的项目工作目录，如把qiime2学习放在目录中。

我们在每次分析开始前，必须先进入工作目录，除非你是一个把什么东西都放在桌面上还很工作更有效率的人。

样本元数据

Sample metadata

在开始分析之前，我们需要阅读样本元数据，以熟悉本研究中使用的样本信息。示例元数据作为Google 表格提供。你可以通过选择，以制表符分隔的文本格式下载该文件。或者，以下命令将作为制表符分隔的文本下载示例元数据，并将其保存在文件。这个文件在本教程中一直被用到。

提示：Keemei是一个用于验证示例元数据的Google Sheets插件。在开始任何分析之前，样本元数据的验证非常重要。尝试按照Keemei网站上的说明安装Keemei，然后验证上面链接的示例元数据电子表格。该电子表格还包括一个带有一些无效数据的表格，以便使用Keemei进行测试。

提示：要了解关于元数据的更多信息，包括如何格式化元数据以便与QIIME 2一起使用，请参阅元数据教程。

下载和导入数据

Obtaining and importing data

下载在本次分析中使用的序列。在本教程中，我们将处理完整的序列数据的一小部分，以便命令能够快速运行(减少等待时间)。

创建子目录并下载实验测序数据

用于输入到QIIME 2的所有数据都以QIIME 2对象的形式出现，其中包含有关数据类型和数据源的信息。因此，我们需要做的第一件事是将这些序列数据文件导入到QIIME 2对象中。

这个QIIME 2对象的语义类型是。 QIIME 对象是包含多样本混合的序列文件，这意味着序列尚未分配给样本（因此包括和文件，其包含与中的每个序列相关联的条形码）。要导入其他格式的序列数据，请参阅导入数据教程。

导入数据：生成qiime2要求的对象格式。time统计计算时间。

输出对象:

: 查看 | 下载

提示：

上面的和由文档中的命令创建的QIIME 2对象和可视化链接。例如，上面的命令创建了单个文件，上面链接了相应的预计算文件（输出结果）。你可以查看预计算的QIIME 2对象和可视化而不需要安装额外的软件（例如，QIIME 2）。

QIIME 1用户：

在QIIME 1中，我们一般建议通过QIIME执行样本拆分（例如，使用或），因为这个步骤还执行序列的质量控制。现在我们将样本拆分和质量控制步骤分开，因此你可以使用混合多样本序列（如我们在此所做的）或拆分后的序列开始QIIME 2分析。

拆分样品

Demultiplexing sequences

为了混合序列进行样本拆分，我们需要知道哪个条形码序列与每个样本相关联。此信息包含在样品元数据文件中。你可以运行以下命令来对序列进行样本拆分（命令指的是这些序列是根据地球微生物组计划标准方法添加的条形码，并且是单端序列）。QIIME 2对象包含样本拆分后的序列。第二个输出文件 () 包括Golay标签错误校正的详细，在本教程中不作讨论 (你可以使用 查看该结果)。

用时1m3.844s

输出结果：

: 查看 | 下载

: 查看 | 下载

在样本拆分之后，生成拆分结果的统计信息非常重要。这允许我们确定每个样本获得多少序列，并且还可以获得序列数据中每个位置处序列质量分布的摘要。

结果统计

输出可视化结果： 查看 | 下载

图1. 样本拆分结果统计结果——样本数据量可视化图表。

主要分为三部分：上部为摘要；中部为样本不同数据量分布频率柱状图，可下载PDF，下部为每个样本的测序量。上方面板还可切换至交互式质量图页面。如下图2。

图2. 交互式质量图查看页面。

同样为三部分：上部为每个位置碱基的质量分布交互式箱线图，鼠标悬停在上面，即可在下面(中部)文字和表格中显示鼠标所在位置碱基质量的详细信息；下部为拆分样本的长度摘要(一般等长测序无差别)。

注：

所有QIIME 2可视化对象（即使用参数指定的文件）将生成一个文件。你可以使用查看这些文件。我们提供了用于查看可视化的第一个命令，但是对于本教程的其余部分，我们将告诉你在运行可视化程序之后查看结果可视化，这意味着你应该在生成的文件上运。

这条命令的显示需要图形界面的支持，如在有图型界面的Linux上，但仅使用SSH登陆方式无法显示图形。

推荐使用 https://view.qiime2.org 网址显示结果

目前命令行方式想要查看结果可能很多使用服务器人员无法实现 (即依赖服务器安装了桌面，本地依赖XShell+XManager或其它ssh终端和图形界面软件)

本地查看可解压，目录中的data目录包括详细的图表文件，主要关注 pdf 和 html 文件，目录结构如下。

qzv文件解压后文件详细，可直接访问打开结果报告式网页。里面的重要结果，全部可以通过此网页进行索引。

序列质控和生成特征表

Sequence quality control and feature table construction

QIIME 2插件多种质量控制方法可选，包括DADA2、Deblur和基于基本质量分数的过滤。在本教程中，我们使用DADA2和Deblur两种方法分别介绍这个步骤。这些步骤是可互相替换的，因此你可以使用自己喜欢的方法。这两种方法的结果将是一个QIIME 2特征表和一个代表性序列对象，对象包含数据集中每个样本中每个唯一序列的计数（频率），对象将中的特征ID与序列对应。

译者注：此步主要有DADA2和Deblur两种方法可选，推荐使用DADA2，2016年发表在Nature Method上，在阴道菌群研究中比OTU聚类结果看到更多细节，详见《扩增子分析还聚OTU就真OUT了》；相较USEARCH的UPARSE算法，目前DADA2方法仅去噪去嵌合，不再按相似度聚类，结果与真实物种的序列更接近。

注意：本节中此次存在两种可选方法时，你将创建具有特定方法名称的对象（例如，使用dada2去噪生成的特性表将被称为）。在创建这些对象之后，你将把两个选项之一的对象重命名为更通用的文件名（例如，）。为对象创建特定名称，然后对其进行重命名的过程仅允许你选择在本步骤中使用的两个选项中之一完成教程，而不必再次关注该选项。需要注意的是，在这个步骤或QIIME 2中的任何步骤中，你给对象或可视化的文件命名并不重要。

QIIME1 用户注意：

QIIME 2对象等价于QIIME 1 OTU或BIOM表，QIIME 2对象等价于QIIME 1代表序列文件。由于DADA2和Deblur产生的“OTU”是通过对唯一序列进行分组而创建的，因此这些OTU相当于来自QIIME 1的100%相似度的OTU，通常称为序列变体。在QIIME 2中，这些OTU比QIIME 1默认的97%相似度聚类的OTU具有更高的分辨率，并且它们具有更高的质量，因为这些质量控制步骤比QIIME 1中实现更好。因此，与QIIME 1相比，可以对样本的多样性和分类组成进行更准确的估计。

方法1. DADA2

Option 1: DADA2

DADA2是用于检测和校正（如果有可能的话）Illumina扩增序列数据的工作流程。正如在插件中实现的，这个质量控制过程将过滤掉在测序数据中鉴定的任何序列（通常存在于标记基因Illumina测序数据中，用于提高扩增子测序质量），并同时过滤嵌合序列。

方法需要两个用于质量过滤的参数：，它去除每个序列的前m个碱基(如引物、标签序列barcode)；，它在位置n截断每个序列。这允许用户去除序列的低质量区域、引物或标签序列等。为了确定要为这两个参数传递什么值，你应该查看上面由生成的文件中的交互质量图选项卡。

读者思考时间：基于上图对拆分样品的统计结果，如何设置和的参数值。

—p-trim-left 截取左端低质量序列，我们看上图中箱线图，左端质量都很高，无低质量区，设置为0；

—p-trunc-len 序列截取长度，也是为了去除右端低质量序列，我们看到大于120以后，质量下降极大，甚至中位数都下降至20以下，需要全部去除，综合考虑决定设置为120。

单端序列去噪, 输入样本拆分后结果；去除左端 0 bp (—p-trim-left，有时用于切除低质量序列、barocde或引物)，序列切成 120 bp 长(—p-trunc-len)；生成代表序列、特征表和去噪过程统计。

下面的步骤计算量较大，有34个样本，26万条序列，计算大约消耗10分钟。

生成三个输出文件：

: dada2计算统计结果。查看 | 下载

: 特征表。查看 | 下载

: 代表序列。 查看 | 下载

对特征表统计进行进行可视化

输出样本统计表：，查看 | 下载

内容为每个样本，输入、过滤、去噪和非嵌合的统计，并支持按列排序，检索和功能，用于样本异常筛选，特征表抽平标准化非常有用。

表格前3行示例如下:

我们的下游分析，将继续使用dada2的结果，需要将它们改名方便继续分析：

方法2. Deblur

Deblur使用序列错误配置文件将错误的序列与从其来源的真实生物序列相关联，从而得到高质量的序列变异数据，主要为两个步骤。首先，应用基于质量分数的初始质量过滤过程，是Bokulich等人2013年发表的质量过滤方法。

详者注：Deblur是本软件作者作为通讯作者2013发表于Nature Methods的重要扩增子代表序列鉴定方法，截止19年8月25号Google学术统计引用1259次，

引文：Bokulich, Nicholas A., et al. “Quality-filtering vastly improves diversity estimates from Illumina amplicon sequencing.” Nature methods 10.1 (2013): 57. https://doi.org/10.1038/nmeth.2276 作者只将自己的方法作为了备选，分析教程中首选了dada2方法。

按测序碱基质量过滤序列

输出对象:

: 序列质量过滤后结果。 查看 | 下载

: 序列质量过滤后结果统计。 查看 | 下载

注意：在Deblur的论文中，作者使用了当时推荐的过滤参数。而这里使用的参数基于最新的经验，效果更好。

接下来，使方法应用于Deblur工作流程。此方法需要一个用于质量过滤的参数，即截断位置n长度的序列的。通常，Deblur开发人员建议将该值设置为质量分数中位数开始下降至低质量区时的长度。在本次数据上，质量图（在质量过滤之前）表明合理的选择是在115至130序列位置范围内。这是一个主观的评估。你可能不采用该建议的一种原因是存在多个批次测序的元分析。在这种情况的元分析中，比较所有批次的序列长度是否相同，以避免人为引入特定的偏差，全局考虑这些是非常重要的。由于我们已经使用修剪长度为120 bp用于分析，并且由于120 bp是基于质量图的结果，这里我们将使用参数。下一个命令可能需要10分钟才能运行完成。

详者注：deblur最大缺点就是慢，本次只分析了33个样品，共177,092条序列。而实际研究中大项目会有成千上万的样本，1亿-10亿条序列，此步分析可能需要几个月甚至根本无法完成，不推荐。

deblur去噪16S过程，输入文件为质控后的序列，设置截取长度参数，生成结果文件有代表序列、特征表、样本统计。

注：在测试服务器上单线程运行时间为5m29s，比作者测试时间快了1倍

输出结果:

过程统计。 查看 | 下载

: 特征表。 查看 | 下载

: 代表序列。 查看 | 下载

注意： 本节中使用的两种命令生成包含汇总统计信息的QIIME 2对象。为了查看这些汇总统计数据，你可以分别使用和命令来分别可视化这两种命令的输出文件。

输出结果:

: 质量过程统计表，同上面提到的统计表类似。 查看 | 下载

示例如下：包括6列，第一列为样本名称，2-6列分别为总输入读长、总保留高读长、截断的读长、截断后太短的读长和超过最大模糊碱基的读长的数量统计。我们通常只关注2，3列数量即可，其它列常用于异常的输助判断。

deblur-stats.qzv: deblur分析统计表，有分析中每个步骤的统计表 查看 | 下载

图3. deblur去噪和鉴定ASV处理过程统计结果

如果你想用此处结果下游分析，可以改名为下游分析的起始名称：

这处演示不运行下面两行代码，前面添加”#”号代表注释，需要运行请自行删除行首的“#”

详者注：记住，以上两种方法只选择一种即可。推荐dada2速度更快一些，步骤也少一些。有精力的条件下，可以两种方法都试试，比较一下两种方法哪个结果更适合自己。其实每种方法都有存在的意义，而且也有适用的范围，要在具体的项目中，结合背景知识分析哪种方法结果更好时才知道。

特征表和特征序列汇总

FeatureTable and FeatureData summaries

在质量筛选步骤完成之后，你将希望探索数据结果。可以使用以下两个命令进行此操作，这两个命令将创建数据的可视化摘要。特性表汇总命令（）将向你提供关于与每个样品和每个特性相关联的序列数量、这些分布的直方图以及一些相关的汇总统计数据的信息。特征表序列表格命令将提供特征ID到序列的映射，并提供链接以针对NCBI nt数据库轻松BLAST每个序列。当你想要了解关于数据集中重要特性的更多信息时，可视化将在本教程的后续分析中非常有用。

输出结果:

table.qzv: 特征表统计。查看 | 下载

rep-seqs.qzv: 代表序列统计，可点击序列跳转NCBI blast查看相近序列的信息。查看 | 下载

图4. 图中展示了特征表的统计结果

上为摘要、中间为样本数据量分布和图，下方为特征出现频率的统计表和图。

图5. 交互式查看每组剩余样本量

右侧还有进一步查看每个特征的频率和在样本中出现的次数

构建进化树用于多样性分析

Generate a tree for phylogenetic diversity analyses

QIIME 2支持几种系统发育多样性度量方法，包括、和。除了每个样本的特征计数（即QIIME2对象）之外，这些度量还需要将特征彼此关联结合有根进化树。此信息将存储在一个QIIME 2对象的有根系统发育对象中。为了生成系统发育树，我们将使用插件中的工作流程。

首先，工作流程使用程序执行对中的序列进行多序列比对，以创建QIIME 2对象。接下来，流程屏蔽（mask或过滤）对齐的的高度可变区(高变区)，这些位置通常被认为会增加系统发育树的噪声。随后，流程应用基于过滤后的比对结果生成系统发育树。FastTree程序创建的是一个无根树，因此在本节的最后一步中，应用根中点法将树的根放置在无根树中最长端到端距离的中点，从而形成有根树。

详者注：多序列比对和建树在分析中是计算量很大的步骤，本测试数据量很小，只用了14秒，实际上千个样本，可能会使用几十分钟，甚至几小时至几天

输出结果文件:

: 多序列比对结果。查看 | 下载

: 过滤去除高变区后的多序列比对结果。查看 | 下载

: 有根树，用于多样性分析。查看 | 下载

: 无根树。查看 | 下载

Alpha和beta多样性分析

Alpha and beta diversity analysis

QIIME 2的多样性分析使用插件，该插件支持计算α和β多样性指数、并应用相关的统计检验以及生成交互式可视化图表。我们将首先应用方法，该方法将(特征表[频率])抽平到用户指定的测序深度，然后计算几种常用的α和β多样性指数，并使用为每个β多样性指数生成主坐标分析（PCoA）图。默认情况下计算的方法有：

划重点：理解下面4种alpha和beta多样性指数的所代表的生物学意义至关重要。

α多样性

香农(Shannon’s)多样性指数（群落丰富度的定量度量，即包括丰富度和均匀度两个层面）

可观测的OTU(Observed OTUs，群落丰富度的定性度量，只包括丰富度）

Faith’s系统发育多样性（包含特征之间的系统发育关系的群落丰富度的定性度量）

均匀度Evenness（或 Pielou’s均匀度；群落均匀度的度量）

β多样性

Jaccard距离（群落差异的定性度量，即只考虑种类，不考虑丰度）

Bray-Curtis距离（群落差异的定量度量，较常用）

非加权UniFrac距离（包含特征之间的系统发育关系的群落差异定性度量）

加权UniFrac距离（包含特征之间的系统发育关系的群落差异定量度量）

需要提供给这个脚本的一个重要参数是，它是指定重采样（即稀疏/稀疏）深度。因为大多数多样指数对不同样本的不同测序深度敏感，所以这个脚本将随机地将每个样本的测序量重新采样至该参数值。例如，提供，则此步骤将对每个样本中的计数进行无放回抽样，从而使得结果表中的每个样本的总计数为500。如果任何样本的总计数小于该值，那么这些样本将从多样性分析中删除。选择这个值很棘手。我们建议你通过查看上面创建的表文件中呈现的信息并选择一个尽可能高的值（因此每个样本保留更多的序列）同时尽可能少地排除样本来进行选择。

读者思考时间：

查看QIIME 2的 对象，尤其是交互式可视化表格。对于采样深度，应该选择什么值呢？根据这个选择，分析中多少个样本将被排除？在命令中，你将分析的总序列是多少条呢？

译者注：下面多样性分析，需要基于重采样/抽平(rarefaction)标准化的特征表，标准化采用无放回重抽样至序列一致，如何设计样品重采样深度参数呢？

如是数据量都很大，选最小的即可。如果有个别数据量非常小，去除最小值再选最小值。比如此分析最小值为917，我们选择1109深度重采样，即保留了大部分样品用于分析，又去除了数据量过低的异常值。本示例为近10年前测序技术的通量水平，454测序时代抽平至1000条即可，现在看来数据量很小。目录一般采用HiSeq2500或NovaSeq6000的 PE250模式测序，数据量都非常大，通常可以采用3万或5万的标准抽平，仍可保留90%以上样本。过低或过高一般结果也会波动较大，不建议放在一起分析。

计算核心多样性

此步计算耗时9秒。在大数据时，可能会计算更多时间。尤其是样本量增加，计算量会随样本平方增长。

输出对象(13个数据文件):

core-metrics-results/faith_pd_vector.qza: Alpha多样性考虑进化的faith指数。 查看 | 下载

core-metrics-results/unweighted_unifrac_distance_matrix.qza: 无权重unifrac距离矩阵。 查看 | 下载

core-metrics-results/bray_curtis_pcoa_results.qza: 基于Bray-Curtis距离PCoA的结果。 查看 | 下载

core-metrics-results/shannon_vector.qza: Alpha多样性香农指数。 查看 | 下载

core-metrics-results/rarefied_table.qza: 等量重采样后的特征表。 查看 | 下载

core-metrics-results/weighted_unifrac_distance_matrix.qza: 有权重的unifrac距离矩阵。 查看 | 下载

core-metrics-results/jaccard_pcoa_results.qza: jaccard距离PCoA结果。 查看 | 下载

core-metrics-results/observed_otus_vector.qza: Alpha多样性observed otus指数。 查看 | 下载

core-metrics-results/weighted_unifrac_pcoa_results.qza: 基于有权重的unifrac距离的PCoA结果。 查看 | 下载

core-metrics-results/jaccard_distance_matrix.qza: jaccard距离矩阵。 查看 | 下载

core-metrics-results/evenness_vector.qza: Alpha多样性均匀度指数。 查看 | 下载

core-metrics-results/bray_curtis_distance_matrix.qza: Bray-Curtis距离矩阵。 查看 | 下载

core-metrics-results/unweighted_unifrac_pcoa_results.qza: 无权重的unifrac距离的PCoA结果。 查看 | 下载

输出对象(4种可视化结果):

core-metrics-results/unweighted_unifrac_emperor.qzv: 无权重的unifrac距离PCoA结果采用emperor可视化。 查看 | 下载

core-metrics-results/jaccard_emperor.qzv: jaccard距离PCoA结果采用emperor可视化。查看 | 下载

core-metrics-results/bray_curtis_emperor.qzv: Bray-Curtis距离PCoA结果采用emperor可视化。查看 | 下载

core-metrics-results/weighted_unifrac_emperor.qzv: 有权重的unifrac距离PCoA结果采用emperor可视化。查看 | 下载

图6. 以weighted_unifrac距离的PCoA结果交互式可视化为例，可用鼠标托动空间查看每个样本的分布位置。

这里，我们将参数设置为1103。这个值是根据样本中的序列数量来选择的，因为它与接下来几个序列计数较高的样本中的序列数量接近，并且它比序列较少的样本中的序列数量高。这将允许我们保留大部分样品。具有较少序列的三个样本将从分析和任何使用这些结果的下游分析中删除。

注意：根据DADA2特征表汇总选择1103的采样深度。如果使用的是Deblur特性表而不是DADA2特性表，则可能需要选择不同的采样深度。应用上一段的逻辑来帮助你选择合理的采样深度。

注意：在许多Illumina测序结果中，你将观察到一些序列计数非常低的例子。你通常希望通过在此阶段采样深度选择更大的值来从分析中剔除它们。

在计算多样性度量之后，我们可以开始在样本元数据的分组信息或属性值背景下探索样本的微生物组成差异。此信息存在于先前下载的示例元数据文件中。

我们将首先测试分类元数据列和alpha多样性数据之间的关系。我们将在这里为系统发育多样性（群体丰富度的度量）和均匀度进行可视化操作。

Alpha多样性组间显著性分析和可视化

输出可视化结果：

core-metrics-results/faith-pd-group-significance.qzv。查看 | 下载

core-metrics-results/evenness-group-significance.qzv。查看 | 下载

图7. 以为例将互探索不同元数据条件下组间差异，可用鼠标选择不同元数据的列名，切换分组方式，探索对应的生物学意义。

问题：哪些分类样本元数据列与微生物群落丰富度的差异密切相关？这些差异在统计学上有显著性吗？

读者思考时间：实验设计中的那一种分组方法，与微生物群体的丰富度差异相关，这些差异显著吗？

详者注：图中可按选择分类方法，查看不同分组下箱线图间的分布与差别。图形下面的表格，详细详述了组间比较的显著性和假阳性率统计。

结果我们会看到本实验设计的分组方式有, , ，只有身体位置各组间差异明显，且下面统计结果也存在很多组间的显著性差异。

在这个数据集中，连续的样本元数据列（例如，）不与α多样性有相关联，所以我们这里不测试这类关联。如果你有兴趣执行这类测试（对于这个数据集或其他数据集），可以使用命令。

接下来，我们将使用PERMANOVA方法（在Anderson 2001年的文章中首次描述）分析分类型元数据的样本组间差异。以下命令将测试一组样本之间的距离，是否比来自其他组(例如，舌头、左手掌和右手掌)的样本彼此更相似，例如来自同一身体部位(例如肠)的样本。如果你用这个命令的参数，它将执行成对检验，结果将允许我们确定哪对特定组（例如，舌头和肠）彼此不同是否显著不同。这个命令运行起来可能很慢，尤其是当使用参数，因为它是基于置换检验的。因此，我们将在元数据的特定列上运行该命令，而不是在其适用的所有元数据列上运行该命令。这里，我们将使用两个示例元数据列将此应用到未加权的UniFrac距离，如下所示。

输出可视化结果:

: 查看 | 下载

: 查看 | 下载

图8. 不同部分组内和组间差异显著性分析，采用箱线图+统计表呈现

问题：受试者之间的关联和微生物组成的差异在统计学上是否显著？身体部位呢？哪些特定的身体部位对彼此有显著的不同？

同样，我们对于这个数据集所拥有的连续样本元数据中没有一个与样本组成相关，因此这里我们不会测试这些关联。如果你对执行这些测试感兴趣，那么可以使用结合和命令组合使用。

最后，排序是在样本元数据分组间探索微生物群落组成差异的流行方法。我们可以使用Emperor工具在示例元数据下探索主坐标分析（PCoA）绘图。虽然我们的命令已经生成了一些Emperor图，但我们希望传递一个可选的参数，这对于探索时间序列数据非常有用。采于的PCoA结果也是一样的，这使得很容易与Emperor生成新的可视化。我们将采用未加权的UniFrac和Bray-Curtis的PCoA结果生成Emperor图，以便所得到的图将包含主坐标1、主坐标2和实验开始以来的天数(days since the experiment start)的轴。我们将使用最后一个轴来探索这些样本是如何随时间变化的。

输出可视化：

: 查看 | 下载

: 查看 | 下载

图9. 探索样本在第1/2主轴和时间上的分布，调整右侧着色方式和颜色方案可方便观察研究的分类或时间序列结果。

问题：Emperor图是否支持我们在这里执行的其他β多样性分析？（提示：对不同实验元数据进行点着色。）

问题：在未加权的UniFrac和Bray-Curtis PCoA图中，你观察到了哪些差异？

Alpha稀疏取线

Alpha rarefaction plotting

在本节中，我们将使用可视化工具来探索α多样性与采样深度的关系。该可视化工具在多个采样深度处计算一个或多个α多样性指数，范围介于1（可选地控制）和最大采样深度提供值之间。在每个采样深度，将生成10个抽样表，并对表中的所有样本计算alpha多样性指数计算。迭代次数（在每个采样深度计算的稀疏表）可以通过来控制。在每个采样深度，将为每个样本绘制平均多样性值，如果提供样本元数据参数，则可以基于元数据对样本进行分组。

用时13S，本步计算量较大。

输出可视化：

alpha-rarefaction.qzv: 稀疏曲线。 查看 | 下载

图10. 查看按身体部位(body site)分组下可观测(observed) otus的稀疏箱线图，注意观察图中变化以及下面对应样本数据的图。

可视化将有两个图。顶部图是α稀疏图(rarefaction plot)，主要用于确定样品的丰度是否已被完全观察或测序。如果图中的线条在沿x轴的某个采样深度处看起来“平坦(level out)”（即斜率接近于零），这表明收集超过该采样深度的附加序列不太可能观测到新特征。如果绘图中的线条没有变平，这可能是因为尚未充分观察样本的丰富度（由于测序的序列太少），或者它可能是在数据中仍然存在许多测序错误（被误认为是新的多样性）。

当通过元数据对样本进行分组时，此可视化中结果底部的绘图结果非常重要。它说明了当特征表被细化到每个采样深度时，每个组中剩余的样本数量。如果给定的采样深度大于样本的总频率(即，针对样本获得的序列数)，则不可能计算采样深度下样本的多样性。在顶部绘图将不可靠，因为它将计算基于相对少的样本。因此，当通过元数据对样本进行分组时，必须查看底部图表，以确定顶部图表中显示的数据是否可靠的。

注意：提供的参数的值应该通过查看上面创建的文件中呈现的“每个样本的测序量”信息来确定。一般来说，选择一个在中位数附近的值似乎很好用。如果得到的稀疏图中的线看起来没有变平，那么你可能希望增加该值。如果由于大于最大采样深度而丢失了许多样本，则减少该值。

问题1：当通过“body-site”列信息对样本进行分组并查看“observed_otus”指数的α稀疏图时，哪些身体部位显示出足够的多样性覆盖（即稀疏曲线趋于平缓）？在这些身体部位似乎存在多少序列变异？

问题2：当通过“body-site”对样本进行分组并查看“observed_otus”指数的α稀疏图时，“右手掌(right palm)”样本的线看起来在40左右变平，但随后跳到大约140。你认为这里发生了什么？（提示：一定要查看顶部和底部的细节。）

译者注答案：问题2左手掌的多样性从突然40跳至140，而对应的样本量从9个下降为3个(由于测序深度不足)。仅有3次生物学重复样本量太少，偶然性太大，导致的结果波动大但可信度不高。问题1很简单，自己看图吧可以想出答案。

物种组成分析

Taxonomic analysis

在这一节中，我们将开始探索样本的物种组成，并将其与样本元数据再次组合。这个过程的第一步是为的序列进行物种注释。我们将使用经过Naive Bayes分类器预训练的，并由插件来完成这项工作。这个分类器是在上训练的，其中序列被修剪到仅包括来自16S区域的250个碱基，该16S区域在该分析中采用V4区域的515F/806R引物扩增并测序。我们将把这个分类器应用到序列中，并且可以生成从序列到物种注释结果关联的可视化。

注意：物种分类器根据你特定的样品制备和测序参数进行训练时表现最好，包括用于扩增的引物和测序序列的长度。因此，一般来说，你应该按照使用的训练特征分类器的说明来训练自己的物种分类器。我们在数据资源页面上提供了一些通用的分类器，包括基于Silva的16S分类器，不过将来我们可能会停止提供这些分类器，而让用户训练他们自己的分类器，这将与他们的序列数据最相关。