#DATA df1 = mtcars mycolors = df1$cyl #Identify the grouping vector library(vegan) m = metaMDS(df1) x = scores(m) #Extract co-ordinates plot(x, col = as.numeric(as.factor(mycolors)))