Original article
Computer generation of random variables with Lindley or Poisson–Lindley distribution via the Lambert W function

https://doi.org/10.1016/j.matcom.2010.09.006Get rights and content

Abstract

We provide procedures to generate random variables with Lindley distribution, and also with Poisson–Lindley or zero-truncated Poisson–Lindley distribution, as simple alternatives to the existing algorithms. Our procedures are based on the fact that the quantile functions of these probability distributions can be expressed in closed form in terms of the Lambert W function. As a consequence, the extreme order statistics from the above distributions can also be computer generated in a straightforward manner.

Introduction

In this note, we mainly deal with the computer generation of random data from the Lindley and Poisson–Lindley distributions. More specifically, we derive closed-form expressions for the quantile functions of both probability distributions, which are related to the Lambert W function. The particular case corresponding to the zero-truncated Poisson–Lindley distribution will also be considered. The explicit expressions obtained allow us to generate random samples from these probability distributions by means of the inverse transform method.

As it is well-known, the inverse transform method is the most simple procedure to generate random data from a given probability distribution T. This method is used when the quantile function of the random variable T can be expressed in closed form (cf. Fishman [9, pp. 149–156] and also Ross [17, Chapters 4 and 5]). We emphasize that the inverse transform method possesses some advantages over other generation methods. In particular, this method allow to generate samples corresponding to a truncated distribution obtained from T, which sometimes is called restricted sampling, and to this end the closed-form expression of the quantile function of T can be used directly (cf. Fishman [9, p. 152]). Another advantage is that the inverse transform method preserves the monotone relationship between the random variable T and the uniform distribution (cf. Fishman [9, p. 152]); this preservation property is the basis for the use of antithetic variables in order to induce negative correlation between two random variables with the same probability distribution (see also Rubinstein and Kroese [18, Chapter 5]). An additional advantage is that random samples corresponding to the extreme order statistics from T can also be generated in a straightforward manner (see Section 6 for further details). Furthermore, from a computational point of view, we have observed that our algorithms based on the inverse transform method to generate samples from the Lindley, Poisson–Lindley and zero-truncated Poisson–Lindley distributions are faster than the previous known algorithms based on the composition method (see Section 7, part (3)).

To be more precise, we first introduce a general definition of quantile function. Let T be an arbitrary random variable with cumulative distribution function FT(t) := P(T  t), tR, where R denotes the set of real numbers. For any u in the interval (0, 1), the uth quantile of T is usually defined as any solution t of the simultaneous inequalities:FT(t)uFT(t),where FT(t) stands for the left limit of FT at t. The solution of inequalities in Eq. (1) does not need to be unique. For example, if T is an integer-valued random variable then the solution of (1) is an interval on the real line for any fixed u  (0, 1). In order to avoid problems arising from the non-uniqueness, the most commonly accepted definition of quantile function of T is the following:QT(u):=inf{tR:FT(t)u},0<u<1.

In particular, the definition given by Eq. (2) implies that the uth quantile of an integer-valued random variable is also an integer. Moreover, from Eq. (2) it is clear that if FT is continuous and strictly increasing then FT has a unique inverse and thereby QT(u)=FT1(u), 0 < u < 1. We refer the reader to Parzen [16] for a review on the quantile function and its applications.

The Lindley distribution was introduced by the British statistician Dennis V. Lindley based on Bayes’ theorem (cf. Lindley [14], [15]). This continuous probability distribution has been largely overlooked in the literature until recently, when Ghitany et al. [11] have studied some of its statistical properties. These authors have shown that this random model can be a better model than the well-known exponential distribution in some particular cases. For simulation purposes, Ghitany et al. [11] consider that the inverse transform method cannot be used to generate data from the Lindley distribution by arguing that its quantile function has not a closed-form expression and then they use the composition method (cf. Fishman [9, pp. 160–165] for more details on this method). More precisely, their algorithm takes advantage of the property that the Lindley distribution can be represented as a mixture of the exponential and gamma distributions.

On the other hand, the Poisson–Lindley distribution is a discrete random variable introduced by Sankaran [19] to model count data, which is closely related to the Lindley distribution. In fact, it is a mixture distribution obtained from the Poisson distribution by assuming that the Poisson parameter is itself a random variable with Lindley distribution. Ghitany and Al-Mutairi [10], Borah and Begum [4] and Borah and Deka Nath [5] have studied some statistical properties of the Poisson–Lindley and other related probability distributions. In particular, Ghitany and Al-Mutairi [10] have developed an algorithm to generate random data from the Poisson–Lindley distribution based on the fact that this random variable is a mixed distribution. In addition, motivated by an application in molecular biology, Ghitany et al. [12] have also considered the so-called zero-truncated Poisson–Lindley distribution and an analogous algorithm to the one given in [10] is described in [12] to generate random samples from this truncated distribution.

In this note, we show that the quantile functions of the Lindley, Poisson–Lindley and zero-truncated Poisson–Lindley distributions can be expressed explicitly in terms of the Lambert W function, which completes the results given in Ghitany et al. [11], [12] and Ghitany and Al-Mutairi [10]. Therefore, the inverse transform method provides a straightforward procedure for the computer generation of the aforementioned random variables. The remainder of this note is organized as follows. In Section 2 we present a brief introduction to the Lambert W function. In Sections 3 The quantile function of the Lindley distribution, 4 The quantile function of the Poisson–Lindley distribution, 5 The quantile function of the zero-truncated Poisson–Lindley distribution, we give closed-form expressions for the quantile functions of the three probability distributions under consideration. Section 6 deals with the problem of the computer generation of random data from the extreme order statistics corresponding to these probability distributions. Finally, a few remarks on computational aspects are given in Section 7.

Section snippets

The Lambert W function

The Lambert W function has attracted a great deal of attention beginning with Lambert in 1758 and Euler in 1779. The name “Lambert W function” has become a standard after its implementation in the computer algebra system Maple in the 1980s and subsequent publication by Corless et al. [7] of a comprehensive survey of the history, theory and applications of this function. Since this seminal paper, a growing number of applications involving the Lambert W function can be found in the scientific

The quantile function of the Lindley distribution

The Lindley distribution with parameter θ, which we denote by X throughout this section, is specified by the following cumulative distribution function:FX(x):=1θ+1+θxθ+1exp(θx),x>0,θ>0.It may be noted that FX is continuous and strictly increasing so that the quantile function of X is QX(u)=FX1(u), 0 < u < 1. In the following result, we give an explicit expression for QX in terms of the Lambert W function.

Theorem 1

For any θ > 0, the quantile function of the Lindley distribution X is

QX(u)=11θ1θW1θ+1exp(

The quantile function of the Poisson–Lindley distribution

Let Y be a random variable having a discrete Poisson–Lindley distribution with parameter θ. Following the definition introduced by Sankaran [19], the probability mass function of Y is given by

P(Y=k):=θ2(θ+2+k)(θ+1)k+3,k=0,1,2,,θ>0,and the corresponding cumulative distribution function isFY(k)=1θ2+3θ+1+θk(θ+1)k+3,k=0,1,2,,θ>0.Let us denote by QY the quantile function of Y defined according to formula (2). Henceforth, log  ) denotes the natural logarithm and ⌈t⌉ stands for the ceiling of a

The quantile function of the zero-truncated Poisson–Lindley distribution

In this section, we establish an analogous result to Theorem 2 for the zero-truncated Poisson–Lindley distribution with parameter θ, which we denote by Z. The probability mass function of Z is the following

P(Z=k):=θ2θ2+3θ+1θ+2+k(θ+1)k,k=1,2,,θ>0,and the cumulative distribution function of Z is given by

FZ(k)=1θ2+3θ+1+θk(θ2+3θ+1)(θ+1)k,k=1,2,,θ>0.Let us denote by QZ the quantile function of Z defined according to formula (2). Then, we state the following.

Theorem 3

For any θ > 0, the quantile function of

Quantile functions of the extreme order statistics

In this final section, we show how the closed-form expressions provided in Theorems Theorem 1, Theorem 2, Theorem 3 also lead to explicit formulae for the quantile functions corresponding to the extreme order statistics from the Lindley, Poisson–Lindley and zero-truncated Poisson–Lindley distributions, as we sketch below. This is a significant advantage of our procedures to generate random data from these probability distributions in contrast to the algorithms proposed in Ghitany et al. [11],

Final remarks

  • (1)

    From Theorems Theorem 1, Theorem 2, Theorem 3, Monte Carlo simulation studies involving the Lindley, Poisson–Lindley and zero-truncated Poisson–Lindley distributions can be carried out in a straightforward manner by means of the inverse transform method and recalling that the Lambert W function is implemented in the most widespread computer algebra packages, as it was mentioned in Section 2.

  • (2)

    In general-purpose programming languages in which the Lambert W function is not usually implemented as a

Acknowledgments

The author expresses his sincere thanks to the anonymous referee for his/her comments, which led to an improvement of the paper. This work has been partially supported by Diputación General de Aragón (Grupo consolidado PDIE), by the research project MTM2009-14394-C02-01 and by FEDER funds.

References (19)

There are more references available in the full text version of this article.

Cited by (0)

View full text