I am working on a community ecology analysis, and have used the metaMDS function from the vegan package to create a set of points to plot. I tried following this previous example Plotting ordiellipse function from vegan package onto NMDS plot created in ggplot2 However, in my case the for loop is not running correctly. These are my nmds scores:
gen_2019_site.scrs <- as.data.frame(scores(hel_gen_2019_nmds, display = "sites"))
gen_2019_site.scrs <- cbind(gen_2019_site.scrs, "Forest Coverage" = landscape_2019$Genera_2019_FQ)
gen_2019_site.scrs
NMDS1 | NMDS2 | NMDS3 | Forest Coverage |
---|---|---|---|
-0.01239620 | -0.35564430 | -0.3568236724 | F3 |
0.36384679 | -0.18408880 | -0.0618107711 | F1 |
-0.04411892 | 0.38898157 | 0.1333310873 | F5 |
-0.30058417 | 0.03967635 | 0.0041468720 | F2 |
0.03221559 | -0.08918320 | 0.3013787345 | F2 |
0.12416181 | -0.16631897 | 0.1091838635 | F1 |
-0.12732362 | 0.26747627 | -0.2164860869 | F1 |
0.02056514 | -0.25946756 | 0.2799632995 | F2 |
0.13859365 | 0.16943544 | -0.2485090677 | F3 |
0.18413252 | -0.13721912 | 0.1782989227 | F2 |
-0.25055749 | 0.30259844 | 0.0008892452 | F4 |
0.39622215 | 0.34536077 | -0.0881237834 | F1 |
0.18078468 | -0.10871097 | 0.0972200679 | F1 |
0.28217778 | -0.09857696 | -0.3025316708 | F4 |
-0.51108610 | -0.18696202 | -0.1046698976 | F2 |
0.19347030 | 0.09482813 | 0.1153315154 | F2 |
0.13955748 | 0.03051796 | 0.2791992925 | F2 |
-0.56016509 | -0.05928792 | 0.1780053536 | F4 |
-0.04069194 | 0.11660629 | -0.2203193485 | F3 |
-0.05838748 | 0.20457427 | 0.0616040447 | F2 |
0.05138256 | 0.01433274 | 0.0466866895 | F2 |
-0.20692179 | -0.39817145 | -0.1074178632 | F2 |
-0.09968087 | 0.04533097 | -0.2597924682 | F5 |
0.19142410 | -0.24604250 | 0.0643757014 | F2 |
-0.08662087 | 0.26995459 | 0.1168699402 | F5 |
and the code I tried based of the previous example was:
gen_2019_NMDS = data.frame(MDS1 = hel_gen_2019_nmds$points[,1], MDS2 = hel_gen_2019_nmds$points[,2],group=gen_2019_site.scrs$`Forest Coverage`)
gen_2019_NMDS.mean=aggregate(gen_2019_NMDS[,1:2],list(group=gen_2019_NMDS$group),mean)
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.gen2019.forest <- data.frame()
for(g in levels(gen_2019_NMDS$group)){
df_ell.gen2019.forest <- rbind(df_ell.gen2019.forest, cbind(as.data.frame(with(gen_2019_NMDS[gen_2019_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))
}
but when I run the for loop to try and get the data for graphing the ellipses, it runs with no errors or warnings however it remains an object with 0 obs of 0 variables. Ideally I would like to plot 5 ellipses based on the Forest Coverage levels like so.
ggplot(gen_2019_site.scrs, aes(x=NMDS1, y=NMDS2))+
geom_point(aes(NMDS1, NMDS2, colour = factor(`Forest Coverage`) ), size = 2)+
geom_path(data=df_ell.gen2019.forest, aes(x=NMDS1, y=NMDS2,colour=group), size=1, linetype=2)
But currently there's no data to plot the geom_path code. Any suggestions of where I have made a misstep in the coding would be greatly appreciated!