Search code examples
rggplot2pca

How to plot PCA using hellinger transformation in ggplot?


I'm trying to ggplot using Hellinger Transformation on my dataset. It works fine for a regular prcomp function but not Hellingers. How can I plot the data from Hellinger transformed data using ggplot?

library(ggfortify)
library(vegan)
df <- iris[1:4]
pca_res <- prcomp(df, scale. = TRUE)
autoplot(pca_res, data = iris, colour = 
'Species',
     loadings = TRUE, loadings.colour = 'blue',
     loadings.label = TRUE, loadings.label.size = 3)

##Hellinger Transformation
df.hell <- decostand(df, method = "hellinger")
df.hell <- rda(df.hell)

ggplot2::autoplot(df.hell)

autoplot(df.hell, data = iris, colour = 
 'Species',
 loadings = TRUE, loadings.colour = 'blue',
 loadings.label = TRUE, loadings.label.size = 3)

Error: Objects of type rda/cca not supported by autoplot.

Error: Objects of type rda/cca not supported by autoplot.

Edit 1: Even if the first plot can be manually computed in ggplot2, what about the rest of the plots like loading, or ellipses etc? base plot allows for overlay when using Hellingers but doesn't seem like ggplot2 would directly allow for it.


Solution

  • prcomp returns an object of class prcomp, which can be plotted with autoplot. As the error message says, rda function returns an object of class "rda" "cca", which cannot be plotted using autoplot. Therefore, you must extract the bits you need manually:

    data.frame(PC = df.hell$CA$u, species = iris$Species) %>% 
        ggplot(aes(x=PC.PC1, y=PC.PC2)) + 
        geom_point(aes(colour=species))
    

    You can find the relevant parts of the object by doing str(df.hell):

    List of 10
     $ colsum        : Named num [1:4] 0.037 0.0746 0.086 0.0854
      ..- attr(*, "names")= chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width"
     $ tot.chi       : num 0.0216
     $ Ybar          : num [1:150, 1:4] 0.0042 0.00511 0.0042 0.00359 0.00363 ...
      ..- attr(*, "scaled:center")= Named num [1:4] 0.656 0.479 0.498 0.267
      .. ..- attr(*, "names")= chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width"
      ..- attr(*, "dimnames")=List of 2
      .. ..$ : NULL
      .. ..$ : chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width"
      ..- attr(*, "METHOD")= chr "PCA"
     $ method        : chr "rda"
     $ call          : language rda(X = df.hell)
     $ pCCA          : NULL
     $ CCA           : NULL
     $ CA            :List of 7
      ..$ eig    : Named num [1:4] 0.0208691 0.0005348 0.0001951 0.0000205
      .. ..- attr(*, "names")= chr [1:4] "PC1" "PC2" "PC3" "PC4"
      ..$ poseig : NULL
      ..$ u      : num [1:150, 1:4] -0.122 -0.11 -0.119 -0.106 -0.123 ...
      .. ..- attr(*, "dimnames")=List of 2
      .. .. ..$ : NULL
      .. .. ..$ : chr [1:4] "PC1" "PC2" "PC3" "PC4"
      ..$ v      : num [1:4, 1:4] -0.241 -0.508 0.589 0.58 0.375 ...
      .. ..- attr(*, "dimnames")=List of 2
      .. .. ..$ : chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width"
      .. .. ..$ : chr [1:4] "PC1" "PC2" "PC3" "PC4"
      ..$ rank   : int 4
      ..$ tot.chi: num 0.0216
      ..$ Xbar   : num [1:150, 1:4] 0.0042 0.00511 0.0042 0.00359 0.00363 ...
      .. ..- attr(*, "scaled:center")= Named num [1:4] 0.656 0.479 0.498 0.267
      .. .. ..- attr(*, "names")= chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width"
      .. ..- attr(*, "dimnames")=List of 2
      .. .. ..$ : NULL
      .. .. ..$ : chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width"
      .. ..- attr(*, "METHOD")= chr "PCA"
     $ inertia       : chr "variance"
     $ regularization: chr "this is a vegan::rda result object"
     - attr(*, "class")= chr [1:2] "rda" "cca"
    

    enter image description here