Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Overlaying clustering results on an ordination

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.

enter image description here

like image 595
Elizabeth Avatar asked Oct 06 '22 15:10

Elizabeth


1 Answers

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:

NMDS plot with cluster membership overlain


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

adding convex hulls for the clusters


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.

like image 122
Gavin Simpson Avatar answered Oct 10 '22 04:10

Gavin Simpson