Exponentiated extended extreme value distribution: Properties, estimation, and applications in applied fields

1.   Introduction

Data analysis plays a critical role in numerous scientific disciplines, including reliability analysis, health sciences, economics, industry, and environmental studies, among others. To effectively conduct data analysis, it is essential to utilize appropriate models and statistical distributions. The development of new distributions has expanded researchers' options of suitable models, allowing them to accurately examine data characteristics and patterns. Additionally, the introduction of new distributions has increased the accuracy of analysis and predictions, enhancing the understanding of data and the ability to forecast future outcomes. Recently, several methods have emerged to facilitate the derivation of new distributions from existing ones, leading to the creation of more adaptable and improved models that better align with real-world data. In this context, we refer the reader to the literature cited in ^{[1,2,3,4,5,6,7,8,9]}. Xu et al. ^[10] proposed a model to evaluate the reliability of multicomponent systems in dynamic environments, accounting for component correlations and shared environmental effects. Wang et al. ^[11] conducted a study that analyzed clustered panel count data, which are commonly found in biomedical studies with multiple observation points and potential correlations within clusters. They introduced two semiparametric models to address these correlations and prevent biased estimations; these models were supported by simulation studies and a real-world application showcasing their efficacy. Xu et al. ^[12] introduced a bivariate Wiener model to capture the degradation patterns of two key performance characteristics of permanent magnet brakes. They also presented an objective Bayesian method for analyzing degradation data with small sample sizes.

Many newly introduced probability distributions include a hazard rate (HR) function that is important in the analysis of lifetime data, particularly in studies related to survival and reliability. For instance, in the field of reliability engineering, numerous real-life datasets exhibit a distinctive bathtub-shaped HR, which indicates varying behaviors across different phases. In the initial phase, known as the infant mortality phase, the HR starts at a high level, indicating a greater likelihood of early failures in the life cycle of the system. Subsequently, as the system enters the normal life phase, the HR stabilizes and remains relatively constant, indicating a consistent HR during this period. Finally, in the wear-out phase, the HR increases, indicating an escalating probability of failures as the system ages and its components deteriorate. Furthermore, the HR function offers valuable insights into the timing and frequency of failures in medical data. It has extensive applications, such as survival analysis, clinical trial analysis, disease progression modeling, and risk assessment ^[13,14].

Cho et al. ^[15] introduced the exponentiated extreme-value (EEV) distribution, which is defined by the cumulative distribution function (CDF) as follows:

Cho et al. ^[15] examined the characteristics of the EEV distribution and provided approximate maximum likelihood estimators (MLEs) for the scale and location parameters by using multiply type-Ⅱ censored samples. The estimators were evaluated based on the mean squared error (MSE) for different censored samples. One major drawback of the work mentioned above is a lack of HR function, rendering it inadequate for many applications, particularly those involving engineering and medical data. Several conventional distributions cannot model the failure rate patterns observed in real-world data with high accuracy, specifically the characteristic bathtub-shaped trend commonly found in engineering and medical data. While the EEV distribution offers some flexibility, it lacks the ability to fully capture these patterns. To address this limitation, we propose a new three-parameter lifetime model called the exponentiated extended extreme-value (EEEV) distribution. By incorporating an additional parameter, the EEEV distribution can effectively model both increasing and bathtub-shaped HRs, making it a more versatile tool for lifetime data analysis.

In this research paper, we propose a novel distribution by adding an additional variable to the CDF of the EEV distribution. We refer to this new distribution as the EEEV distribution, which possesses several advantageous properties. First, there is the flexibility in representing increasing or bathtub-shaped trends in data. Second, it facilitates the realization of the necessary flexibility when analyzing real-world data, particularly in fields such as engineering and medicine. Third, due to its inclusion of only three parameters, this model is easy to implement. Fourth, its parameters can be estimated by using various methods, and, here, simulation results were used to compare the performance of these estimators. Finally, the EEEV distribution consistently outperformed other competing distributions in terms of goodness of fit. Although the EEEV distribution provides flexibility in the modeling of different failure rate patterns, it is not suitable for all types of data, particularly those exhibiting a decreasing trend. Furthermore, the proposed model was only applied to complete datasets, but its applicability to censored datasets, such as type-Ⅱ progressive censoring, generalized hybrid censoring, and adaptive type-Ⅱ progressive censoring schemes, will be explored in future research.

The paper is organized as follows. Section 2 explores the distribution of the proposed EEEVs. Section 3 establishes some of its statistical properties. Section 4 presents four estimation methods for calculating the EEEV distribution parameters. Section 5 shows the simulation results that were obtained by using four different estimating methods. Section 6 discusses the applicability of the EEEV distribution to three real datasets, illustrating its significance and flexibility. Finally, Section 7 provides a summary of the conclusion.

2.   The EEEV distribution

Within this section, we propose a new model that extends the EEV distribution by incorporating an additional variable into the base CDF for the EEV distribution. The additional variable enhances the flexibility of the new distribution in the modeling of various types of data in different fields, consistently providing a better fit than its competitors. The resulting model is referred to as the EEEV distribution model. The CDF and probability density function (PDF) of the EEEV distribution with the parameter vectors denoted by $\varTheta = (\delta, \, \gamma, \, \eta)$ can be expressed as follows:

and

where $\eta, \, \gamma$, and $\delta$ are the shape, location, and scale parameters, respectively.

The survival function and the HR function of the EEEV distribution are expressed as follows:

and

Figures 1 and 2 depict the graphical performance of the PDF and HR function for the EEEV distribution for various parameter options. Figure 1 shows that the PDF of the EEEV distribution can be unimodal-shaped (see Figure 1(a)) or decreasing (see Figure 1(b)). On the other hand, Figure 2 indicates that the HR function can be increasing (see Figure 2(a)) or bathtub-shaped (see Figure 2(b)).

3.   Statistical properties

3.1.   Quantile and mode

To find the $q$th quantile ($x_{q}$) of the EEEV distribution, one can solve the following equation:

Thus, solving Eq (3.1) yields

Let $u = x_{q}\, {\rm e}^{\delta\, x_{q}};$ it is possible to express $x_{q}$ in terms of $u$ when the value of $\delta$ is positive, as follows:

where

We have verified Eq (3.3) and the power-series expansion for $F(w) = ProductLog[w]$ by using Wolfram Mathematica software, which provides $F(w)$ as the main solution for $z$ in $w = z\, {\rm e}^{z}$. We get

Consequently, it is feasible to represent $x_{q}$ in relation to $u$ by using Eq (3.3) given that $x_{q} = \sum_{k = 1}^\infty\, a_{k}\, u^{k},$ where $a_{k} = \frac{(-1)^{k+1}\, k^{k-2}}{(k-1)!}\, \delta^{k-1}$ and the condition of convergence of this sum is $-\log(1-q^{\frac{1}{\eta}}) < {\rm e}^{-(\gamma+1)}$ (see ^[16,17]). Thus, the quantile for the EEEV distribution is as follows:

Applying $q = 0.25, \, 0.5, \, 0.75$ in Eq (3.4), we obtain the first quartile, the median, and the third quartile of the EEEV distribution, respectively. Also, the skewness (SK) and kurtosis (KU) can be determined according to quartiles as follows:

and

The mode is determined by finding the solution to the following nonlinear equation:

In general, Eq (3.5) lacks an analytical solution, which implies that the mode of EEEV distribution must be estimated by using numerical methods.

3.2.   Ordinary moments

Using Eq (2.2), the $s$th moment of the EEEV distribution can be determined as follows:

Apply the binomial series expansion of $\Big(1-{\rm e}^{-\delta\, x\, {\rm e}^{\delta\, x-\gamma }}\Big)^{\eta-1}$, given by

By substituting Eq (3.7) into Eq (3.6), we obtain

By applying the substitution $u = \delta\, x\, {\rm e}^{\delta\, x-\gamma }$ and employing the same method as previously described, it can be shown that $x = \sum_{j = 1}^{\infty}\, a_{j}\, \big(\frac{{\rm e}^{\gamma}}{\delta}\, u\big)^{j}$, where $a_{j} = \frac{(-1)^{j+1}\, j^{j-2}}{(j-1)!}\, \delta^{j-1}$. This sum's convergence requirement is that $\frac{{\rm e}^{\gamma}}{\delta}\, u < \frac{1}{{\rm e}\, \delta}$ (see ^[17]). Therefore, the $s$th moment of the EEEV distribution can be written as

but

where

Returning to Eq (3.9) gives

After solving the integral presented above, the resulting solution is as follows:

Based on the first four ordinary moments of the EEEV distribution, the measures of SK and KU can be derived as follows:

and

Table 1 provides the values of the first four moments, variance, SK, and KU for the EEEV distribution, considering various combinations of $\eta, \, \gamma$, and $\delta$.

According to Table 1, we can see that the KU of the EEEV distribution falls within the range of $(2.46311, 8.25998)$, whereas the SK can exhibit either negative or positive values, making it suitable for modeling skewed data.

3.3.   Order statistics

Assume that $X_{1}, ..., X_{n}$ denote a random sample from EEEV distribution and let $X_{1:n}, ..., X_{n:n}$ be the corresponding order statistics (OS). The PDF of the $k$th-OS is given by

where $B(k, n-k+1)$ is the beta function. Applying the binomial expansion of $\big[1-F(x; \underline{\Theta})\big]^{n-k}$, Eq (3.11) can be expressed as follows:

By substituting Eqs (2.1) and (2.2) into Eq (3.12), we obtain

As a result,

where $f(x; \delta, \gamma, \eta^{'})$ is the PDF of the EEEV distribution with the parameters $\delta, \gamma, \eta' = (k+l)\, \eta$. Utilizing Eq (3.10), the $r$th moment of the $k$th-OS for the EEEV distribution is given by

3.4.   Incomplete moments

The definition of incomplete moments is as follows:

By utilizing the expansion $\Big(1-{\rm e}^{-\delta\, x\, {\rm e}^{\delta\, x-\gamma }}\Big)^{\eta-1}$ and substituting $u = \delta\, x\, {\rm e}^{\delta\, x-\gamma }$, as explained in the previous subsection, we obtain

Upon solving the integral presented above, the resulting solution is as follows:

where $\gamma_{1} = \Gamma(m_{s}+1)-\Gamma \Big(m_{s}+1, (1+i)\, s\, t\, {\rm e}^{(s\, t-\gamma)} \Big), \, \, \, \, \Re(m_{s}) > -1.$

3.5.   Entropy

The Rényi entropy of $X$ for the EEEV distribution is defined as follows:

By substituting Eq (2.2) into Eq (3.16), we get

By employing series expansions for $(1+\delta\, x)^{\rho}$, ${\rm e}^{-\rho\, \delta\, x\, {\rm e}^{\delta\, x-\gamma }}$, and $\Big(1-{\rm e}^{-\delta\, x\, {\rm e}^{\delta\, x-\gamma }}\Big)^{\rho\, (\eta-1)}$, we obtain the Renyi entropy for the EEEV distribution:

where $I(x; \, \delta, \, d) = \int_{0}^{\infty}\, x^{d}\, {\rm e}^{\delta\, d\, x}$, $\, \quad\, d = \rho+j+l$.

4.   Methods of estimation

Within this section, we estimate the unknown parameters $\delta, \, \gamma$, and $\eta$ of the EEEV distribution by utilizing four different estimation methods. The aim is to illustrate the performance variations among different estimators of the EEEV distribution for different combinations of parameters and sample sizes. The estimation methods employed in this study include the use of MLEs, least-squares estimators (LSEs), weighted LSEs (WLSEs), and Cramér-von Mises estimators (CVMEs). The MLE is a consistent method that converges to true parameter values as the sample size is increased. It is efficient and widely applicable, but it can be sensitive to initial parameter choices and computationally challenging for complex models. The LSE offers simplicity and robustness but may lack efficiency compared to the MLE. The WLSE improves upon the LSE by incorporating weights based on data variance, but it requires knowledge of the variance structure. The CVME, based on the empirical distribution function, demonstrates robustness and reduced sensitivity to tail behavior, but it can be computationally complex. The choice of the most suitable method depends on factors such as sample size, parameter values, data characteristics, and computational resources ^[18,19].

4.1.   MLE

Assuming that $x_{1}, \, x_{2}, \ldots, x_{n}$ denote $n$ independent random variables with the EEEV distribution. Then, the log-likelihood function of the EEEV distribution, denoted as $L(\varTheta; \, x)$, can be expressed as follows:

The MLEs $\hat{\delta}, \, \hat{\gamma},$ and $\hat{\eta}$ can by calculated by maximizing Eq (4.1) by using numerical methods, such as NMaximize in the Wolfram Mathematica software. To obtain the normal equations of $L(\varTheta; \, x)$, we can derive the first partial derivatives of $L(\varTheta; \, x)$ with respect to $\delta, \, \gamma,$ and $\eta$, respectively. By setting these derivatives equal to zero, we have the following equations:

Solving Eqs (4.2)–(4.4) by using numerical methods, such as FindRoot in Wolfram Mathematica software, we can obtain the estimators of the EEEV parameters by the MLE method.

4.2.   LSE

Assume that $x_{1}, \, x_{2}, \ldots, x_{n}$ represent the OS of a random sample of size $n$ drawn from the EEEV distribution. Thus, the LSE ^[20] of the EEEV parameters can be obtained by minimizing the following function:

The previous equation's solution can be obtained by using numerical methods, such as NMinimize in the Wolfram Mathematica software.

Similarly, the LSEs $\hat{\delta}, \, \hat{\gamma},$ and $\hat{\eta}$ can also be determined by solving the following nonlinear equations:

where

and $w_{i} = {\rm e}^{-\delta\, x_i \, {\rm e}^{\delta\, x_i-\gamma }}.$

The solution of $I_{s}\, (x_{i:n}|\varTheta)$ for $s = 1, \, 2, \, 3$ can be obtained numerically.

4.3.   WLSE

Let $x_{1}, \, x_{2}, \ldots, x_{n}$ represent the OS of a random sample of size $n$ drawn from the EEEV distribution. Thus, the WLSE ^[20] of the EEEV parameters can be obtained by minimizing the following function:

Moreover, the WLSEs $\hat{\delta}, \, \hat{\gamma},$ and $\hat{\eta}$ can by obtained by solving the following nonlinear equations:

where $I_{s}\, (x_{i:n}|\varTheta), s = 1, 2, 3$ have been defined in Eqs (4.7)–(4.9).

4.4.   CVME

Let $x_{1}, \, x_{2}, \ldots, x_{n}$ represent the OS of a random sample of size $n$ drawn from the EEEV distribution. Thus, the CVME ^[21] of the EEEV parameters can be obtained by minimizing the following function:

Furthermore, the CVMEs $\hat{\delta}, \, \hat{\gamma},$ and $\hat{\eta}$ can by obtained by solving the following nonlinear equations:

where $I_{s}\, (x_{i:n}|\varTheta), s = 1, 2, 3$ have been defined in Eqs (4.7)–(4.9).

5.   Simulation analysis

Within this section, we detail a simulation study that was conducted to illustrate the effectiveness of different estimators by using data generated from the EEEV distribution. Version 12.3 of Wolfram Mathematica was used for all simulated studies and graphical representations. The numerical procedures were executed by following the algorithm below:

Step 1: Generate random samples of varying sizes, $n =\{25, \, 60, \, 100, \, 200, \, 300\}$, from the inverse CDF of the EEEV distribution, or by utilizing Eq (3.4).

Step 2: Use four sets of parameters to generate the simulation results: set 1 $(\delta = 0.05, \gamma = 1.8, \eta = 0.8)$, set 2 $(\delta = 0.1, \gamma = 1.5, \eta = 0.6)$, set 3 $(\delta = 0.25, \gamma = 2.3, \eta = 0.9)$, and set 4 $(\delta = 0.5, \gamma = 2.6, \eta = 1.3)$.

Step 3: Repeat the process 2000 times for each sample.

Step 4: Calculate the average absolute biases (ABs), MSEs, mean relative errors (MREs), and root mean square error (RMSEs) for the four sets. The computational formulas for the AB, MSE, MRE, and RMSE are provided as follows:

and

where $\hat{\varTheta} = (\hat{\delta}, \, \hat{\gamma}, \, \hat{\eta})^{'}.$

The AB, MSE, MRE, and RMSE results were calculated for all sets; the results are shown in Tables 2 and 3. While Figures 3–6 provide graphical representations of the numerical values. Based on the results of the simulation analysis, as displayed in Tables 2 and 3, we can draw several important conclusions:

● As the sample size $n$ increased, a clear trend in bias reduction was observed across all of the estimating techniques being examined.

● It can be deduced that all estimation methods demonstrate the consistency property, meaning that, as the sample size $n$ increases, the estimators tend to approach the true parameter values.

● All of the estimation techniques performed very well on the task of estimating the EEEV parameters.

● The plots depicted in Figures 3–6 illustrate that the AB, MSE, MRE, and RMSE values diminish to zero with increasing sample size, irrespective of the parameter combinations. These graphical representations suggest that the MLE is the most efficient approach for estimating the parameters of the EEEV distribution.

6.   Applications

In this section, we demonstrate the versatility and significance of the EEEV distribution in the modeling of real-life data by employing three datasets from the fields of medicine and engineering. The first dataset, originally presented by Boag ^[22], details the ages (in months) of 18 patients who passed away due to causes other than cancer. The second dataset, introduced by Aarset ^[23], consists of the failure times of 50 electronic devices that were subjected to life tests starting from time zero. The third dataset, provided by Murthy et al. ^[24], revolves around the service times (in thousand hours) of 63 aircraft windshields of a specific model. Figures 7–9 illustrate the total time test (TTT) plots ^[23] for these three real datasets. Based on Figures 7(a) and 8(a), we can conclude that the empirical HR functions for the first and second datasets exhibited bathtub curves. On the other hand, Figure 9(a) suggests that the third dataset has an increasing HR.

In this analysis, we compared our model with six other competitive lifetime models, i.e., the exponentiated Weibull (EW)^[25], power generalized Weibull (PGW)^[26], EEV ^[15], exponentiated Nadarajah-Haghighi (ENH)^[27], alpha logarithmic transformed Weibull (ALTW)^[28], and logistic Nadarajah-Haghighi (LNH)^[29] distributions. To assess the suitability of these competing distributions, we utilized various goodness-of-fit analytical measures, including the log-likelihood (LL), consistent Akaike information criterion (CAIC), Bayesian information criterion (BIC), Hannan-Quinn information criterion (HQIC), and Cramér-von Mises $(W^*)$, Anderson-Darling $(A^*)$, and Kolmogorov-Smirnov (KS) statistics, along with their corresponding P-value. All calculations were performed by using the Wolfram Mathematica software, as well as the graph construction.

Tables 4, 7, and 10 present the estimated EEEV parameters obtained from four different estimation methods, along with the corresponding goodness-of-fit measures for the Boag, Aarset, and aircraft windshield datasets, respectively. By examining the P-values in Tables 4, 7, and 10, it is suggested that the MLE is suitable for estimating the EEEV parameters for the Boag data, while the WLSE is recommended for the Aarset data. As for the aircraft windshield data, it is advised to employ the CVME method to estimate the EEEV parameters. Tables 5, 8, and 11 provide the MLE method estimates for the EEEV parameters, as well as the estimates for all of the compared distributions based on the Boag, Aarset, and aircraft windshield data. Tables 6, 9, and 12 present the comparison statistics for the Boag, Aarset, and aircraft windshield data. These statistics consist of LL, AIC, BIC, CAIC, HQIC, $A^*$, $W^*$, and KS values with corresponding P-value. The values presented in Tables 6, 9, and 12 indicate that the EEEV distribution outperformed other competing models. This is evidenced by the fact that the EEEV distribution had the lowest values across all measures, as well as the highest P-value.

The estimated PDF, survival function, and HR function, along with the TTT plot for the EEEV distribution and all other considered models, are depicted in Figures 5–9 for the three datasets. The findings presented in Tables 5, 8, and 11 indicate that the EEEV distribution is the most suitable model for fitting the three datasets among all of the investigated distribution models. These findings are further supported by the graphical representations in Figures 7–9.

7.   Conclusions

In this paper, we have proposed and examined the EEEV distribution as an extension of the EEV distribution. Its associated HR function can be bathtub-shaped or increasing. Some of its statistical properties have been derived. The estimation of the parameters for the EEEV distribution was performed by using four different estimation methods. The behaviors of these estimators were assessed via simulation. The practical applicability of the EEEV distribution has been illustrated by analyzing three real-life datasets from the fields of medicine and engineering. The analytical measures indicated that our EEEV distribution provided a good fit compared to other competing distributions.

Author contributions

M. G. M. Ghazal: Conceptualization, Data curation, Methodology, Investigation, Software, Writing-review & editing; Yusra A. Tashkandy: Conceptualization, Investigation, Methodology, Project administration, Funding acquisition, Writing-original draft; Oluwafemi Samson Balogun: Conceptualization, Data curation, Writing-original draft, Formal analysis, Investigation; M. E. Bakr: Methodology, Investigation, Funding acquisition, Writing an original draft. All authors have read and approved the final version of the manuscript for publication.

Use of AI tools declaration

The authors declare that they have not used artificial intelligence tools in the creation of this article.

Acknowledgments

This research project was supported by the Researchers Supporting Project (RSP2024R488), King Saud University, Riyadh, Saudi Arabia.

Conflict of interest

There is no conflict of interest declared by the authors.