I have prey consumption data on two species and I am trying to determine the degree to which their functional response curves differ.
Here is a small example of my data (one separate Excel file for each species, the non-native had 6 repetitions for each density instead of 5):
library(frair)
specie1 <- structure(list(density = c(2L, 2L, 2L, 2L, 2L, 4L, 4L, 4L, 4L,
4L, 8L, 8L, 8L, 8L, 8L, 12L, 12L, 12L, 12L, 12L, 16L, 16L, 16L,
16L, 16L, 20L, 20L, 20L, 20L, 20L, 24L, 24L, 24L, 24L, 24L, 28L,
28L, 28L, 28L, 28L, 32L, 32L, 32L, 32L, 32L, 42L, 42L, 42L, 42L,
42L, 54L, 54L, 54L, 54L, 54L), eaten = c(2L, 2L, 2L, 2L, 2L,
4L, 0L, 4L, 4L, 4L, 8L, 8L, 8L, 6L, 8L, 12L, 12L, 12L, 12L, 12L,
16L, 16L, 16L, 9L, 16L, 20L, 20L, 20L, 18L, 20L, 24L, 24L, 24L,
24L, 24L, 28L, 28L, 28L, 28L, 28L, 30L, 22L, 32L, 32L, 32L, 27L,
33L, 36L, 28L, 36L, 37L, 33L, 27L, 42L, 42L)), .Names = c("density",
"eaten"), class = "data.frame", row.names = c(NA, -55L))
specie2 <- structure(list(density = c(2L, 2L, 2L, 2L, 2L, 2L, 4L, 4L, 4L,
4L, 4L, 4L, 8L, 8L, 8L, 8L, 8L, 8L, 12L, 12L, 12L, 12L, 12L,
12L, 16L, 16L, 16L, 16L, 16L, 16L, 20L, 20L, 20L, 20L, 20L, 20L,
24L, 24L, 24L, 24L, 24L, 24L, 28L, 28L, 28L, 28L, 28L, 28L, 32L,
32L, 32L, 32L, 32L, 32L, 42L, 42L, 42L, 42L, 42L, 42L, 54L, 54L,
54L, 54L, 54L, 54L), eaten = c(0L, 0L, 2L, 0L, 0L, 0L, 0L, 0L,
4L, 0L, 3L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 0L, 7L, 0L, 4L, 4L, 4L,
9L, 15L, 8L, 15L, 16L, 0L, 0L, 1L, 4L, 0L, 0L, 13L, 2L, 2L, 4L,
18L, 5L, 10L, 4L, 7L, 14L, 13L, 0L, 5L, 18L, 0L, 18L, 20L, 15L,
21L, 5L, 9L, 15L, 5L, 1L, 2L, 31L, 20L, 7L, 13L, 10L, 1L)), .Names = c("density",
"eaten"), class = "data.frame", row.names = c(NA, -66L))
specie1_fit <- frair_fit(eaten~density, data=specie1, response='rogersII', start=list(a=1.2, h=0.015), fixed=list(T=24))
specie2_fit <- frair_fit(eaten~density, data=specie2, response='rogersII', start=list(a=1.2, h=0.015), fixed=list(T=24))
frair_compare(specie1_fit, specie2_fit)
#FUNCTIONAL RESPONSE COEFFICIENT TEST
#Response: rogersII
#Optimised variables: a,h
#Fixed variables: T
#Original coefficients:
# a h
#specie1_fit 0.22065 0.48170
#specie2_fit 0.01706 0.54308
#Test: specie1_fit - specie2_fit
# Estimate Std. Error z value Pr(z)
#Da 0.20631 NA NA NA
#Dh 0.03813 NA NA NA
#Warning message:
#In bbmle::mle2(minuslogl = fr_nll_difffunc, start = start, fixed = fixed, :
# couldn't invert Hessian
Then, I tried different setting different values for "a" and "h" (1.2 and 0.015 were default values) but the same error message came up.
I'm using R v3.4.2 on Windows 10.
