Bias and MSE of the estimates under the maximum likelihood method.
In this chapter, a new generalization of the Kumaraswamy distribution, namely the gamma-Kumaraswamy distribution is defined and studied. Several distributional properties of the distribution are discussed in this chapter, which includes limiting behavior, mode, quantiles, moments, skewness, kurtosis, Shannon’s entropy, and order statistics. Under the classical method of estimation, the method of maximum likelihood estimation is proposed for the inference of this distribution. We provide the results of an analysis based on two real data sets when applied to the gamma-Kumaraswamy distribution to exhibit the utility of this model.
- gamma-Kumaraswamy distribution
- Renyi’s entropy
- reliability parameter
- stochastic ordering
The generalization of a distribution by mixing it with another distribution over the years has provided a mathematical based way to model a wide variety of random phenomena statistically. These generalized distributions are effective and flexible models to analyze and interpret random durations in a possibly heterogeneous population. In many situations, observed data may be assumed to have come from such a mixture population of two or more distributions.
Two parameter gamma and a two parameter Kumaraswamy are most popular distribution for analyzing any lifetime data. Gamma distribution is a well-known distribution, and it has several desirable properties .
A serious limitation of the gamma distribution, however, is that the distribution function (or survival function) is not available in a closed form if the shape parameter is not an integer, thereby it requires some numerical methods to evaluate these quantities. As a consequence, this distribution is less attractive as compared to Ref. , which has nice tractable distribution function, survival function and hazard function. In this paper, we consider a four parameter gamma-Kumaraswamy distribution. It is observed that it has many properties which are quite similar to those of a gamma distribution, but it has an explicit expression for the distribution function or the survival functions. The major motivation of this chapter is to introduce a new family of distributions, make a comparative study of this family with respect to a Kumaraswamy family and a gamma family and provide the practitioner with an additional option, with a hope that it may have a ‘better fit’ compared to a gamma family or Kumaraswamy family in certain situations. It is noteworthy to note that the gamma-Kumaraswamy distribution is a generalization of Kumaraswamy distribution with the property that it can exhibit various shapes. ( Figure 1 ). This provides more flexibility to the gamma-Kumaraswamy distribution in comparison with Kumaraswamy distribution in modeling different data sets. The property of left-skewness is a rare characteristic as it is not enjoyed by several generalizations of Kumaraswamy distribution. Our proposed model is different from that of Ref. , where the authors have proposed a generalized gamma-generated distribution with an extra positive parameter for any continuous baseline G distribution.
The rest of the paper is organized as follows. In Section 2, we propose the gamma-Kumaraswamy distribution [GK(α, β, a, b)]. In Section 3, we study various properties of the GK(α, β, a, b) including the limiting behavior, transformation, and the mode. In Section 4, the moment generating function, the moments and the mean deviations from the mean and the median, and Renyi’s entropy are studied. In Section 5, we consider the maximum likelihood estimation of the GK(α, β, a, b). In Section 6, we provide an expression for the reliability parameter for two independent GK(α, β, a, b) with different choices for the parameters α and β but for a fixed choice of the two shape parameters of Kumaraswamy distribution. In Section 7, discussion is made for the moment generating function of the r-th order statistic and also the limiting distribution of the sample minimum and the sample maximum for a random sample of size n drawn from GK(α, β, a, b). An application of GK(α, β, a, b) is discussed in Section 8. Certain characterizations of GK(α, β, a, b) are presented in Section 9. In Section 10, some concluding remarks are made.
2. The gamma-Kumaraswamy distribution
We consider the following class of gamma-X class of distributions, for which, the parent model being
where α, β are positive parameters. Also, is the density function [cumulative distribution function] of the random variable X. Furthermore, is the survival function of the associated random variable X.
If X has density Eq. (1), then the random variable has a gamma distribution with parameters α, β. The reverse happens to be true as well. Here, we consider G(.) to be the cdf of a Kumaraswamy distribution with parameters a, b. Then, the cdf of the gamma-Kumaraswamy (hereafter GK) reduces to
where with is the regularized incomplete gamma function. So the density and hazard functions corresponding to Eq. (2) are given, respectively, by
The percentile functions for GK distribution: The p th percentile x p is defined by F(x p ) = p. From Eq. (2), we have . Define , then , where is the inverse of regularized incomplete gamma function. Hence, .
In the density equation (3), a, b, and α are shape parameters and β is the scale parameter. It can be immediately verified that Eq. (3) is a density function. Plots of the GK density and survival rate function for selected parameter values are given in Figures 1 and 2 , respectively.
If X~GK(a, b, α, β), then the survival function of X, S(x) will be
We simulate the GK distribution by solving the nonlinear equation
where u has the uniform (0,1) distribution. Some facts regarding the GK distribution are as follows:
If X~GK(a, b, α, β), then X m ~GK(a, b, α, β), .
Also, we have the following important result: If X~GK(1, b, α, β), then X 1/a ~GK(a, b, α, β), .
The GK distribution does not possess the reproductive property. In other words, if for any two X 1~GK and X 2~GK then the distribution of the sum S = X 1 + X 2 will not be a GK.
The first result provides an important property of the GK distribution for information analysis is that this distribution is closed under power transformation. The latter result is equally important because it provides a simple way to generate random variables following the GK distribution.
3. Properties of GK distribution
The following lemma establishes the relation between GK(α, β, a, b) distribution and gamma distribution.
Lemma 1. (Transformation): If a random variable X follows a gamma distribution with parameters α and β, then follows GK(α, β, a, b) distribution.
Proof. The proof follows immediately by using the transformation technique. W
The limiting behaviors of the GK pdf and its hazard function are given in the following theorem.
Theorem 1. The limits GK density function, f(x), and the hazard function, , are given by
Proof. Straightforward and hence omitted. W
Theorem 2. The mode of the GK distribution is the solution of the equation where
Proof. The derivative of f(x) in Eq. (3) can be written as
The critical values of Eq. (10) are the solutions of W
Next, we discuss the IFR and/or DFR property of the hazard function for the GK distribution. For this, we will consider the result of Lemma 1. According to Lemma 1, if X~GK(a, b, α, β), then ∼Gamma (α, β). In such a case for the random variable Y, the hazard rate function can be written as
Therefore, . If , is decreasing in t and hence r(t) is increasing, thereby and has a IFR. If , then
is increasing in t, so r(t) decreases and hence has a DFR. Now, since X is a one-to-one function of Y, the hazard rate function of X will also follow the exact pattern.
Let X and Y be two random variables. X is said to be stochastically greater than or equal to Y denoted by if for all x in the support set of X.
Theorem 3. Suppose X~GK and Y~GK If β 1 > β 2, a 1 > a 2 and b 1 < b 2. Then , for integer values of a 1 and a 2.
Proof. At first, we note that the incomplete gamma function is an increasing function of x for fixed α. For any real number , and , we have
This implies that . Equivalently, it implies that , and this completes the proof. W
Note: For fractional choices of a 1 and a 2, the reverse of the above inequality will hold.
4. Moments and mean deviations
For any ,
Upper bounds for the -th order moment: Since , for , from Eq. (13), one can write , provided r/a and j/b+k−1 are both integers. Employing successively, the generalized series expansion of , the characteristic function for X~GK will be given by [from Eq. (3)]
If j/a and k 1/b are integers then in Eq. (14), the second and third summations will stop at j/a and k 1/b, respectively.
If we denote the median by T, then the mean deviation from the mean, , and the mean deviation from the median, can be written as
Using the substitution in Eq. (17), we obtain
where we used successively binomial series expansion.
By using Eqs. (2) and (18), the mean deviation from the mean and the mean deviation from the median are, respectively, given by
One useful measure of diversity for a probability model is given by Renyi’s entropy. It is defined as , where and . If a random variable X has a GK distribution, then we have
Next, consider the integral
Now, using successive application of the generalized binomial expansion, we can write
Hence, the integral in Eq. (21) reduces to
Therefore, the expression for the Renyi’s entropy will be
5. Maximum likelihood estimation
In this section, we address the parameter estimation of the GK(α, β, a, b) under the classical set up. Let X 1, X 2, …, X n be a random sample of size n drawn from the density Eq. (3). The log-likelihood function is given by
The derivatives of Eq. (13) with respect to α, β, a, and b are given by
To estimate the model parameters, numerical iterative techniques must be used to solve these equations. We may investigate the global maxima of the log likelihood by setting different starting values for the parameters. The information matrix will be required for interval estimation. The elements of the 4 × 4 total observed information matrix (since expected values are difficult to calculate), (for ), can be obtained from the authors under request, where . The asymptotic distribution of is , under the regularity conditions, where is the expected information matrix, and is the observed information matrix. The multivariate normal distribution can be used to construct approximate confidence intervals for the individual parameters.
5.1. Simulation study
In order to assess the performance of the MLEs, a small simulation study is performed using the statistical software R through the package (stats4), command MLE. The number of Monte Carlo replications was 20,000 For maximizing the log-likelihood function, we use the MaxBFGS subroutine with analytical derivatives. The evaluation of the estimates was performed based on the following quantities for each sample size; the empirical mean squared errors (MSEs) are calculated using the R package from the Monte Carlo replications. The MLEs are determined for each simulated data, say, for , and the biases and MSEs are computed by
for . We consider the sample sizes at n = 100, 200, and 500 and consider different values for the parameters . The empirical results are given in Table 1 . The figures in Table 1 indicate that the estimates are quite stable and, more importantly, are close to the true values for these sample sizes. Furthermore, as the sample size increases, the MSEs decrease as expected.
|Sample size||Actual value||Bias||MSE|
6. Reliability parameter
The reliability parameter R is defined as , where X and Y are independent random variables. For a detailed study on the possible applications of the reliability parameter, an interested reader is suggested to look at Ref. [4, 5]. If X and Y are two continuous and independent random variables with the cdf’s and and their pdf’s and , respectively, then the reliability parameter R can be written as
Theorem 4. Let X~GK(a, b, α 1, β 1) and Y~(a, b, α 2, β 2), then
Using the series expansion for the incomplete gamma function , and using the substitution , Eq. (34) reduces to
Hence the proof. W
7. Order statistics
Here, we derive the general r-th order statistic and the large sample distribution of the sample minimum and the sample maximum based on a random sample of size n from the GK(α, β, a, b) distribution. The corresponding density function of the r-th order statistic, from Eq. (3) will be
Using the series expression for the incomplete gamma function: , the pdf of can be written as
From Eq. (37), it is interesting to note that the pdf of the r-th order statistic can be expressed as an infinite sum of the GK pdf ’s.
Here, we consider two well-known illustrative data sets which are used to show the efficacy of the GK distribution. For details on these two data sets [6, 7], the second data set in Table 2 is from Ref. , and it represents the fatigue life of 6061-T6 aluminum coupons cut parallel with the direction of rolling and oscillated at 18 cycles per second. The GK distribution is fitted to the first data set and compared the result with the Kumaraswamy, gamma-uniform , and beta-Pareto . These results are reported in Table 3 . The results show that gamma-uniform, GK distributions provide adequate fit to the data. Figure 3 displays the empirical and the fitted cumulative distribution functions. This figure supports the results in Table 3 . A close look at Figure 3 indicates that GK distribution provides better fit to the left tail than the gamma-uniform distribution. This is due to the fact that GK distribution can have longer left tail ( Figure 3 ).
In addition, to check the goodness-of-fit of all statistical models, several other goodness-of-fit statistics are used and are computed using computational package Mathematica. The MLEs are computed using N maximize technique as well as the measures of goodness-of-fit statistics including the log-likelihood function evaluated at the MLEs (l), Akaike information criterion (AIC), corrected Akaike information criterion (AICC), consistent Akaike information criterion (CAIC), the Anderson-Darling (A *), the Cramer-von Mises (W *), and the Kolmogrov-Smirnov (K-S) statistics with their p values to compare the fitted models. These statistics are used to evaluate how closely a specific distribution with cdf (2) fits the corresponding empirical distribution for a given data set. The distribution with better fit than the others will be the one having the smallest statistics and largest p value. Alternatively, the distribution for which one obtains the smallest of each of these criteria (i.e., AIC, AICC, K-S, etc.) will be most suitable one. The mathematical equations of those statistics are given by
where denotes the log-likelihood function evaluated at the maximum likelihood estimates, q is the number of parameters, n is the sample size and , the y i ’s being the ordered observations.
Lieblein and Zelen  proposed a five parameter beta generalized Pareto distribution and fitted the data in Table 4 and compared the result with beta-Pareto and other known distributions. The results of fitting beta generalized Pareto and beta-Pareto from Ref.  are reported in Table 4 along with the results of fitting the Pareto (IV) and GK distributions to the data. The KS value from Table 4 indicates that the GK distribution provides the best fit. The fact that GK distribution has the least number of parameters than beta generalized Pareto and beta-Pareto adds an extra advantage over them. Figure 4 displays the empirical and the fitted cumulative distribution functions. This figure supports the results in Table 4 .
|Distribution||Kumaraswamy||Beta-Pareto||Beta generalized Pareto||gamma-Kumaraswamy|
9. Characterization of GK distribution
In this section, we present characterizations of GK distribution in terms of the ratio of two truncated moments. For the previous works done in this direction, we refer the interested readers to Glänzel [11–14] and Hamedani [15–17]. For our characterization results, we employ a theorem due to Ref. , see for further details. The advantage of the characterizations given here is that cdf F need not have a closed form. We present here a corollary as a direct application of the theorem discussed in details in Ref. .
Corollary 1. Let be a continuous random variable and let and for Then X has pdf (3) if and only if the function η defined in Theorem 5 has the form
Proof. Let X has pdf (3), then
Conversely, if η is given as above, then
Now, in view of Theorem 5, X has pdf (3).
Corollary 2. Let be a continuous random variable and let h(x) be as in Proposition 1. Then, X has pdf (3) if and only if there exist functions g and η defined in Theorem 5 satisfying the differential equation
Remarks 1. (a) The general solution of the differential equation in Corollary 1 is
for 0 < x < 1, where D is a constant. One set of appropriate functions is given in Proposition 1 with D = 0
(b) Clearly, there are other triplets of functions (h, g, η) satisfying the conditions of Theorem 5. We presented one such triplet in Proposition 1.
10. Concluding remarks
A special case of the gamma-generated family of distributions, the gamma-Kumaraswamy distribution, is defined and studied. Various properties of the gamma-Kumaraswamy distribution are investigated, including moments, hazard function, and reliability parameter. The new model includes as special sub-models the gamma and Kumaraswamy distribution. Also, we provide various characterizations of the gamma-Kumaraswamy distribution. An application to a real data set shows that the fit of the new model is superior to the fits of its main sub-models. As future work related to this univariate GK model, we will consider the following:
A natural bivariate extension to the model in Eq. (1) would be
In this case, exact evaluation of the normalizing constant would be difficult to obtain, even for a simple analytic expression of a baseline bivariate distribution function, G(x, y). Numerical methods such as Monte Carlo methods of integration might be useful here. We will study and discuss structural properties of such a bivariate GK model.
Extension of the proposed univariate GK model to multivariate GK models and discuss the associated inferential issues. It is noteworthy to mention that classical methods of estimation, such as for example, maximum likelihood method of estimation might not be a good strategy because of the enormous number of model parameters. An appropriate Bayesian inference might be the only remedy. In that case, we will separately study two different cases of estimation: (a) with non-informative priors and (b) with full conditional conjugate priors (Gibbs sampling). Since the GK distribution is in the one parameter exponential family, a reasonable choice for priors for α and β might well be gamma priors with appropriate choice of hyper-parameters. For prior choices of the parameters that are from the baseline G(.) distribution function, a data-driven prior approach will be more suitable.
A discrete analog of the univariate GK model with a possible application in modeling rare events.
Construction of a new class of GK mixture models by adopting Marshall-Olkin method of obtaining new distribution.