|
今天与大家分享一种分析方法，与Mantel test非常相似，常用于评估物种与环境或功能的关联度，被称为普鲁克分析（Procrustes Analysis）,也称为普氏分析。为什么分享呢？因为当我们只聚焦于两个数据集一致性时，普氏分析似乎更为直观，略胜一筹。简要介绍该分析方法的原理和适用范围、R语言实现过程、提取相应结果进行绘图。
普鲁克分析（Procrustes analysis）的简介 普鲁克分析（Procrustes analysis）是一种用于形状分布的分析方法。数学上：通过不断迭代，寻找标准形状（canonical shape），并利用最小二乘法寻找每个样本形状到这个标准形状的仿射变化方式。普氏分析可基于不同多元数据集的排序构型（≥2组），通过平移、旋转、缩放等转换方式，实现最大叠合（maximal superimposition），用于各数据集间的对比分析。排序方法可选择PCA、PCoA、NMDS等。换句话说，普氏分析的作用可以看作是一种对原始数据的预处理，目的是为了获取更好的局部变化模型作为后续模型学习的基础。如下图所示，每一个人脸特征点可以用一种单独的颜色表示；经过归一化变化，人脸的结构越来越明显，即脸部特征簇的位置越来越接近他们的平均位置；经过一系列迭代，尺度和旋转的归一化操作，这些特征簇变得更加紧凑，它们的分布越来越能表达人脸表情的变化。（剔除刚性部分、保留柔性部分） 
对于微生态领域，该方法通常用于分析自相同样本不同数据集之间是否存在相似性关系。例如物种与环境、物种与功能基因等等，可以分析两个数据集之间的关联度/是否存在潜在的一致性。
下面我们将利用R语言vegan包中procrustes和protest两个函数来实现普氏分析，并利用ggplot2包完成对应分析的绘图。1) 加载R包和自带数据集（varespec和varechem），并对数据进行相应的距离转换 # Procrustes Anaylsis ----------
# 评估物种群落结构与环境因子间是否具存在显著的相关性（两个数据集）
# 加载R包
library(vegan)
# 加载数据
data(varespec)
head(varespec)
data(varechem)
head(varechem)
# 需要先对数据计算其样本间的距离
# 一般物种群落使用bray-curtis距离，而环境因子使用欧式距离
spe.dist <- vegdist(varespec) # 默认Bray-Curtis
env.dist <- vegdist(scale(varechem), "euclid")查看数据格式，需要注意：无论分析何种数据，两个数据表格的样品均要一一对应；数据的格式如下，就是普通的宽格式的表格，第一行为物种或环境因子名，第一列为样本名。 
2）由于两个数据集的属性不同，需要分别对两个数据集进行降维分析，并提取前两轴坐标（具有数据集代表性的线性组合）进行比较。 # 在进行普鲁克分析前，需要先对两个数据进行降为分析，这里使用的是NMDS
# 也可以使用PCA或者PCoA，然后进行普鲁克分析，计算显著性
mds.s <- monoMDS(spe.dist)
mds.e <- monoMDS(env.dist)3) 进行普氏分析，用到procrustes函数，我们先大概了解一下该函数。 
4) 通过对该函数中各参数的了解，可知X为目标矩阵也就是降维后的环境（功能基因等）坐，Y为降维后的物种数据的坐标，因为后续普氏分析中旋转和缩放操作是针对Y，将Y匹配给X。 另外，当symmetric=FALSE时，处于"非对称"模式，X和Y的分配值调换后，普氏分析的偏差平方和（M2）也会随之改变。而当symmetric=TRUE时，从而给出更合适比例的对称统计。“对称”模式下，X和Y的分配值调换后，普氏分析的偏差平方和（M2）不会发生改变，但注意旋转仍将是非对称的。# 以对称模式为例进行普氏分析（symmetric = TRUE）
pro.s.e <- procrustes(mds.s,mds.e, symmetric = TRUE)
summary(pro.s.e)可得知：通过最小二乘求得组坐标之间的偏差平方和（M2统计量）为0.65
5) 进一步通过各样本之间的残差来判断物种群落与环境因子的一致性。 plot(pro.s.e, kind = 2)
residuals(pro.s.e)从图1来看，如果样本中物种与环境一致性（相似性）越近，则对应的残差越小，反之物种与环境的相似性越远，则残差越大（三条辅助线对应的位置分别为残差25%、50%和75%）；图2为对应的样本残差值。 
6) 事后检验，对偏差平方和M2统计量进行置换999次的普氏检验。由于置换999次检验会存在细微的误差，我们设定了种子数，避免重复存在差异。 # 普氏分析中M2统计量的显著性检验
set.seed(1)
pro.s.e_t <- protest(mds.s,mds.e, permutations = 999)
pro.s.e_t
提取对应的结果：
# 偏差平方和（M2统计量）
pro.s.e_t$ss
# 对应p值结果
pro.s.e_t$signif999次置换检验后显示p<0.001，结果是非常显著的 
基于ggplot2包ggplot函数将其结果绘制成图，并利用export包导出图片。首先提取降维后的数据轴1和2的坐标，并且提取转换的坐标；然后进行绘制。library(ggplot2)
# 获得x和y轴的坐标及旋转过的坐标
Pro_Y <- cbind(data.frame(pro.s.e$Yrot), data.frame(pro.s.e$X))
Pro_X <- data.frame(pro.s.e$rotation)
# 绘图
ggplot(data=Pro_Y) +
geom_segment(aes(x = X1, y = X2,
xend = (X1 + MDS1)/2, yend = (X2 + MDS2)/2),
# geom_segment 绘制两点间的直线
arrow = arrow(length = unit(0, 'cm')),
color = "#9BBB59", size = 1) +
geom_segment(aes(x = (X1 + MDS1)/2, y = (X2 + MDS2)/2,
xend = MDS1, yend = MDS2),
arrow = arrow(length = unit(0.2, 'cm')),
color = "#957DB1", size = 1) +
geom_point(aes(X1, X2), color = "#9BBB59", size = 3, shape = 16) +
geom_point(aes(MDS1, MDS2), color = "#957DB1", size = 3, shape = 16) +
theme(panel.grid = element_blank(), # 绘制背景
panel.background = element_rect(color = 'black',
fill = 'transparent'),
legend.key = element_rect(fill = 'transparent'),
axis.ticks.length = unit(0.4,"lines"),
axis.ticks = element_line(color='black'),
axis.line = element_line(colour = "black"),
axis.title.x=element_text(colour='black', size=14),
axis.title.y=element_text(colour='black', size=14),
axis.text=element_text(colour='black',size=12)) +
labs(x = 'Dimension 1', y = 'Dimension 2', color = '') +
labs(title="Correlation between community and environment") +
geom_vline(xintercept = 0, color = 'gray', linetype = 2, size = 0.3) +
geom_hline(yintercept = 0, color = 'gray', linetype = 2, size = 0.3) +
geom_abline(intercept = 0, slope = Pro_X[1,2]/Pro_X[1,1], size = 0.3) +
geom_abline(intercept = 0, slope = Pro_X[2,2]/Pro_X[2,1], size = 0.3) +
annotate('text', label = 'Procrustes analysis:\n
M2 = 0.6297, p-value = 0.001',
x = -0.3, y = 0.3, size = 4,hjust = 0) +
theme(plot.title = element_text(size=14,colour = "black",
hjust = 0.5,face = "bold"))
library(export)
graph2ppt(file="Procrustes.ppt", append=T, height=5, width=5)
Forcino(2015)也证实了Jackson(1995)和Peres-Neto和Jackson(2001)的说法，即Protest是一种强大而准确用于度量数据集间相似性的分析方法，甚至优于Mantel test。然而，在某些情况下，Mantel test产生的结果与Protest不同。例如，Manteltest更接近于三个更相似的数据集比较中的两个数据集的评估评分的变化幅度。因此，我们在分析时要根据我们的目标和实际情况进行选择，特别是聚焦于两个数据集的一致性时，Procrustes可能是更优的方案。 参考链接： 1.https://www.cnblogs.com/nsnow/p/4745730.html 2.https://www./article/8301917564/ 3.https://www.jianshu.com/p/ec8faacac230 4.https://www./science/article/pii/S0031018215001546 (Evaluating the effectiveness of the Mantel test and Procrustes randomization test for exploratory ecological similarity among paleocommunitie）
|
