Search code examples
rpredictionveganaccumulate

vegan accumulation curve predictions


I'm having a matrix of plants(rows) and pollinators(columns) and interaction frequencies within (converted to 0 (no interaction) and 1 (interaction/s present) for this analysis). I'm using the vegan package and have produced a species accumulation curve.

accum <- specaccum(mydata[1:47,], method = "random", permutations = 1000)

plot(accum)

I now would like to predict how many new pollinator species I would be likely to find with additional plant sampling but can't figure in what format I have to include "newdata" within the predict command. I have tried empty rows and rows with zeros within the matrix but was not able to get results. This is the code I've used for the prediction:

predictaccum1 <- predict(accum, newdata=mydata[48:94,])

The error message:

Error in UseMethod("predict") : 
no applicable method for 'predict' applied to an object of class "specaccum"

The error message does not change if I specify: interpolation = c("linear") or "spline".

Could anyone help please?


Solution

  • Not perhaps the clearest way of putting this, but the documentation says:

    newdata: Optional data used in prediction interpreted as number of
              sampling units (sites).
    

    It should be a number of sampling units you had. A single number or a vector of numbers will do. However, the predict function cannot extrapolate, but it only interpolates. The nonlinear regression models of fitspecaccum may be able to extrapolate, but should you trust them?

    Here a bit about dangers of extrapolation: the non-linear regression models are conventionally used analysing species accumulation data, but none of these is really firmly based on theory -- they are just some nice non-linear regression models. I know of some models that may have a firmer basis, but we haven't implemented them in vegan, neither plan to do so (but contributions are welcome). However, it is possible to get some idea of problems by subsampling your data and seeing if you can estimate the overall number of species with an extrapolation from your subsample. The following shows how to do this with the BCI data in vegan. These data have 50 sample plot with 225 species. We take subsamples of 25 plots and extrapolate to 50:

    mod <- c("arrhenius", "gleason", "gitay", "lomolino", "asymp", "gompertz", 
          "michaelis-menten", "logis", "weibull")
    extraps <- matrix(NA, 100, length(mod))
    colnames(extraps) <- mod
    for(i in 1:nrow(extraps)) {
       ## use the same accumulation for all nls models 
       m <- specaccum(BCI[sample(50,25),], "exact") 
       for(p in mod) { 
          ## need try because some nls models can fail
          tmp <-  try(predict(fitspecaccum(m, p), newdata=50))
          if(!inherits(tmp, "try-error")) extraps[i,p] <- tmp
       } 
    }
    

    When I tried this, most extrapolation models did not include the correct number of species among their predictions, but all values were either higher than correct richness (from worst: Arrhenius, Gitay, Gleason) or lower than correct richness (from worst: logistic, Gompertz, asymptotic, Michaelis-Menten, Lomolino, Weibull; only these two last included the correct richness in their range).

    In summary: in lack of theory and adequate model, beware extrapolation.