11

通常のプロット関数の代わりに、ggplot2NMDS プロットを作成するために使用しています。パッケージの関数ordiellipse()を使用して、NMDS プロットにグループを表示したいと思います。vegan

サンプルデータ:

library(vegan)
library(ggplot2)
data(dune)
# calculate distance for NMDS
sol <- metaMDS(dune)
# 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)
# draw dispersion ellipses around data points
ordiellipse(sol, MyMeta$amt, display = "sites", kind = "sd", label = T)

# 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))

で作成した NMDS プロットに ordellipse を追加するにはどうすればよいggplot2ですか?

以下のDidzis Elfertsの答えはうまくいきます。ありがとうございました！ただし、次の ordiellipse を で作成された NMDS プロットにプロットすることに興味がありますggplot2。

ordiellipse(sol, MyMeta$amt, display = "sites", kind = "se", conf = 0.95, label = T)

残念ながら、関数がどのようにveganCovEllipse機能するかについて十分に理解していないため、自分でスクリプトを調整することができません。

4

2 に答える 2

13

まず、NMDSデータフレームに列グループを追加しました。

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

2番目のデータフレームには、各グループの平均MDS1およびMDS2値が含まれており、プロットにグループ名を表示するために使用されます。

  NMDS.mean=aggregate(NMDS[,1:2],list(group=group),mean)

データフレームdf_ellには、省略記号を表示するための値が含まれています。パッケージveganCovEllipseに隠されている関数で計算されます。veganこの関数は、NMDS（グループ）の各レベルに適用され、cov.wt共分散行列を計算するための関数も使用します。

  veganCovEllipse<-function (cov, center = c(0, 0), scale = 1, npoints = 100) 
  {
    theta <- (0:npoints) * 2 * pi/npoints
    Circle <- cbind(cos(theta), sin(theta))
    t(center + scale * t(Circle %*% chol(cov)))
  }

  df_ell <- data.frame()
  for(g in levels(NMDS$group)){
    df_ell <- rbind(df_ell, cbind(as.data.frame(with(NMDS[NMDS$group==g,],
                    veganCovEllipse(cov.wt(cbind(MDS1,MDS2),wt=rep(1/length(MDS1),length(MDS1)))$cov,center=c(mean(MDS1),mean(MDS2)))))
                    ,group=g))
  }

これで、楕円が関数でプロットgeom_path()さannotate()れ、グループ名をプロットするために使用されます。

  ggplot(data = NMDS, aes(MDS1, MDS2)) + geom_point(aes(color = group)) +
    geom_path(data=df_ell, aes(x=MDS1, y=MDS2,colour=group), size=1, linetype=2)+
    annotate("text",x=NMDS.mean$MDS1,y=NMDS.mean$MDS2,label=NMDS.mean$group)

楕円プロットのアイデアは、別のスタックオーバーフローの質問から採用されました。

ここに画像の説明を入力してください

UPDATE-両方の場合に機能するソリューション

まず、グループ列でNMDSデータフレームを作成します。

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

次に、関数の結果をordiellipse()オブジェクトとして保存します。

ord<-ordiellipse(sol, MyMeta$amt, display = "sites", 
                   kind = "se", conf = 0.95, label = T)

データフレームdf_ellには、省略記号を表示するための値が含まれています。パッケージveganCovEllipseに隠されている関数で再度計算されます。veganこの関数は、NMDS（グループ）の各レベルに適用され、ordオブジェクト- cov、centerおよびscale各レベルに格納されている引数を使用するようになりました。

df_ell <- data.frame()
for(g in levels(NMDS$group)){
  df_ell <- rbind(df_ell, cbind(as.data.frame(with(NMDS[NMDS$group==g,],
                  veganCovEllipse(ord[[g]]$cov,ord[[g]]$center,ord[[g]]$scale)))
                                ,group=g))
}

プロットは前の例と同じ方法で行われます。使用される楕円オブジェクトの座標の計算に関してordiellipse()は、このソリューションは、この関数に指定したさまざまなパラメーターで機能します。

ggplot(data = NMDS, aes(MDS1, MDS2)) + geom_point(aes(color = group)) +
  geom_path(data=df_ell, aes(x=NMDS1, y=NMDS2,colour=group), size=1, linetype=2)

ここに画像の説明を入力してください

于 2012-12-10T12:36:45.590 に答える