[R] accessing log likelihood of poison model
Peter Dalgaard
p.dalgaard at biostat.ku.dk
Sun Apr 18 10:38:36 CEST 2004
Renaud Lancelot <renaud.lancelot at cirad.fr> writes:
> >>>And that's great; but I need the log likelihood.
> >>>
> >>>Anyone know?
> >>
> >> The deviance will not suffice? sum(dpois(nfix, fitted(freq.mod),
> >> log.p=T))
> >>
> >>should do the trick otherwise.
> > Thank you, that did the trick. I should note that the method
> > doesn't work when I have a noninteger in Y. I'm doing a "log
> > linear" analysis
> > of a cross table, and one of the cells of the cross table has a zero
> > value. I put a .5 in that cell. With the ".5 analysis" the above
> > method doesn't work. So I had to revert to an analysis with a zero
> > cell to make it work. Which has some problems. So I'm still left
> > looking for the log likelihood of the ".5 analysis". I've had success
> > with Stata using this method, but not R.
> > Thanks again.
>
> ?logLik
Right. Forgot about that. It's the same sum(dpois(....,log=T))
calculation, though (after a detour around the AIC), so that too will
generate -Inf if noninteger values are plugged in. Fair enough: The
likelihood is by definition the probability of obtaining the data, and
the Poisson distribution has probability zero of generating fractional
values.
However, the deviance is well-defined in such cases. I'd suspect that
what Stata outputs is not the true log-likelihood, but a number
obtained by "analytic continuation" (wrong term, I know) of the
Poisson density to nonintegers, e.g.
> p <- function(x,lambda) lambda^x*exp(-lambda)/gamma(x+1)
> p(1,1)
[1] 0.3678794
> dpois(1,1)
[1] 0.3678794
> dpois(.5,1)
[1] 0
Warning message:
non-integer x = 0.500000
> p(.5,1)
[1] 0.4151075
--
O__ ---- Peter Dalgaard Blegdamsvej 3
c/ /'_ --- Dept. of Biostatistics 2200 Cph. N
(*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk) FAX: (+45) 35327907
More information about the R-help
mailing list