5分钟讲清楚怎样把微生物多样性的数据用到文章里

      无论是16S、18S、还是ITS，都是最常用的菌群多样性分析的手段。对于新手，如果收到一份不讲“情面”的16S测序结果析报告，很快就会被各种生态学术语、各种指数、各种分析方法弄得晕头转向。满心喜悦完成实验拿到了结果，如何有效运用这些结果发高分SCI呢？今天小编带大家一起来梳理16S测序报告，迅速理清重点关注对象。

1. 扩增子测序流程

    首先我们要对测序过程有个基础的了解，下图为扩增子测序流程，结题报告包括下图的序列聚类与数据分析两部分。后续所有的数据分析都是基于OTU表进行的，这张表中包括我们测的所有序列中，有多少类型序列（OTU聚类或ASV降噪），然后将聚类完的类型序列与物种数据库比对，得到我们测的所有序列中有哪些物种，每个物种有多少，即丰度的概念。

图：扩增子测序流程

2. 文章常用分析内容

      下图是16s测序的分析内容，展现了主要的分析内容，每个主要分析内容下有更细致的分析内容，这些小的分析内容，其实是可有选择的使用，例如物种组成分析中有物种组成柱状图、分类等级树图、GraPhlAn进化树图、Krona物种组成图，这几个分析的内容都是在解决样本中有哪些细菌，细菌的相对丰度是多少。一般文章中常用的是物种组成柱状图。有的老师会认为，文章中要包含所有的分析内容，其实不然，我们只需要结合自己的研究预期目的选择合适的分析内容即可，并不是多多益善。下面我们就按照SCI文章的要求和写作结构顺序，一一跟大家介绍重点分析内容。

3. 测序结果描述

1. 文章在材料方法描述完以后进入结果部分，首先我们要交代一下这次测序的情况。例如我们测得的原始数据有多少，质控后的有效数据有多少，每个样本质控后的有效数据占原始数据的百分比均高于百分之多少（90%），当然所占百分比越高，本次测序的效果越好。

2. 然后我们可以描述一下稀释曲线图，稀释曲线（Rarefaction curve）主要利用各样本在不同测序深度时的微生物Alpha多样性指数构建曲线，以此反映各样本在不同测序数量时的微生物多样性。若多样性指数为sobs（表征实际观测到的物种数目），当曲线趋向平坦时，说明测序数据量合理，更多的数据量只会产生少量新的物种，反之则表明继续测序还可能产生较多新的物种，可以反映样本中绝大多数的微生物多样性信息。

3. 分类学注释，总体描述本次测序结果，例如，本次测序共得到多少门、纲、目、科、属、种、OTU。

图1. 稀释曲线图

4. 物种组成分析

    分类学组成分析，主要关注物种组成柱状图即可，这是测序中最基本的分析，通过该分析可以得到，每个样本中有哪些物种组成，每个物种的相对丰度是多少。需要注意的是，该分析需要在具体的分类水平上看（门、纲、目、科、属），由于二代测序在种水平上的准确率相对较低，所以在种水平上分析的较少，一般常在门和属水平分析。描述可以从两方面描述，一方面看主要的物种组成是什么（有哪些菌，在某一分类的具体水平），不同的比较组之间物种组成是否有差异。另一方面看同一物种的丰度在不同组间的差异，具体描述为升高或降低。从物种组成柱状图中就可以来比较不同组间物种组成是否有差异，但这个差异不具有显著性之说。 

图2. 物种组成柱状图

1. Alpha多样性分析

    Alpha多样性分析，稀释曲线中已经用到了了Alpha多样性指数，这里的分析更加具体，可以对分组样本进行比较，运用统计学T检验的方法，检测每两组之间的指数值是否具有显著性差异。可以选择自己要用的指数进行比较，一般常用丰富度指数Chao值，多样性指数Shannon，具体描述可以是丰富度与多样性指数升高或降低，不同组间指数是否有显著差异。Alpha多样性分析将样本的菌群群整体研究并转换为具体的指数与p值，通过数字说明群落的变化与差异，在文章中是不可或缺的。

图3. Alpha多样性指数图

2. Beta多样性分析

    Beta多样性分析，该分析是十分重要且常见的分析，几乎所有的微生物测序文章中都会有此分析内容，例如我们常见的PCoA、NMDS、PCA分析。通俗讲就是每一个物种(ASV/OTU)在两个样本之间差异，即是反映这两个样本间群落差异的一个维度。而由于群落中物种的数量巨大，样本间群落的差异往往是多维度的，难以进行直接比较。此时，就需要相关算法对多维数据进行降维。降维后的结果就是平面图中样本点的距离，通过距离的远近代表样本群落结构的相似与差异，这里要注意的是群落结构，因为我们是将每个样本中的每个物种都进行了比较，不是只针对具体的某些物种，所以是从整体分析样本之间的相似性与差异性。例如PCoA分析，我们要注意算法的选择，如不考虑物种丰度的unweighted_unifrac，考虑物种丰度的weighted_unifrac，还有常用的bray_curtis，不同的算法会导致不同的结果，还需要各位投稿人结合自己的研究目的进行选择。

图4. PCoA分析

3. 组间差异显著性分析

    该分析需要与Beta分析结合使用，因为PCoA、NMDS等排序分析只是一种探索分析手段，并不是统计检验，对于排序图中呈现的分布规律，我们需要使用检验手段进行验证，常用的分析有Permanova，Anosim。我们在文章中将Beta分析与p值结合来说明，不同的处理是否使样本之间的微生物群落结构产生显著差异。需要注意的是显著性分析的算法选择要与Beta多样性分析相一致。

5. 物种差异分析与标志物种

1. 韦恩图

    物种差异分析是所有老师最关注的分析。先来说韦恩图，该分析可以知道不同的样本间有哪些物种是共有的，哪些是独有的，通过数据说话，如果共有的物种数占所有物种数的比例高，那么样本的相似性就高。如果样本中各自独有的物种数所占比例高，那么样本之间的相似性差，各自的独特性高。

图5. 样本组Venn图

2. 聚类热图

    物种组成热图也是我们常见的分析，将物种丰度转换成颜色色调的冷暖，此分析看起来难以理解，其实只是视觉上的错觉。热图是在具体分类学水平分析所有样本中总丰度排在前几十的物种，一般是在属水平进行分析。因为筛选了丰度，所以热图是对优势物种进行分析。此分析包括两方面内容，一方面可以描述比较具体某一优势物种在各个样本中的丰度高低，另一方面是描述样本的聚类情况，通过聚类树可以知道优势物种组成在样本之间的相似性和差异性。

图6. 双聚类的物种组成热图

3. 差异菌群的筛选

    通过差异菌群筛选可以找出不同的处理组导致哪些菌群发生变化，丰度上具有显著差异的物种，是扩增子测序的核心目的。例如常用的Metagenomeseq和Lefse分析。MetagenomeSeq主要用于两组样本的差异比较。Lefse的分析更加严格，可以进行多组比较，而且是在所有分类学水平同时进行分析。Lefse分析的结果包含了多级物种层级树图和LDA效应值柱状图，文章中放LDA效应值柱状图即可，LDA效应值柱状图需要注意效应值即该图中的横坐标值，文章中一般默认最低值为2，筛选出的物种为该组的标志物种，这些物种的丰度显著高于在其他组中的丰度。

图7. 基于分类等级树的组间差异分类单元图

图8. 标志物种的LDA效应值柱状图

4. 功能预测分析

    功能预测分析，由于扩增子测序是基因片段测序，不像宏基因组将整个基因组实际进行测序，所以扩增子测序只能是预测功能，与宏基因组的功能相比，预测到的信息是有限的。基于具体的功能数据库进行注释，得到功能通路丰度图，通过功能预测可以知道有哪些功能，哪些功能的丰度相对高，也可以为进一步做宏基因功能提供参考，另外也使文章内容更加完善。

图9. 预测的KEGG二级功能通路丰度图

6. 关联网络分析

    关联网络分析，可以用来展示样本与物种之间的分布情况，通过对不同样本间的物种丰度信息进行相关性分析。在微生物生态学的关联网络图中，每个点可以代表群落中的一个ASV/OTU，两个点之间是正相关还是负相关可以通过连线的颜色进行区分。此分析的作用是探究样本中是否存在关系紧密的微生物，发挥特定的功能，并影响整个群落的组成变化，所以文章中是否选择此分析必须要结合老师的研究思路。

图10. 关联网络分析

    这里举个例子，环境因子关联分析，研究土壤的老师大都会收集土壤的各种理化因子做为环境因子数据，常见的分析如RDA和CCA分析，通过该分析可以反映环境因子、样本、菌群三者之间的关系或者两两之间的关系，可得到影响样本分布的重要环境驱动因子，文章中可以描述哪些环境因子显著影响样本的微生物群落。

图11. CCA分析图

总结

    微生物扩增子测序中还有许多其它分析内容，以上是一些常用的基础分析，在文章结构设计时，可以选择其中几个分析进行组合，有些是必须要有的，有些则需要结合实验研究目的进行选择。对于初次接触测序的老师，一定不要被分析图片的美貌所引诱，看起来复杂和精美的图片可能并不能很好的解释我们要讨论的结果，反而会增加分析难度，根据文章目的取材，掌握最基础的分析内容，足以支撑起一篇文章的结构。