A New Generalized Poisson-Lindley Distribution : Applications and Properties

A new generalized Poisson-Lindley distribution is obtained by compounding Poisson distribution with a two parameter generalised Lindley distribution. The new distribution is shown to be unimodal and over–dispersed. This distribution has a tendency to accommodate right tail, and for particular values of parameter, the tail tends to zero at a faster rate. Various properties like cumulative distribution function, generating function, moments etc. are derived. Knowledge about the parameters is obtained through method of moments, maximum likelihood method and EM algorithm. Moreover, an actuarial application to collective risk model is shown by considering the proposed distribution as primary distribution and exponential and Erlang as secondary distributions. The model is applied to real data sets and found to perform better than other competing models.


Introduction
Researchers in the recent past proposed many distributions, most popular being the generalizations of exponential, gamma and Weibull distributions.Lindley (1958) suggested a one parameter distribution to illustrate the difference between fiducial distribution and posterior distribution.Ghitany et al. (2008) studied many statistical properties of the Lindley distribution and showed that the mathematical properties are more flexible than exponential distribution, as such, Lindley distribution is found to be better model than exponential distribution.Many researchers viz.Zakerzadeh and Dolati (2009), Zamani and Ismail (2010), Nadarajah (2011), Bakouch et al. (2012), Elbatal et al. (2013), Hassan (2014), Shanker (2013) have proposed new classes of distributions by modifying the Lindley distribution and discussed various properties of their proposed generalizations.Sankaran (1970) while modelling count data, introduced Poisson-Lindley distribution, by assuming the parameter of the Poisson distribution, λ, to follow a Lindley distribution.Ghitany et al. (2008) proposed different methods of estimation for the discrete Poisson-Lindley distribution.Ghitany et al. (2008) have used the zero truncated Poisson-Lindley distribution to model count data when the data structurally excludes zero counts.Size-biased Poisson-Lindley Distribution and its application were also studied by Ghitany et al. (2008).Mahmoudi et al.(2010) generalized the Poisson-Lindley distribution of Sankaran (1970) and showed that their generalized distribution has more flexibility in analysing count data.Hernández et al. (2010) in the context of insurance have used the Poisson-Lindley to study the collective risk model when number of claims and size of a single claim are incorporated into the model.Another form of discrete Lindley distribution was introduced by Gómez et al. (2011) by discretizing the continuous Lindley distribution.The model is over-dispersed and competitive with the Poisson distribution to fit automobile claim frequency data sets.After revising some of its properties a compound discrete Lindley distribution is obtained.Gómez et al. (2013) analysed the collective risk model by assuming Erlang distribution to the loss data and generalised Lindley distribution to the claim frequency data.In many potential branches like insurance, clinical trials, engineering, biology etc. the variable of interest is a count variable.Though classical models like Poisson, geometric, negative binomial, and their generalizations (see for example Philippou (1983), Gómez (2010Gómez ( , 2011)) 2014)), are available for count data analysis, it is found that these models are not supportive to capture the right tail behaviour of the data set properly.Few attempts have been made by researchers to capture this right tail behaviour (see for example Gomez et al. ( 2013)), but not much work has been done.The importance of right tail behaviour is illustrated in the following example.In insurance industry, the actuary while fixing the premium, assumes a claim distribution which may or may not be realised after a claim is made, therefore it becomes a challenge to the actuary to develop new statistical distributions which take care of this subtle difference.Besides properly fitting to the available data set, the proposed distribution has to capture specific characteristics like shape and tail behaviour.The tail can either tapers out quickly or slowly.When the tail approaches zero slowly, the extreme values can be properly understood.Therefore, besides fitting properly to the data, the shape of the distribution should capture the behaviour of extreme values also with their tendency to approach zero.This critical observation has motivated the authors to search for a discrete distribution which fades away to zero much more slowly/faster than the classical compound Poisson-Lindley distribution by Sankaran (1970).
In this paper we propose a new generalized Poisson-Lindley distribution which is obtained from Poisson distribution when its parameter λ, follows a two parameter Lindley distribution (T PLD(θ, α)) as suggested by Shanker et al. (2013) defined as The paper is structured as follows: In Section 2, a new distribution is proposed.Various properties like shape, quantile function, generating function and moments are presented in Section 3.An algorithm for simulation of random variable is shown in Section 4. Estimation methods like method of moments, maximum likelihood estimation and EM algorithm are discussed in Section 5. Section 6 addresses the actuarial application to auto mobile insurance of the proposed model.In Section 7, applicability of the proposed model is shown, and compared with other competing probability models.Moreover, by using Bayesian methodology, Bonus-Malus premium for automobile insurance product is being computed and shown in same Section 7. for λ > 0 and θ, α > 0. We denote unconditional distribution of by N GPL(θ, α).
iii For α → ∞, pmf (1) tends to negative binomial N B(2, θ 1+θ ).iv (1) can also viewed as a mixture of G( θ 1+θ ) and N B(2, θ 1+θ ) with mixing proportion θ θ+α .For different values of θ and α, the probability function is evaluated and presented in Figure 1.It can be observed that for fixed θ and increasing α, the distribution accommodates more right tail, whereas, for fixed α and increasing θ, the distribution condenses and the right tail approaches to zero at a faster rate.Our proposed model fits appropriately to those data sets where there is large right tail or the tail approaches to zero at a faster rate.Such data sets are quite common in insurance problems and count data example in biology.The cumulative density function of X ∼ N GPL(θ, α) can be given as 3. Properties of the new generalized Poisson-Lindley distribution

Shape of the probability function
It can be seen that The above expression is a decreasing function in x implying unimodalilty.Further, As such P r(X = x) is log-concave and hence the distribution (1) has an increasing failure rate, see Johnson et al.(2005), p. 209.It can also be verified that, (i) For θ > 1, ( 1) is unimodal and has a unique mode at 0.
The above facts are shown in Fig. 1 for selected values of θ, α.

Quantile function
The quantile function of the N GPL(θ, α) is where W (.) denotes the Lambert W function.
Proof: The c.d.f of the distribution is Multiplying both sides of the above equation by − (1+θ) −θ−θ 2 −α αθ αθ log(1 + θ) and after rearranging the terms, we get Using definition of Lambert-W function (W (z)e W (z) = z, where z is a complex number, see Joŕda (2010)), the right hand side of the above expression is the Lambert W function of the Solving the above equation for x, we get The first three quantiles can be obtained by substituting γ = 1 4 , 1 2 and 3 4 in equation (5).

Generating functions
The Probability Generating Function of (1) is given by The Moment Generating Function works out to

Moments, skewness and kurtosis
The first four raw moments of X can easily be obtained by simple computations and are as follows and the moments of X about the mean are The expressions for skewness  The coefficient of variation (C.V) of the distribution comes out to be The index for dispersion is given by As the new distribution comes from a mixture of a Poisson distribution and two parameter Lindley distribution, the variance-to-mean ratio is greater than one (see Karlis (2005) and Sundt and Vernic (2009), p. 66) the proposed distribution is over-dispersed.

Simulation of random variables
Note that, If Y ∼ P (λ) and λ follows a two parameter Lindley distribution, which can also be represented as a mixture of two independent random variable V 1 ∼ exp(θ) and V 2 ∼ gamma(2, θ) with mixture parameter p = θ θ+α , then following steps can be used to generate N GPL(θ, α) random variates Step 2:

Method of moments
Given a random sample x 1 , x 2 , • • • , x n of size n from (1), the moment estimates, α and θ, of α and θ can be obtained by solving the following equations and where m 1 and m 2 are the first and second sample moments.Solving the above equations, we get Theorem 3: For fixed α, the estimator θ of θ is positively biased, i.e.E( θ) > θ.
Theorem 4: For fixed values of α, the moment estimate θ of θ is consistent and asymptotically normal and distributed as where .

Maximum likelihood estimation
Let x 1 , x 2 , • • • , x n be random observations of size n from our proposed generalized Poisson-Lindley distribution.The log-likelihood function for the vector of parameters Θ = (θ, α) can be written as The associated score function is given by U n = ∂ln ∂θ , ∂ln ∂α , where the ML estimator can be obtained by equating ∂ln ∂θ , ∂ln ∂α ( θ, α) = (0, 0) , hence we get The second partial derivatives are given by The above equations can be solved numerically by using open source R software.

EM algorithm for a new generalize Poisson-Lindley distribution
The MLE of α and θ needs to be computed numerically.Newton-Raphson algorithm is one of the standard methods to estimate the parameters.The algorithm requires second order derivatives of the log-likelihood for all iterations which makes it quite complex.This problem can be handled by another estimation method known as Expectation-Maximization (EM) (see Dempster et al. (1977)).The EM algorithm consists of two steps: the E-step and the M-step.E-Step computes the expectation of the unobservable part given the current values of the parameters and M-step maximizes the complete data likelihood and updates the parameters using the conditional expectations obtained in E-step.This procedure can be useful when there are no closed-form expressions for estimating the parameters and the derivatives of the likelihood are complicated.
To start with, a hypothetical complete-data distribution is defined with joint probability function It is straightforward to verify that the E-step of an EM cycle requires the computation of the conditional expectations of (λ|x i ; θ (h) , α (h) ) and ( λ 1+αλ |x i ; α (h) , θ (h) ), say t (h) i and s (h) i respectively, where (θ (h) , α (h) ) is the current estimate of (θ, α).Using, The EM cycle completes with the M-step, involving complete data maximum likelihood over (θ, α), with the missing λ's replaced by their conditional expectations E(λ|x; θ ( h); α ( h)).Thus the EM iterates until defined convergence criterion is satisfied.

Collective risk model
In non-life insurance portfolio, say, automobile insurance, the aggregate loss (S) is a random variable defined as sum of claims incurred in a certain period of time.Let N be number of claims in a certain period which is a random variable and X i be claim severity random variable, which is independent of N , is the size of i th claim.Thus, aggregate loss is defined as S = N i=0 X i .It is well known that the pdf of S is given as f s (x) = ∞ n=0 p n f n * (x), where p n denotes the probability of n claims (primary distribution) and f n * (x) is the n th fold convolution of f (x), the claim amount (secondary distribution).For more detail on classic risk model, see Freifelder (1974), Rolski et al. (1999), Nadarajah andKotz (2006a and2006b) and Klugman et al. (2012) and reference therein.Here, we consider two such situations: In one case, the primary distribution is as defined in Section 1 and claim severity distribution as exponential distribution with parameter (λ).In the second case, the Erlang distribution with parameters r and λ is the secondary distribution.Erlang loss distribution may arise in insurance settings when the individual claim amount is the sum of exponentially distributed claims.
Theorem 5: If we assume a N GPL(θ, α) as primary distribution and an exponential distribution with parameter (λ) as secondary distribution, then the pdf of aggregate loss random variable S = N i=0 X i is given by whereas, Proof: By assuming that the claim severity follows an exponential distribution with parameter λ > 0, since the n th fold convolution of exponential distribution is gamma distribution with parameter n and λ, the n th fold convolution is given by Then the pdf of the random variable S is given by The pdf of the aggregate loss has a jump of size P r(S = 0) at the origin.The no claim probability comes out to be It is well known in the actuarial literature that the mean of the aggregate loss random variable S can be obtained as E(S) = E(X)E(N ) (see eqn. 9.9 of Klugman(2012)).Hence, when X ∼ exp(λ) then Theorem 6: In collective risk model, if the primary distribution follows N GPL(θ, α) and the secondary distribution follows an Erlang(2,λ) distribution, then the pdf of aggregate loss Proof: By assuming that the claim severity follows an Erlang(2,λ) distribution, the n th fold convolution of exponential distribution is gamma distribution with parameter 2n and λ.The n th fold convolution is given by Then, the pdf of the aggregate loss random variable S is given by The pdf of the aggregate loss random variable has a jump of size P r(S = 0) at the origin.The probability of zero claim comes out to be

Application to real data set
In this section, we fit our proposed distribution to 3 data sets so as to illustrate our claim that our proposed model fits well when compared to other competing models.The first data set has a long right tail and approaches to zero slowly and in the other data sets the right tail approaches zero at a faster rate.In one of the illustrations, we propose a method for calculation of automobile insurance premium under Bonus-Malus system.2. In each of these distributions, the parameters are estimated by using the maximum likelihood method.Further, based on the values of log likelihood and chi-square, we observe that New Generalized Poisson-Lindley (N GPL) distribution provides a satisfactorily better fit for the data set compared to other distributions.
It can be seen that the log-likelihood and Chi-square statistic for the N GPL are lower that those of competing models showing that our model satisfactorily fits better to the data set.(19) Using the data presented in Table 3, we have computed the BMPs which are shown in Table 5.Finally, while calculating the premium under BMS, Lemaire (1979) remarked on the problem of overcharges.It can be seen that N GPL model produces a lower penalization as compared to traditional Poisson-Gamma models (as shown in Table 4 and 5

Final comments
This paper provides a new generalized discrete distribution with an infinite and non-negative integer support.It has been shown that the proposed distribution can be considered as an alternative to well known distribution like Poisson, Poisson-Lindley, negative binomial, weighted generalized Poisson Distribution for the dataset possessing a right tail which approach to zero at a slower/faster rate.Further, in the context to actuarial science, close expression of pdf of aggregate loss random variable is obtained by considering our proposed distribution as primary distribution and exponential and Erlang distributions as secondary distributions.Moreover, our proposed model can also be useful in the determination of Bonus-Malus premium for non-life insurance products.

Illustration 1 :
A dataset representing epileptic seizure counts (see Chakraborty (2010)) are used.Comparability of the proposed model with other distributions like Poisson distribution (P), negative binomial (N B), Generalized Poisson-Lindley (GPL) and Weighted Generalized Poisson (WGP) Distribution (Chakraborty (2010)) have been shown in Table

Table 2 :
Distribution of epileptic seizure countsKlugman et al. (2012, pp.664) presented the distribution of automobile insurance policyholders according to number of claims.We reproduce the data set in Table3.It is found that this data set is highly right skewed and over-dispersed (variance is greater than mean).We mentioned in the text that our proposed model can be applied to such overdispersed data sets.We chose chi-square statisticχ 2 = n (O n − E n ) 2 /E n forcomparison purposes and calculated this statistic for other competing models Poisson, Negative Binomial, Poisson-Lindley, Generalized Poisson-Lindley Distribution.

Table 3 :
Gómez and Vázquez (2003)mobile insuranceIt is known that Poisson distribution is not a suitable choice for automoble insurance claims because of the restriction that mean equals variance for Poisson distribution.As such, Negative Binomial distribution is preferred in this case.In many countries Bonus-Malus system (BMS) is operative wherein policyholders with no claims are given bonus whereas policyholders with claims are punished.Therefore, next year premium depends on the history of the policyholder till this year irrespective of size of the claim.Gómez and Vázquez (2003)calculated the premium under BMS by

Table 5 :
BMP using new generalized Poisson-Lindley distribution