I am trying to reproduce the column ("variable" in FactoMineR::PCA
, "species" in vegan::rda
) contribution percentages to axes from the FactoMineR
package in vegan
. The contributions are coded in FactoMiner::PCA
objects:
library(FactoMineR)
library(vegan)
data(dune)
fm <- FactoMineR::PCA(dune, scale.unit = FALSE, graph = FALSE)
head(round(sort(fm$var$contrib[,1], decreasing = TRUE), 3))
# Lolipere Agrostol Eleopalu Planlanc Poaprat Poatriv
# 17.990 16.020 13.866 7.088 6.861 4.850
By looking at the code for FactoMiner::PCA
, I found out that the contributions are calculated as squared axis coordinates divided by the axis eigenvalue and multiplied by 100%:
head(round(sort(100*fm$var$coord[,1]^2/fm$eig[1], decreasing = TRUE), 3))
# Lolipere Agrostol Eleopalu Planlanc Poaprat Poatriv
# 17.990 16.020 13.866 7.088 6.861 4.850
I fail to replicate the calculations above using a vegan::rda
object:
vg <- rda(dune)
head(round(sort(100*scores(vg, choices = 1, display = "sp",
scaling = 0)[,1]^2/vg$CA$eig[1], decreasing = TRUE), 3))
# Lolipere Agrostol Eleopalu Planlanc Poaprat Poatriv
# 0.726 0.646 0.559 0.286 0.277 0.196
I am obviously doing something wrong, and the discrepancy is probably due to the difference how these two packages calculate the coordinates for columns as eigenvalues for axes are quite similar (identical for my actual dataset), but the coordinates are not:
# vegan eigenvalue for axis 1
vg$CA$eig[1]
# PC1
# 24.79532
# FactoMineR eigenvalue for axis 1
fm$eig[1]
# [1] 23.55555
# vegan column coordinates for axis 1
head(round(scores(vg, choices = 1, display = "sp", scaling = 0)[,1], 3))
# Achimill Agrostol Airaprae Alopgeni Anthodor Bellpere
# -0.176 0.400 0.007 0.155 -0.163 -0.097
#FactoMineR, column coordinates for axis 1
head(round(fm$var$coord[,1], 3))
# Achimill Agrostol Airaprae Alopgeni Anthodor Bellpere
# 0.854 -1.943 -0.033 -0.751 0.791 0.472
# Sum of column coordinates for vegan axis 1 to illustrate the difference
sum(scores(vg, choices = 1, display = "sp", scaling = 0)[,1])
# [1] -0.796912
# Sum of column coordinates for FactoMineR axis 1 to illustrate the difference
sum(fm$var$coord[,1])
# [1] 3.867738
How do I calculate column/species percentage contributions to ordination axes using a vegan
rda
object?
vegan
's unscaled (scaling = 0
) column coordinates have equal sums of squares (i.e. 1 for each axis). You can get "column contributions" by simply squaring the unscaled coordinates:
head(sort(round(100*scores(vg, display = "sp", scaling = 0)[,1]^2, 3), decreasing = TRUE))
# Lolipere Agrostol Eleopalu Planlanc Poaprat Poatriv
# 17.990 16.020 13.866 7.088 6.861 4.850
Sum of these "contributions" equals 100% as said above:
sum(100*scores(vg, display = "sp", scaling = 0)[,1]^2)
# [1] 100