Harmonic Mixture-G Family of Distributions: Survival Regression, Simulation by Likelihood, Bootstrap and Bayesian Discussion with MCMC Algorithm

To study the heterogeneous nature of lifetimes of certain mechanical or engineering processes, a mixture model of some suitable lifetime distributions may be more appropriate and appealing than simpler models. In this paper, a new mixture family of the lifetime distributions is introduced via harmonic weighted mean of an underlying distribution and the distribution of the proportional hazard model corresponding to the baseline model. The proposed class of distributions includes the general Marshall-Olkin family of distributions as a special case. Some important properties of the proposed model such as survival function, hazard function, order statistics and some results on stochastic ordering are obtained in a general setting. A special case of this new family is considered by employing Weibull distribution as the parent distribution. We derive several properties of the special distribution such as moments,hazard function survival regression and certain characterizations results. Moreover, we estimate the parameters of the model by using frequentist and Bayesian approaches. For Bayesian analysis, five loss functions, namely the squared error loss function (SELF), weighted squared error loss function (WSELF), modified squared error loss function (MSELF), precautionary loss function (PLF), and K-loss function (KLF) are considered. The beta prior as well as the gamma prior are used to obtain the Bayes estimators and posterior risk of the unknown parameters of the model. Furthermore, credible intervals (CIs) and highest posterior density (HPD) intervals are also obtained. A simulation study is presented via Monte Carlo to investigate the bias and mean square error of the maximum likelihood estimators. For illustrative purposes, two real-life applications of the proposed distribution to Kidney and cancer patients are provided.


Introduction
In recent years, extensive efforts have been made to present new models in the area of distribution theory and related statistical applications. These studies deal mainly with the modeling of various data sources and finding out the probabilistic structure of the model. In connec-tion with the development of new models, it is worthwhile to note that these new models should have the capability for analyzing a wide range of real observations. Undoubtedly, this is also the most basic concern in the development of the new models from the past to the future. Mixture models are quite versatile and thus have been frequently used in different fields of sciences such as reliability, information theory, economic, engineering, agriculture, to name a few. For example, Mendenhall and Hader (1958), while referring to the practical situations encountered by engineers, pointed out that the failure of a system or a device may be divided into two or more different types of causes. Further, Acheson and McElwee (1952) categorized the failures of the electronic tube into gaseous defects, mechanical defects, and normal deterioration of the cathode in order to find the proportion of failure due to a certain cause. Another example is that of an engineering system which consists of different subsystems. These subsystems may be homogeneous or heterogeneous. Heterogeneity nature of such systems can not be captured by a single probability model but it can be captured through mixture models. In the following, three well-known mixture distributions are reviewed. The usual arithmetic mixture distribution is defined by the weighted mean of two probability density functions (PDFs): The geometric mixture is defined by the normalized geometric mean of two PDFs: provided that the normalizing factor k < ∞ on the support of h. This model has been derived in statistics, and in physics in different contexts. Bercher and Vignat (2009) called (2) the generalized escort probability distribution. A more general mixture distribution, called here as the power mean mixture, is defined by the following PDF: where k = [θf α (x) + (1 − θ)g α (x)] 1 α dx.
The α-mixture family (3) contains (1) as a especial case with α = 1 and (2) as a limiting case for α → 0. Asadi, Ebrahimi, Kharazmi, and Soofi (2018) obtained some interesting informational properties for the above mixture models. The above three mixture distributions play important roles in the information theory and shown that they contain optimal information. For more details see Asadi et al. (2018).
Motivated by optimal information properties of these mixture distributions, we propose a new lifetime distribution based on the harmonic mean mixture of two survival functions. The proposed family of distributions is called harmonic mixture G (HM G) distribution. The new HM G distribution includes Marshall-Olkin family of distributions introduced by Marshall and Olkin (1997).
Our main motivation to introduce this new category of distributions is to provide more flexibility for fitting real data sets compare to the other well-known classic statistical distributions. We first derive certain statistical and reliability properties of the HM G distribution in a general setting and then we consider a special case of this model by employing the Weibull distribution as the parent distribution G. It is called HM W distribution. We provide a comprehensive discussion about the survival regression of the new HM W model. Furthermore, we consider the Maximum likelihood, bootstrap estimation and Bayesian procedures to estimate the unknown parameters of the new model for the real application. In addition, the asymptotic confidence intervals and parametric and non-parametric bootstrap confidence intervals are calculated.
In the Bayesian discussion, we consider different types of symmetric and asymmetric loss functions such as squared error loss, weighted squared error, modified squared error, precautionary, K-loss, linear exponential and general entropy loss functions to estimate the four unknown parameters of HM G model. Since all of the parameters are positive, we use gamma and beta prior distributions. Bayesian 95% credible and highest posterior density (HP D) intervals (see Chen and Shao (1999)) are provided for each parameter of the proposed model. In addition, the asymptotic confidence intervals and parametric and non-parametric bootstrap confidence intervals are calculated for comparison with the corresponding Bayesian intervals.
In addition, a simulation study is performed to investigate M LEs of consistency.
The rest of the manuscript is organized as follows. In Section 2, we introduce a new class of distributions called harmonic mixture-G (HM G) distributions and consider the hazard function, quantiles and a discussion about informational properties in a general setting. In Section 3, we consider the Weibull as the parent distribution and introduce generalized harmonic mixture Weibull distribution. This new model is referred to as HM W distribution. Also, in this section, we plot the density function, hazard function and 3D plots of skewness and kurtosis for the HM W distribution. Some characterization results for the HM W distribution are provided in Section 4. In Section 5, the estimation of the parameters of the HM W distribution are obtained via three methods: maximum likelihood, Bayesian and bootstrap estimations. Bayesian analysis results are provided by considering five well-known loss functions. Also, we study the performance of the maximum likelihood estimates of the parameters of HM W distribution via a simulation study. Section 6 is devoted to a discussion about survival regression of HM W model. In Section 7, the superiority of the new model to some competitor statistical models is shown through different criteria of selection model by analyzing a real example. Moreover, Bayesian analysis and associated plots of the posterior samples corresponding to this data are provided in this Section. Also, in this Section we provide numerical analysis for survival regression of HM W distribution via a Kidney data set. Finally, the paper is concluded in Section 8.

Harmonic mixture-G (HM G) family of distributions
Let X be a continuous random variable with survival and hazard functionsḠ(x) and r(x), respectively. The model with hazard function r P h = αr(x); α > 0, is known as proportional hazards (P H) model, popularized by Cox (1972). The survival function of P H model corresponding to the baseline survival functionḠ(x), is given asḠ α (x). This model is widely used in various applications in many fields especially in survival analysis and economics. In this section, first we propose a new model based on the harmonic mixture mean of two survival functionsḠ(x) andḠ α (x) and then study some of its main properties in a general setting.
Lemma 1. If α = 0, the general Marshall-Olkin family of distribution is obtained as The associated CDF and P DF are given, respectively, by and

Hazard rate function
Here, first we obtain the hazard rate function of HM G model and then two associated lemmas are given. The hazard rate is a fundamental tool in reliability modeling for evaluation of the aging process. Knowing the shape of the hazard rate is important in the reliability theory, risk analysis and other disciplines. The concepts of increasing, decreasing, bathtub shaped (first decreasing and then increasing) and upside down bathtub shaped (first increasing and then decreasing) hazard rate functions are very useful in the reliability analysis. The lifetime distributions with these aging properties are designated as the IFR, DFR, BUT and UBT distributions, respectively. The hazard function of the HM G distribution is In fact the hazard rate function of the new model is a weighted version of the baseline hazard . Lemma 2. In view of (6), we have: • If r(x) is increasing and α ≥ 1 then r G (x) is increasing.
• If r(x) is decreasing and 0 < α ≤ 1 then r G (x) is decreasing.
Proof. The proof is straightforward.
In the following lemma we provide a result regarding the stochastic ordering of the hazard function to the compare proposed model and baseline distribution. First, we recall the following definition. The random variable X is said to be less than the random variable Y in the hazard rate order, X ≤ hr Y , if h X (x) ≥ h Y (x), for all x in the union of supports of X and Y , where h X (x)(h Y (x)) is the hazard rate of X(Y ). For more details see Shaked and Shanthikumar (2007).
Lemma 3. Suppose that X is a random variable with the baseline distribution G and X G is harmonic mixture random variable associated with X, and have density function (1).
Proof. The proof is straightforward.

Harmonic mixture-Weibull (HM W ) distribution
In this section, we specialize HM G distribution, which is described in the previous Section, by choosing special case for baseline distribution G. We apply the HM G method to a specific case of baseline distribution, namely to a Weibull distribution and call this special model, four-parameter HM W distribution. Definition 1. A random varible X has HM W (α, β, λ, θ) distribution, if its P DF is given by where, x > 0, α, β, λ > 0 and 0 < θ < 1.
The CDF corresponding to (7) is The survival and hazard rate functions arē respectively. Some plots of PDF and hazard function for the selected parameter values are given in Figures 1 and 2.

Some properties of the HM W distribution
In this section, we obtain some properties of the HM W distribution, such as the moments and order statistics distribution.

Moments
In this subsection, moments and related measures including coefficients of variation, skewness and kurtosis are presented. Tables of values for the first six moments, standard deviation (SD), coefficient of variation (CV ), coefficient of skewness (CS) and coefficient of kurtosis (CK) are also presented. The rth moment of the HM W distribution, denoted by µ r , is The variance, CV , CS, and CK are given by and respectively.
In order to investigate and analyze the amount of skewness and kurtosis of the new model under the three parameters α, λ, β and θ, 3D diagrams are presented in

Order statistics
Order statistics play an important role in probability and statistics. In this subsection, we present the distribution of the ith order statistic from the HM W distribution. The P DF of the ith order statistic from the HM W P DF , f HM W (x), is given by Using the binomial expansion

Characterization results
In this section we establish certain characterizations of the HMW distribution in three directions: (i) based on two truncated moments; (ii) in terms of the hazard function and (iii) based on the conditional expectation of a function of the random variable. These characterizations will be presented in three subsections.

Characterizations based on two truncated moments
This subsection deals with the characterizations of HMW distribution in terms of a simple relationship between two truncated moments. We will employ Theorem 1 of Glänzel (1987) [1] given in the Appendix A. As shown in [2], this characterization is stable in the sense of weak convergence. (7) if and only if the function ξ defined in Theorem 1 is of the form Proof. If X has density (7), then Conversely, if η is of the above form, then Now, according to Theorem 1, X has PDF (7).
Corollary 4.1.1. Suppose X is a continuous random variable. Let q 1 (x) be as in Proposition 4.1.1. Then X has density (7) if and only if there exist functions q 2 and ξ defined in Theorem 1 for which the following first order differential equation holds Corollary 4.1.2. The differential equation in Corollary 4.1.1 has the following general solution where D is a constant. A set of functions satisfying the above differential equation is given in Proposition 4.1.1 with D = 0. Clearly, there are other triplets (q 1 , q 2 , ξ) satisfying the conditions of Theorem 1.

Characterization based on hazard function
The hazard function, h F , of a twice differentiable distribution function, F , satisfies the following differential equation The following proposition establishes a non-trivial characterization HMW distribution based on the hazard function.
Proposition 4.2.1. Suppose X is a continuous random variable. Then, X has density (7) if and only if its hazard function h F (x) satisfies the following first order differential equation Proof. If X has density (7), then clearly the above differential equation holds. Now, if the differential equation holds, then which is the hazard function corresponding to the PDF (7).

Characterizations based on the conditional expectation of a function of the random variable
Hamedani (2013) [3] established the following proposition which can be used to characterize the HMW distribution.

Inference procedure
In this section, we consider estimation of the unknown parameters of the HM W (α, β, λ, θ) distribution via three methods: maximum likelihood method, bootstrap estimation and Bayesian procedure.

Maximum likelihood estimation
Let x 1 , . . . , x n be a random sample from the HM W distribution and ∆ = (α, β, λ, θ) be the vector of parameters. The log-likelihood function is given by The elements of the score vector are given by respectively.

Bootstrap estimation
The parameters of the fitted distribution can be estimated by parametric (resampling from the fitted distribution) or non-parametric (resampling with replacement from the original data set) bootstraps resampling (see Efron and Tibshirani (1994)). These two parametric and nonparametric bootstrap procedures are described below.
Parametric bootstrap procedure: • Estimate ψ (vector of unknown parameters), sayψ , by using the M LE procedure based on a random sample.
• Generate a bootstrap sample {X * 1 , . . . , X * m } usingψ and obtain the bootstrap estimate of ψ, say ψ * , from the bootstrap sample based on the M LE procedure.

• Repeat
Step 2 N BOOT times.
In case of the HM W distribution, the parametric bootstrap estimators (PBs) of α, β, λ and θ, areα P B ,β P B ,λ P B andθ P B , respectively.

Nonparametric bootstrap procedure
• Generate a bootstrap sample {X * 1 , . . . , X * m } , with the replacement from the original data set.
• Obtain the bootstrap estimate of ψ with MLE procedure, say ψ * , by using the bootstrap sample.
• Repeat Step 2 N BOOT times.
In the case of the HM W distribution, the nonparametric bootstrap estimators (NPBs) of α, β, λ and θ, areα N P B ,β N P B ,λ N P B andθ N P B , respectively.

Bayesian inference
Bayesian inference procedure has been used by many statistical researchers, especially researchers in the field of survival analysis and reliability engineering. In this section, a complete sample data is analyzed through Bayesian point of view. We assume that the parameters α, β, λ and θ of HM W distribution have the following independent prior distributions where a, b, e, f , n 0 and n 1 are positive. Hence, the joint prior density function is In the Bayesian estimation, we do not know the actual value of the parameter, which may be adversely affected by loss when we choose an estimator. This loss can be measured by a function of the parameter and the corresponding estimator. For the Bayesian discussion, we consider different types of symmetric and asymmetric loss functions such as squared error loss function (SELF ), weighted squared error loss function (W SELF ), modified squared error loss function (M SELF ), precautionary loss function (P LF ) and K-loss function (KLF ). These loss functions, associated Bayesian estimators and posterior risks are presented in Table 1.
It is clear from the equations (14) and (16) that there are no closed-form expressions for the Bayesian estimators under the five loss functions described in Table 1. Because of intractable integrals associated with joint posterior and marginal posterior distributions, we need to use a software to solve integral equations numerically via MCMC method. The two most popular MCMC methods are: the Metropolis-Hastings algorithm (Metropolis, Rosenbluth, Rosenbluth, Teller, and Teller 1953;Hastings 1970) and the Gibbs sampling (Geman and Geman 1984). Gibbs Sampling is a special case of the Metropolis-Hastings algorithm which generates a Markov chain by sampling from the full set of conditional distributions.
The Gibbs sampling algorithm can be described as follows: Suppose that the general model f (x|ψ) is associated with parameter vector ψ = (ψ 1 , ψ 2 , ..., ψ p ) and observed data x. Thus, the joint posterior distribution is π(ψ 1 , ψ 2 , ..., ψ p |x). We also assume that ψ 0 = (ψ 1 , ψ 2 , ..., ψ p ) is the initial values vector to start Gibbs sampler. The Gibbs sampler draws the values for each iteration in p steps by drawing a new value for each parameter from its full conditional distribution given the most recently drawn values of all other parameters. In symbols, the steps for any iteration, say iteration k, are as follows: • Starting with an initial estimate (ψ , ..., ψ k−1 p , x ; and so on down to • Draw ψ k p from π ψ p |ψ k 1 , ψ k 2 , ..., ψ k p−1 , x . In case of the HM W distribution, by considering the parameter vector Ψ = (α, β, λ, θ) and initial parameter vector Ψ 0 = c(α 0 , β 0 , λ 0 , θ 0 ), the posterior samples are extracted by the above Gibbs sampler where the full conditional distributions are given as Often Bayesian inference requires computing intractable integrals to generate posterior samples. Using Gibbs sampling, one can obtain samples from the joint posterior distribution. In practice, simulations related to Gibbs sampling are conducted through a special software WinBUGS. WinBUGS software was developed in 1997 to simulate data of complex posterior distributions, where analytical or numerical integration techniques cannot be applied. Also, we can use OpenBUGS software, which is an open source version of WinBUGS. Since there are no prior information about the hyper parameters in (12), one can implement the idea of Congdon (2001) and these parameters can be chosen as a = b = c = d = e = f = n 0 = n 1 = 0.0001. Hence, we can use M CM C procedure to extract posterior samples of (14) by means of Gibbs sampling process in OpenBUGS software.

Monte Carlo simulation study
In this subsection, we assess the performance of the M LEs of the parameters with respect to the sample size n for the HM W (α, β, λ, θ) distribution. The assessment of the performance is based on a simulation study using the Monte Carlo method. First, we generate samples of size n from (7). The inversion method is used to generate samples of the HM W distribution generated as the root of where U ∼ U (0, 1) is a uniform variate on the unit interval. Letα,β,λ andθ be the M LEs of the parameters α, β, λ and θ, respectively. We compute the mean square error (M SE) and bias of the M LEs of the parameters α, β, λ and θ, based on the simulation results of N = 1000 independent replications. results are summarized in Table 2 for selected values of n, α, β, λ and θ. From Table 2 the results verify that M SE and bias of the M LEs of the parameters decrease as sample size n increases. Hence, the M LEs of α, β, λ and θ, are consistent estimators.

The LHMW regression model
Let the random variable X follow a Weibull distribution with P DF and CDF , respectively as Inserting (17) and (18) in (2), we have  The HM W (19) is considered as different parametrization of the HM W distribution. Assume that the random variable X follows the HM W distribution, given in (19). We obtain the log-HM W (LHM W ) distribution by applying Y = log(X) transformation and considering the re-parametrization, k = 1/σ and λ = exp(µ). The resulting P DF is where µ ∈ is the location parameter, σ > 0 is the scale parameter, and α > 0 and θ > 0 are the shape parameters. We refer to equation (20) as the LHM W distribution, say Y ∼ LHMW(α, θ, σ, µ). The survival function corresponding to (20) is given by Let Z = (Y − µ)/σ, then the P DF of the standardized random variable is Figure 9 displays the some possible shapes of the LHM W distribution. It is seen that the LHM W distribution could be a good choice to model left skewed lifetime dependent variable with some covariates. α=0. 5,θ=0.5,σ=1,µ=0.4 α=0.5,θ=0.4,σ=1,µ=0.4 α=0.5,θ=0.3,σ=1,µ=0.4 α=0.5,θ=0.2,σ=1,µ=0.4 α=0.5,θ=0.1,σ=1,µ=0.4 Figure 9: The P DF plots of the LHM W distribution Now, using the LHM W density, we introduce a new location-scale regression model where the regression structure is and v i is the explanatory variable vector and z i is the random error following the density in ( Here, we use the M LE method to estimate the unknown parameters of the LHM W regression model. First, we define some useful mathematical notations. Let y i be a dependent variable, distributed as LHM W , given in (19). The dependent variable is defined as y i = min{log(x i ), log(c i )} where log(x i ) and log(c i ) represent the log-lifetime and logcensoring times, respectively. Additionally, we split the observations into the two set. These are F and C which indicate the individuals with log-lifetime and log-censoring, respectively. In view of these definitions, the log-likelihood function of the LHM W regression model is where τ = (α, θ, σ, β ) is the parameter vector, u i = exp(z i ), z i = (y i − v i β)/σ and r is the number of uncensored observations. The MLEs of τ can be obtained by minimizing the negative value of the log-likelihood function in (24). The optim function of R software is used for this purpose.

Residual analysis
To decide the accuracy of the fitted regression model, we analyze the departure from error distribution by means of the residual analysis. Here, two residuals are used. These are martingale and modified deviance residuals. Residual analysis is an important step of the any regression analysis to check the adequacy of the fitted model. The martingale residuals of the LHM W regression model are where u i = exp((y i −x i β)/σ). Using the martingale residuals, the modified deviance residuals of the LHM W regression are where r M i is the martingale residual. The modified deviance residuals are more acceptable and used than martingale residuals. The reason is that the modified deviance residuals are normally distributed with zero men and unit variance once the fitted regression model is suitable and accurate for the given data.

Practical data applications
In this section, we provide two applications for modeling the HM W distribution to real data sets for illustrative purposes. These applications will show the flexibility and usefulness of the HM W distribution. In order to achieve this goal, first we consider the strengths of 1.5 cm glass fibres data set and the parameter estimations are done by means of three methods (maximum likelihood, Bayesian and bootstrap) which are discussed in Section 5. Second we show the performance of of survival regression model LHM W distribution via maximum likelihood method by analyzing Kidney data set which is associated with covariate variables.

Univariate data modeling
The data set obtained from Smith and Naylor (1987) represents the strengths of 1.5 cm glass fibres, measured at the National Physical Laboratory, England. The observations are as follows; Graphical measure: The total time test (T T T ) plot due to Aarset (1987) is an important graphical approach to verify whether the data can be applied to a specific distribution or not. The T T T plot for this data set presented in Figure 10 indicates that the empirical hazard rate functions of the strengths of glass fibres data is increasing. Therefore, the HM W distribution is appropriate to fit these data.

Bootstrap inference for HM W parameters
Here we obtain point and %95 confidence interval (CI) estimation of parameters for the HM W distribution by parametric bootstrap method for the real data set. We provide results of bootstrap estimation based on 1000 bootstrap replicates in Table 4. It is interesting to look at the joint distribution of the bootstrapped values in a scatter plot in order to understand the potential structural correlation between parameters. The corresponding plots of the bootstrap estimation are shown in Figure 11.

M LE estimation and comparison with other models
Here, we fit the HM W distribution to the strengths glass fibres data set and compare it with the Marshall-Olkin Weibull (M OEW ), Log Logistic (LL) and Weibull densities. Table 5 shows the M LEs of parameters, log-likelihood, Akaike information criterion (AIC), Cramervon Mises (W * ), Anderson-Darling (A * ) and p − value (P ) statistics for the data set. The HM W distribution provides the best fit for the data set as it shows the lowest AIC, A * and W * than other considered models. The relative histogram, fitted HM W , M OEW , LL and Weibull P DF s and corresponding empirical and fitted CDFs are plotted in Figure 13 for the current data set. In addition, the Q − Q and P − P plots for the HM W and other fitted distributions are also displayed in Figure 14. These plots support the results in Table 5. We compare the HM W model with a set of competitive models, namely: (i) The log-logistic distribution: The two parameter log-logistic (LL) distribution density function is given by (ii) The Marshall-Olkin extended Weibull distribution (Ghitany, Al-Hussaini, and Al-Jarallah 2005): The three parameter Marshall-Olkin extended Weibull (M OEW ) distribution density function is given by x > 0, α, β, λ > 0,ᾱ = 1 − α.
(iii) The Weibull distribution: The Weibull (W e) distribution with shape parameter α and scale parameter β has density given by

Bayesian estimation under different loss functions
Here, we provide Bayesian estimation results for the parameters of HM W distribution. Bayesian estimators associated with the parameters of HM W distribution are computed based on the single chain of 20, 000 cycles of Gibbs sampler with a conservative burn-in period of the first 5000 iterations. The corresponding Bayesian point and interval estimations and posterior risk are provided in Table 6 for the Strengths of glass fibres data set. Table  7 provides 95% credible and HP D intervals for each parameter of the HM W distribution. The convergence of Gibbs sampler process are verified through graphical inspection (Trace, Autocorrelation and Histogram plots) of the posterior sampled values. It is observed that Gibbs samples of all the parameter estimates achieve a good stationary phase for both data sets. We provide the posterior summary plots for 10, 000 cycles of Gibbs sampler in Figures  15, 16 and 17 .

Data modeling with covariate variables: Kidney data set
In this application, we asses the performance of LHM W regression model via application to a real data set. The data set represents the recurrence times to infection of Kidney patients. The data set is avaliable in survival package of R software. The censoring rate of the data is approximately 23%. The dependent variable, recurrence times to infection y i , is modeled by age x i1 and gender x i2 (0=male, 1=female). The below regression model is fitted where y i is distributed as LHM W .

Parameter estimation
The estimated parameters of the LHM W regression model is obtained by M LE method using the optim function of R software. The estimated parameters, standard errors (SEs), confidence intervals (CIs) as well as the corresponding p-values are listed in Table 8. From   results, we conclude that the estimated regression parameters are statistically significant at 1% level. According to the estimated regression parameters, the recurrence times to infection decreases when the age of patient increases. Additionally, the recurrence times to the infection of female patients are higher than those of the male patients.

Results of residual analysis
The suitability of the fitted LHM W regression model is evaluated via residual analysis. The plot of the modified deviance residuals and its quantile-quantile (Q − Q) plot are displayed in Figure 18 which reveal that the fitted LHM W regression model provides a good fit to the

Conclusion
A new family of lifetime distributions is introduced via harmonic mixture mean and its main properties are derived. One of the interesting and important properties of the proposed family is that it includes the Marshall-Olkin family of distributions, as a special case. An example of this family is proposed by considering the Weibull as the baseline distribution. We provide survival regression model of the special model and a comprehensive discussion about Bayesian estimation of the parameters are studied under various loss functions. Numerical results of a maximum likelihood, Bayesian and bootstrap procedures for a univariate real data set is investigated. Moreover, the associated plots to evaluate the results obtained from the the three methods are provided. The data analysis proves, empirically, that the proposed distribution provides a better fit than its sub-models and other competing distributions for the current data. Finally, the performance of the survival regression model of the special distribution is examined by considering a real example with covariate variables.
This stability theorem makes sure that the convergence of distribution functions is reflected by corresponding convergence of the functions q 1 , q 2 and ξ, respectively. It guarantees, for instance, the 'convergence' of characterization of the Wald distribution to that of the Lévy-Smirnov distribution if α → ∞.
A further consequence of the stability property of Theorem 1 is the application of this theorem to special tasks in statistical practice such as the estimation of the parameters of discrete distributions. For such purpose, the functions q 1 , q 2 and, specially, ξ should be as simple as possible. Since the function triplet is not uniquely determined it is often possible to choose ξ as a linear function. Therefore, it is worth analyzing some special cases which helps to find new characterizations reflecting the relationship between individual continuous univariate distributions and appropriate in other areas of statistics.
In some cases, one can take q 1 (x) ≡ 1, which reduces the condition of Theorem 1 to E [q 2 (X) | X ≥ x] = ξ (x) , x ∈ H. We, however, believe that employing three functions q 1 , q 2 and ξ will enhance the domain of applicability of Theorem 1.