A New Generalized Family of Lifetime Distributions Motivated by Parallel and Series Structures

Any given system can be represented as a parallel arrangement of series structures. Motivated by this fact, a general family of distributions is introduced, by adding two extra parameters to a distribution (called baseline distribution), twice compounding with power series distribution. The new family can allow various hazard rate curves that compete well with other alternatives in ﬁtting real data. We derive formal expressions for its moments, generating function, mean residual lifetime and other reliability functions. Certain characterizations of the new family are presented in terms of the ratio of two truncated moments as well as based on the hazard rate function. The maximum likelihood estimation technique is used to estimate the model parameters and a simulation study is conducted to investigate the performance of the maximum likelihood estimates. Finally, two applications of the model with real data sets are presented to illustrate the usefulness of the proposed distribution.


Introduction
The reliability of the parallel and series systems with a random number of components has extensive literature. Throughout the last two decades, several distributions have been proposed to model the lifetime of systems with an unknown number of components. These new distributions are built by compounding a lifetime distribution and a member of the power series family of distributions in series or parallel arrangement. Let U denote the number of components in a system connected in series or parallel structures and let X i denote the lifetime of the ith component. Then, the lifetime of the system can be represented as Y = min(X 1 , X 2 , . . . , X U ) or Y = max(X 1 , X 2 , . . . , X U ).
Marshall and Olkin [21] proposed adding a parameter to the lifetime distribution through compounding with the geometric distribution. Their work was extended by Adamidis and Loukas [5] when the two-parameter exponentialgeometric (EG) distribution was introduced. Similarly, Chahkandi and Ganjali [9] proposed the exponential power series (EPS) distribution; Morais and Barreto-Souza [22] introduced the Weibull power series (WPS) distribution, which contains the EPS distribution as a special case. Flores et al. [15] proposed the complementary EPS distribution, complementary to the EPS distribution. Let X i,j be independent and identically distributed baseline random variables with the probability distribution function (pdf) g(x; τ τ τ ) and cumulative density function (cdf) G(x, τ τ τ ), where τ τ τ denoted the vector of parameters of the baseline distribution. The survival function of the baseline distribution isḠ(x; τ τ τ ) = 1 − G(x; τ τ τ ). Let Z 1 , Z 2 , . . . , Z U be independent and identically distributed power series random variables (truncated at zero) with parameter θ and the probability mass function (pmf) [24] P (z; θ) = b z θ z B (θ) , for z = 1, 2, . . ., where b z depends only on z and B (θ) = ∑ ∞ z=1 b z θ z < ∞. Let U be a power series random variable with parameter λ and probability mass function for u = 1, 2, . . ., where a u depends only on u and A (λ) = ∑ ∞ u=1 a u λ u < ∞. Assume X i,j and Z 1 , Z 2 , . . . , Z U are independent random variables. Table 1 shows useful quantities of some members of the power series family (truncated at zero) such as the Poisson, geometric, logarithmic series, negative binomial and binomial distributions. Geometric It can be shown that the cdf of X is for x 0, ξ ξ ξ = (τ τ τ , θ, λ) and survival functionF (x; ξ ξ ξ) = 1 − F (x; ξ ξ ξ). The corresponding pdf and hazard rate functions (hrf) are and for x ≥ 0 respectively. We shall see later that the hrf can be constant, decreasing, increasing, bathtub-shaped and upside down bathtub-shaped for different baseline and power series distributions. We shall refer to X as a member of the lifetime power series power series (Lifetime PS 2 ) family of random variables. The baseline distribution is a soecial case for N = U = 1. The WPS distribution of Lemos Morais and Barreto-Souza [22] and EPS distribution of Chahkandi and Ganjali [9] are the special cases for U = 1 and U = α = 1, respectively. On the other hand, CEPS Flores et al. [15] is a special case for Z i = α = 1 and i = 1, 2, . . .. Since the compound of lifetime and power series distributions with minimum structure, compound of lifetime and power series distributions with maximum structure and baseline distribution are the lifetime PS 2 distributions lwith special cases as the limiting distributions.

Mathematical properties
Let X be a lifetime PS 2 random variable with the parameter vector ξ ξ ξ = (τ τ τ , θ, λ). Using the binomial expansion and ua u λ u−1 , the cdf and pdf of X can be expanded as Details of the derivation and coefficients are given in Appendix B. The density function is an infinite linear combination of the baseline survival function. Furthermore, in Appendix C, we introduce and calculate the following useful quantity where g(.) andḠ(.) denoted the pdf and survival function of the baseline distribution, respectively. So, the mathematical quantities, such as the moments, moment generating and mean residual functions, of this family, can be derived. From (3) the sth raw moment of X is for s > 0, where κ (s, 1, k + j) = M (s, 0, k + j) is the (s, 0, k + j)th probability weighted moment (PWM) of the baseline distribution defined by Greenwood et al. [19] as follows The first four moments with different values of the parameters are calculated numerically for the Weibull geometric Poisson distribution (WGP for brevity, see Section 5.3). These values are given in Table 2. The incomplete moments of X can be determined from (4) using the incomplete probability weighted moment as where the lower incomplete kappa function is defined by The moment generating function (M X (t)) can be expressed as Given survival to time x 0 , the residual life is the period from x 0 until the time of failure. From the expansion of the survival function (Appendix B, equation (14)), the mean residual lifetime of the lifetime PS 2 distribution is given by where the upper incomplete kappa function κ x0 (a, b, c) is introduced in Appendix C.

784
A NEW GENERALIZED FAMILY OF LIFETIME DISTRIBUTIONS MOTIVATED The entropy of a random variable is a measure of the variation of uncertainty. Shannon's entropy [29] of the random variable X with pdf f (x) is defined by −E {log [f (x)]}. Hence, the Shannon's entropy for the lifetime PS 2 distribution can be expressed in the form First, we can show that We use an equation of Qi et al. [27] for the logarithm of a power series throughout and the fact that the derivative of a power series is a power series as well (Appendix A, equation (10)), to show that Finally, Table 1 shows the closed-form expressions of I (λ) = 1 for Poisson, geometric, logarithmic and other members of the power series distributions.
Let X 1 , X 2 ,. . . ,X n be independent and identically distributed lifetime PS 2 random variables. We can write the density of the ith order statistic, X i:n as Using binomial expansion and equations (12) and (13) where c ′ k,z and ϕ m,u were introduced in Appendix A and The pdf of X i:n is obtained by inserting pdf and cdf of the lifetime PS 2 distribution in (5). From the expansion of B ′ (θ) = ∑ ∞ r=1 rb r θ r−1 and equation (6), we can rewrite density of X i:n as for x ≥ 0. Equation (7) is the main result of this section. It reveals that the pdf of the lifetime PS 2 order statistics is a linear combination of the baseline survival functions. So, several mathematical quantities of the lifetime PS 2 order statistics such as ordinary, incomplete and factorial moments, mgf, mean deviations and several others can be found. With this, we obtain a formula for the sth moment of the ith order statistic X i:m : Let U be a random variable with standard uniform distribution. Then, the following transformation of U has lifetime PS 2 distribution: whereḠ −1 (.) is the inverse of the baseline survival function. Furthermore, A −1 (.) and B −1 (.) are the inverse functions of A(.) and B(.). The ωth quantile of the lifetime PS 2 distribution is ])) .

Characterizations of PS 2 family of distributions
This section deals with various characterizations of PS 2 distribution. These characterizations are based on: (i) a simple relationship between two truncated moments and (ii) the hazard function. It should be mentioned that for characterization (i) the cdf may not have a closed form.
We present our characterizations (i) − (ii) in two subsections.

Characterizations based on two truncated moments
In this subsection, we present characterizations of PS 2 distribution in terms of the ratio of two truncated moments. This characterization result employs a theorem due to Glänzel [16] see Theorem 1 of Appendix D. Note that the result holds also when the interval H is not closed. As shown in Glänzel [17], this characterization is stable in the sense of weak convergence.
The random variable X has pdf (1) if and only if the function η defined in Theorem 1 has the form Proof. If X has pdf (1), then and finally Conversely, if η is given as above, then and hence Stat., Optim. Inf. Comput. Vol. 7, December 2019 Now, in view of Theorem 1, X has density (1) . Corollary 1. Let X : Ω → R be a continuous random variable and let q 1 (x) be as in Proposition 1 Then, X has pdf (1) if and only if there exist functions q 2 and η defined in Theorem 1 satisfying the differential equation The general solution of the differential equation in Corollary 1 is where D is a constant. Note that a set of functions satisfying the above differential equation is given in Proposition 1 with D = 1 2 . However, it should be also noted that there are other triplets (q 1 , q 2 , η) satisfying the conditions of Theorem 1.

Characterization in terms of the hazard function
The hazard function, h F , of a twice differentiable distribution function, F , satisfies the following first order differential equation It should be mentioned that for many univariate continuous distributions, the above equation is the only differential equation available in terms of the hazard function. In this subsection, we present a non-trivial characterization of PS 2 in terms of the hazard function.
Proposition 2. Let X : Ω → R be a continuous random variable. The random variable X has pdf (1) if and only if its hazard function h F (x) satisfies the following differential equation Proof. If X has pdf (1), then clearly the above differential equation holds. If the differential equation holds, then d dx from which we arrive at the hazard function corresponding to the pdf (1).

Some special cases
In this section, we introduce some special cases of the lifetime PS 2 family of distributions. To illustrate the flexibility of the distributions, graphs of the pdf and hazard rate functions for selected distributions are presented.

Exponential power series power series (EPS 2 ) distribution
The cdf and pdf of the exponential distribution with scale parameter β > 0 are given by G(x; β) = 1 − e −βx and g(x; β) = βe −βx , respectively. Inserting these expressions in (1) results in the EPS 2 pdf

Lindley power series power series (LPS 2 ) distribution
Consider the Lindley distribution with scale parameter β > 0 and cdf and pdf are given by respectively. Inserting these expressions in (1) results in the LPS 2 density function

Weibull power series power series (WPS 2 ) distribution
The WPS 2 distribution is defined from (1) by taking G(x; α, β) = 1 − e −βx α and g(x; α, β) = αβx α−1 e −βx α for the cdf and pdf of the Weibull distribution with parameters α and β. The WPS 2 pdf is given by for x ≥ 0, α > 0 and β > 0. Figures 2 and 3 display the density and hazard rate functions of the Weibull geometric Poisson (WGP) and Weibull geometric binomial (WGB) distributions for selected parameter values.

Log-logistic power series power series (LLPS 2 ) distributions
Consider the log-logistic distribution with scale parameter β and shape parameter α, where the cdf and pdf (for x ≥ 0) are given by G (x; α, β) =

790
A NEW GENERALIZED FAMILY OF LIFETIME DISTRIBUTIONS MOTIVATED distribution defined by twice compounding the log-logistic and power series distribution. The LLPS 2 pdf is given by for x ≥ 0 and α, β > 0.

Generalized half-normal power series power series (GHNPS 2 ) distributions
Generalized half-normal distribution was introduced by Cooray and Ananda [11] with cdf and pdf respectively. The GHNPS 2 pdf is obtained by inserting these expressions in (1) f for x ≥ 0 and α, β > 0.

Estimation of the parameters
In this section, we discuss the estimation of the parameters of the lifetime PS 2 distribution. Suppose X = (X 1 , X 2 , . . . , X n ) is a random sample from the lifetime PS 2 distribution with the vector of the observed values x = (x 1 , x 2 , . . . , x n ) and parameter vector ξ ξ ξ = (τ τ τ , θ, λ). The log-likelihood function based on the observed sample is The associated score function is given by The functions g(.) andḠ(.) are defined in Section 2. The MLE estimate of ξ ξ ξ (ξ ξ ξ) is obtained by solving the nonlinear likelihood equations U n (ξ ξ ξ) = ( ∂ℓ ∂τ τ τ , ∂ℓ ∂θ , ∂ℓ ∂λ ) = 0. These equations cannot be solved analytically and a statistical software can be used to solve them via iterative numerical methods such as the Newton-Raphson, quasi-Newton, or Nelder-Mead procedures such as the Newton-Raphson, quasi-Newton, or Nelder-Mead procedures. In the data examples section 8, the MLEs were obtained by directly maximizing (9) with respect to the parameters. The 'optim' routine in R based on the method "BFGS" was used for numerical maximization.
For interval estimation of (τ τ τ , θ, λ) and hypothesis tests on these parameters, we obtain the observed information matrix since the expected information matrix is very complicated and requires numerical integration. The (p + 2) × (p + 2) observed information matrix J n (ξ ξ ξ), where p is the dimension of the parameter vector τ τ τ , becomes Under conditions that are fulfilled for parameters in the interior of the parameter space but not on the boundary and large n, the distribution of √ n ) . This multivariate normal distribution can be used to construct approximate confidence intervals and confidence regions for the individual parameters and for the hazard rate and survival functions.

Simulation study
We conducted a simulation study to assess the performance of the maximum likelihood estimation procedure for estimating the GHNGP distribution parameters using (8). Samples of sizes 10, 12, 14, . . . , 200 are generated for parameter vector ξ ξ ξ = (1, 2, 0.5, 0.9) from GHNGP distribution.  We repeated the simulation k = 100 times and calculated the MLEs and the bias and mean squared error (MSE) of the parameter estimates. The empirical results are given in Figures 4 and 5 shows that as the sample size increases, the biases, and MSEs of the estimators decrease. The broken lines in Figure 4 correspond to the biases being zero. The following observations can be made: 1. The bias for α is generally positive, but θ and λ have a negative bias; 2. the biases for each parameter decrease to zero as n → ∞; 3. the biases and the mean squared errors are the smallest for the parameter, α; 4. the biases and the mean squared errors are the largest for the parameter, λ; 5. the mean squared errors for each parameter decrease to zero as n → ∞; 6. the biases and mean squared errors for each parameter appear reasonably small for all n > 60.

Real data examples
In this section, we fit the lifetime PS 2 family of distributions to two real data sets by the method of maximum likelihood. First, we give the MLEs and the corresponding standard errors of the model parameters and the values of the Akaike Information Criterion (AIC), Corrected Akaike Information Criterion (AICc), Bayesian Information Criterion (BIC) and Kolmogorov-Smirnov (K-S) statistics. The lower the values of these criteria, the better the fit. Finally, we provide the histograms of the data sets to have a visual comparison of the fitted density functions. We used the total time on test (TTT) transform procedure proposed by Aarset [1] as a tool to identify the hazard behavior of the distribution. The TTT-transform can illustrate the variety of the hazard rate curves for a lifetime distribution. If the empirical TTT-transform is convex and concave, the shape of the corresponding hazard rate function is decreasing and increasing, respectively. If the TTT-transform is convex and concave, the failure rate function will have bathtub shape. Finally, If the TTT-transform is concave and convex, the failure rate function unimodal will be more appropriate.

Data set 1
The first data set consists of the number of successive failures for the air conditioning system of each member in a fleet of 13 Boeing 720 airplanes. The pooled data, yielding a total of 213 observations, were first analyzed by Proschan [26] and further discussed in Tahmasbi and Rezaei [31] and Chahkandi and Ganjali [9]. Table 3 gives some descriptive statistics for the data set. The new three parameters distributions, exponential geometric Poisson (EGP), exponential geometric binomial (EGB), generalized half-normal geometric geometric (GHNGG) distributions given by the following pdfs were fitted: for x ≥ 0 and α, β > 0. By using the geometric stability property [21] and the reparameterization γ = λ−θ 1−λ , it is easy to show that GHNGG is an extended case of generalized half-normal geometric (GHNG) distribution with extended γ parameter space (−∞, 0) ∪ (0, 1). The parameter space of Poisson distribution in compound distributions could be extended to (−∞, 0) ∪ (0, +∞) and the parameter space of negative binomial distribution could be extended to (−1, 0) ∪ (0, +∞). More similar extension of the parameter space may be done to power series parameters (See Table 1). In the EGNB distributions, we assumed m = 2, so every fitted distribution has three parameters. The fit of these distributions were compared to the fit of the odd Weibull (OW) [10], beta exponential (BE) [23] and exponentiated Lomax (ELD) [2] distributions specified by the pdfs: for x ≥ 0 and α, β, a, b and γ are positive value parameters. The MLEs (and the corresponding standard errors in parentheses) of the parameters are reported in Table 4. Additionally, to compare the models, we adopt four criterions: AIC, AICc, BIC, and K-S. The values in this table indicate that the new models give the best fit among the fitted models. In particular, we can see that the largest log-likelihood value, the largest p-value, the smallest AIC value, the smallest AICc value, and the smallest BIC value are obtained for the EGNB distribution. Figure 7a shows the TTT-plot for the data set 1 has a concave and convex shape, which indicates that the hazard rate function has a unimodal shape. The estimated hrf for the data set 1 is plotted in the Figure 8a. So, the EGNB distribution could be a proper model for the data set 1.

Data set 2
The second data set consists of the strength of 1.5 cm glass fibres, measured at National physical laboratory, England [30]. Table 5 gives some descriptive statistics for this data set. The new distributions, WGP, WGB, generalized half-normal geometric Poisson (GHNGP) and generalized half-normal geometric binomial (GHNGB) distributions given by the following pdfs were fitted: for x ≥ 0 and ξ ξ ξ = (α, β, θ, λ), where α, β > 0. The parameter spaces of power series distribution are present in Table 1. The MLEs of the parameters are computed and the goodness-of-fit statistics for these models are compared with fit of the popular OW, beta Weibull (BW) [14] and beta generalized exponential (BGE) [8] distributions with pdfs: respectively, for x ≥ 0 and α, β, a, b > 0. The MLEs, log-likelihood value, the corresponding standard errors, the Kolmogorov-Smirnov statistic, its p-value, AIC, AICc and BIC values are shown in Table 6. For second data set, we can see that the largest log-likelihood value, the largest p-value, the smallest AIC, AICc, and BIC values are obtained for the new family. However, results are very close to those obtained by other members of this family, it shows that the GHNGP distribution gives the best fit respect to all indices. The density graphs for the fit of the distributions for both data sets are shown in Figure 6. The fitted pdfs of the EGNB and GHNGP distributions captures the observed histograms better than others for data sets 1 and 2 previously. Hence, we can say that the new family of distributions provides the best fit for at least two real data sets.  Figure 7b shows that the TTT-plot for the data set 2 has a concave shape, that indicates that the hazard rate function has an increasing shape. Figure 8b plots the estimated hrf for the data set 2. Hence, the GHNGP distribution could be an appropriate model for the fitting of such data.

Concluding remarks
In this paper, we proposed a new family of distributions with two extra parameters, referred to as the lifetime PS 2 family of distributions, by compounding a lifetime and twice power series distributions. This family is the first lifetime family of distributions, motivated by systems having both parallel and series structures. The number of parallel and series components is taken to be a power series random variable and the lifetime of components is taken to be iid continuous positive random variables. This new family of distributions contains all of the lifetime distributions constructed by Marshall and Olkin (1997) method. We obtained some of its mathematical properties, including quantiles, moments, moment generating function, mean residual lifetime, Shannon entropy and order statistics properties. Certain characterizations of the PS 2 family of distributions are presented. The proposed distribution is applied to two real data sets illustrating better fits than the popular lifetime distributions.
where the coefficients a n for n ∈ N satisfy We use an equation of Gradshteyn and Ryzhik [18] for a power series raised to a positive integer n given by where the coefficients c n,k (for k = 1, 2, . . .) are obtained from the recurrence equation (c) Using (11) where c z,0 = (a 1 ) z and for u = 1, 2, . . ., we have [r(u + 1) − u]a u+1 d z,u−r .
(d) We now obtain some expansion for the power series. Let A(λ) be a power series, then Using (11) we have Finally, using Cauchy product of two series where ϕ m,u = u ∑ i=0 (i + 1) a i+1 c ′ m,u−i .
Appendix B: Proofs of (2) and (3) Using binomial expansion and equation (12) we have where j = z − 1 and l k,j = (jb 1 ) −1 j ∑ m=1 [m(j + 1) − j] b m+1 l k,j−m . Then respectively. These functions arise in the mean residual lifetime and incomplete moments of reliability models. It is clear that if κ (a, b, c) exist for a lifetime distribution, then incomplete kappa functions exists for every real value x 0 . One can obtain an expression for the κ (a, b, c) for some distributions. For others, this quantity can be evaluated numerically. We obtain closed-form expressions for κ for exponential, Weibull, Lomax, Lindley and log-logistic distributions, respectively.
(a) Exponential distribution where B(.,.) is the beta function.
(e) Log-logistic distribution Appendix D: Theorem 1 Theorem 1. Let (Ω, F, P) be a given probability space and let H = [d, e] be an interval for some d < e (d = −∞, e = ∞ might as well be allowed) . Let X : Ω → H be a continuous random variable with the distribution function F and let q 1 and q 2 be two real functions defined on H such that is defined with some real function η. Assume that q 1 , q 2 ∈ C 1 (H), η ∈ C 2 (H) and F is twice continuously differentiable and strictly monotone function on the set H. Finally, assume that the equation ηq 1 = q 2 has no real solution in the interior of H. Then F is uniquely determined by the functions q 1 , q 2 and η , particularly