扩增子SCI套路--2组间差异分析
本文转载自“态昌基因”，已获授权。
差异分析一直是高通量测序数据分析的核心部分。在上期的系列文章中，小昌与大家分享了确定微生物群落组间是否存在差异的常见分析手段，有PCA、Pcoa、NMDS等，在实际的微生物群落研究中，我们常常需要进一步找出是哪些菌（群）引起了群落的差异。只有找出核心影响菌（群），我们才能更明确下一步研究方向。
那么，要找到具体的差异菌（群），主要有哪些方法呢？
1. 统计差异分析
统计差异分析最常用的就是Wilcox秩和检验（Wilcoxon rank-sum test），也叫曼-惠特尼U检验（Mann–Whitney U test），是两组独立样本非参数检验的一种方法。其原假设为两组独立样本总体分布无显著差异，通过对两组样本平均秩的研究来实现判断两组总体的分布是否存在差异。根据p值筛出组间差异显著的物种（一般在门、属和OTU三个水平较常见），随后可用柱状图、盒状图或热图直观展示差异菌群的丰度变化[。也可直接在Windows中用STAMP软件进行差异统计分析和直观展示。
分析步骤如下：
a. 根据物种丰度，通过wilcox 秩和检验来计算每个物种在门(phylum)、属(genus)、OTU(species)上组间的p-value；
b. 用FDR(false discovery rate，错误发现率)方法校正后得到q-value；
c. 根据q-value筛选出在两组间显著差异的物种。
一般会得到统计表格：
说明：mean分别为两组样品物种的平均相对丰度，sd分别是两组样本物种相对丰度的标准差。P值为对两组检验原假设为真的概率值，p<0.05表示存在差异，p<0.01表示差异显著，q值为假发现率。
2. 机器学习分类器
机器学习分类器最常用的是随机森林（Random Forest）的方法。随机森林方法可筛选出对分组效果贡献最高的OTU或物种列表。一般筛选出OTU或物种列表后，可进一步作图展示，如柱状图，热图等。
分析步骤如下：
a. 去除测序量过低的样本以及低丰度OTUs或物种；
b. 用随机森林算法计算，得出对分组效果有贡献度的OTU或物种及其得分值列表；
c. 选取得分值高（如大于0.001）的OTU或物种。
部分结果如下：
说明：一般地，选取Mean_decrease_in_accuracy值大于0.05的OTU，作进一步分析；对于组间差异较小的样本，该值可能会降至0.03。
需要注意的是Wilcox检验强调统计学上的差异，随机森林强调对分组准确率的贡献，二者有较高的一致性，但结果还是存在一定的差异，尤其是重要性排序：统计差异分析检验结果是按照差异程度（p值）依次排序，随机森林是按照对分组贡献重要性依次排序。
下面我们看看常见的差异分析的作图吧。
1）通过Wilcox秩和检验分析得到的三组样本间差异物种（Streptococcus），以盒状图展示。
2）通过Wilcox秩和检验分析得到的两组样本间差异物种，并以柱形图展示。
3）通过Stamp软件中的Welch’s t-test分析得出的两组样本间的显著性差异物种，以及该物种在不同组的分布情况。
4）使用随机森林算法或Wilcoxon秩和检验分析挑取各组样本中存在显著差异OTU。 根据挑选出来的差异OTU，根据其在每个样品中的丰度信息，对物种进行聚类，绘制成热图。
说明：图中越接近蓝色表示物种丰度越低，越接近橙红色表示丰度越高。左边的聚类树是根据各物种间的spearman相关性距离进行聚类；上边的聚类树是采用样本间距离算法中最常用的Bray-Curtis算法进行聚类。图中可按样本分组标上不同的颜色，也可根据样本顺序作图而不画出聚类树。
3. LEfSe分析
LEfSe分析即LDA Effect Size分析，综合了统计学上的差异分析和该差异物种对分组结果的影响力得分值，考量面比较全，在文章中出现频率也比较高，一般在属以上水平做分析和展示。更多相关统计细节请阅读《LEfSe分析，你真的懂了么》。
分析步骤如下：
a. 使用non-parametric factorial Kruskal-Wallis(KW) sum-rank test(非参数因子克鲁斯卡尔-沃利斯秩和检验)检测具有显著丰度差异特征，并找到与丰度有显著性差异的类群；
b. 采用线性判别分析(LDA)来估算每个组分(物种)丰度对差异效果影响的大小；
c. 找出对样品分组产生显著性差异影响的群落或物种。
说明：左边的图为统计两个组别当中有显著作用的微生物类群通过LDA分析（线性回归分析）后获得的LDA分值。右边的图为聚类树，节点大小表示丰度，默认从门到属依次向外排列。红色区域和绿色区域表示不同分组，树枝中红色节点表示在红色组别中起到重要作用的微生物类群，绿色节点表示在绿色组别中起到重要作用的微生物类群，黄色节点表示的是在两组中均没有起到重要作用的微生物类群。图中英文字母表示的物种名称在右侧图例中进行展示。在结果 lefse_LDA.xls 表中：第一列为物种名称；第二列为Abundance丰度值，为丰度较高组的丰度；第三列为 Group，展示差异组名；第四列为 LDA 值，第五列为 P值。默认 LDA 值大于 2，P 值小于 0.05，该物种为差异物种，差异组即两组中（或多组中）丰度高的一组。若Group、 LDA、 P 值均为空，则表示该物种在组间无差异。
了解了怎么判断样本组间微生物群落是否存在差异，知道了如何寻找差异菌群，就万事大吉了吗？当然不是，我们还可以进一步找出引起菌群差异的因素，so 下期再续喽！
菌群的差异分析方法相信大家都get到了，最后来点简单的 ，大家是怎么区分他们的呢？
猜你喜欢
写在后面
为鼓励读者交流、快速解决科研困难，我们建立了“宏基因组”专业讨论群，目前己有国内外七十多位PI，七百多名一线科研人员加入。参与讨论，获得专业指导、问题解答，欢迎分享此文至朋友圈，并扫码加主编好友带你入群，务必备注“姓名-单位-研究方向-职务”。技术问题寻求帮助，首先阅读《如何优雅的提问》学习解决问题思路，仍末解决群内讨论，问题不私聊，帮助同行。
学习16S扩增子、宏基因组科研思路和分析实战，关注“宏基因组”