Analysis of Count Data by Transmuted Geometric Distribution

Subrata Chakraborty; Deepesh Bhati

doi:10.2991/jsta.d.191218.001

<Previous Article In Issue

Download article (PDF)

Volume 18, Issue 4, December 2019, Pages 450 - 463

Analysis of Count Data by Transmuted Geometric Distribution

Authors

Subrata Chakraborty¹, Deepesh Bhati²^{, *}

¹Department of Statistics, Dibrugarh University, Assam, India

²Department of Statistics, Central University of Rajasthan, Rajasthan, India

^*Corresponding author. Email: deepesh.bhati@curaj.ac.in

Corresponding Author

Deepesh Bhati

Received 30 November 2017, Accepted 14 September 2018, Available Online 27 December 2019.

DOI: 10.2991/jsta.d.191218.001 How to use a DOI?
Keywords: Transmuted geometric distribution; EM algorithm; Likelihood Ratio test; Rao score's test; Wald's test
Abstract: Transmuted geometric distribution (TGD) was recently introduced and investigated by Chakraborty and Bhati [Stat. Oper. Res. Trans. 40 (2016), 153–176]. This is a flexible extension of geometric distribution having an additional parameter that determines its zero inflation as well as the tail length. In the present article we further study this distribution for some of its reliability, stochastic ordering and parameter estimation properties. In parameter estimation among others we discuss an EM algorithm and the performance of estimators is evaluated through extensive simulation. For assessing the statistical significance of additional parameter α, Likelihood ratio test, the Rao's score tests and the Wald's test are developed and its empirical power via simulation are compared. We have demonstrate two applications of (TGD) in modeling real life count data.
Copyright: © 2019 The Authors. Published by Atlantis Press SARL.
Open Access: This is an open access article distributed under the CC BY-NC 4.0 license (http://creativecommons.org/licenses/by-nc/4.0/).

1. INTRODUCTION

Chakraborty and Bhati [1] recently introduced the transmuted geometric distribution TGDq,α using the quadratic rank transmutation techniques of Shaw and Buckley [2]. It may be noted that though there is a large number of new continuous distribution in statistical literature which are derived using the rank transmutation technique but TGDq,α is the first discrete distribution derived using this technique. Chakraborty and Bhati [1] investigated various distributional properties, showed applicability of TGDq,α in modeling aggregate loss, claim frequency data from automobile insurance and demonstrated the feasibility of TGDq,α as count regression model by considering data from health sector. In the current article, we discussed some additional theoretical and applied aspects of TGDq,α, which are structured as follows. In Section 3, we present various reliability properties and stochastic ordering of TGDq,α. In Section 4, comparative study of maximum likelihood estimator (MLE) obtained numerically and through EM Algorithm are presented through simulation, whereas in Section 5, detailed hypothesis testing is discussed considering three Wald's, Rao's Score and Likelihood Ratio test (LRT) for testing α=0. To illustrate the applicability of TGD models in different disciplines other than those discussed in Chakraborty and Bhati [1], we consider two real data sets and compare them with different family of distributions in Section 6. Finally, some conclusions and comments are presented in Section 7.

2. TRANSMUTED GEOMETRIC DISTRIBUTION (TGD (q, α))

A random variable (rv) X is said to follow Transmuted geometric distribution TGD with two parameters q and α, in short, TGDq,α if its probability mass function (pmf) is given by

py=1−αq y1−q+α1−q2q2y,y=0,1,⋯ .(1)

The corresponding survival function (sf) is written as

F̄Yy=1−αq y+αq2y,y=0,1,⋯ ,(2)

where 0<q<1,−1<α<1. Following distributional characteristics are presented in Chakraborty and Bhati [1].

For α=0, (1) reduces to GDq with pmf py=1−qq y,y=0,1,⋯ ,0<q<1.
For α=−1, (1) reduces to a special case of the Exponentiated Geometric distribution EGD of Chakraborty and Gupta (2015) with power parameter equal to 2. This is the distribution of the maximum of two independent and identical distributed (iid) GDq rvs.
For α=1, (1) reduces to GDq2 with pmf 1−q2q2y, which is the distribution of the minimum of two iid GDq rvs.
For 0<α<1−1<α<0 the TGDq,α distribution with pmf given in (1), the ratio py∕py−1, y=1,2,⋯ , forms a monotone increasing (decreasing) sequence.
TGDq,α is unimodal with a nonzero mode for −1<α<−q2+q−1 provided q>0.414.
The probability generating function (PGF) of TGDq,α is given by
GYz=1−q1−αq1−z−q2z1−qz1−q2z,|q2z|<1.
The rth factorial moment of Y∼TGDq,α is given by
EYr=1−αr!q1−qr+αr!q21−q2r,
where Yr=YY−1…Y−r+1.

3. RELIABILITY PROPERTIES AND STOCHASTIC ORDERING

There are several situations in reliability where continuous time is not a good scale to measure the lifetime, in production we may be interested in how many units are produced by the machine before failure or health insurance companies are interested how long a patient stays in hospital before discharge/death. In such situations, the discrete hazard rate functions can be used to model ageing properties of discrete random lifetimes. The different hazard rate functions of TGD model and associated results are as follows:

3.1. Reliability Properties

3.1.1. Hazard rate function and its classification

The hazard rate function rXx for X∼TGDq,α is given as

rXx=PX=xSXx=1−αqx1−q+α1−q2q2x1−αqx+αq2x,=1−α1−q+αqx1−q21−α+αqx.

The hazard rate function of TGDq,α is plotted in Figure 1 for various values of parameters to investigate the monotonic properties and it is clear that the hazard rate of TGDq,α is increasing for −1<α<0, decreasing when 0<α<1 and constant if α=0 or 1. Also it can be seen that even when α≠1, the hazard rate approaches to a constant as y increases. The smaller the value of q the faster is the rate of stabilization of the hazard rate.

Theorem 1.

The TGDq,α has increasing, decreasing and constant hazard rate for −1<α<0, 0<α<1 and α=0 or 1 respectively.

Proof.

The hazard rate function of TGDq,α is given as

rYy=1−α1−q+αq y1−q21−α+αq y =1−q1−α+αqy+11−α+αq y.

But q1−α+αqy+11−α+αq y is a decreasing (increasing) function of y for −1<α<0 0<α<1. Hence rYy is increasing (decreasing) function of y for −1<α<0 0<α<1. Constant hazard rates are obtained as rYy=1−q for α=0 and rYy=1−q2 for α=1.

Remark.

The hazard rate of TGDq,α clearly obeys rYy≤1−q for −1≤α≤0 and 1−q≤rYy≤1−q2 for 0≤α≤1.

3.1.2. Mean residual life

Kemp (2004) presented various characterization of discrete lifetime distribution among them the mean residual life (MRL) or life expectancy is an important characteristic, for TGD, the closed expression for MRL is given as

LYy=EY−y|Y≥y=1SYy∑j>ySYj=q1+q1−α+αqy+11−q21−α+αq y.(3)

Theorem 2.

The MRL function given in (3) is monotone decreasing (increasing) function of y depending on −1<α<00<α<1.

Proof.

It can be easily be seen that

△LYy=LYy+1−LYy=1−ααqy+11+q1−α1−q y1−α1−qy+1.

For any choice of α∈−1,1 and q∈0,1, the denominator terms 1−α1−q y and 1−α1−qy+1 are always positive. Moreover, since q∈0,1, therefore △LYy<0 for −1<α<0 indicates decreasing MRL, whereas △LYy>0 for 0<α<1 indicates increasing MRL.

3.2. Stochastic Ordering

Many times there is a need of comparing the behavior of one rv with the other. Shaked and Shanthikumar [3] has given many comparisons such as likelihood ratio order ⪯lr, the stochastic order ⪯st, the hazard rate order ⪯hr, the reversed hazard rate order ⪯rh and the expectation order ⪯E having various applications in different context.

Theorem 3.

Let Y be a rv following TGDq,α and X be geometric rv with parameter p. Then Rz=PY=z∕PX=z is an increasing(decreasing) function of z for −1<α<00<α<1 respectively that is X⪯lrYX⪰lrY.

Proof.

Since Rz=1+α1+qqz−1. Thus, we have Rz≤≥Rz+1 for −1<α<00<α<1 for any q∈0,1.

Corollary.

Following results are direct implications of Theorem 3.

X⪯st⪰stY that is, PX≥z≤≥PY≥z for −1<α<00<α<1 respectively and for all z.
Y⪯hr⪰hrX that is, PY=z∕PY≥z≤≥PX=z∕PX≥z for −1<α<00<α<1 respectively and for all z.
X⪯rh⪰rhY that is, PX=z∕PX≤z≤≥PY=z∕PY≤z for −1<α<00<α<1 respectively and for all z.
X⪯E⪰EY that is, EX≤≥EY for −1<α<00<α<1 respectively and for all z.

Theorem 4.

Let Y1 and Y2 be TGDq1,α and TGDq2,α respectively. Then Y2⪯stY1 iff q1≤q2 for 0<α<1.

Proof.

We know that Y2⪯stY1 if PY2≥y≤PY1≥y for all y, hence for TGDq,α with PY≥y=1−αq2y+αq y and it is clearly seen that for 0<α<1

1−αq12y+αq1 y≤1−αq22y+αq2 y∀yiffq1≤q2.

Hence Y2⪯stY1.

4. PARAMETER ESTIMATION AND THEIR COMPARATIVE EVALUATION

Estimates of the parameters q and α of TGDq,α model can be computed by following five methods (i) sample proportion of 1s and 0s method, (ii) sample quantiles, (iii) method of moments and finally (iv) ML method and (v) ML via EM Algorithm. Moreover, in this section we carry out comparative study of ML estimator obtained numerically and via EM Algorithm utilizing initial estimate from one of the first three methods.

4.1. From Sample Proportion of 1's and 0's

If p0,p1 be the known observed proportion of 0's and 1's in the sample, then the parameters q and α can be estimated by solving the equations

p0=1−α1−q+α1−q2andp1=1−αq1−q+αq21−q2.

4.2. From Sample Quantiles

If t1,t2 be two observed points such that FYt1=γ1,FYt2=γ2, then the two parameters q and α can be estimated by solving the simultaneous equations

γ1=1+α−1q t1+1−αq2t1+1andγ2=1+α−1qt2+1−αq2t2+1.

4.3. Methods of Moments

Denoting the first and second observed raw moments by m1 and m2 respectively, the moment estimates can be obtained by

Either solving the following two equations simultaneously
q1−α+q21−q2=m1andq1+q3−αq3q+2+11−q22=m2,
or by the minimization method proposed by Khan et al. [4] by minimizing EY−m12+EY2−m22 with respect to q and α
q1−α+q21−q2−m12+q1+q3−αq3q+2+11−q22−m22.

4.4. ML Method

Let y=y1,y2,⋯ ,yn⊤ be a sample of n observations drawn from TGD distribution, and Θ=q,α⊤ be the parametric vector. The log-likelihood (LL) function for the corresponding sample is

l=logL=nlog1−q+logq∑i=1nyi+∑i=1nlog1−α+αqyi1+q,(4)

and the score function UΘ,y=∂l∂q,∂l∂α⊤ can be obtained by differentiating LL function with respect to q and α as

∂l∂q=−n1−q+1q∑i=1nyi+∑i=1nαq yi+αyi1+qq yi−11−α+α1+qq yi,∂l∂α=∑i=1n1+qq yi−11−α+α1+qq yi.

The MLE Θ̂ of Θ is obtained by solving the nonlinear system of equation UΘ,y=0. Since the likelihood equations have no closed form solution, the estimator q̂ and α̂ of the parameters q and α can be obtained by maximizing LL function using global numerical maximization techniques using optim package of R. Further, the Fishers information matrix is given by

Iyq,α=−E∂2l∂q2−E∂2l∂q∂α−E∂2l∂q∂α−E∂2l∂α2≈−∂2l∂q2−∂2l∂q∂α−∂2l∂q∂α−∂2l∂α2q=q̂,α=α̂,(5)

where q̂ and α̂ are the MLE's of q and α respectively. Moreover elements of Iyq,α may be computed from

∂2l∂q2=−n1−q2−1q2∑i=1nyi−∑i=1nα1+qyi−1yiq yi−2+2αyiq yi−11−α+α1+qqyi−α1+qyiq yi−1+αqyi1−α+α1+qqyi2,∂2l∂q∂α=∑i=1n1+qyiq yi−1+qyi1−α+α1+qqyi−α1+qyiq yi−1+αqyi1+qqyi−11−α+α1+qqyi2,∂2l∂α2=−∑i=1n1+qqyi−121−α+α1+qqyi.

4.5. MLE through EM Algorithm

The Expected Maximization (EM) algorithm is an useful iterative procedure to compute ML estimators in the presence of missing data or assumed to have a missing values. The procedure follows with two steps called Expectation step (E-Step) and Maximization step (M-Step). The E-step concerns with the estimation of those data which are not observed whereas the M-step is a maximization step for more details one may refer Dempster et al. [5]. Moreover, theorem 3.1 of Redner and Walker [6] ensures the consistency and uniqueness of the estimates obtained by EM procedure.

Let the complete-data be constituted with observed set of values y=y1,⋯ ,yn and the hypothetical data set x=x1,⋯ ,xn, where the observations yi's are distributed with rv X defined as

X=1with probability1+α∕20with probability1−α∕2,(6)

and rv Y be defined as

Y=XZ1:2+1−XZ2:2,(7)

where Z1:2∼GDq2, Z2:2∼EGDq,2 (see [7]) and Xi∼Bernoulli1+α2.

Under the formulation, the E-step of an EM cycle requires the expectation of X|Y;Θk, where Θk=qk,αk is the current estimate of Θ (in the kth iteration). Since the conditional distribution of Xi given Yi is

Xi|Yi,Θk∼bernoulli1+αik2,(8)

with

1+αik2=1+αk1−qk2qk2yi1+αk1−qk2qk2yi+1−αk1−qkqkyi2−qk+1qkyi,(9)

where αk is a set of known or estimated parameters at kth step with known initial values. Thus, by the property of the Binomial distribution, the conditional mean is

EXi|Yi,Θk=1+αik2andVXi|Yi,Θk=1+αik21−αik2.(10)

For M-step: The likelihood function of joint pmf of hypothetical complete-data Yi,Xi,i=1,⋯ ,n is given as

L∗Θ;y,x=∏i=1n1+α2xi1−q2q2yixi ⋅∏i=1n1−α21−xi1−qqyi2−1+qqyi1−xi

and the corresponding complete LL function is given as

ln∗Θ;x,y=log1+α2∑i=1nxi+log1−α2∑i=1n1−xi+log1−q2∑i=1nxi +2logq∑i=1nxiyi+∑i=1n1−xiyilogq+log1−q+log2−qyi1+q.(11)

The components of the score function Un∗Θ=∂ln∗∂α,∂ln∗∂q⊤ are given by

∂ln∗∂α=11+α∑i=1nxi−11−α∑i=1n1−xi,(12)

∂ln∗∂q=−2q1−q2∑i=1nxi+2q∑i=1nxiyi +∑i=1n1−xiyiq−11−q−yiq yi−1+yi+1q yi2−q yi1+q.(13)

The EM cycle will completed with the M-step by using the MLE over Θ, that is, Un*Θ^;y,x=0 with the unobserved xis replaced by their conditional expectations given in (10). Hence we obtain the iterative procedure of the EM algorithm as

α̂k+1=1n∑i=1nαik,q̂k+1=∑i=1n1+αik2yi2q̂k+11−q̂k+12∑i=1n1+αik2−∑i=1n1−αik2yiq̂k+1−11−q̂k+1−yiq̂k+1yi−1+yi+1q̂k+1yi2−q̂k+1yi1+q̂k+1,

where q̂k+1 should be determined numerically.

4.5.1 Standard errors of estimates obtained from EM algorithm

In this section, we obtain the standard errors (se) of the estimators from the EM algorithm using result of Louis [8]. Let z=y,x, then the elements of 2×2 observed information matrix IcΘ,z=∂∂ΘUcΘ;z are given by

∂2ln∗∂α2=−11+α2∑i=1nxi−11−α2∑i=1n1−xi,∂2ln∗∂α∂q=∂2ln∗∂q∂α=0,∂2ln∗∂q2=−21+q21−q22∑i=1nxi−2q2∑i=1nxiyi−∑i=1n1−xiyiq2+q2yi−2qyi+q+yi22−q+1q yi2 +yi−1yiq yi−2+yiyi+1q yi−12−q+1q yi+11−q2.

Taking the conditional expectation of IcΘ;z given x, we obtain the 2×2 matrix

lcΘ;y=−EIcΘ;z|y=dij,(14)

where

whereas computation of

lmΘ;y=VUcx;θ|y=mij,(15)

involve the following terms

m11=11+α+11−α2∑i=1nVXi|y,m12=m21=∑i=1n11+α+11−αyiq−2q1−q2+11−q+yiq yi−1+yi+1q yi2−q yi1+qVXi|y,m22=∑i=1nyiq−2q1−q2+11−q+yiq yi−1+yi+1q yi2−q yi1+q2VXi|y.

Finally, the observed information matrix I can be computed as

IΘ̂;y=lcΘ̂;y−lmΘ̂;y,

and IΘ̂;y can be inverted to obtain an estimate of the covariance matrix of the incomplete-data problem. The square roots of the diagonal elements represent the estimates of the standard errors of the parameters.

4.6. Simulation Study to Evaluate EM Algorithm

Here we study the behavior of ML estimators obtained by direct numerical optimization and also through EM algorithm for different finite sample sizes and for different TGDq,α. Observations from TGDq,α are generated using the quantile function provided in Chakraborty and Bhati [1] (see result 4 of Table 1). In the next two subsections, first we investigate the performance of ML estimators q̂,α̂ for various combinations of parameters q,α in subsection (3.6.1) and then evaluate the performance with respect to varying sample size for fixed parameter values in subsection (3.6.2).

			MLE				EM Algorithm
Parameters		n	bias(α̂)	bias(q̂)	se(α̂)	se(q̂)	bias(α̂)	bias(q̂)	se(α̂)	se(q̂)
q = 0.25	α = −0.75	25	−0.5566	0.0099	1.5154	0.1144	−0.0132	0.0232	0.9739	0.1147
		50	−0.2675	0.0101	0.9151	0.0866	−0.0081	0.0122	0.7215	0.0835
		75	−0.1733	0.0050	0.6880	0.0694	−0.0049	0.0073	0.5881	0.0677
		100	−0.1327	0.0053	0.5780	0.0600	−0.0035	0.0052	0.5137	0.0589

q = 0.5	α = −0.75	25	−0.1348	−0.0149	0.5644	0.0859	−0.0031	−0.0058	0.5664	0.0854
		50	−0.0077	−0.0001	0.3960	0.0619	0.0012	−0.0029	0.3888	0.0601
		75	−0.0196	−0.0012	0.3197	0.0498	−0.0006	−0.0011	0.3155	0.0489
		100	−0.0113	−0.0026	0.2765	0.0432	0.0012	−0.0028	0.2730	0.0424

q = 0.75	α = −0.75	25	−0.0411	−0.0003	0.4012	0.0480	0.0060	−0.0035	0.4190	0.0476
		50	−0.0085	−0.0026	0.2766	0.0333	−0.0002	−0.0021	0.2964	0.0333
		75	−0.0011	−0.0018	0.2242	0.0268	0.0012	−0.0020	0.2337	0.0268
		100	−0.0008	−0.0014	0.1909	0.0227	−0.0009	−0.0007	0.1990	0.0229

q = 0.25	α = −0.30	25	−0.3455	0.0269	1.2484	0.1361	−0.0340	0.0260	1.0043	0.1426
		50	−0.3055	0.0095	0.9006	0.1016	−0.0240	0.0165	0.7203	0.1018
		75	−0.0391	0.0290	0.6466	0.0901	−0.0125	0.0109	0.6045	0.0848
		100	−0.0997	0.0123	0.5770	0.0756	−0.0090	0.0069	0.5269	0.0736

q = 0.5	α = −0.30	25	−0.0310	0.0045	0.6288	0.1097	−0.0037	0.0000	0.6840	0.1153
		50	−0.0249	0.0011	0.4672	0.0803	−0.0031	−0.0009	0.4818	0.0818
		75	−0.0253	−0.0002	0.3880	0.0668	−0.0029	−0.0008	0.3908	0.0664
		100	−0.0258	−0.0008	0.3375	0.0580	−0.0036	−0.0004	0.3333	0.0568

q = 0.75	α = −0.30	25	−0.0503	0.0000	0.5182	0.0625	−0.0010	−0.0046	0.6044	0.0689
		50	−0.0432	0.0010	0.3850	0.0459	0.0069	−0.0030	0.4157	0.0482
		75	−0.0020	0.0003	0.3141	0.0369	−0.0038	0.0000	0.3324	0.0381
		100	−0.0009	0.0000	0.2838	0.0330	−0.0026	−0.0005	0.2899	0.0332

Note: EM = Expected Maximization; MLE = maximum likelihood estimator.

Table 1

Average Bias and Average SE computed by method of maximum likelihood and EM Algorithm method.

4.6.1. Performance of estimators for different parametric values

A simulation study consisting of following steps is carried out for each triplet q,α,n, considering q = 0.25, 0.5, 0.75, α = −0.70, −0.30, 0.30, 0.70 and n = 25, 50, 75, 100.

Choose the value q0,α0 for the corresponding elements of the parameter vector Θ=q,α, to specify the TGDq,α;
Choose sample size n;
Generate N independent samples of size n from TGDq,α;
Compute the ML and EM estimate Θ̂n of Θ for each of the N samples;
Compute the average bias, standard error (SE) of the estimate.

In our experiment we have considered the number of replication N=1000. It can be observed from Tables 1 and 2 that as the sample size increases both average bias and average se decreases.

			MLE				EM Algorithm
Parameters		n	bias(α̂)	bias(q̂)	se(α̂)	se(q̂)	bias(α̂)	bias(q̂)	se(α̂)	se(q̂)
q = 0.25	α = 0.30	25	−0.4174	0.0254	1.1108	0.1524	−0.1138	−0.0206	0.7619	0.1547
		50	−0.2518	0.0178	0.8702	0.1281	−0.0667	−0.0158	0.8095	0.1632
		75	−0.1338	0.0193	0.6331	0.1110	−0.0481	−0.0144	0.6889	0.1386
		100	−0.0878	0.0215	0.5479	0.1013	−0.0367	−0.0032	0.6740	0.1353

q = 0.50	α = 0.30	25	−0.2343	0.0226	0.5962	0.1328	−0.0404	−0.0267	0.7990	0.1700
		50	−0.1440	0.0184	0.4884	0.1085	−0.0335	−0.0354	0.6296	0.1349
		75	−0.0611	0.0142	0.4132	0.0926	−0.0319	−0.0336	0.5801	0.1237
		100	−0.0586	0.0125	0.3970	0.0886	−0.0213	−0.0210	0.5013	0.1072

q = 0.75	α = 0.30	25	−0.0594	0.0127	0.5713	0.0829	−0.0143	−0.0326	0.7882	0.1101
		50	−0.0316	0.0097	0.4540	0.0652	−0.0173	−0.0508	0.6689	0.0923
		75	−0.0250	0.0079	0.3969	0.0568	−0.0177	−0.0607	0.5449	0.0759
		100	−0.0081	0.0050	0.3729	0.0522	−0.0107	−0.0224	0.5029	0.0691

q = 0.25	α = 0.75	25	−0.0975	0.0038	0.0240	0.0042	−0.0234	−0.0305	0.6862	0.1166
		50	−0.0696	0.0138	0.0255	0.0027	−0.0189	−0.0220	0.5239	0.0423
		75	−0.0995	0.0046	0.0751	0.0038	−0.0125	−0.0112	0.4443	0.0259
		100	−0.0358	0.0070	0.0338	0.0027	−0.0101	−0.0071	0.4012	0.0125

q = 0.5	α = 0.75	25	−0.1250	0.0288	0.5170	0.1474	−0.0351	−0.0112	0.5214	0.1627
		50	−0.1162	0.0248	0.4238	0.1186	−0.0158	−0.0131	0.4862	0.1456
		75	−0.0641	0.0140	0.3485	0.1000	−0.0093	−0.0583	0.3675	0.1088
		100	−0.0493	0.0125	0.3422	0.0974	−0.0037	−0.1109	0.6810	0.1963

q = 0.75	α = 0.75	25	−0.1542	0.0350	0.5112	0.0966	−0.0191	−0.0176	0.5221	0.1048
		50	−0.1114	0.0178	0.4014	0.0727	−0.0350	−0.0253	0.4722	0.0858
		75	−0.0786	0.0100	0.3595	0.0638	−0.1072	−0.0986	0.3782	0.0676
		100	−0.0455	0.0100	0.3139	0.0566	−0.1168	−0.1178	0.3662	0.0647

Note: EM = Expected Maximization; MLE = maximum likelihood estimator.

Table 2

Average Bias and Average SE computed by method of maximum likelihood and EM Algorithm method.

4.6.2 Performance of estimators for different sample size

In this subsection, we assess the performance of ML estimators of q,α as sample size n, increases by considering n=25,26,…,200, for q=0.25 and α=−0.5. For each n, we generate one thousand samples of size n and obtain MLEs and their standard error. For each repetition we compute average bias and average squared error.

Figures 2 and 3 show behavior of average bias and average standard error of parameter q and α, for fixed q=0.25 and α=−0.5, as one varies sample size n. The horizontal dotted lines from Figure 2 corresponds to zero value and it is clear in Figure 2 that the biases approach to zero with increasing n. Further, in Figure 3, average standard errors for both parameters (q and α) decrease with increase in n. Similar observations were also noted for other parametric values.

Based on our findings it is clear that EM algorithm produces better ML estimators with smaller average bias as compared to the regular ML estimators while w.r.t. standard error there is not much to choose between the two procedures.

5. TESTS OF HYPOTHESIS

The TGDq,α distribution with parameter vector Θ=q,α⊤ reduces to the Geometric distribution with parameter q when α=0. This additional parameter α controls the proportion of zeros of the distribution relative to geometric distribution and also the tail length. Therefore it is of interest to develop test procedure for detecting departure of α from 0. In this section we develop the LRT, the Rao's score test and the Wald's test for testing the null hypothesis H0:α=0 against the alternative hypothesis H1 : α≠0 and numerically study the statistical power of these tests through extensive simulation.

5.1. LRT, Rao's Score Test and Wald's Test

The LRT is based on the difference between the maximum of the likelihood under null and the alternative hypotheses. The LRT statistics is given by −2logL(Θ^∗;y)L(Θ^;y) where Θ̂∗ and Θ̂ are the MLE obtained under the null and alternative hypotheses respectively. The LRT is generally employed to test the significance of the additional parameter which is included to extend a base model.

The Rao's Score test [9] is based on the score vector defined as the first derivative of the LL function w.r.t. the parameters. Rao's score test statistic UI−1U⊤, where U is the score vector and I is the information matrix derived under the null hypothesis. The score vector and the information matrix, obtained by evaluating the derivatives of the LL function, logL is provided in Section 4.4. Note that the scores actually are the slopes of the likelihood functions.

The Wald's test statistics [10] is based on the difference between the maximum of the likelihood estimate value of the parameter under alternative hypothesis and the value of specified null hypothesis. The Wald's test statistic is given in our case by α̂−α0I[22]−1α̂−α0⊤, where I[22]−1 is the 2,2th element of the inverse of the information matrix I, and α̂ is the MLE of α both under alternative hypotheses, whereas α0 is the value of α as per H0. Note that I[22]−1 is an estimate of the variance of α. Therefore in the present case our Wald's statistic reduces to α̂2Vα̂.

All the test statistics follow asymptotically chi-square distribution with “k” degrees of freedom (df), where “k” is the number of parameters specified by the null hypothesis. In the present case the df is just “1”. For well behaved likelihood function all these tests are based on measuring the discrepancy between null and the alternative hypotheses.

5.2. Statistical Power Analysis

Here we present a simulation based study of the statistical power of LR test, Rao's Score test and Wald's test considering 5% level of significance. Since the tests are asymptotic in nature we have considered four different sample sizes, namely n=100,300,500 and 1000. We have generated 1000 replications for each sample size n. The power of these tests are estimated by proportion of rejection in these 1000 replications. The effect size (ES) is a measure of departure from the null hypothesis which in the present case is given by α−0=α is fixed at −0.7,−0.5,−0.3,−0.1,0.1,0.2,0.5,0.7 for our experiments.

The power of the these test for different sample size, ES and for different parametric value q are presented in Tables 3 and 4. Figures 4 and 5 reveal, as expected, that the power of the test increases with the sample size n and ES. Further, for positive ES, all the tests show increase in power with the increase in either or both ES and sample size, Moreover, for negative ES, power increase in a much faster pace. Power for score test is more than LR test for negative ES where as it is other way for positive ES. For positive ES the power of the tests gets closer with increase in sample size. From the over all observation it is clear that the Wald's test is more reliable than both LR test and Score tests.

q = 0.10
n	100			300			500			1000
α	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald
−0.7	0.028	0.051	0.042	0.117	0.066	0.017	0.252	0.133	0.01	0.476	0.351	0.087
−0.5	0.02	0.046	0.041	0.052	0.039	0.029	0.121	0.066	0.024	0.237	0.134	0.015
−0.3	0.022	0.081	0.056	0.024	0.032	0.037	0.059	0.035	0.034	0.101	0.045	0.019
−0.1	0.013	0.11	0.066	0.02	0.069	0.061	0.036	0.045	0.046	0.059	0.04	0.032
0.1	0.015	0.173	0.074	0.021	0.123	0.109	0.042	0.112	0.099	0.058	0.105	0.071
0.3	0.009	0.239	0.091	0.026	0.209	0.143	0.038	0.212	0.161	0.079	0.237	0.134
0.5	0.005	0.28	0.112	0.041	0.351	0.217	0.069	0.394	0.266	0.142	0.431	0.288
0.7	0.033	0.262	0.132	0.025	0.452	0.236	0.074	0.523	0.336	0.197	0.602	0.436

q = 0.30
n	100			300			500			1000
α	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald

−0.7	0.305	0.565	0.127	0.742	0.851	0.741	0.922	0.963	0.927	0.998	0.999	0.999
−0.5	0.137	0.303	0.047	0.412	0.537	0.389	0.619	0.707	0.620	0.917	0.942	0.924
−0.3	0.074	0.177	0.047	0.147	0.219	0.129	0.207	0.274	0.200	0.457	0.519	0.462
−0.1	0.049	0.098	0.064	0.059	0.076	0.065	0.065	0.076	0.053	0.089	0.109	0.085
0.1	0.041	0.072	0.086	0.052	0.055	0.101	0.056	0.054	0.060	0.080	0.071	0.051
0.3	0.034	0.090	0.129	0.082	0.101	0.151	0.153	0.156	0.180	0.296	0.292	0.202
0.5	0.043	0.153	0.172	0.139	0.213	0.271	0.289	0.336	0.347	0.546	0.575	0.455
0.7	0.276	0.181	0.468	0.265	0.300	0.555	0.367	0.424	0.627	0.634	0.642	0.725

q = 0.45
n	100			300			500			1000
α	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald

−0.7	0.470	0.787	0.563	0.933	0.982	0.956	0.989	0.997	0.993	1.000	1.000	1.000
−0.5	0.241	0.540	0.310	0.611	0.792	0.675	0.835	0.909	0.873	0.993	0.996	0.994
−0.3	0.089	0.279	0.125	0.223	0.368	0.280	0.325	0.465	0.391	0.641	0.729	0.699
−0.1	0.058	0.157	0.089	0.076	0.128	0.105	0.071	0.106	0.085	0.090	0.137	0.122
0.1	0.035	0.083	0.078	0.060	0.062	0.096	0.059	0.055	0.071	0.097	0.073	0.051
0.3	0.033	0.062	0.117	0.117	0.083	0.171	0.210	0.163	0.193	0.396	0.313	0.233
0.5	0.055	0.106	0.199	0.224	0.200	0.316	0.417	0.351	0.427	0.700	0.645	0.539
0.7	0.268	0.121	0.468	0.347	0.227	0.639	0.497	0.377	0.711	0.763	0.685	0.805

Note: LR = likelihood ratio.

Table 3

Power of the LR, Score and Walds Test for different sample sizes n, effect size and parametric value q.

q = 0.6
n	100			300			500			1000
α	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald
−0.7	0.628	0.888	0.760	0.985	0.997	0.996	1.000	1.000	1.000	1.000	1.000	1.000
−0.5	0.281	0.630	0.434	0.745	0.875	0.825	0.920	0.962	0.947	0.997	0.998	0.998
−0.3	0.120	0.351	0.213	0.273	0.457	0.364	0.453	0.627	0.550	0.748	0.830	0.802
−0.1	0.045	0.178	0.113	0.070	0.139	0.110	0.081	0.125	0.108	0.109	0.180	0.150
0.1	0.042	0.103	0.108	0.046	0.054	0.097	0.072	0.046	0.078	0.113	0.080	0.057
0.3	0.043	0.077	0.141	0.143	0.082	0.195	0.267	0.165	0.208	0.450	0.358	0.257
0.5	0.064	0.089	0.202	0.252	0.172	0.350	0.481	0.336	0.415	0.784	0.689	0.583
0.7	0.265	0.083	0.485	0.392	0.188	0.667	0.563	0.387	0.741	0.817	0.697	0.863

q = 75
n	100			300			500			1000
α	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald	LR	Score	Wald

−0.7	0.698	0.941	0.852	0.994	1.000	0.998	1.000	1.000	1.000	1.000	1.000	1.000
−0.5	0.336	0.686	0.532	0.798	0.921	0.868	0.956	0.987	0.979	0.998	0.999	0.999
−0.3	0.142	0.374	0.245	0.297	0.512	0.418	0.488	0.669	0.600	0.810	0.877	0.860
−0.1	0.045	0.188	0.131	0.057	0.150	0.119	0.077	0.153	0.130	0.090	0.167	0.141
0.1	0.045	0.112	0.115	0.060	0.056	0.110	0.095	0.071	0.104	0.092	0.057	0.045
0.3	0.044	0.072	0.139	0.145	0.076	0.179	0.299	0.182	0.236	0.470	0.356	0.245
0.5	0.091	0.078	0.242	0.285	0.164	0.355	0.491	0.321	0.427	0.789	0.683	0.586
0.7	0.316	0.089	0.525	0.412	0.176	0.657	0.585	0.369	0.764	0.845	0.698	0.860

Note: LR = likelihood ratio.

Table 4

Power of the LR, Score and Walds Test for different sample sizes n, effect size and parametric value q.

6. DATA ANALYSIS

For the purpose of illustration, in this section, we consider following two data sets with details as follows:

Number of Fires in Greece (NTG)

The data comprise of numbers of fires in district forest of Greece from period 1 July 1998 to 31 August 1998. The observed sample values of size 123 for these data are the following (frequency in parentheses and none when it is equal to one): 0(16),1(13), 2(14), 3(9), 4(11), 5(13), 6(8), 7(4), 8(9), 9(6), 10(3), 11(4), 12(6), 15(4), 16, 20, 43. The data were previously studied by Bakouch et al. [11] and Karlis and Xekalaki [12].
Insurance Claim Count (ICC)

An insurance count data from Belgium in year 1993 is considered [13] and the data is as follows: 0(57178), 1(5617), 2(446), 3(50), 4(8).

The null hypothesis H0:α=0 against H1:α≠0 are examined utilizing the LR, Rao's Score and Wald's test, and the results along with the descriptive statistics and with respective p-values are presented in Table 5. Moreover, p-value (less than 0.05) for Rao's Score and Wald's test reject the null hypothesis at 5% significance level. The suitability of the proposed TGDq,α model with other competitive distributions namely weighted geometric distribution WGDq,α(Bhati and Joshi [14], Negative Binomial NBr,p, Poisson-weighted exponential distribution PWEDα,θ (see [15]) and GDq is carried out and the LL, Akaiki Information Criteria (AIC), Bayesian Information Criteria (BIC) value are computed for five models for both the datasets. Further the Kolmogorov-Smirnov goodness of fit value with bootstrap p- value is also computed. The results in Table 6 such as maximum LL, minimum AIC, BIC and KS value for TGDq,α reveals that the TGDq,α is the best fitted model and could be consider as competitive model for the datasets considered.

Data Set	Mean	Variance	LR test	Score Test	Wald's Test
NTG	5.398	30.045	3.568 (0.059)	41.018 (0.000)	5.423 (0.020)
ICC	0.106	0.115	9.178 (0.002)	8.268 (0.004)	8.118 (0.004)

Notes: ICC = insurance claim count; LR = likelihood ratio; NTG = number of fires in Greece.

Table 5

Descriptive and test statistic with p-value in parentheses for both the datasets.

Dataset ↓	Model →	WGDq,α	NBr,p	PWEDα,θ	TGDq,α	GDq
NTG	Estimates	(0.802, 1.804)	(1.496, 0.223)	(2.470, 0.247)	(0.802, −0.531)	0.839
	LL	−335.552	−334.529	−335.552	−334.415	−337.165
	AIC	675.104	673.058	675.104	672.830	676.330
	BIC	680.728	678.682	680.728	678.454	679.142
	KS	0.547	0.449	0.525	0.431	0.906
	p-value	0.767	0.882	0.847	0.899	0.432

ICC	Estimates	(0.0845, 0.749)	(1.279, 0.924)	(5.870, 10.837)	(0.085, −0.157)	0.095
	LL	−22063.8	−22064.3	−22063.8	−22063.6	−22068.2
	AIC	44131.6	44132.6	44131.6	44131.2	44138.4
	BIC	44149.7	44150.7	44149.7	44149.3	44147.4
	KS	0.091	0.102	0.091	0.074	0.427
	p-value	0.878	0.862	0.882	0.918	0.525

Notes: AIC = Akaiki Information Criteria; BIC = Bayesian Information Criteria; ICC = insurance claim count; LL = log likelihood; NTG = number of fires in Greece.

Table 6

Comparative study of data fitting.

7. CONCLUSION

The current paper investigates some additional property of the TGDq,α distribution with emphasis on the simulation study of the behaviors of the parameter estimation and also power of tests of hypothesis to check statistical significance of the additional parameter. In the parameter estimation we have presented different methods including the EM algorithm implementation of the MLE. A comparative simulation based evaluation of the EM algorithm based MLE against the usual MLE has revealed the superiority of the former in terms of the bias and mean squared errors. We have also presented data modeling examples to showcase the advantage of the TGDq,α over some of the existing distributions from literature. As such it is envisaged that the present contribution will be useful for discrete data analysts.

CONFLICT OF INTEREST

The authors declare they have no conflicts of interest.

AUTHORS' CONTRIBUTIONS

SC contributes in Section 1–3,5,7 whereas DB contributes in Section 4–7.

Funding Statement

This research is not funded by any agency.

ACKNOWLEDGMENT

The authors acknowledge with thanks the comments of two anonymous referees which helped immensely in improving the presentation of the paper.

REFERENCES

1.S. Chakraborty and D. Bhati, Stat. Oper. Res. Trans., Vol. 40, 2016, pp. 153-176.

2.W. Shaw and I. Buckley, The Alchemy of Probability Distributions: Beyond Gram- Charlier Expansions and a Skew-Kurtotic-Normal Distribution from a Rank Transmutation Map, 2007. Research Report, arXiv:0901.0434v1 [q-fin.ST]

3.M. Shaked and J. Shanthikumar, Stochastic Orders and Their Applications, Probability and Mathematical Statistics, Academic Press, Boston, 1994.

4.M.S.A. Khan, A. Khalique, and A.M. Aboummoh, IEEE Trans. Reliab., Vol. 38, 1989, pp. 348-350.

5.A.P. Dempster, N.M. Laird, and D.B. Rubin, J. R. Stat. Soc. Ser. B, Vol. 39, 1977, pp. 1-38.

6.R.A. Redner and H.R. Walker, SIAM Rev., Vol. 26, 1984, pp. 195-239.

7.S. Chakraborty and R.D. Gupta, Exponentiated geometric distribution: another generalization of geometric distribution, Commun. Stat. Theory Methods, Vol. 44, 2015, pp. 1143-1157.

8.T. Louis, J. R. Stat. Soc. Ser. B Methodol., Vol. 44, 1982, pp. 226-233.

9.C.R. Rao, Large sample tests of statistical hypotheses concerning several parameters with applications to problems of estimation, Math. Proc. Camb. Philos. Soc., Vol. 44, 1948, pp. 50-57.

10.A. Wald, Trans. Am. Math. Soc., Vol. 54, 1943, pp. 426-482.

11.H.S. Bakouch, M.A. Jazi, and S. Nadarajah, Statistics, Vol. 48, 2014, pp. 200-240.

12.D. Karlis and E. Xekalaki, E.A. Lipitakis (editor), in Proceedings of the Fifth Hellenic-European Conference on Computer Mathematics and Its Applications, 2001, pp. 872-877.

13.M. Denuit, Astin Bull., Vol. 27, 1997, pp. 229-242.

14.D. Bhati and S. Joshi, Commun. Stat. Theor. Methods, Vol. 47, 2018, pp. 1510-1527.

15.H. Zamani, N. Ismail, and P. Faroughi, J. Math. Stat., Vol. 10, 2014, pp. 148-154.

<Previous Article In Issue

Download article (PDF)

Journal: Journal of Statistical Theory and Applications
Volume-Issue: 18 - 4
Pages: 450 - 463
Publication Date: 2019/12/27
ISSN (Online): 2214-1766
ISSN (Print): 1538-7887
DOI: 10.2991/jsta.d.191218.001 How to use a DOI?
Open Access: This is an open access article distributed under the CC BY-NC 4.0 license (http://creativecommons.org/licenses/by-nc/4.0/).

Cite this article

ris enw bib

TY  - JOUR
AU  - Subrata Chakraborty
AU  - Deepesh Bhati
PY  - 2019
DA  - 2019/12/27
TI  - Analysis of Count Data by Transmuted Geometric Distribution
JO  - Journal of Statistical Theory and Applications
SP  - 450
EP  - 463
VL  - 18
IS  - 4
SN  - 2214-1766
UR  - https://doi.org/10.2991/jsta.d.191218.001
DO  - 10.2991/jsta.d.191218.001
ID  - Chakraborty2019
ER  -

download .riscopy to clipboard