Skip to content

Annotated code for performing nMDS and perMANOVA

    Non-linear Multidimensional Scaling (nMDS)

    We have seen what a Principal Component Analysis does, how it works, and how to implement it in R. Non-linear Multidimensional Scaling (nMDS) is another useful ordination method for reducing the number of dimensions that describe the data. Regardless of the initial number of dimensions, nMDS disposes the data points in a two- or three-dimensional space until it finds the configuration that produces the minimal amount of distortion. Here is how to produce and represent nMDS ordinations in R:

    data(iris)
ii <- iris[,c(1:4)]
i.sp<-iris[,c(5)]
library(vegan)
i.dist <- vegdist(ii, method="euclidean")
i.mds <- metaMDS(i.dist, try=50)
plot(i.mds, display="sites",type="n",
	xlab="Axis 1", ylab="Axis 2",
	main="Iris species")
# add color
iris$Species
iris$col<-c(rep("red",50),rep("green",50),rep("blue",50))
points(i.mds, display="sites",
    pch=21, cex=1,
    bg=iris$col
    )
text(x=-2.5, y=3,paste("Stress=",round(i.mds$stress,4)))
legend("topright",
	legend=unique(iris$Species),
	pch=21,
	pt.bg=unique(iris$col)
	)

    The outcome is very similar to the one produced by PCA, and it suggests that the three Iris species under study are clearly discernible in the multidimensional space defined by sepal length, sepal width, petal length, and petal width. Let’s see how to test this pattern quantitatively.

    perMANOVA

    Permutation-based Multivariate Analysis of Variance, or PerMANOVA, is the multidimensional version of an Analysis of Variance. It answers the same question as an ANOVA does: does at least one study group differ significantly from the others? PerMANOVA can be performed in R with function adonis from package vegan.

    i.obs <- ii
i.spp <- data.frame(species=iris[,5])
head(i.obs)
head(i.spp)
d.manova <- adonis(i.obs ~ species, method = "euclidean", data= i.spp)
d.manova

# Call:
# adonis(formula = i.obs ~ species, data = i.spp, method = "euclidean") 

# Permutation: free
# Number of permutations: 999

# Terms added sequentially (first to last)

           # Df SumsOfSqs MeanSqs F.Model      R2 Pr(>F)    
# species     2    592.07 296.037  487.33 0.86894  0.001 ***
# Residuals 147     89.30   0.607         0.13106           
# Total     149    681.37                 1.00000           
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

    Post-hoc pairwise comparisons:

    if (!require(pairwiseAdonis)) install_github("pmartinezarbizu/pairwiseAdonis/pairwiseAdonis")
library(pairwiseAdonis)

#?pairwise.adonis
iris.adonis.pw <- pairwise.adonis(x = i.obs,
                    factors= i.spp$species,
                    sim.method='euclidean',
                    p.adjust.m='holm')
iris.adonis.pw <- as.data.frame(iris.adonis.pw)
iris.adonis.pw$F.Model <- round(iris.adonis.pw$F.Model, 2)
iris.adonis.pw$R2 <- round(iris.adonis.pw$R2, 2)
iris.adonis.pw

                    # pairs Df SumsOfSqs F.Model   R2 p.value p.adjusted sig
# 1    setosa vs versicolor  1  257.3267  551.00 0.85   0.001      0.003   *
# 2     setosa vs virginica  1  565.1335  943.80 0.91   0.001      0.003   *
# 3 versicolor vs virginica  1   65.6496   86.77 0.47   0.001      0.003   *

    The post-hoc pairwise comparison indicates that the three Iris species differ significantly from each other.

    Did you enjoy this? Consider joining my on-line course “First steps in data analysis with R” and learn data analysis from zero to hero!

    Exit mobile version