NiuGuoxiang的个人博客分享 http://blog.sciencenet.cn/u/NiuGuoxiang

博文

R语言RDA分析全过程梳理及注意事项

已有 14846 次阅读 2020-2-10 21:11 |个人分类:统计学习|系统分类:科研笔记| RDA, 解释率, 层次分割

R语言RDA分析全过程梳理及注意事项

RDA分析不是简单的一行代码,有很多步骤,若舍去这些步骤,RDA分析就没有了意义。通过RDA分析,不仅需要给出哪些因子最重要,还要给出每个因子的解释率,最好有显著性检验。这样结果更明确,更加有说服力,这里梳理一下过程,避免以后分析有所遗漏。RDA分析详细步骤具体看刘尧科学网博客http://blog.sciencenet.cn/blog-3406804-1182489.html 有些语句也用了这里面的。

 

library(vegan)

library(permute)   

library(lattice)

library(ggrepel)  #作图使用此包,和ggplot2一起用,使文字尽可能不重叠

library (ggplot2)

datarda<-read.table("clipboard",header=T) # 物种数据

#datarda-hel<- decostand(datarda, method = 'hellinger') # 物种数据有时需要转化,特别是有0出现的情况,基于转化的RDA分析tb-RDA

env<-read.table("clipboard",header=T)  #环境数据,有时环境数据需要数据转换以满足数据分析要求,但很多文章里土壤pH是不做转换的

rda_result<-rda(datarda~.,env,scale=FALSE)

summary(rda_result)

r2<-RsquareAdj(rda_result)

rda_noadj<-r2$r.squared #原始R2,不能直接使用

rda_adj <- r2$adj.r.squared  #校正后的R2乘以RDA1和2的解释度,才是作图用的实际解释度

anova.cca(rda_result2,step=1000)  #全模型检验,有时文章需要提供

anova.cca(rda_result2,by='axis',step=1000) # 逐轴检验,看每条轴是否显著,原则上不显著的轴要舍去

#anova(rda_result, by = "term", permutations=999) #这里的anova不是单因素方差分析,只是名字一样,气人不?这里该语句用于逐因子检验,会给出每个环境因子的r2以及显著性P,但不能排除共线性的问题,即每个因子的r2和P值依赖于环境因子的顺序,排在前面的环境因子可能解释度最高且显著,后面的因子虽然不显著可能是因为被前面的解释过了。所以不能直接使用这部分结果。

vif.cca(rda_result) #计算方差膨胀系数,有些因子膨胀系数过大(达到几十或几百),需要排除掉一些因子,可进行前向选择保留必要因子

rda_tb_forward_r<-ordiR2step(rda(datarda~1,env,scale= FALSE), scope = formula(rda_result),  direction = 'forward', permutations = 999) ##前向选择,前向选择会保留一些因子,但总体解释的R2基本不变,但有时共线的是核心因子,保留的因子却是在逻辑上无关紧要的,很麻烦

plot(rda_tb_forward_r, scaling = 1, main = '前向选择后,I 型标尺', display = c('wa', 'cn'))

rda_tb_forward_vp <- varpart(phylum_hel, env['pH'], env[c('DOC', 'SOM', 'AP', 'AK', 'NH4')])##前向选择后对结果执行方差分解,可以将环境因子分成不同组,如该语句将pH归为一组,DOC', 'SOM', 'AP', 'AK', 'NH4 归为一组,结果可以使用韦恩图表示出来。

plot(rda_tb_forward_vp, digits = 2, Xnames = c('TC', 'forward_env'), bg = c('blue', 'red'))##韦恩图

anova(rda(phylum_hel,env['pH'],env[c('DOC','SOM','AP','AK','NH4')]),permutations=999) ##同样对执行完前向选择的结果进行检验,这里环境因子的R2和显著性可以使用

rda_tb_forward_vp <- varpart(phylum_hel, env['TC'], env[c('pH', 'DOC', 'SOM', 'AP', 'AK', 'NH4')]) ##剔除的环境因子并不能保证其对解释变量无影响,可以与保留的环境因子一起分析,以查看剔除的环境因子是否对解释变量同样有影响。

plot(rda_tb_forward_vp, digits = 2, Xnames = c('TC', 'forward_env'), bg = c('blue', 'red')) #图

 

之前说了,因为共线性问题,环境因子对解释变量是否有影响以及每个因子的贡献不易判断,所以可以使用层次分割来分析每个因子的贡献以及显著性。当然增强回归树模型同样可以做出每个因子贡献,但只能适用于解释变量为一个的情况。赖江山老师开发的新包rdacca.hp可以结合RDA或CCA分析进行层次分割。具体看http://blog.sciencenet.cn/blog-267448-1200984.html

library(devtools)

library(hier.part)

library(scales)

library(rdacca.hp)

rdacca.hp (datarda,env,pieplot = "tv",type="RDA")   # 每个解释变量占总变化量的比例

rdacca.hp (datarda,env,pieplot = "tev",type="RDA")  #每个解释变量占总解释变化量的比例

permu.hp (datarda,env, type="RDA", permutations=999) #层次分割,可给出所有因子的贡献以及显著性,个人感觉层次分割要比前向选择后显著性检验要简单一些。




http://wap.sciencenet.cn/blog-3366717-1217929.html

上一篇:R语言使用ggplot2对冗余分析(RDA)等约束排序分析结果作图
下一篇:R语言PCA分析过程梳理及注意事项

1 俞洋

发表评论 评论 (0 个评论)

数据加载中...

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2021-10-16 23:30

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部