Polono {VGAM} | R Documentation |
Density, and random generation for the Poisson lognormal distribution.
dpolono(x, meanlog=0, sdlog=1, bigx=Inf, ...) rpolono(n, meanlog=0, sdlog=1)
x |
vector of quantiles. |
n |
number of observations. Must be a positive integer of length 1. |
meanlog, sdlog |
the mean and standard deviation of the normal distribution
(on the log scale).
They match the arguments in
Lognormal .
|
bigx |
Numeric.
This argument is for handling large values of x and/or
when integrate fails.
A first order Taylor series approximation
[Equation (7) of Bulmer (1974)]
is used at values of x that are greater or equal to this argument.
For bigx=10 ,
he showed that the approximation has a relative error less than
0.001 for values of meanlog and
sdlog ``likely to be encountered in practice''. The default value
means that this approximation is not used. Setting something like
bigx=100 may be a good idea.
|
... |
Arguments passed into
integrate .
|
The Poisson lognormal distribution is similar to the negative binomial in that it can be motivated by a Poisson distribution whose mean parameter comes from a right skewed distribution (gamma for the negative binomial and lognormal for the Poisson lognormal distribution).
dpolono
gives the density, and
rpolono
generates random deviates.
By default,
dpolono
involves numerical integration that is performed using
integrate
. Consequently, computations are very
slow and numerical problems may occur
(if so then the use of ...
may be needed).
Alternatively, for extreme values of x
, meanlog
,
sdlog
, etc., the use of bigx
avoids the call to
integrate
; however the answer may be a little
inaccurate.
For the maximum likelihood estimation of the 2 parameters a VGAM
family function called polono
, say, has not been written yet.
T. W. Yee
Bulmer, M. G. (1974) On fitting the Poisson lognormal distribution to species-abundance data. Biometrics, 30, 101–110.
lognormal
,
poissonff
,
negbinomial
.
meanlog = 0.5; sdlog = 0.5 y = 0:19 proby = dpolono(y, m=meanlog, sd=sdlog) sum(proby) # Should be 1 ## Not run: opar = par(no.readonly = TRUE) par(mfrow=c(2,2)) plot(y, proby, type="h", col="blue", ylab="P[Y=y]", log="", main=paste("Poisson lognormal(meanlog=",meanlog,", sdlog=",sdlog,")", sep="")) # More extreme values; use the approximation and plot on a log scale # Notice the kink at bigx. y = 0:190 proby = dpolono(y, m=meanlog, sd=sdlog, bigx=100) sum(proby) # Should be 1 plot(y, proby, type="h", col="blue", ylab="P[Y=y]", log="y", main=paste("Poisson lognormal(meanlog=",meanlog,", sdlog=",sdlog,")")) # Random number generation table(y <- rpolono(n=1000, m=meanlog, sd=sdlog)) hist(y, breaks=((-1):max(y))+0.5, prob=TRUE, border="blue") par(opar) ## End(Not run)