8. 群落结构差异显著性

非参数置换的多变量方差分析(ANONIS)、相似性分析(ANOSIM)和多反应置换法(MRPP)常被用来评估不同采样点之间微生物群落组成差异的显著性。

一般来说，当我们通过降维、聚类等分析后，可能会观察到微生物群落在不同批次、不同处理组、不同采样点等会有不同的聚集、分类的情况（差异），但这些差异是否显著则需上述统计检验后才能得出。

Update

2021/01/29：把 MRPP分析的函数和输出数据整理简单的封装了一下，输出结果为标准数据框格式。

1. MRPP分析

MRPP（Multi Response Permutation Procedure）分析，是一种用于分析高维度数据组间相似性的统计方法，它以距离矩阵为基础，用于判断组间数据是否具有显著性差异。在生态统计中，可以使用MRPP分析（往往同时配合排序分析来使用），查看不同环境的群落组成结构差异是否显著。

若导入数据为OTU丰度表时，需要使用distance参数指定所依据计算的距离类型获得距离矩阵；或者也可以首先使用vegan包中的vegdist()计算样本间距离，然后再将所得距离数据作为mrpp()的输入。之后指定各样本所属的分组信息，且要注意分组信息中的分组名称与OTU丰度表中样本名称要按顺序对应（可以用match()预先排列好二者的顺序），以及置换检验的次数等，进行MRPP分析。

整体水平差异检验

首先在所有分组水平上，检验整体差异，即检验细菌群落结构组成在整体上是否不具备一致性。
参数解释：

otu_table：常见OTUs丰度表
group_name：分组的列名
metadata：分组和样本信息

library(vegan)
single_mrpp <- function(otu_table, group_name, metadata){
  Transp_otu <- t(otu_table)
  group <- as.vector(unlist(metadata[group_name]))
  
  mrpp_all <- mrpp(Transp_otu, group, 
                   permutations = 999, 
                   distance = 'bray')
  
  results <- as.data.frame(cbind(group_name,
                                 mrpp_all$distance,
                                 round(mrpp_all$A, 3),
                                 round(mrpp_all$delta, 3),
                                 round(mrpp_all$`E.delta`, 3),
                                 round(mrpp_all$Pvalue, 3)))
  
  names(results) <- c("Group", "Distance", "Corrected_A", "Observeed_delta", 
                      "Expected_delta", "P_value")
  
  return(results)
  }

mrpp_all <- single_mrpp(otu, "species", metadata)

结果解读：其中，Significance of delta("P_value)，即显著性p值，小于0.05说明差异显著；Chance corrected within-group agreement A(Corrected_A)，即A值，大于0说明组间差异大于组内差异，小于0说明组内差异大于组间差异；observe delta值越小说明组内差异小，expect delta值越大说明组间差异大。

组间比较

通过A值和observe delta值发现组内有差异，则需要进行组间多重比较。
一般通过循环处理，进行小分组间比较。

Notice：

需要将样本的列名设置为"sample"。

pairwise_mrpp <- function(otu_table, group_name, metadata){
  if (! "sample" %in% colnames(metadata)){
    print("The colnames of your sample is not sample")
    stop()
}else{
  group <- unique(as.vector(unlist(metadata[group_name])))
  Transp_otu <- t(otu_table)
  results <- NULL
  for (i in 1:(length(group) - 1)){
    for (j in (i + 1):length(group)){
      group_ij <- subset(metadata, as.vector(unlist(metadata[group_name])) %in% c(group[i], group[j]))
      otu_ij <- Transp_otu[group_ij$sample,  ]
      grouping <- as.vector(unlist(group_ij[group_name]))
      
      mrpp_result_ij <- mrpp(otu_ij, grouping, 
                                 permutations = 999, 
                                 distance = 'bray') 
      mrpp_result_ij
      results
      results <- rbind(results, 
                       data.frame(paste(group[i], group[j], sep = '/'), 
                                  mrpp_result_ij$distance,
                                  round(mrpp_result_ij$A, 3), 
                                  round( mrpp_result_ij$delta, 3), 
                                  round(mrpp_result_ij$`E.delta`, 3), 
                                  round( mrpp_result_ij$Pvalue, 3)))
    }
  }
}
  results <- data.frame(results, stringsAsFactors = FALSE)
  names(results) <- c("Group", "Distance", "Corrected_A", "Observeed_delta", 
                     "Expected_delta", "P_value")
  return(results)
  }

pairwise_results<- pairwise_mrpp(otu, "species", metadata)

2. 相似性分析（ANOSIM）

在生态统计中，可以使用ANOSIM分析（往往同时配合排序分析来使用，对于ANOSIM，更常与NMDS排序分析（ metaMDS() vegan{}）放一起说明问题），查看不同环境的群落组成结构差异是否显著，组间群落差异是否显著不同于组内差异。

整体水平差异检验

# 使用 OTU 丰度表，则需要在计算时指定所依据的距离类型，这里依然使用 Bray-Curtis 距离
anosim_result_otu <- anosim(otu, group_data$group, permutations = 999, distance = 'bray')        

# 首先根据丰度表计算样本距离，在将所的距离数据作为输入
bray_dis <- vegdist(otu, method = 'bray')
anosim_result_bray <- anosim(bray_dis, group_data$group, permutations = 999)   
# 可将 Dissimilarity ranks 拿出来做一个箱形图可视化

ANOSIM分析基于两两样本之间的距离值排序获得的秩，这样任一两两组的比较可以获得三个分类的数据。以箱线图的形式展示组间与组内的秩的分布，横坐标表示所有样品（Between）以及各分组，纵坐标表示距离（本示例使用Bray-Curtis距离）的秩。当Between组相对与其他每个分组的秩较高时，则表明组间差异大于组内差异。同时图的上方标注了R值与p值两个重要统计指标，便于我们直观地对组间差异是否显著不同于各组内的差异进行判断。一般来说，R>0说明组内距离小于组间距离，也即分组是有效的。

组间比较

同样的，检测出差异后两两比较。也推荐使用 OTU 丰度表作为输入数据，每次筛选分组后重新计算样本距离，避免由于样本数减少可能导致的距离变动而造成误差。

image.png

3. 置换多元方差分析（PERMANOVA）

置换多元方差分析（Permutational multivariate analysis of variance，PERMANOVA），又称非参数多因素方差分析（nonparametric multivariate analysis of variance）、或者ADONIS分析，其本质是基于F统计的方差分析，依据距离矩阵对总方差进行分解的非参数多元方差分析方法。使用PERMANOVA可分析不同分组因素对样品差异的解释度，并使用置换检验进行显著性统计。在生态统计中，可以使用PERMANOVA（往往同时配合排序分析来使用，对于PERMANOVA，更常与PCoA排序分析放一起说明问题），查看不同环境的群落组成结构差异是否显著。

整体水平

所有分组间比较，即整体差异.

#距离矩阵，则直接使用现有的距离矩阵进行分析即可
adonis_result_dis <- adonis(dis~group, group_data, permutations = 999) #根据 group$site 这一列样本分组信息进行 PERMANOVA 分析，随机置换检验 999 次

#（2）若是使用 OTU 丰度表，则需要在计算时指定所依据的距离类型，这里依然使用 Bray-Curtis 距离
adonis_result_otu <- adonis(otu~group, group_data, permutations = 999, distance = 'bray')   #这条命令的详情和上述命令所表示的信息一致

# 或者根据丰度表计算样本距离，在将所的距离数据作为输入
bray_dis <- vegdist(otu, method = 'bray')
adonis_result_bray <- adonis(bray_dis~group, group_data, permutations = 999)     #同上所述

对于各项统计内容：Df，自由度，其值=所比较的分组数量-1；SumsOfSqs，即Sums of squares，总方差，又称离差平方和；MeanSqs，即Mean squares，均方（差）；F.Model，F检验值；R2，即Variation (R2)，方差贡献，表示不同分组对样品差异的解释度，即分组方差与总方差的比值，R2越大表示分组对差异的解释度越高；Pr (>F)，显著性p值，小于0.05说明差异显著。

组间比较

利用循环处理，进行小分组间比较。推荐使用 OTU 丰度表作为输入数据，每次筛选分组后重新计算样本距离，避免由于样本数减少可能导致的距离变动而造成误差

4. 相对丰度差异可视化

通过上述成对显著性分析，得出数据后可以用带星号或者其他显著性标记的柱状图展示。