13

我正在努力完成我在 vegan 和 ggplot2 中创建的 NMDS 图,但无法弄清楚如何将 envfit 物种加载向量添加到图中。当我尝试它说“无效的图形状态”。

下面的例子是从另一个问题(Plotting ordiellipse function from vegan package to NMDS plot created in ggplot2)稍微修改的,但它准确地表达了我想要包括的例子,因为我首先使用这个问题来帮助我将 metaMDS 放入 ggplot2 中:

library(vegan)
library(ggplot2)
data(dune)

# calculate distance for NMDS
NMDS.log<-log(dune+1)
sol <- metaMDS(NMDS.log)

# Create meta data for grouping
MyMeta = data.frame(
  sites = c(2,13,4,16,6,1,8,5,17,15,10,11,9,18,3,20,14,19,12,7),
  amt = c("hi", "hi", "hi", "md", "lo", "hi", "hi", "lo", "md", "md", "lo", 
      "lo", "hi", "lo", "hi", "md", "md", "lo", "hi", "lo"),
row.names = "sites")

# plot NMDS using basic plot function and color points by "amt" from MyMeta
plot(sol$points, col = MyMeta$amt)

# same in ggplot2
NMDS = data.frame(MDS1 = sol$points[,1], MDS2 = sol$points[,2])
ggplot(data = NMDS, aes(MDS1, MDS2)) + 
  geom_point(aes(data = MyMeta, color = MyMeta$amt))

#Add species loadings
vec.sp<-envfit(sol$points, NMDS.log, perm=1000)
plot(vec.sp, p.max=0.1, col="blue")
4

4 回答 4

11

(否则很好)接受的答案的问题,这解释了为什么向量在包含的图中都是相同的长度[请注意,接受的答案现在已被编辑,以按照我在下面描述的方式缩放箭头,以避免遇到 Q&A 的用户感到困惑的是,存储在$vectors$arrows返回的对象的组件中的envfit()是拟合向量的方向余弦。这些都是单位长度,因此@Didzis Elferts 图中的箭头都是相同的长度。这与来自的输出不同plot(envfit(sol, NMDS.log)),并且出现是因为我们通过与排序配置(“轴”)的相关性来缩放矢量箭头坐标。这样,与排序配置关系较弱的物种会得到更短的箭头。缩放是通过将方向余弦乘以打印输出表中显示的值来sqrt(r2)完成的。r2将向量添加到现有绘图时,vegan还尝试缩放向量集,以便它们填充可用绘图空间,同时保持箭头的相对长度。这是如何完成的在细节部分讨论?envfit并且需要使用未导出的函数vegan:::ordiArrowMul(result_of_envfit)

这是一个完整的工作示例,它复制了plot.envfit使用ggplot2的行为:

library(vegan)
library(ggplot2)
library(grid)
data(dune)

# calculate distance for NMDS
NMDS.log<-log1p(dune)
set.seed(42)
sol <- metaMDS(NMDS.log)

scrs <- as.data.frame(scores(sol, display = "sites"))
scrs <- cbind(scrs, Group = c("hi","hi","hi","md","lo","hi","hi","lo","md","md",
                              "lo","lo","hi","lo","hi","md","md","lo","hi","lo"))

set.seed(123)
vf <- envfit(sol, NMDS.log, perm = 999)

如果我们在这一点停下来看看vf

> vf

***VECTORS

             NMDS1       NMDS2     r2 Pr(>r)    
Belper -0.78061195 -0.62501598 0.1942  0.174    
Empnig -0.01315693  0.99991344 0.2501  0.054 .  
Junbuf  0.22941001 -0.97332987 0.1397  0.293    
Junart  0.99999981 -0.00062172 0.3647  0.022 *  
Airpra -0.20995196  0.97771170 0.5376  0.002 ** 
Elepal  0.98959723  0.14386566 0.6634  0.001 ***
Rumace -0.87985767 -0.47523728 0.0948  0.429
.... <truncated>

因此,r2数据用于缩放列NMDS1和中的值NMDS2。最终的情节是通过以下方式产生的:

spp.scrs <- as.data.frame(scores(vf, display = "vectors"))
spp.scrs <- cbind(spp.scrs, Species = rownames(spp.scrs))

p <- ggplot(scrs) +
  geom_point(mapping = aes(x = NMDS1, y = NMDS2, colour = Group)) +
  coord_fixed() + ## need aspect ratio of 1!
  geom_segment(data = spp.scrs,
               aes(x = 0, xend = NMDS1, y = 0, yend = NMDS2),
               arrow = arrow(length = unit(0.25, "cm")), colour = "grey") +
  geom_text(data = spp.scrs, aes(x = NMDS1, y = NMDS2, label = Species),
            size = 3)

这会产生:

在此处输入图像描述

于 2013-05-10T15:38:25.300 回答
7

从添加库开始。另外图书馆grid是必要的。

library(ggplot2)
library(vegan)
library(grid)
data(dune)

进行 metaMDS 分析并将结果保存在数据框中。

NMDS.log<-log(dune+1)
sol <- metaMDS(NMDS.log)

NMDS = data.frame(MDS1 = sol$points[,1], MDS2 = sol$points[,2])

添加物种载荷并将它们保存为数据框。箭头余弦的方向存储在 listvectors和 matrixarrows中。r2要获得箭头的坐标,这些方向值应乘以存储在 中的值的平方根vectors$r。更直接的方法是使用scores()@Gavin Simpson 的回答中提供的功能。然后添加包含species名称的新列。

vec.sp<-envfit(sol$points, NMDS.log, perm=1000)
vec.sp.df<-as.data.frame(vec.sp$vectors$arrows*sqrt(vec.sp$vectors$r))
vec.sp.df$species<-rownames(vec.sp.df)

箭头用 加上geom_segment(),物种名称用geom_text()。这两个任务都使用数据框vec.sp.df

ggplot(data = NMDS, aes(MDS1, MDS2)) + 
  geom_point(aes(data = MyMeta, color = MyMeta$amt))+
  geom_segment(data=vec.sp.df,aes(x=0,xend=MDS1,y=0,yend=MDS2),
      arrow = arrow(length = unit(0.5, "cm")),colour="grey",inherit_aes=FALSE) + 
  geom_text(data=vec.sp.df,aes(x=MDS1,y=MDS2,label=species),size=5)+
  coord_fixed()

在此处输入图像描述

于 2013-02-05T16:35:08.877 回答
3

我可以晚点加点东西吗?

Envfit 提供 pvalues,有时您只想绘制重要参数(纯素食者可以在 plot 命令中使用 p.=0.05 为您做一些事情)。我很难用 ggplot2 做到这一点。这是我的解决方案,也许你找到一个更优雅的解决方案?

从上面Didzis的回答开始:

ef<-envfit(sol$points, NMDS.log, perm=1000)
ef.df<-as.data.frame(ef$vectors$arrows*sqrt(ef$vectors$r))
ef.df$species<-rownames(ef.df)

#only significant pvalues
#shortcutting ef$vectors
A <- as.list(ef$vectors)
#creating the dataframe
pvals<-as.data.frame(A$pvals)
arrows<-as.data.frame(A$arrows*sqrt(A$r))
C<-cbind(arrows, pvals)
#subset
Cred<-subset(C,pvals<0.05)
Cred <- cbind(Cred, Species = rownames(Cred))

如上所述,“Cred”现在可以在 geom_segment-argument 中实现。

于 2014-08-21T11:31:31.080 回答
1

简短补充:为了在“箭头长度充分利用绘图区域”中获得plot.envfit功能的完整表示,ggplot2需要应用一个因素。我不知道上面的答案是否故意遗漏了它,因为加文甚至特别提到了它?只需使用提取所需的比例因子arrow_factor <- ordiArrowMul(vf),然后您可以将其应用于两个 NMDS 列,spp.scrs或者您可以手动执行此操作,如

arrow_factor <- ordiArrowMul(vf)
spp.scrs <- as.data.frame(scores(vf, display = "vectors")) * arrow_factor
spp.scrs <- cbind(spp.scrs, Species = rownames(spp.scrs), Pvalues = vf$vectors$pvals, R_squared = vf$vectors$r)

# select significance similarly to `plot(vf, p.max = 0.01)`
spp.scrs <- subset(spp.scrs, Pvalues < 0.01)

# you can also add the arrow factor in here (don't do both!)
ggplot(scrs) +
  geom_point(mapping = aes(x = NMDS1, y = NMDS2, colour = Group)) +
  coord_fixed() + ## need aspect ratio of 1!
  geom_segment(data = spp.scrs,
               aes(x = 0, xend = NMDS1 * arrow_factor, y = 0, yend = NMDS2 * arrow_factor),
               arrow = arrow(length = unit(0.25, "cm")), colour = "grey") +
  geom_text(data = spp.scrs, aes(x = NMDS1 * arrow_factor, y = NMDS2 * arrow_factor, label = Species),
            size = 3)
于 2020-05-20T12:00:30.487 回答