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
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
