「R数据处理」共线性检验和CAP分析

「R数据处理」共线性检验和CAP分析前言:在进行物种多样性分析的过程中,排序分析有时候能够用来表征样方之间物种的β多样性(即群落的组成变化)。排序分析可以分为间接梯度排序(indirect gradient analysis)和直接梯度排序(direct gradient a

欢迎大家来到IT世界,在知识的湖畔探索吧!

前言:在进行物种多样性分析的过程中,排序分析有时候能够用来表征样方之间物种的β多样性(即群落的组成变化)。排序分析可以分为间接梯度排序(indirect gradient analysis)和直接梯度排序(direct gradient analysis)。间接梯度排序又叫非约束性排序,该分析是为了寻求潜在的或在间接的环境梯度下解释物种数据的变化包括PCA,PCoA,NDMS;直接排序又叫约束性排序。它是指在特定的梯度上(环境轴)上探讨物种的变化情况,包括 RDA, CCA, CAP或dbRDA。排序的目的是在一个可视化的低维空间或平面重新排列目标样本,使样本之间的距离最大程度地反映出平面散点图内样本之间的关系信息。

1. 共线性检验

在进行排序分析之前尤其是在直接梯度排序分析前,一般需要检验使用的环境因子之间是否存在共线性问题。因为,环境因子(相当于自变量)之间如果存在严重的共线性会使物种矩阵被过度解读,从而使排序结果不准确。R语言中Hmisc包提供的varclus函数可以基于聚类分析的算法评估目标变量的共线性和冗余性。

代码如下:

setwd(choose.dir()) rm(list = ls()) library(vegan) library(tidyverse) library(Hmisc) #import data data(varespec)# species data data(varechem)# env data #Co-linearity analysis env<-as.matrix(varechem) co_linear <- varclus(env, similarity="spear") # spearman is the default co_linear plot(co_linear)#

欢迎大家来到IT世界,在知识的湖畔探索吧!

「R数据处理」共线性检验和CAP分析

一般来说当spearman ρ2 > 0.7时则表明分支上的变量可能存在严重的共线性,此时应该根据研究的目的去除其中一个变量。

这里经过检验我们发现铁和铝、钙和镁、钾和硫这6个变量两两间存在显著的共线性,因此我们只取共线性变量的其中之一用于接下来的CAP分析。这里我们选择铝、钙和钾这三个变量与其它变量组成环境因子矩阵。

欢迎大家来到IT世界,在知识的湖畔探索吧!#select env %>% data.frame() %>% dplyr::select(-Fe,-Mg,-S) ->env2 env2
「R数据处理」共线性检验和CAP分析

2. CAP分析

CAP分析使用的是vegan包中的capscale函数。该函数最初是作为约束性近似分析的一个变体而开发的(Anderson & Willis 2003),但这些发展使它与dbRDA相似。然而,它抛弃了具有负特征值的虚构维度,排序和显著性检验区域只基于真实维度和正特征值

与RDA和CCA分析不同的是CAP分析必须控制(约束)某些环境因子,评估未被约束的环境因子对物种矩阵的解释量。这里我们控制Baresoil和Humdepth两个变量。

vare_cap <- capscale(varespec ~. + Condition(Humdepth+Baresoil ), env2, dist="bray") vare_cap
「R数据处理」共线性检验和CAP分析

结果表明了约束部分解释了39.21%的物种变化,控制部分解释了25.09%的物种变化,未被解释的部分有41.39%。

此时我们能看到存在一个Imaginary的解释量,这意味着可能存在负引擎,为了避免这种情况我们可以把负值去除。代码如下:

欢迎大家来到IT世界,在知识的湖畔探索吧!#remove negative eigenvalues with additive constant vare_cap2 <- capscale(varespec ~. + Condition(Humdepth+Baresoil ), env2, dist="bray", add =TRUE) vare_cap2
「R数据处理」共线性检验和CAP分析

此时的结果表明了约束部分解释了37.61%的物种变化,控制部分解释了20.49%的物种变化,未被解释的部分有41.91%。并且没有了Imaginary的解释量。

进一步查看每个轴的重要性,代码如下:

summary(vare_cap2)
「R数据处理」共线性检验和CAP分析

结果表明CAP1解释了15.02%,CAP2解释了12.82%。

plot(vare_cap2)

「R数据处理」共线性检验和CAP分析

3. Anova检验

1)模型整体的显著性检验

anova(vare_cap2)
「R数据处理」共线性检验和CAP分析

结果表明我们用来做例子构建的模型并不显著

2)单一变量的显著性

anova(vare_cap2,by = "term")
「R数据处理」共线性检验和CAP分析

结果表明这些变量中只有N和Al显著影响物种矩阵。

3)轴的显著性检验

anova(vare_cap2,by = "axis")
「R数据处理」共线性检验和CAP分析

轴检验表明了所有的轴对物种矩阵的变化解释量都不显著。

4. 画图

一般来说anova检验不显著就没必要继续画图了,应该考虑选择其它的模型,这里为了演示,我们对以上数据进行画图。

代码如下:

#-----------------plot--------------------- scrs<-scores(vare_cap2,display=c("sp","wa","lc","bp","cn")) names(scrs)

#计算箭头变化倍数并计算箭头的长度

multiplier <- vegan:::ordiArrowMul(scrs$biplot) multiplier df_arrows<- scrs$biplot*multiplier#这一步是计算箭头的长度 df_arrows colnames(df_arrows) rownames(df_arrows) #添加样方坐标 site<-scrs$sites %>% as_tibble(rownames = "sample") site group <- env2 %>% as_tibble(rownames = "sample") %>% mutate(group = if_else(Humdepth>2,"deep","superficial")) %>% select(sample,group) site2 <- site %>% left_join(group, by = "sample") df_arrows=as.data.frame(df_arrows) ggplot(data=site2,aes(CAP1,CAP2)) + geom_hline(aes(yintercept=0),colour="#d8d6d6",linetype=5)+ geom_vline(aes(xintercept=0),colour="#d8d6d6",linetype=5)+ geom_point(aes(color=group),size=3.5,shape=16) + geom_segment(data=df_arrows, aes(x = 0, y = 0, xend = CAP1, yend = CAP2), arrow = arrow(angle= 15,length = unit(0.25, "cm")),color="#000000",alpha=0.5)+#画箭头 geom_text(data=as.data.frame(df_arrows*1.1),aes(CAP1, CAP2, label = rownames(df_arrows)), color="#000000",alpha=0.7)+#箭头对应的信息 scale_color_manual(values = c("#2692d6","#8cc1e3"))+#点的颜色 theme_bw() + theme(panel.grid=element_blank())+ labs(x = "CAP1 (15.02%)", y = "CAP2 (12.82%)")
「R数据处理」共线性检验和CAP分析

免责声明:本站所有文章内容,图片,视频等均是来源于用户投稿和互联网及文摘转载整编而成,不代表本站观点,不承担相关法律责任。其著作权各归其原作者或其出版社所有。如发现本站有涉嫌抄袭侵权/违法违规的内容,侵犯到您的权益,请在线联系站长,一经查实,本站将立刻删除。 本文来自网络,若有侵权,请联系删除,如若转载,请注明出处:https://itzsg.com/75618.html

(0)

相关推荐

发表回复

您的邮箱地址不会被公开。 必填项已用 * 标注

联系我们YX

mu99908888

在线咨询: 微信交谈

邮件:itzsgw@126.com

工作时间:时刻准备着!

关注微信