Is there dpois equivalent for zero-inflated Poisson?

Thierry Onkelinx <thierry.onkelinx at inbo.be>
     on Tue, 22 Mar 2016 13:58:09 +0100 writes:

     > dpois(0, lambda) == e^(-lambda)
     > The wikipedia formula is
     > ifelse(x == 0, zero + dpois(0, lambda) * (1-zero), dpois(x, lambda) *
     > (1-zero))

     > or

     > ifelse(x == 0, zero + dpois(x, lambda) * (1-zero), dpois(x, lambda) *
     > (1-zero))

     > so we can move the dpois() out of the ifelse()

     > ifelse(x == 0, zero, 0)  + dpois(x, lambda) * (1-zero)

Nice!  Thank you, Thierry.

Even nicer for symmetry reasons (and much faster) is *not* to use ifelse(),
so we'd get

dzipois <- function(x, lambda, zero)
   (x == 0) * zero +  dpois(x, lambda) * (1-zero)

However, numerically correctly adding the  'log = FALSE'
argument which all good "density" functions have in R --
((and which you *should* use as  log = TRUE  for the
   log-likelihood if you want to become more professional))
is a bit tricky.

The x == 0 case at least needs care for large lambda and/or
small 'zero' (= pi).
All this is related on how to accurately compute  f(L) = log(1 - exp(- L))
which I call   log1mexp(L) in my still-not-submitted paper
available as one of the Rmpfr vignettes at
https://cloud.r-project.org/web/packages/Rmpfr/vignettes/log1mexp-note.pdf

BTW:  There are at least 3-4 R packages which deal with zero-inflated
       poisson in some ways, notably the recommended  'mgcv'
       package which comes bundled with every R.


Martin Maechler,
ETH ZUrich

     > ir. Thierry Onkelinx
     > Instituut voor natuur- en bosonderzoek / Research Institute for Nature and
     > Forest
     > team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
     > Kliniekstraat 25
     > 1070 Anderlecht
     > Belgium

     > To call in the statistician after the experiment is done may be no more
     > than asking him to perform a post-mortem examination: he may be able to say
     > what the experiment died of. ~ Sir Ronald Aylmer Fisher
     > The plural of anecdote is not data. ~ Roger Brinner
     > The combination of some data and an aching desire for an answer does not
     > ensure that a reasonable answer can be extracted from a given body of data.
     > ~ John Tukey

     > 2016-03-22 13:50 GMT+01:00 Matti Viljamaa <mviljamaa at kapsi.fi>:

     >> And why is the first term of ifelse(x == 0, zero, 0) + dpois(x, lambda) /
     >> (1 - zero)
     >>
     >> ifelse(x == 0, zero, 0)
     >>
     >> rather than something corresponding to
     >>
     >> zero+(1-zero)e^{-lambda}
     >>
     >> https://en.wikipedia.org/wiki/Zero-inflated_model#Zero-inflated_Poisson
     >>
     >> On 22 Mar 2016, at 14:25, Matti Viljamaa <mviljamaa at kapsi.fi> wrote:
     >>
     >> Could you clarify what are the parameters and why it?s formulated that way?
     >>
     >> -Matti
     >>
     >> On 22 Mar 2016, at 14:17, Thierry Onkelinx <thierry.onkelinx at inbo.be>
     >> wrote:
     >>
     >> Dear Matti,
     >>
     >> What about this?
     >>
     >> dzeroinflpois <- function(x, lambda, zero){
     >> ifelse(x == 0, zero, 0) + dpois(x, lambda) / (1 - zero)
     >> }
     >> plot(x, dzeroinflpois(x, lambda = 10, zero = 0.2), type = "l")
     >>
     >>
     >>
     >> ir. Thierry Onkelinx
     >> Instituut voor natuur- en bosonderzoek / Research Institute for Nature and
     >> Forest
     >> team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
     >> Kliniekstraat 25
     >> 1070 Anderlecht
     >> Belgium
     >>
     >> To call in the statistician after the experiment is done may be no more
     >> than asking him to perform a post-mortem examination: he may be able to say
     >> what the experiment died of. ~ Sir Ronald Aylmer Fisher
     >> The plural of anecdote is not data. ~ Roger Brinner
     >> The combination of some data and an aching desire for an answer does not
     >> ensure that a reasonable answer can be extracted from a given body of data.
     >> ~ John Tukey
     >>
     >> 2016-03-22 13:04 GMT+01:00 Matti Viljamaa <mviljamaa at kapsi.fi>:
     >>

     >>> I?m doing some optimisation that I first did with normal Poisson (only
     >>> parameter theta was estimated), but now I?m doing the same with a
     >>> zero-inflated Poisson model which
     >>> gives me two estimated parameters theta and p (p is also pi in some
     >>> notation).
     >>>
     >>> My question is, is there something equivalent to dpois that would use
     >>> both of the parameters (or is the p parameter possibly unnecessary)?
     >>>
     >>> I?m calculating the ?fit? of the Poisson model
     >>>
     >>> i.e. like
     >>>
     >>> x = c(0,1,2,3,4,5,6)
     >>> y = c(3062,587,284,103,33,4,2)
     >>> fit1 <- sum(y)*dpois(x, est_theta)
     >>>
     >>> and then comparing fit1 to the real observations.
     >>> [[alternative HTML version deleted]]
     >>>
     >>> ______________________________________________

     >>
     >>
     >>
     >>

     > [[alternative HTML version deleted]]

     > ______________________________________________