0

I have a community matrix (species as columns, samples as rows) from which I would like to generate a species accumulation curve (SAC) using the specaccum() and fitspecaccum() functions in R's vegan package. In order for the resulting SAC and cumulative species richness at sample X to be comparable among regions (I have 1 community matrix per region), I need to have specaccum() choose the same number of sets within each region. My problem is that some regions have a larger number of sets than others. I would like to limit the sample size to the minimum number of sets among regions (in my case, the minimum number of sets is 45, so I would like specaccum() to randomly sample 45 sets, 100 times (set permutations=100) for each region. I would like to sample from the entire data set available for each region. The code below has not worked... it doesn't recognize "subset=45". The vegan package info says "subset" needs to be logical... I don't understand how subset number can be logical, but maybe I am misinterpreting what subset is... Is there another way to do this? Would it be sufficient to run specaccum() for the entire number of sets available for each region and then just truncate the output to 45?

require(vegan)    
pool1<-specaccum(comm.matrix, gamma="jack1", method="random", subet=45, permutations=100)

Any help is much appreciated.

Christine S
  • 25
  • 6
  • Does it help if you write `subset` instead of `subet`? Anyway, the manual page says that `subset` is a `logical expression indicating sites (rows) to keep`. So if you have `n` "sets" and you want to pick 45 of these, you need a logical variable that is `TRUE` for 45 cases, and `FALSE` for others: these 45 with `TRUE` values will be you subset. To have a random mixing of the cases you need, use `sample` to permute the `subset` vector. Please note that `permutations` refers to the number of permutations of cases within one `subset`, i.e., it does not change the `subset` between permutations. – Jari Oksanen Sep 16 '17 at 10:03

1 Answers1

2

Why do you want to limit the function to work in a random sample of 45 cases? Just use the species accumulation up to 45 cases. Taking a random subset of 45 cases gives you the same accumulation, except for the random error of subsampling and throwing away information. If you want to compare your different cases, just compare them at the sample size that suits all cases, that is, at 45 or less. That is the idea of species accumulation models.

The subset is intended for situations where you have (possibly) heterogeneous collection of sampling units, and you want to stratify data. For instance, if you want to see only the species accumulation in the "OldLow" habitat type of the Barro Colorado data, you could do:

data(BCI, BCI.env)
plot(specaccum(BCI, subset = BCI.env$Habitat == "OldLow"))

If you want to have, say, a subset of 30 sample plots of the same data, you could do:

take <- c(rep(TRUE, 30), rep(FALSE, 20))
plot(specaccum(BCI)) # to see it all
# repeat the following to see how taking subset influences
lines(specaccum(BCI, subset = sample(take)), col = "blue")

If you repeat the last line, you see how taking random subset influences the results: the lines are normally within the error bars of all data, but differ from each other due to random error.

Jari Oksanen
  • 3,287
  • 1
  • 11
  • 15