I need to overlay the clusters generated from cutting a dendrogram a given level of similarity onto an ordination result (NMDS). I have been looking through ade4 and vegan without finding any obvious solutions to this problem.
I am current using Primer-e (see screen shot below) but am finding the graphics a bit limited. Any point in the right direction is greatly appreciated.
This is quite easy with vegan and I have a blog post that explains some of this in detail, but not the bit about the clustering.
Here is a quick example, I will assume you can translate this to whatever packages/code you are using.
Load the package and data set
require(vegan)
data(dune)
Compute the dissimilarity matrix and cluster it, cutting the dendrogram to give 3 groups
dij <- vegdist(dune) ## bray curtis dissimilarity
clu <- hclust(dij, method = "average")
grp <- cutree(clu, 3)
Look at grp
R> grp
2 13 4 16 6 1 8 5 17 15 10 11 9 18 3 20 14 19 12 7
1 1 1 2 1 1 1 1 3 2 1 1 1 1 1 2 2 3 1 1
and notice this now gives use the cluster membership (2nd row) for each sample (top row) in the data set.
Next fit the NMDS
set.seed(2) ## setting a seed to make this reproducible
ord <- metaMDS(dune)
In this example I will colour the points according to cluster membership so I need to define a vector of colours, one per cluster
col <- c("red2", "green4", "mediumblue")
I can now use grp and col to produce a vector of colour names for each point (sample) I plot, by indexing into col using grp. E.g.:
R> col[grp]
[1] "red2" "red2" "red2" "green4" "red2"
[6] "red2" "red2" "red2" "mediumblue" "green4"
[11] "red2" "red2" "red2" "red2" "red2"
[16] "green4" "green4" "mediumblue" "red2" "red2"
All the remains is to plot the NMDS ordintion and add the points and a legend. I suppress any plotting in the plot() call so I can have more control over adding the points in the next line. The third line just adds a legend.
plot(ord, type = "n", display = "sites")
points(ord, col = col[grp], bg = col[grp], pch = 21)
legend("topright", legend = paste("Cluster", 1:3),
col = col, pt.bg = col, bty = "n", pch = 21)
The resulting figure should look like this:
Update: To add convex hulls for each cluster of points to the ordination diagram you can make use of the ordihull() function. Continuing the example from above we add the convex hulls as follows
ordihull(ord, groups = grp, display = "sites")
At which point the figure would look like the one below
Note: vegan's higher-level plot() methods are designed purposely to give a quick and dirty display of the ordination and as such don;t accept vectors of colours or plotting characters. instead we expect you to build up your plots flow lower level methods such as the points() one I use here.
