Bivariate Poisson lognormal {poilog}

Bivariate Poisson Lognormal Distribution

Description

Density and random generation for the for the bivariate Poisson lognormal distribution with parameters mu1, mu2, sig1, sig2 and rho.

Usage

dbipoilog(n1, n2, mu1, mu2, sig1, sig2, rho)
rbipoilog(S, mu1, mu2, sig1, sig2, rho, nu1=1, nu2=1, 
          condS=FALSE, keep0=FALSE)

Arguments

n1

vector of observed individuals for each species in sample 1

n2

vector of observed individuals for each species in sample 2
(in the same order as in sample 1)

mu1

mean of lognormal distribution in sample 1

mu2

mean of lognormal distribution in sample 1

sig1

standard deviation of lognormal distribution in sample 1

sig2

standard deviation of lognormal distribution in sample 2

rho

correlation among samples

S

Total number of species in both communities

nu1

sampling intensity sample 1

nu2

sampling intensity sample 2

condS

logical; if TRUE random deviates are conditional on S

keep0

logical; if TRUE species with count 0 in both communities are included in the random deviates

Details

The following is written from the perspective of using the Poisson lognormal distribution to describe community structure (the distribution of species when sampling individuals from a community of several species).

The following assumes knowledge of the Details section of dpoilog.

Consider two communities jointly and assume that the log abundances among species have the binormal distribution with parameters (mu1,sig1,mu2,sig2,rho). If sampling intensities are nu1 = nu2 = 1, samples from the communities will have the bivariate Poisson lognormal distribution

P(N_1 = \code{n1}, N_2 = \code{n2}; \;\code{mu1},\code{sig1},\code{mu2},\code{sig2},\code{rho}) =

q(\code{n1},\code{n2};\;\code{mu1},\code{sig1},\code{mu2},\code{sig2},\code{rho}) =

\int\limits_{-\infty}^{\infty}\int\limits_{-\infty}^{\infty} g_{\code{n1}}(\code{mu1,sig1},u) g_{\code{n2}}(\code{mu2,sig2},v)\,\phi(u,v;\code{rho})\,du dv,

where \phi(u,v; \code{rho}) here denotes the binormal distribution with zero means, unit variances and correlation rho. In the general case with sampling intensities nu1 and nu2, mu1 and mu2 should be replaced by mu1 + ln nu1 and mu2 + ln nu2 respectively. In this case, some species will be missing from both samples. The number of individuals for observed species will then have the truncated distribution

\frac{q(\code{n1},\code{n2};\;\code{mu1},\code{sig1},\code{mu2},\code{sig2},\code{rho})}{1-q(0,0;\;\code{mu1},\code{sig1},\code{mu2},\code{sig2},\code{rho})}

Value

dbipoilog returns the density
rbipoilog returns random deviates

Author(s)

Vidar Grotan [e-mail] and Steinar Engen

References

Engen, S., R. Lande, T. Walla and P. J. DeVries. 2002. Analyzing spatial structure of communities using the two-dimensional Poisson lognormal species abundance model. American Naturalist 160, 60-73.

See Also

bipoilogMLE for maximum likelihood estimation

Examples

### change in density of n2 for two different values of rho (given n1=10)   
barplot(rbind(dbipoilog(n1=rep(10,21),n2=0:20,mu1=0,mu2=0,sig=1,sig2=1,rho=0.0),
              dbipoilog(n1=rep(10,21),n2=0:20,mu1=0,mu2=0,sig=1,sig2=1,rho=0.8)),
              beside=TRUE,space=c(0,0.2),names.arg=0:20,xlab="n2",col=1:2)
legend(35,0.0012,c("rho=0","rho=0.8"),fill=1:2)


### draw random deviates from a community of 50 species 
rbipoilog(S=50,mu1=0,mu2=0,sig1=1,sig2=2,rho=0.7)

### draw random deviates including zeros
rbipoilog(S=50,mu1=0,mu2=0,sig1=1,sig2=2,rho=0.7,keep0=TRUE)

### draw random deviates with sampling intensities nu1=0.5 and nu2=0.7 
rbipoilog(S=50,mu1=0,mu2=0,sig1=1,sig2=2,rho=0.7,nu1=0.5,nu2=0.7)

### draw random deviates conditioned on a certain number of species 
rbipoilog(S=50,mu1=0,mu2=0,sig1=1,sig2=2,rho=0.7,nu1=0.5,nu2=0.7,condS=TRUE)


### how many species are likely to be observed in at least one of the samples
### (given S,mu1,mu2,sig1,sig2,rho)? 
hist(replicate(1000,nrow(rbipoilog(S=50,mu1=0,mu2=0,sig1=1,sig2=2,rho=0.7))),
     main="", xlab = "Number of species observed in at least one of the samples")

### how many individuals are likely to be observed 
### (given S,mu1,mu2,sig1,sig2,rho)? 
hist(replicate(1000,sum(rbipoilog(S=50,mu1=0,mu2=0,sig1=1,sig2=2,rho=0.7))),
     main="", xlab="sum nr of individuals in both samples")