jfq3/ggordiplots

Adding species arrows to PCA

Closed this issue · 4 comments

Hi @jfq3,

I've made a PCA using the rda function of the vegan package. After taht, I decided to use the ggordiplots package to plot my PCA because I needed to plot standard error centroid at 95% confidence.

My code is the following:
myplot <- gg_ordiplot(pca.rda, groups = mydata$Species, kind="se", conf=0.95, pt.size =1)

I would like now to add the arrows that represents the "species" scores of the PCA.

Is that possible?
Thank you,

jfq3 commented

I am confused by your code. Groups shoud be some environmental variable, not species. And PCA should not be used for ordinating species data, at least not without some prior transformation. What, exactly, are you trying to do?

I am confused by your code. Groups shoud be some environmental variable, not species. And PCA should not be used for ordinating species data, at least not without some prior transformation. What, exactly, are you trying to do?

I'm sorry, I used the wrong term; with species I meant variables. Thank you

jfq3 commented
# Load libraries
library(vegan)
#> Loading required package: permute
#> Loading required package: lattice
#> This is vegan 2.5-7
library(ggplot2)
library(ggordiplots)
#> Loading required package: glue

# Load data
data("dune")
data("dune.env")

# Ordinate the species data; since PCA, transform first
dune.hel <- decostand(dune, method = "hellinger")
ord <- rda(dune.hel)

# Capture the ordination plot with ellipses for Management to a variable
plt.ellipse <- gg_ordiplot(ord, groups = dune.env$Management, ellipse = TRUE, plot = FALSE)

# plt.ellipse is a list with items named:
names(plt.ellipse)
#> [1] "df_ord"      "df_mean.ord" "df_ellipse"  "df_hull"     "df_spiders" 
#> [6] "plot"

# The first 5 items are dataframes for making plots

# Create a dataframe of environmental variables.
# If you have more than one, you can simply subset.
# If you have only one, to prevent the name being lost, do:
env <- data.frame(dune.env$A1)
colnames(env) <- "A1"
head(env)
#>    A1
#> 1 2.8
#> 2 3.5
#> 3 4.3
#> 4 4.2
#> 5 6.3
#> 6 4.3

# Capture the envfit plt to a variable
plt.env <- gg_envfit(ord, env = env, groups = dune.env$Management, plot = FALSE)

# Make labels for the ordination axes
ord.labels <- ord_labels(ord = ord)
ord.labels[1:2]
#> [1] "PC1 30.57%" "PC2 19.98%"

# Make the desired plot from data frames extracted from plt.ellipse and plt.env
plt <- ggplot() +
  geom_point(data=plt.ellipse$df_ord, aes(x=x, y=y, color=Group), size = 3) + # plot points
  geom_path(data = plt.ellipse$df_ellipse, aes(x=x, y=y, color=Group)) + # plot ellipses
  geom_segment(data=plt.env$df_arrows, aes(x=0, xend=x, y=0, yend=y),
               arrow=arrow(angle=20, length=unit(0.5, unit="cm")), color="red") + # plot arrow
  geom_text(data=plt.env$df_arrows, aes(x=x, y=y, label=var), color="red", hjust="outward") + # label arrow
  xlab(ord.labels[1]) + # Add axis labels
  ylab(ord.labels[2]) +
  guides(color=guide_legend(title="Management")) # Change legend title

plt

Created on 2022-04-08 by the reprex package (v2.0.1)

jfq3 commented

I wrote ggordiplots so that data frames from the individual types of plots could be recombined in ways like this. Browse the vignettes and code here on GitHub for other examples and methods.