R help - Aug 2013 - bbmle: mle2: initial value in 'vmmin' is not finite

Dear R users

I am trying to estimate parameter values for function using optim in the mle2
function in the bbmle package but I cannot get past the above error message.

I am trying all the recommendations in the optim manual and the answers to this
question so far in the archive so I assume I have a data problem. But I would
like to know your thoughts.

I am using the 64 bit version of R3.0.1 on a PC running windows 7

There are two parameters I am interested in, "a" and "h"
The parameter value I am most interested in is "a"

I can estimate the values of these parameters using non linear regression (nls):

#INPUT
Killed=5
Time= c(29.28, 54.85, 46.22, 38.81, 21.35, 11.03, 28.69, 12.47, 43.28, 12.46,
38.81, 21.35, 81.03, 75.77, 13.01,  9.47,  6.51, 4.64, 12.59, 26.47)  # time
taken to kill first 5 prey
N0 = rep(c(0.25,0.5,0.75,1), each=5) # prey density/Litre
C=Killed/Time # capture rate
fr2<-data.frame(C, N0)
m1<-nls(C~a*N0/(1+(a*h*N0)), start = list(a = 0.2, h= 2), data=fr2)
summary(m1)

#OUTPUT
Formula: C ~ a * N0/(1 + (a * h * N0))
Parameters:
  Estimate Std. Error t value Pr(>|t|)
a   0.2656     0.1431   1.856    0.080 .
h  -1.9604     2.1447  -0.914    0.373
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1
Residual standard error: 0.2215 on 18 degrees of freedom
Number of iterations to convergence: 13
Achieved convergence tolerance: 6.933e-06


I know that prey depletion is a problem in this experiment as total numbers of
prey are limited so I want to use the rogers random predation equation, using
mle2 in the bbmle package to optimise the selection of the best parameters.
For this the data above are used again but in a slightly different format

# load packages
library(bbmle)
library(emdbook)

# specify the model
rogers.pred = function(N0, a, h, P, T) {
  N0 - lambertW(a * h * N0 * exp(-a * (P * T -
  h * N0)))/(a * h)
}

#write the liklihood function
NLL.rogers = function(a, h, T, P) {
  if (a < 0 || h < 0)
  return(NA)
  prop.exp = rogers.pred(Initial, a, h, P, T)/Initial
  -sum(dbinom(Killed, prob = prop.exp, size = Initial,
  log = TRUE))
}


# specfiy the data
P = 0.03125  # there is only one predator, this is the predator density/Litre
N0 = rep(c(0.25,0.5,0.75,1), each=5) # denisty of prey/Litre
Initial = N0
Time= c(29.28, 54.85, 46.22, 38.81, 21.35, 11.03, 28.69, 12.47, 43.28, 12.46,
38.81, 21.35, 81.03, 75.77, 13.01, 9.47, 6.51, 4.64, 12.59, 26.47) # time taken
to kill first 5 prey
Killed=5
#h=0.005      #  h = handling time and could be included here and fixed to
values (0-5) based on observations and the nls fit above





# place the data in a frame
fr2<-data.frame(Killed, P, N0, T, Initial)

# supply the liklihood function with starting values and run the model
FRR.rogers = mle2(NLL.rogers, start=list(a = 0.2, h=0.5), data=fr2)
summary(FRR.rogers)

#OUTPUT
Error in optim(par = 0.1, fn = function (p)  :     # I GET THIS ERROR ALWAYS
  initial value in 'vmmin' is not finite

In addition: Warning message:
In dbinom(Killed, prob = prop.exp, size = Initial, log = TRUE) :    # I GET THIS
ERROR MOST TIMES BUT NOT ALWAYS
  NaNs produced


Here I get this error message repeated no matter how many combinations of
starting values I use.
Realistic values or not.
Whether I include h as fixed or not.

The text above for the rogers equation runs perfectly with the data specified in
an example given by Bolker in emdbook
http://ms.mcmaster.ca/~bolker/emdbook/book.pdf


I would love to hear thoughts on what the problem might be?




Best wishes

Tom


Dr Tom C Cameron
email tom.cameron at emg.umu.se<mailto:tom.cameron at emg.umu.se>

Web site:
http://www.http://www.emg.umu.se/english/about-the-department/staff/cameron-thomas/<http://www.http//www.emg.umu.se/english/about-the-department/staff/cameron-thomas/>
Tel: 09 786 57 81 (ext 5781)
Ecology & Environmental Sciences
KBC-huset
Ume? universitet
SE-901 87 Ume?
Sverige

Thomas Cameron <tom.cameron <at> emg.umu.se> writes:
>  Dear R users 
> I am trying to estimate parameter values for function
> using optim in the mle2 function in the bbmle package but I cannot
> get past the above error message.
 
 [snip]
 
> # load packages
> library(bbmle)
> library(emdbook)
> 
> # specify the model
> rogers.pred = function(N0, a, h, P, T) {
>   N0 - lambertW(a * h * N0 * exp(-a * (P * T -
>   h * N0)))/(a * h)
> }
> 
> #write the liklihood function
> NLL.rogers = function(a, h, T, P) {
>   if (a < 0 || h < 0)
>   return(NA)
>   prop.exp = rogers.pred(Initial, a, h, P, T)/Initial
>   -sum(dbinom(Killed, prob = prop.exp, size = Initial,
>   log = TRUE))
> }
> 
> # specify the data
> P = 0.03125  # there is only one predator, this is the predator
density/Litre
> N0 = rep(c(0.25,0.5,0.75,1), each=5) # density of prey/Litre
> Initial = N0
> Time= c(29.28, 54.85, 46.22, 38.81, 21.35, 11.03, 
>  28.69, 12.47, 43.28, 12.46, 38.81, 21.35, 81.03,
> 75.77, 13.01, 9.47, 6.51, 4.64, 12.59, 26.47) 
> # time taken to kill first 5 prey
> Killed=5
> # place the data in a frame
> fr2<-data.frame(Killed, P, N0, T, Initial)
> 
> # supply the liklihood function with starting values and run the model
> FRR.rogers = mle2(NLL.rogers, start=list(a = 0.2, h=0.5), data=fr2)
> summary(FRR.rogers)
> 
> #OUTPUT
> Error in optim(par = 0.1, fn = function (p)  :     # I GET THIS ERROR
ALWAYS
>   initial value in 'vmmin' is not finite
 
  Short answer: you need N0 to be an integer (i.e. number of prey,
not density of prey per litre) if you're going to use it as the
'size'
parameter in a binomial likelihood.  (I haven't tried that yet to
see if that fixes everything, or if there are problems beyond that
step.)  
   I'm a little bit nervous about using this model for your
experimental design (time how long it takes for the first 5 prey
to be eaten -- the usual statement of the problem is number of
prey eaten in a *fixed* time, and the statistical properties
might be different).  If the predation probability were fixed this
would be Gamma distributed, but I'm not sure how it works out with
depletion ...

  Ben Bolker