[R] simultaneously maximizing two independent log likelihood functions using mle2

Ben Bolker bbolker at gmail.com
Mon Oct 17 14:47:20 CEST 2011


Adam Zeilinger <zeil0006 <at> umn.edu> writes:

> I have a log likelihood function that I was able to optimize using 
> mle2.  I have two years of the data used to fit the function and I would 
> like to fit both years simultaneously to test if the model parameter 
> estimates differ between years, using likelihood ratio tests and AIC.  
> Can anyone give advice on how to do this?
> 
> My likelihood functions are long so I'll use the tadpole predation 
> example from Ben Bolker's book, Ecological Data and Models in R (p. 
> 268-270).
> 
> library(emdbook)
> data(ReedfrogFuncresp)
> attach(ReedfrogFuncresp)
> # Holling Type II Equation
> holling2.pred = function(N0, a, h, P, T) {
>    a * N0 * P * T/(1 + a * h * N0)
> }
> # Negative log likelihood function
> NLL.holling2 = function(a, h, P = 1, T = 1) {
>    -sum(dbinom(Killed, prob = a * T * P/(1 + a * h * Initial),
>    size = Initial, log = TRUE))
> }
> # MLE statement
> FFR.holling2 = mle2(NLL.holling2, start = list(a = 0.012,
>    h = 0.84), data = list(T = 14, P = 3))
> 
> I have my negative log likelihood function setup similarly to the above 
> example.  Again, my goal is to simultaneously estimate parameters from 
> the same function for two years, such that I can test if the parameters 
> from the two years are different.  Perhaps an important difference from 
> the above example is that I am using a multinomial distribution (dmnom) 
> because my data are trinomially distributed.


  The most convenient way to do this would be to use the
'parameters' argument to mle2.  For example,

FFR.holling2.byyear <- mle2(NLL.holling2, start = list(a = 0.012,
    h = 0.84), data = c(list(T = 14, P = 3), ReedfrogFuncresp),  
    parameters=list(a~year, h~year))

-- or you could fit a model where only one or the other of the
parameters varied between years.
   R should automatically construct a set of parameters (in 
this case 'a for year 1', 'delta-a between year 1 and year 2',
'h for year 1', 'delta-h between year 1 and year 2'); it
sets the starting values of the 'delta' variables to zero
by default.  (Things get a little bit tricky if you also 
want to set lower and upper bounds.)  [See section 6.6.1 of
the book for more details.]

   Note that I have included ReedfrogFuncresp in the data list.
I now regret using attach() in the book examples.

  Ben Bolker



More information about the R-help mailing list