R falsely identifying abundance data as presence/absence data (fossil package)

Question

I am attempting to use the spp.est function of the package called "fossil" in RStudio. I have created a matrix called "akimiskibb" of abundance data with species as the columns and sites as the rows. When I try to use the function spp.est, I type this:

spp.est(akimiskibb, rand = 10, abund = TRUE, counter = FALSE, max.est = 'all')

The problem comes in because my abundance data has a lot of zeroes, so I get this error message:

Error in if (max(x) == 1) warning("cannot use incidence data for abundance-based analyses. If the data is incidence based, please run this function again with the option of 'abund=FALSE'") : missing value where TRUE/FALSE needed

This function has worked in the past with matrices with a lot of zeroes (which are also abundance data, not presence/absence). I don't know what I am doing wrong.

Has anyone experienced something similar and found a way around this?

Thank you,

Kayla

Data: matrix format:

 *sp1 sp2 sp3 sp4 sp5 sp6 sp7 sp8 sp9 sp10 sp11 sp12 sp13 sp14 sp15 sp16 sp17
 sample1    0   0   0   0   0   0   0   1   0    0    0    0    0    0    1    
 0    0
 sample 2   0   0   0   1   0   0   1  25   7    0   18   12    0    0    0    
 1    1
 sample3    0   0   0   0   0   0   0   3   0    0    3    1    0    0    0    
 5    4
 sp18 sp19 sp20 sp21 sp22 sp23 sp24 sp25 sp26 sp27 sp28 sp29 sp30 sp31 
 sp32
 sample1     0    0    0    0    0    0    0    0    0    0    0    0    0    
 0    0
 sample 2    1    0    1    0    0    0    0    0    0    0    0    3    2    
 0    3
 sample3     0    0    1    0    0   11    0    0    0    0    0    0    0    
 0    1
 sp33 sp34 sp35 sp36 sp37 sp38 sp39 sp40 sp41 sp42 sp43  X
 sample1     0    0    0    0    0    0    0    0    0    0    0 NA
 sample 2    0    0    3    2    1    0    0    1    8    0    0 NA
 sample3     0    0    0    0    0    0    0    0    0    0    0 NA*

dput:

 *structure(list(sp1 = c(0L, 0L, 0L), sp2 = c(0L, 0L, 0L), sp3 = c(0L, 
 0L, 0L), sp4 = c(0L, 1L, 0L), sp5 = c(0L, 0L, 0L), sp6 = c(0L, 
 0L, 0L), sp7 = c(0L, 1L, 0L), sp8 = c(1L, 25L, 3L), sp9 = c(0L, 
 7L, 0L), sp10 = c(0L, 0L, 0L), sp11 = c(0L, 18L, 3L), sp12 = c(0L, 
 12L, 1L), sp13 = c(0L, 0L, 0L), sp14 = c(0L, 0L, 0L), sp15 = c(1L, 
 0L, 0L), sp16 = c(0L, 1L, 5L), sp17 = c(0L, 1L, 4L), sp18 = c(0L, 
 1L, 0L), sp19 = c(0L, 0L, 0L), sp20 = c(0L, 1L, 1L), sp21 = c(0L, 
 0L, 0L), sp22 = c(0L, 0L, 0L), sp23 = c(0L, 0L, 11L), sp24 = c(0L, 
 0L, 0L), sp25 = c(0L, 0L, 0L), sp26 = c(0L, 0L, 0L), sp27 = c(0L, 
 0L, 0L), sp28 = c(0L, 0L, 0L), sp29 = c(0L, 3L, 0L), sp30 = c(0L, 
 2L, 0L), sp31 = c(0L, 0L, 0L), sp32 = c(0L, 3L, 1L), sp33 = c(0L, 
 0L, 0L), sp34 = c(0L, 0L, 0L), sp35 = c(0L, 3L, 0L), sp36 = c(0L, 
 2L, 0L), sp37 = c(0L, 1L, 0L), sp38 = c(0L, 0L, 0L), sp39 = c(0L, 
 0L, 0L), sp40 = c(0L, 1L, 0L), sp41 = c(0L, 8L, 0L), sp42 = c(0L, 
 0L, 0L), sp43 = c(0L, 0L, 0L), X = c(NA, NA, NA)), .Names = c("sp1", 
 "sp2", "sp3", "sp4", "sp5", "sp6", "sp7", "sp8", "sp9", "sp10", 
 "sp11", "sp12", "sp13", "sp14", "sp15", "sp16", "sp17", "sp18", 
 "sp19", "sp20", "sp21", "sp22", "sp23", "sp24", "sp25", "sp26", 
 "sp27", "sp28", "sp29", "sp30", "sp31", "sp32", "sp33", "sp34", 
 "sp35", "sp36", "sp37", "sp38", "sp39", "sp40", "sp41", "sp42", 
 "sp43", "X"), class = "data.frame", row.names = c("sample1", 
 "sample 2", "sample3"))*

packages used:

fossil (made in R version 3.4.4)

Version of R: R x64 3.4.1

Answer 1

OK, not familiar with this package so thank you for introducing it to me. Always great to see how R is being used in so many disciplines.

I see two issues with your applications. First, according to the fossil documentation for spp.est, your data needs to have your samples as columns and your species as rows. The second issue is the "species" X with the NA values. You need to get rid of these because the function can't handle them.

Here's the code:

library(fossil)
library(tidyverse)

akimiskibb <- as.data.frame(t(akimiskibb)) %>%
 filter(!is.na(sample1) == T)

spp.est(akimiskibb, rand = 10, abund = TRUE, counter = FALSE, max.est = 'all')

# which results in 

     N.obs S.obs S.obs(+95%) S.obs(-95%) Chao1 Chao1(upper) Chao1(lower)
[1,]     1  11.8    25.30995   -1.709948 20.90     28.38331     13.41669
[2,]     2  16.0    25.46770    6.532301 27.15     33.37680     20.92320
[3,]     3  20.0    20.00000   20.000000 26.00     29.24037     22.75963
          ACE ACE(upper) ACE(lower)    Jack1 Jack1(upper) Jack1(lower)
[1,] 14.56286   31.16289  -2.037175 16.82501     36.43129    -2.781268
[2,] 19.19286   30.64000   7.745713 22.21008     35.28448     9.135677
[3,] 22.28571   22.28571  22.285714 25.95082     25.95082    25.950820
attr(,"data.type")
[1] "abundance"
Warning messages:
1: In ests[[k]](b) :
  This data appears to be presence/absence based, but this estimator is for abundance data only
2: In ests[[k]](b) :
  This data appears to be presence/absence based, but this estimator is for abundance data only