1

I am relatively new to R and I am trying to get my head around how to do ordination techniques in R, so that I don't need to use other software. I am trying to get a PCA with environmental factors in the place of species. As I have sites which differ qualitatively (in terms of land use) I wanted to be able to show that difference in the final plot (with different colours). Therefore, I used the method a la Gavin Simpson with the package vegan. So far so good. Here is also the code that I used for that:

with(fish, status)
scl <- -1 ## scaling = -1
colvec <- c("red2", "mediumblue")
plot(pond.pca, type = "n", scaling = scl)
with(fish, points(pond.pca, display = "sites", col = colvec[status], scaling = scl, pch = 21, bg = colvec[status]))
head(with(fish, colvec[status]))
text(pond.pca, display = "species", scaling = scl, cex = 0.8, col = "darkcyan")
with(fish, legend("topright", legend = levels(status), bty = "n", col = colvec, pch = 21, pt.bg = colvec))

The problem arises when I try to put arrows for my environmental variables in the ordination plot. If I use biplot and other functions like ordiplot etc. I ll not be able to keep the different colours for my two types of sites, therefore I don't want to use those. If I use the command here:

plot(envfit(pond.pca, PondEnv38, scaling=-1), add=TRUE, col="black")

I get nice arrows, only the are not aligned (and in some cases are completely opposite) with the environmental variables that I ve given with the code before (line 5). I tried to change the scaling but they just cannot align.

Does anyone know how to deal with that problem?

Any tips would be useful.

Alexander Vogt
  • 17,879
  • 13
  • 52
  • 68
Paris St
  • 33
  • 1
  • 2
  • 6
  • You may think you have implied that a particular package was being used but to someone with no experience in that field it does appear so. Also the data would appear needed, so I suggest you construct a small dataset that illustrates the problem and post it or link to it. – IRTFM Oct 20 '13 at 16:45
  • After waiting two hours, I'm voting to close for lack of clarity. – IRTFM Oct 20 '13 at 19:08
  • @DWin Hmm, the OP states they are using **vegan** and it has been tagged as such. I suppose an explicit `require` would have left no ambiguity, but 2 hours? Seriously? Not everyone hangs out here all the time. Premature. – Gavin Simpson Oct 21 '13 at 01:23
  • At no point in your code do you use any environmental data. Line 5 is just display the site scores coloured according to `status`; what does that have to do with `PondEnv38`? – Gavin Simpson Oct 21 '13 at 01:38
  • Sorry about the late reply, but I just came into my office. I mentioned that I used the package "vegan" although I have to admit I am not clear as would want to, as I am not giving an example dataset. I think the problem with not being able to get the arrows right was that I was using too many files and in the end things got messy. Anyway, I fixed the problem now thanks to the example given by Gavin Simpson below and promise to be more clear in future posts. – Paris St Oct 21 '13 at 09:41

1 Answers1

2

It is not clear what you are doing wrong as you don't provide a reproducible example of the problem and I am having difficulty following your description of what is wrong. Here is a fully worked out example for you to follow that does what you seem to being trying to do.

data(varespec)
data(varechem)

ord <- rda(varespec)

set.seed(1)
(fit <- envfit(ord, varechem, perm = 999))

## make up a fake `status`
status <- factor(rep(c("Class1","Class2"), times = nrow(varespec) / 2))

> head(status)
[1] Class1 Class2 Class1 Class2 Class1 Class2

Now plot

layout(matrix(1:2, ncol = 2))
## auto version
plot(fit, add = FALSE)

## manual version with extra things
colvec <-  c("red","green")
scl <- -1
plot(ord, type = "n", scaling = scl)
points(ord, display = "sites", col = colvec[status], pch = (1:2)[status])
points(ord, display = "species", pch = "+")
plot(fit, add = TRUE, col = "black")
layout(1)

Which gives

enter image description here

And all the arrows seem to be pointing as they would if you plotted the envfit object directly.

Gavin Simpson
  • 170,508
  • 25
  • 396
  • 453