The U Family of Distributions: Properties and Applications


 In this article, we develop a new general family of distributions aimed at unifying some well-established lifetime distributions and offering new work perspectives. A special family member based on the so-called modified Weibull distribution is highlighted and studied. It differs from the competition with a very flexible hazard rate function exhibiting increasing, decreasing, constant, upside-down bathtub and bathtub shapes. This panel of shapes remains rare and particularly desirable for modeling purposes. We provide the main mathematical properties of the special distribution, such as a tractable infinite series expansion of the probability density function, moments of several kinds (raw, incomplete, probability weighted…) with discussions on the skewness and kurtosis. The stochastic ordering structure and stress-strength parameter are also considered, as well as the basics of the order statistics. Then, an emphasis is put on the inferential features of the related model. In particular, the estimation of the model parameters is employed by the maximum likelihood method, with a simulation study to confirm the suitability of the approach. Three practical data sets are then analyzed. It is observed that the proposed model gives better fits than other well-known lifetime models derived from the Weibull model.


Introduction
The Weibull distribution has an advantage over the exponential distribution because of its decreasing and increasing hazard rate function (hrf) depending upon the shape parameter. It has been widely used in many areas, mainly for the analysis of reliability data, time to failure equipment and testing. See, for instance, [26] for its use in the analysis of the failure analysis of capacitors, [10] for its application in vacuum pulse technique, and [5] where related control charts was designed under accelerated hybrid censoring.
Due to its vast applications, numerous researchers have attempted to improve the performance of the Weibull distribution by introducing more parameters, developing various generalized or extended Weibull distributions. The most useful of them are the exponentiated Weibull distribution introduced by [17], additive Weibull distribution studied by [27], generalized Weibull distribution created by [18], generalized modified Weibull (GMW) distribution proposed by [8], modified Weibull (MW) distribution introduced by [24] and new extended Weibull distribution proposed by [19], among others.

The U family of distributions
This section sets the mathematical foundation of the article.

Definition
The following result is at the basis of the newly defined family. α > 0, φ ≥ 0 and θ ∈ R such that θφ ≤ 1, the following function is well defined and has the properties of a cumulative distribution function (cdf ): F (x) = 0 for x < a and F (x) = 1 for x > b.
Therefore, e −R(x) (1+φR(x)) θ is decreasing and e −R(x) (1+φR(x)) θ ≤ lim x→a e −R(x) (1+φR(x)) θ = 1, ending the first part of the proof. Let us now prove that F (x) has the properties of a cdf. By attesting that F (x) is increasing. This ends the proof of Proposition 2.1.
Based on Proposition 2.1, by adopting the same notations, we define the U family of distributions (U for unified) by the cdf given by

1)
F (x) = 0 for x < a and F (x) = 1 for x > b, recalling that α > 0, φ ≥ 0 and θ ∈ R such that θφ ≤ 1 (so θ can be negative). The corresponding probability density function (pdf) is obtained as and f (x) = 0 for x ∈ [a, b]. Also, the corresponding hrf is given by where W (x) denotes the well-known Lambert function satisfying W (x)e W (x) = x. Note that, if the inverse function of R(x) exists and has an analytical expression, one can express analytically Q(y). The interests and motivations behind the U family are described below.

Discussion
As suggested by its name, the U family unifies a plethora of existing distributions defined with a polynomial-exponential function as main term. For instance, when a = 0 and b = +∞, it includes various extensions or generalizations of the Weibull, linear failure rate and Lindley distributions. Some of them are presented below.
The list of the referenced distributions belonging to the U family becomes endless if various values for a and b are considered, such as the distributions defined on a bounded interval or the real line, i.e., a = −∞ and b = +∞ (including various symmetric and skewed distributions. . . ).
Also, the U family can be expressed by considering a baseline cdf denoted by G(x), with corresponding pdf denoted by g(x). Indeed, one can remark that the function R(x) = − log(1 − G(x)), corresponding to the cumulative hrf associated to G(x), satisfies all the conditions of Proposition 2.1. This choice of function is also motivated by [2] in the context of the TX family. Hence, by substituting this expression of R(x) into (2.1), we get The corresponding pdf is obtained as Also, the corresponding hrf is given by and h(x) = 0 for x ∈ [a, b]. The qf of the U family, say Q(y), is involved in the following nonlinear equation: F [Q(y)] = y, i.e., if θφ = 0, G[Q(y)] = y 1/α and, if θφ = 0, To the best of our knowledge, under this quite standard formulation involving a baseline distribution, the U family has no analogue in the literature. Even though some special members are well known, the general definition of the U family opens up unexplored horizons that deserve full attention. In this study, with the utopian idea to find the ultimate flexible lifetime distribution, we use the U family to propose a new promising extension of the Weibull distribution.

Extended modified Weibull distribution
This section is devoted to the presentation of the new distribution, as well as its significant mathematical properties.

Definition
After preliminary investigations, a special member of the U family stands up in several desirable aspects. It constitutes a new four-parameters lifetime distribution defined by the cdf given as (2.2) with a = 0, b = +∞, θ = 1, φ = 1 and R(x) = λx + βx k , i.e., x ≥ 0. (3.1) and F (x) = 0 for x < 0. Here, it is assumed that λ + β > 0, allowing λ = 0 if β > 0 and vice-versa, and k > 0. The definition of F (x) uses the same configuration of the MW distribution, with the activation of the parameters θ, φ and α such as θ = 1, φ = 1 and α > 0. For the purposes of this paper, the corresponding distribution is naturally called the extended MW (EMW) distribution and denoted by EMW(α, λ, β, k) when the parameters need to be specified. Another point of view is to remark that F (x) is defined as (2.2) with the cdf of the MW distribution, i.e., G(x) = 1−e −(λx+βx k ) for x ≥ 0, and G(x) = 0 for x < 0, and θ = 1, φ = 1 and α > 0. The remainder of the study presents the main motivations and interests of the EMW distribution, both from a theoretical side and from an applied perspective.

On the pdf and hrf
First of all, the pdf of the EMW distribution is derived from the differentiation of the F (x) according to x; it is thus given by The corresponding hrf is given by and h(x) = 0 for x < 0. The corresponding qf satisfies the following equation: Due to their analytical complexity, the shape properties of f (x) and h(x) are complicated to analyse in full generality. Figure 1 helps in this regard by plotting these functions for different values of the parameters α, β, λ and k. From Figure 1 (i), we see that the pdf may be almost symmetric, left skewed, right skewed and decreasing with reverse J shape. Also, from Figure 1 (ii), increasing, decreasing, increasing-decreasing, constant, upside-down bathtub and bathtub shapes for the hrf are observed. This panel of shapes is not obtained for the former MW distribution as attested by [23: Figure 2]. These plots reveal the high level of flexibility of these functions, which constitutes a great advantage for the EMW distribution; it can be used as an appropriate statistical model for a great diversity of data sets.

Some properties
Here, we present some important properties of the EMW distribution.

Infinite series expansion
We now give the infinite series expansion of the pdf of the EMW distribution in terms of weighted pdfs of the MW distribution.
Proposition 3.1. For x ≥ 0, we can express the pdf of the EMW distribution as g(x; β(j + 1), λ(j + 1), k) denotes the pdf of the MW distribution with parameters β(j + 1), λ(j + 1) and k, i.e., P r o o f. Under the considered assumptions on the parameters, we have e −(λx+βx k ) (1+λx+βx k ) ∈ (0, 1). Hence, the generalized binomial theorem gives Now, by applying the standard binomial theorem two times in a row, we get This ends the proof of Proposition 3.1.
The obtained expansion for f (x) is useful to provide other expansions for determinant measures of the EMW distribution, such as raw moments, incomplete moments, etc. These points will be considered later.
The result below is similar to Proposition 3.1 but with the special function: P r o o f. The proof is rigorously similar to the one of Proposition 3.1. It is enough to notice that, thanks to (3.1) and (3.2), we have so we need to replace α by α(s + 1) in the first binomial coefficient of the series expansion of f (x). This ends the proof of Proposition 3.2.

Raw moments
We now consider a random variable X following the EMW(α, λ, β, k) distribution, i.e., with the cdf and pdf given by (3.1) and (3.2), respectively. Then, the r-th raw moment of X is defined by µ r = E(X r ), i.e., This integral is complex and no closed form exists. For given parameters, we can compute it by employing a mathematical software. On the other side, an infinite series expansion can be provided. Indeed, by applying Proposition 3.1, we get where µ * r+ +1+m(k−1);j denotes the r + + 1 + m(k − 1)-th raw moment of the MW distribution with parameters β(j + 1), λ(j + 1) and k, which are given in the following result: [24: Theorem 3.1], with the "k" of the result replaced by r + + 1 + m(k − 1), the "α" replaced by λ(j + 1), the "β" replaced by β(j + 1) and the "γ" replaced by k. Note that the expressions of these raw moments depend on the following cases: λ > 0 and β > 0, λ = 0 and β > 0, and λ > 0 and β = 0. For practical purposes, we can exploit µ r by truncating the infinite sums by a large integer and thus obtain an acceptable approximation.
Several important quantities follow from the raw moments. In particular, the mean of X is given by µ = µ 1 . Also, the r-th central moment of X is given by µ r = E [(X − µ 1 ) r ] which can be expressed by a direct application of the standard binomial formula and the raw moments. As a famous special case, the variance of X is given as σ 2 = µ 2 and the coefficient of variation is defined by CV = σ/µ. Also, the skewness and kurtosis of the EMW distribution can be studied by using the r-th general coefficient of X given by C r = µ r /σ r ; the skewness coefficient of X is given by C 3 and the kurtosis coefficient of X is given by C 4 . Table 1 provides some numerical results in order to illustrate the influence of the parameters α, λ, β and k on the moment measures introduced above. The R software is used (see [21]). Table 1. Numerical values for the four first raw moments, variance, skewness, kurtosis and coefficient of variation of the EMW distribution for some parameter values.  Table 1 shows that the moments measures of the EMW distribution are very varying; the mean can be closed to 0 as for the configuration (α, β, λ, k) = (0.5, 6, 0.5, 0.5), or large as for (α, β, λ, k) = (6, 0.5, 0.3, 0.5), the variance can be small as for (α, β, λ, k) = (6, 6, 0.5, 6) and large as for (α, β, λ, k) = (6, 0.5, 0.3, 0.5). Also, the EMW distribution can be left skewed as for (α, β, λ, k) = (0.5, 0.5, 0.5, 6) where C 3 = −0.7541 < 0 or right skewed as for (α, β, λ, k) = (0.5, 6, 0.5, 0.5), among others. The kurtosis of the EMW distribution can be of all states: leptokurtic as for (α, β, λ, k) = (6, 0.5, 6, 0.5) since C 4 = 448.8895 > 3, platykurtic as for (α, β, λ, k) = (6, 0.5, 0.3, 0.5) since C 4 = −4.4606 < 3 and mesokurtic since the case C 4 = 3 is clearly possible. These facts illustrate numerically the great flexibility of the EMW distribution on these aspects.

Probability weighted moments
The notion of raw moments can be generalized by considering the probability weighted moments. They naturally appear in several applied contexts, as the order statistics which will be discussed later. The (r, s)-th moment of X is defined by We can proceed as for the raw moments to evaluate it. For an analytical approach, one can use Proposition 3.2 to express µ o r,s as an infinite series expansion involving raw moments of the MW distribution, i.e., j, ,m µ * r+ +1+m(k−1);j .
Naturally, we rediscover the formula for the r-th raw moment by taking s = 0.

Incomplete moments
The incomplete moments of X are defined as truncated versions of the raw moments of X. For any y ≥ 0, the r-th incomplete moment of X at y is given by µ r (y) = E(X r 1 (0,y) (X)), i.e., Again, no direct solution is available. Only a mathematical software can be employed to evaluate it numerically. However, as for the raw moments, an analytical approach based on the infinite series expansion of Proposition 3.1 can be proposed. Indeed, we have where µ * r+ +1+m(k−1);j (y) denotes the r + + 1 + m(k − 1)-th incomplete moment of the MW distribution with parameters β(j + 1), λ(j + 1) and k. It can be determined in a similar manner to the r + + 1 + m(k − 1)-th raw moment. The incomplete moments are the key tools to define several widely used probabilistic measures and functions, as the mean deviations, residual and reverse residual life functions and Lorenz curve, for instance.

Stochastic ordering
The concept of stochastic ordering gives the mathematical backgrounds on how a random variable can be bigger than another. We refer to [15] for further detail and applications. The EMW distribution enjoys a likelihood ratio order described in the following result. Proposition 3.3. Let X 1 and X 2 be two random variables such that X 1 follows the EMW(α 1 , λ, β, k) distribution and X 2 follows the EMW(α 2 , λ, β, k) distribution. Thus, if α 1 < α 2 , we have X 1 ≤ lr X 2 .
P r o o f. Owing to (3.2), for x ≥ 0, the pdfs of X 1 and X 2 are, respectively, given by Hence, the likelihood order X 1 ≤ lr X 2 is equivalent to say that the ratio function f 1 (x)/f 2 (x) is decreasing. Now, we have This ends the proof of Proposition 3.3.

Stress-strength parameter
Here, we express the stress-strength parameter defined by R = P (X 2 < X 1 ), where X 1 and X 2 are two random variables following the EMW distribution. The interest of R is to measure the performance of a system where its strength is modeled by X 1 and a certain stress is modeled by X 2 . The former work on this topic can be found in [7]. More applications are given in [14], and the references therein. Here, we only consider varying α; we suppose that X 1 follows the EMW(α 1 , λ, β, k) distribution and that X 2 follows the EMW(α 2 , λ, β, k) distribution. Thus, by (3.2) and (3.1), when x ≥ 0, pdf of X 1 and cdf of X 2 are, respectively, given by Hence, by the independence of X 1 and X 2 , we can express R as where f * (x) denotes the pdf of the EMW(α 1 + α 2 , λ, β, k) distribution. So, the last integral is equal to one and This reveals that R have a simple expression. For statistical purposes, the estimation of α 1 and α 2 yields a direct estimation of R by the plug-in method.

Order statistics
In practice, the ordered values of a data set, as the smallest value or largest value, can provide important information on the nature of the related phenomenon. Then, an appropriate theory is given by the order statistics of a random sample. Here, we investigate some distributional properties of the order statistics in the EMW distribution setting. Let X 1 , X 2 , . . . , X n be n independent random variables following the EMW distribution. Then, the pdf of the i-th order statistic, say X i:n , is given by i.e., by (3.1) and (3.2), and f i:n (x) = 0 for x < 0. By taking i = 1 and i = n, we obtain the first and last order statistic, respectively. In particular, one can notice that X n:n = sup(X 1 , X 2 , . . . , X n ) follows the EMW(αn, λ, β, k) distribution.
For more probabilistic investigations, one can give an infinite series expansion of f i:n (x). Indeed, by applying the standard binomial formula and Proposition 3.2, for x ≥ 0, we get is given as in Proposition 3.2 with s = u + i − 1, the same for g(x; β(j + 1), λ(j + 1), k). As an immediate application, we can express the r-th raw moment of X i:n as an infinite series expansion involving raw moments of the MW distribution, i.e., Other measures can be described in a similar manner.

Estimation method
For the statistical purposes of the EMW model, the estimation of unknown parameters is an essential task. In this regard, several methods exist in the literature, but the maximum likelihood estimation method remains the most commonly employed. The resulting estimators, called maximum likelihood estimators (MLEs), enjoy desirable properties. In particular, the normal approximation for the MLEs in large sample distribution theory is easily handled either analytically or numerically, allowing to construct asymptotic confidence intervals and test statistics. In this section, we investigate the MLEs of the parameters of the new distribution, provide a simulation study to validate their numerical efficiency and specify four goodness-of-fit statistics to compare the fitted pdfs to any data set. Given n observed values x 1 , x 2 , . . . , x n from a random variable following the EMW(α, λ, β, k) distribution, the MLEs of the parameters in θ = (α, λ, β, k) are FARRUKH JAMAL ET AL.
By solving the equations of the system ∂ (θ)/∂θ = 0 4 simultaneously, we get the MLEs. One can notice that, thanks to the first equation, the MLEs satisfy the following relation: Some numerical iterative techniques may be used for evaluating the MLEs; the global maxima of the log-likelihood can be investigated by setting different starting values for the parameters. Under standard regularity conditions, the asymptotic distribution behind the random estimator version ofθ − θ is the multivariate normal N 4 (0 4 , I(θ) −1 ), where I(θ) = E 4 [J(θ)] is the expected Fisher information matrix. The approximate multivariate normal distribution given by N 4 (0 4 , J(θ) −1 ), where J(θ) −1 is the inverse of the observed Fisher information matrix evaluated at θ =θ, can be used to set up approximate confidence intervals for the model parameters. The elements of J(θ) can be obtained from the authors upon request.

A simulation study
We now perform a simulation study to investigate the numerical performance of the MLEs of the EMW model parameters. In this regard, we generate N = 3000 samples of sizes n = 50, 150 and 300 from the EMW(α, λ, β, k) distribution with various target values for α, λ, β and k, and we determine the average estimates and mean square errors (MSEs) of the parameters. The results of the simulation are reported in Table 2. We use the R program for the random generation and specially, the optim-CG routine to obtain the MLEs. From this table, as expected, we observe that the average estimates approach the corresponding target values as n increases. For n = 300 only, these estimates are almost equal to the corresponding target values. Also, when n increases, we see that the MSEs decrease to zero, with very small values for n = 300.
We complete the previous numerical study by graphical work. We generate N = 3000 samples of sizes n = 10 to 120 from the EMW distribution with target parameters values α = 5, β = 1.2, λ = 0.2 and k = 1.2 and plot the curves of the average biases for the estimates in Figure 2 and the curves of the MSEs for the estimates in Figure 3.
From these figures, it is clear that the average biases and MSEs decrease as n increases; all the curves tend quickly to the axis y = 0. Also, one can remark that n = 80 seems enough to obtain a precise estimation in our considered setting.

Statistical criteria
To check the goodness-of-fit of the fitted models, we compute four goodness-of-fit statistics, namely the Anderson-Darling (A * ), Cramér-von Mises (W * ) and Kolmogorov-Smirnov (K-S) statistics with their K-S p-values. These statistics can determine how closely a distribution fits the empirical distribution of the data. The distribution with a better fit than the others will be the one having the smallest statistics and largest p-values. Mathematically, these statistics are given  by

THE U FAMILY OF DISTRIBUTIONS: PROPERTIES AND APPLICATIONS
where z i = F * (y i ), F * (x) denotes the cdf of the considered distribution and y 1 , y 2 , . . . , y n are the ordered observations. We also consider the Akaike information criterion (AIC) and Bayesian information criterion (BIC) defined as AIC = −2 log L + 2k, BIC = −2 log L + k log(n), respectively, where − log L denotes the minus estimated log-likelihood function and k is the number of parameters. As usual, the best model is the one having the smallest AIC and BIC.

Applications
In this section, we fit some well-known distributions and the EMW distribution to three data sets coming from environmental science and reliability. Then, we select the best model among them by considering the statistical benchmarks presented above.

Data sets
The considered data sets are presented below, roughly named "Carbon dioxide", "Airborne" and "Rainfall".
Carbon dioxide: The first data set represents the annual mean growth rate of carbon dioxide at Mauna Loa, Hawaii. The measurements are given in parts per million by year (ppm/yr). The data are available at the earth system research laboratory website from the following electronic link: https://www.esrl.noaa.gov/gmd/ccgg/trends/gr.html We analyse these data from the period of 1959 to 2014. Airborne: The second data set is constituted from 46 observations of active repair times for an airborne communication transceiver. The unit of the measurement is the hour. The data can be found in [25: p. 156]. Rainfall: The third data set contains the mean of maximum daily rainfall during the period of 1975 to 2004 at 35 stations in Malaysia. This data set was recently studied by [13].
Elementary statistical descriptions of these three data sets are provided in Table 3.    From this figure, we see that three outliers are present for Airborne, explaining the highly right skewed nature of this data set. Also, from the quantile points of view, we see that Carbon dioxide and Rainfall are similar.
Let us now investigate the nature of the underlying hrf behind the data. In this regard, we present the total time test (TTT) plots of the three data sets in Figure 5. Further detail on the definition and interpretation of the TTT plot can be found in [1].
From this figure, we observe that the underlying hrfs are probably increasing for Carbon dioxide and Rainfall, and probably bathtub for Airborne. Those cases are covered by the EMW model, motivating its use to fit these data sets.

Competitors and parametric estimation
The EMW model is compared with the models corresponding to the following distributions, all having support on (0, +∞): • The transmuted-Weibull (TW) distribution proposed by [3]; it is defined by the following pdf: • The MW distribution introduced by [24]; it is defined by the following pdf: f (x) = (λ + βkx k−1 )e −(λx+βx k ) , λ, β, k > 0.
• The odd log-logistic modified-Weibull (OLLMW) distribution proposed by [22]; it is defined by the following pdf: The maximum likelihood method is used for estimating the parameters of all the models. The MLEs with their standard errors (computed by inverting the observed information matrices) are given in Tables 4, 5 and 6 for the three data sets in order, respectively.  From the obtained estimates of the EMW model, based on our theory and the plug-in approach, we can derive estimation of important measures, such as the mean, variance, skewness and kurtosis. Table 7 indicates these estimates. dioxide and Rainfall. Accordingly, we can conclude that the EMW distribution provides the best fit to the current data sets among the compared distributions.
We complete this numerical study by some graphics. Figures 6, 7 and 8 display the estimated pdfs superimposed on the related histograms and the estimated cdfs compared to the related empirical cdfs for the three data sets in order, respectively.  Visually, the fits of the EMW model seem better to the others. We now put an emphasis on the EMW model fitting in Figures 9, 10 and 11, showing the previous graphics without the other models, and by the consideration of the probability-probability (P-P) plots and the quantilequantile (Q-Q) plots, where nice fits of the considered distribution over the data is characterized by a nice adjustment of the scatter plot by the P-P or Q-Q lines, respectively.  Figure 9. Plots of the estimated pdf superimposed on the histogram, estimated cdf compared to the empirical cdf, P-P plot and Q-Q plot for the EMW model for Carbon dioxide.
As expected, in all the situations, very adequate fits are observed for the EMW distribution.  Figure 11. Plots of the estimated pdf superimposed on the histogram, estimated cdf compared to the empirical cdf, P-P plot and Q-Q plot for the EMW model for Rainfall.