rda冗余分析步骤_R包vegan的基于距离的冗余分析（dbRDA）

R包vegan的基于距离的冗余分析(db-RDA) 常规的冗余分析( RDA )，样方或物种的降维排序实质上反映了欧氏距离。有时我们可能期望关注非欧式距离样方或物种关系的 RDA ，就可使用基于距离的冗余分析( Distance-based redundancy analysis ， db-RDA )来解决。 db-RDA 将主坐标分析( PCoA )计算的样方得分矩阵应用在 RDA 中，其好处是可以基于任意一种距离测度(例如 Bray-curtis 距离等，而不再仅仅局限于欧氏距离)进行 RDA 排序。 因此从概念上来看， db-RDA 是一种回归分析结合主坐标(而非主成分)分析的排序方法，但仍然反映了样方与环境的线性响应关系。 db-RDA 在群落分析中应用广泛。 有关RDA、db-RDA等的基本概念描述，可参考前文 。 前篇已经简介了R包vegan的 常规 RDA 方法。鉴于db-RDA在群落分析中非常重要，本篇单独阐述它。本篇以某16S扩增子测序所得细菌群落数据，简介R包vegan的db-RDA过程。 示例文件、 R 脚本等的百度盘链接： https://pan.baidu.com/s/1itLEWnV_sLpKohHleMwggA

某研究对三种环境中的土壤进行了采样(环境A、B、C，各环境下采集9个土壤样本)，结合二代测序观测了土壤细菌群落物种组成，并对多种土壤理化指标进行了测量。获得两个数据集。

(1)文件“otu_table.txt”为通过16S测序所得的细菌OTU丰度表格，下文统称为“物种”，尽管OTU并不完全等于物种。

(2)文件“env_table.txt”记录了多种土壤理化指标信息，即环境数据。

接下来期望通过db-RDA分析这些数据，用以分析群落物种组成与环境梯度的关系。

vegan包中db-RDA的初步计算

首先展示db-RDA的初步计算过程，对于内容详解见下文。

手动分步计算

按照上述db-RDA原理，不妨先从头走一遍，以便加深理解。

先计算PCoA，再将结果作为RDA的输入，最后被动添加物种得分。方法大致如下。

library(vegan)
#读入物种数据，以细菌 OTU 水平丰度表为例
otu otu  
#读取环境数据
env  
##根据原理一步步计算 db-RDA
#计算样方距离，以 Bray-curtis 距离为例，详情 ?vegdist
dis_bray  
#或者直接使用现有的距离矩阵，这里同样为 Bray-curtis 距离
dis_bray  
#PCoA 排序，这里通过 add = TRUE校正负特征值，详情 ?cmdscale
pcoa  
#提取 PCoA 样方得分(坐标)
pcoa_site  
#db-RDA，环境变量与 PCoA 轴的多元回归
#通过 vegan 包的 RDA 函数 rda() 执行，详情 ?rda
db_rda  
#被动拟合物种得分
v.eig db_rda$CCA$v v.eig db_rda$CA$v  
#先作图展示下，详情 ?plot.cca
#样方展示为点，物种暂且展示为“+”，环境变量为向量
par(mfrow = c(1, 2))
plot(db_rda, type = 'n', display = c('wa', 'cn'), choices = 1:2, scaling = 1, main = 'I型标尺，双序图')
points(db_rda, choices = 1:2, scaling = 1, display = 'wa', pch = 19, col = c(rep('red', 9), rep('orange', 9), rep('green3', 9)), cex = 1)
text(db_rda, choices = 1:2, scaling = 1, display = 'cn', col = 'blue', cex = 0.8)
plot(db_rda, type = 'n', display = c('wa', 'sp', 'cn'), choices = 1:2, scaling = 1, main = 'I型标尺，三序图')
points(db_rda, choices = 1:2, scaling = 1, display = 'sp', pch = 3, col = 'gray', cex = 1)
points(db_rda, choices = 1:2, scaling = 1, display = 'wa', pch = 19, col = c(rep('red', 9), rep('orange', 9), rep('green3', 9)), cex = 1)
text(db_rda, choices = 1:2, scaling = 1, display = 'cn', col = 'blue', cex = 0.8)

capscale()函数

即上述分步过程相比，capscale()是vegan中打包好的db-RDA执行函数，可以直接一步出结果。平常说的“CAP(Canonical analysis of principal coordinates)”就是指它，就是db-RDA。

#或者，capscale() 提供了直接运行的方法，详情 ?capscale
#输入原始数据，指定距离类型，例如样方距离以 Bray-curtis 距离为例
#计算结果中包含物种得分
db_rda  
#或者直接输入已经计算/或读入的距离测度
db_rda #但是这种情况下，物种得分丢失，需要被动添加
v.eig db_rda$CCA$v v.eig db_rda$CA$v  
#先作图展示下
#样方展示为点，物种暂且展示为“+”，环境变量为向量
par(mfrow = c(1, 2))
plot(db_rda, type = 'n', display = c('wa', 'cn'), choices = 1:2, scalin