Harris Extended Power Lomax Distribution: Properties, Inference and Applications

In this work, we present a five-parameter life time distribution called Harris power Lomax (𝐻𝑃𝐿) distribution which is obtained by convoluting the Harris-G distribution and the Power Lomax distribution. When compared to the existing distributions, the new distribution exhibits a very flexible probability functions; which may be increasing, decreasing, J, and reversed J shapes been observed for the probability density and hazard rate functions. The structural properties of the new distribution are studied in detail which includes: moments, incomplete moment, Renyl entropy, order statistics, Bonferroni curve, and Lorenz curve etc. The 𝐻𝑃𝐿 distribution parameters are estimated by using the method of maximum likelihood. Monte Carlo simulation was carried out to investigate the performance of MLEs. Aircraft wind shield data and Glass fibre data applications demonstrate the applicability of the proposed model.


Introduction
In the last decades, the Lomax distribution introduced by Lomax (1954), has been discovered to be very useful in several areas of applications most especially in applied sciences such as applications in flood, queue theory, internet traffic control, life testing, wind speed, sea waves and many others. The use of Lomax distribution is motivated by its heavy tail and simplicity in applications. It provides an alternative to the exponential-type distributions (Rayleigh, Weibull, exponential, Gompertz, pareto.), for further study see Bryson (1974). A random variable is said to be distributed according to Lomax distribution, if it distribution function is given by ( ; , ) = 1 − 1 + , > 0 with , > 0. Thus, is a shape parameter and is a scale parameter. The associated probability density function ( ) corresponding to (1) is given by It should be noted that the of Lomax distribution is naturally a special case of some well-known distributions, and this includes Feller-Pareto, Pareto type II, Pareto type IV, Fisher distribution, and many others. However, the Lomax distribution is limited in applications as a result of some of its limitations, and this includes: Lack of flexibility, heavy tailed features, poor fits etc. when used to model real life data which exhibits non-monotonic, bathtub failure rate. Based on the afore-mentioned reasons various efforts have been made to generalize the Lomax distribution in other to induce flexibility into the distribution and also improve its fits for a better modeling capability. Among these are: the gamma Lomax distribution by Cordeiro and Ortega (2015), Marshall-Olkin Extended Lomax distribution by Ghitany et al. (2007), Al -Zahrani and Sagor (2014) developed and studied the Poisson Lomax distribution, the Exponential Lomax distribution was studied by El-Bassiouny et al. (2015), the Weibull Lomax distribution was studied by Tahir et al. (2015), Abdul-Moniem and Abdel-Hameed (2012) developed and studied the Exponentiated Lomax distribution, the transformed-transformer Lomax distribution was studied by Alzaghal and Hamed (2019), Kilany (2016) studied the weighted Lomax distribution, and the Power Lomax distribution was developed and studied by Rady et al. (2016).
For the purpose of this study, let us give a fair description of the Power Lomax ( ) distribution. A random variable has a power Lomax distribution, if its cumulative distribution function ( ) and survival function is respectively The additional parameter produces a tremendous effect on the flexibility of the of distribution in such that the of distribution is decreasing if ≤ 1 and increasing if > 1.
Using the Harris-G distribution proposed and developed by Aly and Benkherouf (2011), we develop the five-parameter distribution called Harris Power Lomax ( ) distribution. Now, Considering a sequence of . . random variables , , . . . , with cumulative distribution function ( ) ( ), probability density function ( ), hazard rate function ℎ( ) and Cumulative hazard rate function ( ). Let = min ( , , . . . , ), taking has a positive integer valued random variable with probability generating function ( ) ( . , ) for c> 0. The survival function ( ) of is given by Where ̅ ( ) = 1 − ( ). Aly and Benkherouf (2011) discoursed several properties of (6) and also provided the basis for the generalization of the Marshall-Olkin class considering the of the Harris distribution (Harris, 1948) for obtaining new distributions. This is given by Although, Aly and Benkherouf (2011) considered > 0, Harris (1948) restricted to an open interval (0,1). This constraint is a result of the fact that Harris distribution arises from a branching process for which each node may originate new nodes with probability proportional to − log( ). Putting (7) in (6) we obtain The corresponding to (8) is given by The equation (9) above reduces to the one proposed and studied by Marshall and Olkin when = 1.
The chief motivation of this study is based on the advantages of the generalized distribution with respect to having a hazard function that exhibits decreasing, increasing and bathtub shapes as well as the flexibility gain in compounding Harris distribution and Power Lomax distributions in modeling real life data. We develop and study a new distribution called the Harris Power Lomax distribution which possesses these properties with a wide range of applications.
The paper is organized as follows. In section 2, we examine the analytical shapes of the probability, survival and hazard functions under investigation and also the quantile function. Several properties of the proposed distribution are considered in section 3, majorly we obtain an expression for the moments (ordinary and incomplete moments), moment generating function, inequality measures (Bonferroni and Lorenz curves), stress-strength parameter, renyl entropy, and order statistics. In section 4, we evaluate the performance of the MLEs using simulation and two applications of HPL distribution to the Air craft windshield data and the Glass fibre data are carried out to examine its modeling potential.
We provide the concluding remarks in section 5.

Harris Power Lomax Distribution
In this sub-section, we present the distribution and derive some of its properties which include cdf, pdf, expansions of the density, hazard function, reversed hazard function, quantile function and sub-models.
Consider ( ) to be any baseline cdf in the interval (0,1), by putting (4) in (8), we obtain the cdf of distribution given by The corresponding to (10) is given by With , , , , > 0. Thus, , , and are the shape parameters and is a scale parameter. Plots of the pdf of distribution for various values of , , and are given in figures (1) and (2). The plots show that the is right (positive) skewed and can be decreasing (L shape).

Expansion of the Density Function
The expansion of the pdf of distribution is described in the sub-section. For > 0 a real non-integer, we use the series representation given by We can rewrite the density of the HPL distribution as Using (12) in (13), we observe that Then (12) can be re-written in reduce form as Consequently, Here, = and ( ) ( ) denote the exponentiated Power Lomax distribution with parameters with scale parameter ( + 1). This suggests that the distribution can be written as a linear combination of Power Lomax density functions. Hence mathematical properties of the distribution can be obtained from those of the properties.

Survival and Hazard Rate Functions
The hazard function for the distribution will be presented in this sub-section. The Survival and hazard function of the are respectively given by for > 0, > 0, > 0 > > 0 >, > 0. The graph of the Survival function is given in figure 3 for various values of the parameters , , , and   Figure 5. Graph of the hazard function of HPL distribution

Quantile Function and Applications
The quantile function of a distribution is an important tool in describing some important properties of the distribution.
In this section, we present the quantile function of the HPL distribution, as well as some related properties, applications, and functions.
In particular, we can obtain the lower quartile( ), middle quartile or median( ), and the Upper quartile ( ) of distribution by respectively taken = 0.25, 0.5 0.75. The expression for the distribution of quartiles is given below as and The mode of distribution can be obtained by solving the following equation ( ( ; , , , , ) ) = 0.

Skewness and Kurtosis Based on Quantile Function
The moments of distribution offers an empirical approach to measure the skewness and kurtosis of a distribution. However, in some instances the moments of a distribution do not always exist. This is a situation with heavy tail distribution such as the Power Lomax and Lomax distribution. In particular, to measure skewness of distribution, we consider Galton Skewness coefficient defined by The sign on the value of provides information on the direction of the skewness of a distribution. In such that, when > 0, the distribution is right skewed, = 0 for symmetric distribution and < 0 for left skewed distribution. The value measures the tail heaviness of the distribution. In general, the bigger the value of is, the heavier is the tail of the distribution. Table 1

Moments of Distribution
Moments plays an important role in statistical analysis and applications. The most important features of a distribution such as measure of central tendency, skewness, kurtosis, dispersion can be studied using moments.
Proposition 3.0. The moment of a random variable following the distribution, denoted by is Proof: By definition, the moment is given by In order to obtain an expression for the moments, we consider the following lemma: , , , ( + 1) = 1 + Consequently, the expression for the moment of distribution is given by By taking = 1, we obtain the mean of , i.e., = . The variance of is obtained by = ( − ) = -. In addition, one can determine the central moment and cumulant of defined by, respectively, Measures of skewness and kurtosis and kurtosis can also be expressed in a similar manner.

Incomplete Moments of Distribution
Proposition 3.1: Let be a random variable following the distribution and for any ≥ 0, let = if ≤ and 0 otherwise. Then, the incomplete moment of X is given by Other quantities can also be derived in similar manner. Other well-known functions defined with the r-incomplete moment include the lower and upper conditional moments of defined respectively, by,

Moment Generating Function of Distribution
Proposition 3.2: if has the ( ; , , ) distribution, then the moment generating function of X, say ( ) is given as Proof: by definition, the moment generating function of a random variable is given by Using Taylor series expansion, Putting equation (23) in (28), we obtain an expression for the moment generating function of distribution as

Inequality Measures
The Bonferroni and Lorenz curves are widely used measures of income inequality of a given population and have various applications in, insurance, economics, medicine, and reliability.

Proposition 3.3: The Bonferroni curve for the distribution is given by
It should be noted that = equal the first about the origin and is obtained by taking = 1 in equation (23).

Stress-Strength Parameter
Suppose we let and be two continuous and independent random variables, where ~ , , and ~ , , , then the series strength parameter, say Ƨ, is defined as Using the pdf and the cdf of HPL given in equation (10) and (11) in (31), we have Now, by applying change of variables = ( ) and = twice in a row, we get Ƨ = ̅ ̅ + + − 1 1, (2 + + )

Renyi Entropy of Distribution
The Renyi entropy can be used to measures the uncertainty in a distribution as defined by Renyi (1961). The larger the value of entropy obtains the greater the level of uncertainty inherent in the distribution.

Order Statistics of Distribution
Let be a random variable that follow the distribution and, for a random sample of size from , say , , . . . , be the order statistics such that : ≤ : ≤. . . ≤ : , where : ∈ , , . . . , for = 1,2, . . . , . In real life testing the study of order statistics is of importance since they found applications in many systems, mainly those that compose of several components that can fail independently of each other. An expression for the order statistics is given by Then putting equation (10) and (13) By taking = 1, we obtain an expression for the first order statistics of distribution and taking = , we obtain an expression for the or the largest order statistics of distribution.

Maximum Likelihood Estimation Method
The maximum likelihood estimation is used to determine the parameters that maximize the likelihood function of the sample data. Suppose we take a random sample , … , from the distribution, the corresponding likelihood function is given by The analytical solution to the system of nonlinear equations presented in (40), (41), (42), (43), and (44) does not exist. Therefore, a numerical method such as OX program, nlminb function in R program etc. is needed to obtain the solution.
The approximate confidence interval of the parameters = ( , , , , ) can be obtained based on the asymptotic distribution of the MLE's of . Considering a large sample under appropriate regularity conditions, the MLEs for the parameters = ( , ,, , ) have approximate multivariate normal distribution and asymptotic variance-covariance matrix which can be used to approximate the inverse of Fisher information matrix. Then 100 1 − 2 % approximate confidence interval of the parameters = ( , , , , ) are: Where ⁄ is the upper percentile of the normal distribution.

Simulation Results of Distribution
We carry out Monte Carlo simulations to show the asymptotic property of the MLEs for the distribution. We calculate Absolute Biases ( ), Standard Error ( ), and mean squared errors (MSEs) of each parameter for different sample sizes. To obtain the results, the process is replicated N = 1,000 times for n = 50, 100, 150, 200, and 250 for fixed choice of parameters = 0.3, = 1.4, = 1.5, = 0.6, and = 1.5. The estimates of the unknown parameters have been obtained by using BFGS method to minimize the total log-likelihood function. The estimated values of the parameters , , , , and with their corresponding Standard Error ( ), Absolute Bias ( ) and Mean Square Error ( ) are displayed in Table 2. The and for an estimator are defined by It can be noticed from the simulated results in Table 2 of the estimators of all parameters decrease with the increasing sample size.  2.964,4.278, 1.506, 2.190, 3.000, 4.305, 1.568, 2.194, 3.103, 4.376, 1.615, 2.223, 3.114, 4.449, 1.619, 2.224, 3.117, 4.485,1.652, 2.229, 3.166, 4.570, 1.652, 2.300, 3.344, 4.602, 1.757, 2.324, 3.376, 4.663. Some descriptive statistics of these data are presented as follows. The minimum observed value is 0.040, while the maximum value is 4.663. The mean, median and variance are 2.563, 2.385, and 1.239, respectively. Since the value of the mean is greater than that of variance, the aircraft windshield data is under-dispersed Data set 2 is obtained from Smith and Naylor (1987), and consists of the strength of 1.5cm glass fibres measured at the National Physical Laboratory, England.
To demonstrate the tractability of distributions, Anderson-Darling ( * ), Cramer-von Misses ( ), Probability value (P-value), Akaike information criterion (AIC), and Hannan-Quinn information criterion (HQIC) are calculated for each of the model considered. The selection procedure for the model that best fit the data is to choose the model having minimum value of these statistics except for the maximum value of the p-value. In this study, numerical results (of maximum likelihood estimates, and goodness of fit criteria) are calculated by using the goodness the package available in R language. The Adequacy package is also used to construct the Total Time on Test (TTT) plot for the two data sets to determine the shape of their hazard function and the Kernel density plot. Figure 6 shows that the two data set has increasing failure rate and figure 7 shows that the Air windshield data moderately skewed to the left and the glass fibre data also skewed to the same direction, The * , , P-value, AIC, and HQIC are given for the sub-models , , , and competitor models, the Kumaraswamy Generalized Power Lomax ( ) by Nagarjuna, V. B. V. (2021), Kumaraswamy Generalized Lomax ( ) distribution by Shams, T. M. (2013). Tables 3, 4 provide the estimated value of the parameters for data set 1 and 2 respectively. Table 4 and 6 shows the goodness of fit results for distribution and competitor models for data set 1 and 2 respectively. These results show that the distribution is more appropriate than the other models considered since it possess the smallest AIC.BIC, HQIC, * and the largest P-value

(i)
TTT plot for Air Windshield data (ii) TTT plot for fibre glass data       Based on Tables 4 and 6, we can conclude that distribution provides the best fit and can therefore be considered the best model in the class of models considered. Figures 8 and 9 gives more information on the flexibility of distribution over its sub-models and other models considered.

Concluding Remarks
In this work, we have presented a new generalization of the distribution, called the distribution. This generalization is obtained by transforming the three-parameter model using Harris-G distribution as suggested by Aly and Benkherouf (2011). The properties of the proposed distribution are discussed. We obtain the analytical shapes of the density and hazard functions of the HPL distribution. We also consider mean deviations, Bonferroni and Lorenz curves, order statistics and Renyi entropy. Maximum likelihood estimation is discussed within the framework of asymptotic log-likelihood inferences including confidence intervals. The five parameter distribution produced monotonically increasing, decreasing, and inverted bathtub hazard rates. In terms of the statistical significance of the model adequacy, the distribution gives to a better fit than the it's sub-models and other competing models considered in this work