\name{Sibuya}
\alias{rSibuya}
\alias{pSibuya}
\alias{dSibuya}
\alias{rSibuyaR}
\alias{dsumSibuya}
\title{Sibuya Distribution - Sampling and Probabilities}
\description{
The Sibuya distribution \eqn{\mathrm{Sib}(\alpha)}{Sib(alpha)} can be defined by
its Laplace transform
\deqn{1-(1-\exp(-t))^\alpha,\ t\in[0,\infty),}{1-(1-exp(-t))^alpha, t in [0,Inf),}
its distribution function% "R man - FIXME: allow to use \mathbb{N} instead:
\deqn{F(k)=1-(-1)^k{\alpha-1\choose k}=1-\frac{1}{kB(k,1-\alpha)},\ k\in\mathbf{N}}{%
F(k)=1-(-1)^k * choose(alpha-1,k)=1-1/(k * beta(k, 1-alpha)), k in N,}
(where \eqn{B}{beta} denotes the beta function) or its probability mass function
\deqn{p_k={\alpha\choose k}(-1)^{k-1},\ k\in\mathbf{N},
}{p_k = choose(alpha,k)*(-1)^(k-1), k in N,}
where \eqn{\alpha\in(0,1]}{alpha in (0,1]}.
\code{pSibuya} evaluates the distribution function.
\code{dSibuya} evaluates the probability mass function.
\code{rSibuya} generates random variates from
\eqn{\mathrm{Sib}(\alpha)}{Sib(alpha)} with
the algorithm described in Hofert (2011), Proposition 3.2.
\code{dsumSibuya} gives the probability mass function of the
\eqn{n}-fold convolution of Sibuya variables, that is, the sum of \eqn{n}
independent Sibuya random variables,
\eqn{S = \sum_{i=1}^n X_i}{S = sum(i=1..n) X[i]}, where
\eqn{X_i \sim \mathrm{Sib}(\alpha)}{X[i] ~ Sib(alpha)}.
This probability mass function can be shown (see Hofert
(2010, pp. 99)) to be
\deqn{\sum_{j=1}^n{n\choose j}{j\alpha\choose k} (-1)^{k-j},\ k\in\{n,n+1,\dots\}.}{%
sum(j=1..n) choose(n,j) * choose(j*alpha,k) * (-1)^(k-j), k >= n.}
}
\usage{
rSibuya(n, alpha)
dSibuya(x, alpha, log=FALSE)
pSibuya(x, alpha, lower.tail=TRUE, log.p=FALSE)
dsumSibuya(x, n, alpha,
method=c("log", "direct", "Rmpfr", "Rmpfr0", "RmpfrM", "diff", "exp.log"),
log=FALSE)
}
\arguments{
\item{n}{
for \code{rSibuya}: sample size, that is, length of the resulting
vector of random variates.
\cr
for \code{dsumSibuya}: the number \eqn{n} of summands.
}
\item{alpha}{parameter in \eqn{(0,1]}.}
\item{x}{vector of \code{\link{integer}} values (\dQuote{quantiles}) \eqn{x} at which to
compute the probability mass or cumulative probability.}
\item{log, log.p}{\code{\link{logical}}; if TRUE, probabilities p are given as log(p).}
\item{lower.tail}{\code{\link{logical}}; if TRUE (the default), probabilities
are \eqn{P(X \le x)}{P(X <= x)}, otherwise, \eqn{P(X > x)}.}
\item{method}{character string specifying which computational method is to be
applied. Implemented are:
\describe{
\item{\code{"log"}}{evaluates the logarithm of the sum
\deqn{\sum_{j=1}^n {n\choose j}{j\alpha\choose x} (-1)^{x-j}}{%
sum(j=1..n) choose(n,j) * choose(j*alpha,x) * (-1)^(x-j)}
in a numerically stable way;}
\item{\code{"direct"}}{directly evaluates the sum;}
\item{\code{"Rmpfr*"}}{are as \code{method="direct"} but use
high-precision arithmetic; \code{"Rmpfr"} returns
\code{\link{double}}s whereas \code{"RmpfrM"} gives
\code{\link[Rmpfr]{mpfr}} high-precision numbers, each trying to
adapt to high enough precision, whereas \code{"Rmpfr0"} does no
adaption.\cr
For all \code{"Rmpfr*"} methods, \code{alpha} can be set to a
\code{\link[Rmpfr:mpfr-class]{mpfr}} number of specified precision and this
will determine the precision of all parts of the internal
computations.
}
\item{\code{"diff"}}{interprets the sum as a forward difference
and computes it via \code{diff};}
\item{\code{"exp.log"}}{is as \code{method="log"} but without
numerically stable evaluation (not recommended, use with care).}
}
}
}
\details{
The Sibuya distribution has \bold{no} finite moments, that is, specifically
infinite mean and variance.
For documentation and didactical purposes, \code{rSibuyaR} is a pure-\R
implementation of \code{rSibuya}, of course slower than \code{rSibuya}
as the latter is implemented in C.
Note that the sum to evaluate for \code{dsumSibuya} is numerically
highly challenging, even already for small
\eqn{\alpha}{alpha} values (for example, \eqn{n \ge 10}{n >= 10}),
and therefore should be used with care. It may require high-precision
arithmetic which can be accessed with \code{method="Rmpfr"} (and the
\pkg{Rmpfr} package).
}
\value{
\describe{
\item{rSibuya:}{A vector of positive \code{\link{integer}}s of
length \code{n} containing the generated random variates.}
\item{dSibuya, pSibuya:}{a vector of
probabilities of the same length as \code{x}.}
\item{dsumSibuya:}{a vector of probabilities, positive if and only if
\code{x >= n} and of the same length as \code{x} (or \code{n} if
that is longer).}
}
}
\author{Marius Hofert, Martin Maechler}
\references{
Hofert, M. (2010).
\emph{Sampling Nested Archimedean Copulas with Applications to CDO Pricing}.
Suedwestdeutscher Verlag fuer Hochschulschriften AG & Co. KG.
Hofert, M. (2011).
Efficiently sampling nested Archimedean copulas.
\emph{Computational Statistics & Data Analysis} \bold{55}, 57--70.
}
\seealso{
\code{\link{rFJoe}} and \code{\link{rF01Joe}} (where \code{rSibuya} is applied).
}
\examples{
## Sample n random variates from a Sibuya(alpha) distribution and plot a
## histogram
n <- 1000
alpha <- .4
X <- rSibuya(n, alpha)
hist(log(X), prob=TRUE); lines(density(log(X)), col=2, lwd=2)
}
\keyword{distribution}