0

I am conducting an ordination analysis in R and I'm having trouble dealing with the results of the function metaMDS, from vegan. More specifically, I'm not able to get species scores from the metaMDS result. I found a similar question here How to get 'species score' for ordination with metaMDS()?, but the answer for it did not work for me.

My data is available here: https://drive.google.com/file/d/1btKzAWL_fmJ80GjcgMnwX5m_ls8h8vpY/view?usp=sharing

Here is the code I wrote so far:

library(vegan)

mydata <- read.table("nmdsdata.txt", h=T, row.names=1)

dist.f <- vegdist(mydata, method = "bray")
dist.f2  <- stepacross(dist.f, path = "extended")

results<-metaMDS(dist.f2, trymax=500)

I was told to use the function stepacross since the species composition of my sites is quite different and I was not getting converging results in metaMDS without using it. So far so good.

My problem arises when I try to plot the results:

plot(results, type="t")

When I run that line I receive the following message on my console: "species scores not available". I tried to follow the approach recommended on the link I cited earlier, by running the code:

NMDS<-metaMDS(dist.f2,k=2,trymax=500, distance = "bray")
envfit(NMDS, dist.f2)

However, it does not seem to work here. It returns me the sites scores and not the species scores as it does when using the data from the post I commented earlier. The only difference in my code is that I use "bray" and not "euclidean" distance. In addition, I get a warning after running envfit: "Warning message: In complete.cases(X) & complete.cases(env): longer object length is not a multiple of shorter object length".

I need both sites scores and species scores to plot the results. Am I missing something here? I'm new to R, so consider that, please,

Any help would be appreciated. Thanks in advance.

mad10les
  • 152
  • 2
  • 9

1 Answers1

1

If you want to have species scores, you must supply information on species. Dissimilarities do not have any information on species, and therefore they won't work. You have two obvious (and documented) alternatives:

  1. Supply a data matrix to metaMDS and it will calculate the dissimilarities, and if requested, also perform the stepacross:
## bray dissimilarities is the default, noshare triggers stepacross
NMDS <- metaMDS(mydata, noshare = TRUE, autotransform = FALSE, trymax = 500)

This also turns off automatic data transformation as you did not use it in your dissimilarities either.

  1. Add species scores to the result after the analysis when you did not supply information on species:
NMDS <- metaMDS(dist.f2, trymax = 500)
sppscores(NMDS) <- mydata

If you transformed community data, you should use similar transformation for mydata in sppscores.

Jari Oksanen
  • 3,287
  • 1
  • 11
  • 15
  • Thank you, @Jari Oksanen, that works perfectly! I have a question, though. When I plot the results from alternatives 1 and 2, the axes of the plots have different scales. That's probably because the dissimilarities were created by different functions, right? Is that something I should be concerned about? Does that interfere with the interpretations of the results? I'm sorry if I asked something stupid. I'm still a beginner. – mad10les Apr 22 '22 at 01:30
  • Firstly, scaling in NMDS is arbitrary – that is, you can multiply (or divide) results by a constant and the solution is just the same (*relative* differences do not change). However, metaMDS tries some tricks to fix the scaling, but that differs when you have access to the community data and when you don't. With community data you try "half-change" scaling: one axis unit means that community dissimilarity halves. With dissimilarities, the scaling is similar to dissimilarities. See section "5. Scaling of the results" in metaMDS help. – Jari Oksanen Apr 22 '22 at 07:34
  • It's been a few days since I asked the question, but I would enjoy your further help @Jari Oksanen. I followed your instructions and constructed the nmds of alternatives 1 and 2. However, the plots are different (and I'm not talking about the scale as in the previous comment): the positions of both sites and species differ between my two graphs. I used the same commands with the dune dataset and the two plots are identical, only varying in the scale. Any ideas why this is happening? I appreciate any help you can provide. – mad10les May 11 '22 at 15:09
  • In case it helps somebody with the same problem, I think I've found the answer. In alternative 1, the stepacross triggered by noshare has the argument path = "shortest". On the other hand, in alternative 2 the matrix of dissimilarities was created with the argument path of stepacross set to "extended". That caused the differences between the graphs. Somebody correct me if I'm wrong. – mad10les May 11 '22 at 16:02